基于規(guī)劃算法的多層地基一維固結(jié)超孔隙水壓力估算*

劉紅軍，上官士青，王 虎

（中國(guó)海洋大學(xué)1.環(huán)境科學(xué)與工程學(xué)院；2.海洋環(huán)境與生態(tài)教育部重點(diǎn)實(shí)驗(yàn)室，山東 青島 266100）

太沙基于1923年提出了至今仍在工程界廣泛采用的飽和黏土一維固結(jié)方程，并用分離變量法給出了解析解。天然地基的成層性使得太沙基一維固結(jié)理論不能直接應(yīng)用，Schiffman等［1］求得了荷載隨時(shí)間呈線性增長(zhǎng)情況下多層地基固結(jié)的解析解；Gibson等［2］對(duì)大變形固結(jié)進(jìn)行了一系列研究，提出了有限差分求解的孔隙比為控制變量的一維大變形固結(jié)方程；謝康和等［3］獲得了多層地基在任意變荷載下的解析解；吳雄志等［4］考慮了滲流力的影響，對(duì)太沙基一維固結(jié)理論進(jìn)行了修正；Chen等［5］采用了DQM算法提出了多層地基在部分透水界面下且在任意荷載下的解；潘旦光等［6］提出了一種基于直接模態(tài)攝動(dòng)算法的適用于變參數(shù)土層的半解析解；吳健等［7］在Gibson地基中壓縮指數(shù)和滲透指數(shù)不同比值條件下，進(jìn)行了解析解推導(dǎo)和有限元數(shù)值分析；Kim等［8］在一維固結(jié)理論中引用了接觸面邊界關(guān)系以適應(yīng)大應(yīng)變的計(jì)算；Huang等［9］分析提出了一種耦合壓縮量和超孔隙水壓力以及非耦合單純計(jì)算超孔隙水壓力的有限元解法；謝新宇等［10］通過電阻網(wǎng)絡(luò)對(duì)太沙基一維固結(jié)解進(jìn)行了模擬。

有限差分在計(jì)算固結(jié)問題上簡(jiǎn)便實(shí)用，對(duì)于時(shí)間本身就是離散的；同樣采用有限元方法嚴(yán)格的計(jì)算固結(jié)問題，仍需要對(duì)時(shí)間進(jìn)行離散。曾攀等［11］給出了彈塑性多孔介質(zhì)的參變量變分原理（增量形式），可以得到x，z，t均離散的勢(shì)能方程，用于構(gòu)建有限元模型。張洪武等［12］采用二次規(guī)劃（數(shù)學(xué)規(guī)劃的一種）的算法求解了彈塑性biot固結(jié)，給出了時(shí)間卷積全量形式的參變量變分原理，通過時(shí)間卷積不需要迭代。無論通過何種方法考慮時(shí)間離散問題都較為復(fù)雜，因此筆者嘗試不考慮時(shí)間離散，建立一種不同固結(jié)度下超孔隙水壓力分布估算方法，這種情況下模型計(jì)算更為簡(jiǎn)便。

本文基于彈性地基及孔隙比與有效應(yīng)力線性關(guān)系的假設(shè)，以分層土總勢(shì)能公式作為規(guī)劃算法求解的目標(biāo)函數(shù)，并將固結(jié)離散為n個(gè)部分進(jìn)行，通過離散的太沙基解析式描述每個(gè)部分的固結(jié)過程中時(shí)間與超孔隙水壓力的對(duì)數(shù)關(guān)系，建立了數(shù)學(xué)規(guī)劃模型。該模型可依據(jù)平均固結(jié)度（以孔壓定義，下同）求解多層地基各單元超孔隙水壓力分布。該方法計(jì)算時(shí)不需要按時(shí)步進(jìn)行計(jì)算，模型求解簡(jiǎn)單，無需考慮矩陣運(yùn)算等求解中的數(shù)學(xué)問題。

1 計(jì)算思路及固結(jié)過程的離散

