Asymptotic Solutions of Plastic Stress and Displacement at V-Notch Tips Under Anti-Plane Shear
-
摘要:
提出了一种确定幂硬化材料反平面V形切口尖端应力和位移渐近解的主导项和高阶项的有效方法。首先通过在弹塑性理论基本方程中引入V形切口尖端应力场和位移场的渐近级数展开，建立以应力和位移为特征函数的非线性和线性常微分方程组。然后采用插值矩阵法求解常微分方程组，可得到多阶应力特征指数和其相对应的特征函数。该方法具有通用性强、精度高等优点，可处理任意开口角度和应变硬化指数的V形切口。典型算例验证了该方法的准确性和有效性。
Abstract:An efficient method was developed to determine the first- and high-order terms of asymptotic solutions of plastic stress and displacement near V-notch tips under anti-plane shear in power-law hardening materials. Through introduction of the asymptotic series expansions of stress and displacement fields around the V-notch tip into the fundamental equations of the elastoplastic theory, the governing ordinary differential equations (ODEs) with the stress and displacement eigen-functions were established. Then the interpolating matrix method was employed to solve the resulting nonlinear and linear ODEs. Consequently, the high-order stress exponents and the associated eigen-solutions were obtained. The presented method, being capable of dealing with the V-notches with arbitrary opening angles and strain hardening indexes under anti-plane shear, has the advantages of great versatility and high accuracy. Typical examples were given to demonstrate the accuracy and effectiveness of this method.
-
Key words:
- elastoplasticity /
- anti-plane /
- V-notch /
- asymptotic solution /
- singularity
-
引 言
准确获得V形切口尖端附近的应力场和位移场[1]对评价材料和（或）结构的失效具有重要意义。Hutchinson[2]、Rice和Rosengren[3]首先研究了幂硬化材料平面裂纹结构在面内荷载作用下的塑性应力奇异性。他们考虑了应力渐近展开式的第一阶应力主导项，并通过应变协调方程将裂纹尖端附近的塑性应力奇异性问题转化为四阶常微分方程的非线性特征值问题，然后结合Runge-Kutta法和试射法求解方程组，获得著名的“HRR”解。之后，Kuang和Xu[4]，Xia和Wang[5]研究了V形切口尖端附近的奇异应力场。在小尺度屈服的假设下，渐近级数展开式的第一阶主导项可用来描述尖端附近的应力分布。倘若尖端区域出现大规模屈服，则应考虑高阶项的贡献。Li和Wang[6]使用应力函数法对平面Ⅰ型裂纹结构尖端附近的应力场进行了二项渐近分析。Sharma和Aravas[7]以应力和位移为基本变量，将裂纹在平面应变条件下的尖端附近塑性解按两项渐近展开。Xia等[8]、Yuan和Yang[9]计算了平面裂纹在Ⅰ型载荷作用下尖端附近应力场的高阶项。Yang等[10]、Chao和Yang[11]研究了裂纹在面内Ⅱ型载荷作用下尖端附近的应力场。
Rice[12]通过hodograph变换研究了单相均质材料V形切口/裂纹结构在反平面剪切荷载下的应力奇异性问题。Amazigo[13]获得了无限大含裂纹结构的完全塑性解。Yang等[14-15]确定了幂硬化材料和线性硬化材料在反平面剪切作用下的高阶渐近弹塑性裂纹尖端应力场。Yuan和Yang[16]结合hodograph变换和渐近分析方法，获得了反平面裂纹在剪切荷载作用下的高阶渐近弹塑性裂纹尖端应力场。Wang和Kuang[17]采用hodograph变换方法，获得了损伤非线性材料在反平面剪切载荷作用下的应力场和应变场的前二阶解。Zappalorto和Lazzarin[18-20]研究了抛物线切口尖端的弹塑性应力场，其为不同非线性规律下的应力和应变分布提供了统一的解析框架。
除解析方法外，Sharma和Aravas[7]获得了平面应变条件下，裂纹尖端附近的前两阶塑性渐近解。然后，Aravas和Blazo[21]研究了应力和位移渐近展开式的第二项对位移场的影响。Loghin等[22]通过有限元特征分析，获得了反平面剪切作用下应力场的前三个应力项，并确定了V形切口/裂纹结构在面内载荷作用下的奇异应力场[23-24]。
切口尖端附近的高阶应力解在弹塑性V形切口的应力场中起着重要作用，但目前鲜见文献报道任意张角情形反平面撕开型V形切口的弹塑性应力奇异性的高阶项解。本文提出了一种研究幂硬化材料V形切口在应对反平面载荷作用下尖端附近应力奇异性的新方法，获得了V形切口/裂纹尖端附近弹塑性应力奇异性的完全渐近解。
1. 反平面V形切口塑性应力奇异性控制方程
考虑张角为2α的反平面V形切口结构，见图1。z轴垂直于rOθ平面，原点位于V形切口尖端。对于反平面问题，应力分量具有以下特征：
σrr=σθθ=σrθ=σzz=0. (1) 在极坐标系下，反平面切口尖端区域的应力场可假设成渐近级数形式：
σij(r,θ)=N∑k=1Akrsk˜σijk(θ),ij=rz,θz, (2) 式中
${s_k}$ 为应力特征指数，且有${s_k} < {s_{k + 1}}$ ；$ {\tilde \sigma _{ijk}}(\theta ) $ 是与${s_k}$ 相对应的应力特征函数；${A_k}$ 称为幅值系数；N是渐近展开式的截断项数。若${s_k} < 0$ ，该应力项在$r \to 0$ 时产生奇异性，${s_k}$ 也称为应力奇异指数(阶)。注意到应力特征指数
${s_k}$ 及其对应的应力特征函数$ {\tilde \sigma _{ijk}}(\theta ) $ 取决于结构的材料特性、切口张角${\text{2}}\alpha $ 以及楔形边的边界条件，而与外荷载无关。故在不考虑体力的情况，反平面问题的平衡方程为∂σrz∂r+1r∂σθz∂θ+σrzr=0. (3) 将式(2)代入式(3)，消去公因子
$ {A_k}{r^{{s_k}-1}} $ ，得到关于${s_k}$ 和$ {\tilde \sigma _{ijk}}(\theta ) $ 的常微分方程如下：(sk+1)˜σrzk(θ)+˜σ′θzk(θ)=0, (4) 式中
$(\cdot)^{\prime }=\partial (\cdot)/\partial \theta$ 。注意反平面问题有2个应力分量特征函数${\tilde \sigma _{r{\textit{z}}k}}(\theta )$ ，${\tilde \sigma _{\theta {\textit{z}}k}}(\theta )$ ，而常微分方程只有一个，因此需添加新方程。考虑幂硬化材料模型，单轴拉伸试验下应力-应变关系满足Ramberg-Osgood定律[25]：
εεy=σσy+γ(σσy)m, (5) 式中
${\sigma _{\rm{y}}}$ ，${\varepsilon _{\rm{y}}}$ 分别是屈服应力和屈服应变；$\gamma $ 是材料系数；$m$ 是材料的应变硬化指数，这里$m = \infty $ 表示理想的弹塑性模型，$m = 1$ 表示线弹性模型。根据弹塑性形变理论，广义的应力-应变关系为εij=1+νEsij+1−2ν3Eσkkδij+32βσm−1esij,ij=rz,θz, (6) 式中
$\beta {\text{ = }}\gamma {\varepsilon _{\rm{y}}}/\sigma _{\rm{y}}^m$ ，$E$ 是弹性模量，$\nu $ 是材料的Poisson比，${\sigma _{\rm{e}}}$ 是等效应力，${s_{ij}}$ 是应力偏张量。如采用von Mises等效应力，有σe=(32sijsij)1/2, (7) sij = σij−13σkkδij, (8) 式中
${\sigma _{ij}}$ 是应力张量，${\delta _{ij}}$ 为符号张量。将式(2)代入式(8)中得
sij(r,θ)=N∑k=1Akrsk˜σijk(θ). (9) 将式(9)代入式(7)可得von Mises等效应力为
σe=[N∑h=1N∑l=1AhAlrsh+sl⌢σhl(θ)]1/2, (10) 式中
⌢σhl(θ)=3˜σrzh(θ)˜σrzl(θ)+3˜σθzh(θ)˜σθzl(θ). (11) 由此，
${\sigma _{\rm{e}}}$ 被写成σe=A1rs1[⌢σ11(θ)]1/2[1+N∑h=2N∑l=2AhAlA21rsh+sl−2s1⌢σhl(θ)⌢σ11(θ) + N∑h=2AhA1rsh−s1⌢σ1h(θ)+⌢σh1(θ)⌢σ11(θ)]1/2. (12) 对式(6)中的
$\sigma _{\rm{e}}^{m - 1}$ 使用二项式展开，有σm−1e=Am−11r(m−1)s1[⌢σ11(θ)]1/2{1+m−12[N∑h=2N∑l=2AhAlA21rsh+sl−2s1⌢σhl(θ)⌢σ11(θ)+N∑h=2AhA1rsh−s1⌢σ1h(θ)+⌢σh1(θ)⌢σ11(θ)] + (m−1)(m−3)9[N∑h=2N∑l=2AhAlA21rsh+sl−2s1⌢σhl(θ)⌢σ11(θ)+N∑h=2AhA1rsh−s1⌢σ1h(θ)+⌢σh1(θ)⌢σ11(θ)]1/2+⋅⋅⋅}. (13) 由
${s_k} < {s_{k + 1}}$ 可知，当$r \to 0$ 时，${r^{{s_h} + {s_l} - 2{s_1}}}$ 和${r^{{s_h}-{s_1}}}$ 是无穷小量。舍去其小于${r^{{s_h} + {s_l} - 2{s_1}}}$ 的高阶无穷小量，式(13)可近似表达为σm−1e≈Am−11rt11⌢σ∗11(θ)+m−12{N∑h=22AhAm−21rth1⌢σ∗1h(θ)+N∑h=2N∑l=2AhAlAm−21rth1[⌢σ∗hl(θ)+⌢σ+hl(θ)]}, (14) 式中
thl=tlh=(m−3)s1+sh+sl, (15) σ∗hl(θ)=[˜σ11(θ)](m−3)/2˜σhl(θ), (16) σ + hl(θ)=(m−3)[˜σ11(θ)](m−5)/2˜σ1h(θ)˜σ1l(θ). (17) 在广义的应力应变关系(6)中，前2项是弹性应变
$\varepsilon _{ij}^{\rm{e}}$ ，即εeij = 1 + νEsij+1−2ν3Eσkkδij=1E(1 + ν)σij； (18) 第3项是塑性应变
$\varepsilon _{ij}^{\rm{p}}$ ，即εpij = 32βσm−1esij. (19) 对于塑性应变部分，将式(9)和式(14)代入式(19)中，
$\varepsilon _{ij}^{\rm{p}}$ 可表示成εpij = 32βσm−1esij=εoij+ε×ij, (20) 式中
εoij = H111rt111σ∗11(θ)˜sij1(θ)+ N∑k=2H11krt11k[σ∗11(θ)˜sijk(θ)+(m−1)σ∗1k(θ)˜sij1(θ)], (21) ε×ij = N∑h=2N∑k=2(m−1)H1hkrt1hkσ∗1h(θ)˜sijk(θ)+N∑h=2N∑l=2N∑k=1m−12Hhlkrthlk[σ∗hl(θ)+σ+hl(θ)]˜sijk(θ), (22) Hhlk = 32β×Am−31AhAlAk, (23) thlk = thl+sk=(m−3)s1+sh+sl+sk. (24) 由式(21)、(22)可以看出，应力特征指数
$ {s_k} $ 与应变$ \varepsilon _{ij}^{\rm{o}} $ 中$ r $ 的幂次$ {t_{11k}} $ 一一对应，可将切口尖端区域位移场假设成uz(r,θ)=N∑k=1H11krλk+1˜uzk(r,θ)+Ti(r,θ), i = r,θ, (25) 式中
{λk=t11k=(m−1)s1+sk,sk=λk+λ1(1m−1), (26) 其中，
$ {\lambda _k} $ 是位移特征指数，$ {\tilde u_{ik}}(\theta ) $ 是位移特征函数，$ {T_i}(r,\theta ) $ 是位移$ {u_i}(r,\theta ) $ 关于r展开式的后续高阶余项（在接下来的推演中被舍去）。接下来，使用反平面问题的几何方程为
{εrz=∂uz∂r,εθz=12r∂uz∂θ. (27) 显然，塑性应变的两部分
$ \varepsilon _{ij}^{\rm{o}} $ ，$ \varepsilon _{ij}^ \times $ 均应满足几何方程(27)。将式(21)、(22)和式(25)代入几何方程(27)，可获得两个新的方程：{σ∗11(θ)˜σrz1(θ)−(λ1+1)˜uz1(θ)=0,σ∗11(θ)˜σθz1(θ)−˜u′z1(θ)=0,k=1, (28) {f1˜σrzk(θ)+f2˜σθzk(θ)−(λ1+1)˜uzk(θ)=0,g1˜σrzk(θ)+g2˜σθzk(θ)−˜u′zk(θ)=0,k⩾2, (29) 式中
{f1=σ∗11(θ)+3(m−1)[˜σ11(θ)](m−3)/2˜σrz1(θ)˜σrz1(θ),f2=3(m−1)[˜σ11(θ)](m−3)/2˜σrz1(θ)˜σθz1(θ),g1=3(m−1)[˜σ11(θ)](m−3)/2˜σrz1(θ)˜σθz1(θ),g2=σ∗11(θ)+3(m−1)[˜σ11(θ)](m−3)/2˜σθz(θ)˜σθz1(θ). (30) 联立常微分方程(4)和应力位移关系(28)、(29)，获得关于位移和应力分量特征函数的常微分方程组：
{(λ1m+1)˜σrz1(θ)+˜σ′θz1(θ)=0,σ∗11(θ)˜σrz1(θ)−(λ1+1)˜uz1(θ)=0，σ∗11(θ)˜σθz1(θ)−˜u′z1(θ)=0,k = 1, (31) {[λk+λ1(1m−1)+1]˜σrzk(θ)+˜σ′θzk(θ)=0,f1˜σrzk(θ)+f2˜σθzk(θ)−(λk+1)˜uzk(θ)=0，g1˜σrzk(θ)+g2˜σθzk(θ)−˜u′zk(θ)=0,k⩾2. (32) 方程(31)是关于第1阶奇异性的非线性常微分方程组特征值问题，而方程(32)是关于高阶项的线性常微分方程组特征值问题，须依赖于第1阶奇异性的解。
2. 反平面V形切口塑性应力奇异性边界条件和定解问题
考虑一般反平面V形切口两楔形边界类型。若边界
$ {\varGamma _{\text{1}}} $ 和$ {\varGamma _{\text{2}}} $ 上的面力$ {\sigma _{\theta z}}(r,\theta ) = 0 $ ，由式(2)得˜σθzk(θ1)=0, ˜σθzk(θ2)=0. 若边界
$ {\varGamma _{\text{1}}} $ 为固支边界，则有$ {u_z}(r,\theta ) = 0 $ ，由式(25)得˜uzk(θ1)=0. (34) 通过求解方程(31)和(32)和边界条件(33)或(34)，获得反平面V形切口尖端区域第1阶和后续高阶项位移特性指数
$ {\lambda _k} $ 、相应应力场(2个分量)和位移场(1个分量)的特征函数解，再将$ {\lambda _k} $ 代入式(26)求得应力特征指数$ {s_k} $ ，这里分别记为$ s_k^{\rm{o}} $ ，$ \tilde \sigma _{ijk}^{\rm{o}}(\theta ) $ ，$ \tilde u_{ik}^{\rm{o}}(\theta ),k = 1,2, \cdots ,N $ 。然而，仅由应力特征指数
$ s_k^{\rm{o}} $ 组成的序列$ s_{\text{1}}^{\rm{o}} < s_2^{\rm{o}} < s_3^{\rm{o}} < s_4^{\rm{o}} < \cdot \cdot \cdot $ 可能并不完整，这是因为在$ \varepsilon _{ij}^ \times $ 中$ r $ 的某些幂次可能与$ {t_{11k}} $ 相等。为便于讨论，将所有塑性应变中可能出现的$ r $ 的幂次罗列如下：εeij: sl, l=1,2,3,⋅⋅⋅, εoij: t11k=(m−1)s1+sk, k=1,2,3,⋅⋅⋅, ε×ij: thlk=(m−3)s1+sh+sl+sk, h=1,2,3,⋅⋅⋅,l=2,3,⋅⋅⋅,k=1,2,3,⋅⋅⋅. 如果弹性应变
$\varepsilon _{ij}^{\rm{e}}$ 中的$ {s_l} $ 与$ {t_{11k}} $ 相等，由式(26)可得{s#k=(1−m)s1+sl,λ#k=sl,1⩽l<k, (35) 其中，
$ s_k^\# $ 和$ \lambda _k^\# $ 分别是由弹性变形控制的应力和位移特征指数，而与之对应的应力和位移特征函数可由如下线性方程求得{[λ#k+λ1(1m−1)+1]˜σrzk(θ)+˜σ′θzk(θ)=0,f1˜σrzk(θ)+f2˜σθzk(θ)−(λ#k+1)˜uzk(θ)=−Δ#rzl(θ),g1˜σrzk(θ)+g2˜σθzk(θ)−˜u′zk(θ)=−Δ#θzl(θ),k⩾2, (36) 式中
{Δ#rzl(θ) = (1 + ν)˜σrzl(θ),Δ#θzl(θ) = (1 + ν)˜σθzl(θ). (37) 类似的情形也可能发生在塑性应变
$ \varepsilon _{ij}^{\rm{o}} $ 和$ \varepsilon _{ij}^ \times $ 中。例如，若$ \varepsilon _{ij}^ \times $ 中最小的幂次$ {t_{122}} $ 和$ {t_{22{\text{1}}}} $ 与某一个$ {t_{1{\text{1}}k}} $ 相等，我们便可获得可能遗漏的特征解：{s×k=2s2−s1,λ×k=2λ2−λ1, or {s×k=3s2−2s1,λ×k=3λ2−2λ1. (38) 通过式(38)，我们可以看出，
$ s_k^ \times $ 取决于前2阶应力指数$ {s_{\text{1}}} $ 和$ {s_{\text{2}}} $ 。与$ s_k^ \times $ 相对应的应力和位移分量特征函数可由式(39)计算：{[λ#k+λ1(1m−1)+1]˜σrzk(θ)+˜σ′θzk(θ)=0,f1˜σrzk(θ)+f2˜σθzk(θ)−(λ#k+1)˜uzk(θ)=−Δ×rzl(θ),g1˜σrzk(θ)+g2˜σθzk(θ)−˜u′zk(θ)=−Δ×θzl(θ),k⩾3, (39) 式中
{Δ×rzl(θ) = m−12{ [σ∗22(θ) + σ+22(θ)]˜srz1(θ) + [σ∗21(θ) + σ+12(θ)]˜srz2(θ)} ,Δ#θzl(θ) = m−12{ [σ∗22(θ) + σ+22(θ)]˜sθz1(θ) + [σ∗21(θ) + σ+12(θ)]˜sθz2(θ)} ,s×k=2s2−s1, (40) {Δ×rzl(θ) = m−12[σ∗22(θ) + σ+22(θ)]˜srz2(θ),Δ#θzl(θ) = m−12[σ∗22(θ) + σ+22(θ)]˜sθz2(θ),s×k=3s2−2s1. (41) 最后，按照从小到大的顺序排列
$ s_k^{\rm{o}} $ ，$ s_k^{\text{\# }} $ 和$ s_k^ \times $ ，便可获得尖端区域完整的弹塑性应力渐近场特征解序列：sok: so1, so2, so3, so4,⋅⋅⋅,k⩾1, s\# k: (1−m)s1+s1, (1−m)+s2, (1−m)+s3,⋅⋅⋅,k⩾2, s×k: 2s2−s1, 3s2−2s1, ⋅⋅⋅,k⩾3. 下面以应力特征指数为例，详细阐述反平面V形切口弹塑性应力渐近解前4阶的确定过程。
1) 第1阶
基于前文的分析，易于发现
$ s_k^{\text{\# }} $ 和$ s_k^ \times $ 不会出现在第1阶特征项中，因此s1 = so1. (42) 2) 第2阶
特征解序列第2阶
$ {s_{\text{2}}} $ 解有2种可能，即s2 = min (43) 3) 第3阶
从第3阶开始，其特征解的确定取决于第1，2阶的结果：
{s_{\text{3}}}{\text{ = }}\min \{ s_{\text{3}}^{\rm{o}},{\text{ }}(2 - m){s_1},{\text{ 2}}{s_2} - {s_1}\} ,\qquad {\text{ }}{s_2} = s_2^{\rm{o}}， (44) {s_{\text{3}}}{\text{ = }}\min \{ s_{\text{2}}^{\rm{o}},{\text{ }}({\text{1}} - m){s_1}{\text{ + }}{s_2}\} ,\qquad {\text{ }}{s_2} = (2 - m){s_1}{\text{.}} (45) 4) 第4阶
与第3阶类似，我们依次确定第4阶和更高阶特征解：
{s_{\text{4}}}{\text{ = }}\min \{ s_{\text{4}}^{\rm{o}},{\text{ }}({\text{2}} - m){s_1}{\text{, 2}}{s_2} - {s_1}\} ,\qquad {\text{ }}{s_2} = s_2^{\rm{o}},{\text{ }}{s_3} = s_3^{\rm{o}}, (46) {s_{\text{4}}}{\text{ = }}\min \{ s_{\text{3}}^{\rm{o}},{\text{ 2}}{s_2} - {s_1}\} ,\qquad {\text{ }}{s_2} = s_2^{\rm{o}},{\text{ }}{s_3} = (2 - m){s_1}, (47) {s_{\text{4}}}{\text{ = }}\min \{ s_{\text{3}}^{\rm{o}},{\text{ }}({\text{2}} - m){s_1},{\text{ 3}}{s_2} - {\text{2}}{s_1}\} ,\qquad {\text{ }}{s_2} = s_2^{\rm{o}},\;{s_3} =2{s_2} - {s_1}, (48) {s_{\text{4}}}{\text{ = }}\min \{ s_{\text{3}}^{\rm{o}},{\text{ }}({\text{1}} - m){s_1}{\text{ + }}{s_2}{\text{, 2}}{s_2} - {s_1}\} ,\qquad {\text{ }}{s_2} = (2 - m){s_1},{\text{ }}{s_3} = s_2^{\rm{o}}, (49) {s_{\text{4}}}{\text{ = }}\min \{ s_{\text{2}}^{\rm{o}},{\text{ }}({\text{1}} - m){s_1}{\text{ + }}{s_2}\} ,\qquad {\text{ }}{s_2} = (2 - m){s_1},{\text{ }}{s_3} = ({\text{1}} - m){s_1}{\text{ + }}{s_2}{\text{.}} (50) 3. 插值矩阵法求解非线性控制微分方程
反平面V形切口的应力奇异性分析转化为解控制方程(31)、(32)，以及相应的边界条件(33)或(34)。我们采用插值矩阵法[26]求解，须先解非线性方程(31)。首先，区间
$ [{\theta _1},{\text{ }}{\theta _2}] $ 被划分为Q个子段。在每个子段内，应力和位移特征函数$ {\tilde \sigma _{ijk}}(\theta ),{\text{ }}{\tilde u_{ik}}(\theta ) $ 其导函数由插值函数表示，本文采用分段2次函数插值。插值矩阵法解得的$ {\lambda _k},{s_k} $ 与对应特征函数$ {\tilde \sigma _{ijk}}(\theta ),{\text{ }}{\tilde u_{ik}}(\theta ) $ 及其导函数具有相同的计算精度。为提高非线性常微分方程(31)的迭代求解效率，首先令
{\hat \sigma _{ijk}}(\theta ){\text{ = }}\sigma _{{\text{11}}}^ * (\theta ){\tilde \sigma _{ijk}}(\theta ), (51) 那么
{\tilde \sigma _{ijk}}(\theta ){\text{ = }}\frac{{{{\hat \sigma }_{ijk}}(\theta )}}{{\sigma _{{\text{11}}}^ * (\theta )}}, (52) {\tilde \sigma '_{ijk}}(\theta ){{ = [}}\sigma _{{{11}}}^ * (\theta ){]^{ - 2}}{{[}}\sigma _{{{11}}}^ * (\theta )] '{\hat \sigma _{ijk}}(\theta ) + {[\sigma _{{{11}}}^ * (\theta )]^{ - 1}}{\hat \sigma '_{ijk}}(\theta ){\text{.}} (53) 将式(52)代入式(11)可得
{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{\sigma } _{{\text{11}}}}(\theta ){\text{ = }}\frac{{{{\hat \sigma }_{{\text{11}}}}(\theta )}}{{{{[\sigma _{{\text{11}}}^ * (\theta )]}^2}}}, (54) 式中
{\hat \sigma _{{\text{11}}}}(\theta ){\text{ = 3}}{\hat \sigma ^2}_{rz{\text{1}}}(\theta ) + {\text{3}}{\hat \sigma ^2}_{\theta z{\text{1}}}(\theta ){\text{.}} (55) 再将式(51)代入式(53)，有
\sigma _{{\text{11}}}^ * (\theta ){\text{ = [}}{\hat \sigma _{{\text{11}}}}(\theta ){{\text{]}}^{{{(1 - m)/(2m)}}}}, (56) 回代式(52)和(53)，得
{\tilde \sigma _{ijk}}(\theta ){\text{ = [}}{\hat \sigma _{{\text{11}}}}(\theta ){{\text{]}}^{{{(m - 1)/(2m)}}}}{\hat \sigma _{ijk}}(\theta ),\tag{57a} {\tilde \sigma '_{ijk}}(\theta ){\text{ = }}\frac{{1 - m}}{{2m}}{{\text{[}}{\hat \sigma _{{\text{11}}}}(\theta ){\text{]}}^{{{(1 - 3m)/(2m)}}}}[{\hat \sigma _{11}}(\theta )]'{\hat \sigma _{ijk}}(\theta ) + {{\text{[}}{\hat \sigma _{{\text{11}}}}(\theta ){\text{]}}^{{{(1 - m)/(2m)}}}}{\hat \sigma _{ijk}}(\theta ){\text{.}} \tag{57b} 最后，将式(57)代入方程(31)得
\left\{ \begin{gathered} \left(\frac{{{\lambda _1}}}{m} + 1\right){{\tilde \sigma }_{rz1}}(\theta ) + \rho (\theta ){{\hat \sigma }_{\theta z1}}(\theta ) + {{\tilde \sigma '}_{\theta z1}}(\theta ) = 0, \hfill \\ {{\tilde \sigma }_{rz1}}(\theta ) - {{\tilde u}_{r1}}(\theta ) - {\lambda _1}{{\tilde u}_{r1}}(\theta ) = 0, \hfill \\ {{\tilde \sigma }_{\theta z1}}(\theta ) - {{\tilde u'}_{z1}}(\theta ) = 0, \hfill \\ \end{gathered} \right. (58) 式中
\rho (\theta ){\text{ = }}\frac{{{\text{1}} - m}}{{2m}}\frac{{[{{\hat \sigma }_{11}}(\theta )]'}}{{{{\hat \sigma }_{11}}(\theta )}}{\text{.}} (59) 在非线性常微分方程(58)中，非线性项
$ \rho (\theta ) $ 未知。因此在迭代求解中，首先令$ \rho (\theta ){\text{ = 0}} $ ，于是方程(58)变成线性特征值问题，采用插值矩阵法求解获得初始的$ \lambda _{\text{1}}^{(0)} $ ，$ \hat \sigma _{ij{\text{1}}}^{(0)}(\theta ) $ 和$ \tilde u_{ij{\text{1}}}^{(0)}(\theta ) $ 。将初值$ \hat \sigma _{ij{\text{1}}}^{(0)}(\theta ) $ 代入式(59)算得$ {\rho ^{(0)}}(\theta ) $ 作为迭代的初值。然后，使用插值矩阵迭代计算非线性常微分程(58)直至收敛。因为插值矩阵法求得的特征值和特征函数及其导数的精度相同，因此，我们选择相邻两次迭代计算的第1阶位移特征指数解的相对误差小于给定小量$ \varepsilon $ 作为判定迭代收敛的准则：\left| {\frac{{\lambda _{\text{1}}^{(z)} - \lambda _{\text{1}}^{(z - 1)}}}{{\lambda _{\text{1}}^{(z)}}}} \right| < \varepsilon ,\qquad z = 1,2,3,4, \cdot \cdot \cdot , (60) 可取
$\varepsilon {\text{ = 1}}{{\text{0}}^{-1}}$ 。4. 反平面V形切口塑性应力奇异性数值算例
对于一般反平面Ⅲ裂纹结构的应力奇异性控制方程求解，本文采用插值矩阵法先求解非线性方程(58)和相应的边界条件(33)或(34)，通过第3节迭代方法计算获得
$ {\lambda _{\text{1}}} $ ，$ {\hat \sigma _{ij{\text{1}}}}(\theta ) $ 和$ {\tilde u_{ij{\text{1}}}}(\theta ) $ 。然后将获得的$ {\lambda _{\text{1}}} $ ，$ {\hat \sigma _{ij{\text{1}}}}(\theta ) $ 和$ {\tilde u_{ij{\text{1}}}}(\theta ) $ 代入线性方程(32)，并结合边界条件(33)或(34)，再采用插值矩阵法求解获得高阶特征对$ {\lambda _{\text{1}}} $ ，$ {\hat \sigma _{ij{\text{1}}}}(\theta ) $ 和$ {\tilde u_{ijk}}(\theta ),{\text{ }}k > 1 $ 。本节主要讨论典型边界条件下反平面Ⅲ型切口/裂纹结构的前4阶弹塑性应力特征指数随不同分段数Q、幂硬化指数m和切口张角2α的变化规律。算例1 考虑两楔边自由的幂硬化材料反平面Ⅲ型裂纹结构(2α=0o)，见图1。
采用插值矩阵法依次求解控制方程(31)、(32)和相应的边界条件(33)。将
$\theta \in [ - {\text{π}},{\text{ }}{\text{π}}]$ 区间划分Q=20，40，80段，插值矩阵法获得反平面Ⅲ型裂纹尖端的前4阶弹塑性应力特征指数。该问题第1阶应力奇异指数的HRR解是$ {s_1} = - 1/(m + 1) $ 。图2给出了m＝5的反平面Ⅲ型裂纹第1阶应力指数$ {s_1} $ 插值矩阵法迭代计算的收敛情况。
图 2 插值矩阵法迭代计算Ⅲ型裂纹第1阶应力指数$ {s_1} $ 收敛情况Figure 2. First-order stress exponent s1 in the iterative calculation of the mode Ⅲ crack with the interpolating matrix method一般而言，V形切口/裂纹的第1阶弹塑性应力特征指数
$ {s_1} $ 小于0，当$ r \to 0 $ 时，切口尖端应力产生奇异性，因而第1阶应力特征指数和其相对应的角函数是表征尖端区域应力分布的主导项。图2显示，该例中插值矩阵法计算的$ {s_1} $ 在迭代次数z >15时收敛。当迭代次数z=25时，$ \theta \in [ - {\text{π}} , {\text{π}} ] $ 区间划分Q=20，40，80段计算的第1阶应力奇异指数s1分别为−0.166160，−0.166658，−0.166658，单机CPU计算时间均在10 s以内，说明当$ Q \geqslant {\text{40}} $ 时，本文方法计算的第1阶应力奇异指数至少有5位有效数字。Yang和Yuan等[14]通过hodograph变换和渐近分析的方法获得了反平面裂纹在剪切作用下尖端塑性应力奇异性的解析解。其第1阶应力主导项的公式为
{\tilde \sigma _{rz1}}(\theta ) = {\left[\frac{{\delta _1^2 + 1 - (\delta _1^2 - 1){\rm{cos}}(2{\chi _1})}}{2}\right]^{1/[{2m(1 - {\delta _1})}]}}\sin (\theta - {\chi _1}),\tag{61a} {\tilde \sigma _{\theta z1}}(\theta ) = {\left[\frac{{\delta _1^2 + 1 - (\delta _1^2 - 1){\rm{cos}}(2{\chi _1})}}{2}\right]^{1/{[2m(1 - {\delta _1})]}}}\cos (\theta - {\chi _1}), \tag{61b} {u_{z1}}(\theta ) = (m + 1){\left[\frac{{\delta _1^2 + 1 - (\delta _1^2 - 1){\rm{cos}}(2{\chi _1})}}{2}\right]^{1/{[2m(1 - {\delta _1})]}}}\sin (\theta - {\chi _1}), \tag{61c} 式中
{\delta _{\text{1}}}=\frac{{\text{1}}}{{{\text{2}}m}}\left(m - 1 - \sqrt {{{(m - 1)}^2} + 4m} \right) ,\; {\chi _{\text{1}}}{\text{ = }}\frac{{\text{1}}}{{\text{2}}}\left[\theta + {\rm{sin}}{^{ - 1}}\left(\frac{{m - 1}}{{m + 1}}\sin \theta \right)\right] {\text{.}} 当区间
$\theta \in [ - {\text{π}} ,{\text{π}}]$ 的分段数Q=20，40，80时，采用本文方法和文献[14]计算不同幂硬化指数m，反平面Ⅲ型裂纹尖端区域应力特征指数见表1，其中上标“#”代表弹性变形控制项，上标“×”代表非独立相关项。表 1 反平面Ⅲ型裂纹前4阶应力特征指数随硬化指数m的变化Table 1. Stress exponent$ {s_k} $ of the mode Ⅲ crack with various m values$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present (Q=20) −0.249 999 −0.166 666 −0.124 999 −0.090 907 −0.071 426 present (Q =40) −0.249 999 −0.166 667 −0.125 000 −0.090 909 −0.071 429 present (Q =80) −0.249 999 −0.166 667 −0.125 000 −0.090 909 −0.071 429 ref. [14] −0.250 000 −0.166 667 −0.125 000 −0.090 909 −0.071 429 $ {s_{\text{2}}} $ present (Q =20) 0.249999# 0.499998# 0.437490 0.365759 0.314434 present (Q =40) 0.249 999# 0.500 001# 0.435 761 0.363 758 0.312 255 present (Q =80) 0.249 999# 0.500 001# 0.435 666 0.363 646 0.312 133 ref. [14] 0.250 000# 0.500 000# 0.435 660 0.363 636 0.312 094 $ {s_{\text{3}}} $ present (Q =20) 0.574 150 0.501 559 0.624 995# 0.727 256# 0.700 294× present (Q =40) 0.572 948 0.500 086 0.625 000# 0.727 272# 0.695 939× present (Q =80) 0.572 879 0.500 004 0.625 000# 0.727 272# 0.695 695× ref. [14] 0.572 876 0.500 000 0.625 000# 0.727 273# 0.695 617× $ {s_{\text{4}}} $ present (Q =20) 0.749 997# 1.166 662# 0.999 979× 0.822 425× 0.785 686# present (Q =40) 0.749 997# 1.166 669# 0.996 522× 0.818 425× 0.785 719# present (Q =80) 0.749 997# 1.166 669# 0.996 332× 0.818 201× 0.785 719# ref. [14] 0.750 000# 1.166 667# 0.996 320× 0.818 182× 0.785 714# 表1显示当硬化指数m固定时，反平面Ⅲ型裂纹的前4阶应力特征指数于不同的分段数（Q=20, 40, 80）有细微变化。当
$ Q \geqslant {\text{40}} $ 时，应力特征指数保持稳定，不再随Q的增加而发生改变，说明采用插值矩阵法取Q=40计算反平面Ⅲ型裂纹的应力特征指数足够准确。表1还显示对于不同的硬化指数m，$ {s_{\text{1}}} < 0 $ ，$ {s_k} > 0 $ (k=2, 3, 4)，这意味着第1阶应力主导项是奇异的，后续高阶项是非奇异的。随着m的增大，$ {s_{\text{1}}} $ 逐渐变大并趋于零，但高阶应力特征指数$ {s_k} $ (k=2, 3, 4)逐渐减小，因此，反平面Ⅲ型裂纹尖端的塑性奇异性随m的增大而减弱。与文献[14]解析解比较可知，本文方法($ Q \geqslant {\text{40}} $ )获得前4阶应力特征指数的计算精度很高。与此同时，本文方法可获得同精度的位移和应力特征角函数。当m=3, 7, 13时，反平面裂纹的第1阶位移和应力特征角函数见图3。图3显示本文方法计算获得的第1阶位移和应力特征角函数和文献[14]结果吻合很好。本文方法还获得了同精度的高阶项位移和应力特征函数，见图4。图4显示反平面裂纹的第3、4阶特征角函数变化相对于第1、2阶特征角函数更为复杂。
图 3 幂硬化反平面Ⅲ型裂纹第1阶应力和位移特征角函数Figure 3. The 1st-order stress and displacement eigen-functions for the mode Ⅲcrack
图 4 幂硬化反平面Ⅲ型裂纹高阶应力和位移特征角函数，m=7Figure 4. The high-order stress and displacement eigen-functions for the mode Ⅲ crack with m=7算例2 两楔边自由的幂硬化材料反平面V形切口结构，见图1。
采用插值矩阵法依次求解控制方程(31)、(32)和相应的边界条件(33)。将
$\theta \in [ - {\text{π}}, {\text{π}}]$ 区间划分Q=20，40，80段，插值矩阵法获得反平面V形切口($\text{2}\alpha {\text{=30}}^{\circ },{\text{150}}^{\circ }$ )尖端的前4阶弹塑性应力特征指数，如表2、3所示。表 2 反平面V形切口($ {\text{2}}\alpha {\text{ = 3}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化Table 2. Stress exponent$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 3}}{{\text{0}}^ \circ } $ for various m values$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present(Q=20) −0.232204 −0.157527 −0.119464 −0.087783 −0.069413 present(Q=40) −0.232206 −0.157527 −0.119465 −0.087784 −0.069415 present(Q=80) −0.232206 −0.157527 −0.119465 −0.087784 −0.069415 ref. [17] −0.232206 −0.157564 −0.119502 −0.087824 −0.069457 $ {s_{\text{2}}} $ present(Q=20) 0.232204# 0.472581# 0.554001 0.472991 0.413208 present(Q=40) 0.232206# 0.472581# 0.551878 0.470411 0.410316 present(Q=80) 0.232206# 0.472581# 0.551876 0.470409 0.410314 ref. [17] − − 0.551851 0.470288 0.409963 $ {s_{\text{3}}} $ present(Q=20) 0.696612× 0.623794 0.597320# 0.702264# 0.763543# present(Q=40) 0.696618× 0.621968 0.597325# 0.702272# 0.763565# present(Q=80) 0.696618× 0.621972 0.597325# 0.702272# 0.763565# ref. [17] − 0.622051 − − − $ {s_{\text{4}}} $ present(Q=20) 0.700023 1.102689# 1.227466× 1.033765× 0.895829× present(Q=40) 0.698686 1.102689# 1.223221× 1.028606× 0.890047× present(Q=80) 0.698682 1.102689# 1.223217× 1.028602× 0.890043× ref. [17] 0.698604 − − − − 表 3 反平面V形切口($ {\text{2}}\alpha {\text{ = 15}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化Table 3. Stress exponent$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 15}}{{\text{0}}^ \circ } $ for various m values$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present(Q=20) −0.080378 −0.060350 −0.049587 −0.039967 −0.033919 present(Q=40) −0.080383 −0.060356 −0.049564 −0.039974 −0.033925 present(Q=80) −0.080383 −0.060356 −0.049593 −0.039974 −0.033925 ref. [17] −0.080384 −0.060369 −0.049614 −0.040000 −0.033956 $ {s_{\text{2}}} $ present(Q=20) 0.080378# 0.181051# 0.247935# 0.319736# 0.373109#` present(Q=40) 0.080383# 0.181068# 0.247970# 0.319792# 0.373175#` present(Q=80) 0.080383# 0.181068# 0.247965# 0.319792# 0.373175#` ref. [17] − − − − − $ {s_{\text{3}}} $ present(Q=20) 0.241134# 0.422452# 0.545457# 0.679439# 0.780137# present(Q=40) 0.241149# 0.422492# 0.545534# 0.679558# 0.780275# present(Q=80) 0.241149# 0.422492# 0.545523# 0.679558# 0.780275# ref. [17] − − − − − $ {s_{\text{4}}} $ present(Q=20) 1.589836 1.556131 1.520007 1.470114 1.426214 present(Q=40) 1.587821 1.553694 1.516597 1.464221 1.417315 present(Q=80) 1.587820 1.553693 1.516597 1.464220 1.417314 ref. [17] 1.587794 1.553562 1.516332 1.463944 1.416870 表 4 反平面V形切口($ {\text{2}}\alpha {\text{ = 6}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化(Q=40)Table 4. Stress exponent$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 6}}{{\text{0}}^ \circ } $ for various m values (Q=40)$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present −0.209034 −0.144787 −0.111379 −0.083025 −0.066275 ref. [17] −0.209035 −0.144820 −0.111421 −0.083068 −0.066323 $ {s_{\text{2}}} $ present 0.209034# 0.434361# 0.556895# 0.613636 0.545403 ref. [17] − − − 0.613578 0.545247 $ {s_{\text{3}}} $ present 0.627102# 0.775407 0.702317 0.664200# 0.729025# ref. [17] − 0.775369 0.702299 − − $ {s_{\text{4}}} $ present 0.851750 1.013509# 1.225169# 1.310297× 1.157081× ref. [17] 0.851659 − − − − 表 5 反平面V形切口($ {\text{2}}\alpha {\text{ = 9}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化(Q=40)Table 5. Stress exponent$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 9}}{{\text{0}}^ \circ } $ for various m values (Q=40)$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present −0.178394 −0.126571 −0.099142 −0.075360 −0.060999 ref. [17] −0.178395 −0.126599 −0.099179 −0.075403 −0.061047 $ {s_{\text{2}}} $ present 0.178394# 0.379713# 0.495710# 0.602880# 0.670989# ref. [17] − − − − − $ {s_{\text{3}}} $ present 0.535182# 0.885997# 0.899889 0.809746 0.736893 ref. [17] − − 0.899826 0.809670 0.736783 $ {s_{\text{4}}} $ present 1.040992 0.970537 1.090562# 1.281120# 1.402977# ref. [17] 1.040889 0.970465 − − − 表 6 反平面V形切口($ {\text{2}}\alpha {\text{ = 12}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化(Q=40)Table 6. Stress exponents$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 12}}{{\text{0}}^ \circ } $ for various m values (Q=40)$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present −0.137145 −0.099974 −0.080061 −0.060462 −0.051585 ref. [17] −0.137146 −0.100000 −0.080094 −0.062500 −0.051628 $ {s_{\text{2}}} $ present 0.137145# 0.299922# 0.400305# 0.483696# 0.567435# ref. [17] − − − − − $ {s_{\text{3}}} $ present 0.411435# 0.699818# 0.880671# 1.027854# 1.013736 ref. [17] − − − − 1.013494 $ {s_{\text{4}}} $ present 1.279549 1.222471 1.162660 1.082304 1.186455# ref. [17] 1.279434 1.222375 1.162546 1.082129 − Wang和Kuang[17]采用hodograph变换的方法求解了应力平面下损伤非线性切口问题，并获得了具有连续损伤的幂硬化材料V形切口应力特征指数的前2阶解，如表2~6所示。
表2和3显示当硬化指数m固定时，反平面V形切口(
$\text{2}\alpha {\text{=30}}^{\circ },150^{\circ }$ )的前4阶应力特征指数于不同的分段数(Q =20, 40, 80)有细微变化。当$ Q \geqslant 40 $ 时，应力特征指数保持稳定，不再随Q的增加而发生改变，说明采用本文方法计算反平面V形切口的应力特征指数取Q=40已足够。因此后续均取Q=40讨论不同张角$\text{2}\alpha {\text{=60}}^{\circ },\text{}{90}^{\circ },\text{}{120}^{\circ }$ 下反平面切口的前4阶应力特征指数随硬化指数m的变化，如表4~6所示。当Q=40，m=7，$\text{2}\alpha {\text{=60}}^{\circ },{\text{ 90}}^{\circ },\text{ }{120}^{\circ },\text{ }{150}^{\circ }$ 时，反平面V形切口结构的前2阶位移和应力特征角函数分别如图5和6所示。
图 5 反平面切口第1阶应力特征角函数(m=7)Figure 5. The 1st-order stress eigen-functions for the mode Ⅲ V-notch with m=7
图 6 反平面切口第2阶应力特征角函数(m=7)Figure 6. The 2nd-order stress eigen-functions for the mode Ⅲ V-notch with m=7表2~6和图5、6显示：当切口张角2α固定时，反平面V形切口应力特征指数
$ {s_k} $ (k=1, 2, 3, 4)的变化规律和算例1所示反平面Ⅲ型裂纹的变化规律相同，即$ {s_{\text{1}}} $ 随m的增大而增大，并逐渐趋近于0，且所有高阶应力特征指数$ {s_k} $ (k=2, 3, 4)均大于0，这意味着高阶项不会产生尖端的应力奇异性。表2~6和图5、6还显示：当m固定时，V形切口尖端的前4阶应力特征指数
$ {s_k} $ (k=2, 3, 4)随切口张角$ {\text{2}}\alpha $ 的增大而增大，因此V形切口尖端的塑性应力奇异性随张角$ {\text{2}}\alpha $ 的增大逐渐减弱。文献[17]给出了反平面V形切口张角$\text{2}\alpha {\text{=30}}^{\circ },{\text{60}}^{\circ },{\text{90}}^{\circ },\text{}{120}^{\circ },\text{}{150}^{\circ }$ 的前2阶应力特征指数，但未考虑弹性项和遗漏塑性项的贡献，其应力特征指数有漏项，因此其应力特征指数的第2阶不一定是真实的第2阶。而本文方法考虑了弹性项和遗漏塑性项的贡献，给出了前4阶应力特征指数，其计算结果较文献[17]的结果更加完整和准确，因此本文方法可准确计算一般V形切口在反平面剪切作用下的塑性应力奇异指数和相对应的位移和应力特征向量。算例3 一边自由一边固结的幂硬化材料反平面V形切口结构，见图1。
采用插值矩阵法依次求解控制方程(31)、(32)和相应的边界条件(33)、(34)，
$ \theta \in [ - {\text{π}} , {\text{π}} ] $ 区间划分Q=40段，采用插值矩阵法获得反平面V形切口($\text{2}\alpha {\text{=0}}^{\circ },{\text{30}}^{\circ },{\text{90}}^{\circ },\text{}{150}^{\circ }$ )尖端的前4阶弹塑性应力特征指数随硬化指数m的变化规律，见表7，我们未发现其他文献给出该例的应力特征解。当$ {\text{2}}\alpha {\text{ = }}{{\text{0}}^ \circ } $ ，m=3, 7, 13时，反平面Ⅲ型裂纹结构的第1阶位移和应力特征角函数见图7。当m=5，$\text{2}\alpha {\text{=0}}^{\circ },{\text{30}}^{\circ },{\text{90}}^{\circ },\text{}{150}^{\circ }$ 时，反平面V形切口结构的第1阶应力特征角函数见图8。表 7 反平面V形切口前4阶应力特征指数随硬化指数m和张角2α的变化 ( Q =40)Table 7. Stress exponents$ {s_k} $ of the mode Ⅲ V-notch for various m values and$ {\text{2}}\alpha $ (Q=40)$ {s_k} $ $ {\text{2}}\alpha $/(°) m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ 0 −0.316987 −0.194766 −0.140313 −0.098797 −0.076224 30 −0.313242 −0.193473 −0.139669 −0.098487 −0.076043 90 −0.300944 −0.188989 −0.137387 −0.097374 −0.075386 150 −0.274774 −0.178251 −0.131633 −0.094459 −0.073631 $ {s_{\text{2}}} $ 0 −0.049037 −0.026170 −0.017810 −0.012034 −0.009085 30 −0.000702 0.008722 0.009252 0.008153 0.006996 90 0.138739 0.116011 0.095654 0.074849 0.061274 150 0.274774# 0.321572 0.272393 0.220346 0.184745 $ {s_{\text{3}}} $ 0 0.218913× 0.142426× 0.104693× 0.074729× 0.058054× 30 0.311838× 0.210917× 0.158173× 0.114793× 0.090035× 90 0.300944× 0.421011× 0.328695× 0.247072× 0.197934× 150 0.380293 0.534753# 0.658165# 0.535151× 0.443121× $ {s_{\text{4}}} $ 0 0.285422 0.213317 0.169289 0.129282 0.104644 30 0.313242# 0.285697 0.228736 0.176303 0.143628 90 0.638404 0.50008 0.410023 0.323996 0.268649 150 0.824322× 0.821395× 0.676419× 0.623227 0.531023 
图 7 反平面裂纹第1阶位移和应力特征角函数Figure 7. The 1st-order stress and displacement eigen-functions for the mode Ⅲ V-notch
图 8 反平面V形切口第1阶应力特征角函数(m=5)Figure 8. The 1st-order stress and displacement eigen-functions for the mode Ⅲ V-notch (m=5)表7显示：当切口张角2α固定时，反平面V形切口在自由-固结边界下的应力特征指数
$ {s_{\text{1}}} $ 随m的增大而增大，并逐渐趋近于0，与自由-自由边界不一样的是并非所有高阶应力特征指数$ {s_k} $ (k=2, 3, 4)均大于0，对于自由-固结边界下反平面裂纹($ {\text{2}}\alpha {\text{ = }}{{\text{0}}^ \circ } $ )，其第2阶$ {s_{\text{2}}} < 0 $ ，尽管s2接近0，这意味着除第一阶主导项外，第二阶应力特征指数也会产生尖端的应力奇异性，其对于应力场的贡献应考虑。表7还显示当m固定时，V形切口尖端的前4阶应力特征指数$ {s_k} $ (k=1, 2, 3, 4)随切口张角$ {\text{2}}\alpha $ 的增大而增大，即反平面V形切口结构在自由-固结边界下尖端的塑性应力奇异性随张角$ {\text{2}}\alpha $ 的增大逐渐减弱。5. 结 论
本文对于一般反平面V形切口/裂纹结构的奇异应力场问题，建立了一种新的求解方法，并研制了其相应的计算程序，给出了典型边界条件下反平面V形切口/裂纹结构3个算例的计算结果，所得结论如下：
1） 反平面V形切口/裂纹结构尖端附近区域的位移和应力场由渐近级数展开式表达，然后将其代入弹塑性力学平衡方程，可得关于反平面V形切口/裂纹尖端附近区域应力奇异性特征值方程的常微分方程组。采用插值矩阵法计算常微分方程组，可得反平面V形切口/裂纹结构尖端附近区域完整的弹塑性应力特征指数和同精度的特征角函数。
2） 本文方法获得的反平面裂纹尖端塑性第1阶应力奇异解和文献结果吻合较好，且本文方法可获得典型边界条件下不同张角反平面V形切口弹塑性应力解的主导项和高阶项。数值算例结果表明了本文方法求解反平面V形切口/裂纹结构弹塑性应力渐近解的精确性和高效性。
致谢 本文作者衷心感谢安徽建筑大学博士启动基金(2020QDZ208)对本文的资助。
-
图 2 插值矩阵法迭代计算Ⅲ型裂纹第1阶应力指数
$ {s_1} $ 收敛情况Figure 2. First-order stress exponent s1 in the iterative calculation of the mode Ⅲ crack with the interpolating matrix method
图 3 幂硬化反平面Ⅲ型裂纹第1阶应力和位移特征角函数
Figure 3. The 1st-order stress and displacement eigen-functions for the mode Ⅲcrack
图 4 幂硬化反平面Ⅲ型裂纹高阶应力和位移特征角函数，m=7
Figure 4. The high-order stress and displacement eigen-functions for the mode Ⅲ crack with m=7
图 5 反平面切口第1阶应力特征角函数(m=7)
Figure 5. The 1st-order stress eigen-functions for the mode Ⅲ V-notch with m=7
图 6 反平面切口第2阶应力特征角函数(m=7)
Figure 6. The 2nd-order stress eigen-functions for the mode Ⅲ V-notch with m=7
图 7 反平面裂纹第1阶位移和应力特征角函数
Figure 7. The 1st-order stress and displacement eigen-functions for the mode Ⅲ V-notch
图 8 反平面V形切口第1阶应力特征角函数(m=5)
Figure 8. The 1st-order stress and displacement eigen-functions for the mode Ⅲ V-notch (m=5)
表 1 反平面Ⅲ型裂纹前4阶应力特征指数随硬化指数m的变化
Table 1. Stress exponent
$ {s_k} $ of the mode Ⅲ crack with various m values$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present (Q=20) −0.249 999 −0.166 666 −0.124 999 −0.090 907 −0.071 426 present (Q =40) −0.249 999 −0.166 667 −0.125 000 −0.090 909 −0.071 429 present (Q =80) −0.249 999 −0.166 667 −0.125 000 −0.090 909 −0.071 429 ref. [14] −0.250 000 −0.166 667 −0.125 000 −0.090 909 −0.071 429 $ {s_{\text{2}}} $ present (Q =20) 0.249999# 0.499998# 0.437490 0.365759 0.314434 present (Q =40) 0.249 999# 0.500 001# 0.435 761 0.363 758 0.312 255 present (Q =80) 0.249 999# 0.500 001# 0.435 666 0.363 646 0.312 133 ref. [14] 0.250 000# 0.500 000# 0.435 660 0.363 636 0.312 094 $ {s_{\text{3}}} $ present (Q =20) 0.574 150 0.501 559 0.624 995# 0.727 256# 0.700 294× present (Q =40) 0.572 948 0.500 086 0.625 000# 0.727 272# 0.695 939× present (Q =80) 0.572 879 0.500 004 0.625 000# 0.727 272# 0.695 695× ref. [14] 0.572 876 0.500 000 0.625 000# 0.727 273# 0.695 617× $ {s_{\text{4}}} $ present (Q =20) 0.749 997# 1.166 662# 0.999 979× 0.822 425× 0.785 686# present (Q =40) 0.749 997# 1.166 669# 0.996 522× 0.818 425× 0.785 719# present (Q =80) 0.749 997# 1.166 669# 0.996 332× 0.818 201× 0.785 719# ref. [14] 0.750 000# 1.166 667# 0.996 320× 0.818 182× 0.785 714# 
下载: 导出CSV
表 2 反平面V形切口(
$ {\text{2}}\alpha {\text{ = 3}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化Table 2. Stress exponent
$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 3}}{{\text{0}}^ \circ } $ for various m values$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present(Q=20) −0.232204 −0.157527 −0.119464 −0.087783 −0.069413 present(Q=40) −0.232206 −0.157527 −0.119465 −0.087784 −0.069415 present(Q=80) −0.232206 −0.157527 −0.119465 −0.087784 −0.069415 ref. [17] −0.232206 −0.157564 −0.119502 −0.087824 −0.069457 $ {s_{\text{2}}} $ present(Q=20) 0.232204# 0.472581# 0.554001 0.472991 0.413208 present(Q=40) 0.232206# 0.472581# 0.551878 0.470411 0.410316 present(Q=80) 0.232206# 0.472581# 0.551876 0.470409 0.410314 ref. [17] − − 0.551851 0.470288 0.409963 $ {s_{\text{3}}} $ present(Q=20) 0.696612× 0.623794 0.597320# 0.702264# 0.763543# present(Q=40) 0.696618× 0.621968 0.597325# 0.702272# 0.763565# present(Q=80) 0.696618× 0.621972 0.597325# 0.702272# 0.763565# ref. [17] − 0.622051 − − − $ {s_{\text{4}}} $ present(Q=20) 0.700023 1.102689# 1.227466× 1.033765× 0.895829× present(Q=40) 0.698686 1.102689# 1.223221× 1.028606× 0.890047× present(Q=80) 0.698682 1.102689# 1.223217× 1.028602× 0.890043× ref. [17] 0.698604 − − − −
下载: 导出CSV
表 3 反平面V形切口(
$ {\text{2}}\alpha {\text{ = 15}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化Table 3. Stress exponent
$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 15}}{{\text{0}}^ \circ } $ for various m values$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present(Q=20) −0.080378 −0.060350 −0.049587 −0.039967 −0.033919 present(Q=40) −0.080383 −0.060356 −0.049564 −0.039974 −0.033925 present(Q=80) −0.080383 −0.060356 −0.049593 −0.039974 −0.033925 ref. [17] −0.080384 −0.060369 −0.049614 −0.040000 −0.033956 $ {s_{\text{2}}} $ present(Q=20) 0.080378# 0.181051# 0.247935# 0.319736# 0.373109#` present(Q=40) 0.080383# 0.181068# 0.247970# 0.319792# 0.373175#` present(Q=80) 0.080383# 0.181068# 0.247965# 0.319792# 0.373175#` ref. [17] − − − − − $ {s_{\text{3}}} $ present(Q=20) 0.241134# 0.422452# 0.545457# 0.679439# 0.780137# present(Q=40) 0.241149# 0.422492# 0.545534# 0.679558# 0.780275# present(Q=80) 0.241149# 0.422492# 0.545523# 0.679558# 0.780275# ref. [17] − − − − − $ {s_{\text{4}}} $ present(Q=20) 1.589836 1.556131 1.520007 1.470114 1.426214 present(Q=40) 1.587821 1.553694 1.516597 1.464221 1.417315 present(Q=80) 1.587820 1.553693 1.516597 1.464220 1.417314 ref. [17] 1.587794 1.553562 1.516332 1.463944 1.416870
下载: 导出CSV
表 4 反平面V形切口(
$ {\text{2}}\alpha {\text{ = 6}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化(Q=40)Table 4. Stress exponent
$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 6}}{{\text{0}}^ \circ } $ for various m values (Q=40)$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present −0.209034 −0.144787 −0.111379 −0.083025 −0.066275 ref. [17] −0.209035 −0.144820 −0.111421 −0.083068 −0.066323 $ {s_{\text{2}}} $ present 0.209034# 0.434361# 0.556895# 0.613636 0.545403 ref. [17] − − − 0.613578 0.545247 $ {s_{\text{3}}} $ present 0.627102# 0.775407 0.702317 0.664200# 0.729025# ref. [17] − 0.775369 0.702299 − − $ {s_{\text{4}}} $ present 0.851750 1.013509# 1.225169# 1.310297× 1.157081× ref. [17] 0.851659 − − − −
下载: 导出CSV
表 5 反平面V形切口(
$ {\text{2}}\alpha {\text{ = 9}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化(Q=40)Table 5. Stress exponent
$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 9}}{{\text{0}}^ \circ } $ for various m values (Q=40)$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present −0.178394 −0.126571 −0.099142 −0.075360 −0.060999 ref. [17] −0.178395 −0.126599 −0.099179 −0.075403 −0.061047 $ {s_{\text{2}}} $ present 0.178394# 0.379713# 0.495710# 0.602880# 0.670989# ref. [17] − − − − − $ {s_{\text{3}}} $ present 0.535182# 0.885997# 0.899889 0.809746 0.736893 ref. [17] − − 0.899826 0.809670 0.736783 $ {s_{\text{4}}} $ present 1.040992 0.970537 1.090562# 1.281120# 1.402977# ref. [17] 1.040889 0.970465 − − −
下载: 导出CSV
表 6 反平面V形切口(
$ {\text{2}}\alpha {\text{ = 12}}{{\text{0}}^ \circ } $ )前4阶应力特征指数随硬化指数m的变化(Q=40)Table 6. Stress exponents
$ {s_k} $ of the mode Ⅲ V-notch with$ {\text{2}}\alpha {\text{ = 12}}{{\text{0}}^ \circ } $ for various m values (Q=40)$ {s_k} $ method m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ present −0.137145 −0.099974 −0.080061 −0.060462 −0.051585 ref. [17] −0.137146 −0.100000 −0.080094 −0.062500 −0.051628 $ {s_{\text{2}}} $ present 0.137145# 0.299922# 0.400305# 0.483696# 0.567435# ref. [17] − − − − − $ {s_{\text{3}}} $ present 0.411435# 0.699818# 0.880671# 1.027854# 1.013736 ref. [17] − − − − 1.013494 $ {s_{\text{4}}} $ present 1.279549 1.222471 1.162660 1.082304 1.186455# ref. [17] 1.279434 1.222375 1.162546 1.082129 −
下载: 导出CSV
表 7 反平面V形切口前4阶应力特征指数随硬化指数m和张角2α的变化 ( Q =40)
Table 7. Stress exponents
$ {s_k} $ of the mode Ⅲ V-notch for various m values and$ {\text{2}}\alpha $ (Q=40)$ {s_k} $ $ {\text{2}}\alpha $/(°) m=3 m=5 m=7 m=10 m=13 $ {s_{\text{1}}} $ 0 −0.316987 −0.194766 −0.140313 −0.098797 −0.076224 30 −0.313242 −0.193473 −0.139669 −0.098487 −0.076043 90 −0.300944 −0.188989 −0.137387 −0.097374 −0.075386 150 −0.274774 −0.178251 −0.131633 −0.094459 −0.073631 $ {s_{\text{2}}} $ 0 −0.049037 −0.026170 −0.017810 −0.012034 −0.009085 30 −0.000702 0.008722 0.009252 0.008153 0.006996 90 0.138739 0.116011 0.095654 0.074849 0.061274 150 0.274774# 0.321572 0.272393 0.220346 0.184745 $ {s_{\text{3}}} $ 0 0.218913× 0.142426× 0.104693× 0.074729× 0.058054× 30 0.311838× 0.210917× 0.158173× 0.114793× 0.090035× 90 0.300944× 0.421011× 0.328695× 0.247072× 0.197934× 150 0.380293 0.534753# 0.658165# 0.535151× 0.443121× $ {s_{\text{4}}} $ 0 0.285422 0.213317 0.169289 0.129282 0.104644 30 0.313242# 0.285697 0.228736 0.176303 0.143628 90 0.638404 0.50008 0.410023 0.323996 0.268649 150 0.824322× 0.821395× 0.676419× 0.623227 0.531023
下载: 导出CSV
-
[1] 李聪, 牛忠荣, 胡宗军, 等. 求解双材料裂纹结构全域应力场的扩展边界元法[J]. 应用数学和力学, 2019, 40(8): 926-937. (LI Cong, NIU Zhongrong, HU Zongjun, et al. Computation of total stress fields for cracked bi-material structures with the extended boundary element method[J]. Applied Mathematics and Mechanics, 2019, 40(8): 926-937.(in Chinese) [2] HUTCHINSON J W. Singular behavior at the end of a tensile crack in a hardening material[J]. Journal of the Mechanics and Physics of Solids, 1968, 16(1): 13-31. doi: 10.1016/0022-5096(68)90014-8 [3] RICE J R, ROSENGREN G F. Plane strain deformation near a crack tip in a power-law hardening material[J]. Journal of the Mechanics and Physics of Solids, 1968, 16(1): 1-12. doi: 10.1016/0022-5096(68)90013-6 [4] KUANG Z, XU X. Stress and strain fields at the tip of a sharp V-notch in a power-hardening material[J]. International Journal of Fracture, 1987, 35: 39-53. doi: 10.1007/BF00034533 [5] XIA L, WANG T C. Singular behaviour near the tip of a sharp V-notch in a power law hardening material[J]. International Journal of Fracture, 1993, 59(1): 83-93. doi: 10.1007/BF00032219 [6] LI Y C, WANG T C. High-order asymptotic field of tensile plane-strain nonlinear crack problems[J]. Scientia Sinica(Series A) , 1986, 29(9): 941-955. [7] SHARMA S M, ARAVAS N. On the development of variable-separable asymptotic elastoplastic solutions for interfacial cracks[J]. International Journal of Solids & Structures, 1993, 30(5): 695-723. [8] XIA L, WANG T C, SHIH C F. Higher-order analysis of crack tip fields in elastic power-law hardening materials[J]. Journal of the Mechanics and Physics of Solids, 1993, 41(4): 665-687. doi: 10.1016/0022-5096(93)90022-8 [9] YUAN F G, YANG S. Crack-tip fields in elastic-plastic material under plane stress mode Ⅰ loading[J]. International Journal of Fracture, 1997, 85(2): 131-155. doi: 10.1023/A:1007361116709 [10] YANG S, CHAO Y J, SUTTON M A. Higher order asymptotic crack tip fields in a power-law hardening material[J]. Engineering Fracture Mechanics, 1993, 45(1): 1-20. doi: 10.1016/0013-7944(93)90002-A [11] CHAO Y J, YANG S. Higher order crack tip fields and its implication for fracture of solids under mode Ⅱ conditions[J]. Engineering Fracture Mechanics, 1996, 55(5): 777-794. doi: 10.1016/0013-7944(96)00054-9 [12] RICE J R. Stresses due to a sharp notch in a work-hardening elastic-plastic material loaded by longitudinal shear[J]. Journal of Applied Mechanics, 1967, 34(2): 287. doi: 10.1115/1.3607681 [13] AMAZIGO J C. Fully plastic crack in an infinite body under anti-plane shear[J]. International Journal of Solids and Structures, 1974, 10(9): 1003-1015. doi: 10.1016/0020-7683(74)90008-0 [14] YANG S, YUAN F G, CAI X. Higher orderasymptotic elastic-plastic crack-tip fields under antiplane shear[J]. Engineering Fracture Mechanics, 1996, 54(3): 405-422. doi: 10.1016/0013-7944(95)00191-3 [15] YANG S, YUAN F G, CHIANG M Y M. Analytical forms of higher-order asymptotic elastic-plastic crack-tip fields in a linear hardening material under antiplane shear[J]. International Journal of Fracture, 1996, 80(1): 59-71. doi: 10.1007/BF00036480 [16] YUAN F G, YANG S. Analytical solutions of fully plastic crack-tip higher order fields under antiplane shear[J]. International Journal of Fracture, 1994, 69(1): 1-26. [17] WANG T J, KUANG Z B. Higher order asymptotic solutions of V-notch tip fields for damaged nonlinear materials under antiplane shear loading[J]. International Journal of Fracture, 1999, 96(4): 303-329. doi: 10.1023/A:1018657316810 [18] ZAPPALORTO M, LAZZARIN P. Analytical study of the elastic-plastic stress fields ahead of parabolic notches under antiplane shear loading[J]. International Journal of Fracture, 2007, 148(2): 139-154. doi: 10.1007/s10704-008-9185-7 [19] LAZZARIN P, ZAPPALORTO M. Plastic notch stress intensity factors for pointed V-notches under antiplane shear loading[J]. International Journal of Fracture, 2008, 152(1): 1-25. doi: 10.1007/s10704-008-9260-0 [20] ZAPPALORTO M, LAZZARIN P. A unified approach to the analysis of nonlinear stress and strain fields ahead of mode Ⅲ-loaded notches and cracks[J]. International Journal of Solids and Structures, 2010, 47(6): 851-864. doi: 10.1016/j.ijsolstr.2009.11.021 [21] ARAVAS N, BLAZO D H. Higher order terms in asymptotic elastoplastic mode Ⅲ crack tip solutions[J]. Acta Mechanica, 1991, 90(1/4): 139-153. [22] LOGHIN A, ZHANG N, JOSEPH P F. A nonlinear finite element eigenanalysis of antiplane shear including higher order terms[J]. Engineering Fracture Mechanics, 2000, 66(5): 441-454. doi: 10.1016/S0013-7944(00)00031-X [23] ZHANG N, JOSEPH P F. A nonlinear finite element eigenanalysis of singular stress fields in bimaterial wedges for plane strain[J]. International Journal of Fracture, 1998, 94(3): 299-319. [24] LOGHIN A, JOSEPH P F. Asymptotic solutions for mixed mode loading of cracks and wedges in power law hardening materials[J]. Engineering Fracture Mechanics, 2001, 68(14): 1511-1534. doi: 10.1016/S0013-7944(01)00050-9 [25] PATWARDHAN P A, NALAVDE R A, KUJAWSKI D. An estimation of Ramberg-Osgood constants for materials with and without Loder’s strain using yield and ultimate strengths[J]. Procedia Structure Integrity, 2019, 17: 750-757. doi: 10.1016/j.prostr.2019.08.100 [26] NIU Z R, GE D L, CHENG C Z, et al. Evaluation of the stress singularities of plane V-notches in bonded dissimilar materials[J]. Applied Mathematical Modelling, 2009, 33(3): 1776-1792. doi: 10.1016/j.apm.2008.03.007 期刊类型引用(0)
其他类型引用(3)
-
下载:
计量
- 文章访问数: 963
- HTML全文浏览量: 371
- PDF下载量: 46
- 被引次数: 3
渝公网安备50010802005915号