Coherent tunneling in Cu 2 +-and Ag 2 +-doped MgO and CaO : Cu 2 + explored through ab initio calculations

P. García-Fernández,1 C. Sousa,2 J. A. Aramburu,1 M. T. Barriuso,3 and M. Moreno1 1Departamento de Ciencias de la Tierra y Física de la Materia Condensada, Universidad de Cantabria, 39005 Santander, Spain 2Departament de Química Física i Centre Especial de Recerca en Química Teórica, Universitat de Barcelona, 08028 Barcelona, Spain 3Departamento de Física Moderna, Universidad de Cantabria, 39005 Santander, Spain Received 28 July 2005; published 10 October 2005


I. INTRODUCTION
It is well known that the existence of a static or dynamic Jahn-Teller ͑JT͒ effect in many transition metal oxides plays a key role for understanding its properties.So, the colossal magnetoresistance found in some manganites is believed to be a consequence of the interplay between a double exchange mechanism and the JT electron-phonon coupling. 1As another relevant example, the formation of JT polarons has recently been proposed as a possible basis for a superconductivity mechanism in cuprates. 25][6] Apart from reflecting the cubic invariance of the E g e g JT coupling such experimental results stressed the key role played by random strains 7 as responsible for static EPR spectra observed for d 9 or d 7 impurities placed substitutionally in cubic halides. 8espite the large amount of experimental information on JT impurities in sixfold coordination there is a crucial question which has not been solved, why a coherent tunnelling occurs for MgO : Cu 2+ while for Cu 2+ or Ag 2+ in SrO a static JT situation is found. 9,6The present work is aimed to explain this puzzling situation by means of ab initio calculations.A good insight into this relevant issue requires to explore the energy barrier, B, between equivalent minima of the E͑Q ϳ 3z 2 − r 2 ; Q ϳ x 2 − y 2 ͒ surface whose microscopic origin is certainly subtle.Indeed no warping appears in the harmonic approximation for a linear JT coupling. 3,8,10side from the actual value of B there are still doubts about other JT parameters characterizing d 9 impurities in cubic oxides like the JT stabilization energy, E JT .So, Guha and Chase 11 and Reynolds et al. 12 estimated E JT in the range 5500-12 000 cm −1 for Cu 2+ -and Ag 2+ -doped MgO and CaO while from the spin-lattice relaxation time, T 1 , a much lower value E JT Ϸ 300 cm −1 was derived. 13In the case of the excited E g state of MgO : Fe 2+ Hjortsberg et al. 14 estimated E JT Ϸ 1300 cm −1 assuming that the separation between the zero-phonon line and the maximum of the associated broad band is just equal to E JT .This figure is thus similar to that reported by Treguenna-Pigot. 15As regards ab initio calculations E JT values in the 300-500 cm −1 range have recently been reported. 10,16earching to improve our knowledge on JT impurities in oxides, ab initio electronic structure calculations on Cu 2+ and Ag 2+ impurities in cubic oxides as a function of both Q = cos and Q = sin coordinates have been carried out.Particular attention has been addressed to calculate the values of the JT stabilization energy, E JT , and the barrier energy, B, as well as to explore their microscopic origin.Very recently 17 ab initio calculations on NaCl: Rh 2+ and triangular molecules like Cu 3 or Ag 3 have shown that B mainly arises from the anharmonicity of the JT mode.The results obtained in this work have been derived using two different methods in order to reinforce the reliability of present conclusions.

