A phase-field study on the oxidation behavior of Ni considering heat conduction

Chao Wang·Shigang Ai·Daining Fang，3

?

RESEARCH PAPER

A phase-field study on the oxidation behavior of Ni considering heat conduction

Chao Wang1·Shigang Ai2·Daining Fang1，3

Phase-field modeling approach has been used to study the oxidation behavior of pure Ni when considering heat conduction.In this calculation，the dependence of the coefficient of the Cahn-Hilliard equation Lc on the temperature T was considered.To this end，high-temperature oxidation experiments and phase-field modeling for pure Ni were performed in air under atmospheric pressure at 600，700，and 800?C.The oxidation rate was measured by thermogravimetry and Lc at these temperatures was determined via interactive algorithm.With the Lc-T relationship constructed，oxidation behavior of Ni when considering heat conduction was investigated.The influence of temperature boundaries on the oxidation degree，oxide filmthickness，and specific weight gain were discussed.The phase-field modeling approach proposed in this study willgive some highlights of the oxidation resistance analysis and cooling measures design of thermal protection materials.

High-temperature oxidation·Phase-field approach·Heat conduction·Thermogravimetry analysis· Nickel

? Shigang Ai sgai@pku.edu.cn

1State Key Laboratory for Turbulence and Complex Systems，Peking University，Beijing 100871，China

2Institute of Engineering Mechanics，Beijing Jiaotong University，Beijing 100044，China

3Beijing Institute of Technology，Beijing 100081，China

1　Introduction

High-temperature metals and alloys are widely used as structuralcomponentsin aircraft.Because ofthe drastic aerodynamic heating effect，the surface temperature of aircraft is rather high.As an important property，oxidation resistance has been intensively investigated for high-temperature materials［1-4］.Oxygen diffusion in the oxide film plays a key role for the oxidation of metal，since oxygen permeates the oxide film and contacts the metal substance，and then the oxidation reaction continues and finally damages the metal. Therefore，knowledge of the oxygen diffusion in the oxide film and metal substrate is essential in order to evaluate their integrity［5］.

The technique widely used in studying oxidation resistance of thermal materials is thermogravimetric analysis，in which the specimen’s specific weightchangeswith the oxidation time.Haugsrud［6］comprehensively studied the oxidation mechanism of Ni in the temperature range 500-1400?C with differentoxygen pressure by an experimentalapproach. In that study，the oxidation kinetics had been investigated by mass spectroscopy and thermogravimetry，the oxide scales had been characterized by atomic force microscopy（AFM）and scanning electron microscope（SEM）.In addition，in Zhou and Shen’s work［7］，the dynamic oxidation kinetics of pure Ni was conducted in air at 650 and 850?C.Then too，Jin et al.［8］studied isothermal and cyclic oxidation behaviors of pure Ni at 1000?C in air，SEM and transmission electron microscopy（TEM）were used to examine the oxide scales formed on Ni substrate.They found that Y-implantation greatly improved the anti-oxidation ability of Ni both in isothermal and cyclic oxidizing experiments.In Chevalier et al.’s work［9］，high temperature oxidation of Ni was investigated in air under atmospheric pressure in thetemperature range 600-900 C.It was found that oxidation kinetic curves deviate from the parabolic law for temperaturesover800?C.Forexample，Korshunov etal.［10］studied the effect of preliminary severe plastic deformation on structure and durability of Niin oxidation.Furthurstudy was done by Zhou et al.［11］who studied the growth of oxide scales on pure Ni oxidized under tensile and compressive stress at 973 K，though thermogravimetric analysis found that oxidation rate of pure Ni was accelerated by the external stress，especially by compressive stress.

