Improving the Functional Control of Aged Ferroelectrics using Insights from Atomistic Modelling

TL;DR

This paper uses atomistic modeling to understand and control aging in ferroelectric materials, revealing how defect interactions influence their properties and offering methods to optimize their performance.

Contribution

It demonstrates that intrinsic defect interactions can explain aging phenomena in ferroelectrics and proposes a universal approach to improve their electromechanical properties.

Findings

Aging phenomena can be reproduced by defect-dipole interactions in simulations.

Dopant concentration affects hysteretic response.

A method to reduce loss and tune properties of ceramics is identified.

Abstract

We provide a fundamental insight into the microscopic mechanisms of the ageing processes. Using large scale molecular dynamics simulations of the prototypical ferroelectric material PbTiO3, we demonstrate that the experimentally observed ageing phenomena can be reproduced from intrinsic interactions of defect-dipoles related to dopant-vacancy associates, even in the absence of extrinsic effects. We show that variation of the dopant concentration modifies the material's hysteretic response. We identify a universal method to reduce loss and tune the electromechanical properties of inexpensive ceramics for efficient technologies.

Full text

Technologies utilising ferroelectric components are ubiquitous in modern devices, being used from mobile phones, diesel engine drive injectors and sonar to print heads and non-volatile memory Genenko et al. (2015); Glaum et al. (2012); Scott and Paz de Araujo (1989); Catalan et al. (2012). Doping with transition-metals has been shown experimentally to improve electromechanical properties of widely used ferroelectrics. For example, doping of BaTiO3, PbTiO3, Pb(Zr1-xTix)O3 (PZT) and (1-x)Pb(Mg1/3Nb2/3)O3-xPbTiO3 (PMNPT) is used to improve the functional properties and efficiency of these simple and cheap oxides Neumann and Arlt (1987); Arlt and Neumann (1988); Ren and Otsuka (1997, 2000); Ren (2004); Zhang and Ren (2005); Zhang et al. (2008a); Wu et al. (2008); Perez-Delfin et al. (2011). However, the fundamental origin of the electromechanical improvements is not understood and requires full characterisation to enable properties to be directly tuned for purpose and functional lifetimes to be accurately predicted.

Dopant interactions can be classified as intrinsic (bulk/volume) or extrinsic (boundary). Extrinsic coupling is associated with domain wall and grain boundary effects. Defects, including dopants and vacancies, migrate to domain walls and subsequently pin their propagation, resulting in fatigue of the material’s switching properties Carl and Härdtl (1978). Intrinsic effects occur independent of interaction with domain walls, as they arise due to the interaction between defect induced dipoles $\vec{p}_{d}$ and the spontaneous polarisation of the domain surrounding the defect site $\vec{P}_{s}$. Strong evidence from electron paramagnetic spin resonance (ESR) and density functional theory (DFT) calculations has shown dopants/impurities, such as iron Fe3+ and copper Cu2+, in PZT (or Mn2+ in BaTiO3) substitute the B-cations as acceptors, which bind to charge compensating oxygen vacancies $V_{O}^{2-}$ to form thermodynamically stable defect complexes Erhart et al. (2007); Eichel et al. (2008). In Kröger-Vink notation, the divalent dopant-vacancy associates can be written as ($B^{{}^{\prime\prime}}_{Ti}+V^{\bullet\bullet}_{O})^{\times}$, where $B^{{}^{\prime\prime}}_{Ti}$ is an unspecified divalent dopant subsituting a Ti4+ site, $V_{O}^{\bullet\bullet}$ is an oxygen vacancy with a +2e charge relative to the defect free site, ′ identifies a negative charge unit (-e), $\bullet$ represents a positive charge unit (+e) and $\times$ stands for charge neutrality. Density functional theory calculations have shown defect-dipoles $\vec{p}_{d}$ spontaneously form for ($Fe^{{}^{\prime}}_{Ti}+V^{\bullet\bullet}_{O})^{\bullet}$ Meštricć et al. (2005) and ($Cu^{{}^{\prime\prime}}_{Ti}+V^{\bullet\bullet}_{O})^{\times}$ Eichel et al. (2008) associates in PbTiO3 and ($Mn^{{}^{\prime\prime}}_{Ti}+V^{\bullet\bullet}_{O})^{\times}$ in BaTiO3 Zhang and Ren (2005); Nossa et al. (2015), and the energetically favourable oriention is along the polar axis [001]. Group-IIIB and group-VB acceptor substitutes on Ti sites in PbTiO3 have been shown to form immobile clusters of dopant-vacancy associates which have different structures when the associate is aligned parallel or perpendicular to the polar axis Zhang et al. (2008b). ($V_{Pb}^{{}^{\prime\prime}}+V_{O}^{\bullet\bullet})^{\times}$ divacancy complexes in PbTiO3 have been calculated to have a local dipole moment twice the bulk value Cockayne and Burton (2004).

