Available online at www.sciencedirect.comsCIENCE dDIRECToTransactions ofNonferrous MetalsSociety of ChinaScienceTrans. Nonferrous Met. Soc. China 16(2006) s591-s597Presswww.csu.edu.cn/ysxb/Simulation for microstructure evolution of Al-Si alloysin solidification processXU Hong(徐宏), HOU Hua(侯華), ZHANG Guo-wei(張國偉)School of Materials Science and Engineering, North University of China, Taiyuan 030051, ChinaReceived 10 April 2006; accepted 25 April 2006Abstract: The numerical simulation for microstructure evolution of AI-Si alloy in solidification process is carried out with phasefield model. The phase field model, solution algorithm and the program of dendrite growth are introduced. The definition of initialcondition, boundary condition and the stability condition of differential format are all included. The simulation results show that theevolution of dendrite morphology is as follows: the initial circle nucleus transforms to the rectangle one firstly, then its cornersdevelop to the four trunks and from which the secondary side branches are generated and even the third side branches are producedfrom secondary ones. The dendrite tip radius decreases quickly at the initial stage and changes slowly at the late stage, which ismainly due to the fact that more and more side branches appear and grow up. The comparisons of dendrite morphology betweensimulated results and investigations by others are also presented. It is proved that the dendrite morphologies are similar in trunks andarms growth, so the developed phase field program is accurate.Key words: phase field model; numerical simulation; microstructure evolution; dendrite morphology; Al-Si alloy(1983) carried out the numerical solution of phase field1 Introductionequations for the first time. In 1993, KOBAYASHI[10]realized the 2-D simulation for the dendrite growth inEarly in 1966, OLDFIELD[1] proposed that theundercooled pure melts by this model coupled withheat source term in heat transfer equation could beanisotropic influence. WHEELER et al[11] quantitativelyrepresented by the function of nucleation rate and growthsimulated the solidification microstructure of pure Nivelocity of crystal grain, and he attempted to simulate thethrough the thermodynamically consistent phase fieldsolidification microstructure of gray casting iron. Yet themodel proposed by WANG et al[12] and PENROSE andmicro simulation had been developed slowly anFIFE[13]. In 1993, WHEELER et al[14] developed theconfined by the corresponding macro simulation duringbinary alloy phase field model, and predicted the solutethe followed decades. Currently, there mainly exist threemicro segregation and the dendrite morphology. In 1996,kinds of models: Deterministic Model, Stochastic ModelKARMA and RAPPEL[15,16] obtained effectiveand Phase- field Model. In this paper, the authors adoptedGIBBS-THOMSON dependence within certain interfacephase-field model, and in this model the phase fieldwidth, and suggested that the interface width could bevariable φ was introduced to describe the morphology,thicker than capillary length, thus developing a newcurvature and movement of phase interface. Coupledhase field model, which can simulate deep undercoolingwith the thermal field, solute field, flow field and othercondition and shorten the calculation time. PROVATASoutside fields, the phase field model can simulate the realet al[17] successfully simulated the single crystal grainsolidification process.morphology through the self-fitted grid technology.Since 1980s, COLLINS and LEVIN[2], LANGERSTEINBACH et al[18] first proposed the phase field[3] and CAGINALP and FIFE[4] proposed the initialconcept of multi phase system and deduced thephase field model. Then, CAGINALP et al[5-7], FIFEgoverning equations. NANILOV and NESTLER[19]and GILL[8, 9] suggested that when the interface widthpresented adaptive finite element simulatins of dendriticis near zero, PFM would be the sharp interface model,and eutectic中國煤化工mary aloys.and they firstly introduced the anisotropic factor. FIX G JTwo-dimensiorYHCNMHGareusedbyCorresponding author: HOU Hua; Tel: +86-1394153099; E-mail: houhua@263 .net.s592XU Hong , et al/Trans. Nonferrous Met. Soc. China 16(2006)REMIREZ and BECKERMANN[20] to test standardphysical state (liquid or solid state) of the system withtheories for free dendritic growth of alloys.respect to time and space. 中= 1 representssolid phase,The phase field method model becomes the research$=0 represents liquid phase， and 中changeshotspot within modeling material science as soon as itcontinuously from 1 to 0 at solid/liquid interface. Fromhas been proposed, so it is of great potential anFig.2, we can easily see that there is a difusion interfaceimportance[21-24].layer between solid and liquid phase， 中changesgradually from a certain constant at one side of the2 Numerical model of phase field modelinterface to another constant at the other side.(PFM)D↑s2.(0)At present, there are such two models for free| r()interface of dendrite growth as sharp interface model and2. (1)diffusion interface model.Sharp interface model is shown in Fig.l, for anclose space 2, the simulated area was divided into solidstates Q (t) and liquid state 2+() by the sharp interfaceSolid phase), across which the system changes suddenly fromsolid state to liquid state and Rt) is taken as no .Liquid phasethickness. So， temperature and material propertiesparameters of the two sides across the interface are notFig.2 Diffusion interface modelcontinuous.There arex∈Q2_()(1)s2.(1)x∈2+(1)(2| r(tBased on GINZBURG-LANDAU free energytheory[3], the phase field equation could be deducedfrom the free energy F or entropy S, which has beencalled free energy function method or entropy functionmethod. So the expression of phase field model isdifferent from various addresses of free energy. In thispaper, the fundamental formula of phase field modelfrom free energy F and entropy S of the system will be| Liquid phase Xpresented respectively.According to GINZBURG-LANDAU free energyFig.1 Sharp interface modeltheory [3], for a enclosed system with volume s2, F is .In diffusion interface model as shown in Fig.2,F=. [|f(I,c,e)+ =e2(V4)2 +-82(Vc)2 |de(3)solid/liquid interface has a certain thickness, and varioustransmission phenomena are completed at this interfacearea. Solid phase, liquid phase and the interface are takenBased on the minimum energy principle， theas different phases represented by a scalar, which canvariation of LYAPOUNOV function and the linearform a group of differential equations together with otherirreversible dynamics, such equations can be deduced asvariables. Therefore, it is not necessary to track interfacefollows: .and separate liquid and solid phase, whereas, use a scalaraφδFto track the phase in system. Thus, the same numericalat=-Mφ 88(4)method can be adopted in the entire simulated zone, and(5)中(x, 1)=1 can raise the calculation precision. Moreover,Dtbe )compared with the sharp interface model, the diffusioninterface is more close to P(x, t)=0, so it has been widelybc_v.M v°(6)Bt中國煤化工applied in various phase transformation processes.Diffusion interface model is adopted in PFM, theAccordingMHCNMH GAu entropyphase field variable中is introduced, which addresses the theory, for the same system as above, the expression of.XU Hong , et al/Trans. Nonferrous Met. Soc. China 16(2006)s593entropy SisHg=LB|T TB(15)S= j{s(4,c,e)--:2(v) +號82(vc}2 d2(7)where LA and Lg are the latent heats of components ASuch equations can also be givenand B, with the melting point Tm4 and TmB respectively.φ=M'$，δ(8)lc_-v.[M.h'(φ)(Hx-H; )v]+Dtdt)c_V .M'V(9)Mc VcSc )vm c(1-c)besS )V .[M。r(,r)vr]- v[M。v(j2v2c)(16)e)Eqn.(16) is the governing equation of concentration.where F and S are respectively the free energy density. It is assumed that solute diffuse hardly in solid phase,function and entropy density function; c is the alloysolute distribution is homogeneous in liquid everywherecomposition; e is the internal energy; ε and δ are theand solute redistribution only takes place between thegradient coefficients of phase field and solute fieldsolid-liquid interfaces.gradient respectively; M$, M'$ are the phase fieldHere,parameters related to the interface dynamics; M。M'are the phase field parameters related to the soluter(,T)=Vr22(u4 -[川)(17)difusion; M，M' are the phase field parametersrelated to the heat transfer.M?！諺(1- c)[D5s-h(>)Dr](18)For a binary alloy, the model is constructed basedon such conditions: 1) Energy and entropy balancewhere Vm is the molar volume, R is the gas constant, Dsequation; 2) The relationship between the thermal-and DL are the solute diffusion cofficients of solid anddynamics driving force and field flux of all kinds ofliquid phase respectively, δ is a gradient coefficient ofexterior field such as temperature field and concentrationphase field.field follows linear laws; 3) Local entropy is positive.It is assumed that the density, special heat andEqn.(11) is the phase field governing equationconductivity coefficient are constants, and equal solidproposed by KARMA[25]:and liquid thermal conductivities are assumed, then thegoverning equation of temperature could be given as(o)中= (r(wo)-xg()) (r"(W(- o)H, +cH]dt~8T= D,V2τ .[4(1-c)+cL2]Jh():(11)c,0twhere t is time, by defining W(0)= =Woas(0) and t(0)=ToQs(0), the interface anisotropy is introduced, Wo is[- ()nrthe interface thickness, to represents the time of atomicmovement between the solid-liquid interfaces, as(0) iswhere Dr is the thermal diffusivity coefficient, and Cp isthe factor of anisotropy and as(0)= I+ycos(k(0, 0)),the special heat.θ=VφVφ denotes the crossing angle betweenEqns.(11), (16) and (19) constitute the phase fielddendrite principal axis and normal direction of interface,model for a single-phase binary alloy.θ;=45° is its initial value, y and k are the exponent andthe modulus of anisotropy respectively. h($) is a3 Numerical methodsfunction constructed for the latent heat released duringthe process of dendrite growth and g() is a function3.1 Initial condition and boundary conditionconstructed for the free energy of the system, h'(中) andIt is assumed that the original nucleus radius is ro;g'($) are the differential coefficients of h($) and g($)the initial conditions are given asrespectively, c is the concentration of the solute, Q is aWhenx+y≤r,balance coefficient.φ=0, T=To, c=Co(20)In Eqn.(11),Whenx2+y2>心,h(D)=φ .(12)g(φ)=φ(1-9)°(13)φ--1, T=-OT,中國煤化工(21)Hg=L(14)where To isDH.C N M H G undercooledT TM)molten melts, ZERO-NEUMANN boundary condition is.s594XU Hong , et al/Trans. Nonferrous Met. Soc. China 16(2006)selected in calculations of phase field and temperaturegrowth code of the aluminum alloy has been developed,field.and post-process program adopted Visual C++ 6.0 werelso carried out to display dendrite growth process3.2 Stability condition of difference equationdynamically. Systematic environmental is P-IV 2.4G,It is inevitable that some truncation errors exist inmemory is 1 G, hard disc is 80 G, and operating system isthe finite difference method, which is stable if the finalWin2000.solution changes lightly with the small fluctuations of theThe flow chart of the dendrite growth program isinitial condition and boundary condition, whereas theshown in Fig.3. .solution is unstable with the contrary condition.Harmonic analysis method is often adopted to discuss the4 Simulation for dendrite growth process ofstability condition: the error in initial value distribution isAl-Si alloyEi,j, and that of next time is 8, j, deduced similarly inturn for the error in p, p+1 time and so on, we can get theThe molten metal is in thermodynamic unstabilityerror increase factor o(k)= 8f/Epj. The differenceduring solidification; nuclei can nucleate and grow dueequation is stable with |o(k)| ≤1, and is unstable withto the crystallization driving force. As the crystal grows|0(k)| > 1. At the same time, Ot must be taken asto a certain size, the solidified interface would becomeunstable under the solute fluctuation and heat disturbance,Ot≤min(22)and primary branches of the dendrite would be formed.(2a/2x2ij.k Ayz2.j.k 021jk )Similarly，the once branches interface would alsobecome unstable during their growth process andwhere 1/Ort,j,k, 1/0yi,jk， 1/0zij,k are the spaceresulted in the secondary branches. Analogously, thesteps of the discrete unit.third branches even higher branches could appear. SoIn this paper, the even grid of space step(Or=Oy)the dendrite could have complicated side branches,is adopted and the time step is Ot (tOr ")/(5:'). Thewhich represent the most complex interface evolutioncoefficient (ε Ot)/(tOx) is maintained to be positive topattern during solidification.ensure the stability.Fig.4 presents the dendrite morphology evolutionduring its growth process. Here, t is the time step. The3.3 Program flow chartblue color and red color represent the liquid metal andAdopting Visual FORTRAN 6.0, 2-D dendritethe solidified part, respectively. It is assumed that thereBeginningAnisotropyRead parametersPhase fieldRead mesh informationSolute fieldSet initial conditionTemperature fieldSet boundary conditionBoundary conditionEnd calculation?TYesOutput graphic ofmicrostructureCalculate critical parameters: dendrite tip growthspeed, radius, twice branch span etc中國煤化工EndMYHCNMH GFig.3 Flow chart of dendrite growth program.XU Hong , et al/Trans. Nonferrous Met. Soc. China 16(2006)s595(a)b)(C)●(d)(e)十++h)美⑥)@Fig.4 Simulated dendrite growth process of Al alloy at different times: (a) t=100 s; (b) t=500 s; (c) t2 000 s; (d) t=3 000 s;(e)t=5 000 s; () =15000 s; (g) =25000 s; (h) =35000 s;(i)t=45 000 s; () =55000中國煤化工has been embedded a circinal nuclei in the liquid matrixAs shownY.HCNMHGstill a circleand its radius is bigger than r'n.nuclei, and at t=500 s, it grows bigger and begins to.s596XU Hong , et al/Trans. Nonferrous Met. Soc. China 16(2006)transform to rectangle. From t=2 000s to t=3 000 s, thebranches are generated and even the third side branchesfour corners of the rectangle grow along their tipre produced from secondary ones. It seems that thedirections and thus form the four dendrite trunks. Tildendrite tip radius decreases quickly at the initial staget= 5 000 s, some lttle secondary side branches appearand changes slowly at the late stage, which is mainly dueand the trunks change thinner, which is because thatto the fact that more and more side branches appear andthe solid phase in the secondary dendrite come from its .grow up. Similarly, the dendrite grows quickly initiallycorresponding trunk. During the following solidificationand gradually changes slowly.process till t=45 000 s, more and more secondarydendrites appear and grow up, at the same time, the 5 Verification of microstructure simulationdendrite trunk become thinner and thinner, the dendritewith investigations by otherstip radius become smaller and smaller. At t= -55000,65 000, 75 000 s, some third side branches begin toFig.5 presents the comparison of microstructureappear and grow up gradually from the secondary sidebetween simulated results and experimental or simulatedbranches.investigations by others with different methods.Therefore, it can be concluded that, during theFrom the above results, we can see that the dendritesolidification process, the initial circle nuclei transformmorphology is similar in trunks and arms growth, so theto the rectangle ones firstly, and then its corners developdeveloped phase field program is accurate. Of course,to the four trunks and from which the secondary sidethis study is still at the stage of theoretical exploration,a)b)Cd)中國煤化工Fig.5 Verification of microstructure simulation results: (a) Ilustration of simulated do: this paper; (b)Simulated dendrite morphology by CA method by Tsinghua[26]; (C) PhotographsfYHCN M H Gerns taken byGLICKSMAN et al[27]; (d) Experimental result of SCN-Ac organic alloys taken by KURZ and ESAKA[28].XU Hong , et al/Trans. Nonferrous Met. Soc. China 16(2006)s597and it is necessary to be linked with engineering[1] WHEELER A, MURRAY B T, SCHAEFER R J. Computation ofdendrites using a phase field model [J]. Phys D, 1993, 66(10):application further.243- -262.[12] WANG s L, SEKERKA R F, WHEELER A A,et al.6 ConclusionsThermodynamically consistent phase-field models for solidification[凹] Phys D, 1993, 69(10): 189- -200.1) The phase field model, solution algorithm and the[13] PENROSE O, FIFE P C. Therrmodynamically consistent models ofphase-field type for the kinetics of phase transitions []. Plhys D,program of dendrite growth are introduced. The1990, 43(10): 4462.definition of initial condition, boundary condition anI4] WHEELER A A, BOETTINGER w J, MCFADDEN G B.Phase-field model for isothermal phase transitions in binary alloys [J].the stability condition of differential format are allPhys Rev A, 1992, 45(0): 7424- 7439.included.15] KARMA A, RAPPEL W J. Phase field method for computationally2) The dendrite morphology evolution: the initialefficient modeling of solidification with artbitrary interface kinetics[]. Phys Rev E, 1996, 53(4): R3017-3020.circle nuclei transform to the rectangle ones firstly, andKARMA A， RAPPEL W J. Numerical simulation ofthen their corners develop to the four trunks and fromthree- dimensional dendritic growth [J]. Phys Rev Lett, 1996, 77(10):which the secondary side branches generated and eventhe third side branches produced from secondary ones.[17] PROVATAS N, GOLDENFIELD N, DANTZIG J A. Adative aridmethods in solidifcation microstructure modeling [A]. THOMAS BThe dendrite tip radius decreases quickly at the initialG, BECKERMAN c. Modeling of Casting, Proceedings of thestage and change slowly at the late stage, which isEighth International Conference on Welding and AdvancedSolidification Processes VIII, Modeling of Casting, Welding andmainly due to the fact that more and more side branchesAdvanced Solidifcation Processes VIII [C]. Califomia, San Diego: Aappear and grow up.Publication of TMS, 1998: 533- -540.3) The comparisons of dendrite morphology[I18] STEINBACH I, PEZZOLLA F, NESTLER B et al. A phase fieldbetween simulated results and investigations by othersconcept for multiphase system [J]. Phys D, 1996, 94(1): 135-147.19] DANILOV D, NESTLER B. Phase-field simulations of soldificationare also presented. It is proved that the dendritein binary and temnary systems using a finite element method [J].morphologies are similar in trunks and arms growth, soJourmal of Crystal Growth, 2005: e177-e182.the developed phase field program is accurate.20] REMIREZ J C, BECKERMANN C. Examination of binary alloyfree dendritic growth theories with a phase field model []. ActaMaterialia, 2005, 53: 1721-1736.References21] ZHAO Y H, CHEN z, WANG Y X, LU Y L. Atomic-scale computersimulation for early precipitation process of Ni75A110V15 alloy [].Progress in Natural Science, 2004, 14(3): 241-246.ASM Trans.1966. 59(2): 945- -9602] ZHAO Y H, HOU H, XU H, et al. Atomic-scale computer simulation[2] COLLINS B, LEVIN H. Difuie intnodel of dffusion-limitedfor early precipitation process of Ni75AI9V16 alloy []. Joural ofcrystal growth [D]. Phys Rev B, 1985, 319;: 6119 -122Material Science and Technology, 2003, 19(S1): 1719- -19.[3] LANGER J s, GRINSTEIN G MAXENKO G Directions in23]ZHAOYH,CENz,LxL.CComputersimulation ofstrain-induced morphological transformation of coherent precipitates[4] CAGINALP G FIFE P c. Phase field methods for interfacial[D]. J Univ Sci Technol Bejjing (M.M.M.), 2003, 10(4): 55 -60.24] ZHAO Y H, CHEN z, WANG Y x, et al. Computer simulation for[5] CAGINALP G FIFE P c. Dynamics of layered interfaces arising .the coarsening behaviour of coherent precipitates [J]. Rare Metalfrom phase boundaries []. SIAM J of Appl Math, 1988, 48(3):Materials and Engineering, 2004, 33(3): 234- -238.(in Chinese)506- -510.25] KARMA A, PAPPEL w J. Qantitative phase-field modelling of[6] CAGINALP G FIFE P C. High-order phase field models anddendritic growth in two and three dimensions []. Phys Rev, 1998,detailed anisotropy []. Phys Rev B, 1986, 34(7): 4940-4943.[7] CAGINALP G et al. Type models as asymptotic limits of the[26] XUQ Y, LIU B C, LIANG Z J, HWANG W-Sing. Procedings ofphase- field equations [U]. Phys Rev A, 1989, 3911): 5887- -5890.Modeling of Casting and Solidification Processes VI, 2004.[8] FIFE P C, GILL G S. The phase-field description of mushy zones [].27] GLICKSMAN M E, et al. Isothermal dendrite growth-a proposedPhys D, 1989, 35: 267- 275.micro gravity experiment []. Metall Trans, 1988, 19A: 1945-1953.[9] FIFE P C，GILL G S. Phase-transition mechanisms for th28] KURZ W, ESAKA H. Solidification microstructures []. Praktischephase-field model under internal heating []. Phys Rev A, 1991,Metallographie, 1988, 25(5): 207-213.43(2): 843-851.(Edited by PENG Chao-qun)[10] KOBAYASHI R. Modeling and numerical simulation of dendritecrystal growth []_ Phy S D, 1993, 63(): 410 -423.中國煤化工YHCNMH G.