均質黃土壩變形計算參數全局敏感性分析

孫一清 黃 天 王桂智 薛宇陽 徐力群方衛華 馮亞新 黃章鑫 金洪杰 崔衛東

(1.河海大學 水利水電學院，南京 210098;2.江蘇省水利勘測設計研究院有限公司，江蘇 揚州 225127;3.上海市堤防泵閘建設運行中心，上海 200080;4.水利部南京水利水文自動化研究所，南京 210008;5.江蘇省東臺市自然資源和規劃局土地整理中心，江蘇 東臺 224200)

大壩是水資源管理的重要基礎設施，中國有9.8萬多座大壩，其中土壩有8萬多座.土壩具有對各類環境的高度適應性，并可以在壩址區就地取材，建造方便等優點，因而被廣泛采用[1].實際建造過程中，隨著壩體施工高度不斷上升，土壩工程變形情況及其安全性會隨之改變.雖然壩體多采用分層施工的方式，但完工后，大壩整體仍然會因為自重的原因產生相應的位移變化.蓄水期大壩的變形又會在自重與水荷載的共同作用下繼續發展，因此對土壩變形分析是十分必要的.

土壩的主體部分多為當地土料，土料具有粒狀非連續固體顆粒的特性，在受到自重和水壓力作用下會產生變形，通常土體的變形分析采用非線性有限元法進行計算，其中鄧肯-張E-B模型是最常用的材料本構模型之一[2-5].鄧肯-張E-B模型參數通常可以根據已建造工程的經驗類比得到，但這樣得到的參數往往精度不夠，隨著施工技術的發展，后來逐步演變為由室內、室外三軸試驗測得，但測得的參數多與實際巖體參數有較大的偏差[6-7].實踐表明，不管使用何種非線性模型，都必須有足夠精確的參數，否則無法達到理想的效果.可見，如何準確、便捷、合理地獲得土料的鄧肯-張E-B模型參數，減少大壩設計施工過程中的工作量是進行土壩變形分析時需要解決的關鍵問題[8].如今大壩參數大都需要通過算法將參數反演問題轉化為一個優化問題來獲取，但這需要大量的試驗數據來為算法的訓練集提供支撐.

近年來，正交試驗法作為一種可代替全面試驗的便捷、科學的試驗方法被廣泛運用于工程中[9-13]，例如王瑞駿等[10]以公伯峽水電站面板堆石壩作為研究對象，利用正交試驗設計對模型參數進行研究，為流變模型參數的選取提供支撐.張文兵等[11]針對大壩熱-流耦合模型，運用正交試驗設計研究各參數對溫度場的影響，發現滲透率和土體比是模型中的重要參數.李炎隆等[12-13]基于正交試驗法，分別研究了堆石料模型參數對堆石體變形和超高鑲嵌面板堆石壩高模量區模型參數對壩體沉降、面板撓度及順坡向應力的影響.目前，應用到大壩工程的正交試驗法多是不考慮交互作用的，仍屬于局部敏感性分析的范疇，而實際工作中，各因素之間的關聯對于目的指標的影響是不能忽視的;除此之外，由于鄧肯-張E-B模型中參數眾多，標準正交表的列數無法同時容納所有鄧肯-張E-B模型中參數及參數之間的交互作用列，因此有必要對大壩變形影響不顯著的鄧肯-張E-B模型參數進行篩選.本文以某均質黃土壩為例，建立有限元模型對大壩變形進行模擬，并通過Plackett-Burman試驗設計對鄧肯-張E-B模型中參數進行顯著性篩選，再用正交試驗對篩選參數進行全局敏感性分析，研究參數及參數交互作用下對土壩變形的敏感性排序.以期為均質土壩鄧肯-張E-B模型參數選取提供參考.

1 鄧肯-張E-B模型原理

通常均質黃土在數值計算中的應力應變關系是非線性的，鄧肯-張E-B模型即為其中之一[14].

根據變形模量的定義，可得切線模量Et為:

式中:Rf為破壞比;S為剪應力比;Ei為初始剪切模量.

式中:σ1為大主應力;σ3為小主應力.

式中:pa為標準大氣壓;K、n可由試驗確定.

根據摩爾-庫侖破壞準則有:

式中:c為黏聚力;φ為內摩擦角.

將S、Ei分別代入式(1)，得切線模量表達式為

2 敏感性分析方法

鄧肯-張E-B模型共有10個計算參數c、φ、Δφ、Rf、K、Kb、n、m、Kur、nur.實際研究中，并不會選擇全部參數來研究其對大壩變形的敏感性分析，由于壩體處于壓實狀態且為散粒體，故Kur、nur、C不參與討論，除此之外將壩料密度ρ納入考慮范圍.

