# Oscillatory Tunnel Splittings in Spin Systems: A Discrete Wentzel-Kramers-Brillouin Approach

## Abstract

Certain spin Hamiltonians that give rise to tunnel splittings that are viewed in terms of interfering instanton trajectories, are restudied using a discrete WKB method, that is more elementary, and also yields wavefunctions and preexponential factors for the splittings. A novel turning point inside the classically forbidden region is analysed, and a general formula is obtained for the splittings. The result is appled to the Fe<sub>8</sub> system. A previous result for the oscillation of the ground state splitting with external magnetic field is extended to higher levels.

The magnetic properties of the molecular cluster \[(tacn)<sub>6</sub>Fe<sub>8</sub>O<sub>2</sub>(OH)<sub>12</sub>\]<sup>8+</sup> (or just Fe<sub>8</sub> for short) are governed by a Hamiltonian

$$=k_2J_z^2+(k_1k_2)J_x^2g\mu _B𝐉𝐇,$$
(1)   
where $`𝐉`$ is a dimensionless spin operator, $`𝐇`$ is an externally applied magnetic field, $`J=10`$, $`k_10.33`$ K, and $`k_20.22`$ K. The zero-field Hamiltonian has biaxial symmetry with easy, medium, and hard axes along $`z`$, $`y`$, and $`x`$ respectively . Very recently, Wolfgang and Sessoli have seen a new effect in this system, viz., an oscillation in the Landau-Zener transition rate between Zeeman levels, as a function of the applied field along $`\widehat{𝐱}`$. These oscillations reflect an oscillation in the underlying tunneling matrix element between the levels in question, and are in this author’s view, the only unambiguous evidence to date for quantum tunneling of a spin of such a large size in a solid state system. This effect is not seen, e.g., in a closely related Mn<sub>12</sub> cluster.

Oscillations as a function of $`H_x`$ in the ground state tunnel splitting $`\mathrm{\Delta }`$ of the Hamiltonian (1) were in fact predicted earlier, on the basis of an instanton calculation. For $`𝐇\widehat{𝐱}`$, there are two instantons, with a complex action differing by a Berry phase that sweeps through odd multiples of $`\pi `$ as $`H_x`$ is varied, leading to a complete quenching of tunneling. In this view the effect arises from destructive intereference between spin trajectories . A different perspective was provided in Ref. , by noting that $``$ is invariant under a $`180^{}`$ rotation about $`\widehat{𝐱}`$ when $`H_z=0`$. The oscillation is due to a symmetry-allowed crossing of levels with different parity under the rotation.

Although the oscillations with $`H_x`$ are easily seen by direct numerical diagonalization of $``$ (see Fig. 1 of , or Fig. 4 of ), it is of interest to understand these features analytically. Since the spin $`J`$ is large, it is natural to use the semiclassical, or $`J\mathrm{}`$ approximation. To some extent, the instanton method already does this. In this paper, we will try and make further progress using a discrete WKB method , which has several advantages. First, it is very difficult to find the next-to-leading terms in the $`J\mathrm{}`$ asymptotic expressions for various physical quantities (such as $`\mathrm{\Delta }`$) using instantons. Second, Wernsdorfer and Sessoli also see an oscillatory rate in the presence of a dc field along $`\widehat{𝐳}`$ that is such as to align the ground level in one well with an excited level in the other. The splitting is now never perfectly quenched as the symmetry of $``$ is destroyed. We do not (although others may) know how to solve this problem with instantons. Third, with an eye to the future, the method provides wavefunctions in addition to energies, which may be used to calculate matrix elements of various perturbations  to Eq. (1) that are present in the actual physical system, and thus study their influence.

Specifically, we will derive a general result \[see Eq. (19)\] for the tunnel splitting between degenerate pairs of levels in a symmetric problem when oscillations are present. Our result is expressed in terms of two action integrals. In the course of doing this, we will encounter a novel feature that does not arise in previous discrete WKB studies, namely, a turning point in the classically forbidden region! We will apply our result to the hamiltonian (1) for $`H_z=0`$, focussing in detail on quenching fields for ground and excited state splittings. The results for the latter are new. The study of imperfect quenching of tunneling that occurs when $`H_z0`$, because the potential is then asymmetric, is much more involved, and will be published separately.

Let us first briefly review the discrete WKB formalism . The starting point is to write Schrödinger’s equation in the $`J_z`$ basis. Let $`|\psi =E|\psi `$, $`J_z|m=m|m`$, $`m|\psi =C_m`$, $`m||m=w_m`$, and $`m||m^{}=t_{m,m^{}}`$ ($`mm^{}`$). Then we have

