巴西圓盤劈裂二維及三維數(shù)值模擬研究

喬 蘭，劉 建，李慶文?，趙國(guó)彥

1) 北京科技大學(xué)土木與資源工程學(xué)院，北京 100083 2) 中南大學(xué)資源與安全工程學(xué)院，長(zhǎng)沙 410083

巖石、混凝土和陶瓷等脆性材料的抗拉強(qiáng)度遠(yuǎn)低于其抗壓強(qiáng)度及抗剪強(qiáng)度，巖石的宏觀剪切破壞一般源于其內(nèi)部大量張拉裂紋的相互作用及貫通[1?2]，Hoek 和 Martin[3]也曾指出深部巖體的破壞主要是應(yīng)力驅(qū)動(dòng)的張拉破壞，因此巖石類脆性材料的抗拉強(qiáng)度是人們尤其關(guān)注的重要參數(shù).相對(duì)于直接拉伸試驗(yàn)，以巴西劈裂試驗(yàn)為代表的抗拉強(qiáng)度間接測(cè)試方法因其操作簡(jiǎn)單方便而在國(guó)內(nèi)外得到廣泛應(yīng)用.國(guó)際巖石力學(xué)學(xué)會(huì)(ISRM)于1978年將巴西劈裂試驗(yàn)推薦為測(cè)定巖石材料抗拉強(qiáng)度的建議方法[4].與此同時(shí)，自巴西劈裂試驗(yàn)方法于1943年被提出以來(lái)，關(guān)于其有效性的探討一直是國(guó)內(nèi)外學(xué)者的研究熱點(diǎn)[5?7].根據(jù)Griffith強(qiáng)度理論，只有在滿足中心起裂的條件下通過(guò)巴西劈裂試驗(yàn)測(cè)得的巖石抗拉強(qiáng)度才是有效的[5]，并且從破壞結(jié)果來(lái)看其破裂面至少應(yīng)表現(xiàn)為與張拉破壞對(duì)應(yīng)的單一徑向破裂面.但是，F(xiàn)airhurst[5]指出對(duì)于壓拉比較低的脆性材料采用小加載角時(shí)起裂點(diǎn)可能由于應(yīng)力集中效應(yīng)而位于巖樣端部.Erarslan等[8?9]、Komurlu 等[10?11]采用 ISRM 推薦的弧形壓模開(kāi)展巴西劈裂試驗(yàn)，發(fā)現(xiàn)圓盤試樣發(fā)生多破裂壓碎破壞.Jaeger和Hoskins[12]提出采用弧形均布載荷開(kāi)展巴西劈裂試驗(yàn)，通過(guò)增大加載角以減小加載端的應(yīng)力集中程度從而促進(jìn)中心起裂.Erarslan 等[8?9]認(rèn)為加載角為 20°或 30°時(shí)可以獲得可靠的抗拉強(qiáng)度.對(duì)于混凝土巴西圓盤試樣，García等[13]建議加載角應(yīng)大于 12°.Gutiérrez-Moizant等[14]則認(rèn)為開(kāi)展巴西劈裂試驗(yàn)最合適加載角為20°.Bahaaddini等[15]通過(guò)開(kāi)展數(shù)值模擬研究發(fā)現(xiàn)，對(duì)于巖石類材料加載角為30°或更高時(shí)才能確保中心起裂.目前，鮮少有學(xué)者探究壓拉比對(duì)試樣起裂破壞模式的影響.此外，國(guó)內(nèi)外學(xué)者比較關(guān)注圓盤試樣受壓直徑上的起裂點(diǎn)位置及破裂發(fā)展方式，較少有學(xué)者探究巴西劈裂三維破裂面的起裂及擴(kuò)展過(guò)程.由于張拉破裂面沿受壓直徑向圓盤內(nèi)部擴(kuò)展，受限于目前的觀測(cè)手段，以高速攝影儀[16]及數(shù)字圖像相關(guān)技術(shù)[17]為代表的表面變形觀測(cè)技術(shù)便不再適用，而數(shù)值模擬方法便顯示出其獨(dú)特的優(yōu)勢(shì).喻勇[6]、Li和 Wong[7]及鄧華鋒等[18]通過(guò)開(kāi)展三維彈性分析探討圓盤內(nèi)部的應(yīng)力分布規(guī)律及起裂點(diǎn)位置，但未模擬三維破裂面的形成過(guò)程.許學(xué)良[19]運(yùn)用PFC3D開(kāi)展巴西劈裂三維模擬，以細(xì)觀裂紋集聚表征張拉破裂面的擴(kuò)展過(guò)程，模擬結(jié)果中破裂面仍然從加載端開(kāi)始起裂，并且起裂階段初始裂紋全部為張拉裂紋.鑒于此，本文借助FLAC3D進(jìn)行連續(xù)介質(zhì)彈塑性分析，一方面開(kāi)展巴西劈裂二維模擬研究，探究圓盤試樣壓拉比及加載角對(duì)試樣起裂破壞模式的影響；另一方面開(kāi)展三維模擬研究，探究圓盤試樣三維破裂面的形成及擴(kuò)展過(guò)程.

