秦愛芳， 趙 飛， 趙小龍

(上海大學(xué)土木工程系，上海200072)

緩沖材料參數(shù)對核廢料處置庫近場影響的二維有限元分析

秦愛芳， 趙 飛， 趙小龍

(上海大學(xué)土木工程系，上海200072)

以高放射性核廢料地質(zhì)處置的FEBEX原位試驗作為數(shù)值計算模型，利用有限元軟件Code-Bright，通過改變膨潤土初始滲透系數(shù)、初始吸力和初始進氣值，得到熱-水-力(thermo-hydro-mechanical，THM)耦合作用下處置庫關(guān)閉后緩沖層飽和度和吸力的變化規(guī)律，以及以上3個因素對這些性狀影響的敏感程度，研究結(jié)果可為核廢料處置庫緩沖材料的選取提供參考．

膨潤土;核廢料處置庫;熱-水-力耦合;數(shù)值模擬

Abstract:The FEBEX in-situ test for geological disposal of high-level radioactive nuclear waste is used as a calculation model．By changing the initial permeability，the initial suction and the air entry value，variation of saturation，and suction under the coupled thermo-hydro-mechanical(THM)action are analyzed using a Code-Bright program after the closure of nuclear waste repository．By analyzing sensitivity of the three parameters on these traits，valuable reference is presented for the selection of buffer materials of the nuclear waste repository．

Key words:bentonite;nuclear waste repository;coupled thermo-hydro-mechanical(THM);numerical simulation

高放射性核廢料是通過深埋于適當(dāng)?shù)膸r石、鹽或粘土層并使用人工屏障隔離核廢料來達到安全處置的目的，其中作為緩沖/回填材料的膨潤土由于具有較高的膨脹性、極低的滲透性和優(yōu)良的核素吸附等性能，常被用作存放核廢料的銅罐和周圍巖石之間的人工隔離屏障［1］．處置庫投入運營后，一方面庫內(nèi)核廢料衰變放熱導(dǎo)致處置庫周圍緩沖/回填材料產(chǎn)生溫度應(yīng)力;另一方面由于圍巖高壓力水頭的作用使周圍巖石內(nèi)水向膨潤土發(fā)生滲透，從而導(dǎo)致高壓實的、富含蒙脫石的膨潤土的飽和度持續(xù)上升，并產(chǎn)生膨脹應(yīng)力．這樣的隔離屏障在熱量和水的作用下會產(chǎn)生非常復(fù)雜的“熱-水-力(thermo-hydromechanical，THM)”現(xiàn)象，并且它們之間相互關(guān)聯(lián)耦合．

為了確保這個隔離屏障的長期安全性，有必要研究這種耦合現(xiàn)象的作用機理．目前，一些發(fā)達國家已經(jīng)對此進行了許多試驗和理論研究工作．以原位試驗為例，有瑞士的“FEBEX”原位試驗、美國的尤卡山試驗、瑞典的“原型倉庫”試驗、德國高勒本地下試驗室等項目．

1985年，中國開始針對高放核廢料處置展開了一系列實質(zhì)性的研究，如甘肅北山預(yù)選區(qū)花崗巖場的選址和場址評價研究［2-3］．2006年初，國家相關(guān)部門頒布《高放廢物地質(zhì)處置研究開發(fā)規(guī)劃指南》［4］．2007年，國家頒布的《核電中長期發(fā)展規(guī)劃》中明確提出了到2020年前建成處置高放核廢料的地下試驗室、到21世紀(jì)中葉建成高放廢物處置庫的目標(biāo)［5］，而其中所涉及的THM耦合作用機理的研究尚處于起步階段．為此，許多學(xué)者正在對此進行積極而有意義的探索．

