Numerical Modelling of Multilayered Coatings – Latest Developments and Applications

– The remarkable mechanical properties of multilayered coatings, such as super-hardness, excellent resistance against cracking, low wear-rate and high thermal-stability, are due to their unique interfacial structures and deformation mechanisms at the nanometer scale. The multilayered coating process itself is a typical multi-scale phenomenon. The modelling of multilayered coatings has become an important topic in research recently, largely due to the recent progress that has been made in the numerical modelling of materials and structures at different length-scales as well as the improved effectiveness achieved in linking such progress in numerical modelling to enable multi-scale modelling. In this paper, numerical modelling for the analysis of multilayered coatings at individual length-scales: Continuum, Molecular and Nano-scale, is reviewed, along with that at multi-scale modelling. Examples are presented showing numerical models obtained using: the Finite Element Method (FEM), Molecular Dynamics (MD), First-principles calculations and Multi-scale modelling, are presented. Their relative limitations are discussed and challenges to their future work highlighted.


Introduction
Multilayered coatings, especially with transition metal nitrides at nanoscale, are currently subjects to intensive research because of their potential for a wide range of technological applications [1,2] and also because they are, fundamentally, of scientific importance [3,4].Figure 1 shows the representative structures of multilayered coatings (a) TiN/CrN and (b) TiAlN/VN.Indeed, the interest in nitride coatings is twofold: transition metal nitrides, which are well classified as refractory compounds, exhibit unusual physical and mechanical properties such as extreme hardness, high thermal-and chemical-stability, good corrosion resistance, and low, stable friction under even dry contact conditions [5][6][7][8][9].These intriguing properties render them applicable for many uses, for instance, in relation to: engineering tools, magnetic and electric components, superconducting devices, etc., whilst periodic coating of these nitrides in the nano-scale dimension [10,11] onto such tools components and devices can enhance most, if not all, of their functions to a level that none of their individual constituents can ever match [12].A representative example is that the multilayered TiN/VN superlattice, produced by alternating depositions of its two host mononitrides, which reaches a hardness value as high as 56 GPa [13][14][15], this being much higher than that of either of its host coating alone.Moreover, multilayered coatings enable the design of constituted materials with desired properties.For example, of all the transition metal nitrides, TiN is the most commonly used nitride.It is of interest largely because TiN is hard and wear-resistant, has a bright golden color, and can be coated onto many engineering materials.However although TiN/VN multilayered coatings exhibit the highest hardness, their widespread applications are currently limited by their propensity to be oxidized, especially at a high temperature, which degrades mechanical properties of TiN/VN multilayered coatings [16][17][18][19].To avoid this problem, multilayered coatings TiN/CrN or TiN/AlN lead to improved thermal stability and excellent oxidation resistance with no notable loss in hardness [20][21][22][23][24][25][26].Transmission electron microscopy (TEM) studies [27] have attributed the property improvement to a thin Cr 2 O 3 layer that can suppress the formation of oxides and voids.
Apart from the sophisticated applications in tooling industry and other engineered components, due to their excellent mechanical properties, especially the hardness, origin from

