高等巖石力學巖石力學有限元法

高等巖石力學巖石力學有限元法第1頁，共31頁，2023年，2月20日，星期二巖石力學有限元法第一節(jié)概論第二節(jié)施工建造過程的模擬第三節(jié)節(jié)理及不連續(xù)面的模擬第四節(jié)多節(jié)理巖體的模擬第五節(jié)無限域單元及其應(yīng)用第六節(jié)巖體工程中的彈塑性問題第七節(jié)無拉力分析及節(jié)理非線性分析第八節(jié)巖土工程三維非線性有限元程序系列第2頁，共31頁，2023年，2月20日，星期二

第一節(jié)概論一、有限元法在工程中的應(yīng)用

有限元自50年代發(fā)展至今，以成為求解復(fù)雜的巖石力學及巖土工程問題的有力工具，并已為愈來愈多的工程科技人員所熟悉。在求解像彈塑性及流變,動力，非穩(wěn)態(tài)滲流等時間相關(guān)性問題，以及溫度場在求解像彈塑性及流變、動力、非穩(wěn)定滲流等時間相關(guān)性問題，以及溫度場、滲流、應(yīng)力場的耦合問題等復(fù)雜的非線性問題中的效能已使它成為在巖石力學領(lǐng)域中應(yīng)用最廣泛的數(shù)值分析手段。特別是近十余年，在工程應(yīng)用方面已有了較大的進展，并引起廣大工程科技人員的興趣。在巖土工程有關(guān)專業(yè)的大學生和研究生中，有限元已被列為獨立的課程。在本章中將著重討論有限元法在巖石力學中應(yīng)用的有關(guān)問題。涉及有限元法基本原理，方法及求解技術(shù)方面的基礎(chǔ)知識。第3頁，共31頁，2023年，2月20日，星期二有限元法基本方程就數(shù)學概念來說，有限元法是通過變分原理(或加權(quán)余量法)和分區(qū)插值的離散化處理把基本支配方程轉(zhuǎn)化為線性代數(shù)方程，把求解待解域內(nèi)的連續(xù)場函數(shù)轉(zhuǎn)化為求解有限個離散點（節(jié)點）處的場函數(shù)值。顯然，這種離散化的處理是一種近似，因面中有當單元劃分得充分小時，才能保證滿意的求解精度。由于單元充分小，在一個微小的單元內(nèi)，未知場函數(shù)可以采用十分簡單的代數(shù)多項式近似地表述。通?？扇槿缦碌牟逯敌问剑?/p>

(4-1)第4頁，共31頁，2023年，2月20日，星期二下標m表示單元的節(jié)點數(shù)目，下標e表示單元的序號。有限元法即是以所有節(jié)點處的1值作為基本未知量。對于二維的工程物理問題，如熱傳導(dǎo)、滲流等，其基本支配方程為如下形式的準諧方程：

(4-2)邊界條件為：①，邊界在上②邊界在上第5頁，共31頁，2023年，2月20日，星期二當=常量時，上式即退化為波松方程。由變分原理可知，滿足方程式（4-2）及相應(yīng)邊界條件的場函數(shù)，應(yīng)使如下的泛函數(shù)駐值，即：

(4-3)在對被考察域A進行離散化的情況下，此泛函駐值轉(zhuǎn)化為多元函數(shù)的駐值：（4-4）

第6頁，共31頁，2023年，2月20日，星期二對于一個單元體的子域內(nèi)，把方程式（4-1）帶入方程式（4-3)并進行變分后即可得到單元的基本方程：

(4-5)

式中即單元特性矩陣：即為節(jié)點未知量組成的矢量：是由在上的邊界條件給定的a（x,y）有關(guān)的矩陣，故僅對a≠0的二類邊界上的單元才予考慮。對整個求解域A由變分駐值條件式（4-4）可得或?qū)憺?/p>

(4-6)第7頁，共31頁，2023年，2月20日，星期二

式中，[k]為由各單元特性矩陣及按節(jié)點號組集得到的總體特性矩陣：{u}為所有節(jié)點的待求值組成的矢量：{p}為與Q（x,y）有關(guān)的矢量及q(x,y)有關(guān)的矢量。解線性方程式（4-6）即可求得場函數(shù)在各單元節(jié)點處的值。對于彈性力學問題，可通過最小勢能原理或虛功原理導(dǎo)出有限元法的基本方程。有限元法求解彈性力學問題通常以位移作為基本未知量。位移在直角坐標系中沿坐標軸X,Y,Z的分量分別表示為，故單元位移模式這時可寫為

(4-7)式中，l為3×3階單位矩陣；Nt稱為單元的位移插值函數(shù)或形函數(shù)（對二維問題l為2×2階，位移分量）；{}為由單元各節(jié)點位移分量形成的矢量稱為移矢量，利用彈性力學的幾何方程及物理方程可導(dǎo)出單元的應(yīng)變及應(yīng)力表達式：第8頁，共31頁，2023年，2月20日，星期二

