混凝土結構中氯離子擴散分析的精細積分法

楊綠峰,洪 斌,高 欽,曾建聰

(廣西大學工程防災與結構安全教育部重點實驗室,廣西南寧 530004)

鋼筋混凝土結構較長時間暴露于海洋等富含氯離子環境后,氯離子會侵入混凝土中并逐漸到達鋼筋表面,造成鋼筋銹蝕,從而影響鋼筋混凝土結構的使用性能和耐久性。由于鋼筋銹蝕的過程并不持續消耗氯離子,所以分析評估氯離子環境下鋼筋混凝土結構耐久性和使用壽命的關鍵是確定氯離子在混凝土中的分布和擴散規律。吳慶令等[1]通過現場海洋暴露試驗,研究了混凝土中氯離子的擴散特性。Collepardi等[2]用Fick第二擴散定律描述氯離子在混凝土中的擴散過程,得到在一定初始條件和邊界條件下的數學解;Liang等[3]考慮了水泥水化物對氯離子擴散的影響,并對氯離子擴散系數進行了修正;宋子健等[4]研究了溶液成分對混凝土中氯離子擴散系數的影響;彭國軍等[5]考慮骨料形狀,對混凝土氯離子擴散系數進行了數值預測。余紅發等[6]提出了考慮多因素擴散模型的一維解析解。楊綠峰等[7]研究了混凝土時變條件下氯離子擴散的封閉解。

由于實際的鋼筋混凝土結構通常具有復雜的幾何形狀和邊界條件,解析方法在應用中有較多困難,人們通常采用數值方法解決此類問題。在Fick第二定律及其數學模型的基礎上,施養杭等[8]提出了氯離子侵入混凝土計算的有限差分法模型;Funahashi[9]也在氯離子擴散分析的有限差分法方面進行了研究;Sergi等[10]應用最小平方法(leastsquares methods)對氯離子在混凝土中的擴散規律進行了一維數值模擬;Han[11]采用有限元法分析氯離子擴散問題,并推導出具體計算格式;楊綠峰等[12]研究建立了氯離子擴散分析的邊界元法。但是上述氯離子擴散分析的有限元法普遍存在兩個問題,其一是忽視了基于半無限大擴散場的Fick第二定理及其邊界條件對數值分析模型的影響,并造成計算誤差;其二是在用有限元求解氯離子擴散微分方程時,通常采用顯式的中央差分法或隱式的Newmark法、Wilson-θ法等,但計算結果的穩定性和計算精度時常不盡如人意。時域精細積分方法[13]可以將時域積分轉換為矩陣相乘,能夠取得非常高的計算精度;文獻[14]利用三次樣條插值函數模擬積分項中的被積函數,建立了求解非齊次動力方程特解的一種精細數值積分法。

由于氯離子侵入混凝土并發生擴散的過程非常緩慢,直接利用精細數值積分格式有時會帶來明顯的誤差。本文利用誤差函數建立了混凝土中氯離子擴散場的補償長度和補償系數,保證了有限元法在該類問題中的計算精度和計算效率;同時,在一致分布矩陣的基礎上建立了集中分布矩陣,驗證了集中分布矩陣具有更好的計算精度。在此基礎上利用精細積分技術研究建立了混凝土中氯離子擴散分析的精細積分有限元法,克服了普通有限元法分析氯離子擴散時存在的問題。

1 氯離子一維擴散的有限元方程

氯離子在混凝土中的擴散過程可以用Fick第二定律來描述:

其邊界條件和初始條件分別為

式中:D為混凝土中氯離子擴散系數;t為混凝土持續暴露于氯離子環境中的時間;C=C(x,t)表示t時刻混凝土試件內深度為x處氯離子濃度;C0為初始氯離子濃度;Cs為混凝土構件表面氯離子濃度。根據變分原理可知,上述控制方程是泛函的變分極值條件。

將擴散場沿 x方向離散為N個單元,可以建立典型單元(圖1)上t時刻的濃度分布函數:式中:N為形函數矩陣;C為單元節點濃度參數列陣 ,且有 :

圖1 擴散單元

C=(C1,C2)TN=(N1,N2)

式中:C1,C2分別為結點1和2上的氯離子濃度;Ni為插值函數 ,其中 N1=1-,N2=(=x/le為無量綱化單元局部坐標,le為單元長度)。

將式(4)代入式(3),可得單元泛函:式中:Ωe表示單元域;﹒C為C對時間t的一階導數;N′為N對坐標x的一階導數。

根據變分原理 δ Π=0可得混凝土中氯離子擴散分析的有限元方程:

其中

式中:M為氯離子分布矩陣;K為氯離子擴散矩陣。

當單元沿著與擴散方向相正交的另一個方向的尺寸取為單位值時,可求得式(7)為

2 一致分布矩陣和集中分布矩陣

