【研究意义】材料的变形与破坏主要源自材料的微结构，如空位、位错、晶界与微裂纹等^[1]。要想深入研究材料的变形与破坏的特性与机理，还需要开展多尺度分析，把宏观分析与微纳观分析结合起来，在更深层次上找到问题的根源和变形破坏机制^[2]。加载的材料通常会表现出可逆行为，当加载除去后，变形就会消失，并且储存能在变形原子键的弹性能得以释放。然而，对于塑性变形和裂纹开裂情况，这一行为是不可逆的。材料在变形中会有能量损耗，其主要来自两个方面:一是断裂产生新表面时所损耗的，二是位错形核和发射位错时所损耗的。因此，对于存在裂纹缺陷的材料变形行为的检验，特别是当位错的形核和运动是由非弹性变形主导的情况，是非常重要的。【前人研究进展】目前已有许多断裂模拟研究，如高英俊等^[3]研究韧性材料的微裂纹扩展和连通，毛鸿等^[4]研究材料裂纹分叉的机理，刘晓波等^[5]对铝裂纹扩展行为进行分子动力学模拟，Song等^[6]用相场法研究铁电体的裂纹尖端，Abdolahi等^[7]则用相场法模拟铁电体材料的断裂，张跃等^[8]进行脆断微裂纹形核的原位观察。在这些模拟研究方法中，既有单尺度分子动力学模拟，也有复杂的多尺度模拟技术。分子动力学模拟有一定的局限性，即它的时间尺度主要适用于原子振动的时间尺度(10^-14~10^-10 s)，难以拓展到原子扩散的时间尺度(10^-6 s)，而且由于分子动力学模拟强烈依赖于原子势函数的选择，这就要求所加载的应变力速率远高于实际工业生产的应变力速率(10² s^-1)，达到10⁷~10⁹量级^[9]。这使其对微裂纹的扩展模拟与实际存在相当距离。【本研究切入点】晶体相场方法(PFC) ^[10-14]是基于密度泛函理论、继承传统相场模型优势而建立的，它能够揭示晶体学结构特性以及空间尺度为原子尺度、时间尺度为扩散时间尺度下的结构演化。由于该方法能够很好地描述扩散时间尺度上的微结构演化行为^[15-17]，不针对特定材料，因此更适合用于研究微纳米尺度上的裂纹扩展细节^[18-19]。【拟解决的关键问题】应用PFC方法研究不同初始晶向倾角条件下，单轴拉伸应变作用的裂纹生长特征和扩展规律，揭示纳米级裂纹扩展机理及其对材料断裂的影响。

基于PFC模型的液态金属研究如枝晶生长已有许多报道，但更多的是用于固态金属研究。对于固态金属材料，其原子的位置呈规则周期性排列，通过引入周期性相场变量，其局域位置的最大值对应于原子的位置；对于均匀相(液相等)中的原子分布为均匀分布，其值为常量。用周期原子密度函数作为相场变量，要符合上述两方面的要求，其表达式^[11-14]可以写成

式中，等号右边第1项反映的是晶格原子密度的周期排列特征，第2项反映的是均匀相(如液相)的原子密度均匀分布。此时，系统无量纲的自由能函数可以写成^[11]

式中，ρ为局域原子密度，γ为与温度有关的唯象参数，∇²为拉普拉斯算子。

在单模近似下，可以求得公式(2) 的一个稳定特解为

式中， $A=\frac{4}{5}{{\rho }_{0}}+\frac{4}{15}\sqrt{-15r-36{{\rho }_{0}}^{2}}$ ，反映晶相原子密度周期结构的振幅；q=$\sqrt{3}$/2，是在能量函数取极小值时所得；ρ₀为平均密度。

采用保守场Cahn-Hilliard动力学方程^{[15, 18]}描述原子密度随时间的演化，该方程具体如下：

其中t为时间变量。对无量纲动力学演化方程(4) 采用半隐式傅里叶伪谱方法^[20-21]求解，其离散形式为

式中${\hat{\rho }}$_k(t+Δt)为傅里叶空间下t+Δt时刻的原子密度， ${{{\hat{\rho }}}_{{\vec{k}}}}\left( t \right)=\int{{{e}^{-i\vec{k}\vec{\gamma }~}}}{{{\hat{\rho }}}_{{\vec{k}}}}\left( \vec{\gamma },t \right)d\vec{\gamma },\hat{\rho }_{k}^{3}\left( t \right)=\int{{{e}^{-i\vec{k}\vec{\gamma }~}}}\hat{\rho }_{k}^{3}\left( \vec{\gamma },t \right)d\vec{\gamma },\vec{k}$ ，为傅里叶空间的波矢，满足 ${{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\rightharpoonup}$}} {k}}}^{2}}={{\left| k \right|}^{2}}$ 。ρ_i,j^t和ρ_i,j^t+Δt分别为t和t+Δt时刻的原子密度。整理方程(5) 后得：