（4-8）（4-9）應(yīng)用最小勢能原理或虛功原理可以推導(dǎo)出單元剛度矩陣的表達式

（4-10，a）

（4-10，b）

（4-10，c）第9頁，共31頁，2023年，2月20日，星期二

式中，{q}為分布體力；{p}為分布面力，通常為坐標的已知函數(shù)，最簡單也最常見的情況是{q}及{p}均為常量，上式的積分則可簡化。對于一般準諧方程描敘的工程物理問題及彈性力學問題，有限元法的公式表述及分析方法具有完全相同的格式，這是有限元法的主要特點。這種統(tǒng)一的格式有利于實現(xiàn)具有廣泛通用性的計算機程序。二、巖體力學問題的特點

大多數(shù)巖土工程問題如結(jié)構(gòu)-巖體相互作用，巖土邊坡、地下工程等，都涉及無限域或半無限域。處理這類問題通常是在有限的區(qū)域進行離散化。為了使這種離散化不會產(chǎn)生大的誤差，必須取足夠大的計算范圍，并應(yīng)使假定的外邊界條件盡可能接近真實狀態(tài)。第10頁，共31頁，2023年，2月20日，星期二理論分析和計算實踐表明，當由于結(jié)構(gòu)或工程開挖釋放荷載作用于巖體某一部位時，對周圍巖體的應(yīng)力及位移有明顯影響的范圍大約是開挖或結(jié)構(gòu)物與巖體作用面的輪廓尺寸的2.5～3倍，在此范圍之外，影響甚小，可忽略其影響?？紤]到有限元離散誤差和計算誤差，為了保證必要的計算精度，計算范圍應(yīng)取不小于3～4倍，如圖4-1所示。在這種情況下，外邊界可采取兩種方式處理：其一為在距離荷載作用足夠遠的外邊界處位移為零；其二為假定外邊界為受力邊界（包括p=0）。這兩種方式都同實際的無限域不完全一致，因而都有一定誤差。這種誤差隨著計算區(qū)域的減小而增大，并且在靠近外邊界處都遠離外邊界處的誤差大，我們稱此現(xiàn)象為“邊界效應(yīng)”。第11頁，共31頁，2023年，2月20日，星期二三、初始地應(yīng)力場與釋放荷載由于自重及地質(zhì)作用，自然狀態(tài)的巖體處于一定的初始地應(yīng)力狀態(tài)，在結(jié)構(gòu)荷載作用圖4-1有限元離散化的計算范圍第12頁，共31頁，2023年，2月20日，星期二

下，巖體內(nèi)的應(yīng)力為荷載產(chǎn)生的應(yīng)力與初始地應(yīng)力之和。因此，初始地應(yīng)力的大小直接影響到計算結(jié)果。正確決定巖體的初始應(yīng)力狀態(tài)是有限元分析中的一個重要問題，至今未能得到妥善的解決。由于初始地應(yīng)力的存在，巖石的開挖將導(dǎo)致部分巖體卸載，成而使一定范圍內(nèi)的應(yīng)力場發(fā)生變化。正確模擬這種開挖效果也是巖體力學問題的一個重要的特點，把由此所得的位移即作為由于工程開挖長生的巖體位移，由此所得的應(yīng)力場同初始地應(yīng)力場迭加即為開挖后的應(yīng)力場，這種模擬開挖效果的方法可如圖4-3所示?！伴_挖釋放荷載”是按照已知的初始地應(yīng)力由等效原則來確定的。

第13頁，共31頁，2023年，2月20日，星期二第二節(jié)施工建造過程的模擬本節(jié)主要討論如何正確模擬開挖及工程建造的過程一、開挖釋放荷載對于具有已知初始應(yīng)力場的巖體開挖工程，沿預(yù)計開挖線上各點的初始應(yīng)力為已知，在進行離散化的情況下，可假沿開挖邊界上兩相鄰節(jié)點之間的初始應(yīng)力呈線性變化，則對于任一開挖邊界點，開挖引起等效釋放荷載為

(1)第14頁，共31頁，2023年，2月20日，星期二

式中，上標i，i-1及i+1為沿開挖邊界上的有限元網(wǎng)格的節(jié)點號如果坐標X，Y軸正沿主應(yīng)力軸，則上式可簡化為

(2)若原始應(yīng)力場均為應(yīng)力場，寄節(jié)點i、i-1、i+1等各點應(yīng)力場相等，則上式可簡化為

(3)第15頁，共31頁，2023年，2月20日，星期二若X，Y軸同應(yīng)力主軸重合，可進一步簡化為

(4)考慮到存在一個初始應(yīng)力場的情況，開挖后的實際應(yīng)力場為初始應(yīng)力場與開挖釋放應(yīng)力場之和。即以上是在已知開挖邊界上各點初始應(yīng)力的情況下計算各邊界節(jié)點上的等效釋放荷載，邊界節(jié)點上的初始應(yīng)力可按下述方法確定：第16頁，共31頁，2023年，2月20日，星期二1、對于已知均勻初始應(yīng)力場，有等效釋放荷載可由(3)式求得

