数值分析 $1解方程

2024-03-25 17:27:10

S1 求解方程

C1 二分法

1）原理：零点定理： $f\in C[a,b],f(a)f(b)<0\implies \exist c\in[a,b],f(c) = 0$ f∈C[a,b],f(a)f(b)<0⟹∃c∈[a,b],f(c)=0

2）实现：

% 二分法计算f(x)近似解
% 输入：
% 	f：函数句柄
% 	a，b：区间端点
% 	tol：解精度
% 输出：近似解
function xc=binsect(f, a, b, tol)
if sign(f(a))*sign(f(b)) >= 0
	error('f(a)f(b) greater than zero ') % 错误
end
fa = sign(f(a));
fb = sign(f(b));
while ((b-a) >> 1) > tol
	c = (a+b) >> 1
	fc = sign(f(c))
	if fc == 0
		break
	end
	if fc^fa < 0 % 解落在[a,c]
		b = c; 
	else
		a = c;
	end
end
xc = (a+b) >> 1

3）精度：在 $n$ n次二分操作后误差上限为 $\frac{b-a}{2^{n+1}}$ 2n+1b−a

精确位数：若误差小于 $0.5*10^{-p}$ 0.5∗10−p，则解精确到小数点后p位

C2 不动点迭代

1）原理：求解 $f(x)+x$ f(x)+x或 $x-f(x)$ x−f(x)或其它 $x = g(x)$ x=g(x)的不动点

不动点： $g(x) = x$ g(x)=x
不动点迭代： $x_{n+1} = g(x_n)$ xn+1=g(xn)，若 $x_n$ xn收敛，则收敛到不动点

2）实现：

% 不动点计算g(x) = x的近似解
% 输入：
% 	g：函数句柄
% 	x0：初始估计
% 	k: 迭代次数
% 输出：近似解
function xc=fpi(g, x0, k)
x(1) = x0;
for i=1:k
	x(i+1) = g(x(i))l
end
	xc = x(k+1);

3）线性收敛： $\lim\limits_{k \to +\infty} |\frac{x_{k+1}-r}{x_k-r}| = S \lt 1$ k→+∞lim∣xk−rxk+1−r∣=S<1，收敛速度 $S$ S

$g(x)$ g(x)连续可微， $|g'(r) |<1$ ∣g′(r)∣<1,则可在一个足够接近的初始估计下逼近不动点。(局部收敛)

4）终止条件：不动点迭代的步数是难以确定的，且不一定收敛

绝对终止条件： $|x_{k+1} - x_k|<TOL$ ∣xk+1−xk∣<TOL
相对终止条件： $\frac{|x_{k+1} - x_k|}{|x_{k+1}|}<TOL$ ∣xk+1∣∣xk+1−xk∣<TOL
混合终止条件： $\frac{|x_{k+1} - x_k|}{\max\{\theta,|x_{k+1}|\}}<TOL$ max{θ,∣xk+1∣}∣xk+1−xk∣<TOL，通常是解在0附近

C3 精度的极限

1）后向误差： $|f(\hat{x})|$ ∣f(x^)∣；前向误差： $|\hat{x}-x|$ ∣x^−x∣

2）计算机误差

在有限精度的计算机当中，在根的周围会出现多个假根，以及符号不稳定现象
多重根领域切线水平，一般更容易造成精度丢失

3）根的敏感公式： $f(r) = 0,f(r+\Delta r)+\epsilon g(r+\Delta r)=0\implies \epsilon \ll f'(r),\Delta r = -\frac{\epsilon g(r)}{f'(r)}$ f(r)=0,f(r+Δr)+ϵg(r+Δr)=0⟹ϵ≪f′(r),Δr=−f′(r)ϵg(r)

