混凝土板裂紋擴展的態型近場動力學模擬

劉一鳴+黃丹+秦洪遠

摘要： 構建考慮混凝土拉壓異性和宏觀斷裂特征的混凝土類材料非局部態型近場動力學本構模型，并通過引入動態松弛、系統失衡判斷和力邊界等效等算法，構建適于分析混凝土類材料和結構變形破壞過程的態型近場動力學數值模擬體系.通過分組模擬和定量計算，分析算法的收斂性、計算精度和效率等問題；在此基礎上開展含不同角度中心裂紋混凝土板的破壞模擬.

關鍵詞： 混凝土； 裂紋擴展； 態型近場動力學； 非局部模型

中圖分類號： O346.1 文獻標志碼： A

0 引 言

材料和結構的破壞機制及其數值模擬是計算力學研究的經典難題，也是諸多工程領域關注的重點.傳統的有限元等數值方法由于在連續介質框架下求解偏微分方程，所以在分析破壞問題時必須預先知道裂紋的存在與否及其位置和尺寸，計算時又需要判斷裂紋是否擴展及擴展路徑，并重新剖分網格，具有一定的復雜性.邊界元、擴展有限元、非連續有限元等方法及虛擬裂紋閉合等措施可以較好地處理不連續問題，但由于依然在連續性假設框架下進行不連續區域的特殊處理，模擬復雜的諸如三維群裂紋、多尺度動力破壞等問題時，依然有待進一步研究.

近年來，基于非局部積分思想的近場動力學（PeriDynamics，PD）方法[1-2]憑借其不需要求解空間微分方程而在模擬大變形及裂紋擴展、爆炸和沖擊破壞等強不連續力學問題方面的突出優勢[3-5]，成為計算力學和相關領域研究的熱點.針對建筑工程中最廣泛使用的混凝土材料和結構，SILLING[6]，DEMMIE[7]和KILIC等[8-9]采用近場動力學方法模擬簡單混凝土梁和柱的沖擊破壞、失穩等問題，GERSTLE等[10-12]在模擬混凝土梁的基礎上，還進一步分析近場動力學模型參數對計算精度的影響.本課題組也曾根據近場動力學思想建模分析常規混凝土構件的拉、壓和沖擊破壞過程.[13-15]

然而，已有的相關工作主要基于常規“鍵型”近場動力學模型，不論是采用早期的單參數微彈脆性模型[6-9]，還是應用改進后的微極鍵型近場動力學模型[11-12]，均無法真正滿足實際混凝土材料和結構破壞模擬特別是定量計算分析的需要（例如，固有的泊松比限制問題[1，8]，三維模擬時單參數微觀脆性模型中泊松比限定于0.25，而微極模型[11-12]的泊松比則必須小于0.25），對材料的泊松比具有一定的限制.為修正包括泊松比限制在內的鍵型近場動力學模型的缺陷，SILLING等[16-17]進一步發展態型近場動力學理論，近年來已成為關注的熱點[18-19]，但其相關工作主要圍繞各種傳統本構模型和態型近場動力學模型之間的轉化來開展，對于混凝土材料和結構破壞這一工程實際問題的態型近場動力學模擬尚未見諸報道.

本文基于態型近場動力學理論，構建考慮混凝土拉壓異性和宏觀斷裂特征的彈脆性態型近場動力學本構模型，并且通過引入動態松弛、系統失衡判斷、力邊界等效等算法，構建完整的統一求解混凝土材料和結構變形破壞問題的數值體系.通過對混凝土板單軸拉伸問題的定量計算和分組模擬，分析方法的收斂性、計算精度和效率.在此基礎上，應用本文的模型和方法模擬含不同角度中心裂紋混凝土板的裂紋擴展與破壞過程.

1 態型近場動力學理論和模型

1.1 態型近場動力學理論

將空間物質視為由帶質量、有代表性體積的系列物質點x組成，其

僅與其近場范圍H（|ξ|<δ）內的其他物質點相互作用，見圖1.圖1中：x和x′分別表示參考構形中2個物質點位置矢量；u和u′分別表示當前構形中2個物質點位移矢量；y=x+u和y′=x′+u′分別表示當前構形中2個物質點的坐標矢量；ξ=x′-x表示2個物質點的相對位置矢量；η=u′-u表示2個物質點的相對位移矢量.

