低幅值應力波在巖石三維節理面的折反射規律理論研究

劉傳正, 張建經, 崔 鵬

(1.中國科學院 山地災害與地表過程重點實驗室 成都山地災害與環境研究所,成都 610041；2.中國科學院大學 工程科學學院，北京 100049；3.西南交通大學 土木工程學院，成都 610031)

巖體中的結構面很大程度上影響了巖體的力學性質。近些年來，由于地震事故頻發，并引發了大量的崩塌、滑坡災害，地震作用下巖土體失穩機理成為當前十分迫切的研究任務之一。而巖石節理對應力波傳播的影響規律對于研究地震巖質邊坡動力響應與失穩機理具有重要意義。

巖石節理力學性質的研究始終貫穿于巖石力學中的研究中，是巖石工程常常所面對的不可忽視的客觀存在。巖石節理力學研究早期是對二維粗糙節理力學進行研究，而對于三維粗糙節理的研究則相對較為有限，但是仍有大量的研究發現，如Belem等[1]通過對節理面幾何統計分析發現，對于自然形成節理，巖石節理面凸起的傾角和曲率以及節理的粗糙度與曲折度，都會表現出一定的各向異性性質。孫輔庭等[2]運用分形數學理論，提出了一種基于三維均方根抵抗角的節理面粗糙度分析方法，分形參數在方位角上表現出明顯的各向異性性質。Jing等[3]通過巖石節理面剪切實驗發現，在各方向上節理起伏角并不是均勻分布的，且節理的剪切變形剛度系數在各個方向上也不相等；進一步研究發現巖石節理起伏角的各向分布與剪切剛度在方向上分布呈正相關，且具有中心對稱的兩個優勢方向，可用橢圓描述。Jing等[4]在巖石節理表面形態各向異性研究的基礎上，給出了巖石節理三維本構模型。類似的，王光綸等[5]給出了具有非線性彈性變形的巖石節理三維本構方程。而Misra[6]則基于三維彈塑性接觸理論、節理面形態的幾何參數分布函數給出了積分形式的巖石節理三維本構力學模型。以上研究均表明中皆從巖石三維節理面的形態作為主要出發點，有些給出三維節理本構方程，然后都過去研究中可以發現三維節理變形本構方程的剛度矩陣中交叉耦合項非零；由于節理變形方程中交叉耦合剛度系數是描述巖石節理力學性質重要參數，如描述了節理的剪脹、減縮現象，而且節理本構方程中交叉耦合剛度系數不僅僅描述這兩種現象，其對應力波在節理上折反射規律的影響相關研究也十分缺乏。

相比巖石節理力學性質的研究，過去在研究巖石節理應力波傳播規律中對節理本構模型則做出了較大簡化，其很大原因是由于問題求解的復雜性。Schoenberg[7]基于位移非連續模型，給出應力波在線剛度非連續面上的折反射解析解，隨后很多學者對此方法進行了研究和補充，如Myer等[8-15]通過理論、實驗對線剛度DDM(Displacement Discontinue Method)模型進行大量研究，證明了該理論的可靠性和實用性。Zhao等[16]給出了縱波垂直入射非線性變形節理面的解析解；俞縉等[17-18]對非線性節理面的應力波傳播模型進行了改進。Daehnke等[19]給出考慮交叉耦合項的二維節理面DDM模型解析解。Misra等[20]在給出的三維節理本構模型的基礎上，分析了應力波的折反射行為，但是其理論框架仍是在傳統的認識范圍內。綜上述，以往研究沒有對三維粗糙節理的應力波折反射行為進行探討，沒有全面的分析三維節理本構方程中的交叉耦合剛度系數項對三種簡諧波傳播的影響。

本文基于巖石節理的細觀接觸模型，討論了巖石節理一般形式的三維應力-應變本構方程，并以此為基礎分析了簡諧波在節理面上的相互轉化機理，并據此推導了三維節理上應力波折反射的理論解，最后重點討論了節理本構方程中的交叉耦合項對應力波折反射系數的影響。

1 巖石節理的形變剛度討論

1.1 三維粗糙節理接觸變形機理

(1)

圖1 節理凸起接觸示意圖