OPEN ACCESS REVIEW ARTICLE
the nature of transition metal nitrides, multilayered coatings, composed of other functional constituents, are also used in the fields such as optical filters [31][32][33][34][35][36][37][38][39][40], magnetic storage [41,42], sensors and microelectronics [43][44][45][46].For instances, Wang et al. [32] fabricated a TiO 2 /SiO 2 multilayered optical filters with graded refractive index profiles, deposited on glass (BK7) and Si(1 0 0) substrates simultaneously and sequentially, by helicon plasma sputtering at room temperature, and found that the measured transmittance spectrum exhibited a reflectance of 99.8% at a central wavelength of 730 nm and high transmittance over the wavelength region outside of the reflected band, as a result of the suppression of the sidelobes.Furthermore, Chhajed et al. [35] designed and produced a nanostructured multilayered tailored-refractive-index antireflection coatings on a glass substrate by using low-refractive-index silica, and reported that the measured average optical transmittance between 1 000 and 2 000 nm is improved from 92.6 to 99.3% at normal incidence.Moreover, Raut et al. [39] described different types of the Anti-Reflective Coatings (ARCs) in greater detail and illustrated the state-of-the-art fabrication techniques. Louis et al. [40] also gave an overview on the progress in the thin film and surface physics involved in multilayered systems with nanometer scale periodicity, wherein a special attention was given to the development of thin diffusion barrier layers between the materials in the multilayered coatings to enhance the optical contrast and to reduce the inter-diffusion.
From the above brief description, it can be seen that the multilayered coatings could be tailored/designed based on the choice of individual constituents for a particular multilayered systems, and these intriguing advantages make them applicable for many uses.
Improvements in mechanical properties due to introducing multilayered coatings, encompassing post-processing behavior, have attracted many studies from both experimental and theoretical.Experimentally, the studies have covered the effects of: (i) preliminary treatment of the substrates; (ii) various kinds of material combination, as well as the thickness ratio of neighboring layers; (iii) post-processing, and evaluation of the mechanical, thermal, and magnetic properties via bending, buckling, nano-indentation, scratching, tensile test, acidic corrosion, high-temperature test, etc. [47,48].For instance, understanding of the origin of property improvement through multilayered coatings is twofold: on one hand, it is attributed to the interfaces in multilayered coatings.Nordin et al. [49] reported a significant slowdown of corrosion for the TiN/CrN multilayered coatings because the oxidation is found to be almost inversely proportional to the number of interfaces in the multilayered coatings, as is confirmed in other previous works [50,51].On the other hand, it is ascribed to the residual stress present in the multilayered coatings [52].Mendibide et al. [28][29][30] successfully deposited TiN/CrN onto the tool steel and investigated the effect of residual stress [53], finding, by transmission electron microscopy (TEM) and synchrotron analysis, that the wear resistance of TiN/CrN is improved due to the mode shift of crack propagation induced by the fluctuating residual stress field across the multilayered coatings.In addition, they claimed that interfaces in the multilayered coatings also play a crucial role in improving the refractory properties and wear resistance.More detailed information of multilayered coatings could be found in the literatures, such as references [54][55][56][57].Correspondingly, various theoretical models have also been developed to explore the mechanisms for the enhanced chemical and mechanical properties of multilayered coatings.However, a better understanding of the underlying mechanisms of multilayered coatings still requires a more novel and more applicable model.In this paper, recent advances in the numerical modelling of multilayered coatings are critically reviewed with a reference to modelling at length-scales and the numerical methods used.
Evidently, one of the notable characteristics of multilayered coatings is that the thickness of individual layer is in the nanoscale, that is to say, it is generally from several nanometers to tens of nanometers while the thickness of the whole system may be in micrometers to millimeters, inducing an abundance of numerical simulations at different scales.Based on the scale applied, the existing numerical models can be generally classified into three categories: continuum models, molecular dynamics models, and nano-scale models, as illustrated in Figure 2.There are also some coupled models or methods combining two of them, such as the quasi-continuum (QC) method, which is a mixed continuum and atomistic approach for simulating the mechanical response of polycrystalline materials.Importantly, there are two basic requirements for building a theoretical  model.One is that the model has to be ''big'' enough to represent the behaviours/properties of the multilayered coating system studied, and another is that the size of model is limited by the computational power.Consequently, a theoretical model may only be able to consider a thin layer containing the coatings-substrate system, composed of single, duplex or a limited number of phases.

Continuum models
Intuitively, the structure of super-hard multilayered coatings at the nanometer scale is still a continuum [58].Thus, a natural attempt to model multilayered coatings is based on conventional fracture mechanics and continuum mechanics scaled down to the dimension of nano-layered.A representative volume element (RVE) for a constitutive model is shown in Figure 3.It is worth noting that for RVE, the layers of two phases are parallel to each other, and arranged alternatively with perfect bonds, and with a constant interlayer spacing.The x 1 x 2   plane of the coordinate frame is parallel to the layers, and the x 3 axis is perpendicular to the layers.If the thickness of each layer is assumed to be sufficiently small, compared to its in-plane size, the variation of both stress and strain within a phase can be neglected.The thickness h is the modulation period of multilayered coatings.Consequently, the numerical methods applied in these models include mostly, Finite-Element Method (FEM) [59], Finite Difference Method (FDM) and the Boundary-Element Method (BEM) [60,61].Typical examples of a multilayered coating, simulated mostly with FEM, are presented in this section, containing various elastic and elastic-plastic material properties, normal and tangential loading conditions, a hard coating/functionally graded substrate system, etc.Finally, a FEM model with cohesive zone modeling, used for the prediction of the formation and propagation of micro-cracks under indentation, is taken as an example of a general procedure for the FE simulation of multilayered coatings.FEM can be used to simulate the behavior of multilayered coatings under either (nano-) indentation or subjected to a scratch test, mainly concerning: 1. Stress distribution.For instance, Djabella and Arnell [62] used a FEM to investigate the contact stresses during the process of elastic compression on coating/substrate systems consisting of a high modulus surface coating on a relatively low modulus substrate wherein two pressure distributions were imposed to consider those arising during the Herzian indentation of a homogeneous elastic half-space by spherical and cylindrical indenters.They found that the stresses of both the outer surface of the coating and the coating-substrate interface are complicated functions of (a) the ratio of the Young moduli of the coating and substrate and (b) the ratio of the coating thickness to the contact halfwidth; and that the stresses which could be responsible for coating failure cannot simultaneously be minimized by appropriate choice of variable.Further, Stephens et al. [63] analyzed initial yielding behavior due to the indentation and friction process between an elastic cylindrical surface and a hard coating/functionally graded substrate system where a thin hard DLC film deposited on a soft Ti-6Al-4V alloy substrate is considered as a model system with consideration of two functional gradient substrate conditions: (a) a gradient in yield strength and (b) a gradient in elastic modulus.They found that in both cases, appropriate gradients result in significant benefit to the reliability of the coated system, compared to the case of an ungraded substrate.Furthermore, Bouzakis et al. [64] introduced a continuous FEM simulation of the nano-indentation hardness test in order to describe the coating's elastic and plastic deformation and herewith, and thereby to extract, precisely and independently of the indentation load, the coating stress-strain curve.Similarly, Tan and Shen [65] employed FEM to examine the relationship between indentation hardness and overall yield strength of multilayered materials, and found that the indentation-derived strength consistently underestimates the composite yield strength within the scope of the continuum-based approach.Additionally, Zhao et al.
[66] simulated nano-indentation-induced deformation in TiSiN-based multilayered coatings by means of FEM, and quantified stress distributions under moderate indentation loading in the structure with particular emphasis on the relationship between stress concentrations and crack initiation.The results showed that structural layering can be used to modify the stress distribution and to lower the overall stress level within the coating.In the case of radial tensile stresses at the coating/substrate interface, a reduction of ~50% was achieved through layering.2. Hardness measurement.FEM has been widely used to simulate elastic and plastic deformations beneath a pointed indenter, hardness measurement with a sharp pyramidal indenter, such as Berkovich or Vickers, is one of the most complex mechanical problems, due to the 3D dimensional phenomenon associated with elastic and plastic deformations with large strain fields.As is well known, FEM is used to extract the true hardness of thin coatings, from the composite hardness measurements, based on the assumption that the film and substrate behavior can be modeled with a bilinear law.However, Pelletier et al. [67] showed that the bilinear model does not define correctly the strain hardening, resulting in an overestimation of the yield stress or one need to adjust the target modulus at each indentation depth.