1 巴西圓盤應(yīng)力解析解

Hondros[20]推導(dǎo)出平面問(wèn)題中弧形均布載荷作用下巴西圓盤內(nèi)部任意一點(diǎn)應(yīng)力分布的解析解，如圖1所示，則加載直徑上的應(yīng)力分布為（以拉應(yīng)力為正）：

圖1 巴西圓盤弧形均布載荷加載示意圖Fig.1 A disk subjected to diametrically distributed loads in the Brazilian test

式中：σr,σθ,τrθ分別為受壓直徑上的徑向應(yīng)力、周向應(yīng)力及剪應(yīng)力，P為施加的外部載荷，D為巴西圓盤直徑，R為圓盤半徑，t為圓盤厚度，2α為加載角度，r為受壓直徑上一點(diǎn)到圓心的距離.在滿足中心起裂的基礎(chǔ)上，巖石抗拉強(qiáng)度的計(jì)算公式為

當(dāng)接觸加載角小于15°時(shí)，采用近似公式計(jì)算巖石抗拉強(qiáng)度引起的誤差不大于2.3%.此外，ISRM推薦的弧形壓模與巴西圓盤接觸在試樣破壞時(shí)產(chǎn)生的加載角為10°左右[4].

2 二維數(shù)值模擬

二維巴西劈裂數(shù)值模擬采用平面應(yīng)變模型，數(shù)值試樣直徑為 50 mm，厚度為 0.5 mm.模型網(wǎng)格采用六面體立方體單元，單元尺寸為0.5 mm，整個(gè)模型共包含7864個(gè)單元，16130個(gè)節(jié)點(diǎn).加載過(guò)程采用位移控制，邊界條件為：固定模型y方向的變形，在模型頂端和底端施加z方向的速度邊界條件，加載速率為每步 1×10?10m.二維加載示意圖如圖2所示.采用彈塑性應(yīng)變軟化本構(gòu)模型及Mohr-Coulomb強(qiáng)度準(zhǔn)則，參考Iddefjord花崗巖[21]的力學(xué)參數(shù)，數(shù)值模型參數(shù)如表1所示.當(dāng)模型單元發(fā)生屈服后，單元內(nèi)聚力及抗拉強(qiáng)度隨塑性剪應(yīng)變及塑性拉應(yīng)變累計(jì)逐漸線性軟化；當(dāng)塑性剪應(yīng)變及塑性拉應(yīng)變累計(jì)至臨界值后，單元內(nèi)聚力及抗拉強(qiáng)度弱化至殘余強(qiáng)度并保持不變.為探究圓盤試樣壓拉比及加載角對(duì)試樣起裂破壞模式的影響，分別考慮 6 種壓拉比 (λ)及加載角 (2α)：λ=20、18、16、14、12和 10及 2α=4.58°、6.88°、9.18°、11.48°、13.78°和 16.10°.針對(duì)每種壓拉比，保持單軸抗壓強(qiáng)度為200 MPa不變.數(shù)值模擬結(jié)果如表2所示.

圖2 二維巴西劈裂加載示意圖Fig.2 Loading setup of the 2D numerical Brazilian test

表1 數(shù)值模型參數(shù)Table 1 Material properties of the numerical model

表2 不同壓拉比及加載角巴西試樣的破壞模式Table 2 Failure modes of the numerical Brazilian disks with different compression?tension ratios (λ) and contact loading angles (2α)

