数学 - 微分方程数值解 - 第 2 章 常微分方程两点边值问题的差分解法 - 2.2 差分格式

2.2 差分格式

列出几个常用的数值微分公式。

引理 2.2.1

h>0c 为常数

  • 如果 g(x)C2[ch,c+h],则有

(2.2.1)g(c)=12[g(ch)+g(c+h)]h22g(ξ0),ch<ξ0<c+h;

  • 如果 g(x)C2[c,c+h],则有

(2.2.2)g(c)=1h[g(c+h)g(c)]h2g(ξ1),c<ξ1<c+h;

  • 如果 g(x)C2[ch,c],则有

(2.2.3)g(c)=1h[g(c)g(ch)]+h2g(ξ2),ch<ξ2<c;

  • 如果 g(x)C3[ch,c+h],则有

(2.2.4)g(c)=12h[g(c+h)g(ch)]h26g(ξ3),ch<ξ3<c+h;

  • 如果 g(x)C4[ch,c+h],则有

(2.2.5)g(c)=1h2[g(c+h)2g(c)+g(ch)]h212g(4)(ξ4),ch<ξ4<c+h;

  • 如果 g(x)C3[c,c+h],则有

(2.2.6)g(c)=2h[g(c+h)g(c)hg(c)]h3g(ξ5),c<ξ3<c+h;

  • 如果 g(x)C6[ch,c+h],则有

(2.2.7)112[g(ch)+10g(c)+g(c+h)]=1h2[g(c+h)2g(c)+g(ch)]+h4240g(6)(ξ6),ch<ξ6<c+h;

2.2.1 差分格式的建立

(1) 建立差分格式

无特殊说明,均假设所考虑的微分方程定解问题存在具有所需阶数的光滑解。

用有限差分法解两点边值问题的第一步是将求解区间 [a,b] 进行网格剖分。将区间 [a,b]m 等分。

网格步长

h=(ba)/m

网格格点

xi=a+ih,0im

网格

Ωh={xi|0im}

称定义在网格 Ωh 上的函数为网格函数。设 v={vi|0im}Ωh 上的网格函数,记

δxvi12=1h(vivi1),δx2vi=1h(δxvi+12δxvi12),D+vi=1h(vi+1vi),Dvi=1h(vivi1),

现在我们在网格格点上考虑定解问题 (2.1.1)(2.1.2),有

(2.2.8)u(xi)+q(xi)u(xi)=f(xi),1im1

(2.2.9)u(x0)=α,u(xm)=β

定义网格函数 U={Ui|Ui=u(xi),0im},由式 (2.2.5) 可得

u(xi)=1h2[u(xi1)2u(xi)+u(xi+1)]h212u(4)(ξi)=δx2Uih212u(4)(ξi)xi1ξixi

得到

(2.2.10)δx2Ui+q(xi)Ui=f(xi)h212u(4)(ξi),1im1

(2.2.11)U0=α,Um=β

忽略小量项 h212u(4)(ξi),有如下近似等式:

(2.2.12)δx2Ui+q(xi)Uif(xi),1im1

(2.2.13)U0=α,Um=β

在上式中用 ui 代替 Ui,并用 = 代替 ,可得

(2.2.14)δx2ui+q(xi)ui=f(xi),1im1

(2.2.15)u0=α,um=β

将式 (2.2.14) 和式 (2.2.15) 称为求解问题 (2.1.1)(2.1.2)差分格式。建立差分格式的过程称为离散化过程

(2) 局部截断误差

(2.2.14) 和式 (2.2.15) 中共有 m+1 个方程,m+1 个待定量 u0u1um。最后的目标是求出 (u0,u1,,um),然后用求得的 ui 作为 u(xi) 的近似值。

如果我们在 (2.2.14) 式中用精确的 Ui 去代替 ui,则我们只能得到近似等式 (2.2.12),称两边的差为局部截断误差。即下式:

Ri=δx2Ui+q(xi)Uif(xi)

一般来说,Rih 相关,即 Ri=Ri(h),又由 (2.2.10) 可知

Ri=h212u(4)(ξi),xi1<ξi<xi+1

