FLAC3D原理..

1、【精品文檔】如有侵權(quán)，請聯(lián)系網(wǎng)站刪除，僅供學(xué)習(xí)與交流FLAC3D原理.精品文檔.2.2 三維數(shù)值模擬方法及其原理2.2.1 FLAC3D工程分析軟件特點FLAC3D是由美國Itasca Consulting Group, Inc. 為地質(zhì)工程應(yīng)用而開發(fā)的連續(xù)介質(zhì)顯式有限差分計算機軟件。FLAC即Fast Lagrangian Analysis of Continua 的縮寫。該軟件主要適用于模擬計算巖土體材料的力學(xué)行為及巖土材料達(dá)到屈服極限后產(chǎn)生的塑性流動，對大變形情況應(yīng)用效果更好。FLAC3D程序在數(shù)學(xué)上采用的是快速拉格朗日方法，基于顯式差分來獲得模型全部運動方程和本構(gòu)方程的步長解，其本構(gòu)方

2、程由基本應(yīng)力應(yīng)變定義及虎克定律導(dǎo)出，運動平衡方程則直接應(yīng)用了柯西運動方程，該方程由牛頓運動定律導(dǎo)出。計算模型一般是由若干不同形狀的三維單元體組成，也即剖分的空間單元網(wǎng)絡(luò)區(qū)，計算中又將每個單元體進(jìn)一步劃分成由四個節(jié)點構(gòu)成的四面體，四面體的應(yīng)力應(yīng)變只通過四個節(jié)點向其它四面體傳遞，進(jìn)而傳遞到其它單元體。當(dāng)對某一節(jié)點施加荷載后，在某一個微小的時間段內(nèi)，作用于該點的荷載只對周圍的若干節(jié)點（相鄰節(jié)點）有影響。利用運動方程，根據(jù)單元節(jié)點的速度變化和時間，可計算出單元之間的相對位移，進(jìn)而求出單元應(yīng)變，再利用單元模型的本構(gòu)方程，可求出單元應(yīng)力。在計算應(yīng)變過程中，利用高斯積分理論，將三維問題轉(zhuǎn)化為二維問題而使其簡

3、單化。在運動方程中，還充分考慮了巖土體所具有的粘滯性，將其視作阻尼附加于方程中。FLAC3D具有一個功能強大的網(wǎng)格生成器，有12種基本形狀的單元體可供選擇，利用這12種基本單元體，幾乎可以構(gòu)成任何形狀的空間立體模型。FLAC3D主要是為地質(zhì)工程應(yīng)用而開發(fā)的巖土體力學(xué)數(shù)值評價計算程序，自身設(shè)計有九種材料本構(gòu)模型：（1）空模型（Null Model）（2）彈性各向同性材料模型（Elastic, Isotropic Model）（3）彈性各向異性材料模型（Elastic, anisotropic Model）（4）德拉克-普拉格彈塑性材料模型（Drucker-Prager Model）（5）莫爾-庫

4、倫彈塑性材料模型（Mohr-Coulomb Model）（6）應(yīng)變硬化、軟化彈塑性材料模型（Strain-Hardening/Softening Mohr-Coulomb Model）（7）多節(jié)理裂隙材料模型（Ubiquitous-Joint Model）（8）雙曲型應(yīng)變硬化、軟化多節(jié)理裂隙材料模型（Bilinear Strain-Hardening/Softening Ubiquitous-Joint Model）（9）修正的Cam粘土材料模型（Modified Cam-clay Model）除上述本構(gòu)模型之外，F(xiàn)LAC3D還可進(jìn)行動力學(xué)問題、水力學(xué)問題、熱力學(xué)問題等的數(shù)值模擬。在邊界條件及

5、初始條件的考慮上，F(xiàn)LAC3D軟件十分靈活方便，可在數(shù)值計算過程中隨時調(diào)整邊界條件和初始條件。FLAC3D具有強大的后處理功能，用戶可以直接在屏幕上繪制或以文件形式創(chuàng)建或輸出打印多種形式的圖形、文字，用戶還可根據(jù)各自的需要，將若干個變量合并在同一幅圖形中進(jìn)行研究分析。FLAC3D軟件還可對各種開挖工程或施加支護工程等進(jìn)行數(shù)值仿真模擬，軟件自身設(shè)計有錨桿、錨索、襯砌、支架等結(jié)構(gòu)元素，可以直接模擬這些支護于圍巖（土）體的相互作用。FLAC3D擁有可以自行設(shè)計的FISH語言，用戶可根據(jù)自身需求，自己設(shè)計材料的本構(gòu)模型、屈服準(zhǔn)則、支護方案、復(fù)雜形狀的開挖方式等工作。特別注意的是，巖石是一種脆性材料，當(dāng)

