基于比例邊界有限元法計算應力強度因子的不確定量化分析

胡昊文,陳燈紅,王乾峰,胡記磊,駱 歡

(1.防災減災湖北省重點實驗室(三峽大學),湖北 宜昌 443002;2.三峽大學 土木與建筑學院,湖北 宜昌 443002)

裂紋存在于眾多結(jié)構(gòu)中,例如混凝土結(jié)構(gòu)通常是帶裂縫工作的,巖塊表面存在著許多裂隙,等。斷裂問題一直都是計算固體力學的難點和重點。其中,應力強度因子是預測載荷作用下結(jié)構(gòu)中裂紋的產(chǎn)生和擴展的重要參數(shù)。分析斷裂問題時,有限元法(finite element method,FEM)以其成熟的商用性和對不同邊界、不同荷載狀況下的普遍適用性成為求解斷裂問題最主要的數(shù)值方法。但傳統(tǒng)有限元法在處理斷裂問題時會面臨一些挑戰(zhàn)。比如三角形或四邊形單元不能準確地計算裂紋尖端的奇異應力,并且裂紋尖端的網(wǎng)格需要幾個數(shù)量級的加密。一些學者使用邊界元法[1-2](boundary element method,BEM)來處理斷裂問題。邊界元法只需離散結(jié)構(gòu)表面,將問題的維數(shù)減少一維。其作為一種半解析法,邊界元法能夠更準確地計算應力強度因子。對于復雜問題,邊界元法所需的基本解很難解析獲得甚至是不存在,同時積分方程還會產(chǎn)生奇異積分。這些困難限制了邊界元法的發(fā)展。

Song等[3-4]提出的比例邊界有限元法(scaled bountary finite element method,SBFEM),是一種新的半解析數(shù)值方法,它克服了傳統(tǒng)有限元法的一些局限性,并借鑒邊界元法的思想,只在結(jié)構(gòu)表面進行離散,將問題維數(shù)減少一維,輻射邊界條件自然滿足,且不需要基本解。這些優(yōu)勢使其成為學者們所青睞的一種數(shù)值方法。該方法已成功應用于無限域問題[5-8]、斷裂力學[9-12]、波傳播[13-14]、接觸問題[15]等多個領(lǐng)域,關(guān)于SBFEM最新研究詳見相關(guān)綜述[16]。當處理線性斷裂力學時,任何類型的奇異應力都能解析表示,且不需要局部細化和富集函數(shù)[17]。此外,比例邊界有限元法多邊形單元的形函數(shù)包含奇異項,可以準確地表征應力集中系數(shù)[18]、應力強度因子[19]和T-應力[20]。

目前成熟的數(shù)值算法,大都是基于確定性的分析。然而,在工程問題中存在諸多不確定性因素,例如結(jié)構(gòu)的幾何尺寸,材料自身的力學參數(shù),荷載的隨機性等等。若忽略這些因素的影響,可能會導致結(jié)構(gòu)失效,甚至引發(fā)災難性的后果。綜合考慮不確定性因素,使數(shù)值分析更符合實際工程,不確定性量化[21](uncertainty quantification,UQ)分析成為近些年的熱點方向。

不確定量化本質(zhì)上是研究輸入的不確定造成的輸出不確定。針對這一特點,學者們處理不確定量化分析的步驟主要分為三步。第一步:構(gòu)建物理模型。不確定量化的結(jié)果取決于高精度的數(shù)值模型。這一步通常采用成熟的確定性算法,諸如有限元法、邊界元法、比例邊界有限元法等。不確定量化與這些數(shù)值算法結(jié)合,產(chǎn)生了隨機有限元法、隨機邊界元法等經(jīng)典研究方法。第二步:描述不確定性參數(shù)。通常使用隨機變量或隨機場。例如攝動法,將一個隨機函數(shù)在其均值附近進行Taylor展開,然后對其進行合理的截斷;多項式混沌展開,在隨機參數(shù)空間對精確解作一個有限階多項式展開,求解相關(guān)系數(shù)以獲得精確解的全部信息;蒙特卡羅模擬(Monte Carlo simulation,MCS)[22-23],假設隨機變量所遵循的概率分布,通過計算隨機抽樣的樣本響應,獲得統(tǒng)計特征(期望、標準差等)。第三步:將隨機變量或隨機場導入數(shù)值模型中,進行不確定性傳播,得到高階矩、概率密度函數(shù)、累積密度函數(shù)等統(tǒng)計特征,用于評估不確定性參數(shù)對結(jié)構(gòu)響應的影響。