上式表明局部截断误差恰好是建立差分格式时丢掉的小量项,它反映了差分格式 (2.2.14) 对原方程的近似程度。

max1im1|Ri|0,(h0)

当上式成立时,我们称差分格式 (2.2.14)(2.2.15) 与微分方程 (2.1.1)(2.1.2)相容的(相容性)。

差分格式建立后,我们自然会问,差分格式是否存在解(存在性)?差分格式如果存在解,解是否唯一(唯一性)?差分格式的解能否作为微分方程定解问题的近似解(收敛性)?计算过程中的误差对解的影响如何(稳定性)?此外,我们将存在性、唯一性和稳定性统称为适定性。我们可以根据需求在不同的函数空间中提出适当的适定性要求。

2.2.2 差分格式的求解

将式 (2.2.14) 两边同乘以 h2,可得

1h[1h(ui+1ui)1h(uiui1))]×h2+q(xi)uih2=f(xi)h2,1im1

进一步得到

ui1+[q(xi)h2+2]uiui+1=f(xi)h2,1im1

对差分格式表示式 (2.2.14) 和式 (2.2.15),用矩阵可重新表示为

[101q(x1)h2+211q(xm1)h2+2101][u0u1um]=[αh2f(x1)h2f(xm1)β]

上述矩阵为三对角矩阵,可以使用追赶法进行求解。

2.2.3 差分格式解的先验估计式

差分格式解的先验估计式在差分格式的分析中起着关键的作用。

Vh={v|v={vi|0im}Ωh上的网格函数}

V˚h={v|v={vi|0im}Ωh上的网格函数,且v0=vm=0}

vVh,引进如下记号

v=max0im|vi|

v=h(12v02+i=1m1vi2+12vm2)

|v|1=hi=1m(δxvi12)2

v1=v2+|v|12

vvv1 均为 Vh 上的范数,分别成为无穷范数,2 范数和 H1 范数。其中,|v|1Vh 上的半范数,但为 V˚h 上的范数,称为差商的 2 范数。

引理 2.2.2

  • u={ui|0im}Vhv={vi|0im}Vh,则有

(2.2.16)hi=1m1(δx2ui)vi=hi=1m(δxui12)(δxvi12)+D+u0v0Dumvm

  • vV˚h,则有

(2.2.17)hi=1m1(δx2ui)vi=|v|12

  • vV˚h,则有

(2.2.18)vba2|v|1,vba6|v|1

  • vV˚h,则对任意 ε>0

(2.2.19)v2ε|v|12+14εv2,

  • vVh,则

(2.2.20)|v|124h2v2

(2.2.21)v2ε|v|12+(1ε+1ba)v2

根据上述引理中给出的不等式,我们估计差分格式的数值解。

定理 2.2.1

v={vi|0im} 为下述差分格式的解

(2.2.22)δx2vi+q(xi)vi=fi,0im

(2.2.23)v0=0,vm=0

则有

(2.2.24)|v|1ba6f

(2.2.25)v(ba)226f

其中,

f=hi=1m1fi2,f=max1im1|fi|

证明: 将式 (2.2.22) 两边同乘以 hvi,并对 i 求和,得到

(2.2.26)hi=1m1(δx2vi)vi+hi=1m1q(xi)vi2=hi=1m1fivi,

由引理 2.2.2 得到

hi=1m1(δx2ui)vi=|v|12

q(x)0,故

hi=1m1q(xi)vi20

Cauchy-Schwarz 不等式,有

hi=1m1fivihi=1m1fi2hi=1m1vi2=fv

因此有

v12fvfba6|v|1

得到

|v|1ba6f

此外,根据引理 2.2.2 中第三点可以得到

vba2|v|1ba2ba6f(ba)226f

证毕。

2.2.4 差分格式解的存在唯一性

我们给出差分格式解的存在唯一性定理。

定理 2.2.2

差分格式 (2.2.14)(2.2.15) 的解是存在唯一的。

证明:考虑差分格式 (2.2.14)(2.2.15)

δx2ui+q(xi)ui=f(xi),1im1

u0=α,um=β

