帶自然邊界條件多元多項式樣條插值及微分方程數(shù)值解

徐應(yīng)祥

（中山大學(xué)新華學(xué)院，廣東 廣州 510520）

一般情況下，大多的偏微分方程很難求得其精確的解析解，因此，如何求其滿足要求的近似解，也即數(shù)值解，成為科學(xué)和工程計算中的最重要內(nèi)容［1］。求解微分方程的數(shù)值方法有很多種，如較流行的差分法、有限元法、邊界元法、譜方法與擬譜方法、辛方法等［1－2］。近年來，逐漸興起另一種求偏微分方程數(shù)值解的方法：無網(wǎng)格方法［3－4］。這種數(shù)值方法不需要像有限元法一樣需要對微分方程的求解區(qū)域進行剖分，只是基于在求解區(qū)域中選定的某一組散亂分布的點處的信息來確定真解的一種近似，因此其在許多條件下，特別是在高維情形的實現(xiàn)有時要比有限元方法方便。

由于現(xiàn)在對于多元散亂數(shù)據(jù)擬合方法［5－17］的研究也取得了許多成果，而散亂數(shù)據(jù)擬合本身就是根據(jù)一組給定點確定一個合適的逼近函數(shù)，因此將散亂數(shù)據(jù)擬合的方法應(yīng)用于偏微分方程數(shù)值求解中，可能會帶來一些新的可靠的數(shù)值方法。文獻［1－4，15］討論了散亂數(shù)據(jù)擬合方法在求解偏微分方程中的應(yīng)用，特別是文獻［15］在基于徑向基的散亂數(shù)據(jù)擬合方法基礎(chǔ)上，對一般的線性偏微分方程給出了基于徑向基函數(shù)的幾種偏微分方程數(shù)值求解的泛函信息插值方法的全部過程，并對基于徑向基函數(shù)的最小二乘方法、配置法、Galerkin方法及Golberg法的基本思想做了總結(jié)，這些方法實質(zhì)都是無網(wǎng)格方法。

本文從任意d元散亂數(shù)據(jù)出發(fā)，給出一種目標泛函較為簡單、更符合實際的散亂數(shù)據(jù)點的多元多項式自然樣條。用這種自然樣條考慮任意d元散亂數(shù)據(jù)的插值問題，則該問題的自然樣條解與以往不同，其基函數(shù)有良好的邊界條件與簡單的表達式。此外，以這種自然樣條為基礎(chǔ)，可以構(gòu)造出一種新的偏微分方程數(shù)值方法，數(shù)值算例表明這種方法是有效的。

1 散亂數(shù)據(jù)多元自然樣條插值

設(shè)表示d個正整數(shù)集的直積，對d維多重指標p＝（p1，p2，…，pd）∈，記用Rd表示d維歐氏空間，若u（x）是d元函數(shù)，其｜p｜階導(dǎo)數(shù)記為：

令i＝（i1，i2，…，id）∈，記單項式（x－a）i＝以及i?。絠1！i2！…id！，再記dx＝dx1dx2…dxd及求和符號

給定Rd中含有N個兩兩不同散亂數(shù)據(jù)點的集合S＝｛xi｜xi＝（x1i，x2i，…，xdi）∈Rd，i＝1，…，N｝，設(shè)方型域以及Sobolev空間X＝Hs（Ω）＝Ws，2（Ω）（s＝｜p｜），顯然滿足連續(xù)嵌入條件s〉d／2。

記Y＝L2（Ω），令Td：X→Y是一個從X到Y(jié)的線性算子，定義為：

再記BN?X是滿足如（1）式的自然邊界條件的所有u∈X構(gòu)成的子集。

其中，ipk＝（i1，…，ik－1，pk，0，…，0）；pik＝（p1，…，pk－1，ik，0，…，0）；axk＝（a1，…，ak－1，xk，…，xd）；xaxk＝（x1，…，xk－1，ak，xk＋1，…，xd）；ik＝0，1，…，pk－1。

定義A：X→Z是一個插值算子，滿足：Au＝（u（x1），u（x2），…，u（xN））。記Aw＝｛u｜Au＝w，u∈X｝為X中所有滿足插值條件Au＝w的函數(shù)u的集合，考慮如下插值問題［11］。

問題1 求σ（x）∈X，滿足：

問題1稱為任意維散亂點帶自然邊界條件多元樣條插值問題，問題的解稱為散亂點帶自然邊界條件多元插值樣條，簡稱多元自然插值樣條。

定理1 算子T的核空間為：

反之，用數(shù)學(xué)歸納法。

（1）當d＝2時，如果u∈N（T2）且滿足邊界條件（1）式，則T2（u）＝0，即于是有：

（2）假設(shè)當d＝m時，如果u∈N（Tm）且滿足（1）式，則u（x。 那么當d＝m＋1時，如果u∈N（Tm＋1）且滿足邊界條件（1）式，則有）＝0，所以），從而由歸納假設(shè)可知：

