文章编号：1004-0609(2007)10-1567-07

基于概率捕获模型的元胞自动机方法

模拟镁合金枝晶生长过程

付振南, 许庆彦, 熊守美

(清华大学 机械工程系 先进成形制造教育部重点实验室，北京 100084)

摘  要：基于枝晶生长的基本传输过程和元胞自动机(Cellular Automaton，简称CA)模型基本原理，建立符合镁合金密排六方结构特点的凝固过程形核、晶体生长的数学模型。模型耦合宏观温度场和微观组织模拟计算，考虑溶质扩散、曲率过冷和各向异性等重要因素的影响，定义界面单元捕获规则，能够模拟枝晶生长的形态。应用本模型对单个等轴晶生长、等轴晶多晶粒生长及定向凝固中的柱状晶生长进行模拟并与实验结果进行对比，模拟结果与实验结果吻合较好，验证了模型的正确性。

关键词：镁合金；密排六方结构；微观组织模拟；元胞自动机方法；枝晶形貌

中图分类号：TG 244; O 242　　     文献标识码：A

Numerical simulation on dendrite growth process of Mg alloy using cellular automaton method based on probability capturing model

FU Zhen-nan, XU Qing-yan, XIONG Shou-mei

(Department of Mechanical Engineering, Tsinghua University, Key Laboratory for Advanced Materials Processing Technology, Ministry of Education, Beijing 100084, China)

Abstract: A new model for the nucleation and dendritic growth of magnesium alloys during solidification process was developed, which was based on cellular automaton(CA) method and the classical transfer equations during dendritic growth and is suitable for the hexagonal close-packed crystal structure. The model coupled the temperature field calculation and microstructure simulation. The solute diffusion, curvature undercooling and anisotropy were also considered in the model, all of which have significant influence on the dendrite evolution. By setting up the CA-cell capturing rules, the dendritic growth morphology can be modeled. Applications of the model in single equiaxed dendritic growth, multiple equiaxed dendritic growth and columnar grain growth in directional solidification were presented. The microstructure simulation results were compared with the experimental results and they were in good agreement in dendritic morphology, which validates the proposed model in the prediction of dendritic growth of magnesium alloys.

Key words: magnesium alloy; hexagonal close-packed structure; microstructure simulation; cellular automaton method; dendritic morphology

镁合金具有密度小、比强度和比刚度高、阻尼性能好和导热性能好等优点，在汽车、通讯电子和航空航天等领域正得到日益广泛的应用^[1?3]。镁合金铸件的力学性能强烈地依赖于其微观组织(特别是晶粒尺寸)^[4]。如果能通过计算机模拟的方法预测出镁合金铸件的微观组织，就可以据此调整生产工艺，从而得到力学性能优良的铸件。

元胞自动机方法自从20世纪90年代初被提出，后来国内外很多研究者应用CA方法对晶粒域和枝晶域的微观组织进行了数值模拟，但是以往的微观组织模拟研究绝大多数针对立方晶系合金开展，而镁合金属于密排六方晶系，采用立方晶系模型无法体现其各向异性特征，因此不能准确揭示镁合金枝晶生长的组织形貌。由于密排六方晶体结构生长建模和验证的复杂性，因此针对镁合金密排六方结构的微观组织模拟研究还处于起步阶段，目前可以查阅到的文献和资料较少。最近才有德国的Eiken等采用相场模型针对镁合金进行了枝晶生长模拟，引入六方各向异性函数描述密排六方晶体结构，并对AZ31合金的定向凝固及等轴生长进行模拟^[5]。

本文作者基于CA方法及枝晶生长的基本传输过程，建立了镁合金凝固过程中枝晶形貌演化的数学模型，模型耦合了微观溶质场计算，考虑了成分过冷、曲率过冷和各向异性等重要因素。该模型能够模拟镁合金密排六方结构树枝晶的复杂形貌。通过对单个等轴晶、等轴晶多晶粒和定向凝固的模拟结果表明，本模型具有较好的理论和应用价值。

1  镁合金的晶体学特征与枝晶形貌

为了准确模拟镁合金的形核与枝晶生长，必须考虑其晶体学特征。镁合金晶体结构的重要特征是具有密排六方结构，在晶体基面上具有六重对称性，在基面的法向上具有二重对称性。根据实验观察^[5?7]，在非平衡凝固条件下，在具有六重对称性的晶体基面上的生长对镁合金枝晶组织起主导影响。在用冰醋酸腐蚀剂腐蚀后^[8]，AZ91D砂型铸件在光学显微镜偏振光下的枝晶结构如图1所示。图2所示为金属型AZ91D铸件的SEM像。可以看出，枝晶主干与二次臂间夹角约为60?，从而直接验证了镁合金具有密排六方结构特征。AZ91D镁合金初生α-Mg相呈树枝晶形貌，一次枝晶互成60?，沿着{0001}基面的<>方向对称生长。由于镁合金的凝固区间宽，在晶粒长大的过程中，枝晶主干上又可长出<>晶向的二次臂，如果一次枝晶的间距足够大，在二次臂上还有可能会长出三次臂，使得α-Mg枝晶呈现高度的分枝形态。考虑镁合金所具有的实际枝晶形貌特征，假设α-Mg基面上的枝晶结构如图3所示。

