Tony102's Blog - 数学
type: channel
link: https://tony102.com/index.php/tag/%E6%95%B0%E5%AD%A6/
items: 2-1_Finding Roots of Equations (1-Bisection Method)
items: Numerical Analysis: Taylor Series
items: Numerical Analysis: Error
items: 一些在做 3D 中遇到的数学知识
items: Discretization process for State Space Model of Mamba
items: Miscellaneous Math for CG 为计算机图形学准备的数学
items: 排序不等式
items: ARC125B Squares
items: Luogu1286 两数之和

2-1_Finding Roots of Equations (1-Bisection Method)
type: item
date: 2025-03-18T00:29:00+08:00
description: Problem StatementGiven some functions $f(x)$, determine the value(s) of x where $f(x) = 0$Especially when $f(x)$ is not analytically.Root finding algorithmsTwo main types:Bracketing MethodsOpen MethodsFollowing steps are essential for any root find algorithms:Guess: one more more initial values for $x$Iterate: using the scheme to improve the guessCheck the convergenceBracketing methodPut it simply, it a kind of method like binary search. We can test if $f(a)f(b) < 0$ to determine whether we need to discard half of the interval.The approximation relative error is$$ \epsilon_r = \frac{x_{i+1} - x_i}{x_{i+1}} $$When $|\epsilon_r|$ becomes less than a specified threshold $\epsilon_s$ the iterations are stopped.For this algorithms, we are able to set the iteration number $n$. The absolute value of error $|E| = |V_t - V_a|$. Let $L = b- a$ be the length of the interval, so we know the interval length after $n$ iteration his $\frac{L}{2^n}$Therefore,$$ |E| \leq \frac{L}{2^n} $$So we can determine the iteration number $n$ as:$$ n = \frac{\ln \frac{L}{E_{max}}}{\ln 2} $$ImprovementTo improve the root position estimation, we can construct a straight line between point $a, b$, and then use the intersection point with x axis as the initial point.As shown in the figure,$$ x_i = b - f(b) \frac{b-a}{f(b)-f(a)} $$

Numerical Analysis: Taylor Series
type: item
date: 2025-03-18T00:14:00+08:00
description: Pre-requisite KnowledgeRolle's theoremLet $f:[a, b] \to \mathbb{R}$ be differentiable. Assume that $f(a) = f(b)$, then there is at least one point $c \in (a, b)$ where $f'(c) = 0$Mean value theoremLet $f:[a, b] \to \mathbb{R}$ be differentiable. Then there is at least one point $c \in (a, b)$ such that $f'(c) = \frac{f(b) - f(a)}{b - a}$Taylor theoremTaylor's theorem is a fundamental result in calculus that provides a way to approximate a function by its Taylor series, which is an infinite sum of terms calculated from the values of the function's derivatives at a single point.Statement of Taylor's Theorem:For a function $ f(x) $ that is $ n $ times differentiable at a point $ a $, Taylor's theorem states that the function can be approximated near $ a $ by the $ n $th-order Taylor polynomial plus a remainder term.The Taylor polynomial $ P_n(x) $ of $ f(x) $ centered at $ a $ is given by:$$ P_n(x) = f(a) + f'(a)(x - a) + \frac{f''(a)}{2!}(x - a)^2 + \cdots + \frac{f^{(n)}(a)}{n!}(x - a)^n $$The remainder term $ R_n(x) $ represents the error in approximating $ f(x) $ by $ P_n(x) $.Remainder Term:There are several forms of the remainder term, including the Lagrange form and the Peano form.The Lagrange form of the remainder is often used and is given by:$$ R_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!}(x - a)^{n+1} $$where $ \xi$ is some point between $ a $ and $ x $.Therefore, the error is given by:$$ E = f(x) - P_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!} (x - a)^{n+1} $$with the upper bound$$ |E| \leq \max_{\xi \in [a,x]} \frac{|f^{(n+1)}(\xi)|}{(n+1)!} (x-a)^{n+1} $$We can note that：$$ f(x) = P_n(x) + O(h^{n+1}), h = (x-a) $$This is also called the truncation error.Approximating Continuous FunctionsLet $f:[a, b] \to \mathbb{R}$ be differentiable. The Stone-Weierstrass Theorem states that for any small number $\epsilon > 0: \exists P(x)$ a polynomial function such that：$$ \max_{x \in [a, b]} |f(x) - P(x)| \leq \epsilon $$