于是對xm＋1積分可得：

綜上可知，對一切d，如果u∈N（Td）且滿足邊界條件（1）式，則

2 散亂點多元插值自然樣條構(gòu)造

定理2（特征定理）σ∈X是多元插值自然樣條的充分必要條件如下：對一切u∈X，Au＝0成立。

證明 記J（u）＝‖Tu‖2，u∈X，Au＝w。一方面，根據(jù)泛函取極值的必要條件，若泛函J（σ＋εu）在σ處取極值，則其在σ處的變分δJ＝0，而〉，故 有〈Tdσ，Tdu〉＝0。另一方面，若〈Tdσ，Tdu〉＝0，則J（σ＋εu）＝J（σ）＋ε2‖Td（u）‖2≥J（σ），即J（u）在σ處取極小。

綜合上述2個方面，定理得證。

引理1 如果σ∈X是多元插值自然樣條，則必存在系數(shù)di與ki（i＝1，…，N）使得：

其中，為Td的共軛算子，且有〈ki，u〉＝u（xi）。

證明 記由N（A）＝｛u∈X｜Au＝0｝，即是算子A的核空間，并記N（A）⊥為N（A）在X中的正交補空間，由特征定理可知，當u∈N（A）時有Tdσ，u〉＝ 〈Tdσ，Tdu〉＝0，因 此Tdσ∈N（A）⊥。

對?u∈X，作泛函λi（u）＝u（xi）（i＝1，…，N），則由Riesz表示定理知存在ki∈X，使得λi（u）＝〈ki，u〉，因此A（u）＝（〈ki，u〉）i＝1，…，N。故若u∈N（A），則有〈ki，u〉＝0，從而有ki∈N（A）⊥，因此｛ki｝i＝1，…，N組成了N（A）⊥的基，而∈N（A）⊥，所 以 存 在 系 數(shù)di使 得

定理3（構(gòu)造定理） 設(shè)σ（x）是散亂數(shù)據(jù)多元自然樣條，則其具有如下的顯式及緊湊的表達式：

其中，gi（x）為自然樣條基函數(shù)，即

并且有：

其中

（4）式中的系數(shù)cj及di滿足方程組：

其中，Q＝（gj（xi））對稱，而矩陣B＝），C＝（cj）T，D＝（di）T，z＝（w）T，0是零矩陣。

證明 顯然xj∈N（Td）（ji≤pi－1－1，i＝1，…，N），所以由引理1可知：

此結(jié)果即為BTD＝0。

對任意X中的元素u，在a1點Taylor展開，則有：

類似地，由Taylor展開有：

由自然邊界條件（1）式，將上述各式依次代入u中有：

及自然邊界條件（1）式，故有：

對Gi1，由 自 然 邊 界 條 件xd）＝0，可知對任意的x2，…，xd有：

類似地可得c2i2，…，cdkd，形式與c1i1類似，這里不再贅述。再由插值條件Aσ＝w便可得系數(shù)cj及di滿足方程組（6）式。

在方程組（6）式中，再進行簡單驗證可知Q是對稱的，從而其系數(shù)矩陣是對稱的。

3 基于自然樣條的偏微分方程數(shù)值方法

對于給定的偏微分方程，在其求解區(qū)域Ω?Rd中任意選定一組點，那么由前述可知，以多元多項式自然樣條基N）可以作成一個有限維的樣條函數(shù)空間SN。如果在Galerkin方法中取Vn＝SN，則可以得到一種基于多元多項式自然樣條的Galerkin方法。下面以一個具體的含有三元偏微分方程的求解為例來說明這種方法。

其中

令

其中，Γi（i＝1，2，3）表示Ω的邊界，且xi＝bi，aj≤xj≤bj（j＝1，2，3，j≠i），則方程的廣義解u∈H1（Ω）是如（8）式變分問題的解。

其中，（Ω）為H1（Ω）中滿足＝0的所有函數(shù)組成的集合。

在Ω中任意選取一組兩兩不同的點則與這些點相應(yīng)的多元多項式樣條基為gi（x）（i＝1，…，N），由其性質(zhì)可知滿足＝0，因此其張成的函數(shù)空間SN?（Ω）。于是用SN代替變分問題 （8）式中的（Ω），記uN＝），vn依次取為gi代入方程（8）式，則得關(guān)于ci的線性方程組：

解此線性方程組得ci，再代入uN，即得邊值問題（7）式的近似解。

在Ω中隨機取220個散亂點，分布如圖1所示。選擇三元p1＝p2＝p3＝2的自然樣條gi（x）張成的子空間SN進行計算，對計算所得的結(jié)果，列出當x3＝0.08、0.16、0.24時3種情形下的精確解u、近似解uN及誤差r＝u－uN在空間中的截面圖形，如圖2所示。列出這3種情形下的數(shù)值近似解的最大誤差與平均誤差，見表1所列。

