馬 峰，田 羽，張小飛

(1.甘肅省水利水電勘測(cè)設(shè)計(jì)研究院有限責(zé)任公司，甘肅 蘭州 730000；2.中國(guó)電建集團(tuán)北京勘測(cè)設(shè)計(jì)研究院有限公司，北京 100024；3.廣西大學(xué)土木建筑工程學(xué)院，廣西 南寧 530004)

0 引 言

混凝土由于具有良好的可塑性、耐久性及經(jīng)濟(jì)性等優(yōu)點(diǎn)，被廣泛應(yīng)用于水利、土木、道路等工程領(lǐng)域中，其性能的優(yōu)劣很大程度會(huì)影響到工程整體的安全性，因此有關(guān)混凝土力學(xué)性能的研究一直被人們所重視。傳統(tǒng)的混凝土宏觀模型一般在研究過(guò)程中會(huì)忽略混凝土內(nèi)部的復(fù)雜結(jié)構(gòu)，將其視為均質(zhì)且各向同性的連續(xù)性材料，這種模型并不能準(zhǔn)確描述混凝土在受到外荷載時(shí)其內(nèi)部從微裂紋的萌生、擴(kuò)展直至最終產(chǎn)生宏觀裂紋的全過(guò)程。

隨著計(jì)算機(jī)技術(shù)的不斷發(fā)展，基于細(xì)觀層次的數(shù)值模擬成為一種研究混凝土力學(xué)性能的新途徑[1]基于細(xì)觀尺度上對(duì)不同材料組分提出的不同假設(shè)，學(xué)者們提出了不同的細(xì)觀力學(xué)模型及研究方法，目前較為常用的細(xì)觀模型有格構(gòu)模型[2]、M-H模型[3]、隨機(jī)骨料模型[4]、隨機(jī)力學(xué)模型[5]、內(nèi)聚力模型[6]等。而近年來(lái)在有關(guān)混凝土裂紋擴(kuò)展方面的研究中，基于內(nèi)聚力模型的研究方法更加被人們所重點(diǎn)關(guān)注，該模型將混凝土視為由骨料、砂漿以及骨料-砂漿間的粘結(jié)面組成的多相型復(fù)合材料，主要通過(guò)在相鄰實(shí)體單元間插入內(nèi)聚力單元來(lái)模擬混凝土內(nèi)部的損傷開裂，采用牽引力分離準(zhǔn)則來(lái)描述內(nèi)聚力單元的損傷和演化過(guò)程，該模型的優(yōu)點(diǎn)在于不僅可以避免裂紋尖端的應(yīng)力奇異性，同時(shí)支持多條裂紋的擴(kuò)展，更加符合實(shí)際情況。

本文基于MATLAB軟件對(duì)前人隨機(jī)骨料投放算法中單個(gè)骨料面積計(jì)算及骨料間交疊判斷方法加以改進(jìn)，編寫了更加高效簡(jiǎn)便的骨料生成與投放程序，同時(shí)通過(guò)自編C++程序，將零厚度內(nèi)聚力單元批量插入到實(shí)體單元邊界處，并選取合理的本構(gòu)關(guān)系，建立了二級(jí)配混凝土內(nèi)聚力模型，之后以物理實(shí)驗(yàn)結(jié)果中抗拉強(qiáng)度和抗拉彈性模量為指標(biāo)，對(duì)模型在單軸拉伸情況下開展數(shù)值模擬研究，從細(xì)觀角度分析混凝土在單軸受拉時(shí)其內(nèi)部從起始損傷直至破壞的全過(guò)程。從而對(duì)本文提出的混凝土模型生成方法的可行性與有效性加以驗(yàn)證。

1 內(nèi)聚力模型的生成方法

1.1 骨料的生成與投放

在實(shí)際工程中所使用的粗骨料一般是通過(guò)人工破碎生成的，大多數(shù)呈現(xiàn)為凸多面體的形狀，通過(guò)X射線對(duì)混凝土進(jìn)行掃描也可以發(fā)現(xiàn)[7-8]，混凝土在任意截面內(nèi)的骨料形狀大致呈凸多邊形。因此對(duì)于二維的數(shù)值模擬來(lái)講，采用凸多邊形的骨料相較圓形骨料能更加真實(shí)地反應(yīng)混凝土的細(xì)觀結(jié)構(gòu)。