描述不確定性參數(shù)步驟中,蒙特卡羅模擬方法以其思想簡單和廣泛適用成為很多學者首選的方法。此外,蒙特卡羅模擬的收斂速度只與隨機采樣次數(shù)N有關(guān),而與求解問題的維度無關(guān)。所以,在很多高維問題的求解過程中,蒙特卡羅方法往往是最合適的甚至是唯一的選擇。處理工程問題時,單點響應計算耗時嚴重,執(zhí)行大量采樣來確保蒙特卡羅模擬收斂的成本是無法負擔的,直接的蒙特卡羅模擬不得不與其他方法相結(jié)合,以節(jié)約計算成本,提高計算效率。Gu[24]從硬件方面減少硬件資源消耗、提高硬件計算的并行性和加速性能,從而減少蒙特卡羅模擬硬件方面的成本;Fishman[25]優(yōu)化蒙特卡羅模擬采樣策略,一定程度上減少了取樣個數(shù),降低了計算量;Ding等[26]使用模型降階算法加速蒙特卡羅模擬過程,分析結(jié)構(gòu)的高維不確定性。

本文構(gòu)建一種新穎的基于比例邊界有限元法求解斷裂問題的不確定量化分析方法。首先,采用比例邊界有限元法計算應力強度因子(stress intensity factor,SIF);然后,使用蒙特卡羅模擬進行不確定量化分析,并對直接的蒙特卡羅模擬進行改進;最后采用奇異值分解,構(gòu)造低階的子空間,用徑向基函數(shù)對子空間進行近似,通過子空間線性組合獲得新的結(jié)構(gòu)響應,實現(xiàn)基于MCS的快速不確定量化分析。

1 基于比例邊界有限元法的蒙特卡羅模擬

1.1 比例邊界有限元法

Wolf等詳細介紹了比例邊界有限元的推導過程,本文只給出一些必要的公式。在二維域中給出SBFEM所需的局部坐標,如圖1所示。選擇一個點O作為比例中心,對比例中心的唯一要求是S的整個邊界可視。徑向局部坐標ξ,也被稱作比例因子,從比例中心指向邊界上任意一點。η為環(huán)向局部坐標。這種由比例系數(shù)和單元局部坐標定義的轉(zhuǎn)換關(guān)系稱為比例邊界轉(zhuǎn)化,并且這種轉(zhuǎn)化關(guān)系是唯一的。對于有限域,0≤ξ≤1,而對于無限域,1≤ξ≤∞。在(η,ξ)點處的位移可以用分段插值來近似

圖1 比例邊界有限元的坐標系統(tǒng)

{u(ξ,η)}=[N(η)]{u(ξ)}=[N1(η)[I],

N2(η)[I],…]{u(ξ)}

(1)

式中:[N(η)]為在環(huán)向上的形函數(shù);u(ξ)為沿徑向線的位移函數(shù),僅與ξ有關(guān);[I]為一個2×2的單位矩陣。

應變表示為

{ε(ξ,η)}=[B1(η)]{u(ξ)},ξ+

(2)

式中,[B1(η)]和[B2(η)]為節(jié)點應變與位移的關(guān)系。

應力表示為

(3)

式中,[D]是彈性矩陣。

在比例邊界坐標中表示控制微分方程后,可以在徑向ξ上運用伽遼金法或虛功原理,得到基于位移的二維比例邊界有限元法的基本方程

[E0]ξ2{u(ξ)},ξξ+([E0]-[E1]+[E1]T)ξ{u(ξ)},ξ-

(4)

式中,當頻率ω=0,式(4)為靜力平衡方程。矩陣[E0]、[E1]、[E2]和[M0]是組裝邊界(ξ= 1)上計算的單元系數(shù)矩陣獲得,其形式為

(5)

(6)

(7)

(8)

值得說明的是,[E0]、[E1]和[E2]與常規(guī)有限元中靜力剛度陣相似,而[M0]類似于質(zhì)量陣。其中,[E0]、[M0]是正定矩陣,[E2]是半正定矩陣,而[E1]是非對稱。式(4)是Euler-Cauchy方程,可解析求解。

1.2 計算應力強度因子

求解式(4)后,可以獲得節(jié)點位移。通過位移場可以提取應力強度因子。位移場的漸進解表示為

{u(ξ)}=[V11]ξ-[S11]{c}+O{ξ2}

(9)