實際研究中發現，各因素的變化不僅會對目的指標產生影響，而且其相互間的聯動也會對目標指標產生影響，在試驗設計中稱這種聯動為交互作用，有時不同因素相互間的聯動對于目的指標的影響可能會超越因素本身.交互作用主要通過因素聯動的數量來區分:一階交互作用主要發生在任意兩個因素之間，高階交互作用通常發生在3個或者3個以上的因素之間.通常在考慮交互作用的全局敏感性分析中，高階交互作用對于目的指標的顯著性影響微乎其微，一般只研究少數因素的一階交互作用.正交試驗通常也只需考察兩兩因素之間的交互作用，而現有的標準正交表所考察的列數無法同時滿足八個計算參數以及八個參數兩兩之間的交互作用，因此需要使用Plackett-Burman試驗設計篩選出更為顯著的參數，以科學合理地減少正交表所需要的列數，隨后再進行考慮交互作用的正交試驗.

由于極差分析法具有一定的局限性，不能區分試驗結果的差異具體是由什么因素引起[15].本文使用的兩種試驗設計方法的結果均采用方差分析法，具體原理如下:

式中:ST為總離差平方和;Si為第i個因素的離差平方和;Se為誤差的離差平方和;k為因素i在水平j下重要試驗的次數為因素i在水平j下重復試驗的第g次所對應的試驗結果，{|i∈(1，m)，j∈(1，r)，g∈(1，k)}.其中

df為各因素自由度，dfT、dfi、dfe分別表示ST、Si、Se的自由度.

各因素的均方為:

誤差的均方為:

檢驗敏感性統計量F值為:

若Fi≥F0.01，說明該因素非常顯著，記為**;若F0.01＞Fi＞F0.05，說明該因素顯著，記為*;若F0.05≥Fi≥F0.1，說明該因素有影響，記為⊙;若F0.1＞Fi＞F0.2，說明該因素有一定影響，記為Δ;若F0.2≥Fi，說明該因素無影響.

2.1 Plackett-Burman試驗設計方法

Plackett-Burman試驗設計主要被運用于多個因素的篩選，該方法旨在從多個未知其對于目的指標影響強弱的因素中，篩選出對目的指標影響顯著的，以便后續的研究.在多個因素中，每個因素取高低兩個水平，通過試驗設計表格研究因素高低變動對于目的指標變動的影響.通過確定影響顯著的因素，達到篩選因素的目的，相應地減少后續研究中的工作量，因此，Plackett-Burman試驗設計被廣泛應用于實際工作中，起到了很好的效果.

Plackett-Burman試驗設計按規則生成，排列可具有不唯一性，對于N次實驗至多可研究(N-1)個因子，但實際因子應該不多于(N-4)個，保留3個以上虛擬變量用以估計試驗誤差.其試驗次數N為4的倍數，但N為2的冪的時候除外.常用的N為N=12、20、24、28、36等.

2.2 正交試驗設計方法

正交試驗設計屬于試驗設計方法的一種，根據是否考慮因素間的交互作用可劃分為兩類.不考慮交互作用的正交試驗屬于局部敏感性分析方法，而考慮交互作用的正交試驗屬于全局敏感性方法.不管是哪種正交試驗設計，都需正交表的參與，正交表用Ln(rm)表示，試驗的目標因素稱為目標指標.標準正交表如圖1所示.每張正交表都是通過正交原理設計出的標準化的試驗安排表格，具有科學性和均勻分布的特點，可涵蓋所有的試驗情況，如L8(27)、L16(215)等.

圖1 標準正交表

3 工程實例

3.1 位置描述

某水庫位于黃土中山丘陵區，該庫大壩全長575 m，壩頂高程1 557 m，最大壩高96 m，防浪墻頂高程1 558.2 m，水庫正常蓄水位1 554 m.水庫總庫容為4587×104m3.大壩河床段典型橫剖面圖如圖2所示.

圖2 河床段典型橫剖面圖

3.2 水庫大壩有限元模型

計算坐標系規定為:將河流走勢方向定義為X軸，沿河流走勢方向指為正向，取計算模型最左下為X軸零點;將壩軸線方向定義為Y軸，由大壩右岸指向左岸為正向，取右岸壩軸線端為Y軸零點;將垂直方向定義為Z軸，上方為正向，模型底部高程1 266 m為Z軸零點.