(2)

(3)

(4)

1.2 三維節理的本構方程

由接觸面的本構方程，即可以得出節理的整體本構方程，但是需知節理上接觸點在空間的分布情況。如圖2所示，設節理單位面積內的接觸點的分布密度為ρc(與節理上下壁面形態、應力水平有關)，而接觸點在高程上的分布密度函數為H(z)，則ρcH(z)dz則是指單位面積節理上的高程z～z+dz的所有接觸點數目。由此可以得出，此時單位面積節理上接觸點的數目為

圖2 節理表面高程分布

(5)

而對于節理上接觸點的分布，應根據具體的情況選用適當的分布函數來描述，例如：卡方分布χ2、伽馬分布Γ(α,β),F分布等。

圖3 接觸面與作用力

(6)

(7)

式中：Δz=[zmin,zmax]。將節理凸起的接觸變形方程式(3)代入式(7)，得

(8)

(9)

這里剛度系數Kij為此時節理變形的剛度系數。由式(9)可以看出節理的變形剛度系數Kij與節理的形態特征、應力狀態密切關。由于在巖石節理的變形過程中，節理上下兩壁之間的接觸面的分布、形態等是動態變化的，因此變形剛度矩陣Kij是動態變化的。三維節理的變形剛度矩陣為3×3的矩陣，非對角線上的剛度系數一般非零-由式(9)可知，在節理為平直節理時，非對角線上的剛度系數為零；或者對于傾角在傾向上為對稱分布，且當節理當前應力狀態為0時，或僅法向應力存在時節理變形矩陣中的非對角線上的剛度系數為零；對于其他情形，需要根據節理面的形態進行具體的分析。然而，這里并不給直接給出節理本構方程的具體形式，且對于動力問題，應力波在三維非線性本構節理上折反射求解異常復雜，不利于理論解析。然而對于低幅值應力波入射時，應力波動幅值相對節理的強度較小時，可以認為此時節理的變形為線剛度模型，因此這里對巖石節理本構方程為三維線剛度模型時的應力波折反射規律進行分析。

2 三維節理的應力波折反射理論分析

令節理的平均面為XOZ坐標平面，法向為Y軸；設應力波入射線所在平面為XOY平面，則P波、SV波的偏振方向在XOY平面內，而SH波的偏振方向為Z軸方向。節理上下兩壁巖石塊體A，B為線彈性介質。由前面的討論可知，當低幅值應力波入射時，節理本構方程可以采用的三維線剛度本構方程dFi=Kijduj來表示。由線性問題的可疊加性可知，此時在應力波入射與折、反射過程中，節理本構方程可以改寫為

(10)

式中:[δσy,δτyx,δτyz]T為節理上的應力變化分量，Pa；[δuy,δux,δuz]T為節理上下壁面的相對位移變化分量，m；kij，i,j=x,y,z為節理的變形剛度系數，Pa/m。且在本文中將處于矩陣K對角線上的kxx,kyy,kzz稱為主剛度系數，而位于對角線之外剛度系數稱為交叉耦合剛度系數。根據熱力學第二定律可以證明，節理剛度系數矩陣Kij必為實對稱陣，即kij=kji。由波動理論可知，沿同一方向傳播的應力波，P波、SV波與SH波之間的偏振方向是相互正交的，如圖4所示當P波垂直入射時，節理面會產生擾動δuy,δσy；SV波垂直入射時，會產生擾動δux,δτxy；而P波、SV波傾斜入射時，都會在節理面產生擾動δuy,δσy和δux,δτxy。而SH波無論是垂直入射還是傾斜入射僅會產生擾動δuz,δτyz。垂直方向的擾動δuy,δσy會因為節理剛度系數kyy產生垂直方向的變形和應力反饋δuy,δσy；再由波動原理可知，一個平面上無數質點作為振源輻射出垂直擾動會形成平面波陣面，從而形成反射和透射的P波。類似的原理，垂直擾動δuy,δσy會因節理剛度系數kxy反射和透射SV波，因剛度系數kyz反射和透射SH波。同理，剪切擾動δux,δτxy和δuz,δτyz也能夠分別產生P波、SV波與SH波。