式中:{c}為邊界積分常數(shù);O{ξn}為一個高階項,在相同的ξn下趨于零,它包含高頻下慣性力的影響;[S]和[V]為式(4)轉(zhuǎn)化為一階常微分方程求解過程中Hamilton矩陣進行Schur分解時產(chǎn)生的特征值和特征向量[27];[S11]為包含所有非正特征值;[V11]為對應的特征值向量。代入到式(9)得

{u(ξ)}=[ψ(r)]{c(r)}+[ψ(s)]ξ-[S(s)]{c(s)}+

[ψ(T)]ξ{c(T)}+[ψn1]ξ-[Sn1]{cn1}+

O{ξ2}

(10)

位移的導數(shù)求解為

{u(ξ)},ξ=-[ψ(s)][S(s)]ξ-[S[s]]-[I]+

[ψ(T)]{c(T)}-[ψn1][Sn1]ξ-[Sn1]-[I]{cn1}+

O{ξ}

(11)

其中,積分常數(shù){c}表示為

{c}={{cn1};{c(T)};{c(s)};{c(r)}}

(12)

在一個線單元的指定局部坐標ξ處,應力場表示為

(13)

[B2(η)][ψ(s)]),

[B2(η)][ψ(T)])

(14)

奇異點周圍應力的漸近解表明,動力學中奇異應力和T-應力的定義與靜力學相同。圖2給出了基于SBFEM多邊形裂紋模型,其中笛卡爾坐標系和極坐標系的原點位于裂紋尖端。對于均質(zhì)板中的裂紋和各向同性雙材料板中的界面裂紋,應力強度因子統(tǒng)一定義為

(15)

圖2 裂紋區(qū)域的離散(比例中心位于裂紋尖端)

并且

(16)

式中:L為特征長度,應力強度因子的幅值與L無關(guān);ε為兩種材料的振蕩系數(shù),只與材料本身有關(guān)。

(17)

式中,r(θ)為從比例中心到沿角度(θ)徑向線的邊界的距離(圖2)。將ξ的矩陣冪函數(shù)改寫為極坐標形式

(18)

使用式(17),將奇異應力式(13)從局部坐標轉(zhuǎn)換為極坐標形式

(19)

根據(jù)式(15)中應力強度因子的形式,廣義應力強度因子表示為

(20)

通過式(20)可以很容易地計算任意(θ)角度下應力強度因子,并且對于確定裂紋擴展是非常重要的。

1.3 蒙特卡羅模擬

進行蒙特卡羅模擬之前,首先簡要介紹帶有不確定性參數(shù)的控制方程。設整個分析域為Ω,邊界為Γ(Γ=Γu+Γt),假定彈性模量E為不確定性參數(shù),即

{σ(E)}ij,j+{f}i-ρ{u(E)}i,tt-μ{u(E)}i,t=0

(21)

對于靜力學問題略去后兩項的慣性力和阻尼力,邊界上的約束條件為

(22)

(23)

式中,{ε}為應變。則應力為

{σ(E)}ij=[D(E)]ijkl{ε(E)}kl

(24)

使用伽遼金方法,運用分部積分和邊界條件,可以得到帶有隨機參數(shù)的弱形式

(25)

使用比例邊界有限元法進行離散,得到帶有隨機參數(shù)的代數(shù)方程

(26)

式中:g0gggggg為節(jié)點位移向量;{F}為荷載向量;[M]、[C]和[K]分別為全局質(zhì)量矩陣、阻尼矩陣和剛度矩陣。值得注意的是,這里假設彈性模量作為隨機參數(shù)存在于每個單元剛度矩陣中,用隨機場f(x)來描述彈性模量,對于每個單元有

(27)

式中:[B]為位移微分矩陣;J為雅可比矩陣。將式(25)中的[B]轉(zhuǎn)化為相應的系數(shù)矩陣[E],得到帶有不確定參數(shù)的比例邊界有限元方程。(對于靜力學問題,略去質(zhì)量陣和阻尼陣。)隨后基于蒙特卡羅模擬的思想,將隨機變量依次代入系統(tǒng)中進行求解,通過式(28)和(29)計算統(tǒng)計特征,以評估不確定參數(shù)對結(jié)構(gòu)響應的影響,即

(28)

(29)

式中:N為樣本個數(shù);{d(j)}為第j個樣本點的結(jié)構(gòu)響應。蒙特卡羅模擬的收斂速度為O(N-1/2),N越大,蒙特卡羅模擬收斂效果越好,計算成本也就越高。為了降低計算成本,提高計算效率。第2章將使用模型降階法,對直接的MCS進行加速求解。

