郭禎 吳明明 胡勁松

本文對(duì)含齊次邊界條件的KdV方程的初邊值問(wèn)題進(jìn)行了數(shù)值研究. 通過(guò)在時(shí)間層進(jìn)行二階精度的Crank-Nicolson差分離散、在空間層進(jìn)行六階精度的外推組合差分離散，本文建立了一個(gè)具有六階空間精度的兩層非線性差分格式.該格式能夠合理地模擬原問(wèn)題的兩個(gè)守恒量. 然后，本文利用能量方法證明了格式的收斂性和穩(wěn)定性. 數(shù)值算例驗(yàn)證了該方法的有效性.

KdV方程； Crank-Nicolson差分格式； 六階精度； 守恒

O241.82A2023.031006

收稿日期： 2022-06-30

基金項(xiàng)目： 國(guó)家自然科學(xué)基金青年基金（11701481）；四川應(yīng)用基礎(chǔ)研究項(xiàng)目（2019JY0387）

作者簡(jiǎn)介： 郭禎（1997-）， 女， 碩士研究生， 主要研究方向?yàn)橛?jì)算數(shù)學(xué). E-mail： gz2080118252@163.com

通訊作者： 胡勁松. E-mail： hjs888hjs@163.com

A conservative difference scheme with 6-order spatial accuracy for the KdV equation

GUO Zhen， WU Ming-Ming， HU Jin-Song

（School of Science， Xihua University， Chengdu 610039， China）

In this paper， we propose a two-level difference scheme with 6-order spatial accuracy for the initial boundary value problem of KdV equation with homogeneous boundary condition. In this scheme， the Crank-Nicolson differential dispersion with 2-order accuracy is used in the time layer and the discretization of space layer is performed by extrapolating difference combination with 6-order accuracy. This scheme can simulate two conservative properties of the original problem reasonably. Then the convergence and stability of the scheme are proved by using the energy method. Finally， numerical examples verify the performance of the scheme.

KdV equation; Crank-Nicolson difference scheme; 6-order accuracy; Conservation

（2010 MSC 65M60）

1 引 言

KdV方程［1］

ut+αuux+uxxx=0（1）

（其中α為非零實(shí)常數(shù)）是荷蘭數(shù)學(xué)家Korteweg和de Vires在研究淺水中的小振幅長(zhǎng)波運(yùn)動(dòng)時(shí)提出的. 因具有無(wú)窮多個(gè)不變量，該方程在冷等離子體中的離子聲波、巖漿流與管道波、泡沫液體中的聲波及海洋內(nèi)波等諸多研究領(lǐng)域有廣泛應(yīng)用［2］. 作為一類典型的非線性色散波動(dòng)方程，KdV方程少有解析解，這就使得對(duì)其數(shù)值解法的研究具有重要價(jià)值［3-12］.

本文考慮如下KdV方程的初邊值問(wèn)題：

ut+αuux+uxxx=0，

x， t∈xL， xR×0， T（2）

ux， 0=u0x，x∈xL， xR（3）

uxL， t=uxR， t=0，

uxxR， t=0，t∈0， T（4）

郭 禎， 等： KdV方程的一個(gè)六階空間精度守恒差分格式

其中u0（x）是一個(gè)已知函數(shù). 初邊值問(wèn)題（2）～（4）滿足如下守恒律［13］：

Qt=∫xRxLux， tdx=

∫xRxLu0x， tdx=Q0（5）

Et=u2L2=u02L2=E0（6）

其中Q（0）與E（0）均為僅與初始條件有關(guān)的常數(shù). 該方程是奇數(shù)階的，且沒有對(duì)稱性，故其數(shù)值求解有相當(dāng)難度. 在已有的求解方法中，文獻(xiàn)［9］提出了精度為Oτ2+h的兩層非線性差分格式和三層線性差分格式，文獻(xiàn)［10］構(gòu)造了具有二階精度的兩層非線性差分格式.盡管這些差分格式能夠模擬KdV方程的部分守恒律，但理論精度較低. 值得注意的是，為了提高數(shù)值解法的理論精度，文獻(xiàn)［11］引入了中間函數(shù)，將方程轉(zhuǎn)化為含有兩個(gè)未知函數(shù)的方程組，進(jìn)而提出了一個(gè)在空間層具有六階精度的緊致差分格式. 然而，該方法需要求解兩個(gè)未知函數(shù)，計(jì)算量大. 此外，文獻(xiàn)［12］利用正交分解方法建立了一個(gè)在空間層具有六階精度的隱式緊致差分格式，但僅驗(yàn)證了方法的可行性而沒有對(duì)方程的守恒律進(jìn)行模擬.

