王金平, 莊清渠

(華僑大學(xué) 數(shù)學(xué)科學(xué)學(xué)院, 福建 泉州 362021)

五階常微分方程的Petrov-Galerkin譜元法

王金平, 莊清渠

(華僑大學(xué) 數(shù)學(xué)科學(xué)學(xué)院, 福建 泉州 362021)

通過區(qū)間剖分，降低數(shù)值逼近多項式的階數(shù)，構(gòu)造滿足試探函數(shù)空間和檢驗函數(shù)空間的基函數(shù)，使得離散問題所對應(yīng)的線性系統(tǒng)的系數(shù)矩陣是稀疏的，并可以進(jìn)行有效地求解.數(shù)值算例驗證了五階常微分方程的Petrov-Galerkin譜元法的有效性和高精度.

五階常微分方程； Petrov-Galerkin譜元法； 基函數(shù)； 數(shù)值實驗

譜方法作為數(shù)值求解微分方程的方法之一，已被廣泛應(yīng)用于自然科學(xué)和工程技術(shù)問題的數(shù)值計算[1-7].近二十年來，奇數(shù)次高階微分方程基于譜方法的數(shù)值模擬也受到不少科研工作者的關(guān)注[8-13].由于奇次高階微分方程的高階項缺少對稱性，若使用配置點(diǎn)法求解奇次高階微分方程邊值問題，如果不能選取恰當(dāng)?shù)呐渲命c(diǎn)或采用高階多項式進(jìn)行逼近，容易導(dǎo)致條件數(shù)比較高，造成計算不穩(wěn)定[10].吳勝等[13]提出Legendre-Petrov-Galerkin譜元法對三階微分方程進(jìn)行求解.由于將計算區(qū)間進(jìn)行剖分，在每個區(qū)間上不需要使用太高的多項式階數(shù)，將問題分解成一系列子問題后再用Schur補(bǔ)過程進(jìn)行求解，計算穩(wěn)定，且具有高精度.本文研究五階常微分方程的Petrov-Galerkin譜元逼近，主要考慮方程的數(shù)值計算.

1 Petrov-Galerkin譜元法

記Λ=(-1,1)，設(shè)α和β為正常數(shù).考慮五階常微分齊次邊值問題，即

首先，將區(qū)間Λ剖分成K(K≥2)個子區(qū)間

Λk=(ak-1,ak),k=1,2,…,K.

上式中:-1=a0

定義區(qū)間Λk到參考區(qū)間Λ的坐標(biāo)變換為

?x∈Λk.

定義分片多項式空間為

PN;K(Λ)∶={u;u|Λk∈PN(Λk),k=1,2,…,K}，

上式中:PN為Λk上次數(shù)不超過N的全體多項式組成的空間.

用N表示離散參數(shù)(N,K).定義試探函數(shù)空間和檢驗函數(shù)空間為

VN={u;u|Λ∈PN,K(Λ);u∈C1(Λ),u(±1)=ux(±1)=ux，x(1)=0},

WN={u;u|Λ∈PN+1,K(Λ);u∈C2(Λ),u(±1)=ux(±1)=ux，x(±1)=0}.

則問題(1)的Petrov-Galerkin譜元逼近形式為找uN∈VN, 使得

?vN∈WN.

(2)

為了方便表達(dá)，將式(2)記為α(uN,vN)=(f,vN)，記

定義

則可以驗證

當(dāng)k=1,2,…,K-1時，令

而且有

}.

將文獻(xiàn)[13]求解三階微分方程的Petrov-Galerkin譜元方法的逼近過程推廣到式(2)的計算中，詳細(xì)計算過程包括如下4個步驟.

?

?

步驟4 由式(7)，(8)可得

?vN∈WN.

(9)

問題(2)被分解成兩套子區(qū)間問題(6),(7)及Schur補(bǔ)問題(8).因此，可以大大提高計算效率.

(αAk-βBk+Ck)Uk=Fk,k=1,2,…,K.

上式中：

由Legendre多項式的正交性可知

j=i+5,

j=i+4,

j=i+3,

j=i+2,

j=i+1,

j=i,

j=i-1,

j=i-2,

j=i-3,

j=i-4,

j=i-5,

j=i-6,

其他.

j=i+1,

j=i,

j=i-1,

j=i-2,

j=i-3,

其他.

j=i-1,

其他.

2 數(shù)值實驗

例1 問題(1)在x∈(-1,1)上的解析解為

u(x)=(x-1)2sin2mπx.

取計算參數(shù)α=β=m=1.在半log尺度下，h=1/2時，L2-誤差ε及H1-誤差隨多項式階數(shù)N的變化情況，如圖1(a)所示.由圖1(a)可知：誤差隨N呈指數(shù)收斂，這說明對于解析解，數(shù)值解具有譜收斂性質(zhì).在log-log尺度下，N=8時，L2-誤差及H1-誤差隨h的變化情況，如圖1(b)所示.由圖1(b)可知:誤差關(guān)于h呈代數(shù)衰減.

(a) 誤差隨N的變化情況 (b) 誤差隨h的變化情況圖1 誤差分別隨N和h的變化圖Fig.1 Errors as a function of N and h respectively

例2 考慮變系數(shù)五階齊次邊值問題,即

此時，問題的Petrov-Galerkin譜元逼近形式為找uN∈VN, 使得

?vN∈WN.

(10)