從表2橫向分析，接觸加載角越大，巴西圓盤試樣越容易發(fā)生中心起裂，因?yàn)樵龃蠹虞d角可以減小圓盤試樣頂?shù)锥说膽?yīng)力集中程度.以λ=14為例，當(dāng) 2α≤ 9.18°時(shí)，巴西試樣發(fā)生端部起裂；當(dāng) 2α≥ 11.48°時(shí)，巴西試樣發(fā)生中心起裂.沿表2縱向分析，壓拉比越大，巴西圓盤試樣越容易發(fā)生中心起裂.以 2α=11.48°為例，當(dāng)λ≥ 14 時(shí)，試樣發(fā)生中心起裂；而當(dāng)λ≤ 12時(shí)，圓盤試樣發(fā)生端部起裂.較大的壓拉比，一方面表示巖樣的抗拉強(qiáng)度相對(duì)較低，因此試樣更容易發(fā)生張拉破壞；另一方面則表示巖樣的抗壓強(qiáng)度相對(duì)較高，抵抗剪切破壞的能力強(qiáng).因此，巴西圓盤試樣的起裂破壞模式受內(nèi)部因素與外部因素共同控制，總體上，圓盤試樣壓拉比越大，外部加載角越大，試樣破壞越容易產(chǎn)生中心起裂，反之則越容易發(fā)生端部起裂.

圖3（a）～（b）為均勻巴西圓盤試樣的加載破壞過(guò)程，包括峰值載荷(Fpeak)前及峰后階段，接觸加載角為 9.18°，壓拉比分別為λ=20、10.當(dāng)λ=20時(shí)，圓盤試樣發(fā)生中心起裂，隨外載進(jìn)一步增加，徑向張拉裂紋自中心起裂點(diǎn)逐漸向兩端擴(kuò)展，并在加載點(diǎn)附近趨于停止，這與Erarslan等[8?9]通過(guò)實(shí)驗(yàn)觀察到的裂紋擴(kuò)展模式相似，陳沙等[22]也通過(guò)數(shù)值模擬獲得過(guò)類似的裂紋擴(kuò)展過(guò)程.而當(dāng)λ=10時(shí)，兩條起裂裂紋從圓盤受力區(qū)域邊緣向圓盤內(nèi)部擴(kuò)展，初始起裂裂紋主要發(fā)生剪切破壞，向試樣中心擴(kuò)展后則轉(zhuǎn)變?yōu)閺埨茐模罱K形成兩條貫通的劈裂縫，十分明顯這是由于各向同性圓盤試樣及外載關(guān)于受壓直徑的對(duì)稱性引起.圖3（c）～（e）為非均勻巴西圓盤試樣的加載破壞過(guò)程，加載角仍為 9.18°，壓拉比分別為λ=20、10、8.運(yùn)用劉建等[23]提出的巖石介質(zhì)細(xì)觀物理力學(xué)參數(shù)統(tǒng)計(jì)分布模型，采用對(duì)數(shù)正太分布描述細(xì)觀單元彈性模量及單軸抗壓強(qiáng)度的空間變異性[24?25]，試樣均質(zhì)度參數(shù)設(shè)為0.3.當(dāng)λ=20時(shí)，數(shù)值試樣同樣產(chǎn)生中心起裂，受模型非均質(zhì)性的影響起裂點(diǎn)偏離圓盤正中心，同時(shí)在圓盤內(nèi)部受壓直徑周圍分布有離散的張拉破裂單元.這是由于一方面材料的非均質(zhì)性改變了圓盤內(nèi)部的應(yīng)力分布，引起局部應(yīng)力集中；另一方面，可能是由于軟弱介質(zhì)單元的過(guò)早破壞引起.而當(dāng)λ=10、8時(shí)，與均勻條件下的模擬結(jié)果明顯不同的是，圓盤試樣發(fā)生端部起裂后最終僅形成一條劈裂縫；在劈裂縫形成之前，試樣兩端出現(xiàn)剪切破壞單元，而試樣內(nèi)部基本無(wú)破壞單元出現(xiàn)，壓拉比為8的數(shù)值試樣更明顯，這與楊志鵬等[26]關(guān)于頁(yè)巖巴西劈裂的聲發(fā)射定位結(jié)果類似.值得指出的是，無(wú)論是均勻試樣還是非均勻試樣，端部起裂均由剪切破壞引起，而劈裂裂紋進(jìn)一步擴(kuò)展則由張拉破壞驅(qū)動(dòng).