在本文中， 我們對(duì)初邊值問(wèn)題（2）～（4）建立了一個(gè)具有六階空間精度的兩層非線性差分格式，給出其收斂性和穩(wěn)定性的證明. 我們的格式能夠合理地模擬守恒量（5），（6）. 該格式首先在時(shí)間層進(jìn)行二階精度的Crank-Nicolson差分離散，然后在空間層進(jìn)行六階精度的外推組合差分離散. 最后我們通過(guò)數(shù)值算例驗(yàn)證了格式的有效性.

2 差分格式及其守恒性

記剖分區(qū)域xL， xR×0， T，時(shí)間步長(zhǎng)為τ，tn=nτ，0≤n≤N， N=Tτ；xj=xL+jh，0≤j≤J，h=xR-xLJ為空間步長(zhǎng). 我們約定C是與空間步長(zhǎng)和時(shí)間步長(zhǎng)均無(wú)關(guān)一般常數(shù)且C>0，并定義如下記號(hào)

unj≡u(píng)xj， tn，? Unj≈uxj， tn，

Unjx=Unj+1-Unjh， Unjx-=Unj-Unj-1h，

Unjx^=Unj+1-Unj-12h， Unjx¨=Unj+2-Unj-24h，

Unjx…=Unj+3-Unj-36h， Unjt=Un+1j-Unjτ，

Un+1/2j=Un+1j+Unj2，

〈Un，Vn〉=h∑J-1j=1UnjVnj，

Un2=〈Un，Un〉， Un∞=max1≤j≤J-1Unj，

Z0h={U=Uj|U0-i=UJ+i=0， i=0， 1， 2， 3;

j=-3， -2， -1，0，1， …， J+2， J+3}.

對(duì)初邊值問(wèn)題（2）～（4），我們構(gòu)造如下的兩層非線性差分格式

Unjt+α ΨUn+1/2j+ΦUn+1/2j=0，

1≤j≤ J-1，1≤n≤N-1 （7）

U0j=u0xj，j=0， 1， …， J （8）

Un∈Z0h，n=0， 1，? …， N （9）

其中

ΨUn+1/2j=Ψ1Un+1/2j+Ψ2Un+1/2j+

Ψ3Un+1/2j，

Ψ1Un+1/2j=12{Un+1/2jUn+1/2jx^+

Un+1/2j2x^}，

Ψ2Un+1/2j=-15{Un+1/2jUn+1/2jx¨+

Un+1/2j2x¨}，

Ψ3Un+1/2j=130{Un+1/2jUn+1/2jx…+

Un+1/2j2x…}，

ΦUn+1/2j=2615Un+1/2jxx-x^-65Un+1/2jxx-x¨+

715Un+1/2jx^x^x¨.

關(guān)于格式（7）～（9）對(duì)守恒量（5）（6）的模擬效果，我們有如下結(jié)果.

定理2.1 差分格式（7）～（9）關(guān)于以下的離散能量守恒：

Qn=h∑J-1j=1Unj=Qn-1=…=Q0（10）

En=Un2=En-1=…=E0（11）

其中 n=1， 2， …， N.

證明 由條件（9）和分部求和公式［14］，我們有

h∑J-1j=1Un+1/2jUn+1/2jx^=-h∑J-1j=1Un+1/2jx^Un+1/2j.

則h∑J-1j=1Un+1/2jUn+1/2jx^=0. 同理可得

h∑J-1j=1Un+1/2jUn+1/2jx¨=0，

