一類偏微分方程的差分解法

袁毅楓

【摘要】本文提出把一類橢圓區(qū)域上的偏微分方程變換成圓區(qū)域上的偏微分方程，然后在極坐標(biāo)下作網(wǎng)格劃分，使所有網(wǎng)格邊界點與邊界重合，并給出相應(yīng)的差分方程和解法，從而避免了一般差分解法由邊界條件轉(zhuǎn)移問題所產(chǎn)生的誤差.

【關(guān)鍵詞】差分；方程；解

一、問題的提出

考慮方程：

λ22ux2+2uy2=f（x，y），（x，y）∈D，

u|D=α（x，y），

其中λ≠0，D=（x，y）x2y2+y2<1.（1）

對此類方程的一般差分解法是首先分別沿x，y軸方向取步長為h和k，作兩族與坐標(biāo)平行的直線：

xi=ih，i=0，±1，±2，….

yj=jk，j=0，±1，±2，….

兩族直線的交點為網(wǎng)格點.

Dk={（xi，yj）|（xi，yj）∈D，并且（xi+1，yj），（xi-1，yj），（xi，yj+1），（xi，yj-1）∈D}為內(nèi)點.

Dh={（xi，yj）|（xi，yj）∈D-Dh}為網(wǎng)格邊界點.

顯然Dh與D不重合.

在Dh上建立差分格式：

Lhuij=f（xi，yj），（xi，yj）∈Dh，

uij=g（xi，yj），（xi，yj）∈Dh.

（Lh為差分算子）于是出現(xiàn)了邊界條件從D到Dh轉(zhuǎn)移的問題，帶來了誤差并逐層傳播，影響所得結(jié)果的準(zhǔn)確性，所以考慮對這類問題做某種變換，同時采取某種特別的解域劃分，使網(wǎng)格邊界點與邊界重合，并給出這種劃分下相應(yīng)的差分方程和解法.

二、建立差分方程

令x=λt，（1）成為

2ut2+2uy2=f1（t，y），（t，y）∈D1，

u|D1=α1（t，y），

其中D1={（t，y）|t2+y2<1}.（2）

令t=rcosθ，y=rsinθ，

即r=t2+y2，θ=arctanyt，（2）成為

2ur2cos2θ+ur·sin2θr-2uθr·cosθ·sinθr+2ur2·sin2θ+2uθr·sinθ·cosθr+ur·cos2θr-2urθ·sinθ·cosθr+uθ·sinθ·cosθr2+2uθ2·sin2θr2+2urθ·sinθ·cosθr+2urθ·sinθ·cosθr-uθ·sinθ·cosθr2+2uθ2·cos2θr2-uθ·sinθ·cosθr2=1r·rrur+1r2·2uθ2=f（r，θ），（r，θ）∈D2，

u|D2=α2（1，θ）=α2（θ）.

其中D2={（r，θ）|0≤r<1，0≤θ≤2π}.（3）

當(dāng)r=0時，方程（3）的系數(shù)無定義，因此，只有在r>0的情形下其解才有意義，為了定出有定義的解，需要補充在r=0處u有界.

沿徑向r和角度θ方向取等步長Δr和Δθ作網(wǎng)格劃分：

ri=i+12Δr，i=0，1，2，…，I.

θj=jΔθ，j=0，1，2，…，J-1.

Δr=1I+12，Δθ=2πJ.

于是方程（3）的差分方程為

1ri（Δr）2δr（rδruij）+1r2i（Δθ）2δ2θuij=f2（ri，θj），（4）

其中，δ和δθ分別表示r方向和θ方向的中心差分算子.

即，δruij=ui+12，j-ui-12，j；

δθuij=ui，j+12-ui，j-12.

利用循環(huán)邊界條件u（r，θ）=u（r，θ+2π），具體寫出方程（4）在j=0和j=J-1上的差分方程：

j=0，1ri（Δr2）[ri+12ui+1，0-（ri+12+ri-12）ui，0+ri-12ui-1，0]+1r2i（Δθ）2（ui，1-2ui，0+ui，j-1）=f2（ri，0），

i=1，2，…，I-1.（5）

j=J-1，1ri（Δr）2[ri+12ui+1，j-1-（ri+12+ri-12）ui，j-1+ri-12ui-1，j-1]+1r2i（Δθ）2（ui，0-2ui，j-1+ui，j-2）=f2（ri，θj-1），

i=1，2，…，I-1.（6）

當(dāng)j=1，2，…，J-2時，差分方程為

1ri（Δr）2[ri+12ui+1，j-（ri+12+ri-12）ui，j+ri-12ui-1，j]+1r2i（Δθ）2（ui，j+1-2ui，j+ui，j-1）=f2（ri，θj），