計算式(7)中分布矩陣M所用到的形函數與計算擴散矩陣K所用到的形函數一樣時,稱 M為一致分布矩陣。由于求解方程(6)時采用一致分布矩陣M有時會導致計算結果不穩定,這里嘗試通過近似數值積分方法重新計算氯離子分布矩陣M。

考慮數值積分近似計算的梯形公式:式中:l為積分域的長度;k為積分點,取積分域的兩個端點;Fk(Ni,Nj)為被積函數F(Ni,Nj)在積分點k上的值,結合式(7)可以定義被積函數:

根據形函數的性質可以求得分布矩陣:

式(12)中只有對角元素不為零,且等于式(9)中相應行的全部元素的疊加。稱式(12)的分布矩陣為集中分布矩陣。

按照普通有限元法標準步驟,可以將單元分布矩陣和單元擴散矩陣集成為總體分布矩陣和總體擴散矩陣,相應地可將單元擴散方程集成為總體擴散方程。以下第3節、第4節的計算格式都是針對混凝土中氯離子擴散分析的總體有限元方程。

3 精細積分法

由于集中分布矩陣 M存在逆矩陣,所以由式(6)集成的總體有限元方程可以改寫為式中A為定常矩陣,且A=-M-1K。

按照微分方程求解理論,齊次方程(6)的通解為

式中eAt為矩陣At的指數函數。當將時間域[t0,ta]離散為p個等步長的子域,則時間步長Δt和時間域內的第i個離散結點分別為

與之相對應,在ti時刻,結點氯離子濃度向量表示為 Ci,且有:

式中指數矩陣 T=eAΔt。

根據矩陣加法定理有:

其中m為正整數,可選擇m=2N,N為正整數,因此m通常是非常大的正整數。由于Δt是有限大的時間步長,則τ=Δt/m通常非常小。由此可得:

式中I為單元矩陣。

將式(18)代入式(17),有:

將上述分解計算過程重復進行至第 n次時,有:

其中

由此可以看出,當 T0確定后,可以按照迭代計算公式(21)求出 Tn,當 n=N時,代入式(20)求得指數矩陣T,從而避免在計算過程中發生大數吃小數的情況;進而根據迭代計算公式(16),可以求解時間域上各個離散結點的氯離子濃度參數列陣Ci。

以上迭代求解方程的過程將時程積分問題轉化為矩陣相乘問題,簡化了計算,稱之為混凝土中氯離子擴散分析的精細積分有限元法。

4 混凝土中氯離子擴散場的補償長度

式(1)的解析解為

式中erf(s)為誤差函數,且:

由于Fick第二定理所導出的解析模型,要求遠端邊界取在無窮遠處,并滿足遠端邊界條件C(x=+∞)=C0,如式(2)所示。但有限元、邊界元等數值方法所建立的計算模型通常要求求解域為有限大,遠端邊界不可能取在無窮遠處,因此對于氯離子一維擴散分析的有限元、邊界元等數值方法[12],需要定義氯離子擴散場的補償長度L,如圖2所示,以此確定遠端邊界的位置,并保證遠端邊界條件C(x=L)=C0仍然能夠得到滿足。

圖2 混凝土試件長度及補償長度

將式(23)的誤差函數用圖表示,見圖3。假定遠端邊界位于x=y=L處,如圖2中虛線所示,它能夠保證擴散方程的遠端邊界條件成立,由式(24)

圖3 誤差函數y=erf(s)

可得

根據式(25)可得

式中,k=2erf-1(1),即有:

式(23)定義的誤差函數如圖3所示,可以看出,當1.5≤s≤2時誤差函數接近于1。因此,通過比較式(23)和式(27),容易看出式(27)成立時必然要求3≤k≤4。這里,k稱為混凝土的氯離子擴散長度補償系數,簡稱為補償系數。通常情況下建議取k=3.5。

根據式(26)可以確定計算擴散場的補償長度L,從而能夠保證遠端邊界條件的成立。補償長度L表示了混凝土的氯離子擴散問題中,正確計算混凝土內部各點的氯離子濃度所需的擴散場的最小長度,如圖2所示。圖2實線 ABCD表示試件的實際邊界,如果試件的實際長度 l≥L,則按照試件的實際長度建立有限元計算模型。反之,L≥l,那么在建立有限元法等數值方法的計算模型時,需要將圖2所示的混凝土試件的遠端邊界BC虛擬地移動至虛線B′C′處,從而在保證遠端邊界條件仍能成立的同時,方便有限元等數值方法計算模型的建立。

5 算例分析

算例1混凝土試件的尺寸為200mm,混凝土表面的氯離子濃度(氯離子與混凝土的質量百分比,下同)Cs=1.0%,初始氯離子濃度C0=0,氯離子擴散系數D=31.536mm2/a。分別采用基于集中分布矩陣的精細積分有限元法(LCM)和基于協調分布矩陣的精細積分有限元法(CCM)計算擴散時間分別為50a和80a時試件不同深度的氯離子濃度值,并將計算結果和封閉解(CFS)相比較。圖4給出了擴散時間50a時的LCM,CCM和CFS的結果。