1.2 二維常規態型近場動力學本構模型

2 數值實現

2.1 局部阻尼引入

為便于應用近場動力學方法求解準靜力問題，根據準靜力結果驗證模型的可靠型，引入經典力學中的動態松弛法，在運動方程中引入局部阻尼項[8，15]，將物質點運動方程轉化為

2.2 數值離散

對平面問題采用四邊形均勻離散（單位厚度），離散間距為Δx，則物質點的代表面積為VΔx=Δx2，對于物質點xi考慮近場范圍內（xj-xi≤δ）的相互作用

2.3 數值算法

近場動力學方法的顯式動力學實現涉及空間離散方式、時間差分格式及數值積分方法，數值積分方法與空間離散方式相關.本文采用中心差分對時間序列進行離散，

2.4 系統失衡力準則

3 數值算例

為驗證本文模型和算法，首先通過對混凝土板的單向拉伸變形進行定量計算分析，進而分別對含不同角度中心裂紋的混凝土板受拉伸時的裂紋擴展和破壞過程進行分析.根據實際情況，將問題簡化為平面應力問題分析.

對粒子系統的力邊界條件施加問題是粒子類方法中的難點.對于長為L，寬為W的矩形板（見圖2），本文在近場動力學模擬時在模型的加載邊界上增加寬度為δ的加載區域，將外力載荷等效為加載區物質點體力施加于加載區物質點上，使模型邊界上應力等效邊界條件一致.這樣做還可以消除往常近場動力學計算中常見的邊界弱化問題.

3.1 板拉伸定量分析

取200 mm×200 mm混凝土方板，密度ρ=2 400 kg/m3，彈性模量E=30 GPa，泊松比v=0.2，體積模量和剪切模量分別為K=16.67 GPa，G=12.5 GPa.在豎直方向施加均布拉伸載荷σ=3 MPa.模型離散情況見圖2，離散物質點間距Δx=1 mm，模型共包含40 600個物質點，其中邊界物質點600個，近場范圍取δ取3Δx.使用控制變量法，通過改變時間步長、人工阻尼因數C等計算參數，觀察迭代過程中板下邊緣中點y方向位移的收斂變化情況.為方便觀察位移隨時間步變化情況，將收斂容差-設為0，程序根據最大迭代步終止計算.

（1）時間步長取Δt=2×10-7s，改變阻尼系數C時的位移收斂情況見圖3.y方向的位移斂于1.011 8×10-5 mm，與彈性力學解1×10-5 mm相比，誤差為1.18%.

由此可見，當阻尼系數相對較小時，位移值在理論解附近劇烈振蕩且無法收斂；隨著阻尼系數增大，位移振蕩幅值降低，且最終收斂于理論解.當阻尼系數達到1×108 kg/（m3·s）時，位移值呈單調收斂，表現出類似超阻尼現象.阻尼系數繼續增大，位移收斂逐漸變緩.合適地選取阻尼系數可以大幅提高準靜態求解的計算效率.根據本文分組比較與分析，建議將局部阻尼項阻尼系數取為5×107～5×108 kg/（m3·s）.

（2）取阻尼系數C=2×108 kg/（m3·s）不變，改變時間步長Δt進行分析，見圖4.隨著時間步長的逐漸增大，位移的收斂也加快，其中時間步長為3×10-7s時收斂速度最快，但當時間步長取Δt=5×10-7s時，差分迭代出現不收斂.關于態型近場動力學計算中動力學迭代時間步長的理論推導目前還需進一步深入研究和討論.

通過對人工阻尼系數和時間步長尺寸的敏感性分析可知，在基于態型近場動力學的準靜力問題分析中，阻尼系數和時間步長的選取對計算效率影響非常關鍵.若選取合適的阻尼系數和時間步長，則計算將快速收斂（如圖3和4）.此外，依靠最大迭代數終止計算具有很強的不穩定性且無法反映計算結果的收斂情況，因此引入收斂準則顯得極為重要.

根據式（19）的失衡力準則，分別選取不同收斂容差，進行分組計算.分別取C=1×108 kg/（m3·s）和Δt=2×10-7 s，計算結果見表1：當收斂容差-=1×10-4時，計算誤差接近1%；當收斂容差繼續減小時，計算誤差僅有細微變化，而收斂迭代步數不斷增加.綜合考慮計算的精度和效率，收斂容差一般取1×10-4即可.