2 奇異值分解

根據(jù)初始樣本的波動范圍,引入m個隨機變量并計算其響應,將這些樣本點的響應[d(αm)]按列排列,構(gòu)造快照矩陣[Ω][28]

[Ω]=[d(α1),d(α2),…,d(αm)]=

(30)

式中:[Ω]∈Rn×m;n為單個輸入變量下的響應個數(shù);m為變量個數(shù)。使用奇異值分解法(singular value decomposition,SVD)將矩陣[Ω]分解

(31)

式中:r=min(m,n);[U]∈Rn×n、[V]∈Rm×m為正交矩陣;uj和vj分別為[ΩΩT]和[ΩTΩ]的特征向量。[Σ]∈n×m是一個對角矩陣,對角元素σj按從大到小排列。定義[ψj]=uj,[ζj(αi)]=σjvij,式(31)可以改寫為

(32)

式中:[ψj]為正交基;[ζj(αi)]為相應的基底。使用式(32),系統(tǒng)響應可以通過簡化表達的[ψj]和[ζj(αi)]的線性組合來表示,且其規(guī)模比初始模型要小得多。為了實現(xiàn)對任意隨機輸入變量的系統(tǒng)響應的連續(xù)近似,應用徑向基函數(shù)(radial basis function,RBF)來對降階子空間中的基底進行插值

(33)

式中:φ為基函數(shù);χ為其權(quán)系數(shù)。本文選取高斯核函數(shù)作為基函數(shù),如下

(34)

(35)

值得注意的是,根據(jù)樣本點不同,上述方程組是無條件非奇異的。將式(33)代入式(32),式(32)改寫為

(36)

可以看出任何系統(tǒng)響應都可以通過使用這種線性組合來近似,無須通過SBFEM進行完整地計算。(當然,快照矩陣的形成仍然需要完整的計算)。圖3給出了本文算法的流程圖。與直接MCS相比,奇異值分解降低了系統(tǒng)的自由度,并提供了一個低階模型。RBF計算相關(guān)系數(shù),將子空間線性組合獲得新的響應。SVD-RBF提高了計算效率,降低了計算成本,這些優(yōu)勢一定程度上擴展了MCS的適用性。

圖3 算法的流程圖

3 數(shù)值算例

3.1 帶斜裂紋的平面板

圖4考慮一個在y方向,受拉力作用的帶有斜裂紋的平面板。其解析解[30]為

(37)

圖4 帶斜裂紋的平面板模型

式中,P為y方向施加的載荷。a是裂紋一半的長度;β是垂直方向與裂紋之間的夾角。彈性模量E=1 Pa,泊松比v=0.3,矩形板長度L=1 m。

表1給出了三種不同的變異系數(shù)下形成快照矩陣的樣本點數(shù)。并引入R2來評估所提出的方法的準確性,其計算公式為

(38)

表1 不同變異系數(shù)條件下本文算法結(jié)果的精度

為了更直觀地展示算法的精度,圖5給出了不同變異系數(shù)下的算法解和解析解的比較。從圖中可以看到,變異系數(shù)控制著變量波動的范圍(變異系數(shù)越大,數(shù)據(jù)波動越大)。在各種變異系數(shù)下,本文算法與解析解吻合較好。結(jié)合表1,可以得出變異系數(shù)越大,保持算法精度所需的初始樣本點也就越多。初始樣本點的數(shù)目決定著計算效率,所以平衡計算效率與計算精度是值得注意的。

圖5 不同變異系數(shù)的預測值與解析解的比較

隨機分析之前,會進行奇異值分解,此時產(chǎn)生了[Σ]矩陣。該矩陣只有主對角線上存在元素,并且從大到小排列。這些元素稱為奇異值,衰減迅速。奇異值越大包含的系統(tǒng)信息就越多,所以截斷[Σ]矩陣是必要的,以節(jié)省儲存空間,提升計算效率。仍使用算例1的初始樣本,并采取三種方案進行隨機性分析。

方案1:使用SBFEM計算1 000個隨機變量下的KI直接進行蒙特卡羅模擬。該計算需要重復1 000次。

方案2:對表1中的原始快照矩陣進行奇異值分解,截斷[Σ]矩陣,將[Σ]矩陣維數(shù)選擇為行與列的最小值。再進行蒙特卡羅模擬。

方案3:在方案2的基礎上,[Σ]矩陣維數(shù)選擇為σcat=10-1σmax。