Formation and propagation of cracks. For examples,
Aoki et al. [68] built a FE model for the calculation of the energy release rate during the coating cracking, aiming at showing that SiC/C multilayered coatings successfully suppress through-thickness coating cracks.Furthermore, Chen [69] used a cohesive zone model embedded into a FE model to assess nucleation and propagation of delamination triggered by an indentation in a multilayer stack with relatively weak metal/ceramics interfaces.Holmberg et al. [70] developed a three-dimensional FEM model for calculating the first principal stress distribution in the scratch tester contact of a diamond spherical tip with 200 lm radius sliding with increasing load on a 2 lm thick titanium-nitride coated steel surface, with consideration of the elastic, plastic and fracture behavior of the surface.The results showed that the first crack is initiated at the top of the coating from bending and pulling actions and grows down through the coating also.Other mechanical behavior such as milling, cutting, and sliding have also been investigated with FEM.For instances, O ¨zel et al. [71,72] applied FEM for predicting the forces, temperatures, and wear rate during the micro-milling of Ti-6Al-4V alloy with fine-grain uncoated and cBN coated micro-end mills, the results revealing advantages of cBN coatings.Grzesik et al. [73] studied the temperature distribution in the cutting zone in the turning processes with differently coated tools and showed the tool-chip interfacial influence on the temperature distributions, as a consequence of using coated tools.Moreover, Coelho et al. [74] combined experiment with FE simulation to investigate the behavior of TiAlN, TiAlN-nano-coating and AlCrN coated tools in finishcutting conditions, and showed that the temperature at the tool-chip can be stable and as high as 800 °C when using TiAlN-based coated tool.Additionally, the presence of the coating layer creates an additional thermal resistance between workpiece and the tool, slightly lowering temperature at the substrate.Holmberg et al. [75] carried out a micro-scale finite-element method modelling and stress simulation of a 2 lm TiN-coated steel surface and showed a reduction of the generated tensile buckling stresses in front of the sliding tip, when compressive residual stress of 1 GPa was considered in the model.4. Residual Stress.FEM can also be used to study the residual stress caused by the different thermal expansion coefficients between two phases of multilayered coatings during thermal cyclic loading or a manufacturing process.
For instance, Ng and Gan [76] presented finite-element coupled heat transfer and elastic-plastic thermal stress analysis using a general purpose of commercial FEM to simulate the in-situ residual stress generation during plasma spraying of duplex thermal barrier coatings.They also proposed a simple method by post-deposition treatment to effectively reduce the residual stress level.Sarikaya and Celik [77] studied the thermal evolution of MgO-ZrO 2 /NiCrAlY coatings on Ni metal and AlSi alloy substrates, using the finite-element method, and showed that stress concentration was at higher levels in the coatings with AlSi substrate than that when with Ni, and that the larger residual stress was obtained with an increase of the coating thickness.Similarly, Teixeira [52,78] presented a FEM model of the residual stress distribution with a layered metal-ceramic composite coating, based on the difference between the coefficients of thermal expansion of the neighboring layers, obtaining an analytical expression of the residual stress.They also discussed the failure modes of a thin coating under residual stress as: (a) delamination (coating under tensile stress and ''weak'' interface, i.e. low adherence); (b) perpendicular micro-cracking (coating under tensile stress and a ''strong'' interface); and (c) buckling and spalling (coating under compressive stress).Moreover, Hsueh and Lee [79] also calculated the elastic thermal residual stress in thermal-barrier coating, which consisted of a substrate, a graded bond coat, and a top coat, with the origin of different coefficients of thermal expansion [80].Liu et al. [81] investigated the effect of a system composed of Al 2 O 3 and SiC on 316L stainless on the residual stress developed in the composite coatings, and showed that the presence of the graded properties and the compositions within the coating did lead to the reduction of the stress discontinuity at the interfaces between the coating and the substrate whereas the magnitudes of the residual stresses on the coating surface and at the coating/substrate interface were dependent on the Al For developing a generic procedure of FE simulation of multilayered coatings, a parameterized has been developed recently [89].The method which considered key geometrical, material, interfacial and loading variables as model variables in the FE-based parameterized modelling, was developed with a view to improving the efficiency and accuracy of the analysis and design of multi-layered coating-systems.The method D. Yin et al.: Manufacturing Rev. 2014, 1, 8 allows for the parameters to be changed easily during a seriesanalysis.Combined with the capability of the prediction of cracking of the coatings, the developed method/model provides an efficient way for investigating the effects of these parameters on the behavior of multi-layered systems.It was demonstrated with analysis of the coated tool-steel (H11) for: (a) a substrate without being pre-heat-treated and (b) two substrates with a shallow and a deep hardened-case respectively (both are treated by plasma-nitriding).The results showed that the casehardening of a substrate has a significant influence on the performance of the surface system with coating, especially on its load-bearing capacity and on its the resistance to the development of cracks in the coating.
It is worthy of noting that in reference [89], a bilinear cohesive-zone model was employed to simulate the development of cracks in the coating of a tool subjected to indentation, the general procedure being as follows: (a) Definition of the geometry model and meshing scheme.
In this step, model dimensions and meshing scheme are established.Figure 4 shows the structure of a multilayered system indented by a spherical indenter and the corresponding geometry.Based on the assumption of an axi-symmetric problem, only one half of this model is shown.The axi-symmetric model is made up of five blocks including: a central and an outer coating layer; a central and an outer hardened case; and the substrate, where: W1 and W2 indicate the width of the central and the outer parts respectively; and T1, T2 and T3 indicate the thickness of the coating, the hardened case and the substrate, respectively.
(b) Materials definition.Any number of materials can be defined for the simulation.Each material is assigned a name.A region in the model can be associated with a particular material through the assignment of its section properties that refer to the material's name.(c) Definition of the contact and boundary conditions.In this step, the surfaces and interactions between surfaces are created.Boundary conditions, which can be recognized as the interactions between the system and environment, also need to be defined.(d) Definition of cohesive elements.The cohesive elements, combined with the FE, using a traction-separation relation in order to describe the interface as a continuum of coatings and the host, are defined in order to simulate the formation and propagation of cracks under indentation.Several traction laws have, to-date, been applied for this kind of analysis, either in a form of exponential relationship [90,91] or in a bilinear relationship [92][93][94] between the stress and the displacement.The project ''Multi-scale modelling for the design and analysis of multi-layered engineering surface systems (M3-2S)'' was undertaken jointly by 11 international partners, funded through the 7th European Framework Programme, wherein the overall objective of the project was to establish integrated, generic, robust multi-scale modelling techniques, based on atomic FE(nano), crystal plasticity FE (micro) and continuum mechanics FE (macro) modelling.For examples, Leopold et al. [96] developed an advanced adaptive FEM software (AAFEM) for the coating-substrate simulation of superlattice TiN/CrN-coatings; Zhuang et al. [97] developed an user-subroutine to implement a multi-axial, continuum damage mechanics (CDM)-based constitutive model into a commercial FE code, and carried out a series of investigations into the TiN/Cu system; Perucca et al. [98] demonstrated performance of the TiN and TiN/CrN nanoscale multilayered coatings on WC cutting inserts for machining GJL250 cast iron.It was shown that the nano-structured TiN coatings, followed by multilayered coatings, have much better performance, and the FE simulation showed a good agreement between the measured cutting-force and the simulated one; and Karimpour et al. [99] presented an inverse method for  identifying properties of the material constituents in their multilayered coatings, combining FE simulations and indentation data, which is extensible to the analysis of multilayered coatings containing any number of constituents.

