水平集法模拟GH4706合金动态再结晶过程

图1 取样位置示意图及GH4706合金初始组织的EBSD分析

Fig.1 Schematic of sampling position (a) and inverse pole figure (IPF) (b) and grain size distribution (c) of initial microstructure of GH4706 superalloy (ND—normal direction, TD—transverse direction, CD—compression direction)

采用Gleeble-3500伺服液压机进行热变形实验。以10 ℃/s的加热速率将试样加热至目标变形温度，保温180 s以降低流动变形行为的各向异性。基于800 MN大型模锻压机的边界条件，设定变形温度为950~1150 ℃，温度间隔为50 ℃，应变速率分别为0.001、0.01、0.1和1 s^-1，将试样压缩至1.2真应变(本工作均采用真应变，以下简称应变)。在1050 ℃和0.1 s^-1条件下进行应变为0.4~0.7的热压缩实验，随后立即对试样进行水淬处理，同时自动采集真实应力-应变曲线数据。

采用配有电子背散射衍射仪(EBSD)的Sigma 500扫描电子显微镜(SEM)表征微观组织形貌及结构。采用线切割切取EBSD试样，使用400~3000号SiC砂纸依次对试样进行打磨。随后，将试样置于10%HClO₄ + 90%CH₃CH₂OH (体积分数)的混合溶液中进行电解抛光，温度在-30 ℃左右，抛光电压为22 V，抛光时间为20 s，扫描步长为1~3 μm。采用AZtecCrystal软件分析EBSD数据。

2 基于水平集法的动态再结晶模型构建

2.1 计算方法

首先从数学的角度给出水平集法与曲线演化理论之间的关系。在域 $Ω$ 上定义水平集函数(Ψ)，作为子域(即单个晶粒G)的晶界距离函数(Γ)。在每个节点上都计算出Ψ的值，晶粒G内部Ψ > 0，晶粒G外部Ψ < 0，晶粒G的晶界处Ψ = 0。描述水平集和晶界的公式可分别表示为^[25]：

