Importance of elastic finite size effects: neutral defects in ionic compounds

Importance of elastic finite size effects: neutral defects in ionic compounds

P.A. Burr School of Electrical Engineering and Telecommunications, University of New South Wales, Kensington, 2052, NSW, Australia.    M.W.D. Cooper Materials Science and Technology Division, Los Alamos National Laboratory P.O. Box 1663, Los Alamos, NM 87545, USA.
July 31, 2019

Small system sizes are a well known source of error in DFT calculations, yet computational constraints frequently dictate the use of small supercells, often as small as 96 atoms in oxides and compound semiconductors. In ionic compounds, electrostatic finite size effects have been well characterised, but self-interaction of charge neutral defects is often discounted or assumed to follow an asymptotic behaviour and thus easily corrected with linear elastic theory. Here we show that elastic effect are also important in the description of defects in ionic compounds and can lead to qualitatively incorrect conclusions if inadequatly small supercells are used; moreover, the spurious self-interaction does not follow the behaviour predicted by linear elastic theory. Considering the exemplar cases of metal oxides with fluorite structure, we show that numerous previous studies, employing 96-atom supercells, misidentify the ground state structure of (charge neutral) Schottky defects. We show that the error is eliminated by employing larger cells (324, 768 and 1500 atoms), and careful analysis determines that elastic effects, not electrostatic, are responsible. The spurious self-interaction was also observed in non-oxide ionic compounds and irrespective of the computational method used, thereby resolving long standing discrepancies between DFT and force-field methods, previously attributed to the level of theory. The surprising magnitude of the elastic effects are a cautionary tale for defect calculations in ionic materials, particularly when employing computationally expensive methods (e.g. hybrid functionals) or when modelling large defect clusters. We propose two computationally practicable methods to test the magnitude of the elastic self-interaction in any ionic system. In commonly studies oxides, where electrostatic effects would be expected to be dominant, it is the elastic effects that dictate the need for larger supercells — greater than 96 atoms.

I Introduction

Finite size effects have been a known limitation since the beginning of atomic scale simulations of solids Catlow and Mackrodt (1982); Freysoldt et al. (2014). These arise when atomic interactions extend beyond the simulation boundaries. Most modern atomic scale simulations methods adopt periodic boundary conditions (PBC) to represent crystalline matter and introduce defects through the use of supercells Leslie and Gillan (1985). Of primary concern for finite-size effects are long range interactions through elastic (strain) fields and electrostatic (Coulomb) fields. Elastic self-interactions of point defects have been studied since the early days of atomic scale simulations, first for simple elemental metals Eshelby (1955); Hardy and Bullough (1967); Gairola (1978); Hardy (1968), and later for ionic crystals Gillan (1984); Leslie and Gillan (1985). This was dictated by necessity, as the exhisting computational resources limited the size of force-field simulations to tens of atoms — now atoms can be modelled Eckhardt et al. (2013). While the computational power available to atomic scale modellers has increased dramatically, the typical simulation size used for density functional theory (DFT) calculations has not increased accordingly, in favour of ever increasing sophistication in the description of the electronic state. However, much of the old knowledge regarding the importance of accounting for elastic self-interactions seems to be lost or ignored in many recent reports.

In ionic materials, where typically point defects are charged, Coulomb interactions are assumed to be predominant. Consequently, a large body of research has focussed on predicting and countering the electrostatic self-interaction energy in DFT, producing a number of charge correction schemes with increasing degrees of complexity and sophistication Makov and Payne (1995); Lany and Zunger (2009); Freysoldt et al. (2009); Hine et al. (2009); Freysoldt et al. (2011); Taylor and Bruneval (2011); Murphy and Hine (2013); Kumagai and Oba (2014). On the other hand, elastic self-interactions, which are thoroughly accounted for in metals and iono-covalent materials Freysoldt et al. (2014); Leslie and Gillan (1985); Varvenne et al. (2013); Jain et al. (2016), have largely been neglected in strongly ionic compounds, being perceived of secondary importance to electrostatic interactions. Here we show that for charge neutral defects, the elastic interactions are non-negligible and lead to a qualitative change in defect stability.

