郭冬冬，周德亮

(遼寧師范大學數(shù)學學院，遼寧 大連 116029)

?

承壓穩(wěn)定井流計算的單位分解配點法

郭冬冬，周德亮

(遼寧師范大學數(shù)學學院，遼寧 大連 116029)

[摘要]用單位分解配點法求解承壓含水層中地下水向井的穩(wěn)定流動問題，該方法擺脫了背景網(wǎng)格的束縛，是一種真正的無網(wǎng)格方法，根據(jù)具體模型計算后發(fā)現(xiàn)其不僅實施簡單，而且計算精度高。

[關鍵詞]單位分解；配點法；穩(wěn)定井流；無網(wǎng)格

無網(wǎng)格(Meshfree)方法是近幾十年來迅速興起的一種數(shù)值計算方法，它很好的克服了有限差分和有限元法等方法需要網(wǎng)格的缺陷，并利用一組散布在求解域中及其域邊界上的節(jié)點來表示(而非離散)該求解域和其邊界，具有收斂快，后處理方便等優(yōu)點。其中單位分解配點法[1]繼承了無網(wǎng)格方法的性質，一方面可以在子域內(nèi)調(diào)節(jié)節(jié)點的密度，也可以改變某個子域的形狀；另一方面局部逼近空間可以包含非多項式函數(shù),從而很好地局部逼近未知解。因此單位分解配點法具有靈活性且近似解有較高的精度。本文將單位分解配點法用于求解地下水的承壓穩(wěn)定井流動的問題。

1配點法

考慮如下定解問題

L{u(x)}=f(x)，x∈D

(1)

B{u(x)}=g(x)，x∈?D

(2)

其中設D是有界區(qū)域，?D是D的邊界，L是偏微分算子，B是邊界微分算子，f(x)，g(x)是給定的函數(shù)。

(3)

一般情況下中心點取配置節(jié)點，αj是待定系數(shù)，N是配點總數(shù)，φ(‖x-ξj‖2)表示徑向基函數(shù)[2]，常見的徑向基函數(shù)有

高斯徑向基函數(shù)：φ(r)=e-(εr)2(ε>0)；

緊支撐徑向基函數(shù)：φ(r)=(1-εr)4+(4εr+1) (ε>0)；

將(3)代入(1)(2)中，使其分別在內(nèi)節(jié)點，一類邊界節(jié)點和二類邊界節(jié)點上成立，這樣線性方程組便可寫成

Aα=b

(4)

其中α=(α1，α2，…，αN)T是待定的未知向量，

b=(f(x1)，…，f(xN))，g(xNi+1)，g(xNi+N1+1)，g(xN)T

s(x)=φ(x)α=φ(x)A-1u=Ψ(x)u，其中Aα=u，u=[u(x1)…u(xN)]T，

φ(x)=[φ(‖x-x1‖)，…，φ(‖x-xN‖)]，Ψ(x)=φA-1

這里的Ψj(x)具有如下性質

2單位分解配點法

?x∈DI(x)={j|x∈Dj}，I(x)≤K

集合I(x)表示的是包含x的子域的編號，包含x的子域個數(shù)不超過K個，這里常數(shù)K與子域的數(shù)量M無關。在單位分解配點法[1]中,定解問題的解函數(shù)u(x)的近似函數(shù)s(x)可以寫成下面形式：

(2.1)

其中sj(x)是u(x)在子域Dj內(nèi)的近似函數(shù)，wj是定義在Dj上的緊支正定權函數(shù),它可以通過Shepard的方法構造,即表示為

(2.2)

其中

(2.3)

是子域Dj上的緊支徑向基函數(shù).在(2.2)中我們選擇Wendland緊支徑向基函數(shù)

(2.4)

如果對每個節(jié)點Xi∈Dj，等式(2.4)中的函數(shù)sj(x)，j=1，…，M，都具有局部插值性質sj(xi)=u(xi)，那么全局單位分解近似值繼承了局部插值的性質，s(xi)也就可以寫成下面的形式：

利用單位分解配點法解決偏微分方程問題時，需要計算在空間微分算子在內(nèi)部節(jié)點上的值和邊界微分算子在邊界節(jié)點上的值，為此需要計算α階的導數(shù)，應用Leibniz定理，可以得到近似解的α階導數(shù)為

其中α和β是多重階數(shù)。

3對井的處理

考慮平面承壓穩(wěn)定井流混合問題[4]

其中,T為導水系數(shù),h(x，y)為水頭函數(shù),ε(x，y)為越流補給強度，q為補給量，qk為第k口井的開采量，ρk為第k口井的半徑。

通過實際計算可發(fā)現(xiàn)，和傳統(tǒng)數(shù)值方法一樣，如果直接從問題中進行離散，用徑向基函數(shù)配點法求解將會在井口附近產(chǎn)生很大的誤差，使結果不可用。原因是井口相對于求解域非常接近為一個點，從數(shù)學角度看就是一個奇點，目前解決這一問題的常規(guī)方法是將問題轉變?yōu)槿缦碌牡葍r問題

由此求得w,立即得到定解問題的解h=w+v。

圖1　計算域及節(jié)點配置圖

圖2　近似解

圖3　解析解

設ε(x，y)光滑無奇性。要保證數(shù)值求解的精度，函數(shù)vk除了要滿足文[6]給出的前述4個條件

4數(shù)值計算

在計算區(qū)域上按圖1所示布置節(jié)點。單位分解配點法得到的水頭函數(shù)近似解、解析解以及絕對誤差如圖2、圖3和圖4所示，本文中選用的是高斯徑向基函數(shù)，權函數(shù)選用緊支徑向基函數(shù)，參數(shù)取160，與解析解比較可知使用單位分解配點法配求解地下水向井流動的問題時具有靈活性，且近似解有較高的精度，由于該方法是在局部區(qū)域上應用配點法，因而形成的矩陣是稀疏的，從而減少了計算的費用。

圖4　絕對誤差

[1]Safdari-Vaighani A,Heryudono A,Larsson E.A Radial Basis Function Partition of Unity Collocation Method for Convection-Diffusion Equations Arising in Financial Applications[J]. Sci.comput,2014.10.4-7.

[2]Sharan M, Kansa E J, Gupta S. Application of multiquadric method fornumerical solution of elliptic partial differential equations[J]. Applied Mathematics Computation, 1997, 84: 275-303.

[3]Franke C, Schaback R. Solving partial differential equations by collocation using radial basis functions[J]. Applied Mathematics Computation, 1998, 93: 73-82.

[4]孫訥正.地下水流的數(shù)學模型和數(shù)值方法[M]. 北京:地質出版社.1981.7-48,110.

[5]楊天行.非穩(wěn)定流的有限元法中承壓水析奇與潛水逐步析奇的一個方法[J].長春地質學院學報.1978(4): 37-54.

[6]周培德.計算幾何——算法分析與設計[M]. 北京;清華大學出版社.1995.88-92.

[收稿日期]2015-11-16

[作者簡介]郭冬冬(1991-)，女，遼寧大連人，在讀碩士研究生，主攻方向：偏微分數(shù)值解法。

[中圖分類號]P641.139

[文獻標識碼]B

[文章編號]1004-1184(2016)01-0051-02