二维Burgers方程的分裂高阶有限差分方法
The Splitting High-Order Finite Difference Method of Two-Dimensional Burgers Equation
DOI: 10.12677/PM.2021.111004, PDF, HTML, XML,  被引量 下载: 489  浏览: 801  科研立项经费支持
作者: 马佳琪, 王 博:中国民航大学理学院,天津
关键词: Burgers方程有限差分方法收敛性稳定性Burgers Equations Finite Difference Method Convergence Stability
摘要: 本文主要研究了二维Burgers方程的分裂高阶差分方法,利用能量的方法,证明了所提出的差分格式在时间上具有二阶收敛率以及在空间上具有四阶收敛率。数值结果验证了算法的精度和有效性。
Abstract: The study is concerned with the splitting higher-order finite difference method of two-dimensional Burgers equation. By using the energy method, the proposed difference scheme is proved to have the second-order convergence rate in time and the fourth-order convergence rate in space. The accuracy and efficiency of the algorithm are verified by numerical findings.
文章引用:马佳琪, 王博. 二维Burgers方程的分裂高阶有限差分方法[J]. 理论数学, 2021, 11(1): 22-31. https://doi.org/10.12677/PM.2021.111004

1. 引言

在如今这个科技快速发展的情势下,偏微分方程的理论形成和发展在数学及科学技术中应用广泛,Burgers方程作为其重要分支,应用于数学及社会生活的许多领域,它作为最简单的非线性对流扩散的模型,在湍流、传热、传质、大气、水资源污染及连续随机过程等众多领域中发挥着极其重要的作用。对Burgers方程的深入研究有助于其他非线性问题 [1] 的研究,所以求解非线性Burgers方程的数值解十分有意义。

近年来,关于Burgers方程的数值研究有不少的成果 [2] [3] [4]。文献 [3] 研究了一维Burgers方程的两种线性化差分格式,文献 [1] [3] 研究了二维Burgers方程的Crank-Nicolson格式,并得到了时间二阶和空间二阶的收敛速度。本文参考文献 [5] 中的技巧,采用维数分裂的技巧,提出二维Burgers方程的分裂高阶有限差分方法,并对其稳定性以及收敛性进行理论性分析,最后通过数值算例进行验证。

2. 分裂高阶有限差分方法

考虑如下二维Burgers方程

u t + u u x + u u y = V 2 u x 2 + V 2 u y 2 , (a)

其初边值条件为

u ( x , y , 0 ) = u ( x , y ) , ( x , y ) Ω , (b)

u ( x , y , t ) = 0 , ( x , y ) Ω , t [ 0 , T ] . (c)

其中矩形区域 Ω = ( 0 , a ) × ( 0 , b ) ,时间周期为 T > 0

2.1. 记号

为建立差分格式,首先对求解区域

Ω = { ( x , y , t ) | 0 x I , 0 y J , 0 t T } 做网格剖分。取x-和y-方向上的空间步长分别为 h x = a I h y = b J ,时间步长为 τ = T N ,其对应网格节点坐标为 u ( x i , y j , t n ) x i = i h x i = 0 , 1 , , I y j = j h y j = 0 , 1 , , J t n = n τ n = 0 , 1 , , N 。记 Ω h x = { x i | 0 x I } Ω h y = { y j | 0 y J } Ω τ n = { τ n | 0 τ N } Ω h τ = Ω h x × Ω h y × Ω τ n

{ u i j n | 0 i I , 0 j J , 0 n N } 为一个网格函数。下面,我们定义一些差分算子:

δ x u i j n = u i + 1 , j n u i j n h x , u i j n = u i , j + 1 n u i j n h y ,

δ t u i j n = u i j n + 1 u i j n τ ,

δ x ^ u i j n = u i + 1 , j n u i 1 , j n 2 h x , δ y ^ u i j n = u i , j + 1 n u i , j 1 n 2 h y ,

δ x 2 u i j n = u i + 1 , j n 2 u i j n + u i 1 , j n h x 2 , δ y 2 u i j n = u i , j + 1 n 2 u i j n + u i , j 1 n h y 2 .

根据泰勒展开式所得 为二阶,为提升至四阶,引入高阶紧算子,

L x u i j n = 1 12 ( u i + 1 , j n + 10 u i j n + u i 1 , j n ) ,

L y u i j n = 1 12 ( u i , j + 1 n + 10 u i j n + u i , j 1 n ) ,