We consider exemplar cases of metal oxides with fluorite structure (CeO, ThO, UO and actinide oxides) that have been extensively studied in the past due to their engineering applications including solid oxide fuel cells (SOFC), electrolyzer cells, ion conductors, catalysts and nuclear fuel. Specifically we consider the formation of charge-neutral Schottky clusters (), which are known to reduce oxygen mobility Genreith-Schriever et al. (2015); Zacate et al. (2000); Grieshammer et al. (2014, 2016), degrade the electrolytic properties of SOFC Grieshammer et al. (2016); Li et al. (2011); Kim et al. (2012) and govern the distribution and retention of gaseous fission products in nuclear fuel Chroneos (2016); Bertolus et al. (2015); Ichinomiya et al. (2009); Cooper et al. (2014a); Liu et al. (2012). We also extend the study to non-oxide ionic compounds with fluorite structure (CaF) and anti-fluorite structure (BeC, also ionic Herzig and Redinger (1985); Tzeng et al. (1998)) to show that the phenomenon is not limited to oxides.

The formation and migration of charge-neutral clusters in these metal oxides has been extensively investigated with ab-initio simulations Grieshammer et al. (2016); Li et al. (2011); Liu et al. (2012); Shao et al. (2017); Murgida et al. (2014); Wang et al. (2014); Zacherle et al. (2013); Keating et al. (2013); Liu and Sickafus (2011); Nakayama and Martin (2009); Xiao et al. (2011); Gupta et al. (2009); Keating et al. (2012); Thompson and Wolverton (2011); Crocombette (2012); Crocombette et al. (2011); Dorado et al. (2013); Geng et al. (2008); Vathonne et al. (2014); Huang (2014), but nearly all of the studies were carried out with simulation cells containing up to 96 atoms ( supercell) with few instances where larger supercells were used Vathonne et al. (2014); Bradley et al. (2015); Kuganathan et al. (2017); Murphy et al. (2014). Of the three possible configuration that a Schottky cluster may exhibit (see figure 1), the DFT studies consistently report the the SD configuration as the lowest energy cluster. One exception is the publication by Yu et al.  Yu et al. (2009), however their reported Schottky energies appear to be orders of magnitude smaller than all other published work. On the other hand, modelling studies employing empirical force-fields, and thus frequently using very large periodic simulations cells or the Mott-Littleton approach Mott and Littleton (1938), often report the inverse: that SD is more favourable than SD Govers et al. (2007); Aidhy et al. (2009); Cooper et al. (2014b) — although, there is large variation in results as the quality of the potentials precludes the apparent stability of Schottky clusters. It is unclear to what extent the discrepancy between empirical and ab-initio results is due to the different level of theory (DFT v.s. force-fields) or to finite size effects. In this paper we show that DFT and reliable force-field potentials are indeed in agreement if the results are compared across the same supercell sizes. Notably, two studies employed supercell (324 atoms) of ThO Kuganathan et al. (2017); Murphy et al. (2014). Murphy et al. Murphy et al. (2014) report similar trends between DFT and force-field simulations for increasing supercell size, however the discrepancy regarding the most favourable cluster configuration still remained. Although not belonging to the family of fluorite compounds, Bradley et al.  Bradley et al. (2015) also employed a 324-atoms supercells of m-HfO to study defect clusters.

Figure 1: Three possible configurations of the bound SD in the fluorite structure, defined by the arrangement of the anion vacancies around the cation vacancy.

In this letter, we show that charge neutral Schottky clusters interact over long ranges through elastic fields, and that the commonly used 96-atom supercell is inadequately small to capture the correct ground state of Schottky clusters. We then show that the root cause of this finite size effect is spurious interaction between strain fields across PBC, and not electrostatic or electronic effects. Additionally, the elastic finite size effect is not described accurately by linear elastic theory, owing to the complex mode of relaxation. We conclude by proposing two methods to estimate, albeit approximately, the magnitude of the elastic self-interaction for any ionic system.

