唐永超，王同科

（天津師范大學(xué)數(shù)學(xué)科學(xué)學(xué)院，天津 300387）

奇異兩點(diǎn)邊值問題改進(jìn)的梯形公式外推方法

唐永超，王同科

（天津師范大學(xué)數(shù)學(xué)科學(xué)學(xué)院，天津 300387）

研究奇異兩點(diǎn)邊值問題的高精度數(shù)值方法.首先，將奇異兩點(diǎn)邊值問題轉(zhuǎn)化為奇異積分的計(jì)算問題.其次，利用改進(jìn)的復(fù)合梯形公式離散奇異積分，針對(duì)幾種不同情形給出了誤差漸近展開式.再次，由誤差估計(jì)式設(shè)計(jì)了一種改進(jìn)的龍貝格算法，利用該算法可以得到問題的高精度數(shù)值解.最后，通過數(shù)值算例說明了算法的有效性.

奇異兩點(diǎn)邊值問題；復(fù)合梯形積分公式；分?jǐn)?shù)階泰勒展開式；誤差漸近展開；龍貝格算法

1 引言與預(yù)備知識(shí)

奇異兩點(diǎn)邊值問題在現(xiàn)代化學(xué)、物理及相關(guān)學(xué)科中具有廣泛的應(yīng)用，其數(shù)值算法的研究一直受到重視，該問題的數(shù)值算法主要有有限差分方法[1]及其樣條求解[2]或外推[3]、Adomian分解方法及改進(jìn)的ADM方法[4-5]，還有研究利用格林函數(shù)將其轉(zhuǎn)化為積分方程，再利用插值離散積分方程進(jìn)行求解[6-7].本研究考慮奇異兩點(diǎn)邊值問題，設(shè)其源項(xiàng)中不含有未知函數(shù).基本思想是將其轉(zhuǎn)化為奇異積分的計(jì)算問題，再利用改進(jìn)的梯形公式及其外推，得到問題的高精度數(shù)值解.

下面簡要給出奇異積分梯形公式的誤差漸近展開式，詳細(xì)推導(dǎo)見文獻(xiàn)[8].設(shè)函數(shù)f（x）在除端點(diǎn)x= a、b外充分光滑，考慮如下積分：

設(shè)f（x）在x=a和x=b處分別有分?jǐn)?shù)階泰勒展開式：

其中：αi、βi為整數(shù)或分?jǐn)?shù)，且滿足α0＜α1＜…＜αm1和β0＜β1＜…＜βm2；r1（x）和r2（x）充分光滑.這里假定α0＞-1，β0＞-1，以保證積分（1）存在.將區(qū)間[a，b]剖分為n等份，步長h=（b-a）/n，節(jié)點(diǎn)xi=a+ih，i=0，1，…，n.在區(qū)間[a，b]上定義如下改進(jìn)的梯形公式：

當(dāng)α0≥0時(shí)，

當(dāng)β0≥0時(shí)，

記

引理1[8]若f（x）在［a，b］上充分光滑，則有

引理2[8]若f（x）在（a，b］上充分光滑，且在x=a處有分?jǐn)?shù)階泰勒展開式（2），則有

其中：ζ（s）是黎曼zeta函數(shù)，當(dāng)s＞1時(shí)，ζ（s）=，并由解析延拓使其當(dāng)s≠1時(shí)有意義；Bk（k=0， 1，2，…）是伯努利數(shù)，滿足：

引理3[8]若f（x）在［a，b）上充分光滑，且在x=b處有分?jǐn)?shù)階泰勒展開式（3），則有

引理4[8]若f（x）在（a，b）上充分光滑，且在x=a和x=b處分別有分?jǐn)?shù)階泰勒展開式（2）和（3），則有

2 計(jì)算格式與誤差估計(jì)

考慮線性奇異兩點(diǎn)邊值問題：

