牟行洋

(江漢藝術(shù)職業(yè)學(xué)院 管理學(xué)院, 湖北 潛江 433100)

0 引言

長(zhǎng)期以來,半線性橢圓型方程一直受到人們的廣泛關(guān)注.其原因是因?yàn)樵S多數(shù)學(xué)物理問題,如源于非線性源的非線性擴(kuò)散理論[1],熱力學(xué)中的氣體燃燒理論[2-3],量子場(chǎng)論和統(tǒng)計(jì)力學(xué)[4-6]以及星球系的重力平衡理論[2,7]都與它有著極大的淵源.而且,數(shù)學(xué)內(nèi)部的許多分支,如幾何中的Yamabe問題[8]和等周不等式[9],調(diào)和分析中的Hardy-Littlewood-Sobolev不等式[10],Yang-Mills泛函的非極小解存在性[11]以及人口動(dòng)力系統(tǒng)[12]等問題都與之有著深刻的聯(lián)系.

許多學(xué)者研究了二維半線性橢圓方程正解的存在性[13-18],其中劉偉和芮洪興利用二重網(wǎng)格混合元算法求解了二維半線性橢圓問題的混合元解[19],大大降低計(jì)算量,提高計(jì)算效率.麥雄發(fā),胡淑娟和丁方允針對(duì)一類半線性橢圓邊值問題,利用有限元給出了一種一般迭代數(shù)值的方法,收斂定理,誤差估計(jì)[20].因此,對(duì)于這一類方程的研究無論是在理論上,還是在實(shí)際應(yīng)用上都有著非常重要的意義.本文利用有限元方法結(jié)合牛頓法研究了二維半線性橢圓方程的邊值問題,推導(dǎo)出了它的牛頓迭代式,該算法采用Matlab編制程序,具有易編程實(shí)現(xiàn),求解速度快,計(jì)算精度高等優(yōu)點(diǎn),證明該方法的可行性和有效性.

1 數(shù)學(xué)模型

考慮如下的二維半線性橢圓方程數(shù)學(xué)問題:

(1)

其中Ω是R2上的有界閉區(qū)域,邊界?Ω=Ω1∪Ω2分片光滑且Ω1與Ω2互不相交,a(x,y)≥a0>0,f(x,y),φ(x,y),β(x,y)都是關(guān)于(x,y)的二維光滑函數(shù),σ為大于等于零的常數(shù),n是邊界上的外法線方向,g(u)是關(guān)于u的非線性函數(shù).當(dāng)這些函數(shù)都是已知條件求解u時(shí),是二維半線性橢圓方程的數(shù)學(xué)模型.對(duì)于問題(1),在Ω1上給出Dirichlet邊界條件,在Ω2上給出Robin邊界條件,特別地,當(dāng)σ=0時(shí)是Neumann邊界條件.

2 牛頓迭代法

引入Sobolev空間H1(Ω)的子空間V={v∈H1(Ω):v|Ω1=0},在式(1)中第一個(gè)式子的兩端用v(x,y)作內(nèi)積,有

(2)

由泰勒公式有

g(uk)=g(uk-1)+g′(uk-1)(uk-uk-1)+o(uk-uk-1),k=1,2,3,…

(3)

由式(2)和式(3)得迭代形式如下

(g(uk-1)+g′(uk-1)(uk-uk-1),v)=(f(x,y),v),k=1,2,3,…

(4)

對(duì)上式應(yīng)用格林公式,得

(5)

將式(1)中第三個(gè)式子變形為

(6)

將式(6)代入式(5),得

(7)

式(7)稱為式(1)的牛頓迭代格式.

3 數(shù)值算例

在單位矩形Ω={(x,y)|0≤x≤1m,0≤y≤1m}上求解問題

其中

a(x,y)=2x+3y,

真解

u(x,y)=x2y2(1-x)(1-y),

將有界區(qū)域Ω沿x,y方向作均勻三角剖分,剖分的步長(zhǎng)h=0.05,采用Lagrange插值基函數(shù)分4種類型討論,分別是一次、二次和三次,此時(shí)它們的自由度分別為441個(gè)、1681個(gè)、3721個(gè).

情形1 四次偶非線性g(u)=u4.