II. THEORY
Following previous works on NaCl: Rh 2+ , 10,17  2− ͑M = Cu, Ag; B =Mg,Ca,Sr͒ clusters centered around the impurity were carried out using the density functional theory ͑DFT͒ methodology as implemented in the Amsterdam density functional ͑ADF͒ 2000.02code. 18The use of finite clusters for describing the properties of these systems is consistent with the highly localized character of the unpaired electron residing solely in the MO 6 complex formed by the impurity and the six nearest O 2− anions. 4,6,12In the present calculations a generalized gradient approximation ͑GGA͒ was used for the exchange-correlation functional in the Becke-Perdew form 19 and the basis set of the best quality in the ADF code was employed.A point charge electrostatic potential was included to take into account the influence of the rest of the ions in the host lattice.
As it is well known, DFT calculations present limitations in the treatment of open shell electronic states, where spin cannot be neglected, so requiring the use of highly correlated calculation methods. 20,17Taking into account the very subtle effects studied in this work, all DFT results have been complemented by very accurate wave-function-based complete active space perturbation theory of order 2 ͑CASPT2͒ calculations as implemented in the MOLCAS 5.3 program. 21ASPT2 calculations include most of the electronic correlation effects by combining two strategies.In a first step, a correct treatment of the large electron correlation effects in the transition metal d shell is treated in a variational way by constructing an accurate multireference CASSCF wave function for the considered states, with an active space containing nine electrons distributed in all posible ways over 10 orbitals, the five d orbitals of the impurity and a correlating set of five virtual orbitals. 21In the second step, the remaining ͑mostly dynamical͒ electron correlation effects are included by many-body second-order perturbational treatment, CASPT2, taking the CASSCF wave function as the zerothorder wave function, in which all valence electrons are correlated.Because of the very high computational cost of these calculations, most of the calculations were performed on MO 6 B 6 4− clusters with only 13 ions.An accurate embedding was constructed surrounding the cluster by ab initio model potentials ͑AIMPs͒ 22 that take into account the short-range repulsion of the rest of the lattice, and a set of point charges.Moreover, calculations carried out for 25-atom MO 6 B 12 B 6 26+ clusters did not improve significantly the results derived on the 13-atom clusters.
In the study of JT impurities in cubic crystals it is quite useful to consider first the isotropic distortion produced by a hypothetical ion whose electronic density resembles that of closed shell structure. 10,17,23This leads to an impurity-ligand distance, R oct , different from R 0 corresponding to the perfect host lattice.In DFT an electronic density displaying cubic symmetry can easily be built for a d 9 impurity through the ͑3z 2 − r 2 ͒ 1.5 ͑x 2 − y 2 ͒ 1.5 configuration. 10,17,23In a second step, adiabatic calculations are carried out for both as a function of the two Q and Q coordinates.For instance, under a Q = ͱ 12u deformation, axial and equatorial metalligand distances are given by R ax = R oct +2u and R eq = R oct − u, respectively.
By means of these calculations the parameters involved in the usual effective Hamiltonian, H eff , describing the E g e g JT effect within the ͕x 2 − y 2 ,3z 2 − r 2 ͖ basis can be derived and related to the equilibrium geometry at adiabatic minima.H eff is usually written as

͑1͒
Here U and U are 2 ϫ 2 Pauli matrices 7,8 and I means the identity matrix.The first term in ͑1͒ denotes the linear JT coupling while the second and third ones are the elastic and anharmonic contributions, respectively. 3,7,8V 2a is related to the angular frequency of the e g mode, e , by M L e 2 =2V 2a , where M L is the ligand mass.The fourth term depicts the so-called quadratic JT coupling.Recent results 10,17,24 have demonstrated that this contribution plays a minor role for explaining the barrier compared to the anharmonic term.For instance, the relative contribution of the anharmonicity to the barrier is found to be equal to ϳ98% and ϳ70% for CaO:Cu 2+ and SrO : Cu 2+ , respectively. 24For this reason the quadratic JT coupling will not be considered in the subsequent analysis mainly devoted to explore the sign and absolute value of the anharmonic constant, V 3a .
The expansion described in ͑1͒ is reasonable provided 0 Ӷ ͱ 12R oct at adiabatic minima, a condition which is usually followed by JT impurities. 10If only the two first terms in ͑1͒ are considered it is found: 3,7,8 implying that warping in E͑ 0 , ͒ depends on the V 3a term in ͑1͒.When this is taken into account the adiabatic energy of the two solutions of H eff can be approximated by 17 where A Ϸ V 3a 0 3 .Equation ͑3͒ means that adiabatic minima in the ground state display a tetragonal symmetry.If A Ͼ 0 ͑A Ͻ 0͒ then the equilibrium geometry at the three adiabatic minima corresponds to a compressed ͑elongated͒ octahedron.The energy barrier between a compressed ͓3 = ͑2p +1͒ , p = integer͔ and an elongated ͑3 =2p͒ geometry in the ground state is just equal to 2A.Around a minimum, E − ͑ 0 ͒ can be approximated by a harmonic oscillator potential whose associated angular frequency is given by c 2 =9͉A͉ / ͑M L 0 2 ͒.A coherent tunnelling among equivalent adiabatic minima is thus favored if c /2Ͼ ͉B͉. 10 A main goal of the present calculations is just to compute the F = ͉B͉ / ͑ c /2͒ parameter from the previously calculated V 3a and 0 values.