圖4 應力波在三維節理上轉化機理示意圖

Fig.4 Mutual transformation mechanism of stress wave on 3D rock joint

因此對三維巖石節理，當交叉耦合剛度系數不為零時，如圖5所示，P波、SV波或SH波會因剛度系數kxy,kyz,kxy而能夠相互轉化。

這里采用DDM(位移不連續方法)來研究應力波折反射問題，且節理的變形本構方程式(10)與應力連續條件作為求解的邊界條件，入射波、反射波與透射波的波動方程與符號約定如表1所示。

節理的變形[δuy,δux,δuz]T可以表示為

(11)

式中：[δuAy,δuAx,δuAz]T為巖石塊體A上節理的位移

(a) P波入射

(b) SV波入射

(c) SH波入射

變化量；[δuBy,δuBx,δuBz]T為巖石塊體B上節理的位移變化量。而由應力波的傳播方向和質點的偏振方向，可以得出巖石塊體A與塊體B的位移分量為

(12)

(13)

再由幾何變形方程，可以得出塊體A與塊體B中的應變分量εAij與εBij。使用廣義胡克定律，即可以得出塊體A、塊體B中的應力分量σAij與σBij。由DDM模型的假設：節理面上的應力連續條件σAij|y=0=σBij|y=0和變形方程式(10)，并令應力波的反射系數與透射系數為未知量，可以得出求解方程

X·E·C=ein·B

(14)

式中:X為求解方程的系數矩陣，與節理變形本構方程和入射波的頻率、角度有關，具體形式為

(15)

式中：λA,μA為巖石塊體A的拉梅常數;λB,μB為巖石塊體B的拉梅常數。E為指數矩陣

(16)

式中：C為反射應力波與透射應力波的幅值系數，具體形式與入射應力波的類型有關，當入射波為P波時，C為

(17)

式中：簡寫為C=[RP→P,RP→SV,RP→SH,TP→P,TP→SV,TP→SH]T，類似的當入射波為SV波時，折反射系數可以寫為C=[RSV→P,RSV→SV,RSV→SH,TSV→P,TSV→SV,TSV→SH]T；當入射波SH波時，折反射系數可以寫為C=[RSH→P,RSH→SV,RSH→SH,TSH→P,TSH→SV,TSH→SH]T。

表1 波動方程與符號約定

ein·B為式(14)的非齊次項，具體形式與入射波類型有關，當入射波為P波時，ein·B可以寫為

einBP=eiκAd(xsin α1-CAdt)×

(18)

當入射波為SV波時

(19)

當入射波為SH波時

(20)

由穩態應力波場可知，應力波在非連續面上傳播時應該滿足Snell定律。對于三維節理面而言，Snell定律這里可以表示為

ω=κAdCAd=κAsCAs=κBdCBd=κBsCBs

(21)

因此，由以上分析和推導可知，對于三維粗糙節理面，由于剛度系數矩陣中交叉耦合項存在，簡諧P波、SV波與SH波同時相互作用，之間能夠發生能量相互轉換，且SV波與SH波具有相同的反射角β2=γ2與透射角β3=γ3。將式(21)代入式(14)，即可以求出節理上的折反射系數C。

3 三維節理應力波折反射規律分析

三維巖石節理對應力波傳播的影響，與以往研究相比最大的不同點在于：三種類型的簡諧波在節理上可以相互轉化，且本構方程中剛度矩陣中交叉耦合項起到重要的作用。基于上節給出的應力波在三維粗糙節理折反射的理論解，這里針對交叉耦合剛度系數對節理折反射應力波的幅值系數影響進行分析。

3.1 入射角對應力波折反射系數的影響