The oxidation behavior of Ni has been reported in numerous articles，whereas，the oxidation mechanism of Ni is complex and the mechanism of oxygen ingress is not fully clear［6］.Generally，experimental thermogravimetric analysis is time consuming，and for some materials their service temperature is too high to be reached in experiments.Therefore，some scholars developed kinetics models to study oxidation.In Suo et al.’s work［12-14］，an analysis model to characterize the residual stress evolutions during an isothermal oxidation process was developed.Based on this model，effects of creep constants，substrate thickness，and intrinsic strain on the residual stress distribution in the oxide scale and metal substrate were discussed.In their work［15］，a coupled diffusion-reaction-mechanics model of the metallic oxidation was developed and analytical solutions of the concentration in steady state were obtained.More study was completed in work by Dong et al.［16-18］，stresses in the oxide film/metal substrate system at high temperature were theoretically analyzed.In work by Ruan et al.［19，20］，an analysis model to elucidate the residual stress evolutions in an oxide scale/metal substrate system during the isothermal oxidation process was developed.The residual stress in the oxide scale/metal substrate system due to oxidation growth strain and creep deformation were discussed.In work by Chatterjee et al.［21］，the isothermal and cyclic oxidation for Ni-base alloys wasdetermined through a modified Wagner’s oxidation model and the solution of coupled elemental diffusion equations.In Zhou’s work［22］，a continuum thermodynamic model was constructed and the stress-diffusion interaction in the oxidation of Cr-Fe alloys was studied.

The phase-field approach is an effective tool to demonstrate phase change of materials［23-26］.Through this approach，one can get a deep understanding on complicated mechanical-oxidative coupling effects by considering phase evolution and mechanicalmodels together.In Zuo and Zhao’s work［27］，the phase-field approach was adopted to investigate the diffusing of Li-ions and successfully applied in lithium ion batteries.In the work by Ma et al.［28］，a phasefield modelwasdeveloped to analyze the oxidation behaviors of multiphase materials.The evolutions of the porosity and the oxidation stress for the oxidized ZrB2/SiC ceramics with different temperatures and different oxygen partial pressures were predicted.In addition，Yang et al.［29，30］studied the oxidative growth stress and the localized oxidation caused by cracks of iron-based alloy using the phase-field method. In the work by Yang et al.［31］，a phase-field method based on the theory of thermodynamics was developed to simulate the oxidation behavior and oxidation-induced growth stress. With that approach，oxygen induced crack propagation in Ni alloy was discussed.

Generally，the above phase-field oxidation modeling was conducted at a certain temperature，e.g.，700，900?C etc.，but no work has been reported considering the nature of heat conduction，whereas，thermal materials are always used in a thermaldynamic environment.So，heatconduction should be included in phase-field modeling.It should not only include the variation of material parameters at different temperature，butalso the effectoftemperature on the Cahn-Hilliard equation.In the current work，firstly，high temperature oxidation experiments and phase-field modeling for pure Ni was performed in air under atmospheric pressure at 600，700，and 800?C.The oxidation rate was measured by thermogravimetric analysis and the Cahn-Hilliard equation coefficient Lc at these temperatures were determined by iterative calculation.Secondly，based on the phase-field modeling and the experimental results，the L c-T relationship was constructed.Finally，oxidation behavior of Ni considering heat conduction wasinvestigated and the influence oftemperature boundarieson the oxidation degree，oxide film thickness，and specific weight gain were discussed.

2　Experiment

2.1Material

The same batch of pure Ni（99.98%）was used in the present study with the size ofφ10.0 mm×5.0 mm.To investigate the effectoftemperature on the oxidation behaviorofNi，the oxidation experiments were conducted in air under atmospheric pressure at 600，700，and 800?C.For each experiment，weight gains of the Ni samples were recorded at the oxidation times of 0.5，1.0，2.0，4.0，8.0，and 16.0 h.For each temperature，24 samples were tested to plot fully specific weight gain-oxidation time curves，72 Ni cylinder samples were prepared and tested.

2.2Uniaxial tensile test

As the materialparameters of Ni，such as the elastic modulus and Poisson’s ratio，are varying with the annealing process of the material，the mechanical parameters of the Ni material used in this study were tested firstly.Unidirectional tensile experiments were performed on an Instron test system，and four samples were tested.The nominal stress-strain curve isillustrated in Fig.1.By the uniaxial tensile experiments，the elastic modulus and Poisson’s ratio of the Ni metal used in the present study is 199.8 GPa and 0.31.