計算模型的邊界范圍:下邊界從大壩建基面向下取200 m，即從大壩建基面高程1 466 m處再垂直向下取兩倍壩高作為地基，地基厚度約為200 m;上游邊界從大壩上游壩腳處向X軸負坐標方向取1.5倍壩高，約為150 m;下游邊界從大壩下游壩腳處向X軸正坐標方向取1.5倍壩高，約為150 m.

底部約束有限元模型的地基下方，位移約束分別約束地基x方向兩側與y方向兩側.壩體典型斷面有限元網格如圖3所示.

圖3 典型斷面有限元網格

3.3 計算參數

壩體黃土與風化泥巖混合料、排水反濾層、開挖翻填料、開挖壓重等各料區土料的本構模型采用鄧肯-張E-B模型.壩體各料區的計算參數見表1.由于水庫黃土均為濕陷性黃土，施工前會發生深厚覆蓋層滲漏及不均勻變形問題，經處理后的壩址區具備建壩條件.

表1 壩體填筑材料計算參數(鄧肯-張E-B模型)

3.4 分級加載

大壩建造全部共有18級加載完成壩體填筑和蓄水.分級加載及蓄水過程如圖4所示.

圖4 有限元計算分級加載及蓄水過程

3.5 數值模擬結果

計算時，首先加載基巖和未清除的覆蓋層，并在分級加載壩體之前將結點位移初始化為零，僅保留單元應力，從而獲得地基初始應力場.下文所述的位移均是指開始施工填筑以后的位移.

通過數值模擬計算大壩竣工期及運行期壩體典型斷面位移分布如圖5所示.具體數值模擬計算成果見表2.

圖5 典型斷面位移分布

表2 堆石體位移二維計算分析結果匯總 (單位:mm)

整體來看，該壩的變形分布符合一般土石壩的分布特點，竣工期壩體水平位移呈對稱分布，運行期壩體蓄水后，在上游水壓力的緩慢作用下，水平位移分布規律由對稱分布逐漸向下游變形，豎向位移呈現均勻下降趨勢，最大值在壩體2/3處，運行期最大豎向位移略大于竣工期最大豎向位移.竣工期和蓄水期壩體豎向位移分別約占最大拔高的1.01%和1.03%，分析原因可能是當地黃土為濕陷性黃土所致.

4 敏感性分析

4.1 Plackett-Burman試驗設計方法

對φ0、Δφ、Rf、K、Kb、n、m參數進行 Plackett-Burman試驗設計，目的指標分別為豎向位移V、向下游水平位移H1和向上游水平位移H2.7個參數作為7個因素，每個因素考察低(-1)和高(+1)兩個水平，其中低水平為設計參數的80%，高水平為設計參數的120%，見表3.

表3 Plackett-Burman試驗因素水平

為保證篩選試驗的準確性，實驗設置一列誤差列，同時保證足夠的試驗次數便于分析出7個參數的顯著性程度，因此設計出一個需試驗24次的Plackett-Burman試驗方案，具體見表4.

表4 Plackett-Burman試驗方案表

表格中試驗順序是隨機的，目的是減少試驗中由于先后不均勻帶來的誤差干擾，進一步確保篩選的準確性.將φ0、Δφ、Rf、K、Kb、n、m分別填入試驗方案表的第2、3、4、5、6、7和8列;將豎向位移、向下游水平位移和向上游水平位移分別填入試驗方案表的第10、11和12列，表中第9列作為誤差列.

結果分析見表5.可以看出，各參數關于豎向位移V的敏感性程度由高到低依次為Rf＞Kb＞φ0＞K＞n＞Δφ＞m，其中對于豎向位移明顯顯著的是Rf、Kb和φ0，有一定影響的是參數K，而n、Δφ和m3個參數對于豎向位移均為無影響.各參數關于向下 游水平位移H1的敏感性程度由高到低依次為Rf＞φ0＞Kb＞K＞Δφ＞n＞m，其中對于向下游水平位移 明顯顯著的是Rf、φ0和Kb，而K、Δφ、n和m4個參 數對于向下游水平位移均為無影響.各參數關于向上 游水平位移H2的敏感性程度由高到低依次為Rf＞φ0＞Kb＞K＞n＞Δφ＞m，其中對于向上游水平位移 明顯顯著的是Rf、φ0、Kb，顯著的是K，而n、Δφ和m3個參數對于向上游水平位移均為無影響.

表5 Plackett-Burman試驗分析結果

綜合來看，Rf、φ0、Kb和K4個參數對于黃土均質壩變形的影響更加顯著，這與文獻[12]得到的結論一致，計算參數需要重點考慮的參數為Rf、φ0、Kb和K，而參數Δφ、n和m對于變形的影響基本可以忽略不計，因此在考慮一階交互的正交試驗中，只考察Rf、φ0、Kb和K的敏感性.