本质是个线性方程组,且该线性方程组的未知数个数等于方程个数等于 m1,因此只要证明系数矩阵可逆或者证明相应的齐次线性方程组只有零解。

令上式的 f(xi)=0α=0β=0,此时线性方程组右端项为零,由差分解的先验估计式定理 2.2.1u=0,因此齐次线性方程组只有零解,差分格式解存在唯一。

2.2.5 差分格式解的收敛性和稳定性

根据解的先验估计式给出关于差分格式解的收敛性和稳定性一些解释。

(1) 收敛性

定理 2.2.3

{u(x)|axb} 为定解问题 (2.1.1)(2.2.2)的解,{ui|0im} 为差分格式 (2.1.14)(2.1.15) 的解,记

ei=u(xi)ui,0im

则有

eM4(ba)2246h2,M4=maxaxb|u(4)(x)|

证明:

将式 (2.2.10)(2.2.11) 与式 (2.2.14)(2.2.15) 相减得误差方程

δx2ei+q(xi)ei=Ri,1im1

e0=0,em=0

由差分格式先验估计式定理 2.2.1

e(ba)226Ri(ba)226h212maxaxb|u(4)(x)|(ba)2h2246M4

证毕。

定义 2.2.1

如果 e=O(hp),则称差分格式在范数 下是 p 阶收敛的。

由定义知,差分格式 (2.2.14)(2.2.15) 在无穷范数 下是二阶收敛的。

(2) 稳定性

实际计算中,误差不可避免。例如计算 f(xi) 时可能会有一个小的误差 gi,考虑差分格式如下:

(2.2.27)δx2vi+q(xi)vi=f(xi)+gi,1im1

(2.2.28)v0=α,vm=β

v=(v0,v1,,vm) 为上述差分方程的解。

εi=viui,0im

将式 (2.2.27)(2.2.28) 与式 (2.2.14)(2.2.15) 相减,得

(2.2.29)δx2εi+q(xi)εi=gi,1im1

(2.2.30)ε0=0,εm=0

称式 (2.2.29)(2.2.30)摄动方程。它的形式与式 (2.2.14)(2.2.15) 类似。由差分格式先验估计式定理 2.2.1,有

ε(ba)226max1im1|gi|

由上式可以看出当右端项很小时,ε 也很小。由此归纳为如下定理:

定理 2.2.4

差分格式 (2.2.14)(2.2.15) 在下述意义下对右端函数是稳定的:给定差分格式

δx2ui+q(xi)ui=fi,1im1

u0=0,um=0

{ui|0im} 为上述差分格式的解,则有

u(ba)226max1im1|fi|

2.2.6 Richardson 外推法

(1) 示例

设未知量 p 的一个近似式为 p0(h),当 h0 时,若 p0(h)p 之间的存在如下关系

(2.2.31)p0(h)=p+αh2+O(h4)

其中 α 为与 h 无关的非零常数,称上式为 p0(h) 的一个渐进展开式。

p0(h)p=αh2+O(h4)

称上式 p0(h)p 的二阶近似值。

在式 (2.2.31) 中用 h/2 代替 h,得到

(2.2.32)p0(h2)=p+α(h2)2+O((h2)4)

可见 p0(h/2) 也是 p 的一个二阶近似值。现在用 4/3 乘以式 (2.2.32) 得到的式子减去 1/3 乘以式 (2.2.31) 得到的式子,可得

(2.2.33)43p0(h2)13p0(h)=p+O(h4)

p1(h) 作为 p 的又一近似值,可见它比 p0(h/2)p0(h) 具有更高的误差项。

我们把 (2.2.33) 称为 Richardson 外推公式,而把这种从低精度的近似值 p0(h)p0(h/2) 经过线性组合得到高精度近似值的方法称为 Richardson 外推法

(2) 差分方法中的外推法

在差分方法中,外推法也是有效的。记 h 为步长,所得差分格式 (2.2.14)(2.2.15) 的解为 ui(h),0im,我们有下述定理。

定理 2.2.5

给定两点边值问题

(2.2.34)w(x)+q(x)w(x)=112u(4)(x),a<x<b