劉月妙等［6］采用有限差分計算程序 FLAC2D，考慮熱-力耦合，模擬在不同處置室間距的條件下，高放廢物處置庫的溫度場、應(yīng)力場和位移場規(guī)律．張玉軍［7-8］研制出了用來分析飽和孔隙介質(zhì)中THM耦合彈塑性問題的二維有限元程序，并對FEBEX原位試驗進行了數(shù)值模擬，以及通過改變初始水頭、初始滲透系數(shù)和初始飽和度對緩沖層飽和過程的影響進行了二維有限元分析．Millard等［9］對圍巖中核廢料處置庫近場THM耦合進行了數(shù)值模擬分析，得出應(yīng)力對THM耦合影響很大．Lee［10］利用數(shù)值模擬和現(xiàn)場調(diào)查的方法研究了加拿大高放核廢料地質(zhì)處置室圍巖的熱-力學(xué)穩(wěn)定性．Abel等［11］采用Van Genuchten模型，通過試驗和數(shù)值模擬的方法研究了MX-80膨潤土干密度、溫度和含水量與吸力的變化規(guī)律等．

高放核廢料處置庫中，作為緩沖/回填材料的膨潤土的主要的功能是，遇水吸濕引起飽和度的上升，并導(dǎo)致土體發(fā)生膨脹變形(主要發(fā)生在膨潤土靠近周圍巖石部分)，以密封膨潤土在砌置時形成的塊體之間和塊體與周圍巖石之間的施工縫隙，以及周圍巖石因處置庫開挖卸載引起的裂縫，進而直接影響核素的遷移．影響膨潤土飽和度、吸力的諸多因素中，初始滲透系數(shù)、初始吸力及初始進氣值比較關(guān)鍵，因此，本研究選用適當(dāng)?shù)臒?、水、力本?gòu)方程及平衡方程，采用Code-Bright軟件［12］對處置庫近場進行THM數(shù)值模擬，得到處置庫近場膨潤土飽和度和吸力的變化規(guī)律，以及初始滲透系數(shù)、初始水壓及初始進氣值對其影響的敏感程度．

1 熱-水-力耦合模型

1．1 水質(zhì)量、能量及力學(xué)平衡方程［13-14］

本研究的基本控制方程為水質(zhì)量守恒、能量守恒以及靜力平衡方程．

1．1．1 水的質(zhì)量守恒方程

水的質(zhì)量守恒方程為

式中，θl和θg分別表示單位體積液相中水的質(zhì)量和單位體積氣相中水(水蒸氣)的質(zhì)量，Sl和Sg分別表示液相和氣相的飽和度，φ為空隙率，fw為單位體積介質(zhì)中外界提供的水的質(zhì)量，jwl和jwg分別表示液相和氣相中水流動的質(zhì)量．jwl和jwg分別由以下幾部分組成:

式中，ql和qg分別為液相和氣相的Darcy流動(相對于固相)，U·為固相的速度，iwg為氣相中水蒸氣的擴散．

1．1．2 能量守恒方程

能量守恒方程為

式中，Es，El和 Eg分別為固、液、氣三相的比內(nèi)能;ρs，ρl和 ρg分別為固、液、氣三相的密度，其中液相和氣相的能量Elρl和Egρg又可表示為

fE為單位體積介質(zhì)中外界提供的能量;ic為熱傳導(dǎo)能量;jEs，jEl和jEg分別表示三相相對于固定參照系的熱對流能量，具體由以下幾部分組成:

其中iag為氣相中空氣的擴散，θal和 θag分別表示單位體積液相中和氣相中空氣的質(zhì)量．

1．1．3 靜力平衡方程

忽略慣性力的靜力平衡方程為

式中，σ為應(yīng)力狀態(tài)向量;b為體力向量．對于飽和多孔介質(zhì)，力學(xué)本構(gòu)方程引入有效應(yīng)力概念，即

式中，mT=(1，1，1，0，0，0)．對于非飽和多孔介質(zhì)的本構(gòu)方程，通常引入兩個獨立的應(yīng)力狀態(tài)量，即凈法向應(yīng)力σ－pgm以及基質(zhì)吸力s=pg－pl，其中pg為氣壓力，pl為水壓力．

