孫玉平，周德亮

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

目前計算地下水水頭函數(shù)的主流方法有有限元和有限差分法，用它們解決問題時，都要把計算域離散為有限單元或網(wǎng)格，最后求出水頭函數(shù)。徑向基函數(shù)配點(diǎn)法是用徑向基函數(shù)構(gòu)造插值函數(shù)并采用配點(diǎn)法作為離散方案的一種無網(wǎng)格方法，與空間維數(shù)無關(guān)，此方法不需要背景網(wǎng)格，效率高，形式簡單。本文利用的最小二乘配點(diǎn)法，是在徑向基函數(shù)配點(diǎn)法的基礎(chǔ)上，采用最小二乘法求解計算改進(jìn)而來，對計算域進(jìn)行節(jié)點(diǎn)離散，并布置輔助點(diǎn)，這樣計算結(jié)果比徑向基配點(diǎn)法精度更高。

1 最小二乘配點(diǎn)法基本原理

在配點(diǎn)型無網(wǎng)格法中，微分方程只在計算域內(nèi)各節(jié)點(diǎn)處嚴(yán)格滿足，因此可能會產(chǎn)生較大誤差。張雄［1］等人除了節(jié)點(diǎn)又在計算域內(nèi)引入了輔助點(diǎn)。近似函數(shù)仍然只通過節(jié)點(diǎn)構(gòu)造，但要求微分方程在所有節(jié)點(diǎn)和輔助點(diǎn)上滿足。此時方程個數(shù)大于未知數(shù)個數(shù)，需要用最小二乘方法求解，因此將此方法稱為最小二乘配點(diǎn)無網(wǎng)格法。

假設(shè)求解一個一維地下水穩(wěn)定流問題，其對應(yīng)滲流微分方程為:

式中:h為水頭，K為滲透系數(shù)，ε為垂直向補(bǔ)給，Ω為研究區(qū)，g(x)為給定的函數(shù)，Γ為研究區(qū)Ω的邊界。

對研究區(qū) Ω 布置節(jié)點(diǎn) x1，x2，…，xN，其中 x1，x2，…，xN0為內(nèi)部節(jié)點(diǎn)，xN0+1，xN0+2，…，xN0+M，為 M 個輔助點(diǎn)，xN0+M+1，xN0+M+2為一類邊界節(jié)點(diǎn)，總節(jié)點(diǎn)數(shù)為N=N0+M+2。

令所有節(jié)點(diǎn)滿足微分方程(1)，邊界上的節(jié)點(diǎn)滿足邊界條件(2)，

式(3)中有N0+M+2個方程，N0個未知數(shù)，需要用最小二乘法求解。為了保證解的精度，邊界條件必須嚴(yán)格滿足，因此將上式改寫成

式(4)中h1為邊界Γ上的節(jié)點(diǎn)參數(shù)組成的向量，h2為域Ω內(nèi)的N0+M個節(jié)點(diǎn)的參數(shù)組成的向量。上式可改寫為

式(5)對應(yīng)于邊界條件，(6)式對應(yīng)于微分方程，令邊界條件嚴(yán)格滿足，由(5)可解得

將(7)式代入(6)式中，得

式(8)中有N0個未知數(shù)，N0+M個方程，因此沒有通常意義下的解，需要在最小二乘意義下求解。在Matlab中用矩陣除法計算此類方程組，會自動利用最小二乘法求解，因此可以求出水頭h(x)的近似解。而徑向基函數(shù)配點(diǎn)法，只是在研究域內(nèi)不添加輔助點(diǎn)。

2 用最小二乘配點(diǎn)法求解非承壓含水層水流問題

自然界的多孔介質(zhì)很多都是非承壓的，最小二乘配點(diǎn)法可以很好地用來解決非承壓含水層中的水流問題。下面通過此方法在兩種非承壓條件下的應(yīng)用來說明這種方法的優(yōu)越性。

2.1 有垂直補(bǔ)給非承壓含水層一維穩(wěn)定流

如下圖1［2］所示，假設(shè)有一非承壓含水層，底板是水平不透水層，有降雨補(bǔ)給ε=0.004 m/d，含水層兩端有兩條河流流過，兩河相距5 000 m，上、下邊界分別為斷面線Ⅰ、Π，水頭分別為 h1、h2，設(shè) h1＞ h2，h1=280 m、h2=220 m，土壤滲透系數(shù)K=0.785 m/d，在計算區(qū)域內(nèi)滿足

分別使用解析方法、徑向基函數(shù)配點(diǎn)法、最小二乘配點(diǎn)法求解。圖2為不同方法的計算結(jié)果，從圖中可以看出，最小二乘配點(diǎn)法比徑向基函數(shù)配點(diǎn)法有更好的精度，更為接近解析解的結(jié)果。

圖1

圖2

2.2 無垂直補(bǔ)給非承壓含水層一維非穩(wěn)定流

如圖3［2］所示，設(shè)有一水庫，初始狀態(tài)下庫岸地下水位和庫水位齊平，使庫水位由驟降至，這樣庫岸內(nèi)地下水向庫內(nèi)排泄，潛水面形成一向水庫降落的曲線，隨著時間的增加，曲線向岸坡遠(yuǎn)處延伸，直至最終穩(wěn)定下來。在區(qū)域內(nèi)滿足

其中土壤的滲透系數(shù)K=0.864 m/d，給水度μ=0.03。

圖4為解析解、最小二乘配點(diǎn)法和徑向基函數(shù)配點(diǎn)法三種方法的計算結(jié)果。從圖中可以看出，隨著時間的增加，最小二乘配點(diǎn)法的計算結(jié)果比徑向基函數(shù)配點(diǎn)法更精確，穩(wěn)定性更好。

圖3

圖4

3 結(jié)語

上述算例表明:最小二乘配點(diǎn)法處理非承壓含水層穩(wěn)定和非穩(wěn)定流的水流問題是有效的，既能節(jié)省計算量又能保證精度，其近似解和解析解相比的誤差比徑向基函數(shù)配點(diǎn)法更小。此方法同樣適用于二維和三維流問題，方法是類似的。

［1］張雄，劉巖.無網(wǎng)格方法［M］.北京:清華大學(xué)出版社.2004.

［2］王連軍.工程地下水計算［M］.北京:中國水利水電出版社.2004.

［3］張雄，劉小虎，宋康祖，等.Least-square collocation mesh less method，Ⅰnt.J.Num.Meth.Engng，2001，51(9):1089 –1100.

［4］Kansa E J.Multiquadrics-a scattered data approximation scheme with applications to computational fluid dynamics Ⅰ.Ⅰ.Comp.Math.Applic.1990，19(8/9):147 – 161.

［5］Kansa E J.Multiquadrics-a scattered data approximation scheme with applications to computational fluid dynamics Ⅰ.Comp.Math.Applic.1990，19(8/9):127 – 145.

［6］周德亮，李靜.最小二乘配點(diǎn)法解一維非穩(wěn)定流問題［J］.地下水:2012，34(2):20-21.