盧夢凱 張洪武 鄭勇剛

(大連理工大學(xué)，國際計(jì)算力學(xué)研究中心，工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室，運(yùn)載工程與力學(xué)學(xué)部工程力學(xué)系，大連116024)

-固體力學(xué)

應(yīng)變局部化分析的嵌入強(qiáng)間斷多尺度有限元法1)

盧夢凱 張洪武 鄭勇剛2)

(大連理工大學(xué)，國際計(jì)算力學(xué)研究中心，工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室，運(yùn)載工程與力學(xué)學(xué)部工程力學(xué)系，大連116024)

固體材料的應(yīng)變局部化行為是導(dǎo)致結(jié)構(gòu)破壞失效的重要因素之一，開展相關(guān)數(shù)值模擬分析對于結(jié)構(gòu)安全性評估具有重要意義.然而由于材料的非均質(zhì)和多尺度特性，采用傳統(tǒng)數(shù)值方法進(jìn)行求解時(shí)通常需要從最小特征尺度離散求解的結(jié)構(gòu)，這將大幅度增加計(jì)算規(guī)模和成本.針對這一問題，本文提出了一種基于嵌入強(qiáng)間斷模型的多尺度有限元方法.該方法從粗細(xì)兩個(gè)尺度離散求解模型，首先在細(xì)尺度單元上引入嵌入強(qiáng)間斷模型來描述單元間斷特性，所附加的跳躍位移自由度則通過凝聚技術(shù)進(jìn)行消除，從而保持細(xì)尺度單元?jiǎng)偠汝嚲S度不變.其次，提出了一種增強(qiáng)多節(jié)點(diǎn)粗單元技術(shù)，其可根據(jù)局部化帶與粗單元邊界相交情況自適應(yīng)動(dòng)態(tài)地增加粗節(jié)點(diǎn)，新構(gòu)造的增強(qiáng)數(shù)值基函數(shù)可以捕捉細(xì)尺度間斷特性，完成物理信息從細(xì)單元到粗單元的準(zhǔn)確傳遞以及宏觀響應(yīng)的快速分析；再次，在細(xì)尺度解的計(jì)算中，將細(xì)尺度解分解為降尺度解與單胞局部攝動(dòng)解，從而消除彈塑性分析時(shí)單胞內(nèi)部的不平衡力.最后，通過兩個(gè)典型算例分析，并與完全采用細(xì)單元的嵌入有限元結(jié)果進(jìn)行對比，驗(yàn)證了所提出算法的正確性與有效性.

多尺度有限元法，嵌入強(qiáng)間斷模型，應(yīng)變局部化，增強(qiáng)粗單元，彈塑性分析

引言

應(yīng)變局部化作為一種典型的破壞形式廣泛存在于諸多工程材料與結(jié)構(gòu)中，如金屬拉伸頸縮、地基基礎(chǔ)破壞、山體滑坡等.局部化問題的特點(diǎn)通常呈現(xiàn)為很大的塑性應(yīng)變局限在一個(gè)狹小的細(xì)帶內(nèi)，隨著局部化帶的演化擴(kuò)展，最終導(dǎo)致結(jié)構(gòu)的失穩(wěn)破壞.針對應(yīng)變局部化問題的數(shù)值分析研究對預(yù)測材料與結(jié)構(gòu)的承載力和穩(wěn)定性具有重要的指導(dǎo)意義.

目前，已有許多計(jì)算方法用來模擬材料或結(jié)構(gòu)的失效行為[16].模擬應(yīng)變局部化或者斷裂問題的方法必須考慮解的客觀性，即網(wǎng)格無關(guān)性.而基于經(jīng)典損傷或者塑性模型的方法，往往存在網(wǎng)格依賴性問題[7].為了解決這一問題，通常需要在模型中引入正則化機(jī)制.目前廣泛采用的正則化模型主要分為兩大類，其中一類是彌散型的正則化模型，主要包括非局部模型[89]、Cosserat模型[10]、率相關(guān)模型[11]和梯度塑性模型[1213]等.另一類則是離散型的模型，如擴(kuò)展有限元(XFEM)[1,1416]、嵌入有限元(EFEM)[1722]、數(shù)值流形方法[23]等.其中EFEM是基于Simo等[24]所發(fā)展的假設(shè)增強(qiáng)應(yīng)變(AES)框架提出的，其通過引入嵌入強(qiáng)間斷模型來提供耗散機(jī)制[17].另外，EFEM分析時(shí)增加的單元附加自由度可以通過凝聚技術(shù)消除，由此可以保證單元?jiǎng)偠汝嚲S數(shù)不變，便于程序?qū)崿F(xiàn).目前EFEM已被廣泛應(yīng)用于應(yīng)變局部化和斷裂破壞問題的模擬分析中.然而，隨著工程分析要求的日益提高，尤其是針對大規(guī)模、非均質(zhì)問題時(shí)，需要采用非常精細(xì)的網(wǎng)格對空間進(jìn)行離散，這將導(dǎo)致計(jì)算量的大幅增大，從而限制了EFEM的大規(guī)模工程應(yīng)用.

