云南某尾礦壩滲流三維計算分析

安曉宇，程 鵬，樂 陶

（1.天津大學，天津300072；2.交通運輸部天津水運工程科學研究所港口水工建筑技術國家工程實驗室，天津300456；3.中交第一公路勘察設計研究院有限公司，西安710075；4.中鋼集團馬鞍山礦山研究院有限公司，馬鞍山243000）

云南某尾礦壩滲流三維計算分析

安曉宇1，2，程 鵬3，樂 陶4

（1.天津大學，天津300072；2.交通運輸部天津水運工程科學研究所港口水工建筑技術國家工程實驗室，天津300456；3.中交第一公路勘察設計研究院有限公司，西安710075；4.中鋼集團馬鞍山礦山研究院有限公司，馬鞍山243000）

針對云南某尾礦庫蓄水后可能存在的繞壩滲漏問題，在充分分析壩址區工程地質和水文地質特征的基礎上，建立壩址區三維有限元計算模型，計算分析二期滿庫洪水和正常運行工況下地下水的滲流情況。計算結果表明，下游壩坡坡腳處地基滲流溢出處的水力梯度較大，建議在下游坡腳3～5 m范圍內地基土層表面布設厚度為0.8～1.0 cm的透水反濾碎石層；同時，由于壩基滲漏量較大，建議在壩體底部鋪設土工膜，減少兩岸繞滲，降低壩體的浸潤線。

尾礦庫；尾礦壩；滲流；三維有限元模型

尾礦庫是堆存尾礦和澄清尾礦水的專用場地，通常由圍護構筑物即尾礦壩（初期主、副壩及后期尾礦堆壩等）和排水構筑物（排水井、管、斜槽及隧洞等）組成。尾礦壩的安全不僅關系到礦山庫區自身的安全，而且還關系到下游及周邊人民群眾的生命財產安全。在尾礦壩建設之前，需要對壩體的滲流穩定性開展研究，分析壩體浸潤線的埋深狀態等特征，這樣可為壩基的抗滑和防滲設計提供參考。使用三維滲流分析的方法，可以更準確的得到庫區壩址處的滲流情況。

本文針對云南某尾礦壩工程，為了解其自身及其地基的防滲、抗滲性能，以及壩體的總滲漏量，采用三維有限元方法，對該尾礦壩及地基進行數值計算，從而定量分析評價壩體的總滲漏量和滲透穩定性，并對防止滲透變形提出滲控措施，對尾礦壩的設計、施工提供了參考依據。

1 工程概況

云南某尾礦庫位于昆明市西山區?？阪偸R哨村南西約1km處石馬哨溝中游，距450萬t/a選礦廠北東約7.0 km，庫區中心距?？阪偧s7 km，距昆明市35 km。石馬哨溝是螳螂川中游右側的一條支溝，壩址上游有兩條分支溝，主溝總長7.0 km，溝谷縱坡約為i=0.024，流向由西南至東北，壩址上游匯水面積4.0 km2；尾礦壩址距上游分水嶺2.5 km，壩址距出口螳螂川約4.5 km，距昆晉鐵路4.0 km。尾礦壩位于庫區下游溝口處，該庫在壩頂標高+2 055 m位置時，總庫容1 762.02萬m3，有效庫容1 497.72萬m3。

1.1 地質條件

根據勘察報告，壩址處的各巖層滲透性較大。擬建的尾礦庫庫區內的地層主要有：植物層（Qpd）、第四系沖洪積層（Qel+pl）、第四系殘坡積層（Qel+al）、上震旦統東龍潭組（Z2dl）及上震旦統陡山沱組（Z2d）巖土層。場地內分布的地層從上至下主要為：植物層（Qpd）、粉質粘土（Qel+pl）、粉砂（Qel+pl）、粘土（Qel+pl）、粘土（Qel+dl）、全風化泥巖、強風化藻白云巖、強風化石英砂巖、中風化藻白云巖、中風化石英砂巖、強風化白云巖、中風化白云巖、溶洞土等。尾礦壩典型地質橫剖面設計圖如圖1所示。1.2 壩體設計方案

壩體分兩期堆筑，二期壩體采用下游式堆筑法，從排水棱體頂部用碎石土沿一期壩貼坡加高，最終壩高80 m，壩頂寬度6.0 m，壩頂標高為+2 055 m，在壩體外坡+1 995 m、+2 010 m、+2 025 m、+2 040 m標高位置設置4條寬2.5 m的馬道，內坡坡比為1：2.25，內坡透水層沿一期透水層往上堆筑，外坡坡比為1：2.5，河床處清基面高程1 970 m，壩頂長度約664 m。尾礦壩縱剖面設計圖見圖2所示。

