CALANIE: Anisotropic elastic correction to the total energy, to mitigate the effect of periodic boundary conditions,☆☆

https://doi.org/10.1016/j.cpc.2019.107130Get rights and content

Abstract

CALANIE (CALculation of ANIsotropic Elastic energy) computer program evaluates the elastic interaction correction to the total energy of a localized object, for example a defect in a material simulated using an ab initio or molecular statics approach, resulting from the use of periodic boundary conditions. The correction, computed using a fully elastically anisotropic Green’s function formalism, arises from the elastic interaction between a defect and its own periodically translated images. The long-range field of elastic displacements produced by the defect is described in the elastic dipole approximation. Applications of the method are illustrated by two case studies, one involving an ab initio investigation of point defects and vacancy migration in FCC gold, and another a molecular statics simulation of a dislocation loop. We explore the convergence of the method as a function of the simulation cell size, and note the significance of taking into account the elastic correction in the limit where the size of the defect is comparable with the size of the simulation cell.

Program summary

Program Title: CALANIE, version 2.0

Program Files doi: http://dx.doi.org/10.17632/3h6xffk9h6.1

Licensing provisions: Apache License, Version 2.0

Programming language: C/C++

Nature of problem: Periodic boundary conditions (PBCs) are often used in the context of ab initio and molecular statics atomic scale simulations. A localized defect in a crystalline material, simulated using PBCs, interacts elastically with its own periodically translated images, and this gives rise to a systematic error in the computed defect formation and migration energies. Evaluating the correction to the total energy resulting from effects of elastic interaction between a defect and its periodic images, to alleviate the contribution to the total energy arising from PBCs, is an essential aspect of any accurate total energy calculation performed using PBCs.

Solution method: The energy of interaction between a localized defect and its periodically translated images is computed in the linear elasticity approximation. The energy of elastic interaction is expressed analytically in terms of the elastic dipole tensor of the defect and elastic Green’s function. Elements of the dipole tensor are computed as a part of the simulation evaluating the formation energy of the defect. Elastic Green’s function and its first and second derivatives are computed numerically from the elastic constants of the material. The method and the corresponding numerical procedures are implemented in the CALANIE computer program. The program evaluates matrix elements of the elastic dipole tensor of a localized defect and the elastic correction to the total energy arising from the use of periodic boundary conditions.

Restrictions: The approach assumes the validity of the linear elasticity approximation. This limits the accuracy of evaluation of the elastic correction, which becomes less precise if the size of the defect is comparable with the size of the simulation cell.

Unusual features: An open source code, containing full detail of the relevant theoretical concepts, algorithms and numerical implementation.

Introduction

Mechanical deformation, or irradiation by energetic particles, produces a variety of defects in a crystalline material, for example dislocation loops, vacancy clusters, voids, Frenkel pairs, and extended dislocations [1], [2], [3], [4]. Defect structures evolve under the effect of external stress and temperature. Defects migrate, segregate and agglomerate as a result of elastic interaction, mediated by the deformation of the crystal lattice [5], [6], [7], [8], [9]. Evolution of defect structures changes mechanical and physical properties of the material [10].

Electronic and atomic scale simulations are the indispensable numerical simulation tools that help understand the fundamental laws driving microstructure evolution and its effect on mechanical and physical properties of the materials. Ab initio density function theory (DFT) calculations [11], [12] are commonly used for computing the formation and migration energies of small defects. The energy of formation of a defect at equilibrium determines the relative probability of its occurrence, whereas the energy of migration determines the rate of thermal transformation of an already formed defect structure. Molecular dynamics [2], [3] and kinetic Monte Carlo [13], [14], [15] simulations provide information about reaction rates and relaxation pathways characterizing complex configurations of defects.

To avoid considering surface effects, simulations are often performed using periodic boundary conditions (PBCs). Through periodic boundary conditions, a spatially localized defect situated in a simulation cell interacts elastically with an infinite number of its own images situated in periodically translated simulation cells [6], [7], [8], [9]. Since elastic fields effectively have infinite range, and the energy Eel of elastic interaction between any two defects varies as the inverse cube of distance R between the defects EelR3 [8], if a relatively small cell is used in a simulation, the elastic energy of interaction between a defect and its periodic images can be substantial. This can affect the accuracy of calculations performed using PBCs and make the total energy data strongly dependent on the cell size. Although in principle the issue can be partially circumvented using a larger simulation cell, in practice this may not necessarily be a realistic option because of the limitations imposed by the available computational resources or numerical algorithms. For example, in a conventional DFT calculation, the simulation cell size is still limited to a few hundred atoms.