$$\underset{nm}{}t_{m,n}C_n+w_mC_m=EC_m.$$
(2)   
We assume that the matrix $`t_{m,n}`$ is real and symmetric, $`t_{m,n}=t_{n,m}`$. In the present problem, we need matrix elements that are off-diagonal by 1 ($`t_{m,m\pm 1}`$) and by 2 ($`t_{m,m\pm 2}`$). This makes Eq. (2) a recursion relation involving five terms, as opposed to three terms in previous work. The physical idea is to view Eq. (2) as a tight-binding model for an electron hopping on a one-dimensional lattice, and use the approximation of semiclassical electron dynamics. This would be exact if the matrix elements of $``$ were constant with $`m`$, and will be systematically justifiable if they are slowly varying with $`m`$. Formally, the latter means that we can find functions $`w(m)`$, $`t_1(m)`$, and $`t_2(m)`$ of a continuous variable $`m`$, such that on the discrete eigenset of $`J_z`$,

$`w(m)`$ $`=`$ $`w_m,`$ (3)   
$`t_\alpha (m)`$ $`=`$ $`(t_{m,m+\alpha }+t_{m,m\alpha })/2,\alpha =1,2,`$ (4)   
and further, that if $`m/J`$ is regarded as a quantity of order $`J^0`$, then $`\dot{w}(m)dw/dm=O(w(m)/J)`$, with similar restrictions on $`\dot{t}_1(m)`$ and $`\dot{t}_2(m)`$. For Eq. (1), these conditions are met if $`J1`$.

If $`w_m`$, $`t_{m,m\pm 1}`$, and $`t_{m,m\pm 2}`$ were constant, the eigenstates of $``$ would be states with $`C_m=e^{iqm}`$, and $`E=w+2t_1\mathrm{cos}q+2t_2\mathrm{cos}(2q)`$. Now we seek a solution in the form $`C_m=e^{i\mathrm{\Phi }(m)}`$ with $`\mathrm{\Phi }=\mathrm{\Phi }_0+\mathrm{\Phi }_1+\mathrm{\Phi }_2+\mathrm{}`$, where $`\mathrm{\Phi }_n=O(J^{1n})`$, and $`\dot{\mathrm{\Phi }}_n=O(\mathrm{\Phi }_n/J)`$. Then, one can show  that up to terms of order $`J^0`$ in $`\mathrm{\Phi }`$, the solution is given as linear combinations of the form

$$C_m\frac{1}{\sqrt{v(m)}}\mathrm{exp}\left(i^mq(m^{})𝑑m^{}\right),$$
(5)   
where $`q(m)`$ is a local wavevector that obeys the eikonal or Hamilton-Jacobi equation,

$$E=w(m)+2t_1(m)\mathrm{cos}q+2t_2(m)\mathrm{cos}(2q)_{\mathrm{sc}}(q,m),$$
(6)   
and $`v(m)`$ is the associated semiclassical electron velocity, which obeys the transport equation

$$v(m)=_{\mathrm{sc}}/q=2\mathrm{sin}q(m)(t_1(m)+4t_2(m)\mathrm{cos}q(m)).$$
(7)   
To talk of tunneling, we must first understand the classically allowed and forbidden regions in the $`m`$ space. As a function of $`q`$ for fixed $`m`$, the semiclassical Hamiltonian $`_{\mathrm{sc}}(q,m)`$ can be viewed as a band energy curve, and its minimum and maximum values define local band-edge functions $`U_\pm (m)`$. The classically accessible region for any energy $`E`$ is thus defined by $`U_{}(m)EU_+(m)`$. \[The first consequence of having five terms in the recursion relation shows up here. In the three term case, the band edges always occur at $`q=0`$ or $`\pi `$. Now, they can occur at values other than these if $`|t_1(m)/4t_2(m)|<1`$. These functions are sketched in Fig. 1 for Eq. (1) with $`H_z=0`$. The minimum, $`U_{}`$, is attained at $`q=0`$ for $`|m|m^{}`$, and at $`q0`$ for $`|m|<m^{}`$. The curve $`U_{}(m)`$ is smooth at $`m=\pm m^{}`$, and the formula for $`m^{}`$ is unimportant.\] Thus, for the energies $`E_a`$ and $`E_b`$ drawn in Fig. 1, the central region is classically forbidden and allowed, respectively. We will focus on states of the first type in what follows.