多尺度方法作為一種可有效解決此類問題的數(shù)值計(jì)算方法，已經(jīng)成為近年來研究的熱點(diǎn)，如計(jì)算均勻化方法[25]、兩尺度漸近分析方法[26]、多尺度有限元方法[27]等.其中多尺度有限元方法由于其宏微觀尺度之間信息傳遞方便準(zhǔn)確而受到了廣泛關(guān)注.多尺度有限元方法最早由Hou等[27]提出，其主要特點(diǎn)是引入了數(shù)值基函數(shù)的概念，從而可以方便地進(jìn)行升尺度與降尺度計(jì)算.多尺度有限元法首先被用于多孔介質(zhì)滲流分析[2829].Zhang等[30]將該方法推廣到固體力學(xué)領(lǐng)域，并在許多問題中得到了應(yīng)用，如動(dòng)力問題[3132]、熱力電耦合問題[33]、桁架非線性問題[34]、多相多孔介質(zhì)耦合問題[32,35]和裂紋擴(kuò)展問題[36]等，但其在連續(xù)體應(yīng)變局部化的多尺度問題方面，還鮮有報(bào)道.

針對單相固體的應(yīng)變局部化問題，本文提出了一種基于嵌入強(qiáng)間斷模型的多尺度有限元方法.首先介紹了嵌入強(qiáng)間斷模型的基本理論，包括控制方程、物理場描述、弱形式控制方程以及本構(gòu)關(guān)系.其次，詳細(xì)說明了所提出算法的基本思想和計(jì)算框架，主要包括細(xì)單元尺度上間斷特性表征方式、升尺度計(jì)算中增強(qiáng)數(shù)值基函數(shù)構(gòu)造方法、自適應(yīng)多節(jié)點(diǎn)粗單元技術(shù)、降尺度計(jì)算中處理單胞內(nèi)部不平衡力的局部攝動(dòng)方法以及詳細(xì)的算法實(shí)現(xiàn)流程.最后，用所提方法計(jì)算了兩個(gè)局部化算例，與基于精細(xì)網(wǎng)格的EFEM解的對比結(jié)果表明了多尺度有限元方法的正確性與有效性.

1 嵌入強(qiáng)間斷模型的基本理論

1.1 控制方程

本文主要研究二維問題.如圖1所示，考慮一個(gè)二維平面有界區(qū)域Ω，被間斷線Γ劃分為區(qū)域Ω+和Ω-.描述小變形下固體強(qiáng)間斷模型的控制方程可以寫為

式中,?為Nabla算子，σ為應(yīng)力，ρ表示密度，b為體力向量，u為位移，upb為邊界?uΩ上的指定位移，tpb為邊界?tΩ上的指定力，nt為?tΩ的外法向向量，n為Γ的法向向量，[[σ]]=σ+-σ-表示間斷應(yīng)力，σ+和σ-分別為間斷線正向區(qū)域的應(yīng)力和背向區(qū)域的應(yīng)力，[[·]]代表間斷運(yùn)算符，tΓ代表間斷線上的力矢量.

圖1 強(qiáng)間斷模型示意圖Fig.1 Illustration of the strong discontinuity model

1.2 物理場描述

在強(qiáng)間斷模型中，為了描述位移的間斷特性，可將位移場u表示為

由式(6)可得應(yīng)變場

1.3 弱形式控制方程

對控制方程(1)運(yùn)用虛功原理，同時(shí)考慮邊界條件(2)～條件(5)，并運(yùn)用分部積分法與非連續(xù)函數(shù)的散度定理(F代表不連續(xù)函數(shù))，即

可以推得式(1)的弱形式(忽略重力項(xiàng))

1.4 本構(gòu)方程

在嵌入強(qiáng)間斷模型中，連續(xù)區(qū)域內(nèi)的應(yīng)力σ與間斷區(qū)域內(nèi)的力tΓ的演化過程分別由以下兩種本構(gòu)模型確定：

(1)連續(xù)型Drucker-Prager塑性本構(gòu)模型

在ΩΓ區(qū)域，采用屈服函數(shù)與塑性勢函數(shù)相同的關(guān)聯(lián)Drucker-Prager塑性模型，其屈服函數(shù)可以寫作

式中，s是σ的偏量部分；q=-σY-HCp，HC為塑性模量，p為等效塑性應(yīng)變，σY為初始屈服應(yīng)力，β為Drucker-Prager模型的材料參數(shù)，本文取為

