孙兆珑， 任克亮， 纪华， 丁丽宏， 张亚婷， 王婷， 杨佳， 李佶敏

(宁夏大学物理与电子电气工程学院， 银川 750021)

生产生活中广泛使用的金属材料，大多是由不规则形状的众多晶粒组成的多晶体，晶粒之间的边界称为晶粒边界或晶粒间界，简称晶界[1]。在多晶体的形成过程中，晶界区域总会包含一些杂质原子，这些杂质原子会削弱或强化晶界处的结合力,进而影响材料的整体性能[2]。

对于多晶材料来说，晶界结构的强弱与分布情况是影响其材料性能的重要因素之一，针对这一问题，学者们分别从实验、数值模拟等手段对其进行研究分析。路君等[3]通过总结前人实验结果，讨论分析了晶粒大小与各类多晶材料力学性能之间的关系，给出了Hall-Petch公式的适用范围。辛存[4]通过晶体塑性有限元法(crystal plasticity finite element modeling，CPFEM)对Voronoi多晶体模型收拉伸载荷的情况进行了模拟，发现晶界处会产生应力集中、接触压力剧变等现象。戴澍等[5]将晶粒简化为正六边形，利用有限元法对多晶材料晶界处的应力情况进行了模拟分析，给出了晶界上正应力分布情况与晶界和外载荷方向夹角的关系。杨静等[6]通过有限单元法中的内聚力单元模拟多晶材料的晶界属性，对多晶材料的疲劳损伤破坏进行了模拟，预测了裂纹的萌生扩展情况。杨丹[7]用分子动力学方法探究了不同晶粒尺寸、孪晶含量等因素对纳米多晶Al的力学性能与变形机制的影响。卢敏[8]分别通过分子动力学方法、有限元方法对NiTi合金的沿晶和穿晶断裂行为进行了模拟分析，结果表明，晶粒尺寸越小，材料的断裂韧性越好。

在目前较为成熟的研究手段中，通过具体实验能够最客观地获取晶粒结构与各项材料参数之间的关系，但成本较高。用有限元方法模拟时，能够较好地得到材料内部的应力分布情况，但在模拟裂纹扩展时，要在裂纹尖端设置奇异单元、加密网格，求解应力强度因子，确定裂纹扩展方向，并同时要求裂纹顶点和裂纹边与有限元的节点单元边重合，相应的网格必须重新划分，使计算的效率和精度同时降低。同时，多晶材料随机多变的晶粒组织构型也对网格划分造成了一定的难度。利用分子动力学方法对多晶材料进行分析时，可以得出微观尺度下原子的行为和受力、能量物理参数，材料的晶体结构，并能够给出外界载荷作用下微裂纹的演化过程，但分析过程需要消耗大量的计算资源。

Silling[9]于2000年提出了近场动力学(peridynamic)方法。该方法是一种基于非局部思想的无网格方法，通过积分方程求解代替微分方程求解，对解决损伤、断裂等不连续、多尺度问题时具有明显的优势，目前已广泛运用于各类材料的破坏分析中[10-11]。Askari等[12]通过近场动力学微弹性本构模型对多晶硅破坏过程进行了模拟，其结果相对传统数值方法更接近实验现象。伊航[13]基于近场动力学理论对功能梯度材料的破坏过程进行了模拟，分析了材料梯度对试件起裂与完全破坏之间的关系。姜晓伟等[14]基于常规态型近场动力学，引入基于能量的失效判断准则，对复合材料的分层扩展过程进行模拟，其结果能够较好地捕捉到“指甲盖形”的分层前沿，与试验结果达到了较好的一致性。孙杰等[15]通过在传统键型近场动力学本构模型中引入非局部长程作用力强度尺寸效应的核函数修正项，部分消除了传统键型近场动力学模型对泊松比的限制，并对多裂纹脆性板的破坏过程进行模拟分析，给出了初始裂纹数量、位置和方向对结构破坏形式、扩展路径和承载能力的影响规律。

因此，现针对多晶金属材料内部微结构的多样性和损伤破坏过程的复杂性，基于Voronoi图基本思想，在Fortran环境下对多晶材料的晶粒组织结构进行模拟，建立晶粒大小、均匀度可控的晶粒组织模型。结合近场动力学方法在处理不连续、多尺度等问题的优势，对外界载荷作用下多晶材料内部裂纹扩展及结构破坏过程进行模拟，对晶界强度、晶粒大小等因素影响下多晶材料的断裂性能进行研究。

1 多晶材料破坏分析的近场动力学模型

1.1 近场动力学理论