Fig.1 Nominal stress-strain curve of Ni in uniaxial tension

2.3Measurements of oxidation

The oxidation rate was measured by thermogravimetry with a microbalance（weight measurement accuracy is 0.01 mg）. Prior to oxidation，the Ni specimens were polished with SiC-paper with different grades of roughness and then ultrasonically cleaned in anhydrous alcohol.After drying by a hair dryer，the specimens were introduced into a muffle for oxidation.The Ni specimens were placed in several Al2O3crucibles and were introduced into the muffle with a crucible clamp.In this way，the desired temperature can be established before the specimens were exposed to the oxidizing environment.When all of the specimens were introduced in the muffle，the timing began.Every time when the oxidizing time reached the setting value（0.5，1.0，2.0，4.0，8.0，and 16.0 h），a crucible retracted out of the muffle with the crucible clamp.As soon as the specimens cooled to room temperature，weights of the oxide specimens were recorded.

3　Modeling approach

3.1Phase-field approach

The phase-field approach based on the Landau theory is a powerfulapproach in materialsphase transformation simulations.In the modeling ofoxidation，the oxygen concentration c is taken as the field variable.In the process of oxidation，c（i，j，t）isthe degree ofoxidation forthe i，j materialpointat t time in a two-dimensional（2-D）setting.ci，j，tcan be taken as the field variable for the phase-field equation，and it obeys the Cahn-Hilliard equation［32］

where Eelasticrepresents the elastic energy density，F(xiàn) is the free energy density，t is the oxidation time，and Lc is the coefficient of the Cahn-Hilliard equation.The free energy can be calculated by the following equation

where α0，α1，and α2are chemical-free energy density coefficients；kcis the coefficient of gradient energy density.

The elastic energy can be given as

In this model，for the oxidation under constant temperatures cases， is the free-stress strain and εijis the real strain. εfree

ijis due to oxidation strain，and it can be written as

While，for the oxidation in thermal conducting cases， is the combined action of oxidation strain and thermal expansion，which can be written as

Making use of Eqs.（2）-（4）in Eq.（1）yields the Cahn-Hilliard equation as

where ωkis the relaxation factor，γkis the Pilling-Bedworth ratio（PBR），and α is the coefficient of thermal expansion. Forthe NiOstudied in this work，ωkandγkare 0.18 and 1.65，respectively［33］；α2=-1.0×108J/m3and kc=5.0× 107J/m3.The phase-field equation（Eq.5）was calculated by the fast Fourier transformation（FFT）approach，which was detailed in another article of Yang et al.［31］.With all the obtained driving forces，iterations are processed in time space by the finite difference method［29-31］.

Fig.2 The computational model for Ni oxidation

3.2Computational model

Figure 2 illustrates the computational model for the bulk Ni metal oxidation.MN is the surface of the Ni sample and it contacts with external oxygen.In the phase-field modeling for pure Ni，c is defined as the degree of oxidation for Ni material.c=1.0 indicates completely oxidized while c= 0.0 presents completely unoxidized.The thickness of the NiO is along the Y-axis and the oxygen diffusion direction is also in the Y-direction.For the 2-D computational model，the X-axis is infinite，so，periodic boundary conditions were applied on edges MM′and NN′in the X-direction.In the Y-direction，the computational model is a half of a Ni plate，therefore，symmetric boundary conditions were applied on the M′N′edge.The whole mode is divided into 50×50 elements.Each element is square with the element size of 1.0μm×1.0μm.

3.3Heat conduction

In this study the oxidation behaviorof Niis investigated with consideration of heat conduction.For a 2-D setting，the heat transfer equation is as below

At the top edge of the computational model MN，the thermal boundary is: where ù is the heating rate.Two types of thermal boundary conditions，the firstthermalboundary and the second thermal boundary are considered at the bottom edge of the computational modelThe two thermal boundaries can be expressed as Eq.（8）.