該問題同樣可采用第1節(jié)提出的計算過程進(jìn)行快速計算.但由于α(x)，β(x)的存在，弱形式中前兩個積分項一般不可精確計算，因此，采用Gauss-Lobatto數(shù)值積分公式替代.

取α(x)=x，β(x)=sin 10x，γ=1，u(x)=sin3πx進(jìn)行計算.在半log尺度下，h=1/2時，誤差隨多項式階數(shù)N的變化情況，如圖2(a)所示.在log-log尺度下，N=8時，誤差隨h的變化情況，如圖2(b)所示.由圖2可知:誤差隨N呈指數(shù)收斂，誤差隨h呈代數(shù)收斂，收斂性質(zhì)跟常系數(shù)方程的計算結(jié)果類似.從而說明Petrov-Galerkin譜元逼近法對求解一些變系數(shù)五階常微分方程也是有效的.

(a) 誤差隨N的變化情況 (b) 誤差隨h的變化情況圖2 依賴于N和h的誤差變化Fig.2 Errors as a function corresponding to N and h respectively

3 結(jié)束語

利用Petrov-Galerkin譜元法對五階常微分方程進(jìn)行求解.首先，通過區(qū)間剖分，將變分問題轉(zhuǎn)化為一系列子問題；然后，構(gòu)造恰當(dāng)?shù)脑囂胶瘮?shù)和檢驗函數(shù)，由此得到稀疏的線性系統(tǒng)再進(jìn)行求解.數(shù)值算例說明Petrov-Galerkin譜元逼近法對一些常系數(shù)、變系數(shù)的五階常微分方程都是有效的，且具有高精度.

[1] GOTTLIEB D，ORSZAG S A.Numerical analysis of spectral methods: Theory and applications[M].Philadelphia：SIAM,1977：1-179.

[2] BERNARDI C，MADAY Y.Approximations spectrales de probl`emes aux limites elliptiques[M].Berlin：Springer-Verlag,1992：1-248.

[3] GUO Benyu.Spectral methods and their applications[M].Singapore：World Scientific,1998:100-358.

[4] KARNIADAKIS G，SHERWIN S J.Spectral/hp element methods for computational fluid dynamics[M].London：Oxford University Press,2005：187-348.

[5] SHEN Jie，TANG Tao.Spectral and high-order methods with applications[M].Beijing:Science Press of China,2006:183-298.

[6] CANUTO C,HUSSAINI M Y,QUARTERONI A,etal.Spectral methods: Fundamentals in single domains[M].Berlin：Springer-Verlag,2006：401-470.

[7] CANUTO C,HUSSAINI M Y,QUARTERONI A,etal.Spectral methods: Evolution to complex geometries and applications to fluid dynamics[M].Berlin：Springer Verlag,2007：237-357.

[8] MA Heping，SUN Weiwei.A Legendre-Petrov-Galerkin and chebyshev collocation method for third-order differential equations[J].SIAM J NUMER ANAL,2000,38(5):1425-1438.

[9] MA Heping，SUN Weiwei.Optimal error estiamtes of the Legendre-Petrov-Galerkin method for the Korteweg-De Vries equation[J].SIAM J NUMER ANAL,2001,39(4):1380-1394.

[10] SHEN Jie.A new dual-Petrov-Galerkin method for third and higher odd-order differential equations: Application to the KDV equation[J].SIAM J NUMER ANAL,2003,41(5):1595-1619.

[11] SHEN Jie，WANG Lilian.Laguerre and composite Legendre-Laguerre dual-Petrov-Galerkin methods for third-order equations[J].Discrete and Continuous Dynamical Systems Series B,2006,6(6):1381-1402.

[12] ZHENG Chunxiong,WEN Xin,HAN Houde.Numerical solution to a linearized KDV equation on unbounded domain[J].Numerical Methods for Partial Differential Equations,2008,24(2):383-399.

[13] 吳勝,莊清渠.三階微分方程的Legendre-Petrov-Galerkin譜元方法[J].華僑大學(xué)學(xué)報(自然科學(xué)版),2013，34(3):344-348.

(責(zé)任編輯： 陳志賢 英文審校： 黃心中)

Petrov-Galerkin Spectral-Element Method for Solving Fifth-Order Ordinary Differential Equations

WANG Jinping， ZHUANG Qingqu

(School of Mathematical Sciences, Huaqiao University, Quanzhou 362021, China)

The polynomial order in the numerical approximation is reduced by partitioning the interval into several subintervals, and appropriate basis functions of the trial and test spaces are constructed. Which leads to a linear system with sparse coefficient matrix. Then, an efficient computational process is introduced to solve the linear system.Numerical experiment results demonstrate the high accuracy and effectiveness to the Petrov-Galerkin spectral-element method. Keywords：fifth-order ordinary differential equation; Petrov-Galerkin spectral-element method; basis functions; numerical experiments

10.11830/ISSN.1000-5013.201703027

2016-01-17

莊清渠(1980-)，男，副教授，博士，主要從事微分方程數(shù)值解法的研究.E-mail:qqzhuang@hqu.edu.cn.

國家自然科學(xué)基金資助項目(11501224); 華僑大學(xué)中青年教師科研提升資助計劃(ZQN-PY201); 華僑大學(xué)研究生科研創(chuàng)新能力培育計劃項目(1400213008)

O 241.8

A

1000-5013(2017)03-0435-06