1．2 熱本構(gòu)方程

本研究假定熱傳導(dǎo)服從Fourier定律，即

式中，ic為介質(zhì)導(dǎo)熱量;T為介質(zhì)溫度的微分;λ為熱傳導(dǎo)系數(shù)，表征材料導(dǎo)熱性能的物性參數(shù)(λ越大，導(dǎo)熱性能越好)．顯然，λ和三相熱傳導(dǎo)系數(shù)λs，λl，λg，空隙率以及飽和度有關(guān)，具有如下關(guān)系式:

式中，Sl為土的飽和度，λdry，λsat分別為完全干燥和完全飽和狀態(tài)下土的熱傳導(dǎo)系數(shù)．

1．3 流體本構(gòu)方程

本研究假定孔隙中液相、氣相遵循Darcy定律，即

式中，ql，qg分別為液相和氣相的流量;Ki為液相或氣相滲透系數(shù)張量(Ki=kriK/μi，i為l或g)，其中kri為液相或氣相相對滲透系數(shù);g為重力向量;K為介質(zhì)固有滲透系數(shù);Δpl，Δpg分別為水壓力和氣壓力的梯度;ρl，ρg分別為液相和氣相的密度．

與以上計算有關(guān)的參數(shù)表達式如下．

(1)固有滲透系數(shù)K與孔隙率直接相關(guān)，固有滲透系數(shù)為

式中，φ0為參考孔隙比，K0為對應(yīng)于φ0的固有滲透系數(shù)．

(2)液相相對滲透系數(shù)為

式中，A為常數(shù);n為指數(shù)，通常取為2～4;Se為有效飽和度，Se與實際飽和度Sl之間的關(guān)系為

式中，Slr為液相殘余飽和度，Sls為液相最大飽和度．

(3)飽和度與基質(zhì)吸力的關(guān)系．

為了建立起液相飽和度與基質(zhì)吸力之間的聯(lián)系，反映非飽和土滯留水的能力，較常見的 Van Genuchten模型的表達式為

式中，P=P0σ/σ0，P0為參照溫度下的進氣值，σ0為參照溫度下的表面張力;P和 σ為計算溫度下的進氣值和表面張力;β為形狀參數(shù)．

1．4 力學(xué)本構(gòu)方程

目前有關(guān)非飽和土的應(yīng)力應(yīng)變本構(gòu)模型很多，但能較好反映非飽和土相關(guān)力學(xué)參數(shù)的模型有非線性彈性模型和巴塞羅那模型．本研究采用非線性彈性模型［14］，該模型體積應(yīng)變的計算采用以下關(guān)系式:

式中，p 為凈平均應(yīng)力，a1，a2，a3為無量綱參數(shù)，s=Pg－Pl為基質(zhì)吸力．剪切應(yīng)力應(yīng)變關(guān)系為線性，剪切模量為G．本模型假設(shè)基質(zhì)吸力的變化只引起體積應(yīng)變．

2 有限元數(shù)值分析

2．1 計算模型［15］

本算例采用FEBEX原位試驗?zāi)Ｐ停瓼EBEX原位試驗場位于瑞士的阿爾卑斯山脈，海拔1 725 m，其埋深約450 m．基巖為由花崗巖和片麻巖等組成的結(jié)晶質(zhì)巖類．該試驗處置坑道水平的直徑為2．4 m，長18．8 m(見圖1)．水平軸向放置兩個加熱器，其直徑為0．9 m，長為 4．54 m．加熱器周圍填充緩沖/回填材料膨潤土．加熱器最初提供恒定的熱量1 200 W，直到靠近加熱器的膨潤土左邊界達到所需的最高溫度100℃;然后控制加熱器中的熱量使膨潤土左邊界的最高溫度保持100℃恒定不變;3年后(1 095 d)斷掉熱源，6個月冷卻期(1 278 d)后拆除熱源．