圖1 220個散亂點的分布

圖2 3種情形下精確解u、近似解uN及誤差r＝u－uN的曲面

表1 數(shù)值解的誤差

通過此例結(jié)果可以看出，用這樣的基于散亂數(shù)據(jù)的多元樣條自然樣條函數(shù)構(gòu)造的Galerkin方法是有效的，但針對不同的問題還要注意如下的問題：

（1）在邊值問題（7）式中的強加邊界條件是齊次邊界條件，而基函數(shù)gi（x）恰好可能選擇使其滿足相應(yīng)的齊次邊界條件，因此可以直接進行Galerkin近似。如果邊值問題中的邊界條件是非齊次的，那么基函數(shù)不能直接滿足這樣的邊界條件，從而要使近似解也能滿足或近似地滿足給定的邊界條件就需要對邊界條件進行如何處理加以考慮。對于自然邊界條件，則根據(jù)多元函數(shù)的Green公式等適當構(gòu)造泛函，使其在變分過程能自動滿足［16－18］。而強加的非齊次邊界條件可以用補償法將其改造成一個近似的自然邊界條件，或用Lagrange乘子法將其吸收進變分過程，從而達到去掉強加邊界的目的［16－18］。

（2）對于基函數(shù)gi（x），由其構(gòu)造可知，其在方型區(qū)域某些邊界上能夠自動滿足一些齊次邊界條件。因此，在具體問題中，對于相應(yīng)的方型區(qū)域上的齊次邊界條件可以選擇或重新構(gòu)造基函數(shù)使其能滿足齊次邊界條件，從而盡可能減少因要對邊界的處理而要在變分方程中增加一些項。

（3）散亂點的選取可以在區(qū)域內(nèi)部，也可在邊界上，但要盡可能使散亂點均勻分布，從而使近似解更能反應(yīng)方程所反映問題的信息和特征。

［1］ Reutskiy S Y.A boundary method of Trefftz type for PDEs with scattered data［J］.Engineering Analysis with Boundary Elements，2005，29：713－724.

［2］ Bouhamidi A，Jbilou K.Meshless thin plate spline methods for the modified Helmholtz equation［J］.Comput Methods Appl Mech Engrg，2008，197：3733－3741.

［3］ Duan Y.A note on the meshless method using radial basis functions［J］.Computers and Mathematics with Applications，2008，55：66－75.

［4］ Firoozjaee A R，Afshar M H.Discrete least squares meshless method with sampling points for the solution of elliptic partial differential equations［J］.Engineering Analysis with Boundary Elements，2009，33：83－92.

［5］ Kazhdan M，Bolitho M，Hoppe H.Possion surface reconstruction［C］／／Proceeding of Eurogaphics Symposium on Geometry Processing，Cagliari，Italy，2006：61－70.

［6］ Guan L T，Liu Bin.Surface design by natural splines over refined grid points［J］.Journal of Computational and Applied Mathematics，2004，163（1）：107－115.［7］ Lai M J，Schumaker L L.Spline functions over triangulations［M］.London：Cambridge University Press，2007：1－10.

［8］ Lai M J.Multivarariate splines for data fitting and approximation［M］／／Neamtu M，Schumaker L L.Approximation Theory XII，San Antonio 2007.Brentwood：Nashboro Press，2008：210－228.

［9］ Zhou T H，Han D F，Lai M J.Energy minimization method for scattered data Hermit interpolation［J］.Applied Numerical Mathematics，2008，58：646－659.

［10］ Chui C K，Guan L T.Multivariate polynomial natural spline for interpolation of scattered data and other applications［C］／／Conte A，et al.Workship on Computational Geometry.World Scientific，1993：77－98.

［11］ Guan L T.A local basis for bivariate polynomial natural splines of scattered data［C］／／Guangzhou International Symposium of Computational Mathematics，1997：10－20.

［12］ Guan L T.Bivariate polynomial natural spline interpolation algorithms with local basis for scattered Data［J］.Jounal of Computational Analysis and Applications，2003，45I（1）：77－101.

［13］ Laurent P J.Approximation et optimization［M］.Paris：Hermann，1972：50－56.

［14］ Bezhaev A Y，Vasilenko V A.Variational theory of splines［M］.New York：Kluwer Academic／Plenum Publishers，2001：45－60.

［15］ Wu Z.Model，method and theory of scattered data fitting［M］.Beijing：Science Press，2007：75－90.

［16］ 陳艷秋，王家正.三角網(wǎng)格上的Lagrange－Stieltjes型有理插值［J］.合肥工業(yè)大學(xué)學(xué)報：自然科學(xué)版，2013，36（9）：1144－1147.

［17］ 李少華，唐 爍，丁 晨.S－型非均勻有理 TC－雙弧的G2Hermit插值［J］.合肥工業(yè)大學(xué)學(xué)報：自然科學(xué)版，2014，37（12）：1523－1527.