Numerical Analysis: Error
type: item
date: 2025-03-17T13:56:00+08:00
description: DefinitionWe define the error as:$$E = V_t - V_a$$where $V_t$ is the ground truth value and $V_a$ is our approximation valueThe definition of Relative error, which is normalized error of $E$ by $V_t$:$$E_r = \frac{E}{V_t}$$If $\epsilon$ is an approximation for the error, then the relative approximate error $\epsilon_r$ is defined as:$\epsilon_r = \frac{\epsilon}{V_a}$ExampleIf we use $\frac{22}{7}$ as an approximation of $\pi$$E = \pi - \frac{22}{7}$ $E_r = \frac{E}{\pi} \approx -0.04\%$$\epsilon_r = \frac{\pi - \frac{22}{7}}{\frac{22}{7}} \approx -0.04\%$

一些在做 3D 中遇到的数学知识
type: item
date: 2024-06-28T11:35:00+08:00
description: Zero-Level Set"Zero-level set"（零级集）是数学和几何学中的一个术语，特别是在描述隐式曲面或隐式函数时常用到。在三维空间中，隐式函数通常表示为 $f(\mathbf{x}) = 0$，其中 $\mathbf{x} $ 是空间中的点。这个等式定义了一个曲面或者一个空间中的某种几何体，称为零级集。换句话说，零级集是满足隐式函数 $ f(\mathbf{x}) = 0 $ 的所有点构成的集合，它们组成了隐式函数定义的曲面或边界。例如，在有符号距离函数（Signed Distance Function，SDF）中，零级集表示所有距离函数值为零的点的集合，这些点恰好位于表面上。在几何建模和计算机图形学中，零级集常用于表示曲面的边界或界面，因为它提供了一种便捷的方式来隐式表示和处理复杂的几何结构。