初始條件設(shè)定為緩沖材料和坑道回填完畢，此時膨潤土和巖石的初始溫度均為12℃，初始各向應(yīng)力分別為28．0 和 0．5 MPa，初始水壓分別為 －97．2和0．5 MPa，初始?xì)鈮壕鶠?0．1 MPa．巖體是飽和介質(zhì)，膨潤土的初始飽和度是0．46，其中巖體、膨潤土的主要參數(shù)如表1所示．

圖1 FEBEX原位試驗近場模型圖(m)Fig．1 Layout of FEBEX in-situtest(m)

表1 主要計算參數(shù)Table 1 Main computation parameters

考慮試驗布局的對稱性和影響范圍，取處置庫軸對稱模型試驗，即從單個熱源中心向外依次為加熱器、緩沖層(膨潤土)和巖石(見圖2)．模型的邊界條件如下:左邊界在水平方向加約束;上、下邊界在垂直方向加約束;右邊界由于耦合作用影響的效果不明顯，所以徑向應(yīng)力、水壓和溫度取恒定值，分別為28 MPa，0．5 MPa和12℃．本研究重點對以下3個部位進行分析，從左向右依次是:加熱器附近膨潤土(節(jié)點238)、膨潤土中部(節(jié)點239)、靠近巖石的膨潤土(節(jié)點240)，其中共有180個四邊形等參數(shù)單元，273個節(jié)點．

2．2 初始滲透系數(shù)的影響

假定膨潤土的初始滲透系數(shù)分別為2×10－15，2 ×10－13和2 ×10－11m/s，初始進氣值為 50 MPa，初始飽和度為 0．46，初始吸力為 97．3 MPa，計算至1 278 d時加熱器周圍膨潤土中三點的飽和度-時間曲線和吸力-時間曲線如圖3和圖4所示．

圖3 不同初始滲透系數(shù)下飽和度-時間曲線Fig．3 Saturation degree-time curves for different initial permeabilities

圖4 不同初始滲透系數(shù)下吸力-時間曲線Fig．4 Suction-time curves for different initial permeabilities

由圖3可見，當(dāng)初始滲透系數(shù)較大時，3個部位的飽和度均隨時間變化較小;當(dāng)初始滲透系數(shù)較小時，3個部位的飽和度變化差異均較大．隨著初始滲透系數(shù)逐漸減小，膨潤土飽和度達到某一穩(wěn)定值所需的時間將延長．對于節(jié)點238的膨潤土，當(dāng)初始滲透系數(shù)較小時，其飽和度隨時間急劇減小至某一值后達到相對穩(wěn)定，并基本維持不變．?dāng)嚅_熱源后(即1 095 d后)進入冷卻期，飽和度略微回升;對于節(jié)點239的膨潤土，不同初始滲透系數(shù)的條件下，其飽和度隨時間變化不大，且隨初始滲透系數(shù)的減小，達到相對穩(wěn)定后膨潤土的飽和度略有降低;對于節(jié)點240的膨潤土，當(dāng)滲透系數(shù)較小時，飽和度隨時間逐漸增大至一定值后達到相對穩(wěn)定．?dāng)嚅_熱源后(即1 095 d后)進入冷卻期，飽和度略微回落．

由圖4可見，對于節(jié)點238的膨潤土，隨著初始滲透系數(shù)逐漸減小，吸力逐漸增大;對于節(jié)點239的膨潤土，其吸力隨時間先減小后增大至某一值后達到相對穩(wěn)定;對于節(jié)點240的膨潤土，其吸力隨時間逐漸減小并很快達到某一相對穩(wěn)定值．當(dāng)初始滲透系數(shù)較小時，對節(jié)點238的膨潤土吸力影響較大，而初始滲透系數(shù)的大小對其他兩個部位影響不是很明顯．

2．3 初始吸力的影響

假定膨潤土內(nèi)的初始吸力分別為117．3，97．3和 77．3 MPa，對應(yīng)的初始飽和度分別為 0．40，0．46和0．58，初始進氣值為50 MPa，初始滲透系數(shù)為2×10－13m/s，計算至1 278 d時加熱器周圍膨潤土中三點的飽和度-時間曲線如圖5所示．