M x u i j n = 1 6 ( u i + 1 , j n + 4 u i j n + u i 1 , j n ) ,

M y u i j n = 1 6 ( u i , j + 1 n + 4 u i j n + u i , j 1 n ) ,

满足

L x 2 u i j n x 2 = δ x 2 u i j n + ο ( h x 4 ) ,

L y 2 u i j n x 2 = δ y 2 u i j n + ο ( h y 4 ) ,

M x u i j n x = δ x ^ u i j n + ο ( h x 4 ) ,

M y u i j n x = δ y ^ u i j n + ο ( h y 4 ) .

2.2. 差分格式的建立

根据抛物方程的分裂技巧,在时间步 ( t n , t n + 1 ] 下,将方程(a)分裂成两个方程

u t + u u x = V 2 u x 2 , (2.2.1)

u t + u u y = V 2 u y 2 . (2.2.2)

U i j n u i j n 的近似解,对于方程(2.2.1)~(2.2.2),我们提出分裂差分格式:

U i j n + 1 2 U i j n τ + 1 3 ϕ ( U i j n + 1 4 ) = V L x 1 δ x 2 U i j n + 1 4 , (2.2.3)

U i j n + 1 U i j n + 1 2 τ + 1 3 ψ ( U i j n + 3 4 ) = V L y 1 δ y 2 U i j n + 3 4 . (2.2.4)

其中

ϕ ( U i j n + 1 4 ) = U i j n + 1 4 M x 1 δ x ^ U i j n + 1 4 + M x 1 δ x ^ ( U i j n + 1 4 ) 2

ψ ( U i j n + 3 4 ) = U i j n + 3 4 M y 1 δ y ^ U i j n + 3 4 + M y 1 δ y ^ ( U i j n + 3 4 ) 2 .

格式的边界条件为

U ( x , y , t ) = 0 , ( x , y ) Ω , t [ 0 , T ] , (2.2.5)

初值条件为

U i j 0 = u ( x , y ) . (2.2.6)

2.3. 守恒性和收敛性分析

为了方便稳定性和收敛性的证明,我们给出一些内积和范数的定义。令 V h = { u i j | ( x i , y j ) Ω ¯ h , u i j | Ω h = 0 } ,我们定义离散的内积和范数

( u , v ) = i = 1 I 1 j = 1 J 1 u i j v i j h x h y ,

( u , v ) x = i = 0 I 1 j = 1 J 1 u i j v i j h x h y ,

( u , v ) y = i = 1 I 1 j = 0 J 1 u i j v i j h x h y ,

u = ( u , u ) , u x = ( u , u ) x , u y = ( u , u ) y .

为了证明紧格式的离散电荷和能量的守恒性,我们给出以下一些引理。

引理2.3.1 [5] 对于网格函数 u , ν V h ,我们有

( δ x 2 u , ν ) = ( δ x u , δ x ν ) x , ( δ y 2 u , ν ) = ( δ y u , δ y ν ) y .

引理2.3.2 [5] 对于网格函数 μ , ν V h ,我们有

( L x u , ν ) = ( u , L x ν ) , ( M x u , ν ) = ( u , M x ν ) ,

( L x 1 u , ν ) = ( u , L x 1 ν ) , ( M x 1 u , ν ) = ( u , M x 1 ν ) .

引理2.3.3 [5] 对于网格函数 U , V V h ,若 U L 2 V L 2 有界,则存在正常数C,使得

( ϕ ( U ) ϕ ( V ) , U V ) C ( U V 2 + δ x ( U V ) x 2 ) .

引理2.3.4对于对称的正定矩阵 L x 1 M x 1 ,我们有

L x 1 = R 1 T R 1 , M x 1 = R 2 T R 2 .

C 0 U n R 1 U n C 1 U n , C 3 U n R 2 U n C 4 U n

其中 R 1 = C h o l ( L x 1 ) R 2 = C h o l ( M x 1 ) C i ( i = 1 , 2 , 3 , 4 ) 为正常数。

证:因为 L x 1 是一个对称的正定矩阵,故存在矩阵 R 1 ,其中 R 1 为上三角矩阵,满足 L x 1 = R 1 T R 1 ,其中 R 1 = C h o l ( L x 1 )

根据矩阵中范数的定义,我们有

C 0 U n R 1 U n C 1 U n

同样,矩阵 M x 1 R 2 也有类似的结果。

引理2.3.5 [5] 对于网格函数 u V h ,我们有