圖1 尾礦壩沿壩軸線地質剖面圖Fig.1 Geological profile along dam axis of tailing dam

2 滲流有限元計算基本原理

假定滲流從圍巖及尾礦入滲表面到出逸斷面的整個滲流區域內為連續水流，且滲流水流充滿整個壩體、壩基滲流孔隙區的連續介質，采用部分概化的方法，整個滲流區域內各不同滲透性的巖體盡可能按實際地質分布情況模擬［1-2］。各部分壩體材料視為滲流各向同性體；各部分壩基均視為滲流非各向同性體，水平向與垂直向滲透系數有一定差異（一般水平向滲透系數約為垂直向的3～4倍）。根據壩基地層的滲透性特征，建立計算區域內非均質各向異性滲流模型。

對于研究區內的某一單元，在該單元中的滲流滿足下列控制方程［3-4］

式中：θw為體積含水量；kwx、kwy、kwz為土體三個方向的滲透系數，隨含水量而變化；H為總水頭，是孔隙水壓力與位置水頭之和（H=uw+ρwgz）；ρw、g分別為水的密度與重力加速度，當忽略水的可壓縮性時，ρw為常量。

當巖土體不飽和時，體積含水量θw隨孔隙壓力變化的函數關系可通過土水特征曲線（SWCC）來確定。Van Genuchten提出了下列表達式（VG模型）

式中：a、n為試驗擬合參數；Se為有效飽和度；θs、θr分別為土體的飽和與殘余體積含水量；uw為孔隙水壓力，為負值。

當巖土體不飽和時，其滲透系數隨含水率（飽和度）變化，可由式（3）計算

式中：Kr為相對透水率，KS為飽和滲透系數。

VG模型中相對透水率Kr與有效飽和度Se的關系如式（4）所示

對應控制方程式（1）的定解條件可以是

圖2 尾礦壩縱剖面設計圖Fig.2 Longitudinal design of tailing dam

式中：H(x,y,z)表示區域內任一點(x,y,z)總水頭；φ0(x,y,z)和q0(x,y,z)為已知函數；n為邊界Γ2外法線方向。

由控制方程（1）和輔助方程（2）～（4）以及定解條件（5）、（6），共同構成描述飽和-非飽和巖土體穩定-非穩定滲流場的數學模型。該模型可統一考慮飽和-非飽和區，計算得到孔隙水壓力為零的點對應在浸潤面上。該模型既可用于計非穩定滲流場，也可用于計算穩定滲流場：當式（1）的左端項為0時，該方程退化描述穩定滲流控制方程。

3 計算模型與參數

3.1 計算程序簡介

滲流計算采用河海大學巖土工程研究所研制的US3D三維有限元程序。該程序有如下主要特點：程序既可用于飽和巖土體滲流計算也可用于非飽和巖土體滲流計算；既適用于穩定滲流計算也適用于非穩定滲流計算；無需變動網格即可計算得到浸潤線和逸出點；可得到計算區域內的滲流水頭分布、水力梯度和滲流力；可模擬土體飽和度變化和滲透系數的非線性；可根據計算要求計算通過任意給定截面的滲流量。

3.2 三維有限元計算模型

計算中X軸為沿壩軸線方向從左岸指向右岸方向，零點定在左壩肩處；Y軸為垂直壩軸線從上游指向下游方向，零點在壩軸線處；Z軸豎直向上，零點與高程零點相同。模型計算區域確定為X軸方向從左右壩肩各向外側60 m范圍；Y軸方向從壩軸線向上下游分別各取約2倍壩底寬；Z軸方向壩底清基面以下約2倍壩高，即計算模型的長寬高尺寸為：1 200×784×255 m。模型沿直壩軸線劃分41個剖面，剖面間距約為20 m，局部區域進行加密。為了能盡可能模擬地層的變化情況，計算模型共設置16個材料分區。模型中主要采用8結點六面體等參單元，局部用6結點五面體或4結點四面體單元過渡，共劃分23 554個單元，23 305個結點。整體三維有限元網格見圖3。

3.3 滲透系數的確定

圖3 三維有限元網格整體效果圖Fig.3 Three?dimensional numerical calculation model and meshing diagram

某尾礦壩滲控分析中，必須采用真實反映實際地質條件和壩體填筑情況的壩體壩基各介質的滲透系數。在參照云南化工院和武漢勘察研究院出具的《某尾礦庫巖土工程勘察報告》和長沙冶金設計研究總院的《某尾礦庫壩體穩定性研究——尾礦壩固結滲流分析模型研究報告》中尾礦滲透系數數量級基礎上，結合類似工程的基本特點，擬定了頭石山尾礦壩體壩基的滲流分析參數，匯總于表1中。