图1  砂型AZ91D铸态的微观组织

Fig.1  As-cast microstructure of sand mould AZ91D casting (under polarized light)

图2  金属型AZ91D铸件的SEM像

Fig.2  SEM photograph of metal mould AZ91D casting

图3  α-Mg基面枝晶结构示意图

Fig.3  Dendritic structural schematics of basal plane in α-Mg phase

2  镁合金微观组织形成的数理模型

在镁合金的凝固过程中，温度分布和溶质分布是影响枝晶生长和最终形貌的重要因素。假设镁铝合金为二元合金，考虑在一个初始均匀温度为T₀、初始成分为C₀的三维区域Ω内，其枝晶凝固过程的数学描述如图4所示意。这里，n为界面法向矢量，为界面的平均曲率, 曲线Γ代表随时间而变化的固液界面前沿。沿着z方向观察xy平面上的二维枝晶形貌如图5所示。

图4  镁合金枝晶凝固的三维描述图

Fig.4  Three-dimensional domain for dendritic solidification of magnesium alloy

图5  镁合金枝晶凝固的二维描述图

Fig.5  Two-dimensional domain for dendritic solidification of magnesium alloy

2.1  凝固过程的数学描述

根据能量守恒原理，传热过程的控制方程为

液相中的溶质扩散方程如下：

固相中的溶质扩散方程如下：

随着固液界面的推进，固相会向液相中排出溶质，固液界面推进过程中保持溶质守恒^[9]：

对方程(4)进行有限差分，可计算出固液界面单元具有x方向的生长速度v_x、y方向的生长速度v_y以及z方向的生长速度v_z，如下所示：

满足如下关系：

f_L+f_S=1                                         (8)

固液界面的溶质满足局部平衡，如下所示：

                                    (9)

式中  为固液界面的固相溶质浓度。

固液界面上的液相浓度可表示如下：

2.2  形核模型

为了更准确地反映实际情况，形核采用基于高斯分布的连续形核模型^[11]。 该模型假设晶核的出现不需要时间，晶粒密度随过冷度增大以慢?快?慢的趋势连续长大。在给定过冷度下的新生晶粒密度如下：

2.3  凝固前沿的单元捕获规则

根据镁合金在二维基面上具有六重对称性的特征，定义了CA方法中二维邻居单元如图6所示，其中黑色方框为某一CA单元，灰色方框为该CA单元的邻居^[12?13]。假设在模拟区域中每个CA单元共有3种状态：液相、固相以及固液界面状态。在凝固的开始阶段，所有的CA单元都为液相。当温度下降到液相线以下，一些CA单元将会随着形核和晶粒生长过程变成固液界面单元(以后是固相单元)。捕获规则定义如下^[14]：1) 当液相单元发生形核或者液相单元被固相单元捕获时，它变为固液界面单元；2) 固液界面单元随着凝固过程的进行，固相分数的增加值?f_S用下式计算：

3)当固相分数等于1时，固液界面单元变为固相，固相单元有捕获周围邻居单元的能力。每个固相单元捕获周围液态邻居单元的可能性通过概率计算。在图5中，固相单元a用黑色方框表示，其下方的液态邻居单元b用灰色方框表示。邻居单元b被捕获的概率p_c由该单元与晶核(即o点)的连线与该晶粒优先生长方向的夹角ω所决定(由于镁合金在二维基面上有六个等同的优先生长方向，故0≤ω≤π/3)，概率捕获函数如式(16)所示。在计算过程中，每个液态邻居单元都会获取一个随机数r(0≤r≤1)，当r≤p_c时，该单元能够被周围凝固的固相单元捕获并由液态单元转换为固液界面单元，然后随着凝固过程的进行，固液界面单元的固相分数不断增加直到变成固相。上述固相单元捕获液态单元的过程将会一直重复直到所有的单元都转变为固相为止。在邻居单元定义规则和概率捕获函数的双重作用下，一次枝晶生长形状将会呈现类似于图6中的黑色轮廓线所示意的六重对称性特征。

图6?邻居单元的定义

Fig.6  Definition of neighbourhood configuration

3  微观组织模拟结果与分析

3.1  单个等轴晶生长模拟