其中,φ為摩擦角.根據(jù)塑性流動(dòng)理論，可以得到塑性應(yīng)變率與內(nèi)變量的演化方程

式中，λ為塑性乘子.基于加法分解下的應(yīng)變率可以記為彈性應(yīng)變率與塑性應(yīng)變率之和，即

針對連續(xù)型塑性本構(gòu)的積分算法，本文采用經(jīng)典的向后歐拉返回映射算法，具體可參見文獻(xiàn)[37-38].

(2)局部化型Drucker-Prager塑性本構(gòu)

在連續(xù)型塑性本構(gòu)里，塑性乘子λ是正則的，而在間斷Γ區(qū)域，塑性乘子λ表現(xiàn)為奇異性，引入

結(jié)合式(7)～式(9)和式(18)，可以得到如下應(yīng)力率的表達(dá)式

式中，Dsk為四階彈性本構(gòu)張量.結(jié)合應(yīng)力率有界的條件，可進(jìn)一步得到

在局部坐標(biāo)n-t下(如圖1所示)，向量tΓm與可以表示為

式中,tΓn= tΓm·n，tΓt= tΓm·t，T = [n t] 和其中χ=tanφ，且根據(jù)式(17)，可得

將式(26)～式(28)代入式(25)，可得間斷Γ區(qū)域內(nèi)的本構(gòu)關(guān)系

2 基于嵌入強(qiáng)間斷模型的多尺度有限元法

本節(jié)中，我們結(jié)合嵌入強(qiáng)間斷模型提出了一種用于處理單相固體局部化問題的多尺度有限元方法.其主要思想是將幾何區(qū)域離散為粗單元與細(xì)單元，在細(xì)單元上用強(qiáng)間斷模型描述間斷特性，再通過構(gòu)造增強(qiáng)的數(shù)值基函數(shù)，正確地將間斷特性從細(xì)尺度傳遞到粗尺度上，從而進(jìn)行局部化問題的多尺度有限元分析.以下給出了針對單相固體結(jié)合嵌入強(qiáng)間斷模型的多尺度有限元方法基本框架，主要包括細(xì)尺度上的有限元離散、增強(qiáng)數(shù)值基函數(shù)的構(gòu)造、粗尺度上的求解以及細(xì)尺度解的計(jì)算.

2.1 細(xì)尺度上的有限元離散

在細(xì)尺度上采用八節(jié)點(diǎn)等參單元進(jìn)行離散，其位移可由節(jié)點(diǎn)位移插值得到，即

式中，ue與ξe分別為節(jié)點(diǎn)位移與間斷跳躍量，ˉN為標(biāo)準(zhǔn)形函數(shù)，JΓ為間斷形函數(shù)，可寫為

式中，應(yīng)變算子ˉB=?symˉN，協(xié)調(diào)算子G和投射算子G?分別為

圖2 局部化單元內(nèi)角點(diǎn)集合的選取Fig.2 Selection for the corner nodes setin a discontinuous element

對控制方程弱形式(14)和式(15)進(jìn)行Galerkin離散可得

進(jìn)一步對式(36)和式(37)進(jìn)行線性化，并采用Newton-Raphson增量求解策略，可得到以節(jié)點(diǎn)位移增量δue和單元間斷量δξe,t作為未知量的離散控制方程[17,21].

式中各矩陣具體表達(dá)式為

2.2 增強(qiáng)數(shù)值基函數(shù)的構(gòu)造

考慮一個(gè)帶有間斷區(qū)域的多尺度單胞(如圖3(a)所示)，其表示的粗單元在形式上也是一種多節(jié)點(diǎn)粗單元，但與常規(guī)多節(jié)點(diǎn)粗單元中一般采用等間隔布置粗節(jié)點(diǎn)的方式不同，其增加的粗節(jié)點(diǎn)位置與個(gè)數(shù)是根據(jù)間斷ΓE與單胞邊界的交點(diǎn)來確定的.如圖3(a)中的粗節(jié)點(diǎn)9,10和11即為分析過程中隨著間斷區(qū)域演化而動(dòng)態(tài)增加的粗節(jié)點(diǎn).該類粗單元所構(gòu)造的數(shù)值基函數(shù)稱為增強(qiáng)數(shù)值基函數(shù)，其可以正確地傳遞細(xì)尺度單元上所描述的間斷特性，從而保證了局部化問題多尺度分析時(shí)的有效性.

圖3 (a)多節(jié)點(diǎn)粗單元，(b)～(d)分別為1號、9號和5號粗節(jié)點(diǎn)的邊界指定位移分布圖Fig.3(a)a multi-node coarse element,the distributions of the prescribed boundary values for(b)coarse node 1,(c)coarse node 9 and(d)coarse node 5