Ii Methodology

Force-field simulations were performed with the LAMMPS code Plimpton (); Plimpton (1995), using the many-body potential of Cooper, Rushton and Grimes (CRG) Cooper et al. (2014b), as this potential set proved to be reliable and transferable Ghosh et al. (2016); Liu et al. (2016); Cooper et al. (2016); Yuan et al. (2015); Rushton and Chroneos (2014) across a wide range of MO compounds. DFT simulations were carried out with VASP Kresse and Furthmüller (1996), using the PBE exchange-correlation functional Perdew et al. (1996), PAW pseudo-potentials with the maximum number of valence electrons available Kresse and Joubert (1999) and a plane-wave cut-off of . Details of the k-point grids and evidence of convergence to within is provided in the supplementary materialSupplMater (). -ZrO was also initially investigated, but due to the stability of -ZrO it was not possible to retain the cubic symmetry during relaxation of defects in large cells, even by constraining some degrees of freedom.

On-site Coulomb correction terms have been included for CeO and UO, following the majority of the published literature Allen and Watson (2014); Keating et al. (2012); Andersson et al. (2007); Wiktor et al. (2017); Vathonne et al. (2014); Dorado et al. (2013); Dorado and Garcia (2013); Thompson and Wolverton (2011); Crocombette (2012); Crocombette et al. (2011); Dorado et al. (2010); Geng et al. (2008); Iwasawa et al. (2006); Kotani and Yamazaki (1992); Yamazaki and Kotani (1991): Dudarev et al. ’s formalism Dudarev et al. (1998) for CeO with , ; and Liechtenstein’s formalism Liechtenstein et al. (1995) for UO with and . U-ramping Meredig et al. (2010) (for CeO) and occupation matrix control Allen and Watson (2014); Dorado et al. (2010) (for UO) were used to avoid metastable states. Details are provided in the supplementary material SupplMater () together with results obtained without +U to emphasis that the findings are not sensitive to the choice of simulation parameters. UO was described with collinear 1k antiferromagnetic ordering, as this is the best collinear approximation of the true (non-collinear 3k AFM Desgranges et al. (2017); Ikushima et al. (2001)) magnetic ordering of UO Wiktor et al. (2017); Geng et al. (2008); Vathonne et al. (2014); Dorado et al. (2013); Dorado and Garcia (2013); Crocombette (2012); Crocombette et al. (2011); Dorado et al. (2010).

In ionic materials, the formation energy of a defect with charge is conventionally calculated as


where and is the DFT total energies of the defective or pristine supercells, is the chemical potential of all species added or removed to form the defect, is the fermi level of the system normalised by the valence band maximum, is a charge correction term, following a number of possible schemes Makov and Payne (1995); Lany and Zunger (2009); Freysoldt et al. (2009, 2011); Taylor and Bruneval (2011); Hine et al. (2009); Murphy and Hine (2013); Kumagai and Oba (2014), and is the energy due to elastic self-interaction, seldom accounted for in ionic compounds. Since (a) the defects considered here are charge neutral, (b) the composition of the defects is precisely one stoichiometric formula unit (i.e. , where the number of formula units in the supercell), and (c) the elastic self-interaction is the subject of the study, the defect formation energy of a Schottky cluster is simplified to:


in line with previous publications Thompson and Wolverton (2011); Crocombette (2012); Crocombette et al. (2011) (except for the explicit inclusion of elastic self-interaction term). This simplification conveniently removes any dependence of our results from external factors such as chemical potential of reference elements, or apparent bang gap of the material.

The linear elastic theory approximation to was calculated with the aid of the aneto script Varvenne et al. (2013) from the stress tensor of the relaxed simulations and using elastic constants obtained from DFT simulations through lattice perturbation to retain self-consistency. Selected calculations were repeated where atomic relaxation was restricted to atoms within a relaxation radius from the defect centre, and all other atoms were kept fixed at the perfect lattice site.

