.1316 北京科技大学学报 第35卷 进而求得裂纹尖端总的应力强度因子表达式.文中3.3有限元中引入判据 首先计算单个氢原子与基体应力场的作用能，而后 在用有限元分析扩展过程时需要确定氢气团 推导出整个氢气团与基体相互作用.计算氢含量时 的尺寸，在文献[门]中将氢气团简化成直径为x0的 所用扩散通量如下所示： 圆形区域，本文中相应的将缺陷旁边平均浓度最高 j=-器B (12) 的单元视为氢气团.如果单元尺寸足够小并且形状 规则（如矩形单元和正方体体单元），这一简化是可 本文认为这个过程当中没有考虑氢原子之间的 以接受的. 相互作用.即只考虑了应力场而忽略了浓度场对氢 本文建议用梯形法将公式（⑨）进行离散的原因 原子扩散的影响.加上了浓度场项后，式(12)修 在于，本文计算结果显示用Sofronis和McMeek- 正为 ig模型计算的结果在加载初期和应力梯度或浓度 CoD. j=-DNcm-k7TEt· (13) 梯度较高的地方会有较大的震荡，往往会造成在加 载初期裂纹就萌生，而后一直稳定，KH值无法达 式(12)和(13)中，了为氢原子流密度，D为扩散系 到KHC.即使单元划分的足够细密，在加载初期也 数，co为氢含量，k气体常数，T为温度，Emt为 无法避免这一现象的产生.另一方面，由于要模拟 氢原子与I型裂纹之间的交互作用力.最后文献[可] 裂纹的扩展，所以在裂纹扩展的方向上（这个方向 计算出氢气团引起的裂纹尖端应力强度因子的改变 可以在建模时判断出来)的单元应当划分均匀且足 量KH,应用到本文当中为 够小，这就保证了梯形法计算的稳定同时避免了计 3G C Ku=-vvā(e-8)Ha,o}（1④ 算结果的震荡 其中，G和v为剪切弹性模量和泊松比，C为氢气 4有限元模拟裂纹萌生和扩展 团中的平均浓度，cm为晶包中单个氢原子的质量分 本文模拟氢在循环载荷作用下向钢中一个椭 数，取值为1/3，a为晶格参数取值为1.0387，e0为单 圆形空洞扩散和富集，最终引起疲劳裂纹萌生和扩 个原子引起的晶格畸变应变，xo为氢气团直径，H 展.图2所示为有限元模型和边界条件，空洞长轴 为几何参数函数 长度为40um,短轴长度为10m.以空洞中心为 3.2裂纹萌生和扩展判据 原点，椭圆长、短轴方向分别为模型的x轴和y轴， 结合本文要求，将材料中的缺陷视为一个短裂 因载荷对称所以只取第一象限一个边长为500m 纹，通过式(14)可以计算出在材料缺陷边缘富集的 的正方形区域为计算模型.图2中C%为初始氢含 氢原子引起的缺陷处应力强度因子的增量KH这 量，疲劳载荷沿y方向加载，F为疲劳载荷，J为 样就可以建立缺陷处裂纹萌生和扩展的判据.以I 扩散通量，U为位移，t为时间.每一个载荷步中首 型裂纹为例，K为缺陷处应力强度因子，KC为I 先需要计算出节点平均应力om,然后将其作为计 型裂纹断裂韧度，则裂纹萌生和扩展的判据为 算氢原子扩散的应力初始条件，所以用于计算应力 和浓度的模型应当保持一致.虽然可以通过插值或 K1+KH≥KIC, 者外推等方法计算任意位置的应力值，但因为有限 即 元计算的节点处应力结果已经是不精确的，不应该 KH≥KIC-KI. (15) 再降低其精度.对于各项同性均匀材料，在建立模型 时应当先判断出裂纹扩展的方向保证节点在此方向 令 上彼此之间拥有最短距离.另一方面，在此方向上 KHC KIC -K1, (16) 的单元将用于计算氢气团平均浓度，所以单元尺寸 氢致裂纹萌生和扩展的判据为 应当足够小且形状规则.以四边形单元为例，采用 国际单位制，根据扩展判据的推导过程可知单元边 KH≥KHC· (17） 长应当小于1.8×10-5m. 通过式(14)可知KH大小与氢气团平均浓度成 4.1问题简化和边界条件 正比，当氢原子在应力场作用下向缺陷处富集，使 (1)本文有限元模型是二维的，模拟较低的疲 缺陷处浓度越来越高，相应的K红自然也会越来越 劳载荷作用下氢原子向材料中一个椭圆形空洞处 大，当增大到KHC时裂纹将扩展，扩展距离为xo 富集，空洞长轴与加载方向垂直，所以疲劳裂纹将 而后氢将向新的裂尖富集，重复先前的过程 在空洞长轴顶端萌生并沿这个方向扩展.根据裂尖· 1316 · 北 京 科 技 大 学 学 报 第 35 卷 进而求得裂纹尖端总的应力强度因子表达式. 文中 首先计算单个氢原子与基体应力场的作用能，而后 推导出整个氢气团与基体相互作用. 计算氢含量时 所用扩散通量如下所示： j = − c0D kT ∇Eint . (12) 本文认为这个过程当中没有考虑氢原子之间的 相互作用. 即只考虑了应力场而忽略了浓度场对氢 原子扩散的影响. 加上了浓度场项后，式 (12) 修 正为 j = −D∇c0 − c0D kT ∇Eint . (13) 式(12) 和 (13) 中，j 为氢原子流密度，D 为扩散系 数，c0 为氢含量，k 气体常数，T 为温度，Eint 为 氢原子与 I 型裂纹之间的交互作用力. 最后文献 [7] 计算出氢气团引起的裂纹尖端应力强度因子的改变 量 KH，应用到本文当中为 KH = 3G 1 − ν · C¯ cm √ πa ¡ ε 0 22 − ε 0 11¢ H (a, x0). (14) 其中，G 和 ν 为剪切弹性模量和泊松比，C¯ 为氢气 团中的平均浓度，cm 为晶包中单个氢原子的质量分 数，取值为 1/3，a 为晶格参数取值为 1.0387，ε 0 为单 个原子引起的晶格畸变应变，x0 为氢气团直径，H 为几何参数函数. 3.2 裂纹萌生和扩展判据 结合本文要求，将材料中的缺陷视为一个短裂 纹，通过式 (14) 可以计算出在材料缺陷边缘富集的 氢原子引起的缺陷处应力强度因子的增量 KH. 这 样就可以建立缺陷处裂纹萌生和扩展的判据. 以 I 型裂纹为例，KI 为缺陷处应力强度因子，KIC 为 I 型裂纹断裂韧度，则裂纹萌生和扩展的判据为 KI + KH>KIC ， 即 KH>KIC − KI . (15) 令 KHC = KIC − KI ， (16) 氢致裂纹萌生和扩展的判据为 KH>KHC . (17) 通过式 (14) 可知 KH 大小与氢气团平均浓度成 正比，当氢原子在应力场作用下向缺陷处富集，使 缺陷处浓度越来越高，相应的 KH 自然也会越来越 大，当增大到 KHC 时裂纹将扩展，扩展距离为 x0. 而后氢将向新的裂尖富集，重复先前的过程. 3.3 有限元中引入判据 在用有限元分析扩展过程时需要确定氢气团 的尺寸，在文献 [7] 中将氢气团简化成直径为 x0 的 圆形区域，本文中相应的将缺陷旁边平均浓度最高 的单元视为氢气团. 如果单元尺寸足够小并且形状 规则 (如矩形单元和正方体体单元)，这一简化是可 以接受的. 本文建议用梯形法将公式 (9) 进行离散的原因 在于，本文计算结果显示用 Sofronis 和 McMeek￾ing 模型计算的结果在加载初期和应力梯度或浓度 梯度较高的地方会有较大的震荡，往往会造成在加 载初期裂纹就萌生，而后一直稳定，KH 值无法达 到 KHC. 即使单元划分的足够细密，在加载初期也 无法避免这一现象的产生. 另一方面，由于要模拟 裂纹的扩展，所以在裂纹扩展的方向上 (这个方向 可以在建模时判断出来) 的单元应当划分均匀且足 够小，这就保证了梯形法计算的稳定同时避免了计 算结果的震荡. 4 有限元模拟裂纹萌生和扩展 本文模拟氢在循环载荷作用下向钢中一个椭 圆形空洞扩散和富集，最终引起疲劳裂纹萌生和扩 展. 图 2 所示为有限元模型和边界条件，空洞长轴 长度为 40 µm，短轴长度为 10 µm. 以空洞中心为 原点，椭圆长、短轴方向分别为模型的 x 轴和 y 轴， 因载荷对称所以只取第一象限一个边长为 500 µm 的正方形区域为计算模型. 图 2 中 C0 为初始氢含 量，疲劳载荷沿 y 方向加载，F 为疲劳载荷，J 为 扩散通量，U 为位移，t 为时间. 每一个载荷步中首 先需要计算出节点平均应力 σm，然后将其作为计 算氢原子扩散的应力初始条件，所以用于计算应力 和浓度的模型应当保持一致. 虽然可以通过插值或 者外推等方法计算任意位置的应力值，但因为有限 元计算的节点处应力结果已经是不精确的，不应该 再降低其精度. 对于各项同性均匀材料，在建立模型 时应当先判断出裂纹扩展的方向保证节点在此方向 上彼此之间拥有最短距离. 另一方面，在此方向上 的单元将用于计算氢气团平均浓度，所以单元尺寸 应当足够小且形状规则. 以四边形单元为例，采用 国际单位制，根据扩展判据的推导过程可知单元边 长应当小于 1.8×10−5 m. 4.1 问题简化和边界条件 (1) 本文有限元模型是二维的，模拟较低的疲 劳载荷作用下氢原子向材料中一个椭圆形空洞处 富集，空洞长轴与加载方向垂直，所以疲劳裂纹将 在空洞长轴顶端萌生并沿这个方向扩展. 根据裂尖