The next step is to derive a generalization of Herring’s formula for the tunnel splitting $`\mathrm{\Delta }`$ for a pair of levels whose mean energy is $`E`$. Proceeding in exact analogy with Ref.  (see also (c,d) or ), we consider a solution $`C_m`$ to Eq. (2) with energy $`E`$, that is (a) localized in the left well of $`U_{}(m)`$, and decays away from that well everywhere including the region near $`m=0`$, and (b) normalized to unit total probability. The behavior of this solution near the right well need not be specified or examined too closely. Up to an irrelevant over all sign, we find

$$\mathrm{\Delta }=\{\begin{array}{cc}2\left[t_{01}C_0(C_1C_1)+t_{02}C_0(C_2C_2)+t_{1,1}(C_1^2C_1^2)\right],\hfill & \\ integer\text{J},\hfill & \\ 2t_{\frac{1}{2},\frac{1}{2}}\left(C_{\frac{1}{2}}^2C_{\frac{1}{2}}^2\right)+4t_{\frac{3}{2},\frac{1}{2}}\left(C_{\frac{1}{2}}C_{\frac{3}{2}}C_{\frac{1}{2}}C_{\frac{3}{2}}\right),\hfill & \text{half-integer }J\text{.}\hfill \end{array}$$
(8)   
To apply Eq. (8), we must find $`C_m`$ in the central region. In principle the procedure is straightforward, and follows conventional WKB. We first find $`C_m`$ in the allowed region, near $`m_0`$, and then use connection formulas to extend it into the forbidden region. For a three term recursion relation, this is done in Ref. . In the present case, we encounter a new difficulty. To see this, we consider points at which $`v(m)`$ vanishes. At all such points, which may be called turning points, the solution (5) diverges, indicating a breakdown of the WKB approximation. Let us now consider a point strictly inside the classically allowed region in the $`E`$-$`m`$ plane. At such a point $`q`$ is not an extremum of $`_{\mathrm{sc}}`$ for fixed $`m`$, i.e., $`v(m)0`$. It is a simple corollary that the points $`E=U_\pm (m)`$ are turning points, corresponding to $`q=0`$, $`\pi `$, or $`\mathrm{cos}^1(t_1/4t_2)`$. These turning points are of the same physical character as those in conventional WKB, and the $`q=0`$ or $`\pi `$ ones are the only ones that arise with a three-term recursion relation.

For our five-term recursion, however, $`v(m)`$ can also vanish if $`\mathrm{cos}q=t_1/4t_2`$, even though $`EU_\pm (m)`$. To see how this can happen, we solve Eq. (6) to get

$$\mathrm{cos}q(m)=\frac{t_1(m)\pm [t_1^2(m)4t_2(m)f(m)]^{1/2}}{4t_2(m)},$$
(9)   
where $`f(m)=w(m)2t_2(m)E`$. Thus, such a turning point may arise when the discrimiant of the quadratic equation for $`\mathrm{cos}q(m)`$ vanishes. Since, by exclusion, such points must necessarily lie in a classically forbidden region, where $`q(m)`$ is not real, it follows that they can only arise in problems where $`|t_1(m)/4t_2(m)|>1`$ for some $`m`$. This fact and Eq. (9) then imply that at such a point, $`\mathrm{cos}q`$ changes from real to complex, i.e., $`q`$ changes from pure imaginary to complex, and the wavefunction accordingly changes from an exponential decay with one sign to a decay with an oscillating sign.

Since WKB breaks down at the forbidden region turning points, we need connection formulas at these points just as for ordinary ones . We will publish the derivation of these formulas elsewhere, and here we only give the result. Let the discriminant in Eq. (9) vanish at $`m=m_c`$, and let $`\mathrm{cos}q`$ be real for $`m<m_c`$, and complex for $`m>m_c`$. It is convenient to define $`q(m)=i\kappa (m)`$ with $`\kappa >0`$ in the region $`m<m_c`$, and to write $`s(m)=iv(m)`$ everywhere. (This definition renders $`s(m)>0`$ for $`m>m_c`$.) We consider the decaying WKB solution in the region $`m<m_c`$:

$$C_m=\frac{A}{2\sqrt{s(m)}}\mathrm{exp}\left(\underset{m_c}{\overset{m}{}}\kappa (m^{})𝑑m^{}\right),m<m_c,$$
(10)   
where $`A`$ is chosen to be real. For $`m>m_c`$, we must consider linear combinations of the type (5), with two choices for $`q(m)`$ which we write as

$$q_{1,2}(m)=i\kappa (m)\pm \chi (m).$$
(11)   
For the solution to continue decaying, we must still have $`\kappa >0`$, and we also choose $`\chi >0`$. Then, both $`\kappa (m)`$ and $`\chi (m)`$ have a kink at $`m=m_c`$. We further define $`s_{1,2}(m)=iv(q_{1,2}(m))`$ via Eq. (7), so that $`s_2=s_1^{}`$. The WKB solution which connects to (10) is then given by

$$C_m=\mathrm{Re}\frac{A}{\sqrt{s_1(m)}}\mathrm{exp}\left(i\underset{m_c}{\overset{m}{}}q_1(m^{})𝑑m^{}\right),m>m_c.$$
(12)   
Note that this is explicitly real, as the reality of Eq. (2) requires. Also, Eqs. (10) and (12) only hold for $`|m_cm|J^{1/3}`$. The connection formula for the growing solution is similar, but is not needed for our present purpose.

The result (8) is exact, but does not reveal the physically important barrier penetration factor. To remedy this, we substitute Eq. (12) in Eq. (8). We consider a situation as in Fig. 1, with minima in $`U_{}(m)`$ at $`\pm m_0`$, and forbidden region turning points at $`\pm m_1`$. The key, clearly, is to simplify Eq. (12) in the region $`|m|<m_1`$. To this end, we substitute Eq. (11) for $`q`$ in Eq. (9) and separate the real and imaginary parts. This yields

$`\mathrm{cosh}\kappa \mathrm{cos}\chi `$ $`=`$ $`t_1/4t_2,`$ (13)   
$`\mathrm{sinh}\kappa \mathrm{sin}\chi `$ $`=`$ $`(4t_2ft_1^2)^{1/2}/4t_2.`$ (14)   
Using these results, it follows that

$$s_1=8t_2(m)\mathrm{sinh}\kappa (m)\mathrm{sin}\chi (m)\mathrm{sin}q_1(m).$$
(15)   
We now specialize to the case of integer $`J`$; the other is similarly analysed, and yields the same result, Eq. (18) below. For Eq. (8), we only need $`C_m`$ for $`m=0`$, $`\pm 1`$, and $`\pm 2`$. The variations in $`\kappa `$, $`\chi `$, $`t_1`$, $`t_2`$, and $`q_1`$ between these points may be ignored as they are of of order $`J^2`$. Hence, to sufficient accuracy one may write (for $`|m|<2`$),

$$C_m=\mathrm{Re}A_2\frac{e^{i(\mathrm{\Omega }+mq_{10})}}{\sqrt{\mathrm{sin}q_{10}}},$$
(16)   
where $`\mathrm{\Omega }=_{m_1}^0q_1(m^{})𝑑m^{}`$, and $`A_2=(8t_{20}\mathrm{sinh}\kappa _0\mathrm{sin}\chi _0)^{1/2}A`$. The suffix 0 denotes quantities evaluated at $`m=0`$; thus $`q_{10}=q_1(0)`$, $`\kappa _0=\kappa (0)`$, etc. To the same accuracy as Eq. (16) one may write $`t_{01}=t_1(0)`$, and $`t_{02}=t_{1,1}=t_2(0)`$ in Eq. (8). If we use Eq. (13), and write every thing in terms of $`t_{20}`$, $`\kappa _0`$, and $`\chi _0`$, then a certain amount of algebra leads to

$`\mathrm{\Delta }`$ $`=`$ $`{\displaystyle \frac{1}{2}}A^2(e^{2i\mathrm{\Omega }}+e^{2i\mathrm{\Omega }^{}})`$ (17)   
$`=`$ $`A^2\mathrm{exp}\left({\displaystyle \underset{m_1}{\overset{m_1}{}}}\kappa (m^{})𝑑m^{}\right)\mathrm{cos}\left({\displaystyle \underset{m_1}{\overset{m_1}{}}}\chi (m^{})𝑑m^{}\right).`$ (18)   
The cosine factor clearly shows the possibility of oscillations.

The next step is to match the WKB wavefunction (10) in the ordinary decaying region to the wavefunction in the allowed region. It is plain that $`A`$ will contain an additional barrier penetration factor $`\mathrm{exp}(_{m_1}^{m_t}\kappa (m)𝑑m)`$, where $`\pm m_t`$ are the ordinary turning points. Omitting the details of the calculation, which are very much like those in conventional WKB, we find that for the $`n`$th pair of levels, provided $`nJ`$,

$$\mathrm{\Delta }=\frac{2\omega _0}{\pi }g_n\mathrm{exp}\left(\underset{m_t}{\overset{m_t}{}}\kappa (m^{})𝑑m^{}\right)\mathrm{cos}\left(\underset{m_1}{\overset{m_1}{}}\chi (m^{})𝑑m^{}\right),$$
(19)   
where $`\omega _0`$ is the small oscillation frequency in the wells near $`\pm m_0`$, and $`g_n=\sqrt{2\pi }\overline{n}^{\overline{n}}e^{\overline{n}}`$ (with $`\overline{n}=n+\frac{1}{2}`$). It need hardly be said that $`m_t`$ and $`m_1`$ depend on the energy and hence on $`n`$.

Equation (19) is a general formula for the splitting in the presence of interference effects. As opposed to an “exponentially accurate” calculation which gives an asymptotically correct result for $`\mathrm{ln}\mathrm{\Delta }`$ as $`J\mathrm{}`$, it is correct for $`\mathrm{\Delta }`$ itself.

We now apply Eq. (19) to the Hamiltonian Eq. (1). The problem is now merely one of quadrature, so we will focus only on the cosine factor. (The full expression for the non-oscillatory part of $`\mathrm{\Delta }`$ including the exact prefactor is exceedingly lengthy and unilluminating. A partial result for the WKB exponent, or Gamow factor may be found in Ref. .) In doing the quadratures, the first step is to find $`w(m)`$, $`t_\alpha (m)`$, etc. Here, any function that reproduces the first two terms in a series in $`1/J`$ is adequate, since Eq. (5) represents only the two leading terms in $`\mathrm{\Phi }(m)`$. We define $`\overline{J}=J+\frac{1}{2}`$, $`\mu =m/\overline{J}`$, $`H_c=2k_1J/g\mu _B`$, $`h_x=JH_x/\overline{J}H_c`$, $`\lambda =k_2/k_1`$, and measure all energies in units of $`k_1\overline{J}^2`$. Then, $`t_1=h_x(1\mu ^2)^{1/2}`$, $`4t_2=(1\lambda )(1\mu ^2)`$, $`2w=(1+\lambda )(1\mu ^2)`$, and $`\omega _0=2[\lambda (1h_x^2)]^{1/2}/\overline{J}`$. Further, let us denote the argument of the cosine in Eq. (19) by $`\mathrm{\Lambda }`$. The turning point $`\mu _1(E)`$ is given by

$$\mu _1^2(E)=1\frac{h_x^2}{1\lambda }\frac{E}{\lambda }.$$
(20)   
Secondly, from Eq. (14), we see that $`\chi (\mu _1(E)\mu )^{1/2}`$, so that to relative order $`1/J`$, we may write

$$\mathrm{\Lambda }=2\overline{J}\underset{0}{\overset{\mu _{10}}{}}\left(\chi (\mu ,E=0)+E\frac{\chi }{E}|_{E=0}\right)𝑑\mu ,$$
(21)   
where $`\mu _{10}=\mu _1(0)`$, and $`E=(n+\frac{1}{2})\omega _0`$ for the $`n`$th pair of levels. At $`E=0`$, we have $`\mathrm{cos}^2\chi =[(1\mu _{10}^2)/(1\mu ^2)]`$, $`\chi /E=\mathrm{cot}\chi /2(1h_x^2\mu ^2)`$. Doing the integrals, one obtains

$$\mathrm{\Lambda }=\pi J\left(1\frac{H_x}{\sqrt{1\lambda }H_c}\right)n\pi .$$
(22)   
The result for $`n=0`$ is the same as in ), while for $`n0`$ it is new. To order $`1/J`$, the vanishing points for higher pairs are the same as those for the lowest one. It must be remembered, however, that since we demanded $`\chi >0`$, one must have $`\mathrm{\Lambda }>0`$. (Finding $`\mathrm{\Lambda }<0`$ means that the oscillatory forbidden region has disappeared.) Thus the highest-field level crossing is successively eliminated as $`n`$ increases, and (including zero and negative values) there are $`2(Jn)`$ fields in all where $`\mathrm{\Delta }`$ vanishes.

###### Acknowledgements.

This work is supported by the NSF via grant number DMR-9616749.