增強(qiáng)數(shù)值基函數(shù)的構(gòu)造方式同常規(guī)的多節(jié)點(diǎn)粗單元類似.記粗單元對應(yīng)的單胞整體等效剛度陣為為剛度陣組裝算子，Sne為單胞內(nèi)細(xì)單元數(shù)目.則構(gòu)造j號粗節(jié)點(diǎn)的k向自由度對應(yīng)的數(shù)值基函數(shù)Njk時(shí)，需滿足下式

式中,ψjk是Njk在邊界?ΩE上的指定值.本文采用多節(jié)點(diǎn)線性邊界條件來進(jìn)行構(gòu)造，圖3(b)給出了構(gòu)造N1k時(shí)ψ1k的分布，即當(dāng)k為ux時(shí)，邊界上細(xì)節(jié)點(diǎn)的uy都約束為0，ux在線段19和18上施加峰值為1的線性分布位移；當(dāng)k為uy時(shí)，邊界上細(xì)節(jié)點(diǎn)的ux都約束為0，uy在線段19和18上施加峰值為1的線性分布位移.ψ9k和ψ5k的分布由圖3(c)和圖3(d)所示.類似地，可以獲得整個(gè)粗單元下所對應(yīng)邊界條件ψ的分布，并求解得到增強(qiáng)數(shù)值基函數(shù)N.

至此，可以建立起細(xì)單元節(jié)點(diǎn)自由度與粗單元節(jié)點(diǎn)自由度之間的關(guān)系如下

式中，n是單胞內(nèi)的細(xì)節(jié)點(diǎn)數(shù)目，δue和δuE分別為細(xì)單元與粗單元上的節(jié)點(diǎn)增量位移向量.在實(shí)際計(jì)算時(shí)，若出現(xiàn)剪切帶與粗單元邊界相交時(shí)，使擁有該邊界的相鄰粗單元同時(shí)成為增強(qiáng)粗單元，并共享增強(qiáng)節(jié)點(diǎn)，從而保證相鄰粗單元間位移的協(xié)調(diào).

2.3 粗尺度上的求解

在獲得了增強(qiáng)數(shù)值基函數(shù)N之后，就可以用其形成單胞的等效剛度矩陣與右端不平衡力向量，可以分別表示為

其中,A為組裝單元矩陣或向量的集成算子[29].粗尺度下的整體剛度矩陣與右端不平衡力向量可以寫作

其中,Cne代表粗單元數(shù)目.最終，粗尺度下的增量位移解可由求得.

2.4 細(xì)尺度解的計(jì)算

在獲得整體宏觀解δuG后，可以得到每個(gè)單胞的宏觀解δuE，再通過降尺度計(jì)算來獲得降尺度解而在多尺度有限元的非線性分析中，其細(xì)尺度解δuuc可以表示為降尺度解與局部攝動(dòng)解之和，即

3 算法實(shí)現(xiàn)

本節(jié)中，我們給出了所提出算法的具體實(shí)現(xiàn)步驟.

(一)將幾何體離散為細(xì)單元與粗單元，設(shè)置材料參數(shù)、算法參數(shù)等.

(二)計(jì)算每個(gè)單胞的等效剛度陣，計(jì)算數(shù)值基函數(shù).

(三)對時(shí)間步n進(jìn)行循環(huán).

(1)n=n+1.

(2)對迭代步k進(jìn)行循環(huán).

①k=k+1.

②重新計(jì)算那些細(xì)單元?jiǎng)偠汝嚫淖兊膯伟麛?shù)值基函數(shù)，并根據(jù)式(51)和式(52)分別形成新的和

(a)對于 ΩΓ區(qū)域，采用連續(xù)型 Drucker-Prager塑性本構(gòu)關(guān)系，用經(jīng)典的向后歐拉返回映射算法進(jìn)行更新[37-38].

(b)對于Γ區(qū)域，采用局部型Drucker-Prager塑性本構(gòu)關(guān)系，用非標(biāo)準(zhǔn)的返回映射算法進(jìn)行更新[19].

⑤判斷是否滿足Eu＜或k達(dá)到指定的最大迭代數(shù)的條件.如果滿足，則進(jìn)入步驟(四)；否則，返回步驟①.

(四)用下式判斷是否有細(xì)單元進(jìn)入局部化[19]

式中，σ1和σ2分別為第一和第二主應(yīng)力，ˉEL為算法參數(shù).如果有新單元進(jìn)入局部化，則運(yùn)行間斷線跟蹤算法[39]，并動(dòng)態(tài)增加所在粗單元的節(jié)點(diǎn)，計(jì)算增強(qiáng)數(shù)值基函數(shù)；否則，直接進(jìn)入步驟(五).

(五)判斷n是否為最大時(shí)間步.如果是，則進(jìn)入步驟(六)；否則，返回步驟(1).