其中：0＜α＜1，f（x）在區(qū)間[0，1]上非充分光滑或奇異.

對(duì)式（9）在區(qū)間[0，1]上直接積分，得到其積分形式的精確解為：

若式（10）右端可以準(zhǔn)確求出，則自然得到方程（9）的精確解.但通常情況下，需使用數(shù)值積分公式求其在[0，1]區(qū)間某點(diǎn)的數(shù)值解.下面使用改進(jìn)的梯形公式和龍貝格外推算法處理.對(duì)上述問題考慮以下2種情形：

情形1：f（ξ）僅在ξ=0處分?jǐn)?shù)階可導(dǎo)，其分?jǐn)?shù)階泰勒展開式如式（2）所示，則F（ξ）=ξ1-αf（ξ）也在ξ=0處分?jǐn)?shù)階可導(dǎo)，其分?jǐn)?shù)階泰勒展開式為

令h=x/n，由引理2有

式（10）右端第2個(gè)積分為正常積分，令G（ξ）=（ξ1-α-1）f（ξ）及h=（1-x）/n，由引理1有

由式（10）、式（12）和式（13），可以得出方程（9）的數(shù)值解及誤差估計(jì).

情形2：f（ξ）在ξ=0和ξ=a處有分?jǐn)?shù)階導(dǎo)數(shù)，分為2種情形，即0＜a＜x和x＜a＜1.

情形2.1：0＜a＜x.此時(shí)，F(xiàn)（ξ）分別在ξ=0和ξ=a處有分?jǐn)?shù)階導(dǎo)數(shù).F（ξ）在ξ=0處的分?jǐn)?shù)階泰勒展開式為式（11），F(xiàn)（ξ）在ξ=a處的分?jǐn)?shù)階泰勒展開式為

將式（10）的第1個(gè)積分分解，得

令h=a/n，由引理4有

對(duì)于正常積分 ，其計(jì)算公式及誤差估計(jì)仍為式（13）.由式（10）、式（13）和式（14），可以得出方程（9）的數(shù)值解及誤差估計(jì).

情形2.2：x＜a＜1.F（ξ）在ξ=0處的分?jǐn)?shù)階泰勒展開式仍為式（11）.G（ξ）在ξ=a處的分?jǐn)?shù)階泰勒展開式為

令h=x/n，I［F；0，x］的梯形公式及誤差估計(jì)式為式（12）.令h=（a-x）/n，由引理3有

再令h=（1-a）/n，由引理2有

由式（10）、式（12）、式（16）和式（17），可以得出方程（9）的數(shù)值解及誤差估計(jì).

前面給出3種情形的復(fù)合梯形公式的誤差估計(jì)式，在實(shí)際計(jì)算時(shí)，應(yīng)對(duì)梯形公式進(jìn)行外推，以消去誤差低階項(xiàng)，提高計(jì)算精度.下面給出修正的龍貝格外推算法，以式（14）為例，記

則龍貝格外推公式為

公式（19）中，精度控制要求為前后2次外推結(jié)果之差絕對(duì)值小于10-10，則可得方程（9）的高精度計(jì)算結(jié)果.

3 數(shù)值算例

例1 考慮奇異兩點(diǎn)邊值問題

方程（20）的精確解為u（x）=1/（1+x3/2）1/2.該例中，右端項(xiàng)f（x）在x=0處分?jǐn)?shù)階可導(dǎo)，其他點(diǎn)處充分光滑，其精確解u（x）在x=0處僅一階可導(dǎo).利用式（12）、式（13）進(jìn)行計(jì)算，并利用式（18）、式（19）進(jìn)行外推，記所得數(shù)值解為uh（x）.表1給出了一些節(jié)點(diǎn)處的精確解u（x）和數(shù)值解uh（x）及其誤差的絕對(duì)值.

表1 例1中計(jì)算結(jié)果及誤差Tab.1 Results and errors of example 1