Fig.3 Ni samples before（white）and after（black）oxidation at 700?C in air under atmospheric pressure

The temperature at each integral point can be calculated by the central difference approach as where α= τλ/h2，λ=K/（Cρ），τ is the time increment，h is the mesh length，K and C are the conductivity and specific heat of the Ni-NiO composite，respectively.are the conductivity，thermal expansion coefficient，molar mass，and specific heat of Ni，KNiO，αNiO，MNiO，and CNiOare the conductivity，thermalexpansion coefficient，molarmass，and specific heatofNiO.The materialparameters for Niand NiO are listed in Table 1，and both the Ni and NiO are assumed as isotropic.

4　Results and discussion

4.1Oxidation kinetics

Ni samples were oxidized in air under atmospheric pressure at 600，700，and 800?C.Specific weight gain of the Ni specimens was recorded at 0.5，1.0，2.0，4.0，8.0，and 16.0 h.Ni specimens before and after were oxidized at 700?C and are shown in Fig.3，and it is clear that the samples went from silverto black afteroxidation.The variation ofmass gain was measured using a microbalance（weight measurement accuracy is 0.01 mg）.As a general rule，the oxidation kinetics of each specimen is determined by using the thermogravimetric analysis.mgain=△m/S，where the value of mgaindenotes the specific weight gain，△m is the mass gain of measured specimens after oxidation，and S is the surface area of NiO. The oxidation kinetics curves ofpure Niare plotted in Fig.4，illustrating the specific weightgain versustime relationships.

Fig.4 Specific weight gain as a function of time for pure Ni oxidized at 600，700，and 800?C

Fig.5 Specific weight gain as a function of square root of oxidation time for pure Ni oxidized at 600，700，and 800?C

On the basisofthe experimentalresultsin Fig.4，generally，the oxidation kinetics follow the linear-parabolic growth law with the incrementofoxidation time at600，700，and 800?C. The specific weight gain increases sharply at first and then the rate decreases with oxidation time.It is probably due to the NiO scale that acts as a protective layer on the metal substance，which reduces the process of oxidation and diffusion by external oxygen.The corresponding parabolic plots are illustrated in Fig.5.The oxidation rate constant Kpis equal to 1.15×10-6mg2·cm-4·s-1at 600?C.When the temperature is 700 and 800?C，Kpis 2.80×10-6mg2·cm-4·s-1and 1.24×10-5mg2·cm-4·s-1，respectively.Itis noted that the oxidation rate of pure Ni increases with the increment of temperature.

4.2LC-T relationship

Based on the material parameters listed in Table 1，firstly，phase-field simulations were performed at 600，700，and 800?C in air（the molar density of oxygen is 9.375× 10-3mol/L）with various values of Lc（1.0×10-10，1.0× 10-9，1.0×10-8，…，1.0×10-4m5/J·s，...）.Special attention was paid to the specific weight gain of the virtual specimen.Then，the specific weight gain calculated by the phase-field simulations was compared with the experimental results.Ifthe specific weightgainby the phase-field modeling islargerthan the experimentalresultthen a smaller Lc willbe assigned in the nextiteration，and，on the contrary，a larger Lc will adopted.At each temperature，the phase-field calculation will be terminated when the thermogravimetric analysis errors between the phase-field modeling and the experiments are less than 10%（at all the recording times）.

For the 2-D computational model studied in the present paper，the specific weight gain can be calculated by

where ρ and MNiis the density and molar mass of the Ni material，MOis the molar mass of the oxygen atom，Ai，jis the area of the element（i，j），ci，j，tis the oxidation degree of the integralpoint（i，j）attime T，NX and NY are the node numbers of the computational model in the X-and Y-directions，respectively.