非線性規(guī)劃是具有非線性約束條件或目標(biāo)函數(shù)的數(shù)學(xué)規(guī)劃，是運(yùn)籌學(xué)的1個(gè)重要分支，廣泛用于經(jīng)濟(jì)學(xué)中生產(chǎn)力分配、車輛調(diào)度、概率論、圖論、拓?fù)鋵W(xué)等的計(jì)算分析。非線性規(guī)劃研究1個(gè)n元實(shí)函數(shù)在一組等式或不等式的約束條件下的極值問題，即為1個(gè)在可行域中求極值的問題。被求極值的函數(shù)稱為目標(biāo)函數(shù)，形成可行域的關(guān)系式稱為約束條件。在太沙基解的經(jīng)典假設(shè)中，地基土被假設(shè)為彈性，將地基土離散為n個(gè)單元后可以將其視為n個(gè)線性彈簧，這些線性彈簧 （代表骨架）的行為特征符合滲流約束下的最小勢(shì)能原理。對(duì)于固結(jié)問題可嘗試?yán)靡?guī)劃算法對(duì)固結(jié)過程中的要素進(jìn)行分配使整個(gè)系統(tǒng)勢(shì)能最小進(jìn)行求解。

本文數(shù)學(xué)模型的計(jì)算流程見圖1，首先將固結(jié)過程離散，分為n部分進(jìn)行，用太沙基解析解分別表示離散后超孔隙水壓力消散的區(qū)域，將這種情況下太沙基解析解中無實(shí)際物理意義的滲透系數(shù)k定義為滲流能力ki。在滲流能力ki總量為一定值的情況下，分配ki的大小使系統(tǒng)勢(shì)能最小，求得在不同固結(jié)度情況下，各土層單元的超孔隙水壓力分布。

圖1 數(shù)學(xué)規(guī)劃算法流程Fig.1 Flow chart for mathematical optimization algorithm

2 基于最小勢(shì)能原理的控制方程

假設(shè)地基為成層彈性體，在均布荷載作用下壓縮，其系統(tǒng)總體勢(shì)能可表示為

式中：∏為土體總勢(shì)能；F為均布外力；Ei為第i單元的彈性模量（壓縮模量）；si為第i層土體的壓縮量；s總為地表沉降，即所有分層土體壓縮量之和。

基于有效應(yīng)力原理，假設(shè)水體無法壓縮（水沒有彈性勢(shì)能），可以得到土體單元有效應(yīng)力與該單元超孔隙水壓力的關(guān)系：

式中：d為某單元有效應(yīng)力增量；du為某單元超空隙水壓力增量。

側(cè)限條件下孔隙比的變化增量與豎向有效應(yīng)力變化增量的關(guān)系：

式中：e為孔隙比，a為壓縮系數(shù)。

某一單元土層壓縮量與孔隙比的關(guān)系為

代入離散的勢(shì)能方程可得：

式中：li為第i單元厚度；q為均布荷載（方向?yàn)樨?fù)）。

由勢(shì)能方程的表述可見，土層總體的勢(shì)能主要由超孔隙水壓力增量du控制，但實(shí)際上式（7）僅是“完整的”彈性骨架的全量形式，與水相關(guān)的項(xiàng)目由于不能時(shí)間離散，不能完整表述，因此略去。模型主要通過解析解的孔隙水壓力時(shí)間關(guān)系約束彈性骨架的變形。Ei在勢(shì)能方程中可以看做1種權(quán)重項(xiàng)，由于本文擬采用的約束條件為太沙基解析解的一種離散，其中已經(jīng)包含了Ei的影響，因此式（7）中的繼續(xù)使用Ei作為1種權(quán)重將重復(fù)考慮彈性模量的影響，應(yīng)降低Ei對(duì)控制方程的權(quán)重。在此引入模量影響參數(shù)Ai來代替Ei的倒數(shù)，用于對(duì)勢(shì)能方程的修正。為了降低壓縮模量Ei對(duì)控制方程的影響，Ai應(yīng)是大小變化不大且受Ei影響的參數(shù)。

經(jīng)大量試算，Ai可取Ei的調(diào)和平均數(shù)，用于協(xié)調(diào)各單元?jiǎng)菽艿谋壤?。?/p>