圖3 不同壓拉比時(shí)巴西圓盤試樣的破壞過(guò)程.（a）均勻試樣 λ=20；（b） 均勻試樣 λ=10；（c）非均勻試樣 λ=20；（d）非均勻試樣 λ=10；（e）非均勻試樣λ=8Fig.3 Failure processes of the numerical disks with different compression–tension ratios (λ) under the loading angle of 9.18°: (a) homogeneous disk λ =20; (b) homogeneous disk λ = 10; (c) heterogeneous disk λ = 20; (d) heterogeneous disk λ = 10; (e) heterogeneous disk λ = 8

3 三維數(shù)值模擬

三維巴西圓盤數(shù)值模型直徑為50 mm，厚度為25 mm；采用六面體網(wǎng)格進(jìn)行劃分，對(duì)受壓直徑附近的網(wǎng)格進(jìn)行適當(dāng)加密，模型共包含69125個(gè)單元，74360個(gè)節(jié)點(diǎn)，模型網(wǎng)格如圖4所示.采用應(yīng)力邊界條件，在受壓區(qū)域內(nèi)沿模型表面施加法向正應(yīng)力，加載速率為每步500 Pa.三維數(shù)值模擬仍采用表1中的力學(xué)參數(shù)；為模擬不同接觸加載角下巴西圓盤三維破裂面的起裂及擴(kuò)展過(guò)程，考慮以下 8種加載角度：2α=29.7°、25.7°、21.7°、17.9°、14.4°、10.9°、7.6°、4.4°.

圖4 三維巴西圓盤模型網(wǎng)格示意圖Fig.4 Grid model of the 3D numerical Brazilian disk

圖5 為加載角 2α=29.7°、25.7°、14.4°和 4.4°時(shí)巴西圓盤試樣三維破裂面的起裂及擴(kuò)展過(guò)程，由于篇幅限制，其他4種加載角下的數(shù)值模擬結(jié)果沒(méi)有給出，但破裂面發(fā)展規(guī)律基本類似.當(dāng)2α=29.7°時(shí)，初始起裂點(diǎn)位于圓盤端面中心，隨外部載荷增加，破裂面以弧形邊界向試樣內(nèi)部擴(kuò)展，在圓盤厚度方向首先貫通后繼續(xù)沿受壓直徑向兩端擴(kuò)展并在加載點(diǎn)附近趨于停止.值得注意的是，貫通后的破裂面沿受壓直徑方向的擴(kuò)展邊界也并非直線，這與平面應(yīng)變假設(shè)下貫通破裂面的直線擴(kuò)展情形相矛盾.當(dāng) 2α=25.7°、14.4°時(shí)，起裂點(diǎn)仍位于圓盤端面受壓直徑上，但卻偏離端面中心，同時(shí)出現(xiàn)4個(gè)對(duì)稱的起裂點(diǎn)；加載角越小，起裂點(diǎn)越靠近上下加載端，并且破裂面擴(kuò)展邊界線仍為弧形.當(dāng)2α=4.4°時(shí)，起裂點(diǎn)位于圓盤端面加載區(qū)域邊緣，即破裂沒(méi)有在受壓中心面上發(fā)展，這與二維模擬中端部起裂的情況類似.事實(shí)上，8種接觸加載角下圓盤試樣的起裂點(diǎn)均位于圓盤端面，而不是圓盤內(nèi)部中心.喻勇[8]及鄧華鋒等[20]通過(guò)對(duì)線荷載作用下的巴西圓盤試樣進(jìn)行三維彈性分析指出，圓盤橫截面上的應(yīng)力值沿試樣厚度方向是有變化的，越靠近兩端面，水平拉應(yīng)力越大，因此試樣破壞首先從端面開(kāi)始.

圖5 不同加載角時(shí)巴西圓盤試樣三維破裂面的起裂及擴(kuò)展過(guò)程.（a）2α=29.7°；（b）2α=25.7°；（c）2α=14.4°；（d）2α=4.4°Fig.5 Fracture initiation and propagation processes of the 3D Brazilian disks with different contact loading angles: (a) 2α = 29.7°; (b) 2α = 25.7°; (c) 2α =14.4°; (d) 2α = 4.4°

4 解析驗(yàn)證及應(yīng)力分析

4.1 二維模擬結(jié)果分析

