基于摄动方法的多尺度损伤表示理论

豆 丁 推 荐

精 品 文 档

↓

基于摄动方法的多尺度损伤表示理论

李杰

①②*

, 任晓丹

①

① 同济大学土木工程学院, 上海200092; ② 土木工程防灾国家重点实验室, 上海200092 * E-mail: lijie@tongji.edu.cn

收稿日期: 2009-11-04; 接受日期: 2010-02-03 国家自然科学基金资助项目(批准号: 90715033)

摘要 结构和材料的损伤破坏包含多个尺度, 单独一个尺度的分析很难正确地反映结构和材料的非线性行为. 我们从摄动均匀化理论出发, 基于不可逆热力学理论, 建立了联系细观尺度与宏观尺度的多尺度能量积分, 再结合经典连续损伤理论, 建立了基于细观微结构计算宏观连续损伤变量的一般方法体系, 即为多尺度损伤表示理论. 该理论将多尺度分析方法与传统连续损伤力学紧密结合在一起, 在此基础上建立的数值算法既能够从细观和宏观两个尺度上反映整体结构的损伤和破坏过程, 同时数值模拟的计算量也能够控制在合理的范围内. 最后的数值算例表明了理论的正确性和有效性. 关键词 多尺度, 损伤表示, 摄动方法, 本构关系, 数值模拟 PACS: 21.10.Jx, 46.90.+s, 46.15.-x

材料的损伤和破坏问题是固体力学研究的核心. 从最早的弹脆性-破坏理论, 到后来的塑性理论和断裂理论, 人们的认识逐步深入, 形成了从极限状态分析到全过程分析, 从一维理论到三维理论, 从现象学到物理学的发展脉络. 科学思想的汇集, 催生了损伤力学. 20世纪80年代以来, 损伤力学引起了各国学者的广泛兴趣, 得以迅速发展, 逐渐形成了比较完善的连续损伤理论体系[1~9], 损伤力学也已经成为结构非线性分析的标准工具, 并在实际工程中得到越来越广泛的运用. 当然, 连续损伤力学在其发展的过程中也表现出了不可避免的局限性, 其中最典型的就是在连续损伤的框架内无法得到损伤演化方程. 现有的理论, 大多直接采用试验拟合曲线或根据经验假定来描述损伤的演化, 这一方面破坏了理论的严密性和完整性, 另一方面也局限了理论的适用范围.

近年来, 随着多尺度理论和数值模拟技术的蓬勃发展, 为这一问题的解决带来了契机.

对于混凝土、岩石等脆性材料组成的结构, 其损伤和破坏的过程至少包含两个尺度. 在细观尺度上, 材料内部包含的大量微裂缝决定了材料的非线性行为. 细观尺度的模型需要精细考虑微裂缝及其扩展对材料性能的影响. 在宏观尺度上, 我们更关心结构整体的宏观响应, 材料就视为连续体, 微裂缝和孔洞的作用可以平均化为对材料宏观力学性质的影响. 两尺度分析的思想将结构非线性分析问题分解为三个部分, 即细观材料分析、宏观结构模拟和两尺度的耦合分析. 对于细观材料分析, 早期的研究大多基于解析方法, 而近年来, 一系列有效的数值方法相继提出并得到了长足的发展, 也为细观材料分析提供了广阔的发展前景. 对于宏观结构分析, 经典的非线性

引用格式: 李杰, 任晓丹. 基于摄动方法的多尺度损伤表示理论. 中国科学: 物理学 力学 天文学, 2010, 40: 344 ~ 352

中国科学: 物理学 力学 天文学 2010年 第40卷 第3期

有限元理论已经发展得相当成熟, 在非线性有限元列式的基础上结合显式、隐式及弧长等方法, 一般结构的连续非线性分析都可以得到较好的分析结果.

对于结构多尺度分析问题, 经典研究包含两类方法, 即串行多尺度方法和并行多尺度方法. 对于串行方法, 一般首先建立细观代表性单元的数值模型, 然后在细观数值结果的基础上采用多尺度理论估计宏观结构参数, 最后将得到的参数用于宏观结构模拟. 最典型的串联方法就是所谓的“均匀化”方法[10~14], 微结构, 考虑最小周期性单元, 定义为单胞(Unit Cell). 在单胞上定义微观坐标系y=(y1,y2,y3), 单胞所在区域记为?y, 而单胞内部的所有空洞和裂缝定义为内边界CY, 而内边界上的面力表示为p. 建模过程暂不考虑体力的影响.

宏观坐标系x与微观坐标系y之间的换算可以借助尺度参数λ表示如下:

y=

x

λ, (1)