3 基于解析解的時(shí)間約束

已知離散后的單層均勻地基固結(jié)排水太沙基解析解為式（2）。在本模型中，為便于計(jì)算，僅取式中第一項(xiàng)參與計(jì)算，同時(shí)為了統(tǒng)一時(shí)間因數(shù)為0時(shí)固結(jié)度為0，將8/π2也略去，則式其簡(jiǎn)化為

式中：Tv為時(shí)間因數(shù)。由于本文著重于超孔隙水壓力分布的計(jì)算，經(jīng)上述簡(jiǎn)化后，時(shí)間因數(shù)與固結(jié)度的關(guān)系已不準(zhǔn)確，但各單元的時(shí)間因數(shù)比例是統(tǒng)一的，可以用于計(jì)算某一固結(jié)度下的超孔隙水壓力分布。建議采用太沙基解析解的前三項(xiàng)作為約束，以提高計(jì)算精度，特別是對(duì)于固結(jié)度小的單元如果采用式（11）計(jì)算誤差較大。對(duì)于各單元可表述為

或采用式（2）前三項(xiàng)提高計(jì)算精度，表述為

由于式（11）、（12）等時(shí)間約束是用于限定各單元滲流能力ki之間的比例關(guān)系（見式13），公式中ki前的常數(shù)均可省略，即式（2）中ki前的常數(shù)π2/4可以省略。式（12）中常數(shù)項(xiàng)0.868是為了統(tǒng)一時(shí)間因數(shù)為0時(shí)固結(jié)度為0，對(duì)原解析解中的8/π2進(jìn)行的補(bǔ)償。

根據(jù)上述離散方法，每部分固結(jié)均消耗一定的滲透能力ki，令這些ki之和為一定值k總可作為一項(xiàng)約束。本模型中可令k總為任意常數(shù)C，常數(shù)C的大小僅影響時(shí)間t與固結(jié)度的關(guān)系，對(duì)于本文所關(guān)注的某固結(jié)度下各單元超孔隙水壓力分布比例不會(huì)產(chǎn)生影響。C一般取值為1～5，且隨單元數(shù)量變化。根據(jù)大量試算，k總為分組單元數(shù)量、各土層滲透系數(shù)的函數(shù)，求得k總的函數(shù)表達(dá)式可使本模型準(zhǔn)確計(jì)算地基中各單元的固結(jié)曲線，有待進(jìn)一步研究。

式中：為第i部分固結(jié)單元的加權(quán)調(diào)和平均滲透系數(shù)。

當(dāng)?shù)鼗鶠閱我痪|(zhì)土層時(shí)，上式可簡(jiǎn)化為

對(duì)于雙面排水的情況，可通過多次試算的方法確定土層中孔隙水壓力的的最高點(diǎn)，以該點(diǎn)作為界面將1個(gè)雙面排水的問題分割為2個(gè)單面排水進(jìn)行計(jì)算?；虿捎煤?jiǎn)單的0～1選擇約束來處理?？稍O(shè)ci為自上而下的滲透能力；為自上而下滲透的調(diào)和平均滲透系數(shù)；為第i土層單元向下滲流的滲流路徑長(zhǎng)；為第i單元至第n單元的加權(quán)算數(shù)平均壓縮模量。將某一單元向2個(gè)方向滲流所消耗的滲流能力聯(lián)立在一起，按式（15）進(jìn)行計(jì)算，模型在求解系統(tǒng)勢(shì)能的最小值時(shí)將對(duì)滲流方向進(jìn)行選擇，給出各單元變量j的值：

對(duì)于1個(gè)健康的模型，第i單元發(fā)生向上滲流時(shí)ci應(yīng)為0，向下滲流時(shí)ki應(yīng)為0。該規(guī)律與j的取值無關(guān)，j僅用于保證算法選擇的正確性，一般情況下可以將j省略，得式（16）。采用式（16）的計(jì)算結(jié)果中，對(duì)于任何單元的ci或ki必有1個(gè)等于0，符合實(shí)際情況：