By the phase-field modeling，Lc at600，700，and 800?Cis 0.77×10-6，2.1×10-5，and 7.1×10-5m5/J·s，respectively. The comparison of the specific weight gain between phase-field modeling and experiments is illustrated in Fig.6.The errors between the oxidation experimental results and the phase-field simulation are 2.91%，3.24%，2.62%，3.77%，4.45%，and 4.08%at 600?C when oxidation time is 0.5，1.0，2.0，4.0，8.0，and 16.0 h，respectively.They are -3.87%，-6.12%，-4.06%，-3.88%，3.82%，and 9.20% when the oxidation temperature is 700?C.For the phasefield modeling at800?C，the calculation errors are-3.75%，-4.65%，-4.39%，-3.91%，-0.75%，and 8.84%.

Table 1 Materials parameters for Ni and NiO［19，20，34］

Fig.6 Comparison of the specific weight gain between phase-field modeling and experiments

Fig.9 The relationship between phase-field coefficient and temperature forthe oxidation ofpure Nimetalin airunder atmospheric pressure

The contour plot of the degree of oxidation of Ni in air underatmospheric pressure for16 h at600，700，and 800?C is illustrated in Fig.7.In which，the legend value isthe degree of oxidation ci，j，tofNimaterial.ci，j，t=1 presents Nioxidized completely and ci，j，t=0 presents Nicompletely unoxidized. From Fig.7，the thicknesses of the NiO film，approximately，are 6，9，and 14μm when the oxidation temperaturesare 600，700，and 800?C.The microstructures ofscales on specimens at different oxidation temperatures were obtained by examining cross-sections，as shown in Fig.8.By the SEM data，when the oxidation temperatures are 600，700，and 800?C，the thicknesses of the NiO film are 5.7，9.5，and 12.9μm.It can be found that the phase-field modeling agrees well with oxidation experiments.

Based on the results of phase-field simulations and oxidation experiments，the relationship of the phase-field coefficient Lc with temperature T can be expressed as Eq.11.The relationship between Lc and T for the oxidation of pure Niin air under atmospheric pressure is illustrated in Fig.9

Fig.7（Color online）The degree of oxidation for Ni in air under atmospheric pressure for 16 h at a 600?C，b 700?C，and c 800?C

Fig.8 Cross-sections of the NiO film on pure Ni oxidized in air under atmospheric pressure for 16 h at a 600?C，b 700?C，and c 800?C

Fig.10（Color online）The degree of oxidation for pure Ni in air for 16 h with the top edge heated to 800?C and the bottom edge is set as the a first thermal boundary and b the second thermal boundary

Fig.11（Color online）The degree of oxidation for pure Ni in air for 16 h with the top edge heated to 1000?C and the bottom edge is set as the a first thermal boundary and b the second thermal boundary

Fig.12 Comparison on the specific weightgain between two different thermal boundaries.B1 indicates the first thermal boundary and B2 indicates the second thermal boundary

where，a and b are two constants，for the oxidation of pure Ni in air under atmospheric pressure，a=-7.7725×10-4and b=1.01139.

4.3Oxidation behavior of Ni with considering thermal conduction

Thermal materials are working in high-temperature environments facing very serious oxidizing.In order to protect the materials and structures via cooling measures，e.g.，water cooling，air cooling，etc.There are adopted the means to take heat away.For the thermal structures which have cooling measures，the temperature field of the materials is highly heterogeneous.Generally speaking，it has several hot ends and cold ends with heat conducting from the hot ends to the low temperature ends.Based on the L c-T relationship established in Sect.4.2，the phase-field approach can be used to calculate the oxidation behavior of Ni material which has cooling measures.It will give some help for the service performance evaluation of thermal materials and the optimization of the cooling systems for thermal structures.

In this section the oxidation behavior of Ni considering thermalconduction was studied by the phase-field modeling. For the computational model shown in Fig.2，the top edge was the hot end and the thermal boundary was adopted by Eq.7 with˙u=0.5?C/s.When the temperature of the hot end reach the threshold，e.g.，800，1000?C，it will be held constant.The bottom edge was the cold end and two types of thermal boundaries，the first thermal boundary and the second thermal boundary as expressed in Eqs.7 and 8 were adopted.In the two thermalboundariesthe initialtemperature of the virtual model both were 100?C and the oxidation time both were 16 h.When the temperature threshold value of the hot end is 800?C，the degree of oxidation for the two Ni virtual samples is shown in Fig.10.Figure 11 illustrates the oxidation degree of the two Ni virtual samples with the hot end heated to 1000?C.

