本文探讨了常微分方程的数值求解方法，包括显式和隐式Euler法，以及刚性问题的挑战。刚性问题是指系统中存在快速和缓慢变化部分，导致需要极小步长以保持稳定性。为了解决刚性问题，文章提到了隐式格式如BDFs的优势。同时，介绍了Krylov子空间的概念及其在大型矩阵求解中的应用，利用Krylov子空间可以降低计算复杂度，避免直接求逆。
常微分方程的数值求解
微分方程整理成关于倒数
y
′
y'
y′的形式
y
′
=
f
(
x
,
y
)
y'=f(x,y)
y′=f(x,y)给定初值，初值问题
y
(
a
0
)
=
y
0
y(a_0)=y_0
y(a0)=y0
离散空间
x
=
x
0
,
x
1
,
x
2
,
.
.
.
x
k
x=x_0,x_1,x_2,...x_k
x=x0,x1,x2,...xk,求离散点上的
y
=
y
0
,
y
1
,
y
2
,
.
.
.
,
y
k
y=y_0,y_1,y_2,...,y_k
y=y0,y1,y2,...,yk值，从而可以用插值形式构建连续函数
y
(
x
)
y(x)
y(x),注意
y
(
x
k
)
≠
y
k
y(x_k)\neq y_k
y(xk)=yk,前者连续，后者离散
显
式
y
k
+
1
=
y
k
+
h
ϕ
(
x
k
,
y
k
,
h
)
显式 y_{k+1}=y_k+h\phi (x_k,y_k,h)
显式yk+1=yk+hϕ(xk,yk,h)
或
隐
射
式
y
k
+
1
=
y
k
+
h
ϕ
(
x
k
,
y
k
,
y
k
+
1
,
h
)
或隐射式 y_{k+1}=y_k + h\phi (x_k,y_k,y_{k+1},h)
或隐射式yk+1=yk+hϕ(xk,yk,yk+1,h)
ϕ
(
x
,
y
,
h
)
\phi (x,y,h)
ϕ(x,y,h)是与导数
f
(
x
,
y
)
f(x,y)
f(x,y)相关的函数，称为增量函数。
-
误差
整体截断误差：求该点的准确值和离散求解值的差 e k = y ( x k ) − y k ~ e_k = y(x_k)-\tilde{y_k} ek=y(xk)−yk~
局部截断误差： 用 y ( x k ) y(x_k) y(xk)这一准确值带入公式推算 y k + 1 y_{k+1} yk+1与 y ( x k + 1 ) y(x_{k+1}) y(xk+1)准确值产生的误差， T k + 1 = y ( x k + 1 ) − [ y ( x k ) + h ϕ ( x k , y ( x k ) , h ) ] T_{k+1}=y(x_{k+1})-[y(x_k)+h\phi (x_k,y(x_k),h)] Tk+1=y(xk+1)−[y(xk)+hϕ(xk,y(xk),h)] -
阶
如果某种方法的局部阶段误差 T k + 1 = O ( h p + 1 ) T_{k+1} = O(h^{p+1}) Tk+1=O(hp+1)，则该方法具有 p p p 阶精度, 或 p阶方法，阶越高越好
主局部截断误差，主项，同阶方法中，主局部截断误差越小越好 -
稳定性
误差累积和传播，误差是否能控制住
绝对稳定：任意点 x k x_k xk计算数值解 y k y_k yk产生的误差总减小， ∣ ∣ ϵ k + 1 ∣ ∣ ≤ ∣ ∣ ϵ k ∣ ∣ ||\epsilon_{k+1}|| \le ||\epsilon_k|| ∣∣ϵk+1∣∣≤∣∣ϵk∣∣, ϵ k = y k − y k ~ \epsilon_k=y_k-\tilde{y_k} ϵk=yk−yk~
采取试验方程分析方法的绝对稳定性 f ( x , y ) = y ′ = λ y , ( λ 复 数 ) f(x,y)=y'=\lambda y,(\lambda复数) f(x,y)=y′=λy,(λ复数)
稳定性往往和步长 h h h 有关
绝对稳定域包含复平面区域，稳定区间是与复平面实轴的交点
唯一解的Lips条件: ∣ ∣ f ( x k + 1 , y k + 1 ) − f ( x k , y k ) ∣ ∣ ≤ L ∣ ∣ y k + 1 − y k ∣ ∣ , L > 0 ||f(x_{k+1},y_{k+1})-f(x_k,y_k)||\le L||y_{k+1}-y_k||, L>0 ∣∣f(xk+1,yk+1)−f(xk,yk)∣∣≤L∣∣yk+1−yk∣∣,L>0
Euler法
改进Euler法（梯形法）
预测-矫正法
R-K类方法
什么是常微分方程组的刚性（stiffness）？
在常微分系统中包含多个组成部分，包括变化较快的部分以及变化较慢的部分，其特征时间差异巨大，使得在采用经典的Euler法求解时，采用很小的步长保证稳定性 ∣ 1 + λ h ∣ ≤ 1 |1 + \lambda h | \le 1 ∣1+λh∣≤1 ,这使得求解计算步数显著增加，甚至加剧误差的累计，可以说这样的系统具有刚性。
解决刚性问题的方法
-
implicit Euler， BDFs (back difference )
隐式格式，具有无条件稳定性，远好于显示方法
y k + 1 = y k + h f ( x k + 1 , y k + 1 ) y_{k+1}=y_k +hf(x_{k+1},y_{k+1}) yk+1=yk+hf(xk+1,yk+1) -
Gear？
-
VODE？
-
Rosenbrock？
Krylov 子空间
有向量 b b b和矩阵 A A A,则通过 A b Ab Ab生产的向量，有 K m K_m Km = { b , A b , A 2 b , A 3 b , . . . , A m b b,Ab, A^2b,A^3b,...,A^mb b,Ab,A2b,A3b,...,Amb},可以张成Krylov子空间，该子空间的任意向量可以通过Krylov矩阵 K m K_m Km 中列向量的线性组合得到。
有何用处
对一个矩阵求解问题而言 A x = b Ax=b Ax=b x = A − 1 b x=A^{-1}b x=A−1b 当矩阵较大时，A的求逆会特别耗时 O ( n 2 ? 3 ) O(n^{2?3}) O(n2?3),而Krylov子空间方法可以降维处理问题，避免求逆的计算量：根据 Krylov方法 A − 1 b ≈ ∑ ( β i A i b ) , i = 0 , 2 , . . m − 1 A^{-1}b ≈ \sum(\beta_iA^ib),i=0,2,..m-1 A−1b≈∑(βiAib),i=0,2,..m−1
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