圖6為P波入射時，折反射系數RP→P,TP→P,RP→SV,TP→SV,RP→SH和TP→SH隨入射角α1的變化規律。應力波的幅值系數隨入射角的變化曲線與過去DDM(無交叉耦合剛度項)所得的幅值系數類似，但是交叉耦合剛度系數的存在使得曲線有所變化，隨著節理中交叉耦合項的增大，與入射波同型的反射波的反射系數RP→P增大，且與入射波同型的透射波的幅值系數TP→P減小。若無交叉耦合剛度系數存在，則不會生成轉換波SH波；對于轉換波SV波，交叉耦合剛度系數的存在一定程度上改變了SV波折反射系數RP→SV,TP→SV的變化規律。但是總體上來說，隨著交叉耦合剛度系數的增加，轉換波的幅值系數是增大的。類似的規律在SV波入射時同樣成立，見圖7所示。SH波入射時，若無交叉耦合剛度系數存在，則轉換波P波與SV波則不會產生，如圖8所示；隨著交叉耦合剛度系數的增大，同型波反射系數RSH→SH增大，同型波透射系數TSH→SH減小，且轉換波折反射系數增大。

3.2 節理剛度系數對應力折反射系數的影響

如圖9所示，P波垂直入射僅會導致節理面上的質點在y軸方向運動，因此剛度系數kzx對應力波傳播沒有影響，見圖9(c)。隨著剛度系數kxy的增大，與入射波同型波的反射系數RP→P和轉換SV波RP→SV，TP→SV的增大，與入射波同型波的透射系數TP→P減小，此時沒有SH波產生，見圖9(a)。類似的，當變量為剛度系數kyz時，同型波的反射系數RP→P和轉換SH波RP→SH，TP→SH隨著kyz的增大而增加，同型入射波TP→P減小，且沒有SV波形成，見圖9(b)。相似的規律在SV波與SH波垂直入射時同樣存在，見圖10與圖11，這里不再贅述。

4 結 論

本文在討論三維巖石節理本構方程的基礎上，通過合理的簡化，分析了平面諧波在節理面上的作用機理，并給出了平面諧波在三維節理上的折反射理論解。并通過參數分析初步探討了P波、SV波和SH波入射時，反射波與透射波的幅值系數隨入射角以及交叉耦合剛度系數的變化規律。本研究得出以下有意義結論：

(1) 簡諧波P波，SV波或SH波在三維節理上可以相互轉化，由于節理的交叉耦合變形剛度系數的存在，使得P波、SV波與SH波在節理面上可以發生波動能量的相互傳遞和轉化。

(2) 三維粗糙節理本構方程中的交叉耦合剛度系數越大，則越多的入射波能量轉化為轉換波能量；且同型反射波幅值隨交叉耦合剛度系數增大而增大，同型透射波幅值隨之減小。

(a) RP→P

(b) TP→P

(c) RP→SV

(d) TP→SV

(e) RP→SH,TP→SH

(a) RSV→SV

(b) TSV→SV

(c) RSV→P

(d) TSV→P

(e) RSV→SH,TSV→SH

(a) RSH→SH

(b) TSH→SH

(c) RSH→P

(d) TSH→P

(e) RSH→SV

(f) TSH→SV

(a) 變量為kxy,且此時kzy=kzx=0

(b) 變量為kzy,且此時kxy=kzx=0

(c) 變量為kzx,且此時kzy=kxy=0

(a) 變量為kxy,且此時kzy=kzx=0

(b) 變量為kzy,且此時kxy=kzx=0

(c) 變量為kzx,且此時kzy=kxy=0

(a) 變量為kxy,且此時kzy=kzx=0

(b) 變量為kzy,且此時kxy=kzx=0

(c) 變量為kzx,且此時kzy=kxy=0

(3) 相對于以往采用僅存在主剛度系數節理描述應力波折反射規律的研究，節理中交叉耦合剛度系數的存在相當于是對原模型的擴展，交叉耦合剛度系數控制了應力波能量流向。

本文對三維巖石節理應力波傳播規律初步探討和分析，通過變形剛度系數考慮了由于復雜節理表面形態對應力波折反射規律的影響。本研究有利于加深對巖體動力學的理解，對于巖土地震工程和工程物探具有一定的積極意義。然而本文由于采用線剛度模型描述低幅值應力波入射時的力-變形行為，忽略由于較大變形節理面嵌合度變化等現象帶來的非線性力學變化規律，因此不適用于高幅值應力波入射的情形。因此下一步工作應該結合數值計算與實驗工作對高應力幅值應力波入射三維節理的傳播行為進行研究。