圖6為λ=20時(shí)不同加載角下二維巴西圓盤受壓直徑上切向應(yīng)力與徑向應(yīng)力的數(shù)值模擬結(jié)果與理論計(jì)算結(jié)果的對(duì)比示意圖.由圖可知，除在加載端附近數(shù)值模擬結(jié)果與解析解存在較小誤差外，受壓直徑上大部分區(qū)域具有較好的吻合性，說(shuō)明了二維數(shù)值模擬的有效性.需要指出的是，由于二維數(shù)值模型采用正六面體單元逼近圓盤試樣進(jìn)行網(wǎng)格劃分，因此在模型頂?shù)锥藶橹本€邊界而非弧形邊界，并且施加的速度邊界條件與受壓直徑平行而非沿圓盤半徑指向圓盤中心，如圖3所示，這是造成上述誤差的原因.在加載端附近，徑向應(yīng)力與切向應(yīng)力均為壓應(yīng)力，靠近圓盤中心，兩者快速減小，并且切向應(yīng)力由壓應(yīng)力逐漸轉(zhuǎn)變?yōu)槔瓚?yīng)力.隨著接觸加載角增大，受壓直徑上張拉區(qū)域逐漸減小，加載端附近壓縮區(qū)域逐漸增大，徑向與切向壓應(yīng)力也逐漸減小，這說(shuō)明端部應(yīng)力集中程度在減弱.

圖6 不同加載角時(shí)二維巴西圓盤受壓直徑上應(yīng)力分布的數(shù)值模擬結(jié)果與理論計(jì)算結(jié)果對(duì)比示意圖.（a） 2α=4.58°；（b）2α=6.88°；（c）2α=9.18°；（d）2α=11.48°；（e）2α=13.78°；（f）2α=16.10°Fig.6 Comparison of normalized stresses along the compressed diameter between the numerical results and Hondros’ solutions with different contact loading angles: (a) 2α = 4.58°; (b) 2α = 6.88°; (c) 2α = 9.18°; (d) 2α = 11.48°; (e) 2α = 13.78°; (f) 2α = 16.10°

表3為不同加載角及壓拉比下二維巴西圓盤試樣抗拉強(qiáng)度的計(jì)算結(jié)果.由表3可知，對(duì)于發(fā)生中心起裂的圓盤試樣，通過(guò)公式（4）計(jì)算的抗拉強(qiáng)度與圓盤試樣的真實(shí)抗拉強(qiáng)度比較接近，相對(duì)誤差在3%以內(nèi)（括號(hào)內(nèi)數(shù)據(jù)為相對(duì)誤差）；而對(duì)于發(fā)生端部起裂的圓盤試樣來(lái)說(shuō)，其計(jì)算抗拉強(qiáng)度明顯低于試樣的真實(shí)抗拉強(qiáng)度，并且壓拉比及加載角越小，兩者相對(duì)誤差也越大.上述分析表明只有在滿足中心起裂的基礎(chǔ)上依據(jù)峰值荷載計(jì)算出的抗拉強(qiáng)度才是有效的，而端部起裂條件下巴西劈裂實(shí)驗(yàn)會(huì)低估巖石的抗拉強(qiáng)度.

表3 不同加載角及壓拉比下二維巴西圓盤試樣抗拉強(qiáng)度計(jì)算結(jié)果Table 3 Calculated tensile strength of 2D numerical Brazilian disks with different contact loading angles and compression –tension ratios

4.2 三維模擬結(jié)果分析

圖7 為接觸加載角 2α=4.4°、14.4°、25.7°、29.7°時(shí)三維巴西圓盤試樣端面(Y=0.025)及中心剖面(Y=0.0125)受壓直徑上徑向應(yīng)力與切向應(yīng)力的分布示意圖，圖中同時(shí)給出了Hondros解析解的計(jì)算結(jié)果.由圖可知，4種加載角下中心剖面受壓直徑上徑向與切向應(yīng)力的數(shù)值模擬結(jié)果與Hondros解析解的計(jì)算結(jié)果十分吻合，而端面受壓直徑上應(yīng)力分布的數(shù)值模擬結(jié)果與Hondros解析結(jié)果存在一定誤差，這是由于三維圓盤試樣及外載關(guān)于中心剖面對(duì)稱，因此中心剖面更貼近平面應(yīng)變條件.圓盤端面受壓直徑上切向拉應(yīng)力大于中心剖面受壓直徑上的切向拉應(yīng)力，2α=25.7°、29.7°時(shí)尤其明顯.圖8為8種加載角下三維巴西圓盤中心剖面及端面受壓直徑上切向應(yīng)力分布示意圖，中心剖面受壓直徑上最大切向拉應(yīng)力始終位于剖面中心；而在端面受壓直徑上，當(dāng)加載角比較小時(shí)，最大切向拉應(yīng)力靠近圓盤頂?shù)锥耍S著加載角增大，最大切向拉應(yīng)力逐漸向端面圓心移動(dòng)，當(dāng)2α=29.7°時(shí)，最大切向拉應(yīng)力出現(xiàn)在端面圓心處.

