简介概要

显微偏析数值模拟的微观模型和数值方法

来源期刊：中国有色金属学报2004年第3期

论文作者：刘永刚孙瑜邱宗文孙国雄

文章页码：335 - 339

关键词：显微偏析；数值模拟；微观模型；数值方法

Key words：microsegregation; simulation; micro models; numerical method

摘要：研究显微偏析数值模拟中的微观模型和数值方法，在微观模型中考虑了树枝晶粗化、固相溶质逆扩散、枝晶尖端过冷、随温度变化的溶质扩散系数等影响显微偏析形成的动力学因素。数值方法中采用变网格技术跟踪移动界面，通过迭代求解溶质扩散方程和溶质守恒方程计算显微偏析参数。并给出了微观计算流程图，这种数值方法可以方便地与宏观计算、相图计算相耦合。模拟算例表明模拟结果与实验结果吻合较好。

Abstract: The numerical methods and micro models in the simulation on microsegregation were investigated. The main kinetic effects that can influence microsegregation were taken into account in the modeling such as secondary arm coarsening, solid state back diffusion, primary tip undercooling and variable diffusion coefficients according to temperature. The deforming grids technology was applied for tracking the moving interfaces, and the iterative method was applied for resolving diffusion equation and solute balance equation. The flow chart for the micro calculation which can be coupled easily with macro calculation and calculation of phase diagrams were also provided. A example calculation was carried out and the simulation results were in good agreement with experiment results.

中国有色金属学报 2004,(03),335-339 DOI:10.19476/j.ysxb.1004.0609.2004.03.002

显微偏析数值模拟的微观模型和数值方法

刘永刚孙瑜邱宗文孙国雄

东南大学机械系,东南大学机械系,东南大学机械系,东南大学机械系南京210096 ,南京210096 ,南京210096 ,南京210096

摘要：

研究显微偏析数值模拟中的微观模型和数值方法,在微观模型中考虑了树枝晶粗化、固相溶质逆扩散、枝晶尖端过冷、随温度变化的溶质扩散系数等影响显微偏析形成的动力学因素。数值方法中采用变网格技术跟踪移动界面,通过迭代求解溶质扩散方程和溶质守恒方程计算显微偏析参数。并给出了微观计算流程图,这种数值方法可以方便地与宏观计算、相图计算相耦合。模拟算例表明模拟结果与实验结果吻合较好。

关键词：

显微偏析;数值模拟;微观模型;数值方法;

中图分类号： TG111.4

作者简介：刘永刚(1976),男,博士研究生.;

收稿日期：2003-06-11

基金：国家自然科学基金资助项目(59974011);

Micro models and numerical methods of microsegregation simulation

Abstract：

The numerical methods and micro models in the simulation on microsegregation were investigated. The main kinetic effects that can influence microsegregation were taken into account in the modeling such as secondary arm coarsening, solid state back diffusion, primary tip undercooling and variable diffusion coefficients according to temperature. The deforming grids technology was applied for tracking the moving interfaces, and the iterative method was applied for resolving diffusion equation and solute balance equation. The flow chart for the micro calculation which can be coupled easily with macro calculation and calculation of phase diagrams were also provided. A example calculation was carried out and the simulation results were in good agreement with experiment results.

Keyword：

microsegregation; simulation; micro models; numerical method;

Received： 2003-06-11

作为凝固过程微观数值模拟的研究方向之一, 显微偏析数值模拟是为了计算铸件或铸锭的枝晶参数及枝晶偏析(包括二次枝晶臂间距、在二次枝晶臂中的浓度分布、固相分数和共晶的数量), 从而预测材料的组织和性能以及指导进行正确的铸件后续加工过程 ^[1,2] 。近年来已经出现了一些先进的微观模拟模型 ^{[3,4,5,6,7,8,9]} , 在这些模型中考虑了固相扩散、枝晶粗化、过冷等对显微偏析的影响。目前, 预测铸件或铸锭不同位置的枝晶偏析, 一般需要宏/微观计算相耦合, 并与相图计算耦合, 关于这2种耦合的数值方法, 作者已另文论述 ^[10] 。本文中作者主要研究显微偏析模拟的微观模型和数值方法, 为进行显微偏析数值模拟的研究提供模型和数值方法基础。为了验证模型与算法的正确性和适用性, 作者以二元合金为例给出了算例分析, 模拟结果与实验测定结果吻合较好。

1主要的显微偏析模型

1.1研究对象