h∑J-1j=1Un+1/2jUn+1/2jx…=0.

又由

h∑J-1j=1Un+1/2j2x^=0，

h∑J-1j=1Un+1/2j2x¨=0，

h∑J-1j=1Un+1/2j2x…=0，

以及

h∑J-1j=1Un+1/2jxx-x^=0，? h∑J-1j=1Un+1/2jxx-x¨=0，

h∑J-1j=1Un+1/2jx^x^x¨=0，

我們有

h∑J-1j=1ΨUn+1/2j=0， h∑J-1j=1ΦUn+1/2j=0（12）

將（7）式兩端乘以h后對(duì)j從1到J-1求和，注意到（12）式，我們有

h∑J-1j=1Unjt=0（13）

由Qn的定義，我們將（13）式對(duì)n遞推可得（10）式.

另一方面， （7）式對(duì)2Un+1/2取內(nèi)積，由（9）式和分部求和公式有

Un2t+2α〈ΨUn+1/2，Un+1/2〉+

2〈ΦUn+1/2j，Un+1/2〉=0（14）

由于

〈Ψ1Un+1/2，Un+1/2〉=

h∑J-1j=1Un+1/2j2Un+1/2jx^+

h∑J-1j=1Un+1/2j2x^Un+1/2j=

h∑J-1j=1Un+1/2j2Un+1/2jx^-

h∑J-1j=1Un+1/2j2Un+1/2jx^=0，

同理可得

〈Ψ2Un+1/2，Un+1/2〉=0，

〈Ψ3Un+1/2，Un+1/2〉=0，

即

〈ΨUn+1/2，Un+1/2〉=0（15）

又

〈Un+1/2xx-x^，Un+1/2〉=0， 〈Un+1/2xx-x¨，Un+1/2〉=0，

〈Un+1/2x^x^x¨，Un+1/2〉=0，

即

〈ΦUn+1/2，Un+1/2〉=0（16）

由En的定義，將（15）（16）式代入（14）式后對(duì)n遞推可得式（11）. 證畢.

3 差分格式的收斂性和穩(wěn)定性

由Taylor公式，如果函數(shù)ux足夠光滑，則當(dāng)h→0時(shí)我們有

ujx^=dudxx=xj+13！h2d3udx3x=xj+

15！h2d5udx5x=xj+Oh6，

ujx¨=dudxx=xj+43！h2d3udx3x=xj+ 165！h2d5udx5x=xj+Oh6，

ujx…=dudxx=xj+93！h2d3udx3x=xj+

815！h2d5udx5x=xj+Oh6，

ujxx-x^=d3udx3x=xj+14h2d5udx5x=xj+

140h4d7udx7x=xj+Oh6，

ujxx-x¨=d3udx3x=xj+34h2d5udx5x=xj+

23120h4d7udx7x=xj+Oh6，

以及

ujx^x^x¨=d3udx3x=xj+h2d5udx5x=xj+

25h4d7udx7x=xj+Oh6.

于是

32ujx^-35ujx¨+110ujx…=

dudxx=xj+Oh6（17）

2615ujxx-x^-65ujxx-x¨+715ujx^x^x¨=

d3udx3x=xj+Oh6（18）