圖7 不同加載角下三維巴西圓盤受壓直徑上應(yīng)力分布的數(shù)值模擬結(jié)果與理論計(jì)算結(jié)果對(duì)比示意圖.（a）2α=4.4°；（b）2α=14.4°；（c）2α=25.7°；（d）2α=29.7°Fig.7 Comparison of normalized stresses along the compressed diameter between the numerical results and Hondros’ solutions with different contact loading angles: (a) 2α = 4.4°; (b) 2α = 14.4°; (c) 2α = 25.7°; (d) 2α = 29.7°

圖8 不同加載角下三維巴西圓盤試樣中心剖面及端面受壓直徑上切向應(yīng)力分布示意圖.（a）中心剖面（Y=0.0125）切向應(yīng)力；（b）端面（Y=0.025）切向應(yīng)力Fig.8 Normalized tangential stresses along the compressed diameter of the surface and middle section of 3D disks with different contact loading angles:(a) middle section (Y=0.0125); (b) surface (Y=0.025)

圖9 為 2α=29.7°、25.7°、14.4°、7.6°、4.4°時(shí)巴西圓盤軸向受壓平面上(X=0)切向拉應(yīng)力分布云圖.十分明顯，圓盤內(nèi)部切向拉應(yīng)力明顯小于圓盤外部切向拉應(yīng)力，這是由于各質(zhì)點(diǎn)受到周圍質(zhì)點(diǎn)的約束程度不一樣引起的.從圓盤中心剖面(Y=0.0125)越靠近圓盤端面，質(zhì)點(diǎn)受到的水平約束也越小，水平張拉變形越大，因此張拉應(yīng)力就相對(duì)越大，所以三維圓盤試樣的最先起裂點(diǎn)出現(xiàn)在圓盤端面.而在圓盤端面上，除 2α=4.4°外，其他 7 種加載角下，起裂點(diǎn)位置就是端面受壓直徑上最大切向拉應(yīng)力出現(xiàn)的點(diǎn).本文運(yùn)用考慮張拉截?cái)嗟腗ohr-Coulomb強(qiáng)度準(zhǔn)則進(jìn)行模擬分析，當(dāng)圓盤端面受壓直徑上某一點(diǎn)的最大切向拉應(yīng)力大于試樣抗拉強(qiáng)度時(shí)就會(huì)發(fā)生張拉破壞，其本質(zhì)為最大拉應(yīng)力強(qiáng)度準(zhǔn)則.Li和Wong[7]運(yùn)用Stacey[27]張拉應(yīng)變準(zhǔn)則對(duì)三維巴西圓盤起裂點(diǎn)位置進(jìn)行分析，當(dāng)圓盤受壓直徑上某一點(diǎn)的最大張拉應(yīng)變超過(guò)某一限值時(shí)就會(huì)發(fā)生開(kāi)裂，其本質(zhì)為最大拉應(yīng)變強(qiáng)度準(zhǔn)則.喻勇[6]指出三維巴西圓盤的應(yīng)力應(yīng)變分布規(guī)律與材料泊松比有關(guān)，而與材料彈性模量無(wú)關(guān).Li和Wong[7]通過(guò)數(shù)值模擬研究發(fā)現(xiàn)，三維圓盤端面圓心處的拉應(yīng)變隨泊松比線性增大，而端面受壓直徑上的最大拉應(yīng)變隨泊松比呈現(xiàn)非線性增大趨勢(shì)，但是其未探討受壓直徑上最大張拉應(yīng)變位置的變化規(guī)律.圖10 為泊松比μ=0.15、0.20、0.25、0.30、0.35時(shí)最大切向拉應(yīng)力點(diǎn)及最大切向拉應(yīng)變點(diǎn)位置隨接觸加載角的變化情況.由圖可知，隨泊松比減小，受壓直徑上最大切向拉應(yīng)力點(diǎn)及最大切向拉應(yīng)變點(diǎn)位置均會(huì)向端面圓心移動(dòng).上述模擬結(jié)果表明，對(duì)于三維巴西圓盤，其端面起裂點(diǎn)位置還受圓盤材料的泊松比影響，較小的泊松比會(huì)促進(jìn)中心起裂.