(六)計(jì)算結(jié)束，輸出結(jié)果.

在步驟⑤中，收斂判斷準(zhǔn)則Eu定義為

式中,N為細(xì)網(wǎng)格節(jié)點(diǎn)總數(shù)(實(shí)際計(jì)算時(shí)去掉位移為零的自由度).

4 算例分析

本節(jié)用所提出算法來模擬二維平面應(yīng)變下的單相固體應(yīng)變局部化問題，并考慮不同數(shù)目粗節(jié)點(diǎn)的粗單元的計(jì)算效力，所得結(jié)果與采用細(xì)單元的EFEM參考解進(jìn)行對比.此外，兩個(gè)算例均采用直接位移控制算法來控制加載，收斂精度參數(shù)ˉEu和ˉEL分別取為10-6和10-3.

4.1 土柱的豎直壓縮

模型尺寸為0.04m×0.14m，底部節(jié)點(diǎn)全約束，在上表面施加豎直向下位移0.1mm.將整個(gè)模型離散為4×14個(gè)單胞，每個(gè)單胞離散為8×8個(gè)細(xì)單元.采用4種不同節(jié)點(diǎn)的粗單元：CE4,CE8,CE4E和CE8E(CEnE中CE代表粗單元，n代表粗單元節(jié)點(diǎn)個(gè)數(shù)，E代表采用增強(qiáng)數(shù)值基函數(shù)技術(shù))對所提出算法進(jìn)行考察，如圖4所示.表1給出了具體的材料參數(shù)，在分析中采用50個(gè)載荷步進(jìn)行加載.

圖4 土柱壓縮算例示意圖Fig.4 Illustration of the compression test of soil column

表1 土柱材料參數(shù)Table 1 Material parameters for the column

圖5給出了4種類型粗單元的多尺度有限元解與EFEM參考解的作用反力與加載位移變化曲線.從圖中可以看出，采用常規(guī)粗單元CE4與CE8的曲線并沒與呈現(xiàn)出結(jié)構(gòu)軟化現(xiàn)象，而采用CE4E與CE8E的結(jié)果則與參考解吻合很好.該結(jié)果表明增強(qiáng)數(shù)值基函數(shù)可以正確地捕捉細(xì)單元上的間斷特性，從而獲得有效的計(jì)算結(jié)果.圖6給出了CE4E,CE8E和參考解的豎直位移云圖.3種方法得到的結(jié)果吻合良好.另一方面，針對不同方法的CPU花費(fèi)時(shí)間比CE4E/EFEM和CE8E/EFEM分別為0.37和0.42，這也表明了多尺度有限元方法在計(jì)算效率上的優(yōu)勢.需要進(jìn)一步指出的是，對于大規(guī)模工程問題而言，其應(yīng)變局部化區(qū)域相對于整體結(jié)構(gòu)而言比本算例更小，本文所發(fā)展的多尺度有限元方法計(jì)算效率將會(huì)大幅度提高.

圖5 不同方法作用反力與加載位移的結(jié)果比較Fig.5 Results of the reaction force versus the imposed displacement with di ff erent methods

圖6 不同方法豎直位移的云圖Fig.6 Contours of the vertical displacement with di ff erent methods

4.2 基礎(chǔ)的水平加載

本算例模擬了基礎(chǔ)在水平位移作用下的局部化問題.模型尺寸為10m×10m，底部節(jié)點(diǎn)豎直方向位移約束，左側(cè)節(jié)點(diǎn)位移全約束，在右側(cè)剛性塊中心處作用水平方向的位移載荷0.15m.整個(gè)基礎(chǔ)離散為10×10個(gè)單胞，每個(gè)單胞離散為6×6個(gè)細(xì)單元，采用兩種不同節(jié)點(diǎn)的粗單元：CE4E和CE8E，具體見圖7所示.表2列出了具體的材料參數(shù).水平載荷設(shè)置成50個(gè)載荷步進(jìn)行加載.

圖7 基礎(chǔ)水平加載算例示意圖Fig.7 Illustration of the horizontal loading test of the foundation

表2 基礎(chǔ)材料參數(shù)Table 2 Material parameters for the foundation

圖8給出了采用CE4E、CE8E和EFEM作用反力與水平位移的變化曲線.從圖中可以發(fā)現(xiàn)，兩種粗單元CE4E與CE8E的結(jié)果都與EFEM結(jié)果吻合非常好，且CE8E的結(jié)果較CE4E更加接近于EFEM參考解.該結(jié)果表明隨著粗節(jié)點(diǎn)的增加，可以有效地提高解的精度.圖9給出了不同時(shí)刻(A～D分別對應(yīng)圖8不同時(shí)刻)下基于CE8E預(yù)測獲得的局部化帶擴(kuò)展圖.圖10給出了CE4E、CE8E和參考解的水平位移云圖.從圖中局部化帶的位置與形狀分析，可以發(fā)現(xiàn)3種方法得到的結(jié)果吻合良好.另外，CE4E/EFEM和CE8E/EFEM的CPU花費(fèi)時(shí)間比分別為0.40和0.47.