圖6給出了KI在不同變異系數(shù)下的統(tǒng)計特征。從圖中可以看出方案1與方案2較為吻合,方案3具有一定的誤差。這是由于方案3截斷了過多的奇異值,導致降階系統(tǒng)“失真”。為了避免系統(tǒng)失真,保留降階優(yōu)勢,下面的算例都采用方案2的截斷方案。

圖6 不同方案下KI的統(tǒng)計特征

3.2 平面板中心彎曲裂紋

本節(jié)考慮一個中心帶有彎曲裂紋的平面板問題[32]。其模型如圖7所示,半徑和中心裂紋角由r和θ定義。本例中,假設r/w=0.1,取r=4.25 m,上下邊界

圖7 中心彎曲裂紋平面板

受均布荷載作用,計算采用平面應變假定。該問題具有解析解,應力強度因子為

(39)

以應力強度因子KI、KII和節(jié)點位移作為響應構(gòu)造快照矩陣,并對應力強度因子歸一化處理。選取初始裂紋角度作為隨機輸入變量,且服從期望μ= π/6,和標準差為σ=2的正態(tài)分布。設置分布參數(shù)后,使用隨機數(shù)生成器,隨機生成100個隨機角度變量,代入到解析解中進行計算,作為蒙特卡羅模擬的參考解。圖8給出了本文算法與解析解的應力強度因子KI,KII的對比。從圖中的結(jié)果來看,本文的算法與解析解的結(jié)果吻合得很好,表明算法具有較高的精度。值得注意的是,降階算法只采用了35個初始樣本點來形成快照矩陣,即只使用了35個樣本點,獲得了100個該變量下的響應結(jié)果,也說明了該算法在計算效率上具有優(yōu)勢。表2具體給出了該算法與基于SBFEM的蒙特卡羅模擬的統(tǒng)計特征和計算成本。表2中的快照矩陣由應力強度因子和98個節(jié)點位移構(gòu)成,自由度的計算方式為變量個數(shù)乘以模型自由度。從表2可以看出,本文算法降低了系統(tǒng)自由度,加速了傳統(tǒng)的基于SBFEM的蒙特卡羅模擬過程,提高了計算效率。

表2 SBFEM與SVD-RBF進行隨機分析的計算成本

圖8 應力強度因子的算法解與解析解的比較

3.3 帶孔洞的雙裂紋矩形板

本節(jié)考慮一個帶有圓孔的矩形板,設圓孔半徑r為不確定變量,研究其對動態(tài)應力強度因子(dynamic stress intensity factor,DSIF)的影響。矩形板的寬度為2b=30 mm,長度為2h=60 mm,包含一個半徑為r的孔,在圓孔處有兩條與水平方向呈30°夾角的等長度裂紋,如圖9(a)所示。其材料特性分別為:密度ρ=5×10-6kg/mm3、泊松比v=0.3,剪切模量G=76.923 GPa。t=0時,板的上下兩端受到均布的階躍荷載p(t)=PH(t)的作用,其中P是荷載的幅值。采用平面應變假定。

(a) 幾何尺寸

使用多邊形比例邊界有限元法分析時,將該模型劃分為8個多邊形單元,如圖9(b)所示。多邊形C1和C2是包含兩個裂紋的單元,并且將比例中心設置為裂紋的尖端。對于其他單元,都采用八節(jié)點的高階單元進行離散,用正方形來表示每個單元兩端的節(jié)點。網(wǎng)格劃分完畢后,進行隨機性分析。將圓孔半徑設為不確定變量,且滿足均值為μ=3.6,標準差為σ=0.2的正態(tài)分布。使用隨機數(shù)生成器產(chǎn)生50個滿足該分布的隨機變量,導入到模型中進行計算,獲得50個變量下的動態(tài)應力強度因子作為蒙特卡羅模擬的參考解。其中,計算動態(tài)應力強度因子的總時間為20 μs,計算時間步長Δt=0.1 μs。

圖10和11給出了降階算法與蒙特卡羅模擬得出的統(tǒng)計特征的對比。從圖中的結(jié)果來看,本文的算法與蒙特卡羅模擬結(jié)果吻合得很好。值得注意的是,降階算法只采用了17個初始樣本點來形成快照矩陣。意味著降階算法只使用了17個樣本點,獲得了50個蒙特卡羅模擬的結(jié)果。為了使數(shù)據(jù)更具說服性,圖12給出了r=3.75 mm時,算法的結(jié)果與參考解的對比。其中,這里r=3.75 mm是算法預測值,而非初始樣本點。從圖12可以看出,算法預測的結(jié)果與參考解基本保持一致。此外,快照矩陣還使用奇異值分解進行降維,進一步節(jié)省儲存空間。為了更清晰地說明算法的優(yōu)越性,表3給出了基于SBFEM分析斷裂問題的框架下,本文算法與直接的蒙特卡羅模擬的計算成本的對比。總的來說,本文算法提高了傳統(tǒng)蒙特卡羅模擬的計算效率,拓展了蒙特卡羅模擬的適用性。

