王學(xué)濱,劉桐辛,田 鋒,錢(qián)帥帥
(1.遼寧工程技術(shù)大學(xué) 計(jì)算力學(xué)研究所,遼寧 阜新 123000; 2.遼寧工程技術(shù)大學(xué) 力學(xué)與工程學(xué)院,遼寧 阜新 123000)
采礦、土木等工程中巖層周期破斷、斷層周期黏滑和反復(fù)爆破都能產(chǎn)生周期沖擊載荷。周期沖擊載荷是誘發(fā)巖爆、垮塌和片幫等災(zāi)害的重要因素。在周期沖擊載荷作用下,巷道頂板和兩幫會(huì)在極短時(shí)間內(nèi)產(chǎn)生大量相互聯(lián)通的裂縫,進(jìn)而將巷道圍巖切割成眾多塊體,它們會(huì)以猛烈方式脫離巷道圍巖,造成災(zāi)害[1-4]。
輕微巖爆發(fā)生時(shí),可以聽(tīng)見(jiàn)清脆的噼啪、撕裂聲,偶有爆裂聲,巖塊彈射初速度小于1 m/s;而強(qiáng)烈?guī)r爆發(fā)生時(shí),可以聽(tīng)見(jiàn)似炸藥爆破的爆裂聲,聲響強(qiáng)烈,巖塊彈射初速度處于5~10 m/s,極易造成重大危害。目前,關(guān)于巖塊彈射初速度已有不少來(lái)自實(shí)驗(yàn)室和現(xiàn)場(chǎng)的數(shù)據(jù)積累[5-6]。例如,王之東等[5]在單軸壓縮條件下對(duì)帶方孔的3種巖石試樣的能量釋放過(guò)程進(jìn)行了觀測(cè),其中巖塊彈射初速度處于1.6~21.0 m/s;在加拿大和南非的巖爆現(xiàn)場(chǎng)[6],巖塊彈射初速度處于7.65~12.60 m/s。對(duì)巖塊彈射初速度估算的方法主要包括理論方法和數(shù)值模擬研究方法。例如,在前者方面,秦劍鋒等[7]基于巖板的屈曲失穩(wěn)理論估算的巖塊彈射初速度在9 m/s;在后者方面,陳滔等[6]基于能量守恒原理和破壞前后模型的應(yīng)變能差估算的巖塊彈射初速度處于8.16~13.60 m/s??陀^地講,巖爆過(guò)程較為復(fù)雜,受多種因素影響。在理論上想準(zhǔn)確估算巖塊彈射初速度極為困難。而且,在基于連續(xù)方法(例如,有限元法)的數(shù)值模擬中,僅能通過(guò)應(yīng)變能來(lái)估算巖塊彈射初速度,而不能模擬出巖塊從圍巖中脫離的過(guò)程,進(jìn)而難以準(zhǔn)確估算巖塊彈射初速度。
對(duì)于巷道的拉裂機(jī)理,目前已取得了一些重要進(jìn)展。例如,李夕兵等[8]采用PFC2D對(duì)動(dòng)載作用下深部巷道圍巖的動(dòng)力響應(yīng)進(jìn)行了數(shù)值模擬,發(fā)現(xiàn)巷道頂板裂紋數(shù)量增加是由靜態(tài)的拉應(yīng)力與應(yīng)力波到達(dá)頂板后反射的拉應(yīng)力波疊加所引起;ZHU等[9-10]采用Autodyn2D對(duì)爆破誘發(fā)的巷道圍巖破壞過(guò)程進(jìn)行了數(shù)值模擬,發(fā)現(xiàn)從自由邊界反射的拉應(yīng)力波引起了距離自由邊界一定距離的環(huán)向裂紋。
周期沖擊載荷3要素是幅值、頻率和持續(xù)時(shí)間。周期沖擊載荷3要素的影響研究一直受到重視[11]。前人已對(duì)周期沖擊載荷頻率的影響開(kāi)展了一定的研究。載荷頻率的取值范圍在不同文獻(xiàn)中不盡相同[12-13]。例如,閆長(zhǎng)斌等[12]采用FLAC在動(dòng)載作用下對(duì)地下巷道群的頻率影響進(jìn)行了數(shù)值模擬,頻率處于5~100 Hz;陳國(guó)祥等[13]采用FLAC對(duì)半正弦波作用下巷道圍巖的破壞過(guò)程進(jìn)行了數(shù)值模擬,頻率處于5~50 Hz?,F(xiàn)有的研究表明,載荷頻率對(duì)巖樣及巖石結(jié)構(gòu)的動(dòng)力響應(yīng)都有一定的影響[13-18]。例如,陳國(guó)祥等[13]發(fā)現(xiàn)隨著頻率的減小,巷道兩幫圍巖的最大垂直和水平應(yīng)力增加;宮鳳強(qiáng)等[16]在常規(guī)靜載和“預(yù)靜載+擾動(dòng)”條件下對(duì)巖石斷裂特性的頻率影響進(jìn)行了實(shí)驗(yàn),發(fā)現(xiàn)隨著頻率的增加,斷裂韌度呈線性減小的趨勢(shì);SU等[17-18]在真三軸加載條件下對(duì)含孔洞巖樣徑向應(yīng)力梯度的頻率影響進(jìn)行了實(shí)驗(yàn),發(fā)現(xiàn)隨著頻率的增加,巖爆更易發(fā)生。目前,該方面研究主要集中在實(shí)驗(yàn)方面和基于連續(xù)方法的數(shù)值模擬方面,而基于連續(xù)-非連續(xù)方法的研究還十分少見(jiàn)。
鑒于單元畸變和局部自適應(yīng)阻尼可能導(dǎo)致單元彈射出模型時(shí)速度失真,為了準(zhǔn)確模擬彈射現(xiàn)象,在自主開(kāi)發(fā)的拉格朗日元與離散元耦合的連續(xù)-非連續(xù)方法[19]的基礎(chǔ)上,對(duì)彈射單元進(jìn)行了剛化處理,研究了周期沖擊載荷作用下巷道圍巖的變形-開(kāi)裂-運(yùn)動(dòng)過(guò)程,還初步分析了周期載荷頻率對(duì)單元彈射初速度和剪、拉裂縫區(qū)段數(shù)目的影響。
拉格朗日元和離散元耦合的連續(xù)-非連續(xù)方法主要包括4個(gè)計(jì)算模塊[19]:應(yīng)力應(yīng)變模塊、節(jié)點(diǎn)分離模塊、接觸力求解模塊和運(yùn)動(dòng)方程求解模塊。
在應(yīng)力應(yīng)變模塊中,采用了混合離散方法,通過(guò)節(jié)點(diǎn)的速度利用高斯定理求解單元的應(yīng)力和應(yīng)變,可以避免沙漏問(wèn)題。
在節(jié)點(diǎn)分離模塊中,通過(guò)引入虛擬裂縫模型,處理應(yīng)變軟化問(wèn)題。分別選取最大拉應(yīng)力準(zhǔn)則和莫爾-庫(kù)侖準(zhǔn)則作為拉裂和剪裂判據(jù)。當(dāng)節(jié)點(diǎn)分離后,虛擬裂縫產(chǎn)生,虛擬裂縫面之間存在法向及切向黏聚力。通過(guò)同時(shí)引入Ⅰ型和Ⅱ型斷裂能計(jì)算法向及切向黏聚力。當(dāng)法向或切向黏聚力降至0時(shí),虛擬裂縫成為真實(shí)裂縫。
在接觸力求解模塊中,采用了基于空間劃分的接觸檢測(cè)方法和基于勢(shì)的接觸力求解方法,具有接觸檢測(cè)效率較高、無(wú)需對(duì)“角-角”接觸問(wèn)題進(jìn)行特殊處理的優(yōu)勢(shì)。
在運(yùn)動(dòng)方程求解模塊中,采用中心差分方法求解節(jié)點(diǎn)的速度和位移,具有計(jì)算效率高、計(jì)算精度較高的優(yōu)勢(shì)。
網(wǎng)格法(有限元法、有限差分法等)容易出現(xiàn)網(wǎng)格或單元畸變問(wèn)題。本文的連續(xù)-非連續(xù)方法屬于網(wǎng)格法。當(dāng)某單元彈射出模型后,該單元的4個(gè)節(jié)點(diǎn)的速度一般并不相同,而且,各節(jié)點(diǎn)的運(yùn)動(dòng)是相互獨(dú)立的,所以,該單元的變形、運(yùn)動(dòng)規(guī)律將極為復(fù)雜,從而可能導(dǎo)致該單元發(fā)生畸變,這會(huì)使應(yīng)力、應(yīng)變和由應(yīng)力引起的彈性力等計(jì)算錯(cuò)誤,從而會(huì)進(jìn)一步加劇該單元畸變。當(dāng)某單元彈射出模型后,在與其他單元碰撞之前和之后,該單元應(yīng)僅受重力作用,而在水平方向不受力,這樣,該單元的各節(jié)點(diǎn)水平速度vx應(yīng)是常量。一旦該單元發(fā)生畸變,將不能保證這一點(diǎn)。另外,本文方法使用局部自適應(yīng)阻尼(與常見(jiàn)的黏性阻尼不同),由于其方向取決于節(jié)點(diǎn)速度的方向,而大小取決于不平衡力(由應(yīng)力引起的彈性力和重力等構(gòu)成)的大小,這也可能導(dǎo)致vx不是常量。
為此,筆者提出一種處理方法以確保彈射單元各節(jié)點(diǎn)速度的計(jì)算結(jié)果正常,即在與其他單元碰撞之前和之后,保持恒定。
對(duì)于彈射單元,首先,需找到應(yīng)力狀態(tài)較為接近于自由狀態(tài)的臨界時(shí)步數(shù)目??紤]到該單元?jiǎng)偙粡椛涑鰰r(shí)常處于壓縮狀態(tài),因而其剛進(jìn)入完全拉伸狀態(tài)時(shí)的應(yīng)力狀態(tài)較為接近于自由狀態(tài)。為此,將該單元最小主應(yīng)力σ1>0(在本文中,σ3≥σ1,σ3為最大主應(yīng)力)對(duì)應(yīng)的時(shí)步數(shù)目作為臨界時(shí)步數(shù)目。然后,將該單元4個(gè)節(jié)點(diǎn)的速度矢量取平均,獲得平均速度矢量,并賦給這些節(jié)點(diǎn),并對(duì)該單元的應(yīng)力清零,即清除彈性力。各節(jié)點(diǎn)的新速度v0取為
(1)
式中,vi(i=1~4)為該單元各節(jié)點(diǎn)的原速度。
上述處理迫使彈射單元在碰撞之前和之后作剛體平動(dòng)。由于彈射單元的應(yīng)力已被清除,單元畸變不會(huì)進(jìn)一步發(fā)展。而且,由于各節(jié)點(diǎn)只受重力(或?yàn)椴黄胶饬Φ奈ㄒ怀煞?和局部自適應(yīng)阻尼力(僅位于垂直方向上,其大小取決于重力的大小,其方向取決于節(jié)點(diǎn)垂直方向速度vy),各節(jié)點(diǎn)在水平方向上將作勻速直線運(yùn)動(dòng)。這樣,即可避免節(jié)點(diǎn)速度的計(jì)算過(guò)于失真,從而確保巖塊彈射現(xiàn)象的模擬結(jié)果盡量真實(shí)。
本文的計(jì)算模型(以下簡(jiǎn)稱(chēng)模型)建立依據(jù)某大巷[20]所處地質(zhì)條件,巖性中硬。
模型尺寸為40 m×40 m,被剖分成160×160個(gè)正方形單元。壓應(yīng)力波施加后的力學(xué)模型如圖1(a)所示。應(yīng)當(dāng)指出,模型的左、右側(cè)面均為透射邊界,即應(yīng)力波經(jīng)過(guò)時(shí)不會(huì)發(fā)生反射。
在巷道頂板,布置了5個(gè)測(cè)點(diǎn),在巷道左幫,布置了6個(gè)監(jiān)測(cè)單元,具體如圖1(b)所示。計(jì)算在平面應(yīng)變、大變形條件下進(jìn)行。
圖1 力學(xué)模型、測(cè)點(diǎn)和監(jiān)測(cè)單元Fig.1 Mechanical model,monitored points and monitored elements
陳建君等[21]將沖擊載荷簡(jiǎn)化為半正弦波,采用的半正弦壓應(yīng)力波P(N)表達(dá)式為
(2)
式中,Pmax為壓應(yīng)力波幅值,取7.3 MPa;頻率f=ω/π;ω為圓頻率;N為時(shí)步數(shù)目。
計(jì)算可分為3個(gè)過(guò)程:
(1)對(duì)開(kāi)挖前的模型進(jìn)行計(jì)算,直到模型處于靜力平衡狀態(tài)(所用N=12 000)。
(2)開(kāi)挖尺寸為6 m×6 m的巷道(所用N=4 000);此后,繼續(xù)計(jì)算,直到模型處于靜力平衡狀態(tài)(所用N=4 000)。
(3)當(dāng)N=20 000時(shí),在模型的上端面施加豎直向下的周期沖擊載荷(半正弦壓應(yīng)力波),此后,繼續(xù)計(jì)算;當(dāng)N=64 000時(shí)最后1個(gè)壓應(yīng)力波波尾傳入模型,此后一段時(shí)間之內(nèi)計(jì)算仍在繼續(xù)。在最后1個(gè)壓應(yīng)力波在圍巖中逐漸消失的過(guò)程中,由于圍巖應(yīng)力場(chǎng)的略微調(diào)整,一些裂縫仍有可能進(jìn)一步發(fā)展,而且,脫離圍巖的單元仍在運(yùn)動(dòng)。因此,多計(jì)算一段時(shí)間是必要的。
共采用4個(gè)計(jì)算方案。方案1~4的f分別為15,25,35及45 Hz。文獻(xiàn)[22]通過(guò)現(xiàn)場(chǎng)觀測(cè)發(fā)現(xiàn),震動(dòng)波的頻率主要集中于50 Hz以內(nèi)。本文選取的f涵蓋了上述范圍。
3.1.1多個(gè)壓應(yīng)力波沖擊下剪、拉裂縫的時(shí)空分布
以方案1為例簡(jiǎn)單分析多個(gè)壓應(yīng)力波沖擊下剪、拉裂縫的時(shí)空分布。
圖2為方案1的剪裂縫與最大主應(yīng)力σ3的時(shí)空分布規(guī)律,黑色線段代表剪裂縫區(qū)段。圖3為方案1的拉裂縫與最大主應(yīng)力σ3的時(shí)空分布規(guī)律,黑色線段代表拉裂縫區(qū)段。由圖2,3可以發(fā)現(xiàn),首先,第1個(gè)壓應(yīng)力波傳至巷道兩幫后,巷道兩幫產(chǎn)生剪裂縫,并逐漸發(fā)展形成V形坑,其內(nèi)產(chǎn)生少量拉裂縫(圖2(a),(b)和圖3(a),(b))。然后,后續(xù)壓應(yīng)力波陸續(xù)傳入模型,巷道兩幫既有V形坑外若干新的V形坑形成,其內(nèi)仍產(chǎn)生少量拉裂縫,與此同時(shí),巷道頂板產(chǎn)生拉裂縫,并發(fā)展形成層裂結(jié)構(gòu),巷道兩幫脫離圍巖的一些單元彈入巷道(圖2(b),(c)和圖3(b),(c))。文獻(xiàn)[1]采用真三軸試驗(yàn)機(jī)對(duì)含預(yù)制矩形巷道的立方體巖樣進(jìn)行壓縮實(shí)驗(yàn),發(fā)現(xiàn)巷道兩幫出現(xiàn)巖塊彈射現(xiàn)象,巖塊明顯堆積于巷道底板;文獻(xiàn)[23]利用水泥類(lèi)膨脹膠凝材料與水反應(yīng)體積驟增的特點(diǎn)對(duì)巷道圍巖進(jìn)行單次沖擊實(shí)驗(yàn),發(fā)現(xiàn)巷道頂板和兩幫均發(fā)生破壞,且兩幫存在V形坑。上述現(xiàn)象的條件盡管與本文有所不同,但上述現(xiàn)象與本文的結(jié)果(圖2(b),(c)和圖3(b),(c))基本一致。最終,最后1個(gè)壓應(yīng)力波傳出模型后,即在圍巖中消失后,剪、拉裂縫停止發(fā)展。本文方法在處理開(kāi)裂、接觸和摩擦等非線性問(wèn)題時(shí)不可避免存在微小誤差,隨著計(jì)算的進(jìn)行,這種誤差可能被放大,這會(huì)導(dǎo)致剪、拉裂縫的分布不具有嚴(yán)格的對(duì)稱(chēng)性(例如,圖2(c),(d)和圖3(c),(d))。這種現(xiàn)象是非線性數(shù)值模擬方法的共性。應(yīng)當(dāng)指出,2個(gè)單元之間的裂縫稱(chēng)之為1個(gè)裂縫區(qū)段,裂縫區(qū)段的形狀為四邊形。若干裂縫區(qū)段連在一起構(gòu)成裂縫??紤]到單元脫離圍巖后裂縫將變得很大,圖2~3僅顯示了各邊長(zhǎng)度均不大于1個(gè)單元邊長(zhǎng)的裂縫區(qū)段。
3.1.2巷道左幫監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx演化
圖4為方案1的1~6號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx-N曲線。由圖4可以發(fā)現(xiàn),在壓應(yīng)力波傳入模型之后,隨著N的增加,2~6號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx(對(duì)于任一監(jiān)測(cè)單元,彈入巷道后其4個(gè)節(jié)點(diǎn)vx均相同)呈現(xiàn)上升—穩(wěn)定—衰減的變化過(guò)程,而1號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx呈現(xiàn)波動(dòng)—穩(wěn)定的變化過(guò)程。應(yīng)當(dāng)指出,當(dāng)N=20 000時(shí)第1個(gè)壓應(yīng)力波開(kāi)始傳入模型,當(dāng)N=64 000時(shí)最后1個(gè)壓應(yīng)力波波尾傳入模型。
圖4 方案1的1~6號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx-N曲線Fig.4 Evolution of horizontal velocities of lower right corner nodes of monitored elements 1-6 with time steps of Scheme 1
在vx有上升趨勢(shì)的階段,2~6號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx分別上升至最大值,此階段位于壓應(yīng)力波傳入模型之后,vx上升的過(guò)程是監(jiān)測(cè)單元逐漸脫離圍巖的過(guò)程。例如,當(dāng)N=28 120~31 990時(shí),2號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx由0.396 m/s增至峰值6.8 m/s,此時(shí),2號(hào)監(jiān)測(cè)單元脫離巷道圍巖(圖2(b))。期間,vx的最大值就是彈射初速度,分別為6.8,5.8,9.7,9.5和10.3 m/s,上述彈射初速度的平均值,即平均彈射初速度,為8.43 m/s。這位于強(qiáng)烈?guī)r爆[24]發(fā)生時(shí)平均巖塊彈射初速度范圍(5.0~10.0 m/s)之內(nèi)。
在vx穩(wěn)定的階段,2~6號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx保持不變。期間,這些監(jiān)測(cè)單元未與其他單元發(fā)生碰撞,即在水平方向上不受力,因而沒(méi)有產(chǎn)生能量損耗。由于巷道的空間有限,這些監(jiān)測(cè)單元一定會(huì)與其他單元發(fā)生碰撞,因此這些監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx無(wú)法總保持不變。這些監(jiān)測(cè)單元的位置和彈射初速度的不同導(dǎo)致vx發(fā)生變化的時(shí)刻不同,有的位于壓應(yīng)力波傳出模型之前,而有的位于壓應(yīng)力波傳出模型之后。
在vx有衰減趨勢(shì)的階段,2~6號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx呈現(xiàn)總體衰減的趨勢(shì)。此階段存在vx反向現(xiàn)象。首先,當(dāng)N=74 661~74 760時(shí),2號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx由6.8 m/s降至-5.16 m/s。這是因?yàn)樵摫O(jiān)測(cè)單元與其他單元發(fā)生碰撞導(dǎo)致反向并產(chǎn)生能量損耗(圖2(c));然后,該監(jiān)測(cè)單元不斷與其他單元發(fā)生碰撞,導(dǎo)致其右下角節(jié)點(diǎn)vx發(fā)生多次反向現(xiàn)象,并呈衰減趨勢(shì)(圖2(d))。
圖5為剛化處理前后2號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx-N曲線。為了呈現(xiàn)對(duì)彈射單元進(jìn)行剛化處理前后的差異,對(duì)方案1重新進(jìn)行了計(jì)算。由此可以發(fā)現(xiàn),剛化處理前2號(hào)監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx雜亂無(wú)章,這并不是該單元與其他單元碰撞的結(jié)果,而剛化處理后2號(hào)監(jiān)測(cè)單元在未與其他單元碰撞時(shí),可以保持恒定的vx,這較為符合該單元在水平方向上不受力的事實(shí),這在一定程度上說(shuō)明了剛化處理的正確性。
圖5 剛化處理前后2號(hào)監(jiān)測(cè)單元的 右下角節(jié)點(diǎn)vx-N曲線Fig.5 Evolution of horizontal velocities of lower right corner nodes of monitored elements 2 with time steps before and after rigid treatments
3.1.3巷道頂板測(cè)點(diǎn)σ3的演化規(guī)律及開(kāi)裂機(jī)理
以方案1為例,簡(jiǎn)單分析巷道頂板各測(cè)點(diǎn)的σ3隨N的演化規(guī)律。圖6為方案1的1~5號(hào)測(cè)點(diǎn)的σ3-N曲線。由圖6可以發(fā)現(xiàn),在壓應(yīng)力波傳入模型(N=20 000)之前,首先,各測(cè)點(diǎn)的σ3穩(wěn)定在-27 MPa,這對(duì)應(yīng)于巷道開(kāi)挖之前模型的靜力平衡狀態(tài);隨后,受巷道開(kāi)挖的影響,各測(cè)點(diǎn)的σ3先經(jīng)歷震蕩上升、后逐漸衰減至穩(wěn)定的變化過(guò)程。在壓應(yīng)力波傳入模型之后,大部分測(cè)點(diǎn)的σ3呈現(xiàn)近似正弦波動(dòng)上升—衰減—穩(wěn)定的變化過(guò)程。
圖6 方案1的1~5號(hào)測(cè)點(diǎn)的σ3-N曲線Fig.6 Evolution of σ3 of monitored points 1-5 with time steps of Scheme 1
下面以4號(hào)測(cè)點(diǎn)為例,詳細(xì)分析近似正弦波動(dòng)上升階段σ3的演化規(guī)律。在此階段中,在σ3-N曲線上,可隱約觀察到13次波動(dòng),這是因?yàn)榉桨?中施加了13個(gè)壓應(yīng)力波,各壓應(yīng)力波均會(huì)對(duì)該測(cè)點(diǎn)的σ3產(chǎn)生影響。圖7為方案1的4號(hào)測(cè)點(diǎn)的σ3-N曲線。由圖7可以發(fā)現(xiàn):
圖7 方案1的4號(hào)測(cè)點(diǎn)的σ3-N曲線Fig.7 Evolution of σ3 of the monitored point 4 with time steps of Scheme 1
(1)σ3首先呈現(xiàn)有規(guī)律的波動(dòng),然后波動(dòng)幅度突然增大,并伴隨著劇烈震蕩。σ3-N曲線的波峰和波谷均有隨著N的增加而增加的趨勢(shì)。
當(dāng)N=20 460~23 760時(shí),σ3出現(xiàn)第1次波動(dòng)。這表明,第1個(gè)壓應(yīng)力波傳至并逐漸傳過(guò)4號(hào)測(cè)點(diǎn)。
當(dāng)N=20 460~22 130時(shí),σ3由-15.8 MPa減小至-30.94 MPa。這表明第1個(gè)壓應(yīng)力波的峰前部分逐漸傳至該測(cè)點(diǎn)。應(yīng)當(dāng)指出,當(dāng)N=20 800~20 970時(shí),σ3有小幅度震蕩。隨后,隨著N的增加,σ3繼續(xù)下降,σ3-N曲線斜率的絕對(duì)值明顯小于N=20 800之前的。也就是說(shuō),σ3小幅度震蕩前后σ3-N曲線的斜率有所不同,前面曲線更陡,而后面更平緩。σ3小幅震蕩的原因是第1個(gè)壓應(yīng)力波的前緣經(jīng)由巷道頂板表面反射的拉應(yīng)力波傳至該測(cè)點(diǎn)。此后,第1個(gè)壓應(yīng)力波的后續(xù)部分繼續(xù)傳過(guò)該測(cè)點(diǎn),并與反射的拉應(yīng)力波發(fā)生疊加,致使σ3-N曲線的斜率改變。
當(dāng)N=22 131~22 480時(shí),σ3-N曲線出現(xiàn)第1個(gè)波谷,這是因?yàn)榈?個(gè)壓應(yīng)力波的峰值部分剛傳至該測(cè)點(diǎn),與此同時(shí),該測(cè)點(diǎn)仍位于第1個(gè)壓應(yīng)力波反射的拉應(yīng)力波中。
當(dāng)N=22 481~23 760時(shí),第1個(gè)壓應(yīng)力波的峰后部分逐漸傳至該測(cè)點(diǎn),σ3由-30.04 MPa增大至-11 MPa。當(dāng)N=23 760時(shí),第1個(gè)壓應(yīng)力波已完全傳過(guò)該測(cè)點(diǎn),與此同時(shí),上述拉應(yīng)力波尚未完全傳出該測(cè)點(diǎn),σ3-N曲線出現(xiàn)第1個(gè)波峰。第2個(gè)壓應(yīng)力波立即傳至該測(cè)點(diǎn),導(dǎo)致σ3下降。應(yīng)當(dāng)指出,與4號(hào)測(cè)點(diǎn)剛受第1個(gè)壓應(yīng)力波影響時(shí)(N=20 460)的σ3相比,第1個(gè)壓應(yīng)力波剛傳出該測(cè)點(diǎn)時(shí)(N=23 760)的σ3更大,這是因?yàn)榈?個(gè)壓應(yīng)力波反射的拉應(yīng)力波的作用。
當(dāng)N=20 460~46 889時(shí),σ3規(guī)律性波動(dòng),期間,波峰和波谷隨著N的增加均有增加趨勢(shì)。波峰變化的原因同前所述。波谷變化的原因?yàn)椋悍瓷涞亩鄠€(gè)拉應(yīng)力波會(huì)對(duì)傳至該測(cè)點(diǎn)的壓應(yīng)力波有一定的衰減作用,致使壓應(yīng)力波的作用不再?gòu)?qiáng)烈。
當(dāng)N=46 890~47 040時(shí),與此前幾次有規(guī)律的波動(dòng)相比,σ3的波動(dòng)幅度突然增大,震蕩加劇,這與該測(cè)點(diǎn)下方節(jié)點(diǎn)的分離(介質(zhì)的拉裂)有關(guān)。該測(cè)點(diǎn)下方節(jié)點(diǎn)的分離將使σ3向圍巖的深部轉(zhuǎn)移,從而將提升該節(jié)點(diǎn)的σ3。此外,該測(cè)點(diǎn)下方節(jié)點(diǎn)的分離將對(duì)應(yīng)力波的已有傳播產(chǎn)生影響,因而該測(cè)點(diǎn)σ3的震蕩將加劇。
(2)節(jié)點(diǎn)發(fā)生分離后,σ3的震蕩幅度有衰減趨勢(shì)。
當(dāng)N=53 190時(shí),4號(hào)測(cè)點(diǎn)的σ3達(dá)到σt(5 MPa),頂板在此處將發(fā)生拉裂。下面,闡明頂板的拉裂機(jī)理。若σ3的波峰不發(fā)生改變,則頂板不可能被拉裂。σ3的波峰隨著N的增加而增加使頂板拉裂成為可能。當(dāng)σ3極小時(shí),波峰的σ3為負(fù),代表壓縮;當(dāng)σ3極大時(shí),波峰的σ3為正,代表拉伸。前一個(gè)壓應(yīng)力波剛完全傳過(guò)該測(cè)點(diǎn)之時(shí),恰是下一個(gè)壓應(yīng)力波剛傳至該測(cè)點(diǎn)之時(shí)。下一個(gè)壓應(yīng)力波剛傳入該測(cè)點(diǎn)時(shí)的σ3總比前一個(gè)壓應(yīng)力波的大,這說(shuō)明該測(cè)點(diǎn)的最大σ3伴隨著應(yīng)力波的不斷傳入和傳出而越來(lái)越高,在不斷累積,這種累積只能源于反射的多個(gè)拉應(yīng)力波的作用,即反射的拉應(yīng)力波每通過(guò)一次該測(cè)點(diǎn)將其σ3提升一次,直至達(dá)到σt。
此后,σ3的震蕩規(guī)律變得復(fù)雜。與此同時(shí),震蕩幅度有衰減的趨勢(shì),直至穩(wěn)定在零值附近。水平展布的裂縫會(huì)阻斷應(yīng)力波的傳播,并引起應(yīng)力波的反射,從而引起該測(cè)點(diǎn)σ3的劇烈震蕩。經(jīng)由巷道頂板表面反射的拉應(yīng)力波產(chǎn)生的拉伸應(yīng)力應(yīng)在接近垂直方向上。所以,頂板將產(chǎn)生水平展布的拉伸裂縫。與此同時(shí),若干傳入模型的壓應(yīng)力波相當(dāng)于增加了模型所受的垂直方向應(yīng)力。由此,頂板將產(chǎn)生垂直方向展布的拉伸裂縫。
此外,由圖6還可以發(fā)現(xiàn),從整體上看,在近似正弦波動(dòng)上升階段,隨著N的增加,大部分測(cè)點(diǎn)的σ3波動(dòng)上升。越靠近巷道頂板表面,σ3的波動(dòng)幅度越小,其中,1號(hào)測(cè)點(diǎn)的最小,5號(hào)測(cè)點(diǎn)的最大。下面對(duì)其原因進(jìn)行解釋。當(dāng)?shù)?個(gè)壓應(yīng)力波的前緣傳至1號(hào)測(cè)點(diǎn)不久,在巷道頂板表面發(fā)生反射;隨后,由于1號(hào)測(cè)點(diǎn)離巷道頂板表面很近,反射的拉應(yīng)力波立即傳至1號(hào)測(cè)點(diǎn),由于拉應(yīng)力波σ3的絕對(duì)值與此時(shí)傳至1號(hào)測(cè)點(diǎn)的壓應(yīng)力波的σ3相差應(yīng)不大,二者疊加將造成壓應(yīng)力波的能量極大衰減,即對(duì)1號(hào)測(cè)點(diǎn)的影響大幅降低。當(dāng)壓應(yīng)力波的前緣經(jīng)由巷道頂板表面反射的拉應(yīng)力波傳至5號(hào)監(jiān)測(cè)節(jié)點(diǎn)時(shí),由于該拉應(yīng)力波在有阻尼的介質(zhì)中已傳播了一定距離,其能量將有所衰減,傳至該測(cè)點(diǎn)的壓應(yīng)力波的σ3的絕對(duì)值將大于該拉應(yīng)力波的,即該拉應(yīng)力波對(duì)5號(hào)測(cè)點(diǎn)的影響將不如對(duì)1號(hào)測(cè)點(diǎn)的大。
3.1.4頻率的影響
(1)對(duì)剪、拉裂縫區(qū)段數(shù)目的影響。
圖8,9分別為方案1~4的剪裂縫區(qū)段數(shù)目Ns和拉裂縫區(qū)段數(shù)目Nt隨N的演變規(guī)律,統(tǒng)計(jì)的Ns和Nt包括圖2~3中顯示與否的裂縫區(qū)段。
由圖8可以發(fā)現(xiàn),各方案的Ns均呈現(xiàn)恒為0—近似階梯增長(zhǎng)—恒定的變化趨勢(shì)。只標(biāo)注了方案1的Ns-N曲線和Nt-N曲線。由圖9可見(jiàn),各方案的Nt均呈現(xiàn)恒為0—近似階梯增長(zhǎng)—恒定的變化趨勢(shì)。
圖8 方案1~4的Ns-N曲線Fig.8 Evolution of the number of shear crack segments with time steps of Schemes 1-4
圖9 方案1~4的Nt-N曲線Fig.9 Evolution of the number of tensile crack segments with time steps of Schemes 1-4
在恒定為0階段,第1個(gè)壓應(yīng)力波的前緣尚未傳至巷道頂板表面。
在近似階梯增長(zhǎng)階段,Nt時(shí)增時(shí)穩(wěn)。以方案1(f=15 Hz)為例。當(dāng)N=23 200~27 520時(shí),Nt-N曲線存在若干個(gè)增長(zhǎng)階梯,期間,巷道兩幫V形坑內(nèi)產(chǎn)生拉裂縫。當(dāng)N=27 521~46 300時(shí),Nt穩(wěn)定在280附近,期間,V形坑內(nèi)拉裂縫幾乎停止發(fā)展,同時(shí),巷道頂板尚未產(chǎn)生拉裂縫。在N=46 301之后,Nt快速增長(zhǎng),期間,反射的多個(gè)拉應(yīng)力波的累積作用不斷提升頂板某一位置最大主應(yīng)力的波峰,巷道頂板產(chǎn)生拉裂縫。
在恒定階段,方案1~4的Nt分別最終穩(wěn)定在2 681,1 110,741和546。期間,模型中不再有壓應(yīng)力波傳入,拉裂縫停止發(fā)展。由此可見(jiàn),隨著f的增加,Nt的最終穩(wěn)定值減小。
綜上所述,方案1~4的Ns和Nt的最終穩(wěn)定值均隨著f的增加而減小。已有研究表明,應(yīng)力波的頻率越高,在介質(zhì)中傳播衰減程度越大[25]。在本文模型中,采用局部自適應(yīng)阻尼,應(yīng)力波在其中傳播,能量自然也會(huì)衰減。上述Ns和Nt的最終穩(wěn)定值依賴于f的數(shù)值結(jié)果可以在一定程度上得到解釋。
(2)對(duì)平均彈射初速度的影響。
由圖4,10可以計(jì)算得到,方案1~4的平均彈射初速度分別為8.43,8.56,10.14和11.07 m/s。由此可見(jiàn),隨著f的增加,平均彈射初速度增加。應(yīng)當(dāng)指出,方案1~4傳入模型的壓應(yīng)力波分別為13,22,30和33個(gè)。傳入的壓應(yīng)力波數(shù)量越多,則巷道圍巖對(duì)即將彈射單元的推動(dòng)作用越頻繁,這會(huì)導(dǎo)致平均彈射初速度越大。
圖10 方案2~4的監(jiān)測(cè)單元的右下角節(jié)點(diǎn)vx-N曲線Fig.10 Evolution of horizontal velocities of lower right corner nodes of monitored elements with time steps of Schemes 2-4
(1)鑒于單元畸變和局部自適應(yīng)阻尼可能導(dǎo)致單元彈射出模型時(shí)速度失真,為了準(zhǔn)確模擬彈射現(xiàn)象,在原始方法的基礎(chǔ)上對(duì)彈射單元進(jìn)行了剛化處理,確保了脫離模型的單元在碰撞前后作剛體運(yùn)動(dòng)。結(jié)果表明,彈射單元任一節(jié)點(diǎn)的水平速度呈現(xiàn)上升—穩(wěn)定—衰減的變化趨勢(shì),這與不進(jìn)行剛化處理的結(jié)果相比更符合實(shí)際。
(2)在壓應(yīng)力波傳入模型之后,巷道頂板左、右對(duì)稱(chēng)線上大部分測(cè)點(diǎn)的最大主應(yīng)力呈現(xiàn)近似正弦波動(dòng)上升—衰減—穩(wěn)定的變化過(guò)程。在某一節(jié)點(diǎn)分離前,最大主應(yīng)力首先呈現(xiàn)有規(guī)律的波動(dòng);然后,波動(dòng)幅度突然增大,并伴隨著劇烈震蕩。
(3)當(dāng)應(yīng)力波傳入模型之后,頂板某一位置的最大主應(yīng)力的波峰和波谷隨著時(shí)間的增加而有增加的趨勢(shì),反射的多個(gè)拉應(yīng)力波的累積作用不斷提升該位置最大主應(yīng)力的波峰,致使拉裂,此即為巷道頂板拉裂機(jī)理。
(4)周期載荷頻率對(duì)單元彈射初速度和剪、拉裂縫區(qū)段數(shù)目有一定影響。隨著頻率的增加,脫離巷道圍巖的單元的平均彈射初速度增加,剪、拉裂縫區(qū)段數(shù)目減小。