δ x u n x 2 4 h x 2 u n 2 , δ y u n y 2 4 h y 2 u n 2 .

引理2.3.6 [5] (Gronwall不等式)设 w n 0 且存在常数 C 3 , C 4 ,使得 w n C 3 + C 4 τ k = 0 n 1 w k 0 < n τ < T

那么

w n C 3 e C 4 T , 0 < n τ < T .

定理2.3.1 (守恒性)分裂高阶有限差分格式(2.2.3)~(2.2.6)满足离散的电荷守恒,即

Q n + 1 = Q n . (2.3.1)

其中

Q n = i = 1 I 1 j 1 J 1 U i j n h x h y . (2.3.2)

对任意的 n 0 成立。

证:对方程(2.2.3)两端同时乘以 h x h y ,i从1到 I 1 求和,j从1到 J 1 求和:

h x h y i = 1 I 1 j = 1 J 1 U i j n + 1 2 U i j n τ + 1 3 h x h y i = 1 I 1 j = 1 J 1 ϕ ( U i j n + 1 4 ) V h x h y i = 1 I 1 j = 1 J 1 L x 1 δ x 2 U i j n + 1 4 = 0

考虑边界条件(2.2.5),有

U 0 j n = U 1 j n = U i 0 n = U i 1 n = U I 1 , j n = U I j n = U i , J 1 n = U i J n = 0

h x h y i = 1 I 1 j = 1 J 1 L x 1 δ x 2 U i j n + 1 4 = h x h y i = 1 I 1 j = 1 J 1 δ x 2 U i j n + 1 4 = h x h y i = 1 I 1 j = 1 J 1 U i + 1 , j n + 1 4 2 U i j n + 1 4 + U i 1 , j n + 1 4 h x 2 = h x h y i = 1 I 1 j = 1 J 1 [ ( U i + 1 , j n + 1 4 U i j n + 1 4 ) ( U i j n + 1 4 U i 1 , j n + 1 4 ) ] = 0

因为

h x h y i = 1 I 1 j = 1 J 1 ϕ ( U i j n + 1 4 ) = h x h y i = 1 I 1 j = 1 J 1 ( U i j n + 1 4 M x 1 δ x ^ U i j n + 1 4 + M x 1 δ x ^ ( U i j n + 1 4 ) 2 ) = h x h y i = 1 I 1 j = 1 J 1 ( M x 1 U i j n + 1 4 δ x ^ U i j n + 1 4 ) + h x h y i = 1 I 1 j = 1 J 1 δ x ^ ( U i j n + 1 4 ) 2 = h x h y ( i = 1 I 1 j = 1 J 1 U i + 1 , j n + 1 4 ) ( i = 1 I 1 j = 1 J 1 δ x ^ U i j n + 1 4 ) + h x h y i = 1 I 1 j = 1 J 1 δ x ^ ( U i j n + 1 4 ) 2 = 0

所以

h x h y i = 1 I 1 j = 1 J 1 U i j n + 1 2 U i j n τ = 0

h x h y i = 1 I 1 j = 1 J 1 U i j n + 1 2 = h x h y i = 1 I 1 j = 1 J 1 U i j n

Q n + 1 = Q n ,因此(2.2.3)~(2.2.6)满足电荷守恒。

推论 (稳定性)差分格式(2.2.3)~(2.2.6)是无条件稳定的。

证:由定理2.3.1可得:对于任意的 n 0 ,有

U n + 1 = U n = U n 1 = = U 1 = U 0 成立

因此,差分格式(2.2.3)~(2.2.6)是无条件稳定的。

定理2.3.2 (收敛性)假设问题(a)~(c)的精确解ν足够光滑,令u是分裂高阶有限差分格式(2.2.3)~(2.2.6)的数值解,对于固定的 T > 0 ,存在一个独立的 τ h x h y 和正常数C,满足:

max 1 n N ν n u n C ( τ 2 + h x 4 + h y 4 ) .

证:令截断误差

ξ i j n = ν i j n + 1 2 ν i j n τ + 1 3 ϕ ( ν i j n + 1 4 ) V L x 1 δ x 2 ν i j n + 1 4 ,

η i j n = ν i j n + 1 ν i j n + 1 2 τ + 1 3 ψ ( ν i j n + 3 4 ) V L y 1 δ y 2 ν i j n + 3 4 .

由(2.2.3)易得:令 e i j n = ν i j n u i j n