Molecular dynamics
It is well known that the mechanical behavior of a material with a very small volume differs from that which is typically observed for that with a large volume.Analysis of such a phenomenon at the micro-scale has been largely with Molecular Dynamics (MD) Simulation.It is worth noting that there are two basic assumptions made in standard MD simulations [100][101][102]: (1) Molecules of atoms are described as a system of interacting material points whose motion is decided dynamically by the potential, and (2) No mass changes in the system, i.e., the number of atoms in the system remain the same during the process.Similarly, in continuum models, MD has also been used for the analysis of indentation.For instance, Saraev and Miller [103] used MD simulations to elucidate details of plastic deformation and the underlying deformation mechanisms during the nano-indentation of thin copper films with epitaxial nickel coatings, and observed several deformation mechanisms, such as dislocation pile-up at the interface, dislocation cross-slip and movement of misfit dislocations.Similarly, Cao et al. [104] performed MD simulations of nano-indentation, following a nanoscratching process, to investigate the mechanical and tribological properties of Ni/Al multilayer with a semi-coherent interface.The results showed that the indentation hardness of Ni/Al multilayers is larger than that of a pure Ni thin film, and the remarkable strength of Ni/Al multilayered coatings is due to the semicoherent interface which acts as a barrier to glide of dislocations during nanoindentation process.Furthermore, Fang and Wu [105] carried out a MD simulation to study effects of indention deformation, contact and adhesion on Al, Ni, and Al/Ni multilayered films.The results showed that when the indention depth of the sample increased, the maximum load, plastic energy, and adhesion increased.Shan et al. [106] studied the effect of a single twin boundary parallel to the indented surface on nanoindentation of Ag(1 1 1) films, and showed that the twin boundary has little influence on the elastic modulus of films.The load for the initial yield is observably reduced when the twin boundary is very near the indented surface due to the nucleation of the glissile dislocations on the slip plane parallel to the surface, rather than the formation of the tetrahedral sessile lock when nano-indentation on (a) The tension of multilayered coatings [107].For instance, mechanical behavior, dislocation nucleation and development, and evolution of an interface in Cu/Ag bilayer films during equal biaxial tension were studied using MD simulations.The results showed that dislocations are prone to nucleate at the interface of Cu/Ag bilayer film, and then propagate towards the free surface.The deformational defects consist of partial dislocations and intrinsic stacking faults, accompanying some extrinsic stacking faults and stacking fault tetrahedral which just occur at the interface and propagate in the Cu layer.(b) Interface-structure formation [108].Assisted with MD simulations with EAM potential, the formation of the interface of a Cu/Ag bilayer during sputtering was studied.The method provides a notable way to elucidate the mechanism of interfacial formation, which plays a dominant role in controlling the mechanical properties of multilayered coatings.
A MD model of atomic deposition is shown in Figure 6.Four growth cases, including Ag/Cu(1 1 1), Cu/Ag(1 1 1), Ag/Cu(0 0 1) and Cu/Ag(0 0 1), were simulated under the deposition of E i = 8 eV (the kinetic energy for incident atoms) and T = 300 K.The results showed that the first-deposited Ag layer grows on the Cu substrate by the layer mode, while the firs deposited Cu layer grows on the Ag substrate by the island mode, determined by the surface energies of the film and substrate.