模拟程序中所用的AZ91D合金的热物性参数如表1所示，在微观域的中心人为种下一个晶核，除此之外微观计算域的他地方不考虑形核的影响。假设液态初始温度为680 ℃，初始溶质浓度均匀，冷却速度为20 K/s。图7(a)所示为均匀微观温度场条件下，AZ91D合金单个枝晶的三维溶质场模拟结果，图7(b)所示为沿着中心在xoy二维平面上的溶质场模拟结果。对比图7(b)和图1可知，单个等轴晶枝晶形态的模拟结果与实验结果得到的等轴枝晶形态非常相似。镁合金枝晶生长表现出显著的密排六方结构特征：晶体的优先生长方向为<>晶向族，一次枝晶臂之间互成60?角左右，在一次枝晶臂主干和二次枝晶臂之间的夹角同样约为60?。从图7(b)中可以看出，6个一次枝晶臂之间夹角并不正好都是60?，6个一次枝晶臂的长度也不是完全相等，这是由于邻居单元在x方向和y方向的定义并不是完全一致的，从而在邻居单元和概率捕获函数的双重作用下，形成了图6所示大致呈现六重对称性特征的枝晶形态。

表1? 模拟中用的AZ91D合金热物性参数

Table 1  Thermal physical parameters of AZ91D alloy used in simulation

图7  单个枝晶的溶质分布模拟图

Fig.7  Simulated solute distribution scheme of equiaxed dendrite: (a) Three dimension result; (b) Two dimension result

从图7(b)还可以看出，在固液界面前沿和枝晶臂间液相中的溶质浓度要高于远处液相中的溶质浓度。因为在枝晶生长过程中，已凝固的固相单元不断捕获周围的液相单元，由于溶质在液相中的溶解度大于固相中的溶解度，因此在液相转变为固相的过程中向邻近液相排出大量溶质，而固液界面处溶质的排出速度大于溶质在液相中的扩散速度，因此造成了固液界面处的溶质富集。根据枝晶生长动力学，溶质浓度越高，过冷度越小，生长速度越小。在枝晶臂间的液相中，由于溶质无法及时扩散出去，使得浓度越来越高，造成枝晶臂间的溶质富集，因此伸向此处的枝晶臂生长速度减缓或停止。在枝晶生长前沿，溶质扩散畅通，且溶质富集程度相对较小，所以生长速度较快。AZ91D合金单个枝晶在晶体学优先生长方向和溶质场的双重影响下，形成了图7所示的枝晶形貌。

3.2  等轴晶多晶粒生长模拟

等轴晶多晶粒生长模拟可以用于预测晶粒度、枝晶间距等参数。采用本文提出的模型对AM50金属型铸件进行等轴晶多晶粒生长模拟，计算区域为1 mm×1 mm的正方形，划分成500×500个CA单元，每个CA单元尺寸为δx=δy=2 μm，初始温度为730 ℃。计算过程中采用宏观温度场与微观组织模拟耦合来获得固液界面前沿局部温度场的数据。

图8(a)所示为凝固进程为90%时的枝晶组织模拟图，图8(b)所示为AZ91D金属型铸件在金相显微镜偏振光下的照片。从图中可以看出，本模型可以模拟随机晶向的多枝晶生长演化进程。在计算域内随机分布的等轴枝晶已经产生了二次枝晶臂，形成六重对称性形状特征的枝晶形貌。枝晶晶粒互相接触并竞争生长，先形核的晶粒较大且抑制后形核晶粒的生长，在两个或多个晶粒之间的晶粒生长受到空间限制和局部溶质富集的影响，其一次枝晶臂或二次枝晶臂生长速度减慢并最终提前结束生长。如果生长前沿不受阻碍，枝晶臂仍会继续向过冷液相中伸展，因此虽然初始形核的枝晶具有规则的六重对称形状，但最终也会形成不规则的枝晶形貌。采用截线法测得图8(a)与图8(b)的平均晶粒尺寸分别为206 μm和227 μm。由此可见，模拟结果和实验结果在枝晶形貌和平均晶粒尺寸方面吻合较好。

图8  等轴晶多晶粒枝晶组织生长模拟与实验对比

Fig.8  Simulated results of dendritic structure of equiaxed solidification in comparison with experimental result: (a) 90% solidified; (b) Experimental result

3.3  定向凝固模拟

图9所示为AZ31镁合金定向凝固的模拟结果与德国Eiken等实验得到的显微照片结果对比^[15]。微观组织模拟的区域为边长2 mm的正方形，剖分为500×500个单元。模拟程序中的温度场设置成与实验条件一致，温度梯度为竖直方向10 K/mm，自上而下温度依次降低，冷却速度为1.5 K/s。

图9  AZ31镁合金定向凝固的模拟结果与实验对比

Fig.9  Simulations of directional growth of AZ31 magnesium alloy in comparison with experimental result: (a) Simulated result; (b) Experimental result^[15]