有:

ξ i j n = e i j n + 1 2 e i j n τ + 1 3 ϕ ( ν i j n + 1 4 ) 1 3 ϕ ( u i j n + 1 4 ) V L x 1 δ x 2 e i j n + 1 4 .

ξ i j n e i j n + 1 2 + e i j n 作内积:

( ξ n , e n + 1 2 + e n ) = 1 τ ( e n + 1 2 2 e n 2 ) + 1 3 ( ϕ ( ν i j n + 1 4 ) 1 3 ϕ ( u i j n + 1 4 ) , e n + 1 2 + e n ) V ( L x 1 δ x 2 e n + 1 4 , e n + 1 2 + e n ) .

由引理2.3.3:

ϕ ( ν i j n + 1 4 ) ϕ ( u i j n + 1 4 ) C ( e n + 1 2 2 + e n 2 + δ x e n + 1 2 x 2 + δ x e n x 2 ) . (2.3.3)

根据引理2.3.4:

| L x 1 δ x 2 e n + 1 4 , e n + 1 2 + e n | = | R 1 δ x 2 e n + 1 4 , R 1 ( e n + 1 2 + e n ) | = | [ R 1 δ x e n + 1 4 , R 1 δ x ( e n + 1 2 + e n ) x ] | = 1 2 R 1 δ x ( e n + 1 2 + e n ) x 2 1 2 δ x ( e n + 1 2 + e n ) x 2 δ x e n + 1 2 x 2 + δ x e n x 2 , (2.3.4)

又:

( ξ n , e n + 1 2 + e n ) ξ n 2 + 2 e n + 1 2 2 + 2 e n 2 , (2.3.5)

将(2.3.3)~(2.3.5)代入,则有:

1 τ ( e n + 1 2 2 e n 2 ) ξ n 2 + 2 e n + 1 2 2 + 2 e n 2 + C e n + 1 2 2 + C e n 2 + C δ x e n + 1 2 x 2 + C δ x e n x 2 + V δ x e n + 1 2 x 2 + V δ x e n x 2 ,

所以:

e n + 1 2 2 e n 2 τ ξ n 2 + τ ( 2 + C ) e n + 1 2 2 + τ ( 2 + C ) e n 2 + τ ( V + C ) δ x e n + 1 2 x 2 + τ ( V + C ) δ x e n x 2 ,

由引理2.3.5得:

e n + 1 2 2 τ ( 2 + C + 4 h x 2 ( V + C ) ) e n + 1 2 2 e n 2 + τ ( 2 + C + 4 h x 2 ( V + C ) ) e n 2 + τ ξ n 2 ,

即:

[ 1 τ ( 2 + C + 4 h x 2 ( V + C ) ) ] e n + 1 2 2 [ 1 + τ ( 2 + C + 4 h x 2 ( V + C ) ) ] e n 2 + τ ξ n 2

我们假设 τ = o ( h x 2 ) ,有 1 τ ( 2 + C + 4 h x 2 ( V + C ) ) > 0

则:

e n + 1 2 2 1 + τ ( 2 + C + 4 h x 2 ( V + C ) ) 1 τ ( 2 + C + 4 h x 2 ( V + C ) ) e n 2 + τ 1 τ ( 2 + C + 4 h x 2 ( V + C ) ) ξ n 2

P 1 = 2 + C + 4 h x 2 ( V + C ) ,则有:

e n + 1 2 2 1 + τ P 1 1 τ P 1 e n 2 + τ 1 τ P 1 ξ n 2 (2.3.6)

类似地:由方程(2.2.4),令 P 2 = 2 + C + 4 h y 2 ( V + C ) ,有:

e n + 1 2 1 + τ P 2 1 τ P 2 e n + 1 2 2 + τ 1 τ P 2 η n 2 (2.3.7)

由(2.3.6),(2.3.7)得到:

e n + 1 2 ( 1 + P 1 τ ) ( 1 + P 2 τ ) ( 1 P 1 τ ) ( 1 P 2 τ ) e n 2 + τ 1 P 2 τ η n 2 + τ 1 P 1 τ 1 + P 2 τ 1 P 2 τ ξ n 2

e n + 1 2 ( 1 + C 1 τ ) e n 2 + C 2 τ ( ξ n 2 + η n 2 ) ( 1 + C 1 τ ) e n 2 + C 2 τ ( τ 2 + h x 4 + h y 4 ) 2