4 模型求解

通過上述討論，將多層地基的一維固結(jié)問題轉(zhuǎn)化為1個(gè)求極值的數(shù)學(xué)規(guī)劃模型：

該模型是1個(gè)在可行域內(nèi)求極值的純數(shù)學(xué)問題，可采用MATLAB優(yōu)化工具箱（Optimization toolbox）或其他通用數(shù)學(xué)軟件求解。根據(jù)式（15），ui和dui可以表述為ki的函數(shù)，模型中未知數(shù)僅有ki，計(jì)算時(shí)調(diào)整t的大小可獲得不同固結(jié)度下的超孔隙水壓力分布。

5 算例

5.1 算例1

文獻(xiàn)［13］對(duì)層狀土體孔隙水壓力的空間變化規(guī)律進(jìn)行了研究，計(jì)算了滲透系數(shù)不同情況下孔隙水壓力的分布趨勢(shì)，筆者采用本文模型進(jìn)行了相同的計(jì)算，結(jié)果如圖2所示。該例僅假設(shè)上層土體與下層土體滲透系數(shù)（k1與k2）的差異，兩層土體的彈性模量相同。時(shí)間、具體滲透系數(shù)的大小、彈性模量的大小并不影響孔隙水壓力分布情況。由圖2可見，僅滲透系數(shù)變化時(shí)，模型計(jì)算結(jié)果與引證文獻(xiàn)基本一致。

5.2 算例2

圖2 算例1相同固結(jié)度下超孔隙水壓力沿深度分布曲線Fig.2 Excess pore pressure distribution curve of example 1at same degree of consolidation

表1 算例2地基參數(shù)表Table 1 Layered soil parameters of example 2

5.3 算例3

采用文獻(xiàn)［14］中的算例探討本文算法在壓縮模量變化較大的地基中的精度問題。該地基為雙層，各土層厚度及力學(xué)性質(zhì)見表2，分為2種情況計(jì)算。情況1為土層1覆在土層2上；情況2為土層2覆在土層1上。2種情況均為單面透水，均布荷載，每層土體離散為10個(gè)單元參與計(jì)算。計(jì)算結(jié)果如圖4，5所示。

本文算法由于采用了基于離散太沙基一維固結(jié)解析解的約束，計(jì)算結(jié)果主要受到固結(jié)系數(shù)Cv控制。在本例中Cv為一恒定值，如不采用Ai進(jìn)行修正，超孔隙水壓力分布的計(jì)算結(jié)果應(yīng)等同于單層均質(zhì)地基的太沙基解，即類似于函數(shù)sin（πz/2H），體現(xiàn)不出雙層地基中k和mv的變化。通過Ai修正后，本文模型的解基本符合文獻(xiàn)［14］中解的變化規(guī)律，說明本文中Ai的修正是有效的。

圖3 算例2相同固結(jié)度下超孔隙水壓力延深度分布曲線Fig.3 Excess pore pressure distribution curve of example 2at same degree of consolidation

表2 算例3地基參數(shù)表Table 2 Layered soil parameters of example 3

圖4 算例3情況1某固結(jié)度下超孔隙水壓力延深度分布Fig.4 Excess pore pressure distribution curve of example 2case 1at same degree of consolidation

圖5 算例3情況2相同固結(jié)度下超孔隙水壓力延深度分布Fig.5 Excess pore pressure distribution curve of example 2case 2at same degree of consolidation

在計(jì)算情況1的結(jié)果中出現(xiàn)了下部單元超孔隙水壓力小于上部單元的情況（壓縮模量突變處），與實(shí)際情況不符，建議在模型中加入一項(xiàng)姿態(tài)約束dui＋1≥dui，使計(jì)算結(jié)果更符合實(shí)際情況。在情況1中本文算法與文獻(xiàn)［14］的差分解誤差較大，約為15%，表現(xiàn)為孔隙水壓力線性增長(zhǎng)的分界點(diǎn)前移了近一個(gè)單元。將情況1下的地基離散的單元數(shù)量增加到20后，解的精度有所提高，但孔隙水壓力線性增長(zhǎng)的分界點(diǎn)與10單元的計(jì)算結(jié)果變化不大。如果對(duì)計(jì)算精度要求較高，本文方法不適用于壓縮模量變化劇烈的地基。