近场动力学理论的基本思想是在某一空间域R内，将实体模型离散成一系列包含着周围实体模型的体积、密度等信息的物质点。这些物质点分别与近场域Hx(近场域半径δ)内的其他物质点之间存在相互影响，这种影响被称为“键”[11]，键的强度通过力密度f进行描述，如图1所示。

u′、u为相邻物质点位移矢量；x′、x为相邻物质点坐标位置矢量；t为当前时刻图1 物质点间矢量关系Fig.1 Vector relations between points of matter

f=f(u′-u,x′-x,t)

(1)

根据牛顿第二定律可得其运动方程

(2)

在近场动力学方法中，为了描述材料的破坏，需要定义“键”的伸长率s和极限伸长率s0。在式(2)中令x′-x=ξ，u-u′=η，其中ξ表示相邻物质点的相对位置矢量，η表示相邻物质点的相对位移矢量，则伸长率s可表示为

(3)

“键”的极限伸长率s0，其大小与该材料杨氏模量E、断裂时的临界能量释放率Gc以及近场域半径δ的选取有关，具体表达式[16]为

(4)

式(4)中：Gc与材料的断裂韧性KIC存在以下关系：

(5)

在近场动力学分析中，通过引入一个与时间t相关的标量函数μ(ξ,t)来描述键的破坏，其形式为

(6)

当物质点之间键的伸长率s超过极限伸长率s0时，则定义“键”断裂，此后该物质点对之间不再产生相互作用，此时μ(ξ,t)的值为0，反之为1。

在材料破坏分析中，通过求解某一物质点与其相邻物质点间的破坏键的数量与初始键总数的比值来定义物质点x在t时刻的损伤度，其函数形式为

(7)

1.2 多晶材料晶界影响

多晶材料由很多个单晶颗粒与包围单晶颗粒的晶界组成，由于晶界偏聚等现象致使材料在晶粒内部和晶界附近有不同的强度系数，材料在外界载荷作用下将发生晶间或穿晶破坏。为了考虑晶界对多晶材料断裂性能的影响，通过引入晶界强度系数β[12]，并建立关系式为

(8)

对于同一晶粒内的物质点之间的“键”的破坏情况仍用式(6)判别；若两物质点属于不同晶粒，即结合式(8)，采用以下布尔函数[式(9)]判别。

(9)

相应的物质点的损伤判定则定义为

(10)

1.3 基于Voronoi图的多晶材料几何模型

在基于近场动力学方法研究多晶材料的断裂性能时，模型结构中必须包含材料的微观结构特点。Voronoi图能够基于一组可控的特定点将对平面进行分隔，能很好地反映晶体的微观结构，目前被广泛用于复杂问题的处理中[4,18]。用Voronoi方法来建立多晶体微观结构图，在保证晶粒随机分布的同时，使得晶格内的每一点到晶粒中心点位置距离最短，在近场动力学分析时，把同一晶格内的所有物质点的属性设定与中心点相同。

尽管Voronoi图可以在Arcgis、Sketchup等软件内自动生成等，但为了建立基于Voronoi图的多晶材料破坏分析的近场动力学模型，需要在Fortran环境下编写生成Voronoi图的程序结构，具体过程如下。

(1)通过随机抽样并结合几何模型尺寸，确定的一系列随机点坐标位置(xi,yi)作为构建泰森多边形的特征点，同时赋予这些特征点一个不同的颜色作为区分。根据Voronoi图定义，这些点代表着晶粒的中心位置。

(2)在几何模型尺寸范围内以均匀分布的形式生成近场动力学物质点。

(3)计算每个物质点与所有特征点之间的距离，筛选出与当前物质点距离最短的特征点，并赋予相同的颜色。

(4)对求解域内的所有物质点重复步骤(3)，可将多晶材料的实体模型转化为可用近场动力学方法分析的Voronoi多晶材料模型。

通过控制特征点的个数及最小间隔，就可以控制晶体材料中晶粒大小及均匀度，如图2为设置不同特征点个数时所对应晶粒大小变化情况。其中不同颜色用于区别每个独立的晶粒组织，其目的主要表征在颜色不同的色块间存在晶界。

图2 不同晶粒个数的晶体结构图Fig.2 Crystal structure diagram of different number of grains

2 数值计算

为了研究多晶金属材料的断裂性能，以7xxx系高强度铝合金为研究对象，材料参数如表1所示。

表1 材料参数Table 1 Material parameters

其中计算模型结构为长×宽=1 mm×1 mm，厚度0.002 mm的正方形板；在近场动力学分析中,取物质点间距Δx=0.002 mm, 近场域半径δ=3.015Δx, 物质点的数量是250 000，时间步长Δt=0.25 ns，在模型中心预制长度为0.2 mm中心穿透裂纹并在竖直方向施加40 m/s速度载荷后，分别对不同晶界强度系数和不同晶粒大小所对应的模型的断裂过程进行分析。