Iii Results and discussion

The effect of finite PBC was first investigated using the CRG potential. Figure 2 shows the formation energy (filled symbols) of bound Schottky clusters in various actinide oxides and CeO. It is clear that, irrespective of the cation species, the SD defect is the most favourable when simulated in large enough supercells (containing 324 atoms or more), but the SD is the most favourable in the smaller simulation cell containing 96 atoms. The results of calculations before geometry relaxation is also presented in the top panel, showing no cross-over of defect energies.

Figure 2: Defect formation energy (before and after relaxation) from CRG potentials versus supercell size (form to ) for three bound Schottky configurations in actinide oxides and CeO. Values for ThO were taken from Murphy et al. (2014).

To show that the crossover was not a peculiarity of the CRG potential form, DFT calculations were performed with supercells containing 96, 324 and 768 atoms on selected oxides (CeO, ThO, UO) as well as CaF and BeC (Figure 3). BeC has considerably small lattice parameter, therefore a further supercell containing 1500 atoms was also considered. It is evident that DFT and force-field calculations are in agreement, and crossover between SD and SD occurs between the supercells (96 atoms) and the supercells (324 atoms). The trend is predicted for CaF and BeC as well as the oxides. Importantly, including the energy penalty predicted from linear elastic theory (dashed lines) does not correct the trend.

Figure 3: Defect formation energies from DFT as a function of supercell size (labels indicate atoms in supercell). UO simulations were limited to 324 atoms, due to the additional complexity and computational cost of OMC.

It is well known that point defects in ionic materials interact chiefly though their charges, and although the Schottky clusters have no overall charge, they are still comprised of three charged defects () that may individually interact across periodic boundaries. In addition, SD and SD also exhibit effective dipoles since the geometrical centre of the positive charges does not align with that of the negative charges. SD does not have an associated dipole since it is a linear defect with mirror symmetry. Nevertheless, to show that electrostatic interactions alone cannot account for this peculiar finite size effect, all charge-charge, charge-dipole, and dipole-dipole interactions have been evaluated independently (Figure 4). Three point charges (two positive and one negative with no overall charge) have been arranged in a dielectric medium with the same configurations as the three bound Schottky defects, and then replicated in a repeating array of to model the effect of PBC, where the distance between replicas was increased progressively to simulate larger supercells. Dielectric constants, lattice parameter and the magnitude of charges are arbitrary. However, a range of ratios have been modelled, all yielding the same qualitative behaviour.

Figure 4: Total electrostatic energy of SD modelled as point charges arranged in a dielectric medium with periodic boundaries. Insert: energy contribution from self-interaction across supercell boundaries.

Figure 4 shows that as the supercell size increases, the energy contribution from self-interaction across periodic boundaries tend to zero (see insets), therefore the total energy of each system converges toward the internal energy of the point charge triplet, consisting entirely of Coulomb interactions. It is clear that the three point charge configurations (representative of the three Schottky clusters) never cross over at any separation, hence, the electrostatic interactions alone do no account for the change in relative stability of Schottky clusters. This is reassuring, given that CaF shows similar trend to CeO, ThO and UO despite a factor of difference in ionic charges (Figure 3).

Beyond electrostatic effects, electronic effects may also lead to a change in behaviour with increasing supercell size. This has been well documented for charge neutral vacancies in Si Mercer et al. (1998); Puska et al. (1998); Centoni et al. (2005), where small simulation cells (64-atoms or fewer) predict a retention of symmetry, while larger cells correctly identify a reduction of local symmetry to , as observed experimentally. The symmetry reduction, and associated energy reduction, is due to a Jahn-Teller distortion Jahn and Teller (1937), and consequently the effect is strongly sensitive to Brillouin zone sampling as well as supercell size Puska et al. (1998). SiC is another example of covalent material where similar effects have been observed Zywietz et al. (1999). Jahn-Teller effects are not limited to covalently bonded materials, in fact they are known to be important in many oxides Dorado et al. (2010); Elfimov et al. (1999). However, the finite-size effect observed in the current study was consistently reproduced with GGA, GGA+U and force field potentials (with fixed charges on atoms), suggesting that electronic effects cannot be at the heart of the matter.