Nano-scale modelling
At the atomistic scale, the first-principles calculation based on the high-precision density functional theory (DFT) [92], is a sophisticated way to provide insight into the bonding configuration, interfacial structure, and electronic properties, which are difficult to obtain by experiments alone at such a scale, as well as the behavior during the ideal tensile/shear process.
Four steps are usually used for building a reasonable interfacial model for multilayered coatings at the atomic scale [109]: (a) Firstly, the preferred textures of each individual layer for the composited system are made sure, e.g. from experimental observations [110]; (b) Secondly, perform convergence calculations on the surface energy with respect to slab thickness are performed to ensure that individual layers of the multilayered coatings model are thick enough to exhibit bulk-like interiors.
Taking the TiN/VN interface for examples, the surface energy, (r 0 s ) under the condition that the chemical potential of respective element equals its bulk total energy can be expressed as follows [111]: where E total slab , E bulk XN and E bulk X represent the total energies of a formula unit of a slab (TiN or VN), XN bulk, and X bulk (X = Ti or V), respectively, and M X and M N denote the number of X and N atoms and A is the surface area.As illustrated in Figures 7b-7g, there are, in principle, six possible terminations for either TiN or VN.(c) Thirdly, considering the interfacial stacking sequences and various terminations of individual surface, a series of possible candidate models could be constructed.As shown in Figure 8, there are three stacking sequences for TiN/VN multilayered coatings.(d) Finally, the adhesion energy (W ad ), a key quantity in predicting the adhesive properties of an interface can be calculated to determine the most energetically stale surface, from [112,113]: where E TIN , E VN , and E IF are the total energy of isolated TiN, VN slab, and their interface, respectively, and A is the interface area.Two steps were taken to estimate the W ad values.Total energies were first calculated for a series of separations as the two rigid slabs were approaching increasingly closer from a large initial separation.Consequently, the calculated energies behaved like a parabola, passing through a minimum at the equilibrium separation.The unrelaxed W ad was obtained by computing the energy difference between the interface at the equilibrium state and the unrelaxed isolated slabs.Then, each isolated slab as well as the interfacial slab was allowed to optimize fully, yielding an estimation of the relaxed value of W ad .
More importantly, first-principles calculations can be employed to investigate original mechanisms from atomic stacking sequences, bonding properties/interactions between neighboring atoms, etc., at the atomistic scale, due to the fact that the atom is the minimum element for the composited systems.
For instance, the two mostly addressed topics for multilayered coatings are Template effect and Superhard mechanism.
For ''Template effect'', the metastable phase can form in multilayered coatings due to the template effect, which means that the crystal structure of the newly deposited layer is under a strong influent of the underlying multilayer system.In this way, the newly deposited layer may crystalline in ametastable structure by forming  a coherent interface with the rest of the multilayer structure.Such a mechanism often leads to strengthening the coatings and thus is desirable.Taking the representative multilayered coatings TiN/ AlN for example, Chawla et al. [114] combined ab initio methods and the finite element (ABAQUS) to study and predict the equilibrium critical thickness up to which the metastable cubic (c) AlN phase is energetically favoured to the stable wurzite(w) variant in TiN/AlN and CrN/AlN bi-layer system.And Zhang and Veprek [115] selected deliberately a couple of representative transition metal (AlN included) as interfacial layers in TiN-based heterostructures, and investigated the thermal stability using first-principles molecular dynamics.The results showed that one-monolayer-thick pseudomorphically stabilized interfacial layer of B1-type AlN is stable within the units within the whole temperature range considered.Moreover, Zhang, Veprek, and Sheng [115][116][117] have studiedmetastablephasesandspinodaldecompositioninTi 1-x Al x N system by ab initio and thermodynamic modelling compared with the TiN-Si 3 N 4 system and found that metastable fcc-Ti 1-x Al x N coatings can easily undergo spinodal decomposition into coherent fcc-TiN and fcc-AlN, but thereis arelatively large barrier for the formation of the stable hcp-AlN.
Furthermore, a comparison with the TiN-Si 3 N 4 system showed that, due to the much higher de-mixing energy of this system as compared to the TiN-AlN one, spinodal decomposition may occur in that system also for semicoherent TiN and Si 3 N 4 phases.Holec et al. [118] selected deliberately a couple of representative transition metal aluminium nitride (TM-Al-N) thin films (Ti-Al-N included) and addressed the structure and phase stability.Especially for Ti 1-x Al x N, the predicted stability regions of the rock salt cubic structures are x 0.7.Similarly, Stampfl and Freeman [119] also investigated the formation, atomic and electronic structure, and stability of metal-nitride interface systems, (1 0 0) AlN/TiN, as well as AlN/VN and VN/TiN in the rocksalt structure, and found that the layerdependent interaction energy between the adlayer surface and the interface, while typically not taken into account, plays an important role in terms of the formation energy for the initial stages of film growth.
The aforementioned references for elucidating the mechanism of the ''Template Effect'' are also applied for the ''Superhard mechanism'' due to the formation of coherent interface, which remarkably decreases the defects, such as initial cracks, voids, etc., as described by Zhang et al. [120].These authors have studied the intrinsic mechanical properties of two AlN polymorphs (fccand hcp-) and found that, in spite of its higher elastic modulus, the fcc-AlN possesses slightly lower shear strength as compared to the stable hcp-AlN.Moreover, the calculated electronic density of states indicates a similar covalent bond property for both polymorphs.Consequently, the higher hardness-enhancement of fcc-TiN/fcc-AlN heterostructures as compared to fcc-TiN/hcp-AlN nanolaminates is not related to the difference of the intrinsic strength and bonding nature between fcc-and hcp-AlN, but most likely to the formation of semi-coherent interface.