表1 某尾礦壩滲透系數匯總表Tab.1 Permeability coefficient of tailing dam

3.4 容許水力梯度的確定

庫區地表出逸處第四系土體的臨界水力梯度可根據《土工原理與計算》［5-6］中所給出下式進行計算

式中：n為土體孔隙率；Gs為土粒比重；ξ為側壓力系數；φ為土的摩擦角；c為單位土體所受的凝聚力；γw為水的容重。

根據地勘報告資料，地表土體n孔隙率為0.4～0.5，土粒比重為2.72～2.74，側壓力系數取0.5，土體摩擦角為19～25°，凝聚力為30～40 kPa。從而可推算得地表土體臨界水力梯度為11～16。若取容許安全系數為2～3，則地表出逸處土體容許水力梯度為3.7～8。

3.5 計算方案

經過討論，確定共對8種方案進行了計算。計算結果表明，由于兩岸基巖的滲透性很大，如果岸坡段壩體底部不設置土工膜，則壩體內浸潤線較高，對壩體穩定不利。因此，建議岸坡段壩體底部也設置土工膜。所以，將岸坡段設置土工膜作為基本方案進行計算分析。由于篇幅有限，在本文中列出其中3種方案的計算結果（表2）。

表2 某尾礦壩滲流計算分析方案Tab.2 Seepage calculation scheme of tailing dam

4 滲流計算結果及分析

4.1 二期滿庫洪水運行工況計算結果

圖4給出了二期滿庫洪水運行工況計算結果（A1方案）河谷處最大橫剖面（壩0+320）總水頭等值線分布情況?？梢钥闯?，受土工膜隔滲作用［7-8］，浸潤線在壩體迎水面處轉折，浸潤線從上游棱體與土工膜間沿水平排水層達到壩體下游側。表明受壩體迎水面及壩底土工膜包裹隔滲，加上壩體底部水平排水層及下游排水棱體作用，壩體內浸潤線保持在壩基附近。在土工膜和下游排水棱體附近，總水頭等值線分布密集，表明受土工膜隔滲和排水棱體排水作用，相鄰區域水頭損失較大。

圖5為A1方案河谷最大橫剖面（壩0+320）水力梯度等值線圖。從壩體內各點水力梯度分布情況來看，水力梯度最大值發生在河谷處下游坡腳排水棱體附近的地基土中，最大水力梯度達1.12。盡管水力梯度大于0.8的單元不多，但還是應引起注意。

由于河谷段2 020 m以下的下游坡腳處的水力梯度較大，建議在下游坡腳3～5 m范圍內地基土層表面布設透水反濾碎石層，厚度0.8～1.0 cm。排水棱體中的水力梯度較小，最大值近0.31。尾礦砂中的水力梯度大都小于0.3。水力梯度均不大，不會發生滲透破壞。此方案的滲漏量為95 500 m3/d。滲漏量較大，主要是兩邊岸坡滲透系數較大的緣故。

4.2 二期滿庫正常運行工況計算結果

圖6給出了二期滿庫正常運行工況計算結果（A2方案）河谷處最大橫剖面（壩0+320）總水頭等值線分布情況。

與A1方案計算結果（圖4）比較，壩體壩基內滲流規律基本相同。A1方案的上下游最大水頭差約84 m，A2方案的水頭比A1方案水頭減小幅度為1.8%。

圖4 A1方案河谷處最大橫剖面總水頭等值線（單位：m）Fig.4 Total head contour map of the valley maximum horizontal section in A1 scheme

圖5 A1方案河谷最大橫剖面水力梯度等值線Fig.5 Hydraulic gradient contour map of the valley maximum horizontal section in A1 scheme

圖6 A2方案河谷處最大橫剖面總水頭等值線（單位：m）Fig.6 Total head contour map of the valley maximum horizontal section in A2 scheme

圖7 A2方案河谷最大橫剖面水力梯度等值線Fig.7 Hydraulic gradient contour map of the valley maximum horizontal section in A2 scheme

圖8 A2方案河谷處最大橫剖面總水頭等值線（單位：m）Fig.8 Total head contour map of the valley maximum horizontal section in A5 scheme

圖9 A5方案河谷最大橫剖面水力梯度等值線Fig.9 Hydraulic gradient contour map of the valley maximum horizontal section in A5 scheme

圖7為A2方案河谷最大橫剖面（壩0+320）水力梯度等值線圖。河谷斷面下游出口處的水力梯度相比A1方案減小為1.07，幅度約4.5%。和A1方案相比，A2方案的壩體與岸坡內水力梯度略有減小。