The source of the cross-over is instead found in the elastic interactions. Figures 2 shows that prior to geometry relaxation there is no cross-over in stability between SD and SD. This is also shown for DFT calculations in the supplementary material SupplMater2 (). Thus relaxation is hampered in the smallest supercell. Figure 5 depicts the atomic displacements caused by the three bound Schottky clusters in the largest DFT cell of CeO. The boundaries of the smaller supercells are superimposed on the image to highlight that the strain field exceeds the bounds of the supercell (96 atoms). The displacement fields are better quantified in Figure 6, where the atomic displacements are plotted as a function of distance from the cation vacancy. The atomic displacements obtained from the largest supercell, reveal that at a distance of from (i.e. between the and boundaries) O atoms are being displaced by as much as , which is not insignificant. More importantly, when comparing the atomic displacements within the first across different supercell sizes, it is evident that the fingerprint of atomic relaxation in the 96-atom supercell is fundamentally different from that of the larger supercells, i.e. not only the magnitude but also the shape of the strain field is different. On the other hand, the fingerprint of atomic displacements changes only marginally when increasing supercell size further.

Figure 5: (110) slice (three atomic layers) of the 768-atoms supercell of CeO. Dashed lines represents the boundaries of the smaller supercells. Grey and red squares represent the Ce and O vacancies respectively, green arrow represent the displacement magnitude scaled up by a factor of 5. First nearest neighbour atoms exceed the colour map with displacements of up to .
Figure 6: Atomic displacements caused by Schottky clusters in CeO supercells. Red and turquoise crosses represent O and Ce atoms respectively. Black vertical lines represent half of the supercell length, following the same key of Figure 5. Beyond those lines, some atoms are closer to the periodic replica of the defect than the central defect.

The fact that the shape of the displacement fields are fundamentally different between the 96-atom supercell and the larger supercells indicates that atomic relaxation is hampered by artificial restoring forces stemming from the PBC. This frustration of atomic relaxation can only provide a positive contribution to the total energy of the system.

Generally, it is possible to predict the energy contribution arising from elastic self-interaction through methods based on linear elastic theory, which combining the elastic dipole tensor of the simulation with the elastic constants of the material Leslie and Gillan (1985); Varvenne et al. (2013). However, Figure 3 shows that the aneto correction, which has been proven successful in a variety of metallic and covalent systems Varvenne et al. (2013); Burr et al. (2015); Agarwal and Trinkle (2016); Garnier et al. (2014); Murali et al. (2015); Varvenne et al. (2014), cannot counter the finite size effects observed here.

The inability of linear elastic theory to predict the energy contribution in ionic compounds is attributed to the complex relaxation pattern caused by the defect. Figure 5 clearly shows that the pattern of atomic displacements is not consistent with a simple compression (acoustic) wave, where all atoms move coherently towards the vacancies to accomodate the defect volume. The relaxation field more closely resembles of an optical mode, where anions and cations exhibit distinct and more complex displacement patterns. This is evidenced more clearly in Figure 7, where the radial component of the displacement vector is compared against that of an isostructural metallic compound (AlAu). In the metallic compound, the atomic relaxations are nearly exclusively towards the vacancy cluster (negative radial displacements), and the effect diminishes nearly monotonically with distance. In the ionic compound, each shell of anions and cations exhibit both inwards and outward displacements (with the exception of the first oxygen shell). As with any optical mode, this behaviour may only arise in compound materials, explaining why it has not been observed in well studied elemental materials. In addition, Figure 7 shows that the phenomenon is restricted to ionic compounds, where alternating shells of anions and cations move in opposite directions in response to displacement of charge (e.g. neighbouring ions).

Figure 7: Radial displacement, , caused by Schottky clusters in CeO and AlAu. , where is the displacement vector and is the position vector with respect to the cation vacancy.