有關(guān)骨料的投放算法目前已有很多研究，馬懷發(fā)等[9]提出將骨料形狀按圓形考慮并全部投放，之后分別對(duì)每個(gè)骨料在圓周上取若干點(diǎn)并按逆時(shí)針依次連接形成“基骨料”，再以每個(gè)“基骨料”的較長(zhǎng)邊為直徑向外延拓生成多邊形骨料，該方法雖然可以很好地控制骨料的投放位置及顆粒數(shù)，但所延拓的多邊形骨料的粒徑可能會(huì)超出該級(jí)配粒徑的范圍；任志剛等[10]對(duì)該方法提出了改進(jìn)，對(duì)“基骨料”延拓產(chǎn)生的新節(jié)點(diǎn)到圓心的距離不能超出該級(jí)配的粒徑范圍且滿足“凸”形的條件；汪奔等[11]通過(guò)將生成的骨料進(jìn)行適當(dāng)?shù)目s放也達(dá)到了類似的效果。以上學(xué)者的算法雖能滿足骨料的“凸”性及各級(jí)配的占比要求，但缺點(diǎn)是算法均較為復(fù)雜，對(duì)于編程初學(xué)者不太友好?；诖?，本文對(duì)以上學(xué)者的算法加以借鑒和改進(jìn)，編制了一套用于二維隨機(jī)凸多邊形骨料的生成與投放算法，該方法不僅可以滿足粗骨料的占比及級(jí)配分布的合理性，同時(shí)也更加簡(jiǎn)單易行、便于編程實(shí)現(xiàn)。

1.1.1 隨機(jī)凸多邊形骨料的生成

本文將混凝土視為由骨料、砂漿及骨料-砂漿界面過(guò)渡區(qū)組成的三相復(fù)合性材料，其中直徑D<5 mm的骨料視作砂漿考慮。首先利用瓦拉文公式[12]建立空間內(nèi)骨料占比與二維某一截面骨料面積之間的關(guān)系，使模型可應(yīng)用于二維平面狀態(tài)的分析中，平面內(nèi)任意一點(diǎn)處含有骨料顆粒直徑D

(1)

式中，Pc為不同級(jí)配骨料在二維某截面內(nèi)出現(xiàn)的概率；Pk為骨料占混凝土截面總面積的百分比；Dmax為最大骨料的顆粒直徑；D0為某一級(jí)配骨料的顆粒直徑。

對(duì)于任意一個(gè)隨機(jī)凸多邊形骨料來(lái)講，其大小和形狀都是隨機(jī)的，本文通過(guò)在極坐標(biāo)下生成隨機(jī)的極半徑R以及極角alpha的方法確定凸多邊形骨料的每個(gè)頂點(diǎn)坐標(biāo)，之后將其轉(zhuǎn)換為直角坐標(biāo)，其表達(dá)式為

(2)

式中，R為骨料第i個(gè)頂點(diǎn)的極半徑；r為某級(jí)配下骨料的代表粒徑；alpha為骨料第i個(gè)頂點(diǎn)的極角；n為骨料的邊數(shù)，本文取5～12之間的任意整數(shù)；i為骨料的第i個(gè)頂點(diǎn)。

1.1.2 隨機(jī)凸多邊形骨料的投放

在完成單個(gè)骨料算法的編制后，即可開始投放，投放之前需要確定二維混凝土模型的長(zhǎng)和寬、骨料在各個(gè)級(jí)配的代表粒徑及相應(yīng)占比，由瓦拉文公式即可計(jì)算出各個(gè)級(jí)配骨料的占比，并將其按模型的尺寸轉(zhuǎn)化為面積。具體的投放步驟為

(1)根據(jù)極坐標(biāo)的方法生成各個(gè)級(jí)配的骨料庫(kù)，并記錄下每個(gè)骨料的頂點(diǎn)及相應(yīng)外圍節(jié)點(diǎn)的坐標(biāo)和骨料的面積。

