Simulation of crystalline pattern formation by the MPFC method

The Phase Field Crystal model in hyperbolic formulation (modified PFC or MPFC), is investigated as one of the most promising techniques for modeling the formation of crystal patterns. MPFC is a convenient and fundamentally based description linking nanoand meso-scale processes in the evolution of crystal structures. The presented model is a powerful tool for mathematical modeling of the various operations in manufacturing. Among them is the definition of process conditions for the production of metal castings with predetermined properties, the prediction of defects in the crystal structure during casting, the evaluation of quality of special coatings, and others. Our paper presents the structure diagram which was calculated for the one-mode MPFC model and compared to the results of numerical simulation for the fast phase transitions. The diagram is verified by the numerical simulation and also strongly correlates to the previously calculated diagrams. The computations have been performed using software based on the effective parallel computational algorithm.


Introduction
The phase field crystal (PFC) model [1] has been proposed to incorporate the physics naturally embedded on nano-length scales and on diffusive time scales.The PFC-model uses the free energy functional of the Swift-Hohenberg (SH) form that describes pattern formation in Rayleigh-Benard convection [2].The PFC-model may also be expressed as a conserved version of the Swift-Hohenberg equation and as an efficient method for simulating liquid-solid transitions [3,4], colloidal solidification [5], dislocation motion and plasticity [6,7], glass formation [8], epitaxial growth [9,10], grain boundary premelting [11], surface reconstructions [12], and grain boundary energies [13].
In the present paper, numerical calculations based on the Phase Field Crystal (PFC) model are presented which may lead to the solution of many practically important problems in metallurgical industry [14].The determination of the conditions of the technological process for metal castings with predetermined properties, prediction of crystal structure defects during casting, evaluation of the quality of special coatings, and many others are among them can be start to solve beginning from micro-level analysis using PFC-method [14,15].

Free energy of system and MPFC equation
The modified phase field crystal model (MPFC) describes an atomic density field ϕ(x,t), which is conserved variable, and expressed by the sixth order in space and second order in time equation [14][15][16][17]: Here t is the time, τ is the relaxation time of the atomic flux to its stationary state, and  represents the chemical potential.The chemical potential may be obtained from functional of free energy: associated to the computational domain Ω.The Eq. ( 1) includes the chemical potential μ as a variational derivative of the free energy functional F: The function of the homogeneous reference energy density f takes on as: where ε is the dimensionless undercooling ε=(Tc-T)/Tc with Tc -critical temperature of phase transition, Ttemperature of the system, and the term α -the coefficient which controls the metastability of transition.For all our calculations, we state α=0 for the unstable phase transitions modeling with zero potential barrier of the phase transformation.Introduction of the inertial term with τ in Eq. ( 1) changes the class of equation from  5), body centered cubic BCC Eq. ( 6), rods Eq. ( 7), stripes Eq. ( 8).The top-left region for ε>0.6, ϕ<0.4 represents region of coexistence of three structures, where the regions of coexistence of BCC+rods, homogeneous+BCC and homogeneous+rods intersect and match each other exactly.As it can be found on the depicted results of direct computation this region exposed to form irregular structures.parabolic to hyperbolic and plays significant role in the analysis of fast phase transformations [16,17].

The diagram of structures and analytical solutions
The structure diagrams are calculated using the solution of the Maxwell area rule, using the chemical potentials and free energy functionals for each structure.The structure selection is based on the thermodynamic rule of the minimal energy for a virtually existed structure.The structure diagram is constructed in coordinates "dimensionless undercooling -relative atomic density" with the amplitude expansion of the free energy (2) for a given crystalline phase [18][19][20].

Free energy of crystalline states
Using the summation of reciprocal lattice vectors [14] we can make the amplitude expansion and get the dimensionless density profiles indicating the symmetry and properties of a given phase (homogeneous, BCC, rods, stripes).For all structures we choose as a reference the homogeneous phase relative density ϕ: We substitute the amplitude representations of each structure Eqs.(5)(6)(7)(8) to the free energy Eq. ( 2) and get a dependence from the complex atomic density vector amplitude module η, reciprocal lattice vector q and dimensionless atomic density ϕ.Following the minimization procedure, an equilibrium wave vector q=qeq and minimize the amplitude η are obtained numerically [18,19].Finally, we substitute this roots to the free energy (2) and obtain the form of the free energy for each structure.