i=1，2，…，I-1.（7）

令（7）中ri=Δr2，得（r0，θj）處的差分方程：

2（Δr）2（u1，j-u0，j）+4（Δr）2（Δθ）2（u0，j+1-2u0，j-u0，j-1）=f2（r0，θj），

j=0，1，2，…，J-1.（8）

至此建立了I×J個網(wǎng)格點上的差分方程.

在J個網(wǎng)絡(luò)邊界點上u（1，θj）=α2（θj），

j=0，1，2，…，J-1.（9）

至此得到了方程（4）的一個差分格式，并且使得網(wǎng)格邊界點與邊界點完全重合.

三、求解差分方程

把（5）（6）（7）（8）分別整理，并對每處（r0，θj）記為

αi，jui，j-αi+1，jui+1，j-αi-1，jui-1，j-αi，j+1ui，j+1-αi，j-1ui，j-1=-f2（ri，θj），

i=0，1，2，…，I-1；j=0，1，2，…，J-1.（10）

網(wǎng)格邊界點上ui，j由（9）代入.

把D2的內(nèi)點按如圖所示的排列給予序號.

于是可得方程組（10）的矩陣形式

AU=F，（11）

其中A為方程組的系數(shù)矩陣，U為uij依序排列的矩陣，F(xiàn)為f2（ri，θj）依序排列的矩陣.

為了證明所得差分方程是可解的，需要用到下述定理.

極值定理：設(shè)差分算子

Lhuij=αijuij-∑blkulk，（xl，yk）∈G（xi，yj），

其中G（xi，yj）表示格式中除（xi，yi）之外相鄰網(wǎng)格點集合，在Dh∪Dh上，算子的系數(shù)滿足

aij>0，（xi，yj）∈Dh，

blk>0，（xi，yj）∈G（xi，yj），

dij=aij-∑blk≥0.

uij是定義在Dk上的函數(shù)，uij≠常數(shù)，如果Lhuij≤0，（xi，yj）∈Dh，則uij不可能在內(nèi)點取正的最大值：

maxuij≤maxuij.

Dh Dh

如果Lhuij≥0，（xi，yj）∈Dh，則uij不可能在內(nèi)點取得負(fù)的最小值：

minuij≤minuij.

Dh Dh

證明略.

可以驗證差分方程（10）滿足極值定理，由此可得差分格式的穩(wěn)定性和收斂性.另外，由極值定理可知，滿足差分方程（10）的uij除非在所有網(wǎng)格點上均取同一數(shù)值，不會在網(wǎng)格區(qū)域的內(nèi)點取到最大值和最小值.由此立刻可以證明，相應(yīng)于α2=0情形的差分方程，即在邊界點取零值的差分方程，只有全為零的解，這說明相應(yīng)于線性方程組（11）的齊次方程組只有零解，因此，線性代數(shù)方程組（11）必有唯一解，這說明了所得的差分方程恒有唯一解.

由此可以采用直接的方法來求解，但由于在實際計算中為了保證精確度，步長往往取得很小，網(wǎng)格點很多，所得的差分方程是一個高階線性代數(shù)方程，而用直接求解法就會因為受計算機存儲量的限制而產(chǎn)生困難.

由于每個內(nèi)點上的值只依賴周圍相鄰四點上的值，相應(yīng)地，線性方程組的系數(shù)矩陣有大量零元素，即是一個稀疏陣，采用迭代法求解就可以充分利用這個優(yōu)點.為此這里給出如下辦法求解這個差分方程.

首先任意給定在網(wǎng)格區(qū)域內(nèi)點（xi，yi）上的數(shù)量{u（0）ij}及{u（1）i-1}（i=0，1，2，…，I-1）作為解的近似，把這組數(shù)值代入迭代式

u（k+1）ij=ωαi+1αiju（k）i+1，j+αi-1αiju（k+1）i-1，j+αi，j+1αiju（k）i，j+1

+αi，j-1αiju（k+1）i，j-1+（1-ω）u（k）ij，（xi，yi）∈Dk，

u（k+1）ij=α2（θ1），（xi，yi）∈Dk

其中ω是一個常數(shù)，可在1<ω<2中選定.當(dāng)k相當(dāng)大時，{u（k）ij}就給出所要求的近似解.通常對充分大的k，當(dāng)相鄰兩次迭代解{u（k-1）ij}，{u（k）ij}間的誤差maxi，j|u（k）ij-u（k-1）ij|或1N∑i，j|u（k）ij-u（k-1）ij|，其中N為網(wǎng)格點的總數(shù)小于某個預(yù)先給定的適當(dāng)小的控制數(shù)ε>0時，就可以結(jié)束迭代過程，完成差分方程的求解.