2、對于非均勻初始應(yīng)力場，一般先用有限元法計算初始應(yīng)力場，然后把要挖去的區(qū)域內(nèi)各單元改變?yōu)椤翱諉卧庇嬎汩_挖結(jié)果當有限元計算直接給出各單元節(jié)點應(yīng)力時，則開挖邊界上各節(jié)點的應(yīng)力也為已知。利用（2）式可以求出。第17頁，共31頁，2023年，2月20日，星期二二、施工過程模擬為了模擬施工建造過程，在確定離散化網(wǎng)格時，必須要考慮各步施工的情況及結(jié)構(gòu)特征。如下圖為一地下洞室，開挖及支護分上、下兩部分進行。施工步驟如下：

圖4-7洞室斷面及施工步驟第18頁，共31頁，2023年，2月20日，星期二第三節(jié)

節(jié)理及不連續(xù)面的模擬一、平面問題節(jié)理單元——Goodman單元采用特殊的“節(jié)理單元”模擬地質(zhì)不連續(xù)面是由R.EGoodman(1976)最先提出來的。一種無厚度4節(jié)點節(jié)理單元也常被稱作Goodman單元，如下圖所示，它是一段直接接觸的平面，承受沿此界面的切向力及法向力。節(jié)理的應(yīng)力-應(yīng)變關(guān)系對局部坐標（s,n)可寫為也可簡化為第19頁，共31頁，2023年，2月20日，星期二式中，分別為切向應(yīng)力及法向應(yīng)力；，為節(jié)理的切向及法向剛度；，為單元兩側(cè)對應(yīng)點的相對切向及法向位移。

無厚度節(jié)理單元假定位移沿節(jié)理單元長度呈線性變化，則可導(dǎo)出節(jié)理單元對于局部坐標s-t的單元剛度矩陣為第20頁，共31頁，2023年，2月20日，星期二對整體坐標的轉(zhuǎn)換為式中，[T]為坐標變換矩陣，它是具有如下的分塊形式的對角陣：

式中

上述節(jié)理單元的缺點是由于無厚度，計算中可能發(fā)生節(jié)理上下面相互“嵌入”的現(xiàn)象，必須對這種嵌入量做人為的限制，以免導(dǎo)致較大誤差。此外，當節(jié)理的上下面發(fā)生相對轉(zhuǎn)角位移時，也講產(chǎn)生誤差。Goodman曾于1976年對節(jié)理單元做了修正，考慮到相對轉(zhuǎn)角的影響，假定以節(jié)理單元的中點作為計算點，則其應(yīng)力-應(yīng)變關(guān)系表示為

第21頁，共31頁，2023年，2月20日，星期二式中，為節(jié)理中點力矩；為相對轉(zhuǎn)角，以節(jié)理的3—4節(jié)點的平面逆時針為正。仍假定位移沿長度線性變化，可導(dǎo)出節(jié)理單元剛度矩陣為二、變厚度節(jié)理單元下圖所示的一種六節(jié)點變厚度的節(jié)理單元，具有更廣泛的適用性，可用于等厚度、變厚度、平面及曲面的節(jié)理，這種單元實際上是細長四邊形的等參單元，單元的位移函數(shù)及坐標變換關(guān)系可寫為第22頁，共31頁，2023年，2月20日，星期二

這種較窄長的單元，在和方向采取不同階次的插值函數(shù)(

方向為二次，方向為一次），即變厚度節(jié)理單元第23頁，共31頁，2023年，2月20日，星期二對于有一定厚度的夾層，可采用一般平面應(yīng)變問題的本構(gòu)關(guān)系，按一般平面等參單元處理，其單元剛度矩陣為式中det即為雅可比轉(zhuǎn)換。應(yīng)用高斯積分完成單元剛度計算，在方向可取三點積分，在方向可取兩點積分。式中由幾何方程導(dǎo)出，彈性矩陣可取一般各向同性或正交異性彈性矩陣，當采用正交異性彈性矩陣時要注意坐標至整體坐標之間的變換。當厚度很小時，可以作為等厚或無厚度的節(jié)理單元。此時各應(yīng)力分量沿厚度方向為常量，則單元退化為上圖的形式。插值函數(shù)可簡化為第24頁，共31頁，2023年，2月20日，星期二三、三維問題節(jié)理單元

由二維節(jié)理單元可方便的推廣至三維問題，常用三維節(jié)理單元如圖所示單元位移插值函數(shù)可表示為

等厚度三維節(jié)理由單元第25頁，共31頁，2023年，2月20日，星期二對圖a插值函數(shù)為圖b插值函數(shù)為座標變換關(guān)系為第26頁，共31頁，2023年，2月20日，星期二考慮節(jié)理上、下面一對點的相對位移，對整體座標有式中為3×3階單矩陣，為對整體座標的節(jié)點位移，通常需要考慮節(jié)理沿其法向及切向的相對位移及應(yīng)力。為此，取節(jié)理中面任一點的法線方向作為軸建立局部動座標系則平面即為過單元該點的切平面，這時有