From Fig.10，the thicknesses ofthe NiOfilm are approximately 6 and 14μm when the first thermal boundary and the second thermal boundary are adopted on the bottom edge of the computational model.It is clear that by adopting the second thermal boundary the oxidation in the material is more serious than that when the first thermal boundary is adopted. It is because for the first thermal boundary the temperature of the cold end is a constant value；it means that heat in the material will be taken away at the cold end.Whereas，for the second thermal boundary?u/?t=0 heat will accumulate in the material，and it will accelerate the oxidation of materials.

The relationship curves between specific weight gain and oxidation time for pure Ni in air under atmospheric pressure considering different thermal boundaries are illustrated in Fig.12.From Fig.12，it is clear that when the first thermal boundary is considered the specific weight gain is smaller. It also found that the specific weight gain-oxidation time curves have three stages.For both the two types of thermal boundaries，in the first stage the specific weight gain does not increase with time，while in the second stage it rapidly increases with the increment of time.In the third stage，for the second thermal boundary the specific weight gain approximatively follows a parabolic law with oxidation time，whereas for the first thermal boundary the specific weight gain-oxidation time curve almost begins to platform. When the top edge was heated to 800?C and the oxidation time was equal to 16 h，the specific weight gains are 0.932 and 0.509 mg/cm2for the two types of thermal boundaries. When the first thermal boundary was adopted，the specific weight gain has a 45.4%reduction.When the top edge was heated to 1000?C，at16 h，the specific weightgains are 3.037 and 0.940 mg/cm2for the two types of thermal boundaries. When the first thermal boundary was adopted，the specific weight gain has a 69.0%reduction.It means that for thermal materials，when the hot end temperature is higher，the cooling measures are more meaningful.

Fig.13 The specific weight gain-oxidation time curves for the oxidation of pure Ni with different temperatures at the cooling end

Fig.14 The specific weight gain-cold end temperature curves for the oxidation of pure Ni at different oxidation times

The oxidation behaviorofpure Niin airwith differentcold end temperatures was studied to investigate the influence of the temperature of the cooling medium.In the phase-field modeling the top edge ofthe computationalmodelwasheated to 1000?C by Eq.7 with˙u=0.5?C/s and then held，and the cold end of the computational model was set to 400，200，100，and 50?C.The specific weight gain of the Ni samples was studied and the specific weight gain-oxidation time curves with different cold end temperatures are illustrated in Fig.13.The specific weightgain-cold end temperature curves at different oxidation times are shown in Fig.14.As shown in Fig.13，the specific weight gain first rapidly increases with oxidation time.About 1 h later，the growth rate of the specific weight gain gradually reduced with the increment of oxidation time.From Fig.14，it is clear that the specific weight gain is lineally increased with the increment of the temperature of the cold end.

5　Conclusions

In this paper，high-temperature isothermal oxidation behaviors of pure Ni were studied based on oxidation experiments and the phase-field modeling approach.Effects of temperatures on the specific weight gain and the Ni oxide film thickness were investigated.Based on the experiments and phase-field simulations，oxidation kinetics of Ni metal in air under atmospheric pressure at different temperatures were studied，and the L c-T relationship was constructed.Based on the L c-T relationship and the phase-field modeling approach，the oxidation behaviorofthe pure Nimaterialwith consideration of thermal conduction was studied.Through the phase-field simulations，some conclusions are drawn.

（1）For the oxidation of pure Ni in air under atmospheric pressure，the relationship between the coefficient of the Cahn-Hilliard equation and the temperature can be expressed as L c=a t（1-bT）and Lc follows the power exponential growth with the increment of temperature.