例2 考慮兩點(diǎn)邊值問題：

方程（21）的精確解為u（x）=x5/4（1-x）4/3.該例中右端項(xiàng)f（x）在x=0和x=1處均奇異，其精確解u（x）在x=0和x=1處均僅一階可導(dǎo).利用式（12）、式（16）和式（17）進(jìn)行計(jì)算，并利用式（18）和式（19）進(jìn)行外推.表2給出了一些節(jié)點(diǎn)處的精確解u（x）和數(shù)值解uh（x）及其誤差的絕對(duì)值.

表2 例2中計(jì)算結(jié)果及誤差Tab.2 Results and errors of example 2

由表1可知，例1的最大絕對(duì)誤差為2.948 75× 10-13，由表2可知，例2的最大絕對(duì)誤差為5.334 62× 10-14.2個(gè)算例的數(shù)值解和精確解都非常接近，優(yōu)于Adomian分解法和插值方法，格式的計(jì)算效果令人滿意.且由例1和例2的計(jì)算結(jié)果可知，無論方程右端的源項(xiàng)是否奇異，本研究方法均得到了高精度的數(shù)值結(jié)果，因此該方法對(duì)于這類奇異兩點(diǎn)邊值問題具有良好的計(jì)算效果.

[1] CHAWLA M M，KATTI C P.Finite difference methods and their convergence for a class of singular two point boundary value problems[J]. Numerische Mathematik，1982，39：341-350.

[2]韓國強(qiáng).一類奇異方程兩點(diǎn)邊值問題的差分-樣條校正解[J].計(jì)算數(shù)學(xué)，1991，13（2）：187-192.

[3] HAN W M.Asymptotic error expansions for numerical solutions of onedimensional problems with singularities[J].Mathematical and Computer Modeling，1991，15（12）：67-76.

[4] SINGH R，KUMAR J.Solving a class of singular two point boundary value problems using new modified decomposition method[J].ISRN Computational Mathematics，2013，23（3）：1-11.

[5] SINGH R，KUMAR J.An efficient numerical technique for the solution of nonlinear singular boundary value problems[J].Computer Physics Communications，2014，185（4）：1282-1289.

[6]STAKGOLD I.Green’s function and boundary value problems[M].New York：John Wiley and Sons，2011.

[7] CEN Z D.Numerical study for a class of singular two-point boundary value problems using Green’s functions[J].Applied Mathematics and Computation，2006，183（5）：10-16.

[8] WANG T K，LI N，Gao G H.The asymptotic expansion and extrapolation of trapezoidal rule for integrals with fractional order singularities[J]. International Journal of Computer Mathematics，2015，92（3）：579-590.

[9] PHILIP J D，PHILIP R.Methods of numerical integration[M].2nd. ed.San Diego：Academic Press，1984.

（責(zé)任編校 馬新光）

Extrapolation of modified trapezoidal rule for a class of singular two-point boundary value problems

TANG Yongchao，WANG Tongke
（College of Mathematical Science，Tianjin Normal University，Tianjin 300387，China）

High precision numerical method for singular two-point boundary value problems（SBVPs）is studied.Firstly，the SBVPs are transformed into the computation of singular integrals.Secondly，the singular integrals are discretized by modified composite trapezoidal rule and the error asymptotic expansions are obtained for some kind of cases.Thirdly，a modified Romberg extrapolation algorithm is designed to obtain high precision numerical results based on the error expansions.Finally，numerical examples show the effectiveness of the algorithm.

singular two-point boundary value problem；composite trapezoidal integration rule；fractional Taylor's expansion；error asymptotic expansion；Romberg algorithm

1671-1114（2015）04-0005-03

O241.82

A

2015-01-22

國家自然科學(xué)基金資助項(xiàng)目（11471166）.

唐永超（1989—），男，碩士研究生.

王同科（1965—），男，教授，主要從事微分方程數(shù)值解方面的研究.