4）误差放大因子：相对前向误差与相对后向误差之比。 $|\frac{g(r)}{rf'(r)}|$ ∣rf′(r)g(r)∣

C4 牛顿法

1）原理： $x_{k+1} = x_k - \frac{f(x_k)}{f'(x_k)}$ xk+1=xk−f′(xk)f(xk)，使用切线与x轴交点逼近根

2）二次收敛： $\lim\limits_{k \to +\infty} |\frac{x_{k+1}-r}{(x_k-r)^2}| = S \lt +\infty$ k→+∞lim∣(xk−r)2xk+1−r∣=S<+∞

$f$ f二阶连续可导， $f(r) = 0,f'(r)\neq 0$ f(r)=0,f′(r)=0，则牛顿方法局部二次收敛到 $r$ r，迭代误差 $\lim\limits_{k\to+\infty} |\frac{e_{k+1}}{e_k^2}| = \frac{f''(r)}{2f'(r)}$ k→+∞lim∣ek2ek+1∣=2f′(r)f′′(r)
$f$ f在 $[a,b]$ [a,b]上 $m+1$ m+1阶连续可微，且 $r$ r是 $m$ m重根，则牛顿迭代法局部收敛到 $r$ r，迭代误差 $\lim\limits_{k\to+\infty} |\frac{e_{k+1}}{e_k}| = \frac{m-1}{m}$ k→+∞lim∣ekek+1∣=mm−1

3）重根的改进迭代： $x_{k+1} = x_k - \frac{mf(x_k)}{f'(x_k)}$ xk+1=xk−f′(xk)mf(xk) ，局部二次收敛

4）迭代失败：迭代过程中出现 $f'(x_k) = 0$ f′(xk)=0

C5 非求导方法

1）割线替代法：使用差商 $\frac{f(x_i) - f(x_{i+1})}{x_i-x_{i+1}}$ xi−xi+1f(xi)−f(xi+1)替代 $f'(x_i)$ f′(xi)的牛顿迭代法

$x_{i+1} = x_i - \frac{f(x_i)(x_i-x_{i+1})}{f(x_i) - f(x_{i+1})}$ xi+1=xi−f(xi)−f(xi+1)f(xi)(xi−xi+1)
收敛速度介于线性和非线性之间（超线性）

2）试位方法：同二分法类似，属于插值方法

$f(a)f(b) \lt 0,c = a + \frac{f(a)}{f(a)-f(b)}(b-a) = \frac{bf(a) - af(b)}{f(a) - f(b)}$ f(a)f(b)<0,c=a+f(a)−f(b)f(a)(b−a)=f(a)−f(b)bf(a)−af(b)
收敛速度不稳定

3）穆勒方法：计算过 $(x_{k-2},f(x_{k-2})),(x_{k-1},f(x_{k-1})),(x_{k},f(x_{k}))$ (xk−2,f(xk−2)),(xk−1,f(xk−1)),(xk,f(xk))的抛物线 $y=p(x)$ y=p(x)的复根，取接近 $x_k$ xk的为 $x_{k+1}$ xk+1

4）逆二次插值 $IQI$ IQI：计算过 $x_{k-2},x_{k-1},x_k$ xk−2,xk−1,xk的抛物线 $x=p(y)$ x=p(y)与x轴交点 $x_{k+1} = x_k - \frac{r(r-1)(x_k - x_{k-1})+(1-r)s(x_k - x_{k-2})}{(q-1)(r-1)(s-1)},p = \frac{f(x_{k-2})}{f(x_{k-1})},q = \frac{f(x_{k})}{f(x_{k-1})},r = \frac{f(x_{k})}{f(x_{k-2})}$ xk+1=xk−(q−1)(r−1)(s−1)r(r−1)(xk−xk−1)+(1−r)s(xk−xk−2),p=f(xk−1)f(xk−2),q=f(xk−1)f(xk),r=f(xk−2)f(xk)

另一种方式是替换后向误差最大的点

5）Brent方法：在有根区间依次选择使用逆二f次插值，割线方法，二分法，以确保至少使误差减少一半

稳定快速收敛

fzero()函数调用使用的就是Brent

fzero(f, [0 1]) % 寻找[0,1]上的根
fzero(f, 1) % 寻找1附近的根
fzero(f) % 自动寻根
fzero(f, optimset('Display', 'iter')) % 显示迭代情况

码农公寓