由圖5可見，初始吸力不同則初始飽和度不同．對于節(jié)點238的膨潤土，初始吸力的大小對其達到穩(wěn)定后的飽和度影響不大．飽和度隨時間急劇銳減至某一值后達到相對穩(wěn)定，斷開熱源后(即1 095 d后)進入冷卻期，飽和度略微回升;對于節(jié)點239的膨潤土，初始吸力的大小對飽和度變化的影響是相近的，隨著初始吸力逐漸減小，相對穩(wěn)定后所能達到的飽和度逐漸增大;對于節(jié)點240的膨潤土，初始吸力的大小對飽和度變化的影響也是相近的，飽和度隨時間逐漸增大到一定值后達到相對穩(wěn)定．?dāng)嚅_熱源后(即1 095 d后)進入冷卻期，飽和度略微回落，而在相對穩(wěn)定后所能達到的飽和度隨初始吸力的減小而增大．

2．4 初始進氣值的影響

假定膨潤土的初始進氣值分別為20，50和80 MPa，初始吸力為97．3 MPa，初始滲透系數(shù)為2×10－13m/s，計算至1 278 d時加熱器周圍膨潤土中三點的飽和度-時間曲線如圖6所示．

圖5 不同初始吸力下飽和度-時間曲線Fig．5 Saturation degree-time curves for different initial suction

圖6 不同初始進氣值下膨潤土內(nèi)飽和度-時間曲線Fig．6 Saturation degree-time curves for different initial air entry value

由圖6可見，對于節(jié)點238的膨潤土，飽和度隨時間急劇銳減至一定值后達到相對穩(wěn)定．?dāng)嚅_熱源后(即1 095 d后)進入冷卻期，飽和度略微回升;對于節(jié)點239的膨潤土，其飽和度隨時間先增大后略微減小至某一值后達到相對穩(wěn)定;對于節(jié)點240的膨潤土，飽和度隨時間逐漸增大至一定值后達到相對穩(wěn)定，斷掉熱源后(即1 095 d后)進入冷卻期，飽和度略微回落．以上3個點的初始進氣值越小，初始飽和度及飽和度穩(wěn)定后的值也越?。?/p>

3 結(jié)論

(1)高放核廢料處置中膨潤土飽和度受初始滲透系數(shù)、初始吸力和初始進氣值影響較大．

(2)當(dāng)初始滲透系數(shù)較大時，3點的飽和度及吸力變化均不大．當(dāng)初始滲透系數(shù)較小時，靠近加熱器處的膨潤土，其飽和度及吸力受初始滲透系數(shù)影響最大，飽和度隨時間急劇減小至某一值后達到相對穩(wěn)定，進入冷卻期后飽和度略微回升，其吸力值卻逐漸增大;中間處的膨潤土，其飽和度隨時間略有降低隨后相對穩(wěn)定，其吸力隨時間先減小后增大至某一值后達到相對穩(wěn)定;靠近巖石處的膨潤土，其飽和度隨時間逐漸增大至一定值后達到相對穩(wěn)定，進入冷卻期后飽和度略微回落，其吸力隨時間逐漸減小并能很快達到某一相對穩(wěn)定值．

(3)初始吸力越小，初始飽和度越大．靠近熱源及巖石處膨潤土的初始吸力的大小，對飽和度變化影響較大．靠近熱源處的膨潤土，其飽和度隨時間急劇減小到某一值并保持相對穩(wěn)定，進入冷卻期后飽和度略微回升，初始吸力的大小對其達到穩(wěn)定后的飽和度沒有影響;中間處的膨潤土，其飽和度隨時間變化不大，初始吸力的大小對飽和度變化影響很小;靠近巖石處的膨潤土，其飽和度隨時間增大到某一值并保持相對穩(wěn)定，進入冷卻期后飽和度略微回落．隨初始吸力的逐漸減小，相對穩(wěn)定后所能達到的飽和度逐漸增大．