\{\begin{matrix} Ψ_{i} (x, y, t) = \pm d (x, y, Γ) (x \in Ω) \\ Γ (t) = \{x, y \in Ω, Ψ_{i} (x, y, t) = 0\} \end{matrix}

(1)

式中，Ψ_i (x, y, t)表示第i个随时间(t)变化的水平集函数，d为Euclidean距离。图2为通过水平集函数获得的初始晶粒和晶界示意图，其中包括3个水平集距离函数Ψ₁、Ψ₂和Ψ₃。采用这些函数依据前述符号规则共同构建了3个晶粒及其晶界的初始形态。

图2

图2 采用水平集函数定义晶界

Fig.2 Definition of grain boundaries using level-set functions (Ψ represents level-set function; Ψ₁, Ψ₂, and Ψ₃ represent level set functions for Grain 1, Grain 2, and Grain 3, respectively)

采用晶粒的运动方程描述变形过程中材料微观组织的演变过程。假设第i个水平集平面封闭曲线(C_i )随时间变化，则曲线C_i (t)的表达式为^[24]：

C_{i} (t) = \{(x, y) | Ψ_{i} (x, y, t) = a\}

(2)

式中，a为t时刻的水平集函数值，用来判别水平集函数数值与0的关系。对式(2)求全导数，根据复合函数链式规则可得：

\frac{d Ψ_{i}}{d t} = \frac{\partial Ψ_{i}}{\partial t} + \vec{\nabla} Ψ_{i} \cdot \frac{\partial (x, y)}{\partial t} = 0

(3)

\frac{\partial (x, y)}{\partial t} = (\frac{\partial x}{\partial t}, \frac{\partial y}{\partial t}) = V

(4)

式中，V表示晶粒的长大速率。因此，式(3)可表示为：

\frac{\partial Ψ_{i}}{\partial t} = - \vec{\nabla} Ψ_{i} \cdot V

(5)

晶界动力学假设晶界运动受2种驱动力作用：毛细效应，即通过减少系统整体边界表面积从而减少晶界能量；存储能量梯度，即高内部存储能晶粒被其他晶粒消耗，从而实现系统整体能量的减少。V可表示为毛细力驱动项(V_c)和存储能量驱动项(V_e)之和^[19,24]：

V = V_{c} + V_{e}

(6)

V_{c} = - M γ ∆ Ψ_{i} \vec{\nabla} Ψ_{i}

(7)

V_{e} = M δ (\dot{ε}) ⟦E⟧

(8)

式中，M为晶界迁移速率；γ为晶界能； $\dot{ε}$ 为应变速率； $δ (\dot{ε})$ 为由储存能量引起的晶界迁移速率影响因子，该参数为不同应变速率下定义的无量纲常数，并通过软件中的多线性插值方法确定； $⟦E⟧$ 为晶粒在长大过程中晶界界面上的储能跳跃，可表示为^[19]：

⟦E⟧ = τ (ρ_{m} - ρ_{i})

(9)

式中，ρ_m为变形晶粒的位错密度，ρ_i 为第i个等轴晶粒的位错密度；τ为与材料相关的位错线能量，可表示为^[20]：

τ = α μ b^{2}

(10)

式中，α为Taylor因子，通常取值为0.5；µ为材料剪切模量，其值取决于温度；b为Burgers矢量模。

将式(6)~(10)代入式(5)，可得到N个水平集函数中任意一个函数的运动方程，即第i个晶粒的运动方程^[23~25]：

\{\begin{matrix} \begin{array}{l} \frac{\partial Ψ_{i}}{\partial t} - M γ ∆ Ψ_{i} + M δ (\dot{ε}) α μ b^{2} (ρ_{m} - ρ_{i}) \vec{\nabla} Ψ_{i} = 0 \\ (\forall i \in \{1, 2, . . ., N\}) \end{array} \\ Ψ_{i} (x, y, Δ t) |_{Δ t = 0} = Ψ_{i}^{0} (x, y) \end{matrix}

(11)

式中， $Ψ_{i}^{0}$ 表示第i个水平集函数在单位时间步长内的初始化值。

M和γ可分别表示为^[20]：

M = M_{0} e x p (\frac{Q_{m}}{R T})

(12)

γ = \{\begin{matrix} \frac{μ b θ_{m}}{4 π (1 - v)} (θ_{i} \geq 15 °) \\ γ_{m} \frac{θ_{i}}{θ_{m}} (1 - l n (\frac{θ_{i}}{θ_{m}})) (θ_{i} < 15 °) \end{matrix}

(13)

式中，M₀为指数前因子，Q_m为晶界迁移的激活能，R为气体常数，T为热力学温度；υ为泊松比，θ_i 为第i个晶粒与其相邻晶粒之间的取向差，γ_m和θ_m分别为晶界转变为大角度晶界时的边界能和取向差。

使用水平集法模拟金属动态再结晶的过程中还需耦合再结晶形核模型，之后才能对晶粒的形核与长大行为进行仿真模拟。研究^[26~29]表明，变形温度与应变速率对形核速率具有重要影响。本工作使用线性形核模型表征形核速率( $\dot{n}$ )^[26]：

\dot{n} = K_{g} (T, \dot{ε}) Φ ∆ t \{\begin{matrix} Φ = 0 (Φ < ρ_{c r}) \\ Φ = 1 (Φ \geq ρ_{c r}) \end{matrix}

(14)

式中， $K_{g} (T, \dot{ε})$ 为单位时间内达到形核条件的任意晶界表面的形核体积，取决于温度和有效塑性应变速率；Δt为单位时间步长； $K_{g} (T, \dot{ε})$ 为与温度和应变速率相关的任意晶界的形核体积；Φ为形核条件函数；ρ_cr为形核的临界位错密度。形核速率和再结晶分数之间存在Johnson-Mehl-Avrami-Kolmogorov (JMAK)动力学关系^[26]：

F = 1 - e x p (- \dot{n} t^{c})

(15)

式中，F为DRX晶粒体积分数，c为时间指数，取值为1^[26]。

ρ_cr可以通过以下公式计算^[30]：

ρ_{c r} = [\frac{- γ \dot{ε} \frac{S}{M δ (\dot{ε}) τ^{2}}}{l n (1 - \frac{S}{H} ρ_{c r})}]

(16)

式中，H和S分别为硬化参数和软化参数，符合L-J位错密度模型，描述了在热变形过程中临界位错密度的演化^[31]。H和S可表示为：

\frac{d ρ}{d ε} = H - S ρ

(17)

H = H_{0} {(\frac{\dot{ε}}{{\dot{ε}}_{0}})}^{- m} e x p (\frac{- m Q}{R T})

(18)

S = S_{0} {(\frac{\dot{ε}}{{\dot{ε}}_{0}})}^{- m} e x p (\frac{- m Q}{R T})

(19)

式中，ε为应变，ρ为位错密度； ${\dot{ε}}_{0}$ 为应变速率校准常数，在有限元软件中取1；H₀和S₀分别为与温度和应变速率相关的初始硬化参数和初始软化参数；Q为基于修正的JMAK型方程的激活能；m为应变速率敏感性系数，反映了材料随应变速率增加的硬化趋势。m的计算公式为^[32]：

m = {\frac{d l g σ}{d l g \dot{ε}}|}_{T}

(20)

式中，σ为流变应力。

形核的出现取决于位错密度，当一个新核区域的平均位错密度达到式(16)中的ρ_cr时，新核出现并形成新的平均位错密度(ρ₀)。然后，根据DDRX形核机制，位错倾向于集中在晶界区域，临界半径为r的新核会出现在晶界处。因此，新核的中心被限制在晶界周围±L的范围内(L = r)，如图3所示。

图3

图3 晶粒形核过程示意图

Fig.3 Schematic of grain nucleation process (ρ—dislocation density, ρ_cr—critical dislocation density, L—critical nucleation distance)

在此基础上，采用DIGIMU®软件建立了描述全晶粒拓扑结构的耦合模型，该模型的构建基于水平集距离函数和经典位错密度模型。基于变形过程的计算框架，实现了初始多晶生成、晶界描述、多晶变形方法、加工硬化和动态回复模型、位错密度演化模型和形核模型以及网格和时间步长的自动调整，如图4所示。在求解的任意一步，均可得到瞬时晶粒形貌和瞬时F。最后，通过多步迭代运算，实现热变形过程中DRX行为与晶粒形貌动态演化的多场多尺度耦合模拟。DIGIMU软件的初始设置如下：(1) 模拟条件为温度950~1150 ℃、应变速率0.001~1 s^-1和应变1.2；(2) 根据图1，初始平均晶粒尺寸设置为177.67 μm，手动输入(Histogram选项)初始晶粒尺寸分布数据以建立初始的微观组织分布，如图5所示；(3) 使用文本编辑器软件输入多场、多尺度有限元模型的实验数据及参数值，建立初始GH4706合金材料数据库(GH4706.mtc)，在模拟设置中选择并加载；(4) 选择默认的硬化回复定律、阈值静态回复定律和临界位错密度形核起始模型。

图4

图4 水平集法耦合位错密度模型模拟动态再结晶(DRX)过程的流程图

Fig.4 Flowchart of dynamic recrystallization (DRX) simulation using level-set method coupled with dislocation density model (L-J represents Laarsaoui-Jonas, LS represents level-set; Ψ_i —the ith level-set function, Δt—time step, M—migration rate of grain boundaries, γ—grain boundary energy, ΔΨ_i —the ith level-set function per unit, α—Taylor factor, µ—material shear modulus, b—Burgers vector, ρ_m—dislocation density of deformed grains, ρ_i —dislocation density of the ith equiaxed grain, $Ψ_{i}^{0}$ —interfaces of the ith level set function per unit; M₀—pre exponential factor, Q_m—activation energy for grain boundary migration, R—gas constant, T—thermodynamic temperature, V—growth rate of grains, $δ (\dot{ε})$ —factors affecting grain boundary migration rate, $⟦E⟧ —$ jump for energy storage of grains; $\dot{ε}$ —strain rate, τ—dislocation line energy, H—hardening parameter, S—softening parameter; ε—strain, H₀—initial hardening parameter, S₀—initial softening parameter, ${\dot{ε}}_{0}$ —calibration constant for strain rate, Q—activation energy, m—strain rate sensitivity coefficient)

图5

图5 DIGIMU软件中的初始组织模型

Fig.5 Initial microstructure model in DIGIMU software

2.2 位错模型参数识别

基于真实应力-应变曲线，标定了材料数据库GH4706.mtc中的位错模型参数。对式(17)进行积分，从而计算H₀和S₀。ρ可表示为：

ρ = \frac{H}{S} - \frac{C}{S} e x p (- S ε)

(21)

式中，C为常数。考虑初始条件为ε= 0和ρ=ρ₀，C可表示为：

C = H - ρ_{0} S

(22)

在应变足够大时，使用Taylor公式表示应力与位错密度之间的关系^[33]：

σ = α μ b \sqrt[]{ρ}

(23)

将式(21)和(22)代入式(23)得到流变应力的表达式：

σ = {[σ_{s a t}^{2} - (σ_{s a t}^{2} - σ_{0}^{2}) e x p (- S ε)]}^{\frac{1}{2}}

(24)

式中，σ_sat为饱和应力，σ₀为初始应力。

对式(24)求导得到：

2 σ θ = S (σ_{s a t}^{2} - σ_{0}^{2})

(25)

式中，θ为硬化率。

根据式(23)和(24)，H和S的关系可以表示为：

H = \frac{S σ_{s a t}^{2}}{{(α μ b)}^{2}}

(26)

综上所述，可通过绘制2σθ-σ²曲线计算出S。图6a为根据真实应力-应变曲线绘制的硬化率(θ= dσ / dε)曲线。可以看出，硬化率随着应力增加而增加，随后迅速减小，直到达到临界再结晶应力。根据图6a进一步绘制2σθ-σ²曲线，如图6b所示。通过图6a和b中曲线的斜率可分别确定σ_sat和S，然后根据式(26)求得H，根据式(18)~(20)可分别求得H₀、S₀和m，如表1~5所示。从结果来看，m对H₀和S₀具有关键性影响。

图6

图6 硬化率(θ)-流变应力(σ)以及2σθ-σ²关系曲线

Fig.6 Relationship curves of σ-θ (a) and 2σθ-σ² (b) (θ—hardening rate, σ—flow stress, σ_sat—saturated stress)

表1 不同温度和应变速率下的饱和应力(σ_sat) (MPa)

Table 1 σ_sat under different temperatures and strain rates

Strain rate / s^-1	950 ^oC	1000 ^oC	1050 ^oC	1100 ^oC	1150 ^oC
0.001	139.5	89.7	59.4	50.1	45.7
0.01	197.7	145.3	107.5	83.4	63.5
0.1	303.5	217.8	182.1	139.7	109.4
1	395.7	304.6	272.4	224.7	194.1

表2 不同温度和应变速率下的软化参数(S)

Table 2 S under different temperatures and strain rates

Strain rate / s^-1	950 ^oC	1000 ^oC	1050 ^oC	1100 ^oC	1150 ^oC
0.001	23.91	24.52	20.55	24.42	22.54
0.01	23.44	25.41	21.94	23.72	24.53
0.1	19.25	21.91	20.63	21.64	20.43
1	18.32	20.83	19.21	19.53	18.94

表3 不同温度和应变速率下的初始软化参数(S₀)

Table 3 S₀ under different temperatures and strain rates

Strain rate / s^-1	950 ^oC	1000 ^oC	1050 ^oC	1100 ^oC	1150 ^oC
0.001	188256	294331	1620368	2574613	3784939
0.01	283204	490727	3164695	1320442	4827380
0.1	356900	680766	1443620	2161982	3855893
1	583927	754935	2195793	2957684	3757594

表4 不同温度和应变速率下的硬化参数(H)

Table 4 H under different temperatures and strain rates

Strain rate / s^-1	950 ^oC	1000 ^oC	1050 ^oC	1100 ^oC	1150 ^oC
0.001	99.5085	50.8095	17.9823	15.8025	9.64638
0.01	237.0564	135.6558	62.0256	41.4393	18.4382
0.1	459.1146	262.7919	170.1042	105.4578	56.3837
1	553.8599	389.5857	295.9049	184.8933	103.8335

表5 不同温度和应变速率下的应变速率敏感性系数(m)

Table 5 m under different temperatures and strain rates

Strain rate / s^-1	950 ^oC	1000 ^oC	1050 ^oC	1100 ^oC	1150 ^oC
0.001	0.184	0.389	0.446	0.256	0.237
0.01	0.194	0.185	0.269	0.282	0.259
0.1	0.161	0.144	0.217	0.232	0.218
1	0.083	0.267	0.288	0.107	0.113

2.3 使用反向传播人工神经网络和Pareto前沿法识别 $K_{g} (T, \dot{ε})$ 和 $δ (\dot{ε})$ 参数

Pareto原理主要包括以下3个方面^[34]：(1) 非劣支配关系，如果1个解至少在1个维度上优于另1个解，且在其他维度上不劣于后者，则称前者在多目标情况下非劣支配后者；(2) Pareto最优解，即在给定资源条件下，无法仅优化某一个目标而不损害其他目标的解，也就是说，Pareto最优解是在多目标优化问题中无法进一步改进的解；(3) Pareto前沿，即所有Pareto最优解的集合。

与其他可根据数据识别的参数不同， $K_{g} (T, \dot{ε})$ 和 $δ (\dot{ε})$ 参数无法通过拟合获得。2者需要根据EBSD数据与模拟结果的不断比较从而进行逆向调整，得到最优的参数结果。本工作使用Pareto优化方法实现此过程，首先以F和平均晶粒尺寸(D)为优化目标，然后将不同应变下模拟结果与实验结果偏差百分比的平均值作为评价指标，分别建立D的误差函数f₁(K_g, δ)和F的误差函数f₂(K_g, δ)：

f_{1} (K_{g}, δ) = \frac{1}{n} \sum_{i = 1}^{n} \frac{D_{p}^{i} - D_{e}^{i}}{D_{e}^{i}}

(27)

f_{2} (K_{g}, δ) = \frac{1}{n} \sum_{i = 1}^{n} \frac{F_{p}^{i} - F_{e}^{i}}{F_{e}^{i}}

(28)

式中， $δ$ 为晶界迁移速率影响因子， $K_{g}$ 为形核体积，n为实验次数； $D_{p}^{i}$ 为平均晶粒尺寸实验值， $D_{e}^{i}$ 为平均晶粒尺寸模拟值； $F_{p}^{i}$ 为DRX分数实验值， $F_{e}^{i}$ 为DRX分数模拟值。

分别对形核体积和晶界迁移速率影响因子在0~0.1 mm³和0~10的范围内进行取值，形成40个参数组进行模拟，将输出结果与实验值按照式(27)和(28)进行误差计算，然后进行反向传播人工神经网络(back propagation-artificial neural network，BP-ANN)训练。利用MATLAB软件生成函数 $f_{1} (K_{g}, δ)$ 和 $f_{2} (K_{g}, δ)$ 的可执行代码，最后使用遗传算法中的多目标Pareto前沿优化方法进行优化。图7为基于反向传播人工神经网络训练和Pareto前沿方法的优化流程。

图7

图7 基于反向传播人工神经网络(BP-ANN)和Pareto前沿方法的仿真参数识别流程

Fig.7 Simulation parameter identification process based on back propogation artificial neural network (BP-ANN) and Paretn-optimization method (f—objective function, f₁—average grain size objective function, f₂—DRX fraction objective function)

使用MATLAB软件中的遗传算法工具箱，设置Pareto分数0.4、种群数量100、进化代数500、适应度函数偏差1 × 10^-10，随后绘制Pareto前沿图，如图8所示，优化后的Pareto前沿结果标注在虚框中。从虚框中的Pareto解集[x, fval]中选择出最佳参数组合，最终求解出K_g和δ (表6)。使用文本编辑器软件更新K_g和δ，建立参数校准的GH4706合金材料数据库。

图8

图8 多目标优化结果的Pareto前沿图

Fig.8 Pareto sketch map of multi-object optimization results (δ—factors affecting grain boundary migration rate, K_g— nucleation volume)

表6 不同参数条件下的晶界迁移速率影响因子(δ)和形核体积(K_g)

Table 6 δ and K_g under different parameters

Strain rate / s^-1	K_g / mm³					δ
Strain rate / s^-1	950 ^oC	1000 ^oC	1050 ^oC	1100 ^oC	1150 ^oC	δ
0.001	5.2 × 10^-6	3.0 × 10^-6	2.0 × 10^-6	6.0 × 10^-7	2.4 × 10^-7	1.7
0.01	5.0 × 10^-5	3.0 × 10^-5	1.7 × 10^-5	6.0 × 10^-6	2.2 × 10^-6	3.1
0.1	8.6 × 10^-4	5.5 × 10^-4	3.5 × 10^-4	6.5 × 10^-5	1.5 × 10^-5	5.5
1	3.1 × 10^-2	2.0 × 10^-2	7.5 × 10^-3	1.5 × 10^-3	1.6 × 10^-4	7.9

如表6所示，K_g取决于温度和应变速率，δ取决于应变速率。对水平集模型参数与工艺参数之间的关系进行分析(图9)。可以看出，lnK_g与lnZ、δ与 $l n \dot{ε}$ 存在线性关系(其中，Z表示Zener-Hollomon参数)，对其分别拟合后得到如下公式：

l n K_{g} = 0.9919 l n Z - 43.751

(29)

δ = 0.9938 l n \dot{ε} + 7.829

(30)

图9

图9 水平集模型参数与工艺参数之间的关系

Fig.9 Relationships between level-set parameters and process parameters (Z—Zener-Hollomon parameter)

(a) lnK_g-lnZ (b) δ-ln $\dot{ε}$

3 模型验证

为了验证模型的有效性，对组织模型(参数为1050 ℃、0.1 s^-1)的仿真结果与热压缩实验获取的实验结果进行比较，包括DRX分数和晶粒尺寸分布。图10为0.4~0.7应变下GH4706合金的晶粒取向分布(grain orientation spread，GOS)图。可以看出，随着应变从0.4增加至0.7，变形晶粒逐渐被等轴晶粒取代，这说明组织分布表现出强烈的应变依赖性。图11为0.4~0.7应变下GH4706合金晶粒尺寸演变的模拟结果。可以看出，随着应变从0.4增加至0.7，变形晶粒不断转化为等轴晶粒，与图10结果非常吻合。图12为晶粒尺寸分布模拟结果与实验结果的对比。可以看出，随着应变的增加，等轴晶粒的体积分数逐渐增加，变形晶粒占比逐渐减小。在应变为0.4时(图12a)，晶粒尺寸更多集中在60 μm以上。随着变形应变从0.5增加至0.7 (图12b~d)，再结晶晶粒占比不断提高(图13c)，晶粒尺寸的集中区间逐渐向左移动，尺寸在100 μm以上的晶粒越来越少，而尺寸位于0~20 μm区间的晶粒越来越多。随着应变增加，变形晶粒占比降低和新形成的等轴晶粒占比增加，同时平均晶粒尺寸逐渐减小。这说明应变对DRX过程中晶粒的均匀性和细化起到重要作用。如图12a和b所示，在小应变条件下，实际晶粒尺寸分布中会出现少量大尺寸晶粒，这可能是由于应变较小时局部位错累积较少，而模拟中的位错分布设置为平均分布。随着应变增加，实际位错分布趋于均匀，此时模拟和实验观测得到的晶粒尺寸分布吻合度显著提高。

图10

图10 不同应变下的晶粒取向分布(GOS)图

Fig.10 Grain orientation spread (GOS) maps at different strains of 0.4 (a), 0.5 (b), 0.6 (c), and 0.7 (d)

图11

图11 不同应变下的晶粒尺寸模拟图

Fig.11 Simulation maps of grain sizes at different strains of 0.4 (a), 0.5 (b), 0.6 (c), and 0.7 (d)

图12

图12 不同应变下实验及模拟的晶粒尺寸分布

Fig.12 Experimental and simulated grain size distributions at different strains of 0.4 (a), 0.5 (b), 0.6 (c), and 0.7 (d)

图13

图13 平均晶粒尺寸和DRX晶粒体积分数随应变演化的模拟与实验结果对比及偏差分析

Fig.13 Comparisons of average grain size (a) and DRX grain volume fraction (c) with different strains acquired from experiment and simulation; and relative deviation (Δ) and standard deviation (ξ) of average grain size (b) and DRX grain volume fraction (d) (Δ_avg—average deviation)

为了比较3种本构模型的预测能力及稳定性，引入相对偏差(Δ)的计算公式：

Δ_{i} = |\frac{E_{i} - P_{i}}{E_{i}}|

(31)

Δ_{a v g} = \frac{1}{x'} \sum_{i = 1}^{x'} |δ_{i}|

(32)

式中，Δ_i 是第i个数据点的相对偏差，E_i 是第i个数据点的实验值，P_i是第i个数据点的预测值，x'是总数据点数，Δ_avg是相对偏差的平均值，能够有效反映模型的预测能力和稳定性。为了对模型预测的精度和稳定性进行更直观的量化评价，引入标准偏差(ξ)作为所有数据中单个点离散度的衡量指标，其表达式如式(33)所示^[35]。ξ越小，表明相对偏差与Δ_avg越接近，反之亦然。越小的ξ表明模型的预测能力越好，预测稳定性越强。ξ可表示为：

ξ = \sqrt[]{\frac{1}{x' - 1} \sum_{i = 1}^{x'} {(Δ_{i} - Δ_{a v g})}^{2}}

(33)

图13为D和F随应变的变化关系，并将模拟与实验结果在0.4~0.7应变范围内进行了比较。图13a和c中的实验和模拟结果表明，随着应变的增加，D和F都在达到ρ_cr (应变0.2~0.3)后迅速变化，然后趋于平稳。D和F的模拟结果均与实验结果吻合较好。图13b和d表明，F的标准偏差(ξ= 5.077)明显大于平均晶粒尺寸的标准偏差(ξ= 3.163)，这说明对于D和F 2种随应变演化的参数来说，模型的预测结果对于后者更加敏感。表7和8列出了不同应变下实验数据与模拟结果的相对偏差和标准偏差。结果表明，2种偏差均在10%以内，这说明，经过Pareto前沿法优化参数后所建立的基于水平集法的动态再结晶模型可有效模拟GH4706合金在热压缩变形过程中的DRX组织演化过程，并可较好地实现不同应变条件下微观组织形貌的清晰可视化，这对于有效调控GH4706合金在热变形过程中的微观组织具有重要意义。

表7 平均晶粒尺寸模拟结果与实验结果的相对偏差

Table 7 Relative deviations between simulated and experimental results of average grain size

Strain	Average grain size / μm		Δ / %	Δ_avg / %	ξ / %
Strain	Experimental	Simulated	Δ / %	Δ_avg / %	ξ / %
0.4	97.1	99.6	2.5	7.75	3.16
0.5	69.9	76.8	9.8
0.6	54.4	60.2	10.6
0.7	39.3	36.1	8.1

表8 DRX晶粒体积分数模拟结果与实验结果的相对偏差

Table 8 Relative deviations between simulated and experimental results of DRX grain volume fraction

Strain	DRX fraction / %		Δ / %	Δ_avg / %	ξ / %
Strain	Experimental	Simulated	Δ / %	Δ_avg / %	ξ / %
0.4	6.9	8.1	14.8	9.75	5.08
0.5	22.4	19.6	12.5
0.6	28.9	32.1	10.3
0.7	42.5	43.1	1.4

4 结论

(1) 通过引入形核和位错等物理机制，根据热压缩实验数据构建了基于水平集法的GH4706合金动态再结晶模型。对于无法拟合获得的参数，使用Pareto多目标优化方法通过使实验与模拟结果偏差值达到最小的方法进行逆向识别，实现了GH4706合金在0.4~0.7应变下的动态再结晶过程的模拟。

(2) 模拟结果与实验结果的对比表明，D和F的模拟结果与实验数据的偏差均在10%以内。这验证了基于水平集法的GH4706合金动态再结晶建模方法的有效性。

(3) 基于水平集法的GH4706合金动态再结晶模拟能够实现不同应变下微观组织形貌演变的清晰可视化。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

Zhang

, Zeng

L R

, Zhao

D Q

, et al.

Comparison study of microstructure and mechanical properties of standard and direct-aging heat treated superalloy Inconel 706

[J]. Mater. Sci. Eng., 2022, A839: 142836

[2]

J H

, Lv

X D

, Dong

J X

, et al.

Research progress of wrought superalloys in China

[J]. Acta Metall. Sin., 2019, 55: 1115

DOI

Wrought superalloys are high temperature alloys produced by casting-forging-hot rolling-cold drawing, including disc, plate, bar, wire, tape, pipe etc. These products are widely used in aviation, aerospace, energy, petrochemical, nuclear power and other industrial fields. In this paper, domestic progress of wrought superalloys in recent ten years was reviewed, including advances in fabrication process, research in new alloys (GH4169G, GH4169D, GH4065 and GH4068 alloy et al.) and new techniques (deforming of FGH4096 alloy, nitriding of NGH5011 alloy and 3D printing of In718 alloy et al.).

杜金辉, 吕旭东, 董建新等.

国内变形高温合金研制进展

[J]. 金属学报, 2019, 55: 1115

DOI

变形高温合金是指通过铸造-变形工艺生产的高温合金，包括盘、板、棒、丝、带、管等产品，该类产品广泛用于航空、航天、能源、石化、核电等工业领域。本文介绍了国内变形高温合金近10年的最新进展，分别从变形高温合金的制备工艺流程，GH4169G、GH4169D、GH4065、GH4068等新合金的研制，以及FGH4096的变形化、NGH5011的氮化、In718合金的3D打印等新技术3个方面展开论述。

[3]

Huang

. Microstructure control and mechanical properties optimization of GH4706 wrought superalloy [D]. Shenyang: Northeastern University, 2014

黄烁

. 变形高温合金GH4706组织控制与力学性能优化 [D]. 沈阳: 东北大学, 2014

[4]

Yang

W H

, deBarbadillo

J J

, Morita

, et al.

A freckle criterion for the solidification of superalloys with a tilted solidification front

[J]. JOM, 2004, 56(9): 56

[5]

Fesland

J P

, Petit

Manufacturing alloy 706 forgings

[A]. Loria E A. Superalloys 718, 625, 706 and Various Derivatives [M]. Warrendale: TMS, 1994: 229

[6]

Sun

Y F

, Xiang

, Zhang

, et al.

Microstructures and mechanical properties of GH4169 superalloy manufactured by selective laser melting

[J]. Acta Metall. Sin., 2024, DOI: 10.11900/0412.1961.2024.00075

孙勇飞, 向超, 张涛等.

选区激光熔化GH4169高温合金的微观组织和力学性能

[J]. 金属学报, 2024, DOI: 10.11900/0412.1961.2024.00075

[7]

, Wang

, Li

, et al.

Effect of the Zener-Hollomon parameter on the dynamic recrystallization kinetics of Mg-Zn-Zr-Yb magnesium alloy

[J]. Comput. Mater. Sci., 2019, 166: 221

[8]

Meng

, Sugiyama

, Yanagimoto

Microstructure of Cr-V-Mo steel processed by recrystallization and partial melting and its effect on mechanical properties

[J]. Mater. Trans., 2014, 55: 921

[9]

F L

, Fu

, Bai

Y R

, et al.

Effects of initial grain size and strengthening phase on thermal deformation and recrystallization behavior of GH4096 superalloy

[J]. Acta Metall. Sin., 2023, 59: 855

DOI [本文引用: 1]

GH4096 alloy were used for disks and shafts of advanced gas turbine engines owing to its excellent properties such as resistance to creep, fatigue, and corrosion as well as microstructure stability up to about 700oC. In this study, GH4096, a hard-to-deform disk superalloy, was processed through an advanced cast and wrought route to avoid the expensive power metallurgy (P/M) route. Many types of full-scale disk forgings possessing homogeneous fine-grained microstructures were successfully carried out, and the ultrasonic inspectability was comparative to that of the alloy produced by the P/M route. The effects of the initial grain size and strengthening phase on hot deformation behavior and dynamic recrystallization (DRX) were studied by OM, SEM, EBSD, and TEM under different deformation parameters. The results showed that as the initial grain size decreased within the temperature range of 1050-1120oC, the flow peak stresses decreased and the fractions of DRX increased. With an increase in the initial grain size, the thermal deformation temperature required for complete dynamic recrystallization decreased, and also the critical strain of dynamic recrystallization decreased. The initial grain size and the strain did not affect the recrystallized grain size when deformed at a sub-solvus temperature. The thermal deformation constitutive equations related to the initial grain sizes were established and the activation energies of thermal deformation related to the original grain sizes were calculated. The effect of γ' phase size on the thermal deformation behavior in as-cast microstructure was studied. In the sub-solvus temperature range, the thermal deformation resistance could be effectively reduced with the increase in the size of γ' phase, the critical strain of DRX was decreased, and the DRX fraction was also increased. The dynamic recrystallization mechanisms related to the γ' phase and initial grain size were also discussed. DRX nucleation takes place at the sub-grains near original grain boundaries for samples with larger initial grain size deformed at sub-solvus temperature. For samples with fine initial grain size, the interface slip of incoherent γ' phase is the significant dynamic softening mechanism during the sub-solvus temperature deformation. For as-cast samples, the main dynamic softening mechanism is original grain boundary bowing out DRX nucleation and coarse second-phase-induced DRX nucleation.

李福林, 付锐, 白云瑞等.

初始晶粒尺寸和强化相对GH4096高温合金热变形行为和再结晶的影响

[J]. 金属学报, 2023, 59: 855

DOI [本文引用: 1]

以镍基变形高温合金GH4096为研究对象，进行了不同组织状态下的热压缩实验，分析了该合金热变形流变行为特点，并采用OM、SEM、EBSD和TEM等手段研究了热压缩过程中初始晶粒尺寸(晶界面积和数量)和强化相对动态再结晶行为及组织演变的影响。结果表明：在1050~1120℃温度范围内，随着初始晶粒尺寸的减小，峰值流变应力降低，动态再结晶分数增加，实现完全动态再结晶所需热变形温度降低，发生动态再结晶的临界应变降低；在亚固溶温度下热变形时，动态再结晶晶粒尺寸与初始晶粒尺寸和变形量没有相关性；计算了与初始晶粒尺寸相关的热变形激活能，并建立了与热变形前的初始晶粒尺寸相关的热变形本构方程；研究了铸态组织中强化相γ'相形态对热变形行为的影响，在亚固溶温度下热变形时，随着强化相γ'相尺寸的增加，可一定程度降低热变形峰值流变应力，降低动态再结晶的临界应变，提高动态再结晶分数。分析讨论了GH4096合金不同组织状态下的主导动态软化机制。

[10]

Avrami

Kinetics of phase change. II Transformation-time relations for random distribution of nuclei

[J]. J. Chem. Phys., 1940, 8: 212

[11]

Sellars

C M

The physical metallurgy of hot working

[A]. Hot Working and Forming Processes: Proceedings of an International Conference on Hot Working and Forming Processes [M]. London: The Metals Society, 1980

[12]

Jonas

J J

, Quelennec

, Jiang

, et al.

The Avrami kinetics of dynamic recrystallization

[J]. Acta Mater., 2009, 57: 2748

[13]

Rollett

A D

, Raabe

A hybrid model for mesoscopic simulation of recrystallization

[J]. Comput. Mater. Sci., 2001, 21: 69

[14]

Zhu

H J

, Chen

, Zhang

H M

, et al.

Review on modeling and simulation of microstructure evolution during dynamic recrystallization using cellular automaton method

[J]. Sci. China Technol. Sci., 2020, 63: 357

[15]

Q H

, Zhang

L W

, et al.

Cellular automaton modeling of dynamic recrystallization of nimonic 80A superalloy based on inhomogeneous distribution of dislocations inside grains

[J]. J. Mater. Eng. Perform., 2018, 27: 4955

[16]

Han

Z Q

. Simulation of dynamic recrystallization based on Monte Carlo method [D]. Ji'nan: Shandong University, 2007

韩振强

. 基于Monte Carlo方法的金属动态再结晶组织模拟 [D]. 济南: 山东大学, 2007

[17]

Zhou

G W

, Li

Z H

, Li

D Y

, et al.

A polycrystal plasticity based discontinuous dynamic recrystallization simulation method and its application to copper

[J]. Int. J. Plast., 2017, 91: 48

[18]

Hallberg

Approaches to modeling of recrystallization

[J]. Metals, 2011, 1: 16

[本文引用: 4]

[19]

Bernacki

, Logé

R E

, Coupez

Level set framework for the finite-element modelling of recrystallization and grain growth in polycrystalline materials

[J]. Scr. Mater., 2011, 64: 525

[20]

Scholtes

, Boulais-Sinou

, Settefrati

, et al.

3D level set modeling of static recrystallization considering stored energy fields

[J]. Comput. Mater. Sci., 2016, 122: 57

[21]

Cruz-Fabiano

A L

, Logé

, Bernacki

Assessment of simplified 2D grain growth models from numerical experiments based on a level set framework

[J]. Comput. Mater. Sci., 2014, 92: 305

[22]

Jin

, Bozzolo

, Rollett

A D

, et al.

2D finite element modeling of misorientation dependent anisotropic grain growth in polycrystalline materials: Level set versus multi-phase-field method

[J]. Comput. Mater. Sci., 2015, 104: 108

[23]

Bernacki

Kinetic equations and level-set approach for simulating solid-state microstructure evolutions at the mesoscopic scale: State of the art, imitations, and prospects

[J]. Prog. Mater. Sci., 2024, 142: 101224

[24]

Agnoli

, Bernacki

, Logé

, et al.

Selective Growth of low stored energy grains during δ sub-solvus annealing in the inconel 718 nickel-based superalloy

[J]. Metall. Mater. Trans., 2015, 46A: 4405

[25]

Gao

. Research on topology optimization for multiscale design of structure-material based on parametric level set [D]. Wuhan: Huazhong University of Science and Technology, 2019

高杰

. 基于参数化水平集的结构/材料多尺度拓扑优化设计研究 [D]. 武汉: 华中科技大学, 2019

[26]

Peczak

, Luton

M J

The effect of nucleation models on dynamic recrystallization I. Homogeneous stored energy distribution

[J]. Philos. Mag., 1993, 68B: 115

[本文引用: 4]

[27]

Cram

D G

, Zurob

H S

, Brechet

Y J M

, et al.

Modelling discontinuous dynamic recrystallization using a physically based model for nucleation

[J]. Acta Mater., 2009, 57: 5218

[28]

Guo

Q M

, Li

D F

, Peng

H J

, et al.

Nucleation mechanisms of dynamic recrystallization in Inconel 625 superalloy deformed with different strain rates

[J]. Rare Met., 2012, 31: 215

[29]

Rios

P R

, Siciliano Jr

F S

, Sandim

H R Z

, et al.

Nucleation and growth during recrystallization

[J]. Mater. Res., 2005, 8: 225

[30]

Beltran

, Huang

, Logé

R E

A mean field model of dynamic and post-dynamic recrystallization predicting kinetics, grain size and flow stress

[J]. Comput. Mater. Sci., 2015, 102: 293

[31]

J C

, Wu

X D

, Liao

, et al.

Simulation of low proportion of dynamic recrystallization in 7055 aluminum alloy

[J]. Trans. Nonferrous Met. Soc. China, 2021, 31: 1902

[32]

Quan

G Z

, Zhang

K K

, An

, et al.

Analysis of dynamic recrystallization behaviors in resistance heating compressions of heat-resistant alloy by multi-field and multi-scale coupling method

[J]. Comput. Mater. Sci., 2018, 149: 73

[33]

Reyes

L A

, Garza

, Delgado

, et al.

Cellular automata modeling for rotary friction welding of Inconel 718

[J]. Mater. Manuf. Processes, 2022, 37: 877

[34]

W Z

. Research on the modeling and optimization algorithms for the slab design and the hot rolling production planning of steel plate [D]. Chongqing: Chongqing University, 2019

呼万哲

. 中厚板坯料设计及其热轧生产计划建模与优化算法研究 [D]. 重庆: 重庆大学, 2019

[35]

Quan

G Z

, Zhang

, Lei

, et al.

Characterization of flow behaviors by a PSO-BP integrated model for a medium carbon alloy steel

[J]. Materials, 2023, 16: 2982