本研究晶体相用二维三角点阵相表示，其中平衡相为条状相、三角相和液相(图 1)。选取原子密度参量为ρ₀=0.49，温度参量γ=-1.0(图 1中B点)，缺口处的参数设置为ρ₁=0.79，γ=-1.0(图 1中A点)。应用公式(3) 设置单晶体结构。计算模拟区域为1 024Δx×512Δy,Δx=Δy=π/3，其原子排列方向与y轴夹角θ为0°、5°、10°。在样品的中心位置，设置一半径r=8的圆形缺口作为初始裂口。由于不涉及材料的物性参数，模拟所用的参数均已无量纲化处理，并将连续空间离散为正方格子，计算时采用周期性边界条件。设置的3组样品参数如表 1所示。

图 1 (Fig. 1)

图 1 单模近似得出的二维相图(阴影部分代表两相共存区)[7] Fig.1 The two-dimensional phase diagram obtained by the single mode approximation(The shadow parts represent the two phase coexistence region)[7]

表 1(Table 1)

表 1 3组不同晶向倾角的裂口样品参数 Table 1 3 groups of the split smaple parameters of different crystal orientations Sample γ ρ0 ρ1 θ ${\dot{\varepsilon }}$ A -1.0 0.49 0.79 0° 6×10-6 B -1.0 0.49 0.79 5° 6×10-6 C -1.0 0.49 0.79 10° 6×10-6
	表 1 3组不同晶向倾角的裂口样品参数 Table 1 3 groups of the split smaple parameters of different crystal orientations

为观察裂纹的萌生和扩展情况，对设置的含有裂口的样品施加拉应变作用。首先，经过5×10⁴步弛豫，得到函缺口的稳态样品(图 2，其放大区域给出实验的坐标体系)，然后再沿y轴方向(图 3所示)施加拉应变。在变形过程中，x方向空间步长保持不变，y方向空间步长随着应变速率在每一时间步长下都有增量d=${\dot{\varepsilon }}$nΔt，其中，${\dot{\varepsilon }}$为无量纲的应变速率，n为施加拉应变的时间步长数，Δt为时间步长，Δy′=Δy(1+n${\dot{\varepsilon }}$Δt)。

图 2 (Fig. 2)

图 2 初始样品 Fig.2 Initial sample

图 3 (Fig. 3)

图 3 y轴施加拉应变的示意图 Fig.3 Schematic diagram of tensile strain applied to y axis

当外加应变量达到临界应变值0.138时，裂口开始起裂(图 4a)，裂纹开裂方向分别向左上[121]方向和右上[211]方向，与y轴的夹角为60°，裂纹长度约为4～5个原子距离，两端裂尖的裂口附近都出现位错结构(图 4a中的放大区域所示)。当外加应变量增加到0.162时，裂纹已经由与y轴夹角60°的向左上和右上扩展，转变为向左边水平方向和右边水平方向扩展(图 4b)，在裂尖的放大细节图上可以观察到裂纹尖端仍然出现位错，但左边位错方向与图 4a相比发生转向，顺时针旋转60°，位错的柏氏矢量沿[211]方向。外加应变量继续加大，裂纹仍然沿左边水平方向和右边水平方向扩展，裂纹长度增加，左边位错矢量方向不变(图 4c)。到扩展后期，外加应变量达到0.210，裂纹总体仍然沿左边水平方向和右边水平方向扩展(图 4d)，裂纹长度约为模拟区域的三分之二，但尖端可以看到较明显的锯齿状，右边裂尖处的位错与图 4a的相比逆时针转动60°。表明裂纹尖端位错滑移方向改变，沿[121]和[211]方向交替变化，形成锯齿边缘结构(图 4d中的放大区域所示)。整个裂纹的形状呈解理状，裂尖前端没有出现塑性变形区，也没有空洞形成和位错发射现象，表明裂纹的扩展属于脆性解理扩展。

图 4 (Fig. 4)

左图为裂纹扩展，右图为对应的应力分布 Left is the crack propagation, right is the stress distribution of the corresponding 图 4 样品A的裂纹扩展演化 Fig.4 Crack propagation evolution of sample A