圖9 不同加載角時(shí)三維巴西圓盤試樣軸向受壓平面上（X=0）切向應(yīng)力分布云圖.（a）2α=29.7°；（b）2α=25.7°；（c）2α=14.4°；（d）2α=7.6°；（e）2α=4.4°Fig.9 Contour plots of normalized tangential stresses (σxx) on the compressed middle section plane (X = 0) of the 3D disks with different contact loading angles: (a) 2α = 29.7°; (b) 2α = 25.7°; (c) 2α = 14.4°; (d) 2α = 7.6°; (e) 2α = 4.4°

圖10 不同泊松比時(shí)最大切向拉應(yīng)力點(diǎn)及最大切向拉應(yīng)變點(diǎn)位置隨接觸加載角的變化情況.（a）最大切向拉應(yīng)力點(diǎn)；（b）最大切向拉應(yīng)變點(diǎn)Fig.10 Change in the position of the maximum tangential tensile stress point and strain point with the contact loading angles under different Poisson’s ratios: (a) the maximum tangential tensile stress point; (b) the maximum tangential tensile strain point

表4為9種加載角下三維巴西圓盤試樣抗拉強(qiáng)度計(jì)算結(jié)果，除上述8種加載角外增加2α=33.7°這一情況，表中Distance表示起裂點(diǎn)至端面圓心的距離.當(dāng) 2α=4.4°時(shí)，圓盤試樣發(fā)生端部起裂，峰值載荷為5427 N，計(jì)算抗拉強(qiáng)度最小，而相對(duì)誤差最大為65.5%；隨著加載角增大，起裂點(diǎn)逐漸向圓盤端面圓心移動(dòng)，峰值載荷逐漸增大，計(jì)算抗拉強(qiáng)度逐漸接近試樣的真實(shí)抗拉強(qiáng)度.但是，當(dāng)2α=29.7°及33.7°時(shí)，圓盤試樣發(fā)生端面中心起裂，此時(shí)計(jì)算抗拉強(qiáng)度分別為 7.17 MPa、7.16 MPa，仍低于試樣的真實(shí)抗拉強(qiáng)度，相對(duì)誤差均為10.4%.上述計(jì)算及分析結(jié)果表明，即便在保證中心起裂的條件下，由于三維效應(yīng)巴西劈裂試驗(yàn)可能仍會(huì)低估巖石試樣的抗拉強(qiáng)度.但從工程角度而言，其結(jié)果是偏于保守的，有利于工程穩(wěn)定.

表4 不同加載角下三維巴西圓盤試樣抗拉強(qiáng)度計(jì)算結(jié)果Table 4 Calculated tensile strength of the 3D numerical Brazilian disks with different loading angles

5 討論

本文運(yùn)用Mohr-Coulomb強(qiáng)度準(zhǔn)則進(jìn)行模擬，由于圓盤破裂面上發(fā)生張拉破壞，因此其實(shí)質(zhì)為最大拉應(yīng)力強(qiáng)度準(zhǔn)則.當(dāng)圓盤受壓直徑上最大切向拉應(yīng)力位于圓盤圓心時(shí)，隨外載增加圓盤最終就會(huì)出現(xiàn)中心起裂.由于直接對(duì)公式（2）關(guān)于r求極值十分困難，因此采用數(shù)值枚舉法分析公式（2）在不同加載角下的最大值的取值位置，r/R、α取值間距分別為 0.001、0.1，結(jié)果當(dāng)α≤ 30°時(shí)公式（2）的最大值均在r/R=0即圓心處取得.這是二維數(shù)值模擬中試樣出現(xiàn)中心起裂的內(nèi)在原因.

Griffith強(qiáng)度準(zhǔn)則的實(shí)質(zhì)為當(dāng)一點(diǎn)的Griffith等效拉應(yīng)力等于其抗拉強(qiáng)度時(shí)就會(huì)發(fā)生破壞.對(duì)于線荷載作用下的平面巴西圓盤，其受壓直徑上的應(yīng)力分布為