3.2 含中心裂紋板的裂紋擴展模擬

在圖2所示的混凝土板中心位置設置20 mm×1 mm的初始裂紋，并在豎直方向施加均布拉伸載荷，計算混凝土材料的抗壓強度和抗拉強度.對含不同角度中心裂紋的混凝土板進行分級準靜力加載，觀察裂紋開展情況，并確定其開裂載荷.時間步長取為Δt=2×10-7s，局部阻尼系數取C=2×108 kg/（m3·s），收斂容差取為-=1×10-4.

在上述平面模型中，近場范圍為3倍物質點間距即δ=3 mm，裂紋寬度小于近場范圍，在形成物質點對時斷開所有跨越裂紋線的物質點對形成初始裂紋.

對于典型的平面I-II復合型裂紋擴展問題，根據脆性材料的最大拉應力準則可確定起裂角應滿足sin θ0+（3cos θ0-1）cot β=0，其中θ0為起裂角，β為裂紋與加載方向夾角.當β分別為30°，45°，60°和90°時，理論計算所得起裂角θ0分別為60°，53.1°，43.2°和0°.

逐級加載過程如下：在初始加載時，以0.030 MPa為載荷增量進行施加，每次確定系統平衡后施加下一增量步；在系統出現損傷后，退回一個增量步將載荷增量降低為0.003 MPa進行加載，模擬裂紋的穩定擴展，直至最后失穩擴展，系統無法平衡，不滿足式（19）.

圖5a～5d分別顯示含45°初始裂紋混凝土板中從裂尖出現損傷、演化成裂紋起裂直至最后裂紋失穩擴展直至貫穿的全過程.根據分級加載算法可以確定，含45°初始裂紋混凝土板的最終破壞載荷為0.441 MPa，即載荷增加到0.441 MPa時，系統無法達到平衡狀態，無法承受下一載荷增量，隨著動力學方程的繼續迭代，裂紋自然的快速失穩擴展.

其他含不同初始角度裂紋混凝土板的破壞過程與含45°初始裂紋板類似，均由彈性變形、損傷累積、起裂和穩定擴展、失穩擴展等階段組成，主要區別在于最終的破壞載荷以及裂紋擴展路徑，詳見圖6.由圖6可知，盡管初始裂紋的傾角不同，裂紋最終基本會沿著垂直于加載方向快速擴展，與試驗結果一致.各種情況下裂紋的起裂角也與根據最大周向拉應力計算的理論結果較為吻合.此外，由模擬結果可知，當初始裂紋與加載方向夾角β增大時，I型裂紋所占的比例增加，板的承載能力也隨之下降，最終破壞載荷減小.這一結果也符合實際情況.

4 結 論

（1）考慮混凝土類材料的拉壓異性和宏觀斷裂

特征，建立混凝土新的態型近場動力學本構模型，并通過引入動態松弛、系統失衡判斷、力邊界等效等系列算法，構建完整的適合于求解混凝土材料和結構變形破壞過程的態型近場動力學數值體系.

（2）通過對典型混凝土板的準靜態定量變形分析，驗證本文模型和算法的收斂性、計算精度和計算效率，并討論局部阻尼、時間步長和收斂容差等計算變量的選取對態型近場動力學計算精度和效率的影響，給出合理取值范圍建議.

（3）通過分級加載實現混凝土類材料和結構從彈性變形、損傷累積、裂紋穩定擴展以及失穩破壞全過程的連續模擬與定量分析，并成功開展含不同角度中心裂紋混凝土板的破壞過程模擬，模擬所得起裂角和裂紋擴展路徑與試驗結果和理論值吻合，對結構破壞載荷的定量計算結果與變化規律符合實際情況.

參考文獻：

[1] SILLING S A. Reformulation of elasticity theory for discontinuities and long-range forces[J]. Journal of the Mechanics and Physics of Solids， 2000， 48（1）：175-209. DOI： 10.1016/S0022-5096（99）00029-0.

[2] 黃丹， 章青， 喬丕忠， 等. 近場動力學方法及其應用[J]. 力學進展， 2010， 40（4）： 448-459.