Ageing is simply defined as the change in a material’s properties over time. It has been proposed that in aged ferroelectrics, defect-dipoles produced from dopant-vacancy associates will slowly rotate to align in parallel with the domain symmetry to minimise its energy state Lambeck and Jonker (1986); Ren (2004); Kimmel et al. (2012). The co-alignment and subsequent correlated behaviour of these aged defect-dipoles has been proposed to create a macroscopically measurable internal bias, which in turn has been conjectured to be responsible for experimentally observed ageing phenomena, including a 10-40 fold increase in piezoelectric coefficients, shifts in the hysteresis along the electric field axis and pinched/double hysteresis loops typically associated with antiferroelectrics Genenko et al. (2015); Ren (2004); Zhang and Ren (2005); Zhang et al. (2008a).

In this letter, we use large scale classical molecular dynamics to model ageing arising from defect-dipoles of dopant-vacancy associates in tetragonal bulk lead titanate (PbTiO3). We show that all the experimentally observed large signal effects (P-E and S-E hysteresis) of aged prototype perovskite ferroelectrics; pinched and double hysteresis, shifted hysteresis and a large recoverable electromechanical response can be reproduced from intrinsic effects alone and we identify the microscopic mechanisms of each case.

We study ideal and aliovalent-doped bulk PbTiO3 using classical molecular dynamics (MD) as implemented in the DL$\_$POLY code Todorov et al. (2006). We use the adiabatic core-shell interatomic potentials derived in Gindele et al Gindele et al. (2015) that reproduces the properties of bulk and thin films of PbTiO3 in excellent agreement with DFT calculations Gindele et al. (2015); Chapman et al. (2017). The prototype PbTiO3 has been chosen as it has a single ferroelectric phase, which reduces competing effects and because it is a parent compound for two of the most widely used ferroelectric materials in industry (PZT/PMNPT).

In this study we investigate volume effects, therefore, three-dimension periodic boundary conditions are implemented to mimic an infinite crystal, devoid of surfaces, interfaces and grain boundaries. We choose a moderate supercell constructed from $12\times 12\times 12$ unit cells, approximately 125 nm3, corresponding to 8,640 atoms (for the ideal bulk). This system size is large enough for ensemble sampling but sufficiently small to prevent the formation of 90*∘* domain walls. We use the Smooth Particle Mesh Ewald (SPME) summation for the calculation of Coulomb interactions. Coupling between strain and polarisation is enabled using the constant-stress Nosé-Hoover ($N\sigma T$) ensemble with thermostat and barostat relaxation times of 0.01 ps and 0.1 ps, respectively. A 0.2 fs timestep is used in all instances. Initial calculations were run at 100 K to prevent diffusion of the vacancies Morozov and Damjanovic (2010); Smyth (1994) and to allow the correct characterisation of each effect. The temperature dependence, for the range from 50 K to 400 K, is then investigated.

We calculate polarisation - electric field (P-E) hysteresis using a quasistatic approach. Starting at 0 kV/mm, the electric field is cycled between the limits $\pm 150$ kV/mm in 16.7 kV/mm intervals. For each field strength the system is restarted using the coordinates, velocities and forces from the previous calculation and equilibrated for 4 ps to enable the system to equilibrate following the E-field impulse. This is followed by an 8 ps production run over which statistics are collected (total of 12 ps per iteration). We calculate the local polarisation by considering conventional Ti-centred unit cells as implemented in references Sepliarsky and Cohen (2011); Gindele et al. (2015); Chapman et al. (2017). Further details are provided in the Supplementary Information.