2.1 不同晶界强度系数时的断裂特性

为了研究晶界强度对多晶材料断裂性能的影响，在计算模型内随机生成大小形状不同的600个晶粒后(其中平均晶粒直径d=46 μm)，对晶界强度系数β分别取1、2、4和0.5时多晶材料的断裂过程进行模拟研究，得到了如图3～图6所示的损伤云图。

从图3～图6所示的损伤云图中，通过观察物质点的损伤度变化情况，可较为清晰的反映裂纹扩展的动态过程，当β=1时，晶界强度与晶内强度一致，可认为模型是均质材料，晶界的存在不产生影响。此时随着加载，模型中裂纹呈对称形式沿水平方向两端扩展，直到2 950时间步完全贯通模型；当β=2、β=4时，随着加载，裂纹的扩展会受到晶界的阻碍，在晶界处会产生小幅偏转及分叉，由于晶界分布的随机性导致裂纹不再呈现对称性，表现为穿晶破坏。与β=1相比，如果取相同加载时间，裂纹在β=2、β=4时扩展尺寸较短，结构完全破坏所需时间较长；当β=0.5时，材料在晶界处强度较低，随着加载，裂纹会沿晶界扩展，表现为晶间破坏，与β≥1情形相比，裂纹扩展速度最快，材料完全破坏用时短。

图3 β=1时裂纹扩展过程Fig.3 Damage cloud of the model when β=1

图4 β=2时裂纹扩展过程Fig.4 Damage cloud of the model when β=2

图5 β =4时裂纹扩展过程Fig.5 Damage cloud of the model when β=4

图6 β =0.5时裂纹扩展过程Fig.6 Damage cloud of the model when β=0.5

在近场动力学方法中，通过判定某物质点与其相邻物质点间“键”的破坏来模拟裂纹扩展。同时，也可以用每个物质点的损伤度进行累计来反映模型的破坏情况。如图7所示，在晶界强度系数取不同值时，加载初期，材料损伤度均为水平直线，表示材料还未出现新的损伤。随着加载的继续，损伤度开始增加，对应裂纹开始扩展，扩展速率较慢，当继续加载，裂纹进入快速扩展期。晶界强度系数不同，裂纹起裂时间不同。当β<1，裂纹开始扩展时间较早，β≥1时裂纹开始扩展时间较晚。

图7 不同晶界强度系数累计损伤度随时间变化Fig.7 Cumulative damage curves of different grain boundary strength coefficients over time

2.2 不同晶粒大小时的断裂特性

为了研究晶粒大小对多晶材料断裂性能的影响，在计算模型中，分别对晶粒平均直径取65、46、38、33 μm 4种情况的断裂性能进行了研究。由于实验中7xxx系铝合金多以准解理断裂为主[19]，即脆性穿晶断裂，本文选取其晶界强度系数β=4，对比分析了4种晶粒大小对应的损伤特性。

如图8所示，随着晶粒平均直径的减小，晶界成分在模型中占比提高，对裂纹扩展的阻碍作用也在增强,在相同速度载荷下，裂纹完全贯穿模型所需时间增加,即材料的抗断裂性能越好。

图8 不同晶粒尺寸累计损伤度随时间变化Fig.8 Cumulative damage curves of different grain sizes over time

在相同速度载荷作用下，通过求解不同晶粒大小所对应的材料完全破坏时间步，建立完全破坏时间步与晶粒平均直径d倒数的平方的关系曲线如图9所示，结果发现晶体材料的抗断裂性能与晶粒平均直径倒数的平方呈线性关系，且晶粒细化有助于提高多晶材料的力学性能。

图9 晶粒尺寸对模型完全破坏时间的影响Fig.9 Effect of grain size on complete failure time of the model

3 结论

以近场动力学理论为基础，引入Voronoi图方法建立多晶材料的晶体结构图，以Voronoi图的特征点的个数来控制晶体的大小，通过晶界强度系数来区别晶粒内部和晶界处的强度关系，建立了描述晶体结构破坏的近场动力学模型。对7xxx系铝合金的拉伸破坏过程进行模拟，得到以下结论。

(1)通过控制晶界强度系数β，该模型能够模拟多晶材料穿晶断裂或晶间断裂两种不同的断裂模式，同时模拟结果表明材料在穿晶断裂时的抗断裂性能要优于晶间断裂的抗断裂性能。

(2)通过对不同晶粒大小模型的裂纹扩展模拟发现，在穿晶断裂模式下，完全破坏时间步与晶粒平均直径d倒数的平方呈正相关线性关系，即晶粒细化有助于提高材料的抗断裂性能。