A2方案為正常運行工況，其上游水位比A1方案洪水運行工況的水位低1.5 m。計算結果顯示：當上游水位略有下降時，壩體滲漏量略有下降；從A1方案的95 500 m3/d減小到92 400 m3/d，減小幅度3.2%。

4.3 二期滿庫洪水+降雨運行工況計算結果

圖8給出了二期滿庫洪水+降雨運行工況計算結果（A5方案）河谷處最大橫剖面（壩0+320）總水頭等值線分布情況。圖9為A5方案河谷最大橫剖面（壩0+320）水力梯度等值線圖。與二期滿庫洪水運行工況（A1方案）計算結果（圖4）比較，壩體壩基內滲流規律基本相近。水力梯度比A1方案的計算結果略大，其中，河谷段地基土層中的最大水力梯度增大8.0%。

比較A5和A1方案，結果顯示，總體上降雨作用對壩體滲流場影響不大，壩坡及岸坡水力梯度均不大，大部分區域小于0.5。

5 結論

根據前述計算結果，得出如下結論與建議：（1）采用擬定的基本計算參數和各種工況條件下，某尾礦庫滲流場符合一般規律。壩體和地基中大部分區域水力梯度小于0.5，遠小于其容許水力梯度；但在河谷段高程2 020 m以下的下游壩坡坡腳處地基滲流溢出處的水力梯度較大，建議在下游坡腳3～5 m范圍內地基土層表面布設厚度為0.8～1.0 cm的透水反濾碎石層。（2）尾礦庫總滲漏量主要取決于地基（風化石英砂巖及白云巖）的滲透性大小，壩基滲漏量較大，為了減小總滲漏量，需在壩體底部鋪設土工膜，這樣可以有效的減少兩岸的繞滲，降低壩體的浸潤線。

［1］Van Genuchren.A closed?form equation for predicting the hydraulic conductivity of unsaturated soil［J］.Soil Sci.Am.J.，1980，44 （5）：892-898.

［2］齊清蘭，張力霆.尾礦庫滲流場的數值模擬與工程應用［M］.北京：中國水利水電出版社，2011.

［3］胡云進.裂隙巖體非飽和滲流分析及其工程應用［M］.杭州：浙江大學出版社，2009.

［4］駱祖江，張弘，李會中，等.烏東德水電站壩址區地下水滲流三維非穩定流數值模擬［J］.巖石力學與工程學報，2011，30（2）：341-347. LUO Z J，ZHANG H，LI H Z，et al.Three?dimensional Transient Flow Numerical Simulation of Ground Water Seepage in Dam Site of Wudongde Hydropower Station［J］.Chinese Journal of Rock Mechanics and Engineering，2011，30（2）：341-347.

［5］殷宗澤.土工原理［M］.北京：中國水利水電出版社，2007.

［6］錢家歡，殷宗澤.土工原理與計算：二版［M］.北京：中國水利水電出版社，1996.

［7］SL/T231-98，聚乙烯（PE）土工膜防滲工程技術規范［S］.

［8］GB/T17643-2011，國家標準：土工合成材料聚乙烯土工膜［S］.

Analysis of 3D seepage of a tailing reservoir in Yunnan

AN Xiao?yu1,2,CHENG Peng3，YUE Tao4
（1.Tianjin University,Tianjin 300072,China;2.Tianjin Research Institute for Water Transport Engineering,National Engineering Laboratory for Port Hydraulic Construction Technology,Tianjin 300456,China;3.CCCC First Highway Consultants Co.,LTD.,Xi′an 710075,China;4.Sinosteel Maanshan Institute of Mining Research Co.,LTD., Maanshan 243000,China）

To investigate the leakage problem below a tailing dam in Yunnan province,a three dimensional fi?nite element model was developed according to the hydrogeological and engineering geological characteristics.The effects of groundwater seepage under full?reservoir and normal condition were studied and analyzed respectively. The calculated results show that the foundation′s seepage overflow point at downstream slope toe presents higher hy?draulic gradient.Thus the pavement of permeable filter gravel layer(0.8-1.0 cm)above 3-5m range of foundation soil at downstream slope toe was proposed.It is suggested that the geotechnical membrane should be laid at the bot?tom of the dam to reduce the seepage around both sides and lower the soaking line of the dam.

tailing reservoir;tailing dam;seepage;three dimensional finite element model

TU 46;O 242.21

A

1005-8443（2016）06-0641-05

2016-06-01；

：2016-07-11

安曉宇（1988-），男，天津市人，工程師，主要從事環境巖土工程研究。

Biography：AN Xiao?yu（1988-），male，engineer.