A dopant-vacancy concentration $n_{d}=100(N_{Ti}^{ideal}-N_{B^{{}^{\prime\prime}}})/N_{Ti}^{ideal}$ is introduced into the supercell initially containing $N_{Ti}^{ideal}$ Ti atoms, by randomly selecting a total of $N_{B^{{}^{\prime\prime}}}$ Ti atoms to be replaced with generic divalent dopants $B_{Ti}^{{}^{\prime\prime}}$. Each dopant is coordinated by six nearest neighbouring oxygen-sites from which a charge compensating oxygen vacancy, $V_{O}^{\bullet\bullet}$, can be introduced. This configuration mimics ($B^{{}^{\prime\prime}}_{Ti}+V^{\bullet\bullet}_{O})^{\times}$ dopant-vacancy associates observed from ESR experiments (Figure 1a). In experiments it is observed that the properties of an aged sample can be removed by heating above the Curie temperature for a long period and then rapidly quenching. It has been hypothesised that during this ‘un-ageing’ process in the cubic phase of the prototype ferroelectric, each orientation of the defect-dipole is equally probable such that vacancies will thermally hop between the neighbouring oxygen site adjacent to the dopant and eventually 1/6 defect-dipoles will populate each of the six possible directions Ren (2004). These are then frozen when quenched into the ferroelectric phase. Ferroelectrics can then be intentionally aged again by applying a bias field for a significantly long period. It is believed this causes defect-dipoles to align. Even in the absence of an ageing field, defect-dipoles in a sample left for a long period will align with the spontaneous polarisation of the domain Zhang et al. (2008a). When constructing the supercell for a particular simulation, the choice of which oxygen is removed neighbouring the dopant depends on the aged/unaged condition:

(1) Unaged condition. To simulate unaged tetragonal PbTiO3 we assign $N_{B^{{}^{\prime\prime}}}/6$ defect-dipoles along each of the six possible orientations causing the total moment to cancel, Figure 1b.

(2) Aged condition. To simulate an aged PbTiO3 sample, each $V_{O}^{\bullet\bullet}$ is selected to situate on the oxygen-site along the ageing direction (defined below) relative to its associated dopant. For these simulations we arbitrarily choose the ageing direction along $+\hat{x}$ (see Figure 1a). This initialises all defect dipoles $\vec{p}_{d}$ as parallel, polarised along [$\bar{1}00$] as shown in Figure 1c.

The ageing direction is defined relative to the driving field for the hysteresis characterisation. If the defect dipoles are co-aligned with the driving field we label this as aged($\parallel$) (Figure 1d), whereas perpendicular alignments are labelled aged($\perp$) (see Figure 1e). The strain is calculated as $\Delta\epsilon=(c_{0}-c)/c$ where $c_{0}$ is the relaxed lattice constant (parallel to the drive field orientation) under no applied field. Further details of the computational methodology are provided in the Supplementary Information.

Firstly we discuss results obtained at 100 K, to observe the ideal behaviour without thermal diffusion or hopping. The calculated P-E hysteresis of PbTiO3 in response to an external driving field is shown in Figure 2a for the defect-free bulk, unaged and the two aged conditions. In all instances, the response is highly non-linear, typical of ferroelectrics. For the ideal bulk case a symmetric, square loop indicative of a hard ferroelectric is observed. We note our bulk coercive field $E_{c}^{int}$ corresponds to the material’s intrinsic coercive field, which greatly exceeds those measured experimentally for Pb-based ferroelectrics Glaum et al. (2012). This is because our model excludes grain boundaries, surfaces and domain walls which would all act as nucleation sites, which lower the energy barrier for reversal in physical samples. Our result of 130 kV/mm matches other MD models Zeng and Cohen (2011) and is in excellent agreement with the intrinsic coercive field of 150 kV/mm calculated using density functional perturbation theory Sai et al. (2002).