A possible way forward is to introduce an elastic correction to the calculated formation energy. A first order correction, assuming the linear elasticity approximation, can be derived using the elastic dipole tensor formalism [6], [7], [8], [9], which only requires knowing the elements Pij of elastic dipole tensor of the defect and the elastic constants tensor Cijkl of the material. This information can be readily derived from the same DFT or molecular statics calculation.

An elastic dipole tensor fully defines the elastic field produced by a defect in a material [8], [16], [17]. The elastic strain field associated with a localized defect, or even a large but still localized agglomerate of defects, can be expressed in an explicit analytical form using the notion of the dipole tensor. From the dipole tensor it is also possible to evaluate the relaxation volume tensor of a defect or an ensemble of defects [17], [18]. By considering all the defects in a certain volume element of the material as a compound object characterized by its dipole tensor, it is possible to formulate a continuum model spanning the spatial scale many orders of magnitude larger than an atomistic simulation [17]. In addition, the notion of the dipole tensor enables treating interactions between defects. A dipole tensor can be defined for an arbitrarily large configuration of defects, for example the entire defect structure created in a collision cascade simulation can be described by a dipole tensor, enabling extending the treatment to a macroscopic scale [17], [18].

In previous studies, we derived analytical equations for treating the elastic fields of defects in a simulation cell using periodic boundary conditions [8]. We have also derived equations for evaluating the elastic correction to the energy of a localized defect [9], and implemented them in our program CALANIE. It is appropriate to make this code, suitable for evaluating the elastic correction to the total energy, and for calculating the elastic dipole tensor of a defect in a simulation cell, available as an open source computer program. Full numerical and algorithmic aspects of the code are described below.

In what follows we review the fundamental theory and explain the meaning of various equations. We also discuss the numerical implementation of the method, followed by the details of the compilation procedure, and the format of input and output files. We give two examples illustrating applications of the code. The first example involves ab initio calculations of point defects and vacancy migration in FCC gold. This example illustrates the applicability of CALANIE to both equilibrium and non-equilibrium configurations. The second example illustrates molecular statics calculations of mesoscopic size dislocation loops. We investigate the numerical convergence of elements of the dipole tensor and formation energy of a defect as a function of the simulation box size, and the significance of applying elastic correction to the formation energy in the limit where the simulation cell is relatively small.

Section snippets

Elastic dipole tensor

In a continuum elasticity theory, the elastic strain energy of a defect in an infinite medium is defined as a volume integral over the entire space: ED=12Vσij(r)ϵij(r)dV,where ϵij and σij are the coordinate-dependent elastic strain and stress fields. Assuming the validity of the linear elasticity approximation, we write σij=Cijklϵkl, where Cijkl is the rank four elastic constant tensor. The above equation now acquires the form ED=12VCijklϵkl(r)ϵij(r)dV.In the presence of infinitesimal

Algorithm

It is not practically feasible to compute EDDtotal and EDDcorr by summing up an infinite number of terms in the series. Provided that we include the same number of terms in both series, the sum of them, EDD, converges in the limit where the cutoff distance is sufficiently large [8]. Calculating EDDtotal is relatively straightforward, since we can evaluate the second derivative of elastic Green’s function numerically. The calculation of EDDcorr is somewhat more involved as it requires evaluating

Compilation of the program

CALANIE is a code written in C++. It can be compiled using any modern C++ compiler, including Intel and GNU compilers. No linking to external libraries is required. The code can be compiled in two different ways, for two different purposes. The first one is for general type ab initio calculations. Using g++, one can compile CALANIE using the following command line

$ g++ -DABINITIO -DSTRESSeV -o calanie CALANIE_2.0.cpp

or

$ g++ -DABINITIO -DSTRESSGPa -o calanie CALANIE_2.0.cpp

Option -DABINITIO

Input and output

CALANIE uses two input files. They are input_data and input_elastic. These files need to be located in the same directory in order to execute CALANIE. Both are ASCII files.

When we use option -DABINITIO, in the input_data file we need to specify the translation vectors, the linear scaling factor, and the residual stresses in the perfect cell and in the cell containing a defect. They should be specified using the following format

box_ref_11 ???

 box_ref_12 ???

 box_ref_13 ???

 box_ref_21 ???

 box_ref_22

Ab initio calculations: Point defects in FCC gold

