第 1章各向同性多孔弹性本构模型 多孔材料是自然界中常见的材料类型,例如沙子、砂岩、页岩、金属泡沫、骨头等,其结构都是由固体材料和固体材料之间的孔隙构成的。本书将研究如下一种多孔材料:它由固体材料与孔隙两部分构成,其中固体材料有一定的支撑作用,固体材料之外的部分称为孔隙,孔隙中包含连通的孔隙与不连通的孔隙两类。连通的孔隙中充满某一种特定的孔隙流体,即假设该多孔材料中的流体是饱和的,如图 1.1所示,但对不连通的孔隙部分没有类似要求。 图 1.1多孔充液弹性本构模型微观示意图 可见,多孔材料是一种非均匀又含有孔隙流体的材料,它的力学响应非常复杂,而固体变形往往也和流体扩散现象耦合在一起。为了避免这种微观响应的复杂性,Biot与 Willis [3–5,14]在 1940—1950年代介绍了一种多孔充液 ①弹性本构模型,将其固体部分和被流体浸润的不连通孔隙联合起来看作一种弹性连续介质,称为“固体骨架”(solid skeleton), ①本书中“充液”实际上指充有流体,其中孔隙流体可以是液体或气体。为遵循汉语习惯而统 称为“充液”。 或简称“固体”。进一步地,将固体骨架和充满流体的连通孔隙看作多孔弹性介质,称为“ Biot介质”,并假设该材料具有弹性响应。以后如不特殊声明,“孔隙”一词专指连通的孔隙。 为便于读者理解,本章将从最简单的各向同性多孔弹性本构模型开始介绍。对于更一般的各向异性多孔弹性本构模型的介绍见第 3章。为方便不熟悉张量记号的读者阅读,本书分别使用张量记号与分量记号书写,表达同一个公式。本书只考虑准静态问题。 1.1各向同性多孔弹性本构模型的建立 在广义胡克定律中,应变只与局部应力成线性关系。 Biot为了处理固体变形与孔隙流体渗流耦合在一起的多孔材料,在广义胡克定律的基础上构造了多孔弹性本构模型。这个模型在经典的广义胡克定律的基础上,引入了两个新的流体场变量:在应力一侧加入孔隙流体压力 (或简称 “孔隙压力” ) p、在应变一侧加入孔隙流体体积分数变化量 γ,这是一对功共轭的场变量。基于这种框架,多孔介质中的微观细节就可以忽略,并被看作一种包含了孔隙流体的弹性介质。如非特殊说明,本章所述的应力型和应变型变量 {ξij ,ρij , p, γ}等均表示增量,即本构方程为增量型本构。变量 ξij ,ρij, p, γ都从初始状态算起,起始值均为零。 1.1.1本构方程中的应变-应力关系 ζij(σij,p) 众所周知,在各向同性的弹性本构关系 (广义胡克定律 )中,由于线弹性假设,应变与应力之间有如下的本构关系 ρij (ξij): 1 η ζ = σ . (σ : Σ)Σ, 2G 2G(1 + η) (1.1) 1 η ρij = ξij . ξkkαij 2G 2G(1 + η) 式中, G为胡克弹性材料的剪切模量, η为泊松比, αij称作“克罗内克 α” (Kronecker delta),定义为 1,i = j αij = (1.2)0,i =.j . .. .. 式中,αij为二阶张量 Σ的分量。式 (1.1)右端第二项中 σ : Σ = ξijαij = ξkk,符号“ :”表示双点积,即张量 σ的分量 ξij与张量 Σ的分量 αij两 33 两乘积 (共 9项)之和 ξijαij,称为“指标 i的缩并与指标 j的缩 ∑∑ i=1 j=1 并”。在对两张量进行双点积时,规定在双点积符号“ :”临近的两对指标 ij中,位于前面的两个指标 i进行缩并,同时位于后面的两个指标 j进行缩并。凡是重复出现两次的下标称作“哑标”。这里用到了爱因斯坦求和约定,这是一种针对求和的简写约定:对于一个多项式内的某一项,若该单项式中有两个相同的下标,则该项表示对该下标所有可取的范围求和;若有多组相同的两个下标,则表示对所有下标可取的范围依次求和。如果不是需要被缩并的哑标,一般情况下应避免出现相同指标的情况。例如: ∑ 3 ξkk = ξiii=1 在本章中约定:所有使用拉丁字母的下标 (如 i,k等)取值范围为 {1, 2, 3},而使用希腊字母的下标 (如应,Σ等)取值范围为 {1, 2}。本构关系式 (1.1)常见于弹性力学或材料力学课本,见文献 [15]。 应变张量 ζ可以分解为球形张量 (简称“球量” ) ζ球与偏斜张量 (简称“偏量”) ζ偏两部分: ζ = ζ球 + ζ偏 =1 ρkkΣ + ζ偏 , = ρ球 3 偏 (1.3) ρij ij + ρij 式中, ρ球 =1 ρkkαij ,ρ偏 = ρij . 1 ρkkαij (1.4) ij ij 33 同样应力张量 σ也可以分解为球量 σ球与偏量 σ偏两部分: σ = σ球 + σ偏 =1 ξkkΣ + σ偏 , ξij = ξij 球 +3 ξij 偏 (1.5) 式中, ξij 球 = 13 ξkkαij,ξij 偏 = ξij . 13 ξkkαij (1.6) 则各向同性弹性本构关系式 (1.1)可表示为 ζ球 =1 σ球 , ζ偏 =1 σ偏 (1.7)3K 2G 式中,第一项表达了弹性胡克定律的球量部分,即体积变形部分;第二项表达了偏斜部分,即剪切变形部分。注意偏量部分不仅含有代表剪切的非对角项 (i .11,ρ2233},三项之和 ρkk等于 = j),而且也含有对角项 {ρ偏偏 ,ρ偏偏 零,表示对应的体积变形为零。K为胡克弹性材料的体积模量: 2G(1 + η) K = (1.8) 3(1 . 2η) 而如前文所述,由于多孔弹性材料在应力与应变侧分别加入了 p与 γ,必须分别定义 ρij (ξkl,p)与 γ(ξkl,p)。这里出现在圆括号外的下标 i, j可任选 1, 2, 3之一;但出现在圆括号里作为自变量的下标则表示全部分量。有时为了简便,括号内外一律记作 i, j或 k, l。 对于 ρij (ξkl,p),由线弹性假设,可以要求它们满足 ① 1 η 1 ζ = σ . (σ : Σ)Σ + pΣ, 2G 2G(1 + η)3H (1.9) 1 η 1 ρij = ξij . ξkkαij + pαij 2G 2G(1 + η)3H 式 (1.9)中的前两项来自于式 (1.1),而新增加的第三项来自于如下假设:孔隙流体不抗剪,因而孔隙流体的压力变化只对应变球量产生影响。由式 (1.9)可知,该球量影响可写作 ρll =1 ξkk +1 p (1.10) 3KH 式中, K是多孔充液介质 (即 Biot介质 )当孔隙流体压力保持常数不变 (增量 p =0)时的体积模量,K与同一条件 (p =0)下的剪切模量 G、泊 ①这里隐含一个假设:对于 Biot介质,当应力或孔隙流体压力保持常数不变时,其体积模量 K、剪切模量 G和泊松比 ζ等弹性材料常数与保持不变的应力和孔隙流体常数值无关,称为“胡克弹性材料常数”。 松比 η保持和胡克弹性材料中一样的关系 ①: 1 τξkk 2G(1 + η) K == (1.11) 3 τρll3(1 . 2η) p=0 而 H是一个新的材料常数,它表示当介质整体外载球应力不变 (ξkk =0)时,体积应变增长随孔隙压力增长的比例关系,即 1 τρll (1.12) Hτp νkk=0 该常数写作 H是为了与早期文献 [3,7,10]一致。 式 (1.9)可被进一步整理得 ρij = 1 2G Eij . η 2G(1 + η) Ekkαij (1.13) 式中, Eij = ξij + αpαij (1.14) 即 Biot有效应力,α为 Biot有效应力系数: α = K (1.15) H 在大多数有关多孔弹性本构模型的书籍文献中, α比 H更为常用。读者可将式 (1.14)和式 (1.15)代入式 (1.13)自行验证其正确性。可见,基于有效应力的多孔弹性本构方程 (式 (1.13))形式上与广义胡克定律 (式 (1.1))一样。 基于 α的定义,式 (1.9)也可写作 1 ηα(1 . 2η) ρij = ξij . ξkkαij + pαij (1.16) 2G 2G(1 + η)2G(1 + η) 由式 (1.13)也容易推导出 ξij (ρij,p): 2Gη ξij =2Gρij + ρkkαij . αpαij (1.17) 1 . 2η 本书符号规定:应力 ξij拉伸为正,应变 ρij伸长为正,压力 p受压为正。 ①本书采用在偏导运算符的竖线右下角书写条件的方式表示在对应条件的前提下求偏导数。例如式 (1.11)表示当 p =0时,求 νkk对于 ρll的偏导数。 1.1.2折算体积分数〈的引入 在 1.1节的开头,已经介绍了 Biot在广义胡克定律的基础上添加了孔隙流体压力 p和一个孔隙流体体积分数变化量 γ。但这不免给读者留下疑惑:究竟什么是 γ,为什么要引入它,又为什么说它们是功共轭的? 另一方面,由图 1.1可知,多孔材料被分割为三个部分:固体材料、充满流体的连通孔隙和不连通的孔隙。为何要将孔隙区分为连通与不连通的两部分? 实际上,这是由孔隙流体的客观性质决定的:考察孔隙流体的渗流过程时,需要对其列出流体质量守恒方程 (见 1.4.1节基本方程中的式 (1.92))。显然,只有连通孔隙中的孔隙流体才应当建立质量守恒方程。假设孔隙流体是可压缩的。为了方便起见,把由质量按初始密度 ν0计算的体积称为“折算体积”。定义 m为流入变形前每单位初始体积 Biot介质中的孔隙流体质量 ①,并定义 γ = m/ν0为 Biot介质每单位初始体积接受来自周围或外界的流体折算体积 ②。故有 m = m1 . m0 = ν0(γ1 . γ0)= ν0γ (1.18) γ = γ1 . γ0 (1.19) 式 (1.18)中带有下标 0的量表示初始值,带有下标 1的量表示变形后的值,而不带下标的量依前述定义表示增量。 m0表示 Biot介质在变形前每单位体积所含的孔隙流体质量, γ0 = m0/ν0表示初始状态的孔隙体积分数或体积比, γ1 = m1/ν0。应当注意的是,在这种定义规则下,不带下标的量 (即增量 )在变形前的初值为零,例如 m0表示初始孔隙流体质量,而增量 m在变形前为零。这也是为了与其他场变量 (如应力、应变等)的增量标记保持一致。在本节中,至此共有 {m, γ, ν}及将在式 (1.22)和式 (1.23)中定义的 v四个物理量采用了这类定义方式。在初始状态, m0,γ0,ν0,v0都不是零,但 m, γ, ν, v表示增量都等于零。 ①读者应注意区分本节定义的孔隙流体质量 m、第 3章中定义的各向异性材料常数二阶张量 m,以及横观各向同性中的材料常数标量 m。 ②在文献 [10]中 Detournay与 Cheng称 β为“流体含量” (fluid content)。由于 β = m/θ0, θ0为常数,因此 β实际上代表质量,但具有体积的量纲。 如前所述,在微元整体 .中,除了充满流体的连通孔隙部分 .p之外,固体材料部分与不连通的孔隙应当被混在一起看待。包含了所有固体材料和不连通的孔隙部分的并集,用 .s表示,即在本章开头所述的固体骨架,并有 .=.s ∞ .p,V = Vs + Vp (1.20) 式中,V = |.| ,Vp =.p ,Vs = |.s|表示对应部分的体积。关于固体骨架变形更详细的讨论见 3.1节中的相关分析。基于这样的分割,可以定义材料的 (初始 )孔隙比 φ0(也称“ (初始 )孔隙体积分数”): Vp0 Vs0 φ0 = , 1 . φ0 = (1.21) V0 V0 式中, V0,Vp0,Vs0表示 V, Vp,Vs的初值。孔隙比 φ0是一个材料常数,它的值就定义为介质在变形前的初始孔隙比。随着 Biot介质的变形, Vp/V可能会发生变化,但是在线性弹性理论中, .Vp = Vp . Vp0以及 Vp/V . Vp0/V0都是一阶小量[12]。 本书分析的是小变形、小位移的线弹性理论,在大部分场合中,将孔隙比替换为初始孔隙比对相应公式都只会带来高阶小量的改变。为略去高阶小量,凡是在这种场合中,本书都用初始孔隙比 φ0作为 Vp/V的替代。 另一方面,也需要考察连通的孔隙部分的体积分数增量,将其定义为 ① .Vp v = (1.22) V0 并按照式 (1.18)处所描述的规则,记 Vp0 Vp v0 = ,v1 = (1.23) V0 V0 ①在文献 [7]中 Rice与 Cleary定义了如下的表观体积分数 (apparent fluid volume fraction): Vp vRC = V0这实际上与式 (1.23)中使用的 v1为同一个物理量。本书为了统一记号,使用不带下标的物理量表示增量,分别使用下标 0和 1表示变形前与变形后的值,故与上述定义有所区别。也请读者注意,书中定义的 v与本书作者发表的文章 [16]中的 v定义不同。更多细节见第 78页的脚注①。 因而依然满足 v = v1 . v0,也同样按照约定记密度的变化量 ν = ν1 . ν0, v0与 v1各表示 Biot介质每单位初始体积变形前的初始孔隙体积与变形后的最终孔隙体积。若孔隙增大, v为正,否则为负。但应注意到,增量 v由两部分构成: ①由孔隙压力变化导致的孔隙内流体本身体积膨胀或收缩; ②由其他微元的流体流入或流出导致的孔隙体积分数变化。从式 (1.18)中也可以得到印证。因而有 m = m1 . m0 = ν1v1 . ν0v0 . ν0(v1 . v0)+ v0(ν1 . ν0) = ν0v + v0ν (1.24) 式中的“ .”实际上略去了一个二阶小量 νv =(ν1 . ν0)(v1 . v0)。尽管对比定义式 (1.21)1和式 (1.23)1即可知 v0 = φ0,但本书为了方便读者区分材料常数孔隙比 φ0与孔隙体积分数增量 v,还是使用了不同的符号表述这两个物理量。 若进一步假设孔隙流体在一定范围内的密度是随压力线性变化的,即 νν1 . ν0 p = == Cf p (1.25)ν0 ν0 Kf 式中,Kf为孔隙流体的体积模量,Cf = 1/Kf为体积模量的倒数。则将式 (1.18)和式 (1.25)代入式 (1.24),可得 m = ν0γ . ν0(v + v0Cf p) (1.26a) γ = v + φ0Cf p (1.26b) 式 (1.26a)的物理意义如下:类似于式 (1.23)下所说的关于孔隙体积变化分两部分组成,从 Biot介质微元周围或外部流入微元每单位初始体积中的流体质量 m也由两部分组成。第一部分 (右端第一项 )为填满增大的孔隙体积 v = v1 . v0(由 v0增大到 v1)所需,而第二部分 (右端第二项)为提高流体密度 ν = ν1 . ν0(由 ν0增大到 ν1)所需。式 (1.26b)的物理意义同式 (1.26a),只不过用“折算体积”来代替质量。 作为 Biot介质的本构关系,若以应力型变量 ξkl与 p为自变量,前文中只有应变型变量 ρij(ξkl,p)的表达式 (1.9)和式 (1.16),还缺少 γ(ξkl,p)。其实 (ρij,γ)和 (ξkl,p)是两对互为函数的变量: ρij = ρij(ξkl,p) (1.27a) γ = γ(ξkl,p) (1.27b) { ξkl = ξkl(ρij ,γ) (1.28a) p = p(ρij ,γ) (1.28b) {(1.27b)下面讨论如何建立式。 与式 (1.9)之前关于圆括号内外的下标记法类似,这里出现于左端的下标 i, j, k, l表示任选1, 2, 3之一;但出现在右端的圆括号内作为自变量的下标则表示全部分量。式 (1.27a)就代表式 (1.9)或式(1.16)。如果同时还有式 (1.27b),就可以由式 (1.27a)和式 (1.27b)解出式 (1.28a)和式 (1.28b)。 作为弹性介质的特性,存在应变能函数 W (ρij ,γ)与余能函数 W*(ξij ,p)使得 τW (ρij ,γ) ξij = (1.29a)τρij τW (ρij ,γ) p = (1.29b) τγ ρij = τW*(ξij,p) (1.30a) . . . .. . . . . τξij τW*(ξij,p) γ = (1.30b) τp 那么利用微积分中函数的二阶偏导数与求导先后次序无关的定理,可由式 (1.29)和式 (1.30)得到 ① τξij(ρij ,γ) τp(ρij ,γ) = (1.31) τγ τρij τρij(ξij ,p) τγ(ξij ,p) = (1.32) τp τξij ①对应变分量 ρij或应力分量 νij求导,当 i .j时,在连续介质力学中有两种方法。一种 = 是张量分析中的方法,将 ρij与 ρji(或 νij与 νji)看作独立的变量。另一种是工程记法:对剪应力,考虑 νij = νji,两者不加区分但改记为 πij,πji与 πij是同一个变量;而对剪应变,也考虑 ρij = ρji,但改用 7ij = 7ji = ρij + ρji,7ji与 7ij是同一个变量。有些细节详见文献 [17],本书采用张量分析中的方法。 对各向同性材料,将式 (1.9)代入式 (1.32),可得 τγ τρij 1 == αij (1.33) τξij τp 3H 故对于线弹性情况,必有 (利用 ξij αij = ξkk) γ =1 ξkk + CCHp (1.34) 3H 式中第二项的系数 CCH是一个新的常数: CCH τγ (1.35) τp νkk=0 常数 CCH在文献 [10]中也写作 1/R ′。①在本书中写作 CCH是为了与第 3章中的记号统一,其中下标 CH表示该常数是 Cheng[18]定义的。由式 (1.35)可知,常数 CCH表示材料在不受力的状态下,每提升单位孔隙压力后单位初始体积材料能够吸收孔隙流体 (按折算体积计算)的能力。 两对方程式 (1.29)和式 (1.30)、式 (1.31)和式 (1.32)之间,属于自变量与函数的转换关系,这种关系常出现在连续介质力学中,即勒让德变换 (Legendre transformation,见文献 [17])。 现在来计算 Biot介质的变形功。把应变型变量 (ρij ,γ)从初始状态 (ρij =0,γ = 0)到最终状态 (ρij ,γ = γ1 . γ0)的加载过程分割为许多小段 (dρij , dγ)。对于 Biot介质每单位初始体积 (V0 = 1),每一小段的变形元功为 dW = σ : dζ + p dγ = ξij dρij + p dγ (1.36) 因为弹性材料的元功转化为应变能的增量 dW,所以式 (1.29)恰好是式 (1.36)构成全微分的条件,它的物理意义是:在应变空间 (ρij,γ)中从 ①在两篇关于多孔弹性本构的经典文献 [10]和文献 [7]中都定义了材料常数 R,但应当注意不要把这两个常数搞混,它们之间满足 11 φ0 = . (1.37) RRC RDC Kf 式中,RRC即文献 [7]中的 R;RDC即文献 [10]中的 R ′,也是本书中的 1 ;φ0为材料的初始CCH 孔隙比;而 Kf为孔隙流体的体积模量,见式 (1.25)。