The spurious interaction results in a change in shape of relaxation near the defect core, which cannot be captured by the dipole tensor: the dipole tensor may be obtained from (a) the (anisotropic) stress tensor and/or strain tensors on the cell, or (b) by convergent summation of atomic displacements and/or restoring forces on atoms Nazarov et al. (2016). In both cases the dipole tensor may only capture a truncation of the strain field, and not a change in shape of the core of the field (i.e. the “optical relaxation”), which is illustrated particularly clearly in Figure 6 by comparison between the and supercells for which there is a notable change in relaxation displacements for atoms even within of the defect. These core atoms are away from the defect replica in the supercell.

As further evidence that the source of the self-interaction energy is the inhibition of elastic relaxation, we have performed calculations where only atoms within a given radius of the defect centre were allowed to relax while all other atoms were “clamped” at perfect lattice positions, see Figure 8. Provided that the relaxation radius is less then half the supercell length, , this constraint ensures that no strain is transmitted across the PBC. Thus the relaxation field, although artificially truncated at the given radius, is entirely due to the defect and not its periodic replicas.

Figure 8: Formation energies of Schottky defects in CaF as a function of relaxation radius, whereby all atoms beyond that radius are forcefully kept fixed in their perfect lattice sites. Lower panels provide an enlarged view of the crossover region. Vertical black line represent the largest radius that may be accommodated in the supercell without overlap . First and last point in each plot represent a single point calculation and a complete relaxation, respectively.

When the relaxation radius is less then half the supercell length (i.e. the region of atoms that are allowed to relax fits entirely within the supercell), SD appears less favourable than SD, consistent with the findings obtained from very large supercell simulations. Only when the relaxation radius is larger than half the cell length (i.e. atoms within the sphere are responding to the strain field of the defect and its replicas), does the apparent cross-over in stability manifest in the 96-atom supercell. This is a further confirmation that the phenomenon is due to the interaction of strain fields across PBC in the smaller supercells.

Constraining atomic positions as a function of distance from the defect allows one to isolate the relaxation strain component of the defect energy. And by iterating that procedure over increasing relaxation radii, it is evident that the behaviour is distinctly not monotonic, with two clear steps at the and and nearest neighbour. This is at odds with the behaviour expected from a continuum elastic medium.

It is concluded that if inadequately small supercells are used, 96 atoms or less for fluorite structure oxides, the spurious self interaction is beyond the reach of linear elastic theory, and thus cannot be corrected for without further calculations. It is not always possible to increase the simulation size, especially when performing calculations with computationally expensive methods, such as hybrid functionals, ab-initio molecular dynamics or time-dependent DFT. Thus, here we propose two computationally practicable methods to ascertain whether finite size effects are significant in the system of interest:

  1. One option is to use a lower level of theory method to perform the convergence test. In the current work, we have shown that CRG force-field potential was as predictive as DFT in the analysis of the spurious self-interaction. Thus, if reliable force-field potentials are for the system of interest, these may be used to test supercell convergence. Similarly, for hybrid calculations, one could use LDA or GGA functionals to test the supercell size convergence. Note however, that this is a two-stage approach: in the first instance one must test that within the same supercell size the two methods yield similar trends. Quantitative agreement is not expected, but any qualitative agreement between the two methods observed in the small simulation size is expected to be preserved in larger simulations cells (as shown in the current work). Note that this method cannot be extended to electronic self-interaction, as these are highly sensitive to the description of the electronic states Centoni et al. (2005).

  2. The alternative approach is to use a range-dependent constraint method with , as illustrated in Figure 8, to test whether there is qualitative agreement between the fully relaxed simulations and those where the strain fields were not allowed to transfer across the PBC. Incidentally, this method also significantly reduces the computational cost of energy minimisation on account of the reduced internal degrees of freedom. The size of the discrepancy may be considered as an indicator of the extent to which strain extends beyond the PBCs.

Iv Conclusions

