数值分析

说明

数值分析

数值积分

\[\int_a^bf(x)dx \approx \sum_{k = 0}^{n}A_kf(x_k) \]

当所有\(A_k > 0\)时，求积公式稳定

---

等分区间的时候，有newton-cotes插值积分公式
这玩意就是等间距插值，然后积分多项式，得到的系数
高阶插值数值不稳定，一般都是用低阶的多项式插值

比如辛普森积分的表： \(\frac{1}{6}\) \(\frac{4}{6}\) \(\frac{1}{6}\)

\[\int_a^bf(x)dx \approx \frac{(b - a)}{6} * (f(a) + 4 * f(\frac{a + b}{2}) + f(b)) \]

多项式插值

任何一个连续函数都可以展开成多项式

拉格朗日插值

\[f(x) = \frac{(x - x_1) * (x - x_2)}{(x_0 - x_1) * (x_0 - x_2)} * f(x_0) + \frac{(x - x_0) * (x - x_2)}{(x_1 - x_0) * (x_1 - x_2)} * f(x_1) + \frac{(x - x_0) * (x - x_1)}{(x_2 - x_0) * (x_2 - x_1)} * f(x_2) \]

差商和牛顿插值

定义差商为 ：\(f[x_i,x_j] = \frac{f(x_i) - f(x_j)}{x_i - x_j}\)
\(f[x_i,x_j，x_k] = \frac{f[x_i,x_j] - f[x_j,x_k]}{x_i - x_k}\)
有\(f[x_i,x_j,x_k] = f[x_k,x_j,x_i]\)(也就是变量的顺序无所谓),并且差商可以表示为f(x_i),f(x_j)...的线性组合
有牛顿插值多项式

\[f(x) = f(x_0) + f[x_0,x_1](x - x_0) + f[x_0,x_1,x_2](x - x_0)(x - x1) + ... + f[x_0,x_1,...,x_n](x - x_0)(x - x_1)...(x - x_{n - 1}) \]

---

等距节点的插值
我们有\(f[x_k,x_{k+1},...,x_{k + m}] = \frac{1}{m!h^m} * \Delta^mf_k\)
然后带入牛顿插值公式，有

\[f(x_0 + th) = c(t,0) * f_0 + c(t,1) * \Delta f_0 + c(n,2) * \Delta f_0 + ... \]

来个栗子
f(x) = x^2

f(2....6) = {4,9,16,25,36}
差分三次
4 9 16 25 36
5 7 9 11
2 2 2
那么

\[f(2 + 1 * t) = c(t,0) * 4 + c(t,1) * 5 + c(t,2) *2 = (t + 2)^2 \]

数值求根

主要解决两个问题
1.什么样的求根算法会收敛
2.收敛速度
3.初始值如何确定

收敛性定理

牛顿迭代法

从几何上看就是求切点与x轴的交点，再从该点继续迭代
假设有f(x),做泰勒展开，有