圖8 不同方法作用反力與加載位移的結(jié)果比較Fig.8 Results of the reaction force versus the imposed displacement with di ff erent methods

圖9 CE8E獲得的局部化帶擴(kuò)展圖(A～D對應(yīng)于圖8中的不同時(shí)刻，此圖中變形放大了5倍)Fig.9 Illustrations of localization band propagation with the CE8E(A～D correspond to di ff erent times in Fig.8,the deformation is scaled by 5 for clarity)

圖10 不同方法水平位移的云圖Fig.10 Contours of the horizontal displacement with di ff erent methods

5 結(jié)論

本文提出了一種針對單相固體應(yīng)變局部化分析的多尺度有限元方法.首先，在細(xì)單元內(nèi)引入嵌入強(qiáng)間斷模型，并通過凝聚技術(shù)在單元級別上消除間斷自由度，使得整體多尺度算法流程可以保持不變且易于程序?qū)崿F(xiàn)；其次，對于具有局部化細(xì)單元的單胞提出了一種增強(qiáng)的多節(jié)點(diǎn)粗單元，其可以動(dòng)態(tài)地根據(jù)間斷位置來增加粗節(jié)點(diǎn)個(gè)數(shù)，構(gòu)造的增強(qiáng)數(shù)值基函數(shù)可以正確捕捉細(xì)尺度上的間斷特性；然后，在細(xì)尺度解的計(jì)算中，應(yīng)用局部攝動(dòng)法來消除單胞內(nèi)部的不平衡力.最后，通過兩個(gè)應(yīng)變局部化算例，驗(yàn)證了所提出算法的正確性與有效性.此外，該方法可以較容易地推廣用于處理多場多相強(qiáng)間斷問題.

1 Belytschko T,Black T.Elastic crack growth in finit elements with minimal remeshing.International Journal for Numerical Methods in Engineering,1999,45(5)：601-620

2 Schrefle BA,Secchi S,Simoni L.On adaptive refinemen techniques in multi-fielproblems including cohesive fracture.Computer Methods in Applied Mechanics and Engineering,2006,195(4-6)：444-461

3 Miehe C,Mauthe S.Phase fielmodeling of fracture in multiphysics problems.Part III.Crack driving forces in hydro-poroelasticity and hydraulic fracturing of fluid-saturate porous media.Computer Methods in Applied Mechanics and Engineering,2016,304：619-655

4 Chen Z,Hu W,Shen L,et al.An evaluation of the MPM for simulating dynamic failure with damage di ff usion.Engineering Fracture Mechanics,2002,69(17)：1873-1890

5 Yang PF,Gan Y,Zhang X,et al.Improved decohesion modeling with the material point method for simulating crack evolution.International Journal of Fracture,2014,186(1-2)：177-184

6 Boyce BL,Kramer SLB,Fang HE,et al.The Sandia Fracture Challenge：blind round robin predictions of ductile tearing.International Journal of Fracture,2014,186(1-2)：5-68

7 Bazant ZP,Belytschko T.Wave propagation in a strain-softening bar：Exact solution.Journal of Engineering Mechanics,1985,111：381-389

8 Pijaudier-Cabot G,Bazant ZP.Nonlocal damage theory.Journal of Mechanical Engineering,1987,113：1512-1533

9 Wol ffC,Richart N,Molinari JF.A non-local continuum damage approach to model dynamic crack branching.International Journal for Numerical Methods in Engineering,2015,101(12)：933-949

10李錫夔,唐洪祥.壓力相關(guān)彈塑性Cosserat連續(xù)體模型與應(yīng)變局部化有限元模擬.巖石力學(xué)與工程學(xué)報(bào),2005,24(9)：1497-1505(Li Xikui,Tang Hongxiang.Pressure-dependent elasto-plastic Cosserat continuum model and finit element simulation of strain localization.Chinese Journal of Rock Mechanics and Engineering,2005,24(9)：247-253(in Chinese))

11 Needleman A.Material rate dependence and mesh sensitivity in localization problems.Computer Methods in Applied Mechanics and Engineering,1988,67(1)：69-85

12 de Borst R,M¨uhlhaus HB.Gradient-dependent plasticity：Formulation and algorithmic aspects.International Journal for Numerical Methods in Engineering,1992,35(3)：521-539