Discretization process for State Space Model of Mamba
type: item
date: 2024-05-31T14:35:00+08:00
description: DefinitionIn Mamba, we have the formula:$$ h'(t)=Ah(t) + Bx(t) $$Here, $A, B$ represent the parametric matrix, and $h$ represents the hidden state, and $x$ represents the input.It is basically from control theory, which is originally:$x'(t) = Ax(t) + Bu(t)$Here we inherit the initial statement of the state space theory, where $x$ represent the hidden state, and $u$ represents the input at this time step.MotivationBy observation, we will realize that:$x'(t) = Ax(t)$which reminds us the $(e^x)' = e^x$.This is where the discretization process from.Move the $Ax(t)$ to the LHS,$$ x'(t) - Ax(t) = Bu(t) $$Let it be:$e^{-At} x'(t) - e^{-At}x(t) = e^{-At} Bu(t)$Note that $F(t) = e^{-At}x(t)$, then $F'(t) = -Ae^{-At}x(t) + e^{-At} x'(t)$Also we know that $F'(t) = e^{-At}Bu(t)$Through the integration, we can get a formular which contains both $u$ and $x$:$$ \forall \lambda \in (-\infty, +\infty), F(t) = F(\lambda) + \int_{\lambda}^{t} F'(\tau) \mathrm{d}\tau $$Specificly, for the convinience, let $\lambda = 0$, substitute $F(t) = e^{-At}x(t)$ and $F'(t)$ into the above equation$$ e^{-At}x(t) = x(0) + \int_{0}^{t} -Ae^{-A\tau}x(\tau) + e^{-A\tau}x(\tau) \mathrm{d} \tau $$It also equals to:$$ e^{-At}x(t) = x(0) + \int_{0}^{t} e^{-A\tau} Bu(\tau) \mathrm{d} \tau $$Divide both sides by $e^{-At}$$$ x(t) = e^{At}x(0) + e^{At}\int_{0}^{t} e^{A\tau} Bu(\tau) \mathrm{d} \tau $$Generally,$$ x(t_k) = e^{At_k} x(0) + e^{At_k} \int_{0}^{t_k} e^{-A\tau} Bu(\tau) \mathrm{d} \tau $$To get the form like $x(t_{k+1}) = ?? x(t_k)$$x(t_{k+1}) = e^{A(t_k + (t_{k+1}-t_k))}x(0) + e^{A(t_k + (t_{k+1} - t_k))} \int_{0}^{t_{k+1}} e^{-A\tau} Bu(\tau) \mathrm{d} \tau$Get the simplication,$$ x(t_{k+1}) = e^{A(t_{k+1} - t_k)} [e^{At_k}x(0) + e^{At_k} \int_{0}^{t_k} -e^{A\tau} Bu(\tau) \mathrm{d}\tau] + e^{At_{k+1}} \int_{t_k}^{t_{k+1}} e^{-A\tau} bu(\tau) \mathrm{d} \tau $$Notice that, the term in $[]$ equals to $x(t_k)$. Hence, it is equivalent to:$$ x(t_{k+1}) = e^{t_{k+1} - t_k} x(t_k) + \int_{t_k}^{t_{k…