Figure 2a shows the hysteresis of an aged single domain simulated sample, with a defect concentration of 1.38%, in response to a driving field perpendicular to the direction in which the material was aged. Interestingly, when the system is equilibrated with no applied field the spontaneous polarisation $\vec{P}$ reorientates parallel to the ageing direction (See Supplementary Figure 2 and cartoon schematic in Figure 1d). This shows the internal bias created from the defect-dipoles is sufficient to overcome the switching barrier Kimmel et al. (2012). This observation provides direct evidence supporting the work of Zhang et al Zhang et al. (2008a) who observed that non-switching defect-dipoles from (Mg${}^{{}^{\prime\prime}}_{Ti}-$V${}_{O}^{\bullet\bullet})^{\times}$ associates in BaTiO3 create restoring forces that promote reversible domain switching. Under the application of the perpendicular driving field there is an almost linear response until 67 kV/mm ($\approx E_{c}^{int}/2$), at which point the field strength is sufficient to switch the polarisation parallel to the drive field. As the electric field decreases to zero, the polarisation again reorientates along the ageing axis such that no remnant polarisation $P_{r}$ remains in the poling direction. Thus in our work, the iconic double-hysteresis indicative of aged ferroelectrics is observed without the requirement of either domain walls or grain boundaries Ren (2004); Genenko et al. (2015).

When poling parallel to the ageing orientation (Figure 1e) the system exhibits a shifted hysteresis curve along the electric-field axis as shown in Figure 2a. Such an effect is well documented in the literature when there is a preferred orientation of the defect dipoles in the poling direction Rojac et al. (2016); Carl and Härdtl (1978); Lambeck and Jonker (1986).

In the unaged simulation we observe a symmetric square E-P hysteresis loop (Figure 2a). The computed coercive field of the unaged PbTiO3 is reduced relative to the ideal bulk value by 35% ($0.65E_{C}^{int}$). The reduction of the coercive field from $E_{c}^{int}$ occurs because the dopant-vacancy associates break local symmetry creating localised areas where the activation energy for nucleation of reverse domains is reduced Vopsaroiu et al. (2010).

The electrostrain (S-E) of each condition is shown in Figure 2b. Symmetric butterfly S-E curves are observed for both the ideal ($n_{d}=0$%) and unaged ($n_{d}=$1.38%) simulations, in excellent agreement with unaged ferroelectrics measured by experiment Glaum et al. (2012); Balke et al. (2009); Zhang et al. (2004). In order to test the validity of results in comparison to experiment, we calculate the d33 piezoelectric tensor coefficient of bulk PbTiO3 using a fluctuation-perturbation theory approach Garcia and Vanderbilt (1998) and second derivative matrices using the GULP package Gale (1997) (see Supplementary Information for further details). We find $d_{33}=47\pm 1$ pC/N from the fluctuations, and 49 pC/N from the derivatives, which agree with each other and agree well with values measured on polycrystalline PbTiO3 films (52-65 pC/N) Kighelman et al. (2002).

Ageing parallel to the poling field is shown to induce an asymmetric S-E hysteresis (Figure 2b). Instances of large asymmetric S-E loops have been experimentally reported in a range of ferroelectric materials Tan et al. (2014). In Shi et al. (2015), the authors report a strain difference of 0.15% in Li doped (Bi0.5Na0.4K0.1)0.98Ce0.2TiO3 ceramics which they propose is due to alignment of (Li${}^{{}^{\prime\prime\prime}}_{Ti}$-V${}_{O}^{\bullet\bullet})^{{}^{\prime}}$ associates. Using our prototypical system, we provide evidence that the asymmetry is likely to arise from an excess orientation along the poling direction and is a general feature of ageing that may be exploited for technological applications.

It has been observed that ageing of doped BaTiO3 is capable of producing a large recoverable non-linear electric field induced strain of 0.75%; far greater than those measured in PZT or PMNPT Ren (2004); Zhang and Ren (2005). It is argued that a strong restoring force from aligned defect-dipoles enables polar axis rotation parallel to the defect-dipoles enabling reversible switching of 90*∘* domains and could lead to the realisation of strain values of 6% in PbTiO3. In Figure 2b, we show that ageing perpendicular to a poling field in PbTiO3 leads to a large recoverable strain in excess of 4.5% (black-diamonds). This large non-linear strain arises from the reorientation of the polar axis along the ageing direction due to the internal bias from the defect dipoles at subswitching fields (see Supplementary Figure 2. $c\to a$, $\Delta\epsilon=(c_{0}-a)/a$). Observing the switching behaviour, we find that, in this instance, the 90*∘* switching occurs via near-homogeneous polarisation rotation over a small field range rather than nucleation and growth of 90*∘* domains - a switching mechanism predicted in bulk PbTiO3 Zeng and Cohen (2011) and BaTiO3 Kimmel et al. (2012). Therefore, we further show the volume effect of the dopant-vacancy associates to be the fundamental cause of this ageing phenomenon and a domain wall mechanism is not required for the full reproduction of experimental observations. Our ageing results are in excellent agreement with a complementary bond valence model study of ageing in ideal BaTiO3 using fixed dipoles introduced into the crystal structure Liu and Cohen (2017).