圖4 擴散50a時氯離子濃度分布

從圖4可以看出,LCM和CCM的結果都與CFS比較接近,能夠給出正確的氯離子濃度分布規律,其中LCM的結果更為準確。為了便于比較LCM和CCM,這里通過和CFS比較計算兩種方法的誤差,見圖5,可以看出,LCM的計算結果精度優于CCM。

圖5 LCM和CCM的相對誤差比較

算例2制作一批混凝土試件,尺寸為150mm×150mm×200mm,養護28d后將5個面用環氧樹脂封閉,只留下1個150mm×150mm的面暴露于海水中,混凝土的氯離子擴散系數D=9.252×10-12m2/s=291.77mm2/a,試件中氯離子初始濃度C0=0,試件表面氯離子濃度Cs=0.565%,利用精細積分有限元法(PIFEM)計算暴露時間分別為50a和80a時試件內部的氯離子分布濃度,并與封閉解進行比較。在精細積分有限元法的計算模型中,取時間步長Δt=0.1a。

首先按照試件實際長度建立精細積分有限元法離散模型,計算試件暴露50 a和80a時的氯離子濃度分布,并將計算結果與CFS相比較,如圖6(a)所示,可以看出此時PIFEM計算結果與CFS之間的誤差相當大。

圖6 算例2計算結果

按照本文提出的擴散場補償長度公式(26)確定有限元法離散模型(圖2),并重新利用PIFEM計算。當暴露時間分別取50a和80a時,根據式(26),可以求得擴散場補償長度分別為544mm和688mm,二者都大于試件的實際長度,所以應根據擴散場補償長度建立有限元離散模型。將計算結果同CFS相比較,如圖6(b)所示,可以看出,按照本文提出的擴散場補償長度理論建立有限元法計算模型,計算結果和封閉解吻合很好,能夠取得很高的計算精度。

6 結 論

a.基于混凝土中氯離子擴散場的補償長度及其表達式,建立了有限元法離散模型,克服了傳統數值方法離散模型中半無限大假設導致的計算誤差,從而保證計算結果的正確性。

b.將一致分布矩陣轉換為集中分布矩陣,并與精細積分技術相結合,研究了混凝土中氯離子擴散分析的有限元法控制方程的求解方法,不僅保證了計算結果的穩定性,而且能夠取得很好的計算精度。

[1]吳慶令,余紅發,王甲春,等.現場海洋區域環境中混凝土的Cl-擴散特性[J].河海大學學報:自然科學版,2009,37(4):410-414.

[2]COLLEPARDI M,MARCIALIS A,TURRIZZANIR.The kinetics of penetration of chloride ions into the concrete[J].Cemento,1970(4):157-164.

[3]LIANG M T,WANG K L,LIANG C H.Service life prediction of reinforced concrete structures[J].Cement and Concrete Research,1999,29:1411-1418.

[4]宋子健,蔣林華.溶液成分對混凝土中氯離子擴散系數影響的研究綜述[J].水利水電科技進展,2011,31(5):82-87.

[5]彭國軍,鄭建軍,周穎瓊.考慮骨料形狀時混凝土氯離子擴散系數預測的數值方法[J].水利水電科技進展,2009,29(6):13-16.

[6]余紅發,孫偉,鄢良慧,等.混凝土使用壽命預測方法的研究[J].硅酸鹽學報,2002,30(6):686-701.

[7]楊綠峰,李冉.混凝土中氯離子二維擴散規律的解析研究[J].水利水電科技進展,2009,29(3):20-23.

[8]施養杭,羅剛.含多種因素的氯離子侵入混凝土的有限差分計算模型[J].工業建筑,2004,34(5):7-10.

[9]FUNANSHI M.Prediction corrosion-free service life of a concrete structure in a chloride environment[J].ACI Materials Journal,1990,87(6):581-587.

[10]SERGI H,PAGE C.Diffusion of chloride and hydroxyl ions in cementitious materials exposed to a saline environment[J].Magazine of Concrete Research,1992,44:63-69.

[11]HAN S H.Influence of diffusion coefficient on chloride ion penetration of concrete structure[J].Construction and Building Materials,2007,21:370-378.

[12]楊綠峰,陳正,王炎炎,等.混凝土中氯離子二維擴散分析的邊界元法[J].硅酸鹽學報,2009(7):1110-1117.

[13]ZHONG W X,WILLIAMS F W.A precise time step integration method[J].Journal of Mechanical Science,1994,208:427-430.

[14]富明慧,劉祚秋,林敬華.一種廣義精細積分法[J].力學學報,2007,39(5):672-677.