A comparison of different approaches to enforce lattice symmetry in two‐dimensional crystals

The phase‐field crystal (PFC) model allows for the resolution of atomic‐scale structures on diffusive time scales. It is based on an approximation of the two‐particle direct correlation function in the free energy, which provides the symmetry of the lattice structure. Various approaches have been proposed to model common lattices. We here only focus on two‐dimensional crystals, review the different approaches and demonstrate that one, which is based on affine linear deformations, violates rotational invariance. The consequences of this physical inconsistency are shown numerically using a pseudospectral algorithm for the higher order equations.


Introduction
We consider phase-field methods working on the molecular scale, termed phase-field crystal (PFC) models as introduced in [1,2]. These are phase field models with a periodic order parameter, which allow the resolution of atomic-scale structure and defects that evolve on inherently diffusive time-scales. Their local state is characterized by time averaged particle densities that might be related to molecular dynamics under appropriate conditions [3]. The PFC model of [1,2] can be derived from classical density functional theory via expanding the two-particle direct correlation function, c 2 (k), up to the 4th order term, a procedure that yields a single peak in the direct correlation function, as detailed in [4,5]. PFC models have many of the complex physical mechanism that need to be built into traditional phase field models such as strain relaxation, elasticity, and arbitrary grain orientation. Within PFC models they arise holistically from the form of the free energy. The dimensionless free energy measured relative to a reference liquid can be expressed in terms of reduced temperature and reduced density ψ, leading to a Swift-Hohenberg type expression: where all quantities are dimensionless. In 2D, the vs. ψ phase diagram of the PFC model includes stability regions for the homogeneous liquid, for a crystalline phase of triangular structure, and for a striped phase, and appropriate coexistence regions in between. We are here only interested in the crystalline phase. The model contains a single model parameter that incorporates a combination of the expansion coefficients of the direct correlation function, which in turn can be related to the compressibility of the reference liquid, the bulk modulus of the crystal, and the interatomic distance. Fitting the phase diagram to that of real matter, at least within a limited temperature/density window was already done in [6] and has been extended to higher order expansions of the two-particle direct correlation function [7,8]. In order to model other lattice structures a multi-mode PFC (MMPFC) model was introduced in [9] extending approaches in [6,10,11]. Assuming that the two-particle direct correlation function has N peaks and is isotropic the respective free energy reads where Λ, b i , τ and ξ are phenomenological constants and the N peaks of the two-particle direct correlation function are located at wave numbers Q i . [9] demonstrate that a system with three competing length scales (N = 3) can order into all five Bravais lattices in 2D. Another approach to achieve the same goal is considered within the extended PFC (XPFC) model [10,12]. Here the twoparticle direct correlation function c 2 (k) is approximated in k-space as the envelope of a set of Gaussians and with peaks chosen by the primary k-vectors defining the crystal structure. The free energy reads with t and v appropriate parameters.
Also more general fully nonlinear approaches can be used to realize different lattice structures, see [11,13,14]. The nonlinear modification of the PFC energy reads The terms g 2n ψ 2 ∇ 2n ψ 2 and s 2n |∇ψ| 2n couple density waves with different k-vectors and the relative orientation of the k-vectors becomes important [11]. These models, maybe with exception of the last, due to its increased complexity, have been intensively used to analyse the influence of lattice symmetry on both the structure and dynamics of crystalline and defected systems in 2D materials, see e.g. [15][16][17][18]. The dynamics is thereby based on the assumption of overdamped conservative dynamics, which read The equation of motion, eq. (5), can also be deduced from dynamic density functional theory [5].
In the crystalline state ψ can be expanded in terms of its Fourier components A q and the reciprocal lattice vectors q: ψ(x) = ψ 0 + q A q e iq·x , where ψ 0 is the average rescaled density. In 2D, q = nk 1 + mk 2 , where n and m are integers and k 1,2 are the principle reciprocal lattice vectors. With Q i , in the MMPFC model, chosen such that the magnitudes of the critical wave vectors correspond to the three shortest wave numbers of the target lattice we obtain |q| = k 1 m 2 + µ 2 n 2 + 2mnµ cos θ, where µ = k 2 /k 1 and θ the angle between k 1 and k 2 . A square lattice can, e.g., be obtained by Q i=0,1,2 = 1, √ 2, 2, leading to two possible sets of density wave triples: two with |q| = 1 and one with |q| = √ 2, or two with |q| = √ 2, and one with with |q| = 2. For appropriate parameter sets to realize different lattice symmetries we refer to [9]. For a square symmetry in the XPFC model a minimum of two peaks is needed, c 2 (k) = max(c 2,0 (k), c 2,1 (k)) and c 2, The effect of temperature on the elastic properties is seen in the width of the peaks and modeled by ξ i and A i is a Debye-Waller factor controlling the height of the peaks. Within the NLPFC model g 2n and s 2n can be chosen to lead to square ordering. That is, the energy is minimized if the k-vectors become perpendicular [11].
All these approaches require higher order derivatives or even nonlinear higher order terms in the free energy, which makes their numerical solution difficult. Thus the question arises, if there are simpler models to obtain a square lattice. An alternative route to obtain the five Bravais lattices in 2D has been outlined in [19], based on an idea in [20]. It considers an affine lattice deformation (ALD) using hexagonal symmetry as the reference structure [21]. The general form with transformation matrix A reads For the original PFC model [1,2] with a crystalline phase of triangular structure, we have A = I, with I the identity matrix and a square lattice can, e.g., be obtained by A square = 1 −1/ ∂x1∂x2 , see [21]. The PFC energy for this approach reads with the equation of motion as eq. (5). This approach has recently been used in the context of battery electrodes [22]. We will in the following concentrate only on square lattices and compare the different modeling approaches. The results for the other Bravais lattices in 2D will be similar.