6、外荷載達(dá)到巖石強度后，材料發(fā)生斷裂破壞，產(chǎn)生弱化現(xiàn)象，應(yīng)屬于彈塑性體。在FLAC3D中，一般對于彈塑性材料，判斷其破壞與否的基本準(zhǔn)則有兩個，即Drucker-Prager準(zhǔn)則和Mohr-Coulomb準(zhǔn)則。根據(jù)室內(nèi)巖石力學(xué)性質(zhì)試驗結(jié)果，其典型應(yīng)力應(yīng)變曲線反映出巖體破壞包絡(luò)線符合莫爾庫倫屈服準(zhǔn)則，故本次建立的本構(gòu)力學(xué)模型選擇莫爾庫倫彈塑性材料模型為宜。2.2.2 FLAC3D分析計算原理計算所采用的數(shù)學(xué)模型是根據(jù)彈塑性理論的基本原理（應(yīng)變定義、運動定律、能量守衡定律、平衡方程及理想材料的連續(xù)性方程等）而建立的。2.2.2.1 基本約定在數(shù)學(xué)及數(shù)值模型的表達(dá)式中，符號有一定的約定含義，一般A表示張

7、量，Aij表示張量A的（i,j）分量，a表示矢量，ai表示矢量a的i分量，a,i表示a對xi的偏導(dǎo)數(shù)。xi，ui，vi和dvi/dt，（i=1,3）分別表示一點的位置矢量分量、位移矢量分量、速度矢量分量和加速度矢量分量。2.2.2.2 數(shù)學(xué)模型（一）柯西（Cauchy）應(yīng)力張量與柯西公式對于一個具有體積V的封閉曲面s的物體，在其上取一表面元素Ds，這個表面元素的單位外法向矢量為n，在某一時刻t，在表面元素對于連續(xù)介質(zhì)中一點，作用著對稱的應(yīng)力張量sij，根據(jù)Ds上作用有力DP，則極限稱為表面力。若用ti表示T的分量，則在三維直角坐標(biāo)系中可有關(guān)系式 （1）這個關(guān)系式稱為柯西公式，其中，sij稱為柯

8、西應(yīng)力張量。（二）應(yīng)變速率和旋轉(zhuǎn)速率如果介質(zhì)質(zhì)點具有運動速度矢量v，則在一個無限小的時間dt內(nèi)，介質(zhì)會產(chǎn)生一個由vidt決定的無限小應(yīng)變，對應(yīng)的應(yīng)變速率分量xij為 （2）而其旋轉(zhuǎn)速率分量wij為 （3）（三）運動及平衡方程根據(jù)牛頓運動定律與柯西應(yīng)力原理，如果質(zhì)點作用著應(yīng)力sij與體力bi，且具有速度vi，則在無限小時間段dt內(nèi)，它們之間的關(guān)系為 （4）式中，r為質(zhì)點密度。（4）式稱為柯西運動方程。當(dāng)質(zhì)點的加速度為零時，上式變?yōu)殪o力平衡方程 （5）（四）本構(gòu)方程上述（4）式與（5）式組成的方程組中含有9個方程，15個未知量，其中12個是應(yīng)力與應(yīng)變速率分量，3個是速度分量。其余6個關(guān)系式則由本構(gòu)

9、方程提供，本構(gòu)方程一般具有如下形式 （6）式中，為應(yīng)力變化速率，H表示一個特定的函數(shù)關(guān)系，k為與荷載歷史有關(guān)的參數(shù)。2.2.2.3 數(shù)值模型FLAC3D的數(shù)值剖分網(wǎng)格在計算中是按照四面體進(jìn)行的，四面體的節(jié)點也既是網(wǎng)格剖分的節(jié)點，因此，每個計算單元有4個面和4個節(jié)點（見圖2-1）。單元面單元節(jié)點圖2-1 四面體單元示意圖（一）空間微分的有限差分逼近對于一個計算單元，若內(nèi)部各質(zhì)點速度為一連續(xù)的矢量場v，則根據(jù)高斯（Gauss）積分原理有 （7）式中，n為外法向單位矢量場。由于單元體的應(yīng)變速率是連續(xù)的，因此可以近似認(rèn)為速度是線性變化的，則（7）式可用下面的求和公式近似逼近 （8）式中，V為單元體體積

10、，為單元體某一面的面積，f為單元體面數(shù)，f=1,4，為面平均速度的i分量。由于速度場是線性的，則有 （9）式中，l為單元體節(jié)點數(shù)，l=1,4；vil為l節(jié)點的i速度分量。將（9）式代入（8）式可得 （10）根據(jù)正交原理有 （11）則將（10）式兩邊除以V，并將（11）式代入可得 （12）因此有 （13）（二）運動方程的節(jié)點公式根據(jù)前面對質(zhì)點運動方程的討論，對于連續(xù)介質(zhì)，當(dāng)處于平衡狀態(tài)時，其平衡方程為 （14） （15）介質(zhì)可由若干個作用著體力B，各自產(chǎn)生一定變形的四面體組成，設(shè)節(jié)點力為fn，n=1,4，利用虛功原理，假定單元體節(jié)點具有速度dvn，內(nèi)部具有變形速率dx，則由節(jié)點力fn和體力B所做