Elastic correction can be readily applied in the context of a calculation of formation and migration energies of point defects. We have applied CALANIE to improve the quality of ab initio defect energies in FCC gold, which were partially described in a study by Hofmann et al. [30]. The calculations were performed for a vacancy and several self-interstitial atom (SIA) defects, where the latter included a 100 dumbbell, an octahedral site interstitial, a 110 crowdion, and a 110 dumbbell.

All

Conclusion

In this study, we presented a summary of the fundamental theory, algorithms and numerical implementation of computer program CALANIE, intended for the evaluation of anisotropic elastic interaction energy correction associated with the use of periodic boundary conditions (PBCs). The theory is based on the linear elasticity approximation. The elastic interaction of a defect with its periodic images is approximated and evaluated using the elastic dipole and elastic Green’s function formalism. The

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

This work has been carried out within the framework of the EUROFusion Consortium and has received funding from the Euratom research and training programmes 2014–2018 and 2019–2020 under grant agreement No. 633053 and from the RCUK Energy Programme [grant No. EP/T012250/1]. To obtain further information on the data and models underlying the paper please contact [email protected]. The views and opinions expressed herein do not necessarily reflect those of the European Commission.

References (45)

  • SandA.E. et al.

    J. Nucl. Mater.

    (2014)
  • ClouetE. et al.

    Acta Mater.

    (2008)
  • JourdanT.

    J. Mech. Phys. Solids

    (2019)
  • KresseG. et al.

    Comput. Mater. Sci.

    (1996)
  • MillsG. et al.

    Surf. Sci.

    (1995)
  • PlimptonS.

    J. Comput. Phys.

    (1995)
  • LiJ.

    Modelling Simulation Mater. Sci. Eng.

    (2003)
  • CaiW. et al.

    Imperfections in Crystalline Solids

    (2016)
  • SandA.E. et al.

    Europhys. Lett.

    (2016)
  • ByggmästarJ. et al.

    J. Phys.: Condens. Matter

    (2019)
  • DudarevS.L. et al.

    Phys. Rev. B

    (2010)
  • VarvenneC. et al.

    Phys. Rev. B

    (2013)
  • VarvenneC. et al.

    Phys. Rev. B

    (2017)
  • DudarevS.L. et al.

    Phys. Rev. Mater.

    (2018)
  • MaP.-W. et al.

    Phys. Rev. Mater.

    (2019)
  • Durrand-CharreM.

    Microstructure of Steels and Cast Irons

    (2003)
  • KohnW. et al.

    Phys. Rev.

    (1965)
  • HohenbergP. et al.

    Phys. Rev.

    (1964)
  • MasonD.R. et al.

    J. Phys.: Condens. Matter

    (2014)
  • SubramanianG. et al.

    Phys. Rev. B

    (2013)
  • MartínezE. et al.
  • LeibfriedG. et al.

    Point Defects in Metals

    (1978)
  • Cited by (16)

    • Materials for in-vessel components

      2022, Fusion Engineering and Design
      Citation Excerpt :

      Still, knowing the volumes of defects is essential for modeling radiation effects on the macroscopic scale these volumes determine the magnitude of radiation-induced swelling or contraction [66] of materials, resulting from the accumulation of radiation defects. Relaxation volumes of a vacancy and several high symmetry configurations of self-interstitial atom defects in some selected bcc metals, computed using density functional theory [67–69], are summarized in Table5. The results are still sensitive to the choice of exchange-correlation functional, but the margin of error is very small in comparison with experimental measurements.

    • Elastic dipole tensors and relaxation volumes of point defects in concentrated random magnetic Fe-Cr alloys

      2021, Computational Materials Science
      Citation Excerpt :

      Since most of the results for Fe-Cr alloys were obtained using fixed volume simulations cell, defined by the lattice parameter of pure Fe a = 2.831 Å, all the results for Cr given in this study were also computed assuming this lattice parameter. The computed formation energies include the correction term resulting from the periodic boundary conditions and the requirement of vanishing average strain [14,45]. The formation energies of approximately 300 vacancy and 400 dumbbell configurations derived from fixed cell volume DFT simulations are shown, as functions of Cr content and the number of Cr atoms in the local environment of a defect, in Fig. 4a–d.

    View all citing articles on Scopus

    This paper and its associated computer program are available via the Computer Physics Communication homepage on ScienceDirect (http://www.sciencedirect.com/science/journal/00104655)

    ☆☆

    The review of this paper was arranged by Prof. D.P. Landau.

    View full text