(2)按照骨料粒徑的大小在指定區(qū)域內(nèi)逐個(gè)投放骨料，并判斷新生成的骨料與之前的骨料是否產(chǎn)生重疊，若重疊則重新投放新的骨料，反之記錄下相應(yīng)的骨料頂點(diǎn)坐標(biāo)，并將該骨料的面積累加至已完成投放的骨料面積中。

(3)判斷各個(gè)級(jí)配已投放的骨料面積是否達(dá)到或接近所計(jì)算的值，由于每次投放很難達(dá)到完全精確，故設(shè)置了3%的容錯(cuò)率，若未達(dá)到預(yù)期值，則繼續(xù)投放，反之則結(jié)束該級(jí)配骨料的投放，進(jìn)行下一級(jí)配，直至完成所有級(jí)配骨料的投放。

本文隨機(jī)凸多邊形骨料的生成與投放的流程見圖1。

圖1 骨料生效與投放流程

本文借鑒馬懷發(fā)[9]、任志剛[10]、汪奔[11]等學(xué)者有關(guān)隨機(jī)骨料投放算法的研究，并對(duì)算法加以改進(jìn)。改進(jìn)有：①對(duì)單個(gè)骨料面積計(jì)算的方法進(jìn)行了改進(jìn)，采用MATLAB中的polyarea(x，y)函數(shù)，將單個(gè)骨料的頂點(diǎn)坐標(biāo)按逆時(shí)針?lè)较蛞来未鎯?chǔ)在所創(chuàng)建的x和y向量中，然后利用該函數(shù)即可返回該骨料的面積，該方法相較馬懷發(fā)、任志剛等將多邊形按頂點(diǎn)分割成三角形累加計(jì)算的方法更加簡(jiǎn)便且更加高效，也節(jié)省了編制相應(yīng)代碼的時(shí)間。②對(duì)骨料的交疊判斷采取了不同的算法，如圖2所示，首先對(duì)每個(gè)骨料頂點(diǎn)坐標(biāo)擴(kuò)大1.1倍并按逆時(shí)針?lè)较蛞来螌⑾噜弮牲c(diǎn)之間均分成8等分，形成一圈骨料外圍節(jié)點(diǎn)，然后利用inpolygon(x，y，xv，yv)函數(shù)判斷每個(gè)新生成的骨料外圍節(jié)點(diǎn)坐標(biāo)集合(x，y)是否在每個(gè)已生成的多邊形骨料外圍節(jié)點(diǎn)坐標(biāo)集合(xv，yv)的邊界內(nèi)，如果函數(shù)返回值為1，則骨料頂點(diǎn)產(chǎn)生交疊，投放失敗，如果不在邊界內(nèi)，則表明該骨料投放成功，并記錄下相應(yīng)的頂點(diǎn)坐標(biāo)。該方法運(yùn)行效率高且代碼易于編制，僅需記錄下各骨料的頂點(diǎn)及外圍節(jié)點(diǎn)坐標(biāo)即可完成交疊判斷。③取消將投放骨料顆粒個(gè)數(shù)作為最終某一級(jí)配下所有骨料投放成功的判定條件，對(duì)于隨機(jī)多邊形骨料來(lái)講，每個(gè)骨料的面積都是不等的，如果將骨料的個(gè)數(shù)列為判定條件，不僅會(huì)增加判定次數(shù)，降低運(yùn)行效率，同時(shí)也變相加入了人為的干擾因素。

圖2 骨料交疊判斷示意

1.2 內(nèi)聚力單元的批量插入

雖然通常情況下混凝土內(nèi)部的開裂主要產(chǎn)生于骨料-砂漿交界面、水泥砂漿內(nèi)部，但仍有部分骨料可能會(huì)產(chǎn)生斷裂，為了更加精確的描述混凝土內(nèi)部的損傷開裂，將骨料內(nèi)部也視為潛在開裂區(qū)。本文通過(guò)自編的C++算法對(duì)inp文件進(jìn)行批處理，實(shí)現(xiàn)了零厚度內(nèi)聚力單元精確且高效的批量插入。具體的算法步驟為