13 Peerlings RHJ,de Borst R,Brekelmans WAM,et al.Gradient enhanced damage for quasi-brittle materials.International Journal for Numerical Methods in Engineering,1996,39(19)：3391-3403

14 Dolbow J,Mo¨es N,Belytschko T.Discontinuous enrichment in fi nite elements with a partition of unity method.Finite Elements in Analysis and Design,2000,36(3-4)：235-260

15莊茁,成斌斌.發(fā)展基于CB殼單元的擴(kuò)展有限元模擬三維任意擴(kuò)展裂紋.工程力學(xué),2012,29(6)：12-21(Zhuang Zhuo,Cheng Binbin.Development of X-FEM on CB shell element for simulating 3D arbitrary crack growth.Engineering Mechanics,2012,29(6)：12-21(in Chinese))

16石路楊,余天堂.多裂紋擴(kuò)展的擴(kuò)展有限元法分析.巖土力學(xué),2014,35(1)：263-272(Shi Luyang,Yu Tiantang.Analysis of multiple crack growth using extended finit element method.Rock and Soil Mechanics,2014,35(1)：263-272(in Chinese))

17 Simo JC,Oliver J,Armero F.An analysis of strong discontinuities induced by strain-softening in rate-independent inelastic solids.Computational Mechanics,1993,12(5)：277-296

18 Armero F.On the characterization of localized solutions in inelastic solids：an analysis of wave propagation in a softening bar.Computer Methods in Applied Mechanics and Engineering,2001,191(3-5)：181-213

19 Borja RI,Regueiro RA.Strain localization in frictional materials exhibiting displacement jumps.Computer Methods in Applied Mechanics and Engineering,2001,190(20-21)：2555-2580

20 Benkemoun N,Gelet R,Roubin E,et al.Poroelastic two-phase material modeling：theoretical formulation and embedded finit element method implementation.International Journal for Numerical and Analytical Methods in Geomechanics,2015,39(12)：1255-1275

21 Linder C,Armero F.Finite elements with embedded strong discontinuities for the modeling of failure in solids.International Journal for Numerical Methods in Engineering,2007,72(12)：1391-1433

22鄭利濤,胡志強(qiáng),唐洪祥.強(qiáng)間斷分析方法在土工結(jié)構(gòu)物漸進(jìn)破壞過程中的應(yīng)用.巖土力學(xué),2012,33(9)：2771-2780(Zheng Litao,Hu Zhiqiang,Tang Hongxiang.Application of strong discontinuity analysis to progressive failure process of geotechnical structures.Rock and Soil Mechanics,2012，33(9)：2771-2780(in Chinese))

23徐棟棟,鄭宏,楊永濤等.多裂紋擴(kuò)展的數(shù)值流形法.力學(xué)學(xué)報(bào),2015,47(3)：471-481(Xu Dongdong,Zheng Hong,Yang Yongtao et al.Multiple crack propagation based on the numerical manifold method.Chinese Journal of Theoretical and Applied Mechanics,2015,47(3)：471-481(in Chinese))

24 Simo JC,Rifai MS.A class of mixed assumed strain methods and the method of incompatible modes.International Journal for Numerical Methods in Engineering,1990,29(8)：1595-1638

25 Kouznetsova V,Geers MGD,Brekelmans WAM.Multi-scale constitutive modelling of heterogeneous materials with a gradientenhanced computational homogenization scheme. International Journal for Numerical Methods in Engineering,2002,54(8)：1235-1260

26 Yang ZQ,Cui JZ,Sun Y,et al.Multiscale computation for transient heat conduction problem with radiation boundary condition in porous materials.Finite Elements in Analysis and Design,2015,102：7-18

27 Hou TY,Wu XH.A multiscale finit element method for elliptic problemsincompositematerialsandporousmedia.JournalofComputational Physics,1997,134(1)：169-189

28 Dostert P,Efendiev Y,Hou TY.Multiscale finit element methods for stochastic porous media fl w equations and application to uncertainty quantification Computer Methods in Applied Mechanics and Engineering,2008,197(43-44)：3445-3455

29 Efendiev Y,Hou TY.Multiscale Finite Element Methods.[s.l.]：Springer Science&Business Media,2009

30 Zhang HW,Wu JK,Fu ZD.Extended multiscale finitelement method for elasto-plastic analysis of 2D periodic lattice truss materials.Computational Mechanics,2010,45(6)：623-635

31卓小翔,劉輝,楚錫華等.非均質(zhì)材料動(dòng)力分析的廣義多尺度有限元法.力學(xué)學(xué)報(bào),2016,48(2)：378-386(Zhuo Xiaoxiang,Liu Hui,Chu Xihua,et al.A generalized multiscale finit element method for dynamic analysis of heterogeneous material.Chinese Journal ofTheoretical and Applied Mechanics,2016,48(2)：378-386(in Chinese))