(2.2.35)w(x)=0,w(b)=0

设上述问题具有光滑解 u(x),则

max0im|u(xi)[43u2i(h2)13ui(h)]|=O(h4)

其中,h 为步长,xi=a+ih

2.2.7 紧差分格式

(1) 紧差分格式的构造

考虑函数 u(x+h) 的展开

u(x+h)=u(x)+u(x)h+u(x)2h2+u(x)6h3+

D 为微分算子,上式可形式写成

u(x+h)=ehDu(x)

由此可得

u(x+h)2u(x)+u(xh)h2=ehD2+ehDh2u(x)=2(D22!+h2D44!+h4D66!+)u(x)=u(x)+h212u(4)(x)+h4360u(6)+

从而得

u(x)=u(x+h)2u(x)+u(xh)h2h212u(4)(x)h4360u(6)=δx2u(xi)h212u(4)(xi)h4360u(6)(xi)+O(h6)

对上式两边求导可得

u(4)(xi)=δx2u(xi)h212u(6)(xi)+O(h4)

由以上两式得到

112[u(xi1)+10u(xi)+u(xi+1)]=δx2u(xi)+h4240u(6)(xi)+O(h6)

w={wi|0im}Vh,定义算子

Awi={112(wi1+10wi+wi+1),i=1,,m1wi,i=0,m

在节点 xi 处考虑微分方程

u(xi)+q(xi)u(xi)=f(xi),i=1,,m1

用算子 A 作用到上式,得

Au(xi)+A[q(xi)u(xi)]=Af(xi),1im1

由于

Au(xi)=δx2u(xi)+h4240u(6)(xi)+O(h6)

因此有

δx2Ui+A[q(xi)Ui]=Af(xi)+h4240u(6)(xi)+O(h6),1im1

忽略小量项得到

δx2Ui+A[q(xi)Ui]Af(xi)

ui 代替 Ui,结合边界条件,得到如下差分格式

(2.2.36)δx2ui+Aq(xi)ui=Af(xi),i=1,,m1

(2.2.37)u0=α,um=β

可以证明差分格式 (2.2.36)(2.2.37) 是唯一可解的,并且上述差分格式在无穷范数下是四阶收敛的且是稳定的。

xi1xixi+1 三点所构造的差分格式中差分格式 (2.2.36)(2.2.37) 的精度达到最高阶 O(h4),称其为紧差分格式

(2) 紧差分格式的求解

对式 (2.2.36) 进行算子展开,得到

ui12ui+ui+1h2+112[q(xi1)ui1+10q(xi)ui+q(xi+1)ui+1]=112[f(xi1)+10f(xi)+f(xi+1)]

关于 ui 进行系数合并

(1+h212q(xi1))ui1+(2+10h212q(xi))ui+(1+h212q(xi+1))ui+1=h212[f(xi1)+10f(xi)+f(xi+1)],i=1,,m1

根据上式以及边界条件可把差分格式写成线性方程组,与之前线性方程组方程形式一致,仍然采用追赶法求解。

posted on   Black_x  阅读(1465)  评论(0编辑  收藏  举报

编辑推荐:
· AI与.NET技术实操系列:基于图像分类模型对图像进行分类
· go语言实现终端里的倒计时
· 如何编写易于单元测试的代码
· 10年+ .NET Coder 心语,封装的思维:从隐藏、稳定开始理解其本质意义
· .NET Core 中如何实现缓存的预热?
阅读排行:
· 分享一个免费、快速、无限量使用的满血 DeepSeek R1 模型,支持深度思考和联网搜索!
· 基于 Docker 搭建 FRP 内网穿透开源项目(很简单哒)
· 25岁的心里话
· ollama系列01:轻松3步本地部署deepseek,普通电脑可用
· 按钮权限的设计及实现
< 2025年3月 >
23 24 25 26 27 28 1
2 3 4 5 6 7 8
9 10 11 12 13 14 15
16 17 18 19 20 21 22
23 24 25 26 27 28 29
30 31 1 2 3 4 5

统计

点击右上角即可分享
微信分享提示