显微偏析计算主要在二次枝晶臂间的微观尺度内进行, 如图1所示的0<x<x(t)部分, t代表时间。这个区域的大小x(t)为二次枝晶臂间距的一半, 当凝固进行时是不断向右扩展的, 这就是树枝晶的粗化现象, 与此同时, 固/液界面s(t)也在不断移动。显微偏析模拟的输出参数就是凝固过程中树枝晶二次臂间距、二次臂内的溶质分布和共晶相体积分数的变化等。在图1中采用的平板状几何模型近似树枝晶形状, 尽管是对枝晶部分的简单近似, 但有研究结果表明, 除了极少数情况下, 如冷却速度的突然变化或者存在很强的强制对流, 这种几何模型对于大多数合金系统和凝固过程都可以有效描述近似枝晶形貌 ^[11] 。

图1 显微偏析模型的枝晶形态

Fig.1 Dendritic morphology in microsegregation model

1.2主要的数学物理模型

粗化方程采用Roósz等提出的动力学模型 ^[12] :

[λ(t)]³=∫ $_{0}^{t}$ BMd t (1)

式中 λ(t)为二次枝晶臂间距的一半; B为一形状因子, 对于特定合金有特定的数值, 对于本文所选的Al-4.5%Cu合金, B取13.125 ^[13,14] ; 粗化参数M的表达式为

$Μ = \frac{γ D_{L} Τ}{m (1 - k) Δ Η c_{L}^{*}} (2)$

式中 ΔH为熔融焓; γ为固液界面能; D_L为溶质在液相中的扩散系数; T为温度; m为液相线斜率; k为平衡溶质分配系数; c^*_L为固液界面上液相中的溶质浓度。

固相中的溶质再分配符合菲克第二定律:

$\frac{\partial c_{S}}{\partial t} = \frac{\partial}{\partial x} (D_{S} \frac{\partial c_{S}}{\partial x}) (3)$

式中 c_S为固相溶质浓度; D_S为溶质在固相中的扩散系数。

考虑在二次臂内区域的封闭系统, 溶质分配应满足以下守恒控制方程:

c^*_L[x(t)-x_S]+∫ $_{0}^{x_{S}}$ c_Sd x=x(T)c₀ (4)

式中 x(t)为体积元的长度; x_S为固相长度, c₀为合金名义成分。

枝晶尖端过冷采用KGT模型迭代计算 ^[15] , 该模型是Kurz, Giovanola和Trivedi在研究枝晶尖端溶质平衡的基础上建立的, 定量地描述了枝晶尖端半径r、枝晶尖端生长速度v_g和尖端过冷度ΔT之间的关系。主要包含以下各式:

$\begin{array}{l} Δ Τ = m \cdot c_{0} \cdot [1 - \frac{1}{1 - Ω (1 - k)}] (5) \\ Ω = \frac{c_{L}^{*} - c_{0}}{c_{L}^{*} (1 - k)} = Ι_{v} (Ρ_{e}) (6) \\ Ρ_{e} = \frac{r \cdot v_{g}}{2 D_{L}} (7) \\ r = 2 π [\frac{Γ}{m \cdot G_{c} \cdot ξ_{c} (Ρ_{e}) - G}]^{1 / 2} (8) \end{array}$

式中 Ω为枝晶尖端的溶质过饱和度; P_e为溶质的Péclet数; ξ_c为与P_e有关的常数 ^[15] ; G_c为浓度梯度; G为温度梯度; I_v为Ivantsov函数, 可表示为

I_v(P_e)=P_e·exp(P_e)·E_i(P_e) (9)

其中 E_i为指数积分函数,

$E_{i} (Ρ_{e}) = \int_{Ρ_{e}}^{\infty} \frac{\exp (- z)}{z} \cdot d z, z$ 为积分变量; Г为Gibbs-Thomson系数, $Γ = \frac{γ \cdot Τ}{Δ Η}$ 。

2微观模型的数值方法

宏观计算出铸件各结点的冷却曲线以后, 输入显微偏析微观模型中计算各结点的显微偏析参数, 其流程如图2所示。

在微观计算中, 对于每一结点, 首先读取宏观计算的冷却曲线, 计算树枝晶尖端过冷, 在KGT模型中需要用到温度梯度G和凝固速度或生长速度(固/液界面前进速度)v_g, 而由冷却曲线仅仅可以得到该结点的温度和冷却速度v_c。温度梯度需要由相邻结点的温度给出, 已知冷却速度和温度梯度后, 生长速度可以由下式得出 ^[16] :

$v_{g} = \frac{v_{c}}{G} (10)$

得到枝晶尖端过冷度以后, 将局部凝固时间t_f离散为一定数目的时间段。在每一时间步长计算开始时, 首先从热力学计算相图读取界面平衡信息, 求解粗化方程获得此时二次臂间距的值, 然后就可以使用Scheil方程计算出此时的固相体积分数, 并可以由得到的固相体积分数通过下式估算此时固相长度和固相长度在此时间步长内的增加量Δx_est:

$f_{s} = \frac{x_{S} (t)}{x (t)} (11)$