即在周期性或者统计均匀性微结构假定的基础上, 采用摄动理论计算均匀化材料参数的方法. 这类方法一般只用于线性结构多尺度分析, 算法效率较高; 对于非线性体系, 由于在加载过程中材料的性质发生了不可逆变化, 这类方法很难直接应用. 对于并行算法, 一般需要同时建立宏观结构数值模型和细观数值模型, 以多尺度理论作为纽带在宏观和细观模型之间实时传递数据. 这类算法可用于结构线性或非线性多尺度分析, 宏观和细观的耦合效应被实时考虑, 算法精度较高, 当然所耗费计算量亦很大. 常用的并联多尺度方法有MAAD方法[15]、准连续体方法[16]和大分子动力学方法[17].

在前述研究的基础上, 本文试图基于能量的表达式建立一种能够有效考虑材料非线性的串行多尺度模型. 具体思路是: 首先从材料的细观模拟出发, 基于摄动多尺度理论和连续损伤理论框架, 通过多尺度分析建立损伤演化的多尺度表示理论, 即通过数值模型, 计算得到材料的宏观损伤演化规律, 然后结合宏观连续损伤理论, 对结构进行非线性分析和模拟. 较之传统的串行与并行多尺度研究, 这一研究一方面可望为连续损伤理论中的损伤演化规律建立严格的理论基础, 另一方面又大大降低了材料多尺度非线性模拟的计算量, 可以方便地应用于实际结构的模拟.

1 摄动均匀化理论

考虑具有周期性微结构材料的二维弹性结构分

析问题, 如图1所示. 首先在整体结构上定义宏观坐标系x=(x1,x2,x3), 结构所在区域记为?, 边界记为Γ. 位移边界条件ui=ui作用在边界Γu上, 而面力t作用在边界Γt上. 另一方面, 由于材料具有周期性

对于定义在宏观结构域内并包含细观结构的函数, 可以用两个尺度的坐标表示如下:

Φλ(x)=Φ(x,y), (2)

上标λ 表示宏观函数中的细观结构尺度, 即两尺度函数可由尺度参数变换到统一尺度. 后面分析中将以尺度参数为小参数对状态量做摄动展开并截断, 为了保证摄动方法的精度, 尺度参数λ应该取足够小; 换言之, 单胞较之宏观结构要足够小, 即所谓“宏观足够小”. 那么结合式(1), 函数对空间坐标的偏导数可以表示为

?Φλ(x)?Φ(x,y)?Φ?x==+1?Φ, (3) i?xi?xiλ?yi

对于图1中的宏观材料结构, 建立考虑细观结构的平

衡方程如下: ?σλ

ij

?x=0, (4)

j

其中σλij

包含微空洞的复合材料结构中任意点的应力. 平衡方程对应的边界条件为 uλi=ui, 在边界Γu上, (5)

σλijnj=ti, 在边界Γt, (6)

弹性应力应变关系表示为

σλCλij=ijklεkl, (7)

其中应变定义为

ελ1??uλk

?uλl?

kl=2?+

(8) ??xl?x?,k?

而uλi是复合材料结构任意一点的位移.

上述方程中的场函数uλ, ελ和σλ中均完整的考虑了微观结构的信息, 针对上述方程的数值求解必须对微观结构进行精细的建模, 如此巨大的计算

345

李杰等: 基于摄动方法的多尺度损伤表示理论

图1 宏观结构和微观结构

量是一般数值模拟不能够承受的. 同时空间尺度的差异也极容易导致微分方程数值求解的奇异和发散. 为了减小数值模拟的计算量, 同时为了避免数值奇异, 就需要从两个尺度上分别求解方程.

考虑到位移场uλ包含的尺度参数λ为小参数, 那么可以对位移场做下述摄动展开:

[n]

述应力的摄动展开:

λσij=

∞

[]

1

λnσijn, (13) ∑n=?

[n]

为 其中应力各阶摄动函数σij

[n][n]

σij=Cijklεkl. (14)

将应力的摄动展开代入平衡方程(4)并整理,

u(x)=∑λu(x,y), (9) 可得

λ∞

n

[n]

n=0

其中u(x,y)为u的第n阶摄动函数. 考虑应变的定义式(8)以及函数的空间偏导数(3), 可得下述应变

摄动展开:

λ

[?1][n][n+1]∞??σij??σij1?σijn

++λ??=0, (15) ∑????xyλ2?yjn=?1jj??

由于λ的任意性, (15)式恒成立的条件是λ的各阶系

数等于0, 因此可得下述各阶平衡方程: λ?ulλ?∞n[n]1??ukλ[?1] εkl=?+?=∑λεkl, (10) ??σij

2??xl?xk?n=?1

=0,??y?jλ[n]

的各阶摄动函数εkl(x,y)为 其中应变εkl (16) ?[n]

[n+1]?σij??σij[?1]y