Electrostatic self-interactions have long been thought to be the dominant source of finite size effects in ionic materials. Here we have shown that even in ionic compounds, elastic self-interactions are non-negligible, causing qualitative and quantitative changes to the energy and structure of defects. This finite size effect is the source of numerous inaccurate reports in the ab-initio literature regarding the defect stability of many important functional oxides including CeO, ThO and UO. The magnitude of the elastic self-interaciton is such that it warrants the need for simulations sizes larger than the widely used 96-atoms supercell for cubic oxide compounds. In addition, we show that the spurious self-interaction cannot be countered by simple linear elastic theory approaches.

We considered the exemplar cases of charge-neutral Schottky clusters in fluorite-structured oxides and perform DFT and force-field simulations in supercells of increasing size up to 1500 atoms. Supercells of 96 atoms were not sufficiently large to capture the most stable neutral Schottky clusters compared to larger cells. This behaviour is observed also for non-oxide ionic compounds with related structures (prototypical fluorite CaF, and antifluorite BeC), but was not observed in iso-structural metallic compounds. This phenomenon is also insensitive to the level of theory used (EAM, DFT, DFT+U).

We provide evidence that the finite size effect is not due to electrostatic self-interaction (accounting for charge-charge, charge-dipole and dipole-dipole interactions), nor electronic effects (owing to the presence of the self-interaction in force-field calculations with fixed charges) instead it is caused by spurious interaction of strain fields across PBC. This is confirmed through the use of range-constrains, where selected atoms were fixed to the perfect lattice sites, thereby artificially truncating the strain field before the PBC. In these simulations the correct order of defect energies was restored. Surprisingly, the spurious self-interaction energy cannot be accounted for, and therefor corrected with, linear elastic theory. The failure of linear elastic theory is attributed, through careful analysis of the displacement fields, to the complex relaxation patterns observed in the ionic compounds, akin to optical modes.

Our findings also resolve a long-standing discrepancy between DFT and force-field simulations regarding the ground state of neutral defect clusters in actinide oxides. This had typically been attributed to the simplified atomic interactions of the empirical force-field methods, but we show that good agreement is in fact obtained when comparing across the same supercell sizes. The fact that the phenomenon is observed at all levels of theory can be exploited to the practitioners’ advantage, using computationally simpler methods to test the degree of elastic self-interaction in the system of interest.

Although we have focussed primarily on bound Schottky defects, the findings are likely relevant to other neutral clusters and their migration pathways, such as and in yttria-stabilised-zirconia and gadolinia-doped ceria. Rather than just being considered a spurious size effect for dilute limit calculations, the elastic effects discussed here also apply to real materials with high concentrations of Schottky defects, providing an insight into the interaction between two or more Schottky clusters that come within 10–16 of each other. This may be an important factor in the nucleation of voids in nuclear fuels and SOCF.

The current work should serve as a cautionary tale for future simulations of all ionic systems, especially those where computational requirements dictate the use of small supercell sizes. Whilst proving that 96-atom supercells are inadequate for point defect analysis in fluorite-structured oxides, every solid state system would have different supercell requirements. In the current paper we provide two computationally efficient methods to ascertain whether elastic finite-size effects are significant for the system of interest: 1) comparison with lower level of theory (as exemplified here using force-field CRG potentials), 2) comparison with range-constrained simulations (where selected atoms have been fixed to the perfect lattice sites to hamper propagation of strain fields).

The authors would like to acknowledge S. T. Murphy, N. Kuganathan, M. J. D. Rushton and R. W. Grimes for the useful discussion. M. W. D. Cooper is funded by the US Department of Energy, Office of Nuclear Energy, Nuclear Energy Advanced Modeling Simulation (NEAMS) program. Los Alamos National Laboratory, an affirmative action/equal opportunity employer, is operated by Los Alamos National Security, LLC, for the National Nuclear Security Administration of the U.S. Department of Energy under Contract No. DE-AC52- 06NA25396. P. A. Burr acknowledges the Australian Nuclear Science and Technology Organisation and the Tyree Foundation for their charitable financial support.


Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description