III. RESULTS AND DISCUSSION
Consistently with the size of the employed clusters in all present calculations more than 98% of the unpaired electron density is found to reside within the complex formed by the impurity and ligands.Tables I and II gather the parameters for Cu 2+ -and Ag 2+ -doped cubic oxides obtained through CASPT2 calculations, which according to previous results 20,17 is the most reliable method for the present systems.For comparison, computed values are also given for NaCl: Rh 2+ displaying a static JT effect. 25For the same reason, Table I also includes the DFT values of R oct , 0 , E JT , and e , although the poor self-consistent field convergence in the computed E ± ͑ , ͒ curves are prevented in a few cases to obtain reliable values for the parameters collected in Tables I and II.
CASPT2 and DFT values of R oct and 0 given in Table I are reasonably similar although the latter method leads in general to smaller values for E JT .A similar pattern was previously obtained through Hartree-Fock calculations. 26ence, our results support a E JT value in the range 300-500 cm −1 for MgO : Cu 2+ .A similar value ͑E JT = 394 cm −1 ͒ has been obtained for the excited 5 E g state of MgO : Fe 2+ .The much larger estimate of E JT by Hjortsberg et al. 14 can be ascribed to the neglect of the contribution of the symmetric a 1g mode to the ⌬E Ϸ 1100 cm −1 separation between the electric dipole origin and the maximum of the associated broad band.For transition-metal impurities at octahedral sites ⌬E Ϸ S a a + S e e , 27 where E JT ϵ S e e and S e ͑S a ͒ is the Huang-Rhys factor of the e g ͑a 1g ͒ mode.Present DFT calculations give S a a = 360 cm −1 thus supporting a E JT value again not far from 500 cm −1 rather than close to 1500 cm −1 . 14Moreover, the value V 2a = 1.85 eV/ Å 2 derived for CaO : Cu 2+ leads to e = 250 cm −1 which is not far from the Raman peak estimated at ϳ300 cm −1 . 11t us now focus on calculated parameters in Table II whose main features are the following: ͑i͒ The values of ͉B͉ obtained for Cu 2+ -and Ag 2+ -doped MgO and CaO : Cu 2+ are smaller than 20 cm −1 and thus at least 25 times smaller than the figure derived for NaCl: Rh 2+ .͑ii͒ The F parameter is found to be smaller than 0.4 for Cu 2+ -and Ag 2+ -doped MgO and CaO : Cu 2+ , while it is equal to 6.4 for NaCl: Rh 2+ .͑iii͒ On passing from a MgO to a SrO host lattice F increases by a factor equal to 20 and 7 for Cu 2+ and Ag 2+ impurities, respectively.These trends strongly support that MgO is a more suitable host lattice than SrO for observing dynamic EPR spectra of d 9 ions in agreement with experimental facts. 4,6,9͑iv͒ Similarly to the results obtained for NaCl: Rh 2+ , 10 the barrier B is found to be dominated by the anharmonic contribution for Cu 2+ and Ag 2+ in cubic oxides.However, the sign of B and V 3a in the latter cases is in general opposite to that encountered in NaCl: Rh 2+ where the elongated RhCl 6 4− complex has been demonstrated to be significantly decoupled from the rest of the host lattice. 10͑v͒ The evolution of E JT / e in Table II exhibits a similar trend as the F parameter.However, direct information about the barrier is only embodied in the last parameter.
A microscopic explanation of the main trends of B and F requires a detailed analysis of the equilibrium geometry displayed in Table I.First, it can be noticed that the calculated R oct distances follow rather closely the value of R 0 , the differences between R oct and R 0 being less than 5%.Therefore, according to the theory of elasticity, the second neighbors A 2+ ͑A = Mg, Ca or Sr͒ in ͗100͘ directions ͑Fig.1͒ are expected to be shifted by less than 1%.This argument justifies a posteriori our previous assumption to keep A 2+ ions at the host lattice positions.It is worth stressing that results in Table I are qualitatively different from those for impurities in halide lattices.In these systems R oct is essentially determined by the relative size of impurity and ligand and only slightly dependent on the host lattice.This idea can be quantified by means of the parameter f = ⌬R oct / ⌬R 0 , the ratio of the variations of R 0 and R oct induced by a change of host lattice.For Cr 3+ in elpasolites 28 or Fe 3+ in fluoroperovskites 29 f has been determined to be ഛ0.1 while a value f = 0.3 reflects the behavior of Mn 2+ impurities in cubic fluoroperovskites. 27For Cu 2+ and Ag 2+ in cubic oxides f values derived from Table I Values of coupling constants V 1e , V 2a , and V 3a , energy barrier B, and F = ͉B͉ / ͑ c /2͒ and E JT / ប e ratios obtained by means of CASPT2 calculations.Units of V i␥ are eV/ Å i while B is in cm −1 .If B Ͼ 0 the minimum corresponds to a compressed configuration and the top of the barrier to an elongated one.The opposite happens for B Ͻ 0. lie in the 0.6-0.8range.This strongly suggests that the interaction of the ligand ion with the nearest A 2+ cation of the host lattice ͑Fig.1͒ dominates over the impurity-ligand interaction.In a ball-and-spring model ͑Fig.2͒ the energy of a given M 2+ -O 2− -A 2+ bond in Fig. 1 is given by 27 where x M , x L , and x A mean the position of the impurity, the ligand and the A 2+ ion lying in the same ͗100͘ direction, respectively, while r 0 and r 0 Ј are parameters describing the two springs.If the impurity and the A 2+ ion are fixed the energy changes induced by variations of the ligand around the equilibrium position, x L0 , can be approximately described by