11、的外功率與sij所做的內(nèi)功率應(yīng)當(dāng)相等。外功率E可用下式表示 （16）而內(nèi)功率I為 （17）由（13）式，對于恒定應(yīng)變速率單元體可有 （18）由于應(yīng)力張量是對稱的，定義矢量Tl （19）可得 （20）將（15）式代入（16）式得 （21）式中，Eb和EI分別是體力rbi和慣性力產(chǎn)生的外功率，對于單元體內(nèi)恒定的體力rbi，Eb可寫成 （22）而EI可以寫成 （23）如前所述，單元體內(nèi)部速度場以線性變化，為方便描述，選取單元體質(zhì)心為原點，、為坐標(biāo)軸，可有 （24）式中，Nn(n=1,4)是具有如下形式的線性函數(shù) （25）、（n=1,4）則由下面方程決定 （26）dnj為克羅內(nèi)克爾（Kronecker

12、）記號，根據(jù)質(zhì)心定義可得 （27）用克萊默（Cramers）法則解得，則 （28）而 （29）則 （30）在穩(wěn)定狀態(tài)下，內(nèi)功率與外功率必定相等，則有 （31）上式最后一項，如果單元體內(nèi)r恒定，則根據(jù)質(zhì)心定義有 （32）若將rV/4看作是假定的節(jié)點質(zhì)量mn，則上式變?yōu)?（33）故 （34）根據(jù)牛頓定律有 （35）式中，F(xiàn)i<l>為節(jié)點l所受力的i分量，M<l>為節(jié)點l的質(zhì)量，為節(jié)點l在Fi<l>作用下產(chǎn)生的加速度，nn為包含在全部連續(xù)介質(zhì)中單元節(jié)點總數(shù)。設(shè)質(zhì)點質(zhì)量，不平衡力F<l>為 （36）（三）時間微分的顯式有限差分對于上面（35）式，可寫成

13、（37）FLAC3D計算中，對于上式的逼近是在一個足夠小Dt時間段內(nèi)，采用中心差分格式，即上式變?yōu)?（38）同理，節(jié)點的坐標(biāo)差分公式與位移差分公式分別為 （39） （40）（四）本構(gòu)方程的增量形式FLAC3D中，假定在Dt時間段內(nèi)，速度是不變的，本構(gòu)方程的增量形式可表示為 （41）記 （42）則Deij為Dt時間段內(nèi)應(yīng)變增量。應(yīng)力增量由下式確定 （43）式中，由下式定義 （44）而 （45）2.2.2.4 本構(gòu)模型（一）彈塑性理論的增量關(guān)系一般來講，彈塑性體破壞的判斷準(zhǔn)則具有如下形式 （46）式中，f為一特定的函數(shù)關(guān)系，為n維應(yīng)力矢量i分量。彈塑性體的變形為彈性變形與塑性變形的總和，即 （47

14、）彈性應(yīng)變與應(yīng)力之間的關(guān)系為 （48）其中，Si為線性函數(shù)。塑性變形可用下式表示 （49）式中，l為一常數(shù)，g為的某一函數(shù)關(guān)系。將（47）、（49）式代入（48）式得 （50）對于新應(yīng)力，仍有 （51）當(dāng)為線性函數(shù)時，（51）式可寫成 （52）其中，。由于，再將（50）式代入（52）式可得 （53）若記則有 （54）因此有 （55）則新應(yīng)力為 （56）（二）莫爾庫倫（Mohr-Coulomb）模型本構(gòu)關(guān)系本次計算的本構(gòu)模型選擇的是莫爾庫倫彈塑性材料模型，彈塑性體可產(chǎn)生彈性及塑性兩部分變形，根據(jù)虎克定律（Hookes law），應(yīng)力應(yīng)變的關(guān)系為 （57）其中，E為剛度矩陣。對于彈性變形，應(yīng)力增量

15、可由下式確定 （58）式中，a1、a2可由剪切模量與體積模量得出（58）式還可寫成 （59）（三）破壞準(zhǔn)則與流動法則若，則莫爾庫倫破壞準(zhǔn)則具有形式f s = 0，其中 （60）對于拉張破壞，有f t = 0的形式，其中 （61）上兩式中，f為內(nèi)摩擦角，c為內(nèi)聚力，s t為抗拉強度， （62）對于材料的抗拉強度s t，其值不可能超過s3（見圖2-2），st的最大值為 （63）對于塑性流動狀態(tài)，設(shè)g s與g t分別為剪切破壞和拉張破壞所對應(yīng)的函數(shù)關(guān)系，有 （64） （65）式中，y為膨脹角， （66）f t = 0s3s1= 0=As3=s1=Bf s = 0s tc / tan fC圖2-2 莫爾庫倫破壞準(zhǔn)則示意圖流動準(zhǔn)則可定義為函數(shù)關(guān)系，在s 1、s 3平面內(nèi)將和以上部分分為兩個區(qū)域（見圖2-3），其函數(shù)表達(dá)式如下1區(qū)2區(qū)圖2-3 莫爾庫倫模型塑性流分區(qū) （67）其中，a p與s p由下式定義 （68）從圖2-3中可見，若某點的應(yīng)力狀態(tài)在1區(qū)時，則發(fā)生剪切破壞，在2區(qū)時