开始阶段，在y向的拉应变作用下，裂口发生应力集中；当外加拉应变量达到临界应变值0.138时，裂口向两边发射滑移位错，这是典型的局部塑性变形的表现。由图 5a可见，距裂口约2个原子距离的左右两边各出现一个位错，左右两边位错矢量方向相反，柏氏矢量几乎沿水平方向，且此时没有出现与裂口相连接的裂纹。当外加应变量达到0.162时(图 5b)，左右两边已经形成明显的裂纹并与裂口相连接，由图 5b的放大细节图可知，左边的裂尖有一个位错，其柏氏矢量方向与水平方向夹角约60°。此时位错滑移受阻，裂尖应力集中，裂纹前进转向，裂纹边缘呈锯齿状，即左边裂纹扩展方式为扩展-转向-扩展。总体来说，裂纹向左下延伸与水平方向呈约30°，即[211]方向；右端裂纹向右边水平方向延伸。当外加应变量达到0.186时(图 5c)，左边裂纹向左下方扩展，扩展方式改为解理扩展，方向为[211]；右边裂纹扩展方式仍然为扩展-钝化-扩展，裂纹形状呈锯齿状，总体上裂纹扩展方向水平向右；右边裂尖也有一个位错，其柏氏矢量方向为向上右偏60°，方向为[101]。到扩展后期(图 5d)，左边裂纹位错滑行阻力较小，形成解理裂纹，解理裂纹扩展方向为[121]，其裂尖位错方向与图 5c相同，依旧为向下水平右偏60°，裂纹右边仍然为锯齿状扩展，裂尖的刃位错的半原子面垂直向上。

图 5 (Fig. 5)

左图为裂纹扩展，右图为对应的应力分布 Left is the crack propagation, right is the stress distribution of the corresponding 图 5 样品B的裂纹扩展演化 Fig.5 Crack propagation evolution of sample B

由图 6a可见，由于裂口应力集中，当外加应变量达到临界应变值0.138时，裂口向左右两侧各发射一个位错，位错滑移方向分别为方向为向上右偏沿[211]方向和向下左偏沿[211]方向，位错在滑移过程中沿滑移线诱发生成一系列孤立的空位，空位发展成空洞，然后空洞长大连通，形成裂纹。在空位前端，左边的刃位错滑移方向沿[121]方向，右边的刃位错滑移方向沿[121]方向(图 6a，b)。当外加拉应变量继续增加时，靠近裂口的空位发展成空洞，空洞相互连接，形成一条与水平方向夹角30°的倾斜裂纹(图 6b)。由图 6b右下角放大图可见，在裂尖处出现应力集中，在裂尖生成新的位错。随着外加拉应变量的进一步增加，在裂口处向左下和右上延伸的过程中，裂尖应力集中，引起裂尖钝化并使裂纹扩展方向转向；远离裂口处的空位则相互连接成新的裂纹，方向与原主裂纹扩展方向相同。在拉应变作用下，远离主裂纹的空位末端处的位错沿新的方向滑移，出现裂纹分叉(图 6c，d)，形成二次裂纹扩展。主裂纹钝化转向的尖端也出现位错(图 6d)，由此可见主裂纹两端钝化转向后开始呈解理状生长，而远离裂口的空位形成的新裂纹也开始钝化转向，该新裂纹的裂尖处也出现位错，位错的方向为左侧的位错向上右偏60°，右侧的位错向下左偏60°。该二次裂纹与主裂纹的扩展方向夹角约为60°。整个区域形成几条不连通的裂纹，表明裂纹主要是韧性断裂模式扩展。

图 6 (Fig. 6)

左图为裂纹扩展，右图为对应的应力分布 Left is the crack propagation, right is the stress distribution of the corresponding 图 6 样品C的裂纹扩展演化 Fig.6 Crack propagation evolution of sample C

本研究采用晶体相场法，研究具有不同初始晶向倾角的样品在施加拉应变作用下的裂纹尖端扩展行为，得出如下主要结论：

1) 晶向角为0°、5°的样品，裂纹主要是解理脆性断裂模式扩展；晶向角为10°的样品，裂纹主要是韧性断裂模式扩展。随着晶向角的变化，裂纹从脆性向韧性裂纹扩展模式转化。

2) 裂纹作解理脆性扩展时，在裂尖处存在一对位错锁结构，由于位错的柏氏矢量垂直于裂纹扩展方向，这导致在裂尖处存在应力应变集中，原子键呈解理方式断开，裂纹沿光滑直线扩展；裂纹作韧性断裂扩展时，在裂尖前方不远处发生位错滑移运动，且在位错滑移轨道上产生空洞，空洞生长扩大并连通，形成裂纹，此时裂纹边缘呈锯齿形状。