Multiscale models
As the result of the confluence of parallel computing power, experimental capabilities to characterize structure-property relations down to the atomic level, and theories that admit multiple length scales, there has been rapid growth of activities in multiscale modelling.In order to make the computations tractable, multiscale models generally make use of a coarse-fine decomposition.An atomistic simulation method, such as MD, is used in a small subregion of the domain in which it is crucial to capture the individual atomistic dynamics accurately.A continuum simulation is used in all other regions of the domain in which the deformation is considered to be homogeneous and smooth.Since the continuum region is usually chosen to be much larger than the atomistic region, the overall domain of interest can be considerably large.The aforementioned QC methods, based on a combined MD-FE technique, can be applied to a multiscale modelling simulation, where MD is only used in localized regions in which the atomic-scale dynamics are important, and a continuum simulation (FE) used everywhere else.For instance, Sun et al. [121,122] have applied QC methods to simulate the nanometric cutting of crystal copper.Their study demonstrated that multi-scale simulation is feasible, not withstanding that there is still more work needed to be done to make multi-scale simulation more practical.The QC method is a directly coupled model between two scales (MD and FE).The unique characteristics of such a directly coupled method is that there is a hand shaking (HS) region to connect various scales wherein the message exchange, such as the velocity of particles (atoms or molecules), temperature, etc., is always a big challenge.At the same time, there is another way to utilize the advantages of various single scales, that is to say, the usage of one single scale could provide parameters/evidence for further investigation in another single scale.For instance, both FEM and MD simulation for nano-indentation showed that the calculated mechanical properties intrinsically depend on the interface contact conditions of the film/substrate [123].Especially, investigations into the ''template effect'' from first-principles calculations have ascribed the super-hardness of multilayered coatings as the formation of coherent interfaces between the neighboring individual layers.Aiming at evaluating the residual stress caused by lattice mismatch between neighboring layers, a approach based on the representative volume element (RVE), as shown in Figure 3, incorporating with first-principles calculations, is proposed from continuum mechanics principles [124].At the atomic scale, the evidence for the good bonding between the multilayered coatings has been presented by first-principles calculation, as shown in Figure 9.It can be seen in these figures that, at first, the charge distributions of TiN, CrN show a spherical symmetry.Similar charge distributions can also be observed between the interfacial Ti-N bonds and the Ti-N bonds deeper to TiN: the majority of the charge being located on both Ti and N atoms, with slight but visible distortions directed towards each other.Not surprisingly, for TiN/CrN multilayered coatings, charges that are accumulated on the interfacial N come predominantly from adjacent Ti and Cr atoms, as shown in the right panel of Figure 9.It was concluded that the interfacial bonding Ti-N (or Cr-N) is mainly ionic, yet maintains a small degree of covalency, which indicates the well bonded interfaces of the TiN/CrN multilayered coatings.After these, material parameters, such as Young's Modulus (E), Poisson's ratio (v), lattice constant (a), have also been calculated.In the following, the constitutive relationship of each phase can be expressed as: ðfor the a-phase ðCrNÞÞ; ð4Þ where De * is mismatch strain induced by the lattice mismatch between neighbouring layers.Based on the assumption that the in-plane strain and the out-of-plane stress components in both TiN and CrN phases should be identical to the corresponding values in the RVE to meet the compatibility and equilibrium conditions [125], the analytical solution of residual stress can be obtained and the relationship between various parameters can also be described.The whole process presents another way to take advantage of modelling at various length-scales either in a top-down approach or in a bottom-up one.