6 結(jié)論

（1）以最小勢(shì)能原理作為目標(biāo)函數(shù)，太沙基解析式的一維離散作為約束的數(shù)學(xué)規(guī)劃模型用于計(jì)算多層地基一維固結(jié)中某一固結(jié)度下超孔隙水壓力分布是可行有效的，且不需按時(shí)步計(jì)算。其計(jì)算結(jié)果主要受土層固結(jié)系數(shù)Cv影響，并通過目標(biāo)函數(shù)中的參數(shù)A修正，修正后計(jì)算結(jié)果與差分算法計(jì)算結(jié)果基本吻合。

（2）規(guī)劃算法在計(jì)算中對(duì)壓縮模量進(jìn)行了算術(shù)平均處理，如果對(duì)計(jì)算要求精度較高，該方法不適用于各土層壓縮模量變化劇烈的地基。

（3）根據(jù)大量試算，新定義的滲流能力的總和k總為分組單元數(shù)量、各土層滲透系數(shù)等參數(shù)的函數(shù)，準(zhǔn)確確定k總的表達(dá)式可使該模型計(jì)算固結(jié)度于時(shí)間的關(guān)系，有待進(jìn)一步研究。

［1］Schiffman R L，Stein J R.One-dimensional consolidation of layered systems［J］.J Soil Mech and Found，1970，964：1499-1504.

［2］Gibson R E，Schiffman R L，Cargill K W.Theory of one-dimensional consolidation of saturated clays.Finite non-linear consolidation of thick homogeneous layers［J］.Canadian Geotechnical Journal，1981，18（2）：280-293.

［3］謝康和，潘秋元.變荷載下任意層地基一維固結(jié)理論［J］.巖土工程學(xué)報(bào)，1995，17（5）：80-85.

［4］吳雄志.考慮滲流力時(shí)對(duì)太沙基一維固結(jié)理論的修正［J］.巖土力學(xué)，2004，25（8）：1279-1282.

［5］Chen R P，Zhou W H，Wang H Z，et al.One-dimensional nonlinear consolidation of multi-layered soil by differential quadrature method［J］.Computers and Geotechnics，2005，32（5）：358-369.

［6］潘旦光，樓夢(mèng)麟.變參數(shù)土層一維固結(jié)的半解析解［J］.工程力學(xué)，2009，26（1）：58-63.

［7］吳健，謝新字，朱向榮.飽和土體一維復(fù)雜非線性固結(jié)特性研究［J］.巖土力學(xué)，2010，31（1）：81-86.

［8］Kim H J，Jose L M.Numerical analysis of one-dimensional consolidation in layered clay using interface boundary relations in terms of infinitesimal strain［J］.International Journal of Geomechanics，2011，11（1）：72-77.

［9］Huang J，Griffiths D V.One-dimensional consolidation theories for layered soil and coupled and uncoupled solutions by the finite-element method［J］.Geotechnique，2010，60（9）：709-713.

［10］謝新宇，王龍，劉開富.飽和土一維固結(jié)的電阻網(wǎng)絡(luò)模擬［J］.巖土力學(xué)，2011，32（11）：3455-3460.

［11］曾攀，錢令希.巖土中彈塑性滲透固結(jié)問題的參變量變分原理［J］.力學(xué)學(xué)報(bào)，1991，23（4）：484-490.

［12］張洪武，鐘萬勰.土體固結(jié)彈塑性分析的參數(shù)二次規(guī)劃理論及有限元解［J］.巖土力學(xué)，1995，16（1）：35-45.

［13］欒茂田，錢令希.層狀飽和土體一維固結(jié)分析［J］.巖土力學(xué)，1992，13（4）：45-56.

［14］Pyrah I C.One-dimensional consolidation of layered soils ［J］.Geotechnique，1996，46（3）：555-560.