4.2 考慮交互作用的正交試驗

在上述Plackett-Burman試驗設計的基礎上，對Rf、φ0、Kb和K4個參數進行考慮交互作用的正交試驗.目的指標分別為豎向位移、向下游水平位移和向上游水平位移.Rf、φ0、Kb和K4個參數分別作為4個因素，每個參數考察低(-1)和高(+1)兩個水平，其中低水平為設計參數的80%，高水平為設計參 數的120%.具體因素和各因素水平見表6.

表6 正交試驗因素水平參數取值

由于不僅要考慮上述4個參數的敏感性，還需要考慮Rf×φ0、Rf×Kb、Rf×K、φ0×Kb、φ0×K以及Kb×K的敏感性.因此選擇能夠容納4因素及其交互作用的表頭設計.L16(215)表的交互作用見表7.

表7 L 16(215)表頭設計

選擇正交表時，由于需考慮4個參數，選擇4因子表頭，將Rf、φ0、Kb和K分別填入表頭中第1、2、4和8列，將Rf×φ0、Rf×Kb、φ0×Kb、Rf×K、φ0×K以及Kb×K分別填入表頭中3、5、6、9、10和12列，因此選用能容納12因素2水平的正交表，其中第7和11列可作為誤差列.

除了考慮因素水平外，還需要保證對精度的要求，試驗次數越多，試驗精度越高.故本次考慮交互作用的正交試驗選擇L16(215)的正交表.具體試驗方案 及有限元計算結果見表8.

表8 正交試驗方案及變形結果

對上表中數據及結果進行方差分析，結果見表9.關于豎向位移V的敏感性程度由高到低依次為Rf＞Kb＞φ0＞Rf×Kb＞φ0×Kb＞K＞Rf×φ0＞Rf×K＞Kb×K＞φ0×K，其中Rf、Kb和φ0明顯顯著，Rf×Kb、φ0×Kb、K、Rf×φ0和Rf×K顯著;關于向下游水平位移H1的敏感性程度由高到低依次為Rf＞φ0＞Rf×φ0＞Kb＞Rf×Kb＞φ0×Kb＞K＞Rf×K＞φ0×K＞Kb×K，其中Rf、φ0、Rf×φ0、Kb、Rf×Kb和φ0×Kb明顯顯著;關于向上游水平位移H2的敏感性程度由高到低依次為Rf＞φ0＞Rf×φ0＞Kb＞Rf×Kb＞φ0×Kb＞K＞Rf×K＞φ0×K＞Kb×K，其中Rf、φ0、Rf×φ0、Kb、Rf×Kb、φ0×Kb和K明顯顯著.綜合來看Rf、Kb和φ0是對大壩變形敏感性最強的3個參數，三者對于大壩變形的敏感性均強于參數K，考慮其交互作用時，Rf×φ0、Rf×Kb以及φ0×Kb之間的交互作用對于大壩變形的敏感性是不可忽視的.

表9 正交試驗方差分析

5 結 論

為研究鄧肯-張E-B模型參數對于大壩變形的影響，本文以某均質土壩作為研究對象，通過建立靜力分析模型，分析大壩運行期和竣工期的基本變形規律，并借助Plackett-Burman試驗設計和正交試驗設計對鄧肯-張E-B模型參數進行敏感性分析.得到如下主要結論:

1)由于壩址區主要存在黃土濕陷、深厚覆蓋層滲漏及不均勻變形問題，經處理后的壩址區具備建壩條件.土壩的變形分布規律滿足一般土石壩的分布特點.

2)通過Plackett-Burman試驗設計對7個鄧肯-張E-B模型參數進行研究，可以發現，從單個因素對于土壩變形的影響來看，7個計算參數中Rf、φ0、Kb和K4個參數是對于大壩變形明顯顯著的參數.

3)使用正交試驗設計對Rf、φ0、Kb和K這4個參數以及4個參數之間的一階交互作用進行全局敏感性分析，分析結果表明，Rf、Kb和φ0是對大壩變形敏感性最強的3個參數，三者對于大壩變形的敏感性均強于參數K，而Rf×φ0、Rf×Kb以及φ0×Kb之間的交互作用對于大壩變形的敏感性是不可忽視的，甚至其交互作用的敏感性大于單個參數的敏感性.

4)綜上可知，對均質黃土土壩鄧肯-張E-B模型參數的全局敏感性分析是十分必要的，本文研究內容可為土壩變形計算參數的選擇提供一定的參考，但是否可以推廣到其它級配及物理特性的土壩，需要進一步驗證.