Discussion
Multi-scale modelling can synthesizes atomistic and continuum methods to provide a more natural and coherent description of material phenomenon than that obtained from the fragmented perspective for a limited individual scale.It has experienced some success with metals, but the usage for transition nitride multilayered coatings has been lacking.Especially, there is essentially a request for research related to multiscale modelling to deal with defect/boundary/cracks.The challenges on the future development of numerical modelling for multilayered coatings may be from the following aspects: 1. Continuum models.Over the past few decades, continuum models have dominated materials modelling research.This approach of predicting material deformations and failures by implicitly averaging atom-scale dynamics and defect evolutions over time and space, however, is valid only for large enough systems that include a substantial number of defects.As a result, numerous experimental observations of multilayered coatings cannot be readily explained within the continuum mechanics framework: dislocation patterns in fatigue and creep, surface roughening and crack nucleation in fatigue, inherent inhomogeneity of plastic deformation, the statistical nature of brittle failure, plastic flow localization in shear bands, and the effects of size, geometry and stress state on yield properties [126], etc.Moreover, the accuracy of the models needs to be validated for each case dealt with.2. Molecular dynamics models.A MD model could only explore the phenomenon occurring at nanometric scale.Its applications to the analysis of multilayered coatings are limited by two factors: (a) Availability of the values of potentials for the atomic interactions to be examined.As well known, a key factor that determines the accuracy of a MD simulation significantly is the potential used to describe the interaction between atoms.A great effort has been made to develop various interatomic potentials fitted to different first principles, empirical and semiempirical models, but only a few potentials for relatively simple material systems have been developed until now.Typically, multilayered coatings consist of at least two elements, and fewer potentials are available due to the increased complexity of interatomic action.On the other hand, simulations will often give acceptable results only if the potential was fitted for the problem under consideration, or one that is very similar.In other words, interfaces of multilayers play an important role in their properties.Therefore, more attention should be paid to choice and development of potentials with proper descriptions of interfacial properties before these multilayered systems could be simulated with MD accurately.(b) A typical MD simulation is still restricted to a very small system consisting of several million atoms or less and timescales in the order of picoseconds.A general acceleration methodology would have revolutionary implications that would relax current limitations on the use of MD. 3. First principles calculations.Many of the first principles calculations are confined to zero temperature, wherein the basic quantity is a Hamiltonian for the system expressed, in terms of the appropriate degrees of freedom.In principle, Hamiltonian methods can be extended to equilibrium at finite temperatures by using the free energy.Inherently non-equilibrium situations are, however, fundamentally different and a general approach is not yet within grasp.At the same time, the limitation to the computational scale, i.e. the maximum number of atoms which can be included in a numerical model for first principles, is limited.Generally, the maximum number of atoms used in a model does not exceed 100.The concept of multi-scale modelling embodies a comprehensive description of a material, which requires understanding over both multiple time and length scales.Although the modelling for crystalline materials has coupled the nano-and intermediate-scale in a highly effective manner, the challenges remaining are the for longer time and length scales, largely due to the high requirement on computational effort associated with high-resolution calculations over such scales.In addition, multi-scale models also require that descriptions at all levels be consistent with each other.The requirement for consistency often demands the models for different process or pieces of a system to be coupled and solved concurrently.Consequently, multiscale modelling, on one hand, would still need more effective coupling methods as an enabler.On the other hand, it does not satisfactorily address issues such as disparate time scales in the two regions, and it provides only a rather simplified treatment of the interfaces between the atomistic and continuum regimes.
For the aforementioned unique properties of multilayered coatings: the thickness of the individual layers of a coating would be at the nano-scale while the entire thickness of the coating is at the macro-or meso-scale.Single-scale methods such as first-principles calculations, Molecular Dynamics and the Finite Element Method will have difficulty in analyzing such multilayered coatings alone, due to the limitations to the time and length scales that each method is confined to.With further development of accurate inter-atomic potentials for the material-components in multilayered coatings and further development of time-integration algorithms used in the simulation, the MD and Multiscale modelling (MD coupled) will become a prominent tool for elucidating complex physical phenomenon.Although multi-scale modelling is important for the design and analysis of multilayered coatings, effort for the development of the modelling at individual length-scales could not be neglected, that is to say, efficiency and effectiveness of multiscale modelling not only concern how the models for individual length-scales are coupled/linked but also depend on how individual modelling methods and techniques for a particular length-scale could advance.