\[f(x) = \frac{f(x_0)}{0!} + \frac{f'(x_0)(x - x_0)}{1!} + \frac{f''(x_0)(x - x_0)^2}{2!} + ... \]

用一阶泰勒展开来近似，有

\[f(x) \approx f(x_0) + f'(x_0)(x - x_0) = 0 \]

解上述近似方程，得到

\[x = x_0 - \frac{f(x_0)}{f'(x_0)} \]

二分

二分k次后得有根区间\([a_k,b_k]\),令\(x' = \frac{a_k + b_k}{2}\)
那么,

\[误差 = |x - x'| \leq \frac{b_k - a_k}{2} = \frac{b - a}{2^{k + 1}} \]

通过\(\frac{b - a}{2^{k+1}} < \varepsilon\),求得\(k > \frac{Ln(b - a) - Ln(2\varepsilon)}{Ln(2)}\)

函数最佳逼近

有时候需要插值函数在整个区间上和多项式比较逼近
定义\(f(x)\)在区间\([a,b]\)上的p范数为

\[||f(x)||_p = (\int_a^b{|f(x)|}^pdx)^{\frac{1}{p}} \]

对于不同的范数，有不同的最佳逼近理论，我们主要关心\(||f(x)||_2\)和\(||f(x)||_{\infty}\)下的最佳逼近

各种定义

内积：

\[(f,g) = \int_a^b{p(x)f(x)g(x)dx} \]

同欧式空间中的内积一样，有
1.对称性
2.非负性 \((f,f) \geq 0\) ,当且仅当\(f(x) = 0,(f,f) = 0\)
3.齐次性 \((\alpha f,g) = \alpha (f,g)\)
4.分配率 \((f + g,h) = (f,g) + (f,h)\)
还可以定义平方度量

\[||f||_2 = (f,f)^{\frac{1}{2}} \]

有了内积，就可以求两个函数的夹角

\[\theta = arccos\frac{(f,g)}{||f||_2,||g||_2} \]

---

正交函数族
若函数族\(\{\phi(x)\}_0^{\inf} = \{\phi_0(x),\phi_1(x),...,\phi_n(x),...\}\)
满足关系

\[（\phi_i,\phi_j） = \int_a^b{p(x)\phi_i(x)\phi_j(x)dx} = \begin{cases} 0, \qquad \qquad i\neq j \\ A_i > 0, \quad i = j \end{cases} \]

---

只要给定区间和权函数\(p(x)\),就可以对线性无关的函数族\(\{1,x,x^2,...,x^n,...\}\)做施密特正交化变换得到正交多项式序列 $ {\phi(x)}_0^{\inf} $,并具有如下性质 ：

1.\(\phi(x)\)是最高次项系数为1的n次多项式
2.任何次数不超过n的多项式集合均可以表示为\(\{\phi_0(x),\phi_1(x),...,\phi_n(x),...\}\)的线性组合
3.当\(i \neq j\)时，\((\phi_i(x),\phi_j(x)) = 0\),且\(\phi_j(x)\)与任意一次数小于j的多项式正交
4.\(\{\phi_0(x),\phi_1(x),...,\phi_n(x),...\}\)存在递推关系
5.设\(\{\phi(x)\}_0^{\inf}\)是在\([a,b]\)上带权\(p(x)\)的正交多项式序列，则\(\phi_n(x)(n \geq 1)\)的n个根都是区间\((a,b)\)上的单根

---

常用的几个正交多项式
1.Legendre正交多项式
区间为\([-1,1]\),权函数\(p(x) = 1\)时，由\(\{1,x,x^2,...,x^n,...\}\)做施密特正交化变换得到的多项式称为Legendre正交多项式
2.Chebyshev多项式
3.第二类Chebyshev多项式
4.Laguerre多项式
5.Hermite多项式

最佳平方逼近

定理 ：对于任意的连续函数，存在最佳平方逼近多项式
定理 ：对于任意的连续函数，正交函数族有唯一最佳平方逼近

---

最佳平方逼近多项式的求法
设\(\{\phi_i(x)\}\)是标准正交函数族，对于任意连续函数f(x)
最佳平方逼近多项式\(f'(x) = c_0 * \phi_0(x) + c_1 * \phi_1(x) + ... + c_n * \phi_n(x)\)
其中，\(c_i = (f,\phi_i) = \int_a^b{p(x) * f(x) * \phi_k(x)dx}\)

最佳一致逼近

各种定义
在范数(sup为最小上确界)

\[||f||_{\inf} = sup_{a \leq x \leq b}|f(x)| \]

的意义下定义两个函数的距离：

\[d(f,g) = ||f - g||_{\inf}|f(x) - g(x)| \]

---

一致收敛
若一个函数序列\(\{ f_n(x) \}_{n = 1}^{\inf}\)在如上定义的距离的意义下满足

\[lim_{n \rightarrow \inf}||f - f_n|| = 0 \]

则称\(f_n(x)\)在\([a,b]\)上一致收敛于\(f(x)\)

---

定理 ：对于任意的连续函数f(x)存在唯一的最佳一致逼近多项式

---

如何求最佳一致逼近多项式
假设有原函数f(x)和逼近函数f'(x)
记\(偏差 = sup_{a \leq x \leq b}|f(x) - f'(x)|\)
若存在x',使得f'(x) - f(x) = 偏差，那么x'被称为正偏差点
若存在x',使得f'(x) - f(x) = -偏差，那么x'被称为负偏差点
定理 ：f'(x)为f(x)的最佳一致逼近多项式的充分必要条件是
在闭区间[a,b]上存在 f'(x)的阶 + 2 个交错点组
数值算法Remes algorith算法就是基于交错点组的算法

矩阵运算

矩阵计算库

常见的定理

计算矩阵乘法和计算矩阵的逆复杂度同阶，或者说，计算二者的难度一样
通过构造某些常数倍大小的矩阵，即可通过计算矩阵乘法得到逆矩阵，或者反之

LU分解

去掉LUP中选主元的部分即可，这样可能导致除0或者数值不稳定的情况出现

LUP分解

数值稳定性好，速度比高斯消元快3倍，在对同一个矩阵做多次解方程的时候，能减少运算量
也比高斯消元好写...不知道为什么没有流行起来？
LUP分解的复杂度为\(O(n^3)\)，分解好后，解方程只需要\(O(n^2)\)
值得注意的是，L + U - I = LUP分解之后的A矩阵

---

LUP求逆矩阵：
求\(AX = I_n\)
设\(X_i\)为第i列，
那么求逆等价于求n次 \(AX_i = I_i\)

---

LUP求行列式:
\(det(a) * det(P) = det(L) * det(U)\)
注意P是一个置换矩阵，需要求一下逆序对判断行列式的正负，或者在LUP的过程里面维护一下逆序对

---

LUP分解就是找出三个\(n * n\)的矩阵L、U、P满足：\(PA = LU\),其中
L是一个单位下三角矩阵（对角线全是1）
U是一个上三角矩阵
P是一个置换矩阵
可以证明对于任意的非奇异矩阵，都有这样一个分解

---

现在为了求解\(Ax = B\)
\(PAx = PB\)，即\(LUx = PB\)
问题就等效到求两个三角方程组，

\[ \begin{bmatrix} 1 & & \\ l_{21} & 1 & \\ l_{31} * & l_{32} & 1 \\ \end{bmatrix} * \begin{bmatrix} y_1 \\ y_2 \\ y_3 \\ \end{bmatrix} = \begin{bmatrix} b_1 \\ b_2 \\ b_3 \\ \end{bmatrix} \]

显然有\(y_i = b_i- \sum_{j = 1}^{i - 1}{l_{ij} * y_j}\)

\[ \begin{bmatrix} u_{11} & u_{12} & u_{13} \\ & u_{22} & u_{23} \\ & & u_{33} \\ \end{bmatrix} * \begin{bmatrix} y_1 \\ y_2 \\ y_3 \\ \end{bmatrix} = \begin{bmatrix} b_1 \\ b_2 \\ b_3 \\ \end{bmatrix} \]

显然有\(y_i = (b_i- \sum_{j = i + 1}^{n}{u_{ij} * y_j}) / u_{ii}\)

对称正定矩阵与最小二乘逼近

什么是最小二乘逼近
平面上有n个散点\((x_i,y_i)\)

找一条直线去拟合这些点，最小化\(\sum_{i = 1}^{n}{l_i^2}\)
可以找任意多项式去拟合散点，就是平方度量下的最佳逼近

---

设A为一个矩阵
当且仅当\(A = A^T\),且对于任意的\(x \in R^n\)，都有\(x^TAx > 0\)
则A为对称正定矩阵（同样有负定矩阵,即<0）（\(x^TAx \geq 0\)称为半正定矩阵）

---

如何判断是不是正定矩阵？
1、求出A的所有特征值。若A的特征值均为正数，则A是正定的；若A的特征值均为负数，则A为负定的。

2、计算A的各阶主子式。若A的各阶主子式均大于零，则A是正定的；若A的各阶主子式中，奇数阶主子式为负，偶数阶为正，则A为负定的。

---

对称正定矩阵都是非奇异矩阵，可以对他们进行LU分解，而不必担心出现除0的情况
并且对称正定矩阵的每一个主子式都是对称正定的

---

如果A是一个对称正定矩阵，\(A_k\)是A的\(k*k\)主子式，则A关于\(A_k\)的Schur补也是对称正定的

最小二乘逼近

就是对点进行拟合，直接用高次多项式很容易出现过拟合的问题，所以
最小二乘逼近，就是在平方度量下（具体见数值分析），用m次多项式拟合n个点，其中\(m < n\)

---

逼近的正规描述
给定m个点\({(x_1,y_1),(x_2,y_2)...(x_m,y_m)}\),我们希望找出一个函数\(F(x)\)
近似误差 \(\eta_i = F(x_i) - y_i\) 在平方度量下最小，即\(\eta_1^2 + \eta2^2 ... + \eta_m^2\)最小
函数\(F(x)\)的形式依赖于我们所遇到的问题，在此我们假定它是一个线性加权和

\[F(x) = \sum_{j = 1}^{n}{c_j * f_j(x)} \]

我们选择\(f_j(x) = x^{j - 1}\)
也就是说，\(F(x)\)是一个多项式

---

设

\[Ac = \begin{bmatrix} f_1(x_1) & f_2(x_1) & ... & f_n(x_1) \\ f_1(x_2) & f_2(x_2) & ... & f_n(x_2) \\ ... & ... & ... & ... \\ f_1(x_m) & f_2(x_m) & ... & f_n(x_m) \\ \end{bmatrix} \begin{bmatrix} c_1 \\ c_2 \\ ... \\ c_n \\ \end{bmatrix} = \begin{bmatrix} F(x_1)\\ F(x_2)\\ ... \\ F(x_m)\\ \end{bmatrix} \]

那么显然，\(\eta = Ac- y\)
我们的目标是最小化\(||\eta|| = (\sum_{i = 1}^{m}{\eta_i^2})^{\frac{1}{2}}\)
然后我们求\(\frac{d||\eta||^2}{dc_k} = 0\)
然后问题被转化为\((Ac - y)^TA = 0\)
等效于\(c = ((A^TA)^{-1}AT)y = A^+y\)
其中\(A^+\)被称为矩阵A的伪逆矩阵

例子教学

假定有5个点
\((x_0,y_0) = (-1,2)\)
\((x_1,y_1) = (1,1)\)
\((x_2,y_2) = (2,1)\)
\((x_3,y_3) = (3,0)\)
\((x_4,y_4) = (5,3)\)
我们希望用一个一元二次多项式\(F(x) = c_0 + c_1x + c_2x^2\)来拟合
构造：

\[A = \begin{bmatrix} 1 & x_0 & x_0^2 \\ 1 & x_1 & x_1^2 \\ 1 & x_2 & x_2^2 \\ 1 & x_3 & x_3^2 \\ 1 & x_4 & x_4^2 \\ \end{bmatrix} \]

我们有\(c = ((A^TA)^{-1}AT)y = A^+y\)，即

\[c = \begin{bmatrix} 1.200 \\ -0.757 \\ 0.214 \\ \end{bmatrix} \]

那么最终的多项式为\(F(x) = 1.200 - 0.757x + 0.214x^2\)