(4)以上3個點的初始進氣值越小，初始飽和度及飽和度穩(wěn)定后的值也越?。?/p>

［1］ 孫文靜，孫德安，閆威．側(cè)限應(yīng)力狀態(tài)下非飽和膨潤土的變形特性［J］．上海大學(xué)學(xué)報:自然科學(xué)版，2010，16(1):105-110．

［2］ 徐國慶．2000～2040年我國高放廢物深部地質(zhì)處置研究初探［J］．鈾礦地質(zhì)，2002，18(3):160-167．

［3］ 王駒，范顯華，徐國慶，等．中國高放廢物地質(zhì)處置十年進展［M］．北京:原子能出版社，2004:13-26．

［4］ 國防科學(xué)技術(shù)工業(yè)委員會，科學(xué)技術(shù)部，國家環(huán)境保護總局．高放廢物地質(zhì)處置研究開發(fā)規(guī)劃指南［Z］．2006．

［5］ 王駒．高放廢物深地質(zhì)處置:回顧與展望［J］．鈾礦地質(zhì)，2009，25(2):71-77．

［6］ 劉月妙，王駒，蔡美峰，等．熱-力耦合條件下高放廢物處置室間距研究［J］．鈾礦地質(zhì)，2009，25(6):373-378．

［7］ 張玉軍．熱-水-應(yīng)力耦合模型及FEBEX原位試驗二維有限元分析［J］．巖土工程學(xué)報，2007，29(3):313-318．

［8］ 張玉軍．三個因素對緩沖層飽和過程影響的二維有限元分析［J］．巖土力學(xué)，2006，27(6):853-857．

［9］ MILLARD A， REJEB A， CHIJIMATSU M， et al．Numerical study of the THM effects on the near-field safety of a hypothetical nuclear waste repository—BMT1 of the DECOVALEX Ⅲ project．Part 2:Effects of THM coupling in continuous and homogeneous rocks［J］．International Journal of Rock Mechanics＆ Mining Sciences，2005，42:731-744．

［10］ LEE C F．Thermo-mechanical stability of rock around a nuclear waste disposal vault［C］ ∥ International Conference on Computational Methods in Structural and Geotechnical Engineering．1994:82-90．

［11］ ABEL C J，MARíA V V，GóMEZ-ESPINA R，et al．Adaptation of the Van Genuchten expression to the effects of temperature and density for compacted bentonites［J］．Applied Clay Science，2009，42:575-582．

［12］ OLIVELLA S， GENS A， CARRERA J． Numerical formulation for a simulator(Code-Bright)for the coupled analysis of saline media［J］．Engineering Computations，1996，13:87-112．

［13］ OLIVELLA S， CARRERA J， GENS A， et al．Nonisothermal multiphase flow of brine and gas through saline media［J］．Transport in Porous Media，1994，15:271-293．

［14］ RUTQVIST J，BORGESSON L，CHIJIMATSU M，et al．Thermohydromechanics of partially saturated geological media;governing equations and formulation of four finite element models［J］．International Journal of Rock Mechanics＆ Mining Sciences，2001，38:105-127．

［15］ GENS A， GARCIA-MOLINA J， OLIVELLA S， et al．Analysis of a full scale in situ test simulating repository conditions［J］．International Journal for Numerical and Analytical Methods in Geomechanics，1998，22:515-548．

Two-Dimensional FEM Analysis of Near Field Influence of Buffer Material Parameters on High Level Radioactive Nuclear Waste Repository

QIN Ai-fang， ZHAO Fei， ZHAO Xiao-long
(Department of Civil Engineering，Shanghai University，Shanghai 200072，China)

TU 443

A

1007-2861(2012)05-0539-06

10．3969/j．issn．1007-2861．2012．05．018

2011-00-00

秦愛芳(1966～)，女，副教授，博士，研究方向為巖土力學(xué)與工程．E-mail:qinaifang@shu．edu．cn