式中：r1，r2分別為受壓直徑上一點(diǎn)到上下加載點(diǎn)處的距離，滿足r1+r2=2R.將上述兩式代入Griffith強(qiáng)度準(zhǔn)則，則Griffith等效拉應(yīng)力為

在圓心處，σG取最小值，為 2P/πDt；而在加載點(diǎn)處，σG為無(wú)窮大.因此，根據(jù)Griffith強(qiáng)度準(zhǔn)則加載點(diǎn)處最先起裂.喻勇[6]通過(guò)數(shù)值模擬研究同樣發(fā)現(xiàn)數(shù)值試樣圓心處等效拉應(yīng)力最小，而上下加載點(diǎn)處等效拉應(yīng)力最大.如果假設(shè)圓心處首先起裂，將r1=r2=R代入上式，則

這就是通常的抗拉強(qiáng)度的計(jì)算公式.對(duì)于弧形均布荷載作用下的平面巴西圓盤，受壓直徑上的Griffith等效拉應(yīng)力為

同樣采用數(shù)值枚舉法分析其在不同加載角下受壓直徑上最大等效拉應(yīng)力點(diǎn)的位置變化特征.結(jié)果如圖11所示，圖中同時(shí)給出了三維巴西圓盤端面受壓直徑上最大等效拉應(yīng)力位置隨加載角度的變化情況.對(duì)于公式（9），當(dāng) 2α≥18.2°時(shí)最大等效拉應(yīng)力點(diǎn)位于圓心處.而對(duì)于三維巴西圓盤來(lái)說(shuō)，以泊松比μ=0.20為例，由于三維效應(yīng)，當(dāng)2α≥29.7°時(shí)，最大等效拉應(yīng)力點(diǎn)才移動(dòng)到端面圓心處.當(dāng)發(fā)生中心破壞時(shí)，將r/R=0代入上式，則

圖11 最大Griffith等效拉應(yīng)力位置隨加載角度的變化情況Fig.11 Change in the position of the maximum Griffith equivalent tensile stresses with the contact loading angles

上式即為根據(jù)Griffith強(qiáng)度準(zhǔn)則推出的弧形均布荷載作用下平面巴西圓盤的抗拉強(qiáng)度計(jì)算公式.對(duì)于近似公式，當(dāng) 2α=20°相對(duì)誤差為 4.0%，而當(dāng) 2α=30°相對(duì)誤差達(dá)到 9.0%.需要指出的是，公式（4）是根據(jù)最大拉應(yīng)力強(qiáng)度準(zhǔn)則導(dǎo)出的.隨著加載角度減小，公式（4）及（10）逐漸逼近公式（8）.由于一般巖石的泊松比介于0.15～0.25，通過(guò)三維巴西圓盤端面受壓直徑上最大切向拉應(yīng)力及最大Griffith等效拉應(yīng)力位置隨加載角度的變化趨勢(shì)（圖10（a）及圖11），本文建議加載角度至少為 30°，因?yàn)榇藭r(shí)無(wú)論在最大拉應(yīng)力強(qiáng)度準(zhǔn)則還是在Griffith強(qiáng)度準(zhǔn)則下，端面初始破壞點(diǎn)均位于圓心處.

6 結(jié)論

（1）二維數(shù)值模擬結(jié)果表明，接觸加載角及壓拉比越大，巴西圓盤試樣越容易發(fā)生中心起裂；端部起裂由剪切破壞引起，而劈裂裂紋進(jìn)一步擴(kuò)展則由張拉破壞驅(qū)動(dòng)；端部起裂條件下巴西劈裂實(shí)驗(yàn)會(huì)低估巖石的抗拉強(qiáng)度.

（2）三維數(shù)值模擬結(jié)果表明，初始起裂點(diǎn)位于三維圓盤端面，隨加載角增大其沿受壓直徑逐漸向端面圓心移動(dòng)；當(dāng)圓盤發(fā)生端面中心起裂時(shí)，三維破裂面以弧形邊界向試樣內(nèi)部發(fā)散擴(kuò)展；對(duì)于三維巴西圓盤，其端面起裂點(diǎn)位置還受圓盤材料的泊松比影響，較小的泊松比會(huì)促進(jìn)中心起裂.

（3）無(wú)論圓盤試樣發(fā)生中心起裂還是端部起裂，由于三維效應(yīng)巴西劈裂試驗(yàn)可能都會(huì)低估巖石的抗拉強(qiáng)度.