εkl=εkl(u[0]),+=0, n≥?1.??x (11) y?jj?xy[n]εkl=εkl(u[n])+εkl(u[n+1] ),n≥0,

下面我们求解上述得到的摄动方程组, 求解过

其中应变算符定义为

程中只考虑位移场uλ的0阶和1阶摄动项, 对高阶摄

?x1??uk?ul?x

+动项实施截断. ?ε(u)=εkl(u)=??,2xx??k???l

(12) 考虑(9)~(16)式中的第1阶项, 可得下述平衡 ?

1??uk?ul??yy方程: εεuu()()==+??,kl?2??yl?yk?[?1]y??σij?Cijklεkl(u[0])

==0, (17) 将应变摄动展开(10)代入应力应变关系式(7), 可得下?yj?yj

346

中国科学: 物理学 力学 天文学 2010年 第40卷 第3期

[0]xyx

? =Cijklεkl(v[0])+Cijklεkl

?χpq(y)??εpq(v)

yx[0] u[0](x,y)=v[0](x), (18) ?=Cijpq+Cijklεkl(24)?χpq(y)??εpq(v).

可知位移场uλ的零阶摄动函数u[0]表示宏观结构的定义平均化算符 平均化位移场. 1

=d?, (25) 考虑(9)~(16)式中的第2阶项并结合(18)式, 可得Vy∫?y

下述平衡方程:

根据Bakhvalov和Panasenko[10]的讨论, 上述方程的

解与细观坐标y无关, 只是宏观坐标x的函数, 表 示为

=

1

λyxy

Cijklεkl(v[0])+Cijklεkl(v[0])+Cijklεkl(u[1])

{}

{}?σ[?1]σ[0]ij?ij? ?x+=?Cy

ijklεklj?yj?xj

?(v[0])?? +? ?y?j?Cxijklεkl(v[0])+Cyijklεkl(u[1])??=0.(19)考虑到方程(19)是u[1]的线性方程, 它的解可以表示

为下述分离变量的形式[10]:

u[1]

(x,y)=v[1]

(x)+χ(y):εx

??v[0]

(x)??, (20)

其中χ是微观坐标y的3阶张量函数.

在推导上述两个基本解的过程中需用到周期性

条件, 即要求单胞在宏观上具有周期对称性, 当然, 后来的研究中又将这种周期对称性弱化为统计均匀性[14]. 在实际应用过程中, 上述周期性条件要求单胞单元要足够大, 包含足够的微观结构信息, 代表微观结构的行为, 即所谓“微观足够大”.

将(20)式代入平衡方程(19), 可得下述关于χ的方程

?

?y{C+Cy

ijpq

ijklεkl

?(y)?j

?χpq?}=0, (21)

再结合单胞的内边界条件

σ[0]

ij

nj=pi, (22)

有

{C+Cy

ijpqijklεkl??χpq(y)??}nj=pi, (23)

其中n表示单胞内边界的外法向量. 结合式(21)与(23)

即可求解出3阶张量函数χ.

考虑应力摄动展开, 并截断了高阶摄动项, 可得

应力近似表示如下:

σλ=1[?1][0][1]

ijλσij+σij+λσij+…≈1

[?1][0]λσij+σij其中Vy表示单胞的体积. (24)式两边作用平均化算

符, 有 σλij=σij

=Cyijpq+Cijklεkl??χx[0]pq(y)??εpq(v)

={Cyijpq+Cijklεkl

??χx[0]pq(y)??}εpq(v).

(26)

(26)式给出了平均化应力应变关系, 显然其系数张量

Cy

ijpq=Cijpq+Cijklεkl

??χpq(y)??

={Iy

klpq+εkl

??χpq(y)??}Cijkl,

(27)

就是基于微观结构得到的宏观平均化弹性张量. 其中Iklpq为4阶单位张量.

至此我们得到了微缺陷作用下材料的宏观均匀化弹性张量, 从理论上建立了材料多尺度分析的基本框架. 但是, 3阶张量χ的求解和积分都非常的复

杂, 需要耗费大量的计算量. 为了建立更为实用的多尺度材料模型, 本文试图以连续损伤理论作为支撑, 建立的多尺度损伤表示理论. 2 连续损伤理论框架

连续损伤力学的发展也经历由直观到抽象的过程. 最初损伤大都直观的定义为损伤材料与未损材

料的刚度比[1]或者面积比[2]. 后来为了同时考虑塑性和损伤等多种非线性因素的影响, 研究者将不同非线性机制的作用统一定义为能量, 然后在不可逆热力学的基础上建立了完整的理论体系[5,7,8].

根据应变等效假定, 有效应力可以定义为作用

在材料未损伤面积上的应力(图2)

σeij≡Cijklεkl, (28) 其中Cijkl为材料初始弹性刚度. 为了考虑塑性应变的影响, 定义如下应变分解: 347