e 0 = 0 ,由Gronwall不等式:

max 1 n N ν n u n C ( τ 2 + h x 4 + h y 4 ) .

3. 数值算例

通过算例验证数值结果,考虑二维Burgers方程

u t + u u x + u u y = V 2 u x 2 + V 2 u y 2

区域为 [ 0 , 1 ] × [ 0 , 1 ] ,精确解为

u ( x , y , t ) = exp ( 2 π 2 t ) sin ( π x ) cos ( π y )

初始条件和边界条件有精确解直接得到。

表1表2列出了分裂高阶有限差分格式具有空间四阶的收敛率、时间二阶的收敛率,这个结果和我们的理论分析是一致的。

Table 1. When t = 10 6 , space step error and convergence order, where τ = 10 7 h x = h y = h

表1. 当 t = 10 6 时,空间步长的误差和收敛阶,其中 τ = 10 7 h x = h y = h

Table 2. When t = 0.1 , space step error and convergence order, where h x = h y = 0.1

表2. 当 t = 0.1 时,时间步长的误差和收敛阶,其中 h x = h y = 0.1

我们考虑下方的精确解

u ( x , y , t ) = exp ( k π 2 t ) sin ( π x ) cos ( π y )

表3可以看出,当k的值发生变化时,时间步长和空间步长分别为 h x = h y = 0.2 τ = 10 7 ,我们提出的差分格式的电荷的误差足够小,并且始终为正值,因此差分格式满足电荷守恒,也验证了定理2.3.1的分析结果。

Table 3. When h x = h y = 0.2 τ = 10 7 , space error and convergence order of different exact solutions

表3. 当 h x = h y = 0.2 τ = 10 7 时,不同精确解的时间步长的误差和收敛阶

图1图2分别展示了 τ = 10 7 , h x = h y = 0.05 时,差分格式的数值解和精确解的图像。可以看出,精确解和数值解的三维立体图是完全一致的,而且曲面的几何特性也是完全相同的。而且我们发现,当空间步长 h x , h y 越小时,网格剖分越密集,Burgers方程的解模拟的更加精密。

Figure 1. When τ = 10 7 , h x = h y = 0.002 , numerical solutions of difference schemes

图1. τ = 10 7 h x = h y = 0.002 时,差分格式的数值解

Figure 2. When τ = 10 7 , h x = h y = 0.002 , the exact solutions of difference schemes

图2. τ = 10 7 h x = h y = 0.002 时,差分格式的精确解

4. 结论

通过抛物型方程的高精度分裂方法 [6] [7],我们提出了二维Burgers方程的高阶分裂有限差分方法。借用能量法,我们证明了所提格式是守恒的,并且满足时间二阶、空间四阶的收敛性。最后,数值算例和理论结果一致。

基金项目

大学生创新创业训练计划项目(项目编号202010059074)。

参考文献

[1] 陈景良, 邓定文. 非线性延迟波动方程的两类差分格式[J]. 理论数学, 2020, 10(5): 508-517.
[2] 黄展鹏. 解Burgers方程的修正局部Crank-Nicolson方法[D]: [硕士学位论文]. 乌鲁木齐: 新疆大学, 2009.
[3] Sun, H. and Sun, Z.Z. (2015) On Two Linearized Difference Schemes for Burgers Equation. International Journal of Computational Methods, 92, 1160-1179.
https://doi.org/10.1080/00207160.2014.927059
[4] Xu, P.P. and Sun, Z.Z. (2009) A Second Order Accurate Difference Scheme for the Two-Dimensional Burgers System. Numerical Methods for Partial Differential Equations, 25, 172-194.
https://doi.org/10.1002/num.20335
[5] Wang, B., Sun, T.J. and Liang, D. (2019) The Conservative and Fourth-Order Compact Finite Difference Schemes for Regularized Long Wave Equation. Journal of Computational and Applied Mathematics, 356, 98-117.
https://doi.org/10.1016/j.cam.2019.01.036
[6] Wang, B., Liang, D. and Sun, T.J. (2017) The Conservative Split-ting High-Order Compact Finite Difference Scheme for Two-Dimensional Schrodinger Equations. International Journal of Computational Methods, 15, Article ID: 1750079.
https://doi.org/10.1142/S0219876217500797
[7] 陈宁, 顾海明. 抛物型方程的高精度交替方向法[J]. 青岛科技大学学报(自然科学版), 2009, 30(1): 3-5.