(1)通過(guò)上述算法建立隨機(jī)骨料模型，并進(jìn)行網(wǎng)格劃分，本文采用自由劃分網(wǎng)格法，對(duì)骨料和砂漿區(qū)域分別創(chuàng)建集合，然后導(dǎo)出該模型的inp文件。

(2)讀取inp文件中所有的節(jié)點(diǎn)及單元信息，并按照骨料與砂漿建立相應(yīng)集合進(jìn)行存儲(chǔ)，之后分別對(duì)各單元集合內(nèi)的節(jié)點(diǎn)進(jìn)行遍歷，記錄下每個(gè)節(jié)點(diǎn)出現(xiàn)的次數(shù)，每個(gè)節(jié)點(diǎn)重復(fù)出現(xiàn)了幾次就表明需要分裂幾次，并根據(jù)原節(jié)點(diǎn)信息在相同的位置生成新的節(jié)點(diǎn)，如圖3所示，節(jié)點(diǎn)3由于被6個(gè)單元所公用，因此需要將其分裂5次。

圖3 算法示意

(3)將新產(chǎn)生的節(jié)點(diǎn)重新分配給相應(yīng)的實(shí)體單元，并在各實(shí)體單元間創(chuàng)建零厚度的內(nèi)聚力單元，并分別創(chuàng)建相應(yīng)的集合，如圖3所示，綠色為骨料-砂漿界面單元，橙色為骨料內(nèi)部單元，紫色為砂漿內(nèi)部單元。

(4)將各集合內(nèi)新的節(jié)點(diǎn)、單元信息重新覆蓋之前的inp文件，并導(dǎo)入ABAQUS，完成零厚度內(nèi)聚力單元的插入。

以建立的二級(jí)配混凝土模型(尺寸為150 mm×150 mm，骨料占比約為46.9%)為例，根據(jù)本文的骨料生成與投放算法及內(nèi)聚力單元的批量插入算法所建模型如圖4所示，分別按集合顯示了各單元分布情況。圖中模型共37 614個(gè)節(jié)點(diǎn)，由18 657個(gè)三角形平面應(yīng)力單元(CPS3)以及12 538個(gè)內(nèi)聚力單元組成，其中采用本文C++程序插設(shè)內(nèi)聚力單元僅需8 min左右，說(shuō)明本文編制算法具有高效性和可行性。

圖4 二級(jí)配混凝土模型

2 混凝土單軸受拉破壞細(xì)觀模擬

2.1 模型的生成

為驗(yàn)證本研究提出的內(nèi)聚力模型生成算法的可行性和有效性，利用王仲華等[13]以灰?guī)r作為粗、細(xì)骨料時(shí)所開展的二級(jí)配混凝土(28 d齡期)單軸拉伸試驗(yàn)結(jié)果作為數(shù)值模擬的對(duì)象，試驗(yàn)所得混凝土的抗拉強(qiáng)度平均值為2.36 MPa，抗拉彈性模量為37.7 GPa。

本試驗(yàn)采用開發(fā)的骨料生成與投放算法及內(nèi)聚力單元的批量插入算法，生成150 mm×150 mm的二級(jí)配混凝土內(nèi)聚力模型，選用全局尺寸為2 mm的三角形網(wǎng)格進(jìn)行單元?jiǎng)澐?，荷載施加方式為位移勻速加載，對(duì)試件的底邊施加豎向的位移約束，同時(shí)最中間節(jié)點(diǎn)施加水平位移約束，圖5為模型荷載施加示意。

圖5 二級(jí)配混凝土細(xì)觀模型

2.2 本構(gòu)關(guān)系的選取