（2）An active cooling measure（the first thermal boundary adopted）will significantly reduce the oxidation of thermal materials at high temperature.For the thermal structures，the hotend temperature ishigher，and the cooling measures are more meaningful.

（3）For pure Ni oxidation the specific weight gain is initially rapid and then decreases to become parabolic with oxidation time at different temperatures.The oxidation rate of pure Ni increases with an increment of temperature.

Acknowledgments The projectwassupported by the Beijing Jiaotong University（Grant C15JB00080）.

1.Zhang，H.，Lopez，H.F.:Role of nanoceria coatings on the high temperature oxidation of stainless steels:an atomistic simulation. Comput.Mater.Sci.79，634-638（2013）

2.Ai，S.G.，F(xiàn)u，H.L.，He，R.J.，et al.:Multi-scale modeling of thermal expansion coefficients of C/C composites at high temperature. Mater.Des.82，181-188（2015）

3.Panicaud，B.，Grosseau-Poussard，J.L.，Dinhut，J.F.:General approach on the growth strain versusviscoplastic relaxation during oxidation of metals.Comput.Mater.Sci.42，286-294（2008）

4.Zumpicchiat，G.，Pascal，S.，Tupin，M.，et al.:Finite element modelling of the oxidation kinetics of zircaloy-4 with a controlled metal-oxide interface and the influence of growth stress.Corros. Sci.100，209-221（2015）

5.Heo，T.W.，Chen，L.Q.，Wood，B.C.:Phase-field modeling of diffusionalphase behaviorsofsolid surfaces:a case study ofphaseseparating LiXFePO4electrode particles.Comput.Mater.Sci.108，323-332（2015）

6.Haugsrud，R.:On the high-temperature oxidation ofnickel.Corros. Sci.45，211-235（2003）

7.Zhou，X.W.，Shen，Y.F.:A comparative study of pure nickel and the Ni-CeO2nanocrystalline coatings:microstructural evolution， oxidation behavior，and thermodynamic stability.J.Mater.Sci.49，3755-3774（2014）

8.Jin，H.M.，Zhang，J.F.，Yan，K.，et al.:Rare earth effects on high temperature oxidation of pure nickel at 1000?C.Process Nat.Sci. 14，373-376（2004）

9.Chevalier，S.，Desserrey，F(xiàn).，Larpin，J.P.:Oxygen transport during the high temperature oxidation of pure nickel.Oxid.Met.64，219-234（2005）

10.Korshunov，L.G.，Chernenko，N.L.:Effect of preliminary severe plastic deformation on structure and durability of nickel subjected to oxidation.Plast.Metal Metallogr.116，512-518（2015）

11.Zhou，C.H.，Ma，H.T.，Wang，L.:Comparative study of oxidation kinetics for pure nickel oxidized under tensile and compressive stress.Corros.Sci.52，210-215（2010）

12.Suo，Y.H.，Yang，X.X.，Shen，S.P.:Residual stress analysis due to chemomechanical coupled effect，intrinsic strain and creep deformation during oxidation.Oxid.Met.84，413-427（2015）

13.Suo，Y.H.，Yang，X.X.，Shen，S.P.:Reaction-diffusion-stress coupling effect in inelastic oxide scale during oxidation.Oxid.Met. 83，507-519（2015）

14.Suo，Y.H.，Shen，S.P.:General approach on chemistry and stress coupling effects during oxidation.J.Appl.Phys.114，164905（2013）

15.Suo，Y.H.，Shen，S.P.:Coupling diffusion-reaction-mechanics model for oxidation.Acta Mech.226，3375-3386（2015）

16.Dong，X.L.，F(xiàn)eng，X.，Hwang，K.C.:Oxidation stressevolution and relaxation ofoxide film/metalsubstrate system.J.Appl.Phys.112，023502（2012）

17.Dong，X.L.，F(xiàn)ang，X.F.，F(xiàn)eng，X.，et al.:Diffusion and stress coupling effect during oxidation at high temperature.J.Am.Ceram. Soc.96，44-46（2013）