͑5͒
Calculated values of the k / kЈ ratio for the present systems are collected in Table III.The value k / kЈ = 2.5 found for NaCl: Rh 2+ indicates a significant decoupling of the RhCl 6 4− complex from the ions of the host lattice.By contrast the vibrational coupling of CuO 6 10− units ͑where active electrons reside͒ to the host lattice is, as expected, much stronger, the k / kЈ values lying around 0.4.In fact, for cases like MgO:Cu 2+ ͑but not for NaCl: Rh 2+ ͒ the impurity and the host cation have the same nominal charge.
B and F values are greatly sensitive to 0 which in turn depends on the ratio V 1e / V 2a ͓Eq.͑2͔͒.Although, as expected, both V 1e and V 2a decrease on passing from MgO to SrO, the ratio V 1e / V 2a and thus 0 experience a significant increase.This fact is quite consistent with the well-known increase of the Huang-Rhys factors S a and S e when R 0 increases, 27,30 explored first in Ref. 31.For instance S e 2 is proportional to V 1e / e 3/2 . 27A similar behavior has also been observed 32,27 for the Stokes shift ͑depending on 2S a a +2S e e ͒ and is thus consistent with the increase of E JT on going from MgO to SrO.
Although the important increase of 0 along the MgO → CaO → SrO host lattices favors a transition from a dynamic to a static regime it is also crucial to consider the influence of the anharmonic constant V 3a on both B and F. As shown in Table II, on passing from NaCl: Rh 2+ to MgO:Cu 2+ V 2a increases by a factor of 3 while ͉V 3a ͉ decreases by a factor of 2. Moreover V 3a Ͻ 0 for the former system while V 3a Ͼ 0 for the latter.According to Eq. ͑5͒ V 3a is proportional to the difference between the two anharmonic constants g and gЈ that are negative for ordinary bond types.The negative sign of both V 3a and B parameters found for NaCl: Rh 2+ is consistent with the relevant decoupling of RhCl 6 4− from the rest of the lattice. 10In addition, the anharmonicity of the Rhu Cl bond contributes to the higher stability of the elongated over the compressed form of the RhCl 6 octahedron.On the contrary, if the relative importance of kЈ increases, V 3a is expected to be softened and can even experience a change of sign leading to positive B values.The small k / kЈ ratio in Table III for SrO : Cu 2+ is accompanied by a rather large, positive V 3a and a sizeable B = + 98 cm −1 .Hence, the interpretation of our theoretical result based on the ball-and-spring model supports the existence of a compressed equilibrium geometry as demonstrated in the static EPR spectrum of SrO : Cu 2+ . 9Similar observations hold for CaO:Cu 2+ , for which the same equilibrium geometry has been inferred from the Raman spectra corresponding to vibronic levels. 11assing from first-row to second-row transition-metal ions in general more covalent bonds are formed.So, it can be expected that the replacement of the 3d Cu 2+ ion for the 4d ion Ag 2+ increases the relative importance of the impurityligand constant k.We obtain a k / kЈ ratio of 0.75 for SrO:Ag 2+ , which is the largest value found in all six systems studied here ͑Table III͒.In agreement with this, we find a very small V 3a value, close to zero, i.e., at the border of the inversion of the relative stability of elongated and compressed forms of the AgO 6 10− complex.Nevertheless, B is still found to be positive ͑71 cm −1 ͒ ͑Table II͒, whereas experimental EPR data on this system indicate however an elongated octahedral geometry. 6