本文模型中，骨料、砂漿單元采用三角形平面應(yīng)力實(shí)體單元(CPS3)劃分網(wǎng)格，為提高計(jì)算效率，采用線彈性的本構(gòu)關(guān)系。骨料-砂漿界面單元、骨料內(nèi)部單元、砂漿內(nèi)部單元采用內(nèi)聚力單元(COH2D4)模擬，內(nèi)聚力單元描述裂縫尖端的斷裂方法為：假設(shè)在真實(shí)裂紋尖端之前存在一個(gè)“斷裂過(guò)程區(qū)”，該區(qū)域裂尖的張開位移與應(yīng)力的關(guān)系采用牽引力分離準(zhǔn)則來(lái)表示。目前常用的牽引力分離準(zhǔn)則有雙線型、梯型、指數(shù)型[14]等，雙線型牽引力分離準(zhǔn)則適用于描述脆性材料的斷裂，因此本文采用該準(zhǔn)則來(lái)描述混凝土的斷裂。

目前ABAQUS中共有4種損傷起裂準(zhǔn)則：最大名義應(yīng)力準(zhǔn)則、最大名義應(yīng)變準(zhǔn)則、平方名義應(yīng)力準(zhǔn)則和平方名義應(yīng)變準(zhǔn)則?；炷猎趩屋S受壓條件下的破壞可能受到法向和切向應(yīng)力的共同作用，同時(shí)存在I型(張拉)和II型(剪切)兩種斷裂模式，因此本文考慮了兩種裂紋模式的組合，采用了平方名義應(yīng)力準(zhǔn)則來(lái)判定單元的損傷起裂，即

(3)

采取了基于能量演化準(zhǔn)則中的Benzeggagh-Kenane(B-K)混合斷裂準(zhǔn)則來(lái)描述單元的損傷演化，該準(zhǔn)則適用于包含I型、II型裂紋的混合型斷裂，表達(dá)式如下：

(4)

2.3 參數(shù)校調(diào)

一般混凝土內(nèi)部各相組分的力學(xué)性能關(guān)系為：骨料>水泥砂漿>骨料-砂漿交界面[15]，因此本著這一原則對(duì)砂漿內(nèi)部及骨料-砂漿粘結(jié)面單元的參數(shù)進(jìn)行校調(diào)，通過(guò)不斷修正單元材料參數(shù)直至與物理試驗(yàn)結(jié)果相吻合。為使結(jié)果更加精確，共建立3個(gè)150 mm×150 mm的二級(jí)配混凝土試件，求解后取其結(jié)果的平均值，最終確定的材料參數(shù)見表1。

表1 材料參數(shù)

對(duì)3個(gè)150 mm×150 mm的二級(jí)配混凝土模擬試件分別進(jìn)行求解，得到圖6的單軸拉伸應(yīng)力應(yīng)變曲線，具體計(jì)算結(jié)果為：二級(jí)配混凝土中3個(gè)試件的抗拉強(qiáng)度最大值分別為2.26、2.31、2.35 MPa，抗拉強(qiáng)度平均峰值應(yīng)力為2.29 MPa；抗拉彈性模量為35.38、35.12、35.11 GPa，平均值為35.20 GPa，與王仲華[13]的試驗(yàn)數(shù)據(jù)對(duì)比，抗拉強(qiáng)度最大誤差不超過(guò)2.8%，抗拉彈性模量最大誤差不超過(guò)6.6%。

圖6 單軸拉伸應(yīng)力應(yīng)變曲線

2.4 單軸受拉破壞全過(guò)程模擬

以試件1為例，對(duì)混凝土在單軸受拉時(shí)的全過(guò)程進(jìn)行分析，圖7給出了試件的單軸受壓應(yīng)力-應(yīng)變曲線，A、B、C、D等4個(gè)點(diǎn)分別對(duì)應(yīng)圖8中4個(gè)具有代表性的損傷演化階段。為了清楚的顯示不同類型單元的損傷破壞，用紅色代表骨料-砂漿粘結(jié)面單元，綠色代表砂漿內(nèi)部單元，黃色代表骨料內(nèi)部單元，c表示試件在受拉方向所放大的倍數(shù)。

圖7 試件1單軸拉伸應(yīng)力應(yīng)變曲線

圖8 試件1損傷演化示意