18.Dong，X.L.，F(xiàn)eng，X.，Hwang，K.C.:Stress-diffusion interaction during oxidation at high temperature.Chem.Phys.Lett.614，95-98（2014）

19.Ruan，J.L.，Pei，Y.M.，F(xiàn)ang，D.N.:Residual stress analysis in the oxide scale/metal substrate system due to oxidation growth strain and creep deformation.Acta Mech.223，2597-2607（2012）

20.Ruan，J.L.，Pei，Y.M.，F(xiàn)ang，D.N.:On the elastic and creep stress analysis modeling in the oxide scale/metal substrate system due to oxidation growth strain.Corros.Sci.66，315-323（2013）

21.Chatterjee，A.，Srikanth，S.，Sanyal，S.，et al.:Kinetic modeling of high temperature oxidation of Ni-base alloys.Comput.Mater.Sci. 50，811-819（2011）

22.Zhou，H.G.，Qu，J.M.，Cherkaoui，M.:Finite element analysis of oxidation induced metal depletion at oxide-metal interface.Comput.Mater.Sci.48，842-847（2010）

23.Hong，L.，Liang，L.，Bhattacharyya，S.，et al.:Anisotropic Li intercalation in a LixFePO4nano-particle:a spectral smoothed boundary phase-field model.Phys.Chem.Chem.Phys.18，9537-9543（2016）

24.Momeni，K.，Levitas，V.I.:A phase-field approach to solid-solid phase transformations via intermediate interfacial phases under stress tensor.Int.J.Solids Struct.71，39-56（2015）

25.Deckelnick，K.，Elliott，C.M.，Styles，V.:Double obstacle phase field approach to an inverse problem for a discontinuous diffusion coefficient.Inverse Prob.32，045008（2016）

26.Wu，X.H.，Wang，G.，Zhao，L.Z.，et al.:Phase field simulation of dendrite growth in binary Ni-Cu alloy under the applied temperature gradient.Comput.Mater.Sci.117，286-293（2016）

27.Zuo，P.，Zhao，Y.P.:A phase field model coupling lithium diffusion and stress evolution with crack propagation and application in lithium ion batteries.Phys.Chem.Chem.Phys.17，287-297（2015）

28.Ma，Y.J.，Yao，X.F.，Hao，W.F.，et al.:Oxidation mechanism of ZrB2/SiC ceramics based on phase-field model.Compos.Sci. Technol.72，1196-1202（2012）

29.Yang，F(xiàn).，Liu，B.，F(xiàn)ang，D.N.:Interplay between fracture and diffusion behaviors:modeling and phase field computation.Comput. Mater.Sci.50，2554-2560（2011）

30.Yang，F(xiàn).，F(xiàn)ang，D.N.，Liu，B.:A theoretical model and phase field simulation on the evolution of interface roughness in the oxidation process.Modeling and simulation in materials science and engineering 20，015001（2012）

31.Yang，F(xiàn).，Liu，B.，F(xiàn)ang，D.N.:Analysis on high-temperature oxidation and growth stressofiron-based alloy using phase field method. Appl.Math.Mech.32，757-764（2011）

32.Rokkam，S.，El-Azab，A.，Millett，P.，et al.:Phase field modeling of void nucleation and growth in irradiated metals.Model.Simul. Mater.Sci.17，064002（2009）

33.Huntz，A.M.:Stresses in NiO，Cr2O3 and Al2O3 oxide scales. Mater.Sci.Eng.A 201，211-228（1995）

34.Huang，B.Y.，Li，C.G.，Shi，L.K.，et al.:China Materials Engineering Canon:Inorganic Non-ferrous Metal Materials Engineering，358-361.Chemical Industry Press，Beijing（2006）.［in Chinese］

15 March 2016/Revised:17 May 2016/Accepted:24 May 2016/Published online:11 August 2016

?The Chinese Society of Theoretical and Applied Mechanics；Institute of Mechanics，Chinese Academy of Sciences and Springer-Verlag Berlin Heidelberg 2016