然后求解溶质扩散方程和溶质守恒方程。求解溶质扩散方程时使用的固相溶质浓度为该结点当前时刻通过求解溶质扩散方程得到的各微观结点的溶质浓度, 液相溶质浓度为相图液相线给出的液相成分。通过这个求解溶质分布的过程又可以得到一个具体计算出的固相长度增加量Δx_cal, 当这2个值不符合收敛条件时, 用计算出的Δx_cal替代Δx_est, 再次求解扩散方程与守恒方程, 直到迭代收敛为止。值得注意的是, 由于在微观计算中, 时间步长是由等分局部凝固时间t_f预先给定的, 所以求解扩散方程必须采用有限差分的全隐式求解方法, 因为这个预先给定的时间步长未必满足显式差分的求解稳定性条件。

图2 显微偏析数值模拟微观数值方法

Fig.2 Scheme of numerical method of micro model for microsegregation simulation

重复上述过程, 直至该结点的温度到达合金的共晶温度, 即所有的微观时间步长循环完毕后, 认为剩余液相瞬时转变为共晶体, 然后进入下一结点的计算, 当所有宏观节点都运算完毕以后, 就可以得到整个铸件或铸锭的显微偏析参数。

在微观计算中, 需要解决的关键问题是处理好2个移动的界面: 一个是由于树枝晶粗化引起的控制容积的扩展, 可以由粗化方程得到; 另一个是在凝固过程中固/液界面的向前移动。关于移动的固/液界面的跟踪, 由于在每一微观计算时间步长内, 固相长度增加是不同的, 采用简单的固定网格、固定时间步长技术计算是存在困难的。本文采用定时间步长、递增节点数目的变网格技术跟踪移动的界面, 即在每一微观时间步长内增加一定的固相长度, 也就增加了一个固相结点, 参与本次的迭代运算。离散化方法如图3所示, 可见, 每一时间步长凝固出的固相长度被分为2部分, 分别属于不同的结点参加固相扩散计算。由于固相内溶质浓度是递增的, 这样就避免了由此结点代表一个整时间段形成固相长度所带来的误差。结点所在位置为每一微观时间步长后界面所在的位置, 从而达到显式跟踪自由边界的目的。

图3 微观计算离散化方法

Fig.3 Numerical dispersed method in micro calculation

3模拟结果的验证与讨论

以水冷铜型内冷却的Al-4.5%Cu合金铸锭为例进行了仿真计算, 铸锭尺寸为d 50 mm×100 mm, 水冷温度为300 K, 浇注温度980 K, 铸锭上表面直接和空气接触, 下表面绝热。考虑到圆柱形的旋转对称性, 取其纵截面作为二维模拟计算区域。计算结果表明, 所采用的模拟算法有很好的计算效率, 整个模拟过程在Pentium PC-500 MHz、内存128 MB微机上的运行时间大约为50 min。

3.1冷却曲线

温度场计算结果的准确性是显微偏析计算的基础和保障, 为了验证宏观温度场计算的结果, 测定了铸锭内若干点的冷却曲线, 图4所示为实验测定的冷却曲线与模拟结果的对照。测定点距离圆心18 mm, 图中细线为模拟结果, 方形散点曲线为实验测定的冷却曲线。由图4可以看出, 计算结果与实验测定结果吻合较好。

图4 实验测定冷却曲线与计算结果的对比

Fig.4 Comparison of cool curves between calculated results and experimental results

3.2二次臂间距与共晶相分数

表1给出了不同位置处的二次臂间距、共晶相分数的计算结果与实验测定结果, 此处得到的实验数据均是统计结果。图像采集与数据的测定均采用自主开发的商业金相检测软件CIMAS, 数据的可靠性是经过验证的。在金相分析过程中, 由图像直接测定的结果实际上并不是共晶相分数, 而是CuAl₂相的体积分数, 需要把测定结果换算成共晶相的分数再和模拟结果进行对比。由表1可以看出, 仿真计算获得了与实验一致的结果, 随着离开圆心距离的减小, 树枝晶二次臂间距在增加, 共晶相数量在增加, 因为模拟过程并没有考虑实际铸锭凝固过程中共晶温度以下固溶体析出的第二相数量, 所以测定的共晶分数大于模拟的结果, 这与理论分析也是一致的。

表1 不同位置处二次臂间距与共晶相分数的实验结果与模拟结果的对比

Table 1 Comparison between calculated resultsand experimental results of secondary dendritearm spacing and volume fraction ofeutectic for different locations

Distance from center/ mm	Secondary dendrite arm spacing/μm		Volume fraction of eutectic%
Distance from center/ mm	Experimental result	Calculated result	Experimental result	Calculated result
19	26	26.3	7.91	7.64
15	26	26.5	8.39	8.17
11	28	26.6	9.08	8.75
8	29	27.0	10.20	9.15