IV. FINAL REMARKS
The present results, together with those recently obtained on layered perovskites, 23 stress the importance of ab initio calculations for gaining a better insight into the subtle phenomena related to the JT effect.This is particularly true as regards the barrier energy, B, on which the experimental information is certainly scarce.This parameter plays a key role for the existence of coherent tunnelling, whose relevance for quantum information processing has recently been pointed out. 33It is worth noting that static situations can experimentally be found even if F Ӷ 1.It was first pointed out by Ham 7 that defects can destroy the coherence and lock at T = 0 K the wave function in the lowest potential well.In this sense  remote defects leading to energy differences among the three wells of only ϳ1-10 cm −1 can suffice for producing a static situation.
As a main conclusion, it has been shown in this work that both the sign of V 3a and the nature of the ground state can be modified simply by considering the elastic coupling of the complex to the rest of the lattice.Another relevant result concerns the variation of the average impurity-ligand distance along the series of cubic oxides.While for transitionmetal impurities in fluorides or chlorides with sixfold coordination the impurity-ligand distance is usually close to the sum of ionic radii, 27,29 for present cases it follows the host lattic parameter.This different behavior can be related to the bigger deformability of the O 2− ion in comparison to more rigid ions like F − or Cl − . 34As the electronic properties associated with the impurity essentially depend on the actual value of the impurity-ligand distance it can be expected 35 that the charge transfer transitions undergo a strong redshift on passing from MgO : Cu 2+ to SrO : Cu 2+ .Further work along this line is underway.

FIG. 2 .
FIG. 2. ͑Color online͒ Model representing three atoms ͑M = impurity metal ion, L = oxygen ligand, and A = next neighbor in the ͗100͘ direction͒ connected by two springs with parameters ͑k , g͒ and ͑kЈ , gЈ͒, respectively.x i represents the position of atom i.

TABLE I .
Calculated values of R oct and 0 ͑in Å units͒ and E JT and ប e ͑in cm −1 units͒.CASPT2 results are shown in normal characters while DFT ones are shown in italic characters.The experimental metal-ligand distances R 0 ͑in Å units͒ for perfect host lattices are given for comparison.

TABLE III .
Values of harmonic force constants k and kЈ ͑in eV/ Å 2 ͒ from DFT calculations.