表3 SBFEM與SVD-RBF進行隨機分析的計算成本

圖10 無量綱動應力強度因子KI的統(tǒng)計特征

圖11 無量綱動應力強度因子KII的統(tǒng)計特征

圖12 本文算法解與參考解的動應力強度因子對比

3.4 帶有裂紋的重力壩-地基系統(tǒng)的地震響應分析

本節(jié)研究巖石地基對大壩動力斷裂的影響。大壩的幾何尺寸如圖13所示。大壩材料為混凝土,其力學參數(shù)為:彈性模量Ec=31.0 GPa,泊松比vc=0.15,密度ρc=2 643 kg/m3;下部是巖石地基,部分力學參數(shù)vr=0.25,ρr=2 643 kg/m3。依照文獻[33],選取初始裂紋的高程為59.25 m,初始長度為1 m,輸入的地震動為Koyna地震波。設定巖基的彈性模量為不確定性參數(shù),其波動范圍為[10～15],步長為0.1選取50個Er作為蒙特卡羅模擬的結(jié)果。整個計算采用平面應變假定。

圖13 Koyna大壩幾何模型 (m)

采用多邊形單元進行離散,劃分93個多邊形單元458節(jié)點,其網(wǎng)格如圖14所示。其中截斷邊界采用11個三節(jié)點線單元離散,相似中心選取為壩基的中點處,采用基于連分式算法的高階透射邊界模擬壩基無限域。選取動態(tài)應力強度因子KI作為響應,計算時間為10 s,步長Δt=0.01 s。

圖14 SBFEM網(wǎng)格

圖15給出了本文算法和直接蒙特卡羅模擬獲得的統(tǒng)計特征的對比圖。從圖15可以看到,兩種方法得到的統(tǒng)計特征基本一致。值得注意的是,算法只使用了11個初始樣本點形成快照矩陣,且采用方案2的策略對快照矩陣進行降階處理。大大節(jié)省了儲存空間,提高了計算效率。為了進一步說明算法的準確性,圖16給出了Er=12 GPa算法的結(jié)果,與已有文獻結(jié)果[33]對比,這里的Er=12 GPa是算法預測得到結(jié)果而非初始樣本點。從圖16可以看出,算法解與參考解吻合得很好。KI的峰值發(fā)生在4 s處,在彈性模量Er的波動范圍內(nèi),峰值處KI從8 641 527.9 Pa·m0.5波動到7 023 890.08 Pa·m0.5相對波動范圍為18%,可見巖基參數(shù)的不確定性對動態(tài)應力強度因子有較大影響,因此研究不確定參數(shù)對結(jié)構(gòu)響應的影響是很有意義的。

圖15 歸一化動應力強度因子KI的統(tǒng)計特征

圖16 Er=12 GPa時本文算法解與參考解的結(jié)果對比

4 結(jié) 論

本文構(gòu)建了一種新的使用比例邊界有限元法計算應力強度因子的蒙特卡羅模擬框架,并得到如下結(jié)論:

(1) 采用多邊形比例邊界有限元法分析斷裂問題時,將比例中心設置在裂紋尖端,使裂紋尖端不需要網(wǎng)格細化,可直接提取應力強度因子,為不確定量化提供了高精度的計算模型。

(2) 采用蒙特卡羅模擬計算響應的統(tǒng)計特征,量化不確定參數(shù)對結(jié)構(gòu)響應的影響,并使用模型降階法對傳統(tǒng)的蒙特卡羅模擬進行改進,直接降低物理問題的自由度,快速獲得新的結(jié)構(gòu)響應,并討論了奇異值矩陣的合理截斷。同時,考慮了不同變異系數(shù)對結(jié)構(gòu)響應的影響,對大的變異系數(shù)(0.4)的結(jié)果也有較高的精度。

(3) 數(shù)值算例表明,對于不同的變量和不同的物理模型,使用SVD-RBF加速MCS具有良好的精度。

并且能顯著降低計算成本,提升了傳統(tǒng)蒙特卡羅模擬的適用性。