为了再现柱状枝晶的竞争生长特征，假设在计算区域底端的长度方向间隔120 μm形成一个初始晶核，共有15个晶核并具有不同的枝晶生长晶向。微观组织模拟的枝晶生长形态如图9(a)所示，晶向角度不同的各晶核同时开始捕获周围的液态单元，不断向上生长，在主枝晶干方向又生出二次枝晶臂，密排六方晶系镁合金的枝晶与立方晶系枝晶最大的区别是枝晶干与枝晶臂之间的夹角是60?而不是90?。由于枝晶之间的竞争生长，最后只有部分枝晶臂发展充分，而其它枝晶由于受到抑制不能充分生长，存在于底部区域，无法进一步长大。图9(b)所示为德国Eiken等的实验结果，通过对比可以看出，模拟结果与文献报道的实验结果在枝晶形貌上比较相似。

REFERENCES

[1] Aghion E, Bronfin B. Magnesium alloys development towards the 21st century[C]//Switzerland: Trans Tech Publications. Materials Science Forum, 2000, 350/351: 19?28.

[2] Mordike B L, Ebert T. Magnesium properties- applications-potential[J]. Mater Sci Eng A, 2001, A 302: 37?45.

[3] Friedrich H, Schumann S. Research for a “new age of magnesium” in the automotive industry[J]. Journal of Materials Processing Technology, 2001, 117: 276?281.

[4] 赵浩峰, 刘红梅, 张永忠, 范晋平, 蔚晓嘉, 郝兴明. 镁合金晶粒细化的研究[J]. 轻合金加工技术, 2004, 32(2): 10?13.
ZHAO Hao-feng, LIU Hong-mei, ZHANG Yong-zhong, FAN Jin-ping, WEI Xiao-jia, HAO Xing-ming. Study about refining grain of magnesium alloy[J]. Journal of Light Metals, 2004, 32(2): 10?13.

[5] Janin E, Bernd B, Ingo S. Phase field simulations of microstructure evolution during solidification of magnesium-based alloys[C]//Modeling of Casting, Welding and Advanced Solidification Processes-Ⅺ, TMS, 2006: 489?493.

[6] Dahle A K, Lee Y C, Nave M D, Schaffer P L, StJohn D H. Development of the as-cast microstructure in magnesium aluminium alloys[J]. Journal of Light Metals, 2001(1): 61?72.

[7] 胡汉起. 金属凝固原理[M]. 第2版. 北京: 机械工业出版社, 2000: 128?134.
HU H Q. Principle of solidification[M]. 2nd Edition. Beijing: China Machine Press, 2000: 128?134.

[8] Maltais A, Dube D, Fiset M, Laroche G, Turgeon S. Improvements in the metallography of as-cast AZ91 alloy[J]. Materials Characterization, 2004, 52: 103?119.

[9] Nastac L. Numerical modeling of solidification morphologies and segregation patterns in cast dendritic alloys[J]. Acta Mater, 1999, 47(17): 4253.

[10] 梁作俭, 许庆彦, 柳百成, 等. 单相合金凝固过程微观组织的三维数值模拟[J]. 金属学报, 2004, 40(4): 439?444.
LIANG Zuo-jian, XU Qing-yan, LI Jia-rong, et al . Three dimensional numerical simulation of solidification microstructure of single phase alloy[J]. Acta Metallurgica Sin, 2004, 40(4): 439?444.

[11] Thévoz Ph, Desbiolles J L, Rappaz M. Modeling of equiaxed microstructure formation in castings[J]. Metall Trans A, 1989, 20A: 311?322.

[12] Rappaz M, Gandin Ch A. Probabilistic modeling of microstructure formation in solidification processes[J]. Acta Metall, 1993, 2(41): 350.

[13] Nastac L, Stefanescu D M. Stochastic modeling of microstructure formation in solidification processes. Modeling Simul[J]. Mate Sci Eng, 1997(5): 391?420.

[14] LIU Ying, XU Qing-yan, LIU Bai-cheng. A modified cellular automaton method for the modeling of the dendritic morphology of binary alloys[J]. Tsinghua Science and Technology, 2006, 11(3): 153?159.

[15] Bernd B, Janin E, Munekazu O, Gerald K, Martin F, Rainer S F, Ingo S, Andreas B P. Controlling microstructure in magnesium alloys: A combined thermodynamic, experimental and simulation approach[J]. Advanced Engineering Materials, 2006, 4(8): 241?247.

基金项目：国家重点基础研究发展计划资助项目(2006CB605208)；清华－东洋镁合金成形技术研究开发中心项目

收稿日期：2007-03-26；修订日期：2007-07-28

通讯作者：熊守美, 教授, 博士; 电话：010-62773793; E-mail: smxiong@tsinghua.edu.cn

(编辑 陈爱华)