HUANG D， ZHANG Q， QIAO P Z， et al. A review on peridynamics（PD） method and its application[J]. Advance in Mechanics， 2010， 40（4）： 448-459.

[3] LAI X， REN B， FAN H， et al. Peridynamics simulations of geomaterial fragmentation by impulse loads[J]. International Journal for Numerical and Analytical Methods in Geomechanics， 2015， 39（12）：1304-1330. DOI： 10.1002/nag.2356.

[4] FAN H， BERGEL G L， LI S. A hybrid peridynamics-SPH simulation of soil fragmentation by blast loads of buried explosive[J].International Journal of Impact Engineering， 2016（87）： 14-27. DOI： 10.1016/j.ijimpeng.2015.08.006.

[5] SILLING S A， ASKARI A. Peridynamic model for fatigue cracking， SAND2014-18590[R]. Sandia Report. New Mexico， 2014.

[6] SILLING S A， DEMMIE P， WARREN T L. Peridynamic simulation of high-rate material failure， SAND2007-3464C[R]. Proceedings of 2007 Applied Mechanics and Materials conference. Austin TX， 2007.

[7] DEMMIE P N， SILLING S A. An approach to modeling extreme loading of structure using peridynamics[J]. Journal of Mechanics of Materials and Structures， 2007， 2（10）： 1921-1945.

[8] KILIC B. Peridynamic theory for progressive failure prediction in homogeneous and heterogeneous materials[D].Tucson： The University of Arizona， 2008.

[9] KILIC B， MADENCI E. Structure stability and failure analysis using peridynamic theory[J]. International of Journal of Non-linear Mechanics， 2009， 44（8）： 845-854.

[10] GERSTLE W， SAU N， SILLING S A. Peridynamic modeling of plain and reinforced concrete structures[C]// Proceedings of 18th International Conference Structural Mechanics in Reactor Technology（SMiRT18）. Beijing， 2005.

[11] GERSTLE W， SAU N， AGUILERA E. Micropolar peridynamic constitutive model for concrete[C]// Proceedings of 19th International Conference Structural Mechanics in Reactor Technology （SMiRT19）， Toronto，2007： 12-17.

[12] GERSTLE W H， SAU N， SAKHAVAND N. On peridynamic computational simulation of concrete structures[J]. ACI Special Publication， 2009， 65： 245-264.

[13] 沈峰， 章青， 黃丹， 等. 基于近場動力學理論的混凝土軸拉破壞過程模擬[J]. 計算力學學報， 2013， 30（S1）： 79-83. DOI： 10.7511/jslx2013z017.

SHEN F， ZHANG Q， HUANG D， et al. Damage and failure process of concrete structure under uni-axialtension based on peridynamics modeling[J]. Chinese Journal of Computational Mechanics， 2013， 30（S1）： 79-83. DOI： 10.7511/jslx2013z017.

[14] HUANG D， ZHANG Q， QIAO P Z. Damage and progressive failure of concrete structures using non-local peridynamic modeling[J].Science China： Technological Science， 2011， 54（3）： 591-596. DOI： 10.1007/s11431-011-4306-3.

[15] HUANG D， LU G， QIAO P. An improved peridynamic approach for quasi-static elastic deformation and brittle fracture analysis[J]. International Journal of Mechanical Sciences， 2015（94/95）： 111-122. DOI： 10.1016/j.ijmecsci.2015.02.018.

[16] SILLING S A， EPTON M， WECKNER O， et al. Peridynamic states and constitutive modeling[J]. Journal of Elasticity， 2007， 88（2）： 151-184. DOI： 10.1007/s10659-007-9125-1.

[17] SILLING S A. Linearized theory of peridynamic states[J]. Journal of Elasticity， 2010， 99（1）： 85-111. DOI： 10.1007/s10659-009-9234-0.

[18] LE Q V， CHAN W K， SCHWARTZ J. A two-dimensional ordinary， state-based peridynamic model for linearly elastic solids[J]. International Journal for Numerical Methods in Engineering 2014， 98（8）： 547-561. DOI： 10.1002/nme.4642.

[19] WU C T. Kinematic constraints in the state-based peridynamics with mixed local/nonlocal gradient approximations[J]. Computational Mechanics， 2014， 54（5）： 1255-1267. DOI： 10.1007/s00466-014-1055-8.