(1)微裂紋成核階段。此過(guò)程發(fā)生于應(yīng)力上升階段，該階段混凝土試件已經(jīng)出現(xiàn)了塑性變形，由于骨料與砂漿單元在力學(xué)性能上差異較大，在其交界面處容易產(chǎn)生應(yīng)力集中，導(dǎo)致骨料-砂漿粘結(jié)面上已經(jīng)產(chǎn)生損傷，在骨料表面生成了一些與受拉方向接近垂直的微裂紋。

(2)微裂紋擴(kuò)展階段。隨著混凝土在受拉方向的應(yīng)變不斷增大，這些微裂紋向相鄰骨料的方向逐漸擴(kuò)展，即擴(kuò)展至了砂漿內(nèi)部單元區(qū)域，由圖8b可以發(fā)現(xiàn)，在微裂紋的擴(kuò)展過(guò)程中，相鄰粒徑較大的骨料之間更容易產(chǎn)生連接貫通，這可能是由于隨著骨料粒徑的增大，其表面潛在的微裂紋會(huì)越多，因此相鄰兩顆大粒徑骨料之間的砂漿內(nèi)部單元更傾向于連接貫通，有逐漸形成裂紋貫通帶的趨勢(shì)。

(3)微裂紋連接貫通階段。此時(shí)微裂紋逐漸從骨料表面及其附近區(qū)域擴(kuò)展貫通，在混凝土內(nèi)部大致形成了一條主裂紋帶，由圖8c可知，此時(shí)的骨料-砂漿粘結(jié)面單元及砂漿內(nèi)部單元仍然存在，即代表宏觀裂紋還沒(méi)有真正產(chǎn)生，在主裂紋帶逐漸形成的過(guò)程中，部分區(qū)域的微裂紋由于主裂紋帶的擴(kuò)展而受到限制，甚至?xí)S著主裂紋帶的繼續(xù)增加而逐漸呈現(xiàn)閉合的趨勢(shì)。

(4)宏觀裂紋形成階段：該階段大部分骨料-砂漿粘結(jié)面單元及砂漿內(nèi)部單元已經(jīng)完全損傷而被刪除，最終形成了一條肉眼可見的宏觀裂紋帶，裂紋帶近似垂直于受拉方向，但這條裂紋帶并非是一條連續(xù)的曲線，這可能是由于粗骨料的形狀影響了裂紋的擴(kuò)展路徑所導(dǎo)致的。在整個(gè)單軸拉伸模擬過(guò)程中，骨料內(nèi)部并未發(fā)現(xiàn)黃色單元，表明骨料內(nèi)部單元沒(méi)有產(chǎn)生損傷。

3 結(jié) 論

(1)基于MATLAB對(duì)前人隨機(jī)骨料投放算法進(jìn)行改進(jìn)，優(yōu)化了單個(gè)骨料面積計(jì)算及骨料間交疊判斷的方法，編寫了更加高效簡(jiǎn)便的骨料生成與投放算法，之后采用C++編程，實(shí)現(xiàn)在相鄰實(shí)體單元間快速、成批次地插入內(nèi)聚力單元，并根據(jù)單元所處的部位進(jìn)行集合劃分，完成內(nèi)聚力模型的建立。

(2)以王仲華所做的物理試驗(yàn)結(jié)果中抗拉強(qiáng)度和抗拉彈性模量為指標(biāo)，對(duì)本文建立的內(nèi)聚力模型進(jìn)行校調(diào)，與試驗(yàn)數(shù)據(jù)對(duì)比，抗拉強(qiáng)度最大誤差不超過(guò)2.8%，抗拉彈性模量最大誤差不超過(guò)6.6%。

(3)從單軸拉伸模擬結(jié)果來(lái)看，混凝土在均勻受拉時(shí)大致經(jīng)歷了4個(gè)階段：微裂紋成核階段、微裂紋擴(kuò)展階段、微裂紋連接貫通階段和宏觀裂紋形成階段，這與實(shí)際中混凝土單軸受拉全過(guò)程大致相似，從而驗(yàn)證了本研究提出的內(nèi)聚力模型生成算法的可行性和有效性。