The results of an investigation into the effect of temperature on the ageing phenomenon in PbTiO3 is shown in Figure 3a for $n_{d}=1.38$%. As the temperature increases, a decrease in the effective coercive field and saturation polarisation are observed, corresponding to a narrowing of the double hysteresis as indicated by the trend arrows. This is analogous to the behaviour known for the square loop of the ideal prototype. Near room temperature under poling fields comparable to $E_{c}$, vacancy hopping becomes thermally activated causing limited events whereby a subset of defect dipoles reorientate. This was observed by tracking the displacement of each oxygen atom relative to its initial position. This reorientation can create asymmetric loops as seen at 400 K ($0.67T_{c}$), clearly demonstrating that at high temperatures/large fields the defect-dipoles can readily realign, elucidating the microscopic mechanism for aged $\to$ unaged transitions. No defect-dipoles were observed to switch below 300 K (0.5Tc). At 300 K a single hopping event was observed (1/24 vacancies) and two (1/12 vacancies) at 400 K, over the full hysteresis. We note that due to the relatively short simulation times these hopping frequencies will be under-sampled for accurate statistics and will be an interesting subject for future investigation.

Defect concentrations close to and above 1.38% are shown to form closed double hysteresis loops described previously (Fig. 3b). As the concentration is increased the enclosed area of the hysteresis loops decrease due to the increased strength of the internal bias, which lowers the barrier for the reorientation of the polar axis. For intermediate defect concentrations (0.78% in this model), we find pinched hysteresis loops are produced. This form of P-E loop is the most common large signal observation noted in experimental studies of aged ferroelectrics Rojac et al. (2016); Glaum et al. (2012); Morozov and Damjanovic (2008). We find that the work dissipated (area enclosed by the P-E loop) decreases with the dopant level (Fig. 3c). Increased defect concentrations start to pinch the square loop which, upon further increases, leads to a closed double hysteresis and gradual reduction of area. Thus, the dissipated energy losses, effective coercive fields and hysteretic behaviour of ferroelectric materials can be controlled by varying the applied fields and dopant levels. We note that in our study we are limited by the constraint of zero total dipole moment in our unaged simulation cell, which restricts the number of dopants $N_{B^{{}^{\prime\prime}}}$ to factors of six. Thus, the concentrations identifying pinching and double hysteresis are, in fact, upper bounds.

In conclusion, we use molecular dynamics to model ageing in boundary-free single domain doped PbTiO3. We show that all the large-signal characteristics of ageing: pinched/double hysteresis, hysteresis shifts and large recoverable non-linear strains, can be reproduced from intrinsic effects of defect-dipoles from dopant-vacancy associates alone, resulting from the net defect dipole orientation with respect to the poling field. Varying the concentration of dopants was found to modify the material’s hysteretic response, suggesting a mechanism for tuning ferroelectric and electromechanical properties for enhanced device performance. This work identifies and clarifies the microscopic mechanisms involved the ageing phenomena and suggests practical methods to inexpensively improve functional performance of ferroelectric ceramic based technologies.

Funding was provided by the EPSRC (EP/G036675/1) via the Centre for Doctoral Training in Molecular Modelling and Materials Science at University College London and the National Measurement Office of the UK Department of Business Innovation and Skills. Computer services on Archer were provided via membership of the UK’s HPC Materials Chemistry Consortium funded by EPSRC (EP/L000202). We acknowledge the use of the UCL facilities LEGION and GRACE, and computational resources at the London Centre for Nanotechnology. REC acknowledges support of the US Office of Naval Research, the ERC Advanced grant ToMCaT, and the Carnegie Institution for Science.