32張洪武,盧夢凱,鄭勇剛.非均質(zhì)飽和多孔介質(zhì)彈塑性動(dòng)力分析的廣義耦合擴(kuò)展多尺度有限元法.計(jì)算力學(xué)學(xué)報(bào),2016,33(4)：454-461(Zhang Hongwu,Lu Mengkai,Zheng Yonggang.General coupling extended multi-scale finit element method for the elastoplastic dynamic analysis of heterogeneous saturated porous media.Chinese Journal of Computational Mechanics,2016,33(4)：454-461(in Chinese))

33 Lv J,Yang K,Zhang HW,et al.A hierarchical multiscale approach for predicting thermo-electro-mechanical behavior of heterogeneous piezoelectric smart materials.Computational Materials Science,2014,87：88-99

34 Liu H,L¨u J.An equivalent continuum multiscale formulation for 2D geometrical nonlinear analysis of lattice truss structure.Composite Structures,2017,160：335-348

35 Zhang HW,Lu MK,Zheng YG,et al.General coupling extended multiscale FEM for elasto-plastic consolidation analysis of heterogeneous saturated porous media.International Journal for Numerical and Analytical Methods in Geomechanics,2015,39(1)：63-95

36 Lu MK,Zhang HW,Zheng YG,et al.A multiscale finit element method with embedded strong discontinuity model for the simulation of cohesive cracks in solids.Computer Methods in Applied Mechanics and Engineering,2016,311：576-598

37 Simo JC,Hughes TJR.Computational Inelasticity.[s.l.]：Springer Science&Business Media,2006

38 de Souza Neto EA,Peric D,Owen DRJ.Computational Methods for Plasticity.New Jersey：John Wiley&Sons,2011

39 ParvanehSM,FosterCD.Onnumericalaspectsofdi ff erentupdating schedules for tracking fracture path in strain localization modeling.Engineering Fracture Mechanics,2016;152：26-57

EMBEDDED STRONG DISCONTINUITY MODEL BASED MULTISCALE FINITE ELEMENT METHOD FOR STRAIN LOCALIZATION ANALYSIS1)

Lu Mengkai Zhang Hongwu Zheng Yonggang2)
(International Research Center for Computational Mechanics，State Key Laboratory of Structural Analysis for Industrial Equipment，Department of Engineering Mechanics，Dalian University of Technology，Dalian 116024，China)

Strain localization is a common factor that may lead to the failure of solid structure and its numerical analysis becomes an important aspect for the structural safety evaluation.Due to the heterogeneity and multiscale nature,however,traditional numerical methods need to resolve the structure at the fin scale to obtain reasonable results,which increases drastically the computational scale and cost.To solve this problem,an embedded strong discontinuity model based multiscale finit element method is proposed here.In this method,both coarse and fin scale elements are used to represent the structure.The embedded strong discontinuity model is firs introduced into the fin element to describe the discontinuity and the corresponding additional displacement jump degree of freedom on the elemental level can be eliminated with the condensation technique,which keeps the dimensions of the sti ff ness matrix unchanged.Then,an enhanced multi-node coarse element technique is proposed,which can adaptively insert coarse nodes according to the intersection between the discontinuity line and coarse element boundary and thus guarantees the proper transformation ofinformation between the fin and coarse elements.The problem can then be e ff ectively solved on the coarse scale level.Moreover,a solution decomposition technique,in which the fin scale solution is decomposed into the downscaling and local perturbation solutions,is adopted to eliminate the unbalance forces within the unit cell in the elasto-plastic analysis.Finally,two representative examples are presented to demonstrate the accuracy and e ff ectiveness of the proposed method through the comparisons with the results of the embedded finit element method.

multiscale finit element method,embedded strong discontinuity model,strain localization,enhanced coarse element,elasto-plastic analysis

O242.21,O344.3

：A

10.6052/0459-1879-16-397

2016–12–27 收稿，2017–03–12 錄用,2017–03–13 網(wǎng)絡(luò)版發(fā)表.

1)國家自然科學(xué)基金(11232003,11672062)和中央高校基本科研業(yè)務(wù)費(fèi)(DUT14YQ217)資助項(xiàng)目.

2)鄭勇剛，教授，主要研究方向：多尺度與多場耦合計(jì)算力學(xué).E-mail：zhengyg@dlut.edu.cn

盧夢凱，張洪武，鄭勇剛.應(yīng)變局部化分析的嵌入強(qiáng)間斷多尺度有限元法.力學(xué)學(xué)報(bào),2017,49(3)：649-658

Lu Mengkai,Zhang Hongwu,Zheng Yonggang.Embedded strong discontinuity model based multiscale finit element method for strain localization analysis.Chinese Journal of Theoretical and Applied Mechanics,2017,49(3)：649-658