Construction of structure diagram
For the determined above free energies, we construct the diagram of the coexistence of 3D structures in terms of the expressions Fh (for homogeneous phase), FBCC, Fstripes, Frods dependent from driving force ε and dimensionless atomic density ϕ.With this aim, we write the Maxwell area rule [13,14,[18][19][20][21] and solve the equations for the entire range of control parameter values 0 ≤ ε ≤ 1.The example of the system of equations showed below representing the conditions for the coexistence of the homogeneous and BCC phases: Fig. 1 shows the domains of existence and coexistence of homogeneous (liquid or metastable melt) phase, the crystalline BCC, rods and stripes structures.We calculate the coexistence separately for the each pair of the structures and then select the existence regions for each structure by its minimal free energy as a criteria for selection for a given region of parameters.

Numerical calculations of structure formations by the MPFC method
The numerical formulation of the MPFC equation is based on the isogeometric analysis (IGA) [22,23] for the spatial discretization.Such method allows us to perform the C 2continuous functions for the discretization of the sixthorder differential equation in primal form.The integration in time is executed by the generalized-α method which was suggested by Chung and Hulbert [24].Thus an implicit finite element method (FEM) for solving the problem of modified PFC was performed.Within this procedure, a finite-element regular hexagonal mesh is used for simulations.The details of the algorithm and its software implementation are given in [25] where, in particular, was shown that the numerical solution of the MPFC problem is consistent to the diagram obtained by Jaatinen and Ala-Nissila in [26].The computational efficiency of such algorithm was shown in [27,28].In [29,30], the analysis of solution's accuracy and physical stability was presented using this algorithm.The present numerical simulations have been provided using the special software [25] based on the PETIGA and PETSc libraries [31].For the three dimensional computational domain the spatial mesh has been generated by C 2 -continuous quadratic elements.We set an initial distribution with a small spatial gradient of the atomic density.Such an initial distribution has no influence on the modeled crystal pattern in our simulations.The boundary of the computational domain has as periodic conditions.

Results and discussions
The initial state of the atomic density has been characterized by the average value of ϕ.Then, the system is quenched with a strength given by the driving force parameter ε.These two values define a position in the phase diagram (see Fig. 1).The solid lines in Fig. 1 represent the structure diagram calculated for three dimensional structures given by the static PFC-equations with free energy (2) and chemical potential (3) in onemode approximation.As it can be seen in the Fig. 1, the present modeling gives the following set of structures: (1) stripes, (2) mixture of stripes and rods, (3) rods, (4) mixture of BCC and rods, (5) BCC structure, (6) mixture of BCC and liquid, (7) multiple patterns highly sensitive to values of undercooling and atomic density, (5) liquid.All of these modeled structures are perfectly consistent with the regions of the predicted structure diagram.Also is important to note that our diagram correlates with the diagram Jaatinen et al. [26] that also was verified by the presented method [25].The coincidence of the classical PFC model predictions and numerical dynamical calculations shows the viability of the presently developed algorithm for the numerical solution of the modified three dimensional PFC-equation.The most complicated region of the diagram is the coexistence of three structures (the top-left region for ε>0.6, ϕ<0.4, where homogeneous-BCC-rods regions of coexistence intersect) in which dynamical MPFC system has a tendency to evaluate to irregularly mixed structures.These mixed structures may be obtained in modeling very often that confirms the idea about evolution of the system through a series of the metastable states.
The modeling of the real systems could face problems with the determination of the crucial parameters, which have to be overcomed by measuring the structure factor of the liquid and crystal state.The necessary data my be provided by the X-Ray diffraction measurements which gives us the reference parameters of the liquid state (density), equilibrium wave number and structure factor for both phases.The perspective system for the calculation of the nature of transition may be the pure Iron with the Homogeneous-BCC-FCC-BCC phase transformations [18][19][20][21].Moreover there is an attempt to use the molecular dynamics data for such calculations [13,21].
The evolution away the mixed and defected structures can be provided by introducing the noise term.The noise is affected on the structure with the "shake-up" effect [30,32].Indeed, locally metastable and structurally mixed states have to be correlated in time and in space, due to the overlapping of the fluctuating atomic configurations [33,34].In this case, a colored noise should be added to the MPFC equation.Previously, in the analysis of spinodal decomposition, it was notable shown that the consistently introduced colored noise influenced the pattern formation in a great extent [34].

Conclusions
Structure diagrams were calculated for the one-mode PFC model and compared with the numerical simulation of the full MPFC equation for the fast phase transitions.The computations were performed using parallel algorithms based on isogeometric analysis, which is known extension of the finite element method.The evolution of three dimensional crystal structures to their steady equilibrium state was simulated for various atomic densities and temperatures using the MPFC-equation and compared with the stationary solutions for regions of structural coexistence.The latter is obtained from amplitude expansion of free energies of the PFC-equations.The obtained diagram verified by the direct simulation and also it is strongly correlate to the diagram Jaatinen et.al. [26].
The proposed method makes possibilities to extend the description of the fast dynamical phase transitions by the calculation for metastable structural boundaries.The MPFC formalism allows us to model the bunch of the real liquid-crystalline and liquid-smectic phase transitions even for two-dimensional systems [15,18,19].The extension of the method with the two-mode approximation creates the way to model the complicated three dimensional crystalline structures (such as crystals with FCC-and HCP-lattices).Further studies on the simulation of crystalline patterns by the MPFC method will allow developing an effective software tool for research and develop of materials in industry, including manufacturing and processing technologies.