f(x,y)=-4xy2(1-x)(1-y)+2x2y2(1-y)-2y2(2x+3y)(1-x)(1-y)-4xy2(2x+3y)(1-y)-

6yx2(1-x)(1-y)+3x2y2(1-x)-2x2(2x+3y)(1-x)(1-y)-4yx2(2x+3y)(1-x)+u4.

將g(u)代入式(7),得牛頓迭代式(8).取初值U0=0,根據(jù)已知的邊界條件進(jìn)行迭代計(jì)算,此處僅討論Dirichlet邊界條件,Neumann和 Robin邊界條件討論方法類似進(jìn)行.在不同的基函數(shù)下,數(shù)值解剖分圖見圖1,真解與數(shù)值解的誤差剖分圖見圖2-4,對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間計(jì)算結(jié)果見表1.

表1 在不同基函數(shù)下對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間的計(jì)算結(jié)果

(8)

情形2 五次奇非線性g(u)=u5.

f(x,y)=-4xy2(1-x)(1-y)+2x2y2(1-y)-2y2(2x+3y)(1-x)(1-y)-4xy2(2x+3y)(1-y)-

6yx2(1-x)(1-y)+3x2y2(1-x)-2x2(2x+3y)(1-x)(1-y)-4yx2(2x+3y)(1-x)+u5.

將g(u)代入式(7),得牛頓迭代式(9).取初值U0=0,根據(jù)已知的邊界條件進(jìn)行迭代計(jì)算,在不同的基函數(shù)下, 數(shù)值解剖分圖見圖5,真解與數(shù)值解的誤差剖分圖見圖6-8,對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間計(jì)算結(jié)果見表2.

表2 在不同基函數(shù)下對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間的計(jì)算結(jié)果

(9)

情形3 奇有界g(u)=sinu.

f(x,y)=-4xy2(1-x)(1-y)+2x2y2(1-y)-2y2(2x+3y)(1-x)(1-y)-

4xy2(2x+3y)(1-y)-6yx2(1-x)(1-y)+3x2y2(1-x)-

2x2(2x+3y)(1-x)(1-y)-4yx2(2x+3y)(1-x)+sinu.

將g(u)代入式(7),得牛頓迭代式(10).取初值U0=0,根據(jù)已知的邊界條件進(jìn)行迭代計(jì)算,在不同的基函數(shù)下, 數(shù)值解剖分圖見圖9,真解與數(shù)值解的誤差剖分圖見圖10-12,對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間計(jì)算結(jié)果見表3.

(10)

表3 在不同基函數(shù)下對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間的計(jì)算結(jié)果

情形4 偶有界g(u)=cosu.

f(x,y)=-4xy2(1-x)(1-y)+2x2y2(1-y)-2y2(2x+3y)(1-x)(1-y)-4xy2(2x+3y)(1-y)-

6yx2(1-x)(1-y)+3x2y2(1-x)-2x2(2x+3y)(1-x)(1-y)-4yx2(2x+3y)(1-x)+cosu.

將g(u)代入式(7),得牛頓迭代式(11).取初值U0=0,根據(jù)已知的邊界條件進(jìn)行迭代計(jì)算,在不同的基函數(shù)下,數(shù)值解剖分圖見圖13,真解與數(shù)值解的誤差剖分圖見圖14-16,對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間計(jì)算結(jié)果見表4.

(11)

表4 在不同基函數(shù)下對(duì)應(yīng)的迭代次數(shù)、L2誤差、最大誤差和CPU運(yùn)行時(shí)間的計(jì)算結(jié)果

4 結(jié)論

利用有限元結(jié)合牛頓法求解了一類二維半線性橢圓方程邊值問題,體現(xiàn)了此方法的可行性和有效性.無論是四次偶非線性情形還是五次奇非線性情形,在相同的基函數(shù)下對(duì)應(yīng)誤差的大小變化不大,但是,隨著插值基函數(shù)次數(shù)的增加,兩類誤差越來越小,且代數(shù)精度較高.奇有界情形兩類誤差的大小變化不大,代數(shù)精度較高且相對(duì)較穩(wěn)定;偶有界情形兩類誤差的大小雖然變化不大,但代數(shù)精度相對(duì)較低.