格式（7）～（9）的截?cái)嗾`差定義為

rnj=unjt+Φun+1/2j+αΨun+1/2j（19）

則由Taylor公式以及（17）（18）式可知

rnj =Oτ2+h6（20）

引理3.1 U∈Z0h，恒有 Ux^2≤Ux2， Ux¨2≤Ux^2及Ux…2≤Ux^2.

證明 由于Ux2=Ux-2，由Cauchy-Schwarz不等式有

Ux^2=14h∑J-1j=1Ujx+Ujx-·

Ujx+Ujx-≤Ux2，

Ux¨2=14h∑J-1j=1Uj+1x^+Uj-1x^·

Uj+1x^+Uj-1x^≤Ux^2，

Ux…2=19h∑J-1j=1Uj+2x^+Ujx^+Uj-2x^·

Uj+2x^+Ujx^+Uj-2x^≤

h9·3·∑J-1j=1Uj+2x^2+Ujx^2+

Uj-2x^2≤Ux^2.

證畢.

定理3.2 設(shè)u（x， t）足夠光滑且u0∈H2.則格式（7）～（9）的解Un以范數(shù) · 收斂到原問(wèn)題（2）～（4）的解，收斂階為Oτ2+h6.

證明 記enj=unj-Unj. 由式（19）減去式（7）得

rnj=enjt+Φen+1/2j+αΨun+1/2j-

αΨUn+1/2j（21）

類似（16）式，有〈Φen+1/2， en+1/2〉=0. 則（21）式對(duì)2en+1/2取內(nèi)積可得

〈rn， 2en+1/2〉=en2t+2α〈Ψun+1/2- ΨUn+1/2， en+1/2〉（22）

因?yàn)閡x， t在有界閉區(qū)域xL， xR×0， T內(nèi)具有連續(xù)三階偏導(dǎo)數(shù)uxxx，所以

uL∞≤C， uxL∞≤C（23）

再由微分中值定理有

un+1/2jx=uxj+1， tn+tn+12-uxj， tn+tn+12h=

xuξj，? tn+tn+12，

其中 xj≤ξj≤xj+1. 故

un+1/2x∞≤C（24）

又由

Ψ1un+1/2j-Ψ1Un+1/2j=12un+1/2j（en+1/2j）x^+

12en+1/2j（un+1/2j）x^+（en+1/2jun+1/2j）x^-

Ψ1en+1/2j

及〈Ψ1en+1/2， en+1/2〉=0，結(jié)合（23）（24）式以及引理3.1有

〈Ψ1un+1/2-Ψ1Un+1/2， en+1/2〉=

12h∑J-1j=1［un+1/2j（en+1/2j）x^+en+1/2j（un+1/2j）x^］en+1/2j-

2h∑J-1j=1en+1/2jun+1/2j（en+1/2j）x^=

12h∑J-1j=1en+1/2j（un+1/2j）x^en+1/2j-

h∑J-1j=112h（un+1/2jen+1/2j+1en+1/2j-un+1/2jen+1/2jen+1/2j-1）=

12h∑J-1j=1en+1/2j（un+1/2j）x^en+1/2j+h2∑J-1j=1en+1/2j+1en+1/2j

（un+1/2j）x≤C（en+12+en2）（25）

同理可得

〈Ψ2un+1/2-Ψ2Un+1/2， en+1/2〉≤

C（en+12+en2）（26）

〈Ψ3un+1/2-Ψ3Un+1/2， en+1/2〉≤

C（en+12+en2）（27）

又

〈rn， 2en+1/2〉≤rn2+12（en+12+en2）（28）

將（25）～（28）式代入（22）式，整理得

en+12-en2≤τrn2+

Cτen+12+en2（29）

將（29）式由0到n-1遞推求和得

en2≤e02+τ∑n-1l=0rl2+Cτ∑nl=0el2（30）

由（20）式有

τ∑n-1l=0rl2≤nτmax0≤l≤n-1rl2≤ T·O（τ2+h6）2.

再由e02=0知式等價(jià)于

en2≤O（τ2+h6）2+Cτ∑nl=0el2.

由離散的Gronwall不等式［14］，有en≤O（τ2+h6）. 證畢.

類似定理3.2的證明，我們有

定理3.3 在定理3.2的條件下，差分格式（7）～（9）的解Un以范數(shù) · 關(guān)于初值無(wú)條件穩(wěn)定.

4 數(shù)值算例

當(dāng)參數(shù)α=6時(shí)，KdV方程的孤波解［15］為u（x， t）=Asech2kx-ω? t，其中A=r2，k=r2，ω=r3/22（r是常數(shù)）. 由此可知KdV方程的物理邊界（漸近邊界）條件滿足當(dāng) x →∞時(shí)，有u（x，t）→0，t>0. 因此，我們只需取-xL0，xR0，問(wèn)題（2）～（4）與KdV方程的Cauchy問(wèn)題就是一致的. 在數(shù)值計(jì)算時(shí)我們?nèi)ˇ?6，r=0.5. 問(wèn)題（2）～（4）的初值函數(shù)可取為 u0x=14sech224x. 固定xL=-20，xR=40，T=32. 為驗(yàn)證差分格式（7）～（9）式的數(shù)值解在不同范數(shù)下的理論精度均為O（τ2+h6），定義

orderl2=log2enh，τ/e8nh2，τ8，

orderl∞=log2enh，τ∞/e8nh2，τ8∞.

對(duì)于τ和h的不同取值，數(shù)值解在不同時(shí)刻的誤差及差分格式（7）～（9）精度的檢驗(yàn)參見表1，2，差分格式（7）～（9）對(duì)不變量（5）（6）的數(shù)值模擬見表3.

從數(shù)值結(jié)果可以看出，差分格式（7）～（9）明顯具有6階空間精度，并且能夠合理地模擬守恒性質(zhì)（5）（6），因而是有效的.

參考文獻(xiàn)：

［1］ Korteveg D J， de Vries G. On the change of long waves advancing in a rectangular canal and on a new type of long stationary waves ［J］. Phil Mag， 1895， 39： 422.

［2］ Crighton D G. Applications of KdV ［J］. Acta Appl Math， 1995， 39： 39.

［3］ Khattak A J. A comparative study of numerical solutions of a class of KdV equation ［J］. Appl Math Comput， 2008， 199： 425.

［4］ Benkhaldoun F， Seaid M. New finite-volume relaxation methods for the third-order differential equations ［J］. Commun Comput Phys， 2008， 4： 820.

［5］ Dutykh D， Marx C， Francesco F. Geometric numerical schemes for the KdV equation ［J］. Comp Math Math Phys+， 2013， 53： 221.

［6］ Zabusky N J， Martin D K. Interaction of "solitons" in a collisionless plasma and the recurrence of initial states ［J］. Phys Rev Lett， 1965， 15： 240.

［7］ CourtèsC， Frédéric L， Frédéric R. Error estimates of finite difference schemes for the Korteweg–de Vries equation ［J］. IMA J Numer Anal， 2020， 40： 628.

［8］ Wang H P， Wang Y S， Hu Y Y. An explicit scheme for the KdV equation ［J］. Chin Phys Lett， 2008， 25： 2335.

［9］ Shen J， Wang X， Sun Z.The conservation and convergence of two finite difference schemes for KdV equations with initial and boundary value conditions ［J］. Numer Math： Theory Me， 2020， 13： 253.

［10］ Wang X， Sun Z. A second order convergent difference scheme for the initial-boundary value problem of Korteweg-de Vires equation ［J］. Numer Meth Part D E， 2021， 37： 2873.

［11］ 趙修成， 黃浪揚(yáng). KdV方程的一個(gè)緊致差分格式 ［J］. 工程數(shù)學(xué)學(xué)報(bào)， 2015， 32： 876.

［12］ Zhang X， Zhang P. A reduced high-order compact finite difference scheme based on proper orthogonal decomposition technique for KdV equation ［J］. Appl Math Comput， 2018， 339： 535.

［13］ Miura R M， Gardner C S， Kruskal M D. Korteweg-de Vries equation and generalizations （II） ：? existence of conservation laws and constants of motion ［J］. J Math Phys， 1968， 9： 1204.

［14］ Zhou Y L. Applications of discrete functional analysis to the finite difference method ［M］. Beijing： International Academic Publishers， 1990.

［15］ Dehghan M， Shokri A. A numerical method for KdV equation using collocation and radial basis functions ［J］. Nonlinear Dynam， 2007， 50： 111.

引用本文格式：

中 文：? 郭禎， 吳明明， 胡勁松. KdV方程的一個(gè)六階空間精度守恒差分格式［J］. 四川大學(xué)學(xué)報(bào)：? 自然科學(xué)版， 2023， 60：? 031006.

英 文：? Guo Z， Wu M M， Hu J S. A conservative difference scheme with 6-order spatial accuracy for the KdV equation ［J］. J Sichuan Univ：? Nat Sci Ed， 2023， 60：? 031006.