Numerical approach
The equation of motion eq. (5) can be solved via a pseudo-spectral algorithm [23,24] using periodic boundary conditions, which provides a relatively simple implementation for all approaches. We here consider a time discretization which simply takes all linear terms implicitly and all non-linear terms explicitly. Only NLPFC looses much numerical stability and for efficient simulation more elaborate methods should be used [25]. Examples for other efficient discretisations, based on finite elements, can e.g. be found in [26,27].

Physical consistency
Unfortunately the ALDPFC model, eq. (6), is not invariant under rotation. To see this, we consider transformed differential operator in a rotated coordinate system: x and x are coordinates in the original and rotated frame of reference and R the associated rotation matrix. For general ALD A differs from ∇ t x A t A∇ x . This can be easily checked by taking A Square and R = 0 1 −1 0 . In two dimensions only A ∈ SO(2) commutes with R and leads to rotational invariance. But this transformations does not allow to modify the crystalline symmetry. This discrepancy can also be seen in the simulations in Fig. 1, where the particle density peaks are elongated in the [1,-1]-direction, see Fig. 2 for details. The elongation is independent on the local lattice structure.

ALDPFC
MMPFC XPFC NLPFC The rotational variance of the ALDPFC only stabilize square ordering for a single orientation. The other modeling approaches do not show any preferred direction. All density peaks are ordered in squares independent on orientation, see Fig. 1. But the shape of the density peaks differs, see Fig 2. NLPFC has quite broad peaks with a symmetry between density minima and maxima, both ordered in squares. NLPFC only includes the first wave vectors describing a square phase. The triangular phase is inhibited by the nonlinear extensions. XPFC and MMPFC include the first two resp. three wave vectors describing the square ordering. This increases the number of considered wave vectors in the correlation function and thus leads to sharper density peaks and a broken symmetry between density minima and maxima.

ALDPFC
MMPFC XPFC NLPFC Fig. 2: Crystalline symmetry for the differnet modelling approaches. The pictures are taken form the same simulation as Fig. 1 and peaks in the density wave are dark. The affine transformation leads to elongated peaks in the density wave and rotation lead to deformed square ordering. MMPFC, XPFC and NLPFC leads to square ordering, but the shape of the peaks depend on the used model.
Violation of rotational invariance and with it Galilean invariance as a basic principle of Newtonian mechanics, makes models based on eq. (6) inappropriate. On the other side, the MMPFC, XPFC and NLPFC models fulfill rotational invariance. The density peaks are isotropic and the square lattice is achieved for all orientations. Even if not compared in detail, also the structure of the grain boundaries agree between the MMPFC, XPFC and NLPFC models and differ from the ALDPFC approach. www.gamm-proceedings.com

Conclusion
Various extensions of the PFC model allow to enforce different lattice symmetries. Even if the approach based on an ALD [19,20,22] looks attractive at a first glance, it turns out to be physically inconsistent, as it violates rotational symmetry. All other approaches the MMPFC, XPFC and also the NLPFC model allow to reproduce lattice symmetries consistently. At least within a pseudo-spectral algorithm with periodic boundary conditions the appearing higher order derivatives in these models can be treated relatively easily. Our simple discretization of the NLPFC needed two order of magnitude smaller time steps and more elaborated discretization methods are needed for efficient simulations e.g. [25].