Miscellaneous Math for CG 为计算机图形学准备的数学
type: item
date: 2024-04-18T19:26:00+08:00
description: VectorCross Product叉积，简单回顾一下$$ ||\vec{a} \times \vec{b}|| = ||a||||b|| \sin \phi $$而 $\vec{a} \times \vec{b}$ 的方向可以用右手系来判断。在笛卡尔坐标系中，其结果可以用：$$ \vec{a} \times \vec{b} = (y_az_b - z_ay_b, z_ax_b - x_az_b, x_ay_b - y_ax_b) $$Orthonormal BasesConstructing a Basis from a Single Vector利用单个向量构造正交基首先取该向量的方向向量：$$ \vec{w} = \frac{\vec{a}}{||\vec{a}||} $$之后选取任意一个与 $\vec{w}$ 不共线的向量 $\vec{t}$，利用叉积来构造与 $\vec{w}$ 正交的单位向量$$ \vec{u} = \frac{\vec{t} \times \vec{w}}{||\vec{t} \times \vec{w}||} $$Constructing a Basis from a Single Vector将 $\vec{u}$ 构造中的 $\vec{t}$ 替换成第二个向量 $\vec{b}$ 即可。Curves and Surfaces2D Implicit Curves一种常见的描述二维曲线的方法是隐等式（Implicit Equation）。二维曲线的一般形式是：$$ f(x,y) = 0 $$其中， 点 $(x,y)$ 是 $f(x,y)=0$ 的根，表示该点在曲线上。反之则不在。例如，一个圆可以表示为 $f(x,y) = (x-x_c)^2 + (y - y_c)^2 -r^2$ 。因为在曲线上的点并非显式地表达出来，而是要解这个等式，所以成为 Implicit Equation。此外，我们还可以用向量来压缩我们的表达式。如果我们有 $\vec{c} = (x_c, y_c)$ 和 $\vec{p} = (x,y)$，那么以 $(x_c, y_c)$ 为圆心 $r$ 为半径的原则可以用这些满足以下条件的位置向量来表示：$$ (\vec{p} - \vec{c}) \cdot (\vec{p} - \vec{c}) - r^2 = 0 $$由点积的性质可知，该式也等价于：$$ \vert \vert \vec{p} - \vec{c} || -r^2 = 0 $$2D Gradient这一节我们讨论2D梯度假如我们认为 $f(x,y)$ 是一块高地且高度$=f(x,y)$，那么梯度向量指向的是最大向上斜率的方向。梯度向量 $\nabla f(x,y)$ 由下式给出：$$ \nabla f(x,y) = (\frac{\partial f}{\partial x}, \frac{\partial f}{\partial y}) $$在一个在隐曲线 $f(x,y) = 0$ 上的点测量的梯度向量是与 tangent vector 在该点垂直的一个向量。这个垂直的向量通常称为 normal vector （法向量）。因为梯度向量指向的是最大向上斜率的方向，这说明了梯度向量指向了 $f(x,y) > 0$ 的区域。3D Implicit Surfaces类似于 2D 的情形，我们可以简单给出 3D Implicit Equation 来描述一个 3D 平面$$ f(x,y,z) = 0 $$同理，可以用向量来表示。有 $\vec{p} =(x,y,z)$ ，上式可以重写为：$$ f(\vec{p}) = 0 $$Surface Normal to an Implicit SurfaceSurface Normal，即法线。其目的是为了减小计算量，是一个垂直于平面的向量。任意一个在平面上的点都有可能有不同的法向量。类似地，梯度为2D隐曲线指出了一条法线，在点 $\vec{p}$ 处的3D平面发现由隐函数的梯度给出：$$ \vec{n} = \nabla f(\vec{p}) = (\frac{\partial f(\vec{p}}{\partial x}, \frac{\partial f(\vec{p})}{\partial y}, \frac{\partial f(\vec{p})}{\partial z}) $$Implicit Planes上面提到了法线。我们是否可以用法线来定义一个平面？假如我们有一个无限平面，穿过点 $\vec{a}$，且在该点有法线 $\vec{n}$。则该平面的 Implicit Equation 可以表达为：$$ (\vec{p} - \vec{a}) \cdot \ve…

排序不等式
type: item
date: 2024-01-11T10:48:00+08:00
description: 定义对于两个序列 $a_1 \geq a_2 \geq \dots \geq a_n$ 和 $b_1 \geq b_2 \geq \dots \geq b_n$，则有不等式：$\sum\limits_{i=1}^{n} a_ib_i\geq \sum\limits_{i=1}^{n} a_ib_{\pi(i)} \geq \sum\limits_{i=1}^{n} a_1b_{n-i+1}$其中 $\pi(i)$ 是任意一个 $n$ 阶排列。证明假如排列 $\pi(i)$ 使 $\sum\limits_{i=1}^{n} a_{i}b_{\pi(i)}$ 最大，那么 $\pi(i)$ 是一个确定的排列即 $\forall i, \pi(i) = i$。假设该结论不成立，那么 $\pi(i)$ 为一个乱序 $n$ 阶排列且 $\exists i, \pi(i) \neq i$，设 $i$ 是排列 $\pi$ 的第一个出现 $\pi(i) \neq i$ 的位置。$i$ 是第一个这样的位置，那么 $\pi(i) = j > i$，并且一定存在一个 $k >i$ 有 $\pi(k) = i$。因为 $i<j$ 那么 $b_i \geq b_j$，同时有：因为 $i < k$，那么 $a_i \geq a_k$。因此：$$ (a_k - a_j)(b_j - b_i) \geq 0 \Longrightarrow a_ib_i + a_jb_j \geq a_ib_j + a_jb_i $$所以通过使 $\pi(i) = j, \pi(k) = i$ 变成 $\pi(i) = i, \pi(j) = j$， $\sum\limits_{i=1}^{n} a_ib_{\pi(i)}$ 单调不减。证毕。

ARC125B Squares
type: item
date: 2021-10-30T11:02:00+08:00
description: LinkDescription要求： $x^2 - y = z^2$, 满足 $1 \leq x, y \leq N, N \leq 10^{12}$Sol移项可得：$x^2 - z^2 = y$ ，利用平法差公式：$(x+z)(x-z) = y$令 $p = (x+z), q = (x-z)$，因为 $x, y, z \in Z$，所以 $p \geq q$，现在要求 $pq = y$，可以枚举 $q$，然后 $O(1)$ 计算合法的 $p$ 的个数。注意范围：$q \leq p \leq \frac{N}{q}, p \equiv q \mod {2}$Code#include <bits/stdc++.h> #define int long long using namespace std; signed main() { int n, ans = 0; cin >> n; const int mod = 998244353; for (int i = 1, l, r; i * i <= n; ++ i) { l = min(i, n / i), r = max(i, n / i); ans = (ans + ((r - l) / 2ll + 1) % mod) % mod; } printf("%lld\n", ans); return 0; }

Luogu1286 两数之和
type: item
date: 2020-10-17T08:31:00+08:00
description: LinkSol先对给你的这个数组排序，我们设排完序后的数组为$a[1],a[2],\dots,a[n]$，给你的数组升序排序后是$b[1],b[2],\dots,b[n * (n-1) /2]$考虑$b[1]$是怎么被$a[]$相加而得的，因为$b[1]$是$b$中最小的元素，那么$b[1]=a[1]+a[2]$第二小的$b[2]=a[1] + a[3]$由$a[2]$开始往后相加的元素一定在$b[x] = a[1]+a[k]$开始，这样的$b[x]$我们只需要$O(n)$的枚举找到就可以了做法就出来了每次找到$a[k]$这个元素，把$a[1]+a[2], a[1]+a[3],\dots,a[1]+a[k]$这些元素找出来，再把他们在$b$中去掉，再从$a[2]$开始，以此类推时间复杂度$O(n^3)$Code#include <bits/stdc++.h> using namespace std; #define int long long const int SIZE = 1e7 + 5; int n, m, cnt; int a[SIZE], b[SIZE]; std::multiset < int > s; std::vector < vector < int > > ans; inline int read() { char ch = getchar(); int f = 1, x = 0; while (ch < '0' || ch > '9') { if (ch == '-') f = -1; ch = getchar(); } while (ch >= '0' && ch <= '9') { x = (x << 1) + (x << 3) + ch - '0'; ch = getchar(); } return x * f; } inline void solve(int p) { if ((a[1] + a[2] + a[p]) & 1) return; b[1] = (a[1] + a[2] + a[p]) / 2 - a[p]; b[2] = a[1] - b[1], b[3] = a[2] - b[1]; if (b[1] <= 0) return; s.clear(); for (int i = 3; i <= m; ++ i) { if (i != p) s.insert(a[i]); } for (int i = 4; i <= n; ++ i) { b[i] = *s.begin() - b[1]; s.erase(s.begin()); for (int j = 2; j < i; ++ j) { auto it = s.find(b[i] + b[j]); if (it == s.end()) return; s.erase(it); } } std::vector < int > tmp; for (int i = 1; i <= n; ++ i) tmp.push_back(b[i]); ans.push_back(tmp); return; } inline int cmp(vector < int > &a, vector < int > &b) { for (int i = 0; i < n; i++) { if (a[i] == b[i]) continue; return a[i] > b[i]; } } signed main() { while (std::cin >> n) { m = n * (n - 1) / 2; for (int i = 1; i <= m; ++ i) a[i] = read(); std::sort(a + 1, a + m + 1); for (int i = 3; i <= m; ++ i) { if (i == 3 || a[i] != a[i - 1]) solve(i); } if (!ans.size()) { ans.clear(); puts("Impossible"); continue; } std::sort(ans.begin(), ans.end(), cmp); for (int i = 0; i < ans.size(); ++ i) { for (int j = 0; j < n; ++ j) printf("%d ", ans[i][j]); puts(""); } ans.clear(); } return 0; }