Conclusions and remarks
Multilayered coatings are currently a subject of intensive research because of their potential for a wide range of technological applications and also due to their scientific importance.In this paper, development of the following numerical methods and models have been reviewed: (i) Continuum models; (ii) Molecular Dynamics simulations; and (iii) Nano-scale (first-principles calculations) modelling, for which examples of typical applications have been given.
In continuum models, FEM has been applied mainly for investigating the stress distribution, formation and propagation of crack under the (nano-) indentation or scratch test, formation of residual stresses, etc.Similarly, Molecular Dynamics (MD) simulations have also been applied mainly to the analysis of nano-indentations.Applications of first principles calculations are mainly presented on two kinds of mechanisms: ''template effect'' and ''super-hardness''.Multiscale modelling, such as that with QC methods, has also been discussed.
In summary, multi-scale modelling is unambiguously a particularly vibrant field of research related to engineering and sciences.Much of the excitement is being fuelled by increasingly sophisticated experiments that are capable of probing mater at the level of individual atoms or molecules.Interpretation of their experimental results require molecular and multiscale models to be connected with experimental observations.It is felt that more effort on the numerical modelling of multilayered coatings is needed.For example, when 400 nano-layers coated onto a substrate material needs to be modelled, there has to be a more effective method and more efficient way to achieve it.Developing accurate descriptions of inter-atomic potentials for multi-material elements, improving time-integration algorithms used in numerical simulations and developing more effective and efficient coupling methods for multiscale modelling are urgent tasks to be undertaken in order to develop a powerful numerical tool for high-efficiency performance prediction and process design for multi-layered coatings on engineering components, machining and forming tools.

Figure 3 .
Figure 3. Schematic illustration of an RVE of a multilayered coating.

Figure 2 .
Figure 2. Schematic numerical simulation process from Macro-to Micro-to Nano-scale.

4 D.
Yin et al.: Manufacturing Rev. 2014, 1, 8 (e) Validation of the FE model.The FE model is validated through a series of indentation tests.A comparison of the experimental results with finite-element results is shown in Figure5.It should be mentioned that the adhesion energy values used in Cohesive-zone modelling is critical and that these can be obtained either from calculations based on experimental data or by means of computations based on Molecular Dynamics or First principles calculation based on high-precision density functional theory (DFT)[95].(f) FE analysis.After validation of the model has been done, the model can be used to analyse cases of interest.The analysis can be re-activated for a different set of input-variables, simply by re-setting values of these variables in an input file, without need to modify the model through using a CAD/CAE software.

Figure 4 .
Figure 4. Geometric model of a multi-layered system including the coating, the hardened case and the substrate.

Figure 5 .
Figure 5.A comparison between (a) Experimental results and (b) Finite Element ones.

Figure 6 .
Figure 6.Schematic diagram of the MD model of atomic deposition.

Figure 8 .
Figure 8. Schematic plots of (a) OT, (b) SL, and (c) TL stacking sequences.Upper part shows top view and lower part side view.The interfaces between N1-terminated VN and Ti-terminated TiN are given only.For top views, only the Ti layer nearest to the interface is shown for clarity.The view direction for bottom panels is along [1 " 1 0] and the location of interface is marked with horizontal dotted lines.The top and bottom portions of the interface have been omitted.

Figure 9 .
Figure 9. Contour maps of charge density (left) and its difference (right) for TiN/CrN multilayered coatings, taken along the (1 " 1 0) plane.The difference in charge density represents redistribution of charge in the interface relative to the isolated system.The interface is indicated by a horizontal line and atoms that intersect the contour plane are labeled.The upper scale denotes the magnitude of charge in the left panel and lower scale that in the right panel scale.
[85]larly, Toparli et al.[83]determined the thermal residual stresses in WC-Co/Ni-Al coating layers deposited on 316L stainless steel substrates, developed during and after thermal cycling, showing that the calculated tensile stresses were higher than the compressive stresses.Recently, Far et al.[84]developed a FEM model to evaluate the stress induced by the thermal cycle in a typical plasma-sprayed thermal-barrier coating system, taking into account the effect of thermal and mechanical properties, morphology of the top-coat/bond-coat interface, and oxidation on the local stresses that are responsible for the micro-crack nucleation during cooling, especially near to the metal/ceramic interface.5.Other applications.Apart from the applications described above, FEM has also been used for the following particular applications.For examples, Pan et al.[85]compared the sliding contacts of gradient layer, sandwich layer and ordinary multilayered coatings, subjected to a similar surface load and under elastic-plastic deformations, and found that the positive gradient layer coating can improve the stress and strain fields near the layer-substrate inter- [86]3 and SiC coating thickness.Using FEM, Sayman et al.[82]carried out a transient thermal analysis in WC-Co/Cr-Ni multilayered coatings deposited on a 316L steel substrate during the cooling process from 1 273 K to room temperature.face,and in the substrate, meanwhile, the sandwich layer coating can greatly reduce the maximum shear stress on the interface of the layer-substrate and the stress at the crack tip.Conversely, the stress and strain distributions of the ordinary layer have not evidently been improved, compared to those of the single layer.Their work provided a useful reference for the selection and design of coatings.Furthermore, Bemporad et al.[86]designed the thickness and Ti distribution in multilayered Ti/TiN coatings with FEM-based optimization leading to a significant increase (45%) in adhesion between the layers and the WC-Co interlayer (Ti-6Al-4V substrate) and load bearing capacity of the coated system, compared to monolayered TiN, without reduction in superficial hardness and in the load bearing capacity.Lakkaraju et al.[87] modeled multilayer Cr/CrN thin films and optimized them to have effective ''load support'' by the films on stiff A2 steel and compliant 2024-Al substrates.The