format: text

text:
# Perturbative Formulation and Non-adiabatic Corrections in Adiabatic Quantum Computing Schemes

## Abstract

Adiabatic limit is the presumption of the adiabatic geometric quantum computation and of the adiabatic quantum algorithm. But in reality, the variation speed of the Hamiltonian is finite. Here we develop a general formulation of adiabatic quantum computing, which accurately describes the evolution of the quantum state in a perturbative way, in which the adiabatic limit is the zeroth-order approximation. As an application of this formulation, non-adiabatic correction or error is estimated for several physical implementations of the adiabatic geometric gates. A quantum computing process consisting of many adiabatic gate operations is considered, for which the total non-adiabatic error is found to be about the sum of those of all the gates. This is a useful constraint on the computational power. The formalism is also briefly applied to the adiabatic quantum algorithm.

Recently a considerable amount of attention has been paid to the idea of using geometric phases accumulated by an adiabatically time-dependent Hamiltonian to realize quantum gates lit1 ; lit2 ; lit3 . Construction of universal gates by geometric quantum teleportation was studied, with the analysis of errors from imperfect control ellinas . On the other hand, adiabatic evolution is also the basis of the so-called adiabatic quantum algorithms farhi , for which the speed and the overall time have been analyzed vandam .

For these quantum computing schemes to work, it was supposed that the adiabatic limit is retained. However, in practice, and particularly in the case of quantum computation, where the advantage lies in speedup and the operation time should be shorter than the decoherence time, the evolution is required to be completed in a finite period of time. Therefore, it is important to know the full picture of the evolution of the quantum state and the non-adiabatic correction, which gives rise to error if the adiabatic limit is necessary for the designed quantum computing scheme. Here we develop a general formulation of adiabatic quantum computing, applicable to the previously proposed quantum computing schemes and to any slowly varying Hamiltonian. As an adiabatic perturbation theory, it accurately describes the quantum evolution in a perturbative way, in which the adiabatic limit is the zeroth-order approximation. As an application, an examination is made on the non-adiabatic errors in several previously proposed implementations of the adiabatic geometric gates. We also investigate the non-adiabatic error in an entire quantum computing process consisting of many adiabatic gates, which has not been considered previously. Finally we briefly discuss the adiabatic quantum algorithm, noting that such an algorithm can still be implemented even if the non-adiabatic correction is not vanishingly small.

If the evolution of a time-dependent Hamiltonian is sufficiently slow, the adiabatic theorem tells that in the adiabatic limit and under such conditions as continuity, non-crossing and differentiability, an instantaneous eigenstate at an initial time evolves to a state close to the corresponding instantaneous eigenstate at a later time book .

In general, using the instantaneous eigenstate $`|\varphi _n(t)`$, one can always expand the state of the system $`|\psi (t)`$ as

$$|\psi (t)=\underset{n}{}a_n(t)|\varphi _n(t)\mathrm{exp}[i\eta _n(t)],$$
(1)
where $`\eta _n(t)=\frac{i}{\mathrm{}}_0^tE_n(\tau )𝑑\tau `$ is the dynamic phase. Then the Schrödinger equation $`i\mathrm{}_t|\psi (t)=H(t)|\psi `$ leads to

$$_ta_n(t)=\underset{m}{}a_m(t)\varphi _n(t)|_t\varphi _m(t)\mathrm{exp}[i\eta _m(t)i\eta _n(t)],$$
(2)
which, together with the initial condition $`a_n(0)\varphi _n(0)|\psi (0)`$, determines $`|\psi (t)`$.

First suppose $`|\psi (0)`$ is a non-degenerate eigenstate $`|\varphi _n(0)`$. Then in the adiabatic limit, one obtains $`|\psi (t)|\varphi _n(t)\mathrm{exp}[i\gamma _n(t)+i\eta _n(t)]`$, where $`\gamma _n(t)=_0^t\varphi _n|_\tau \varphi _n𝑑\tau =\varphi _n|_\mu \varphi _n𝑑x^\mu `$ is the geometric or Berry phase berry .

However, when non-adiabatic correction is considered, the exact state should be the solution of Eq. (2). For a slowly varying $`H(t)`$, using a perturbative approach, one can obtain

$$\begin{array}{c}U(t)|\varphi _n(0)=\mathrm{exp}[i\gamma _n(t)+i\eta _n(t)][|\varphi _n(t)+\\ \mathrm{}_{mn}\frac{|\varphi _m(t)\varphi _m(t)|_t\varphi _n(t)}{E_m(t)E_n(t)}]+\mathrm{}.\end{array}$$
(3)
In general, as in quantum computing, $`|\psi (0)=_na_n(0)|\varphi _n(0)`$ is a superposition of different eigenstates. Then linearity of quantum evolution implies

$$|\psi (t)=\underset{n}{}a_n(0)U(t)|\varphi _n(0),$$
(4)
where each $`U(t)|\varphi _n(0)`$ is as given in Eq. (3). Therefore, $`a_n(t)exp[i\eta _n(t)]=_m\varphi _n(t)|U(t)|\varphi _m(0)a_m(0)`$. From (3), one obtains $`\varphi _n(t)|U(t)|\varphi _n(0)\varphi _n(t)|U^{(0)}(t)|\varphi _n(0)=\mathrm{exp}[i\gamma _n(t)+i\eta _n(t)]`$, while for $`nm`$, $`\varphi _n(t)|U(t)|\varphi _m(0)\varphi _n(t)|U^{(1)}(t)|\varphi _m(0)=\mathrm{}\mathrm{exp}[i\gamma _m(t)+i\eta _m(t)]\varphi _n|_t\varphi _m/[E_n(t)E_m(t)]`$. Here $`U^{(k)}`$ refer $`k`$th order term. Since $`E_n(t)E_m(t)`$, $`\eta _n(t)\eta _m(t)`$. If one implements an all-geometric gate, in which the instantaneous basis states are $`|\varphi _n(t)`$ and $`|\varphi _m(t)`$, the difference between $`\eta _n(t)`$ and $`\eta _m(t)`$ needs to cancelled out by using a certain method lit1 .

In the presence of degeneracy of eigenstates, denote the eigenstates as $`|\varphi _{\alpha _n}^n(0)`$, where $`n`$ labels the energy levels, while $`\alpha _n`$ labels the different eigenstates in the subspace $`n`$. As the generalization of Eq. (3), we obtain

$$\begin{array}{c}U(t)|\varphi _{\alpha _n}^n(0)=\mathrm{exp}[i\eta _n(t)][|\chi _{\alpha _n}^n(t)+\\ \mathrm{}_{mn}_{\beta _m}\frac{|\chi _{\beta _m}^m(t)\chi _{\beta _m}^m(t)|_t\chi _{\alpha _n}^n(t)}{E_m(t)E_n(t)}]+\mathrm{},\end{array}$$
(5)
with $`|\chi _{\alpha _n}^n(t)=_{\beta _n}V_{\alpha _n\beta _n}^n(t)|\varphi _{\beta _n}^n(t)`$, $`V^n(t)=P\mathrm{exp}_0^tA^n(\tau )𝑑\tau `$, where $`A_{\alpha _n\beta _n}^n\chi _{\beta _n}^n(\tau )|_t\chi _{\alpha _n}^n(\tau )`$ is the connection in the subspace $`n`$. $`V^n`$ may be called non-abelian geometric phase or Wilczek-Zee (WZ) phase wz . In the zeroth order, $`U(t)`$ is block-diagonal, each block being a WZ phase in the subspace of a set of degenerate eigenstates. In the adiabatic limit, as a unitary transformation, a non-abelian geometric phase, i.e. the first term in (5), may be used to realize a quantum gate lit2 .

With the existence of degeneracy of eigenstates, a general superposition state can be written as

$$|\psi (t)=\underset{n}{}\underset{\alpha _n}{}a_{\alpha _n}^n(t)|\chi _{\alpha _n}^n(t)\mathrm{exp}[i\eta _n(t)].$$
(6)
By choosing an appropriate basis for each degenerate subspace, the initial state can always be expanded in such a way that its projection in each degenerate subspace is a single eigenstate $`|\varphi _{\beta _n}^n(0)`$, i.e. $`|\psi (0)=_na_{\beta _n}^n(0)|\chi _{\beta _n}^n(0)`$, with $`|\chi _{\beta _n}^n(0)=|\varphi _{\beta _n}^n(0)`$. Therefore

$$|\psi (t)=\underset{n}{}a_{\beta _n}^n(0)U(t)|\chi _{\beta _n}^n(0).$$
(7)
Thus $`a_{\alpha _n}^n(t)exp[i\eta _n(t)]=_m_{\beta _m}\chi _{\alpha _n}^n(t)|U(t)|\chi _{\beta _m}^m(0)a_{\beta _m}^m(0)`$. From (5), $`\chi _{\alpha _n}^n(t)|U(t)|\chi _{\beta _n}^n(0)\chi _{\alpha _n}^n(t)|U^{(0)}(t)|\chi _{\beta _n}^n(0)=\mathrm{exp}[i\eta _n(t)]\delta _{\alpha _n\beta _n}`$, while for $`nm`$, $`\chi _{\alpha _n}^n(t)|U(t)|\chi _{\beta _m}^m(0)\chi _{\alpha _n}^n(t)|U^{(1)}(t)|\chi _{\beta _m}^m(0)=\mathrm{}\mathrm{exp}[i\eta _m(t)]\chi _{\alpha _n}^n(t)|_t\chi _{\beta _m}^m(t)/[E_n(t)E_m(t)]`$.

Through this formulation, it becomes clear that the adiabatic quantum computing is based on $`\varphi _n(t)|U^{(0)}(t)|\varphi _n(0)`$ or $`\chi _{\alpha _n}^n(t)|U^{(0)}(t)|\chi _{\beta _n}^n(0)`$, with higher-order terms neglected. Besides, while the previous proposals of adiabatic geometric gates are based on closed paths, there is nothing in principle against using open paths, as far as the corresponding geometric phases can be detected wu . Another noteworthy point, which was not pointed out before, is that when the qubits under a gate operation is entangled with other qubits, the linearity of quantum evolution guarantees that the gate operation is still given by Eq. (3) or (5), where the eigenstates are those of this concerned gate; one may include in the coefficients $`a_n(0)`$ or $`a_{\alpha _n}^n(0)`$ the states of the other qubits projected in the same branch as the eigenstates of the gated qubits. This is crucial for the possibility that different adiabatic geometric gates can be networked.

There is a significant difference in the uses of (abelian) Berry phase and WZ phase to realize a quantum gate, under adiabatic limit. For a Berry phase gate, it is necessary to have $`d`$ non-degenerate states, where $`d`$ is the Hilbert space dimension of the gate. For a WZ phase gate, one intentionally restrict the gate in a single degenerate eigenspace. A quantum gate based on WZ phase is more advantageous than that based on Berry phase, on the aspect that for the former, in the adiabatic limit, the state is always an instantaneous eigenstate of the Hamiltonian, hence there is no dynamical phase difference between the basis states, and it is more stable against environmental perturbation.

The non-adiabatic correction or error at time $`t`$ is $`ϵ(t)[U(t)U^{(0)}(t)]|\psi (0)=_{k=1}^{\mathrm{}}U^{(k)}(t)|\psi (0)U^{(1)}(t)|\psi (0)`$. The adiabatic limit means $`|ϵ(T)|1`$. For a Berry phase gate, with $`|\psi (0)=_na_n(0)|\varphi _n(0)`$,

$$ϵ(t)\underset{n}{}a_n(0)\underset{mn}{}\varphi _m(t)|U^{(1)}(t)|\varphi _n(0)|\varphi _n(0).$$
(8)
For a WZ phase gate geometric gate at $`E_n`$, with $`|\psi (0)=_{\alpha _n}a_{\alpha _n}(0)|\varphi _{\alpha _n}^n(0)`$,

$$ϵ(t)\underset{\alpha _n}{}a_{\alpha _n}(0)\underset{mn}{}\underset{\beta _m}{}\chi _{\beta _m}^m(t)|U^{(1)}(t)|\chi _{\alpha _n}^n(0)|\varphi _{\alpha _n}^n(0).$$
(9)
Note that the first order correction at time $`t`$ is determined only by eigenvalues, eigenstates and their time derivatives at $`t`$, hence is history-independent. This simplifies the analysis. The time derivatives do depend on the details of time-dependence. However, since only the path is specified cen , without the fine control of the dynamics, numerically it suffices to obtain the order of magnitude. The first-order correction is $`\mathrm{}/\mathrm{\Delta }T`$, where $`T`$ is the time duration of the gate operation, $`\mathrm{\Delta }`$ is the minimum energy gap with other eigenstates. It is the presumption of “slow variation” or perturbative approach that $`\mathrm{}/\mathrm{\Delta }T<1`$. The $`k`$-th order correction is $`(\mathrm{}/\mathrm{\Delta }T)^k`$.

As applications, we now apply the above results to several physical implementations previously proposed. The first proposal, based on Berry phase, uses NMR lit1 . The Hamiltonian, in the rotating frame, is $`H(t)=𝐑(t)𝐈`$, where $`𝐑(R_x,R_y,R_z)=\mathrm{}(\omega _1\mathrm{cos}\varphi ,\omega _1\mathrm{sin}\varphi ,\omega _0\omega )`$, $`𝐈=\frac{1}{2}(\sigma _x,\sigma _y,\sigma _z)`$, $`\omega _0`$ is proportional to the static magnetic field in $`z`$ direction, $`\omega _1`$ is proportional to the RF magnetic field in $`xy`$ plane, $`\omega `$ is its angular frequency, $`\varphi `$ is its initial phase. The instantaneous eigenstates is $`|1(t)=1/\sqrt{R}[(R_xiR_y)/\sqrt{RR_z}|+\sqrt{RR_z}|]`$, with eigenvalue $`R/2`$, and $`|0(t)=1/\sqrt{R}[(R_xiR_y)/\sqrt{R+R_z}|\sqrt{R+R_z}|]`$ with eigenvalue $`R/2`$, where $`R|𝐑|`$. From this, one obtains, for $`n=0,1`$, $`U_{nn}(0,t)=\mathrm{exp}[i\gamma _n(t)+i\eta _n(t)]`$. $`\eta _0(t)=Rt/2`$, $`\eta _1(t)=Rt/2`$. The Berry phase $`\gamma _n(t)`$ is, in the case of a cycle path $`C`$, $`\mathrm{\Omega }(C)/2`$ for $`n=0`$ and $`\mathrm{\Omega }(C)`$ for $`n=1`$, where $`\mathrm{\Omega }(C)`$ is the solid angel that $`C`$ subtends at $`R=0`$. It is straightforward to write down $`U_{01}(0,t)`$ and $`U_{10}(0,t)`$. For a gate operation of a period $`T`$, the order of magnitude of these two matrix elements, and thus the non-adiabatic correction, is $`\mathrm{}/RT=1/\sqrt{\omega _1^2+(\omega _0\omega )^2}T`$. The two-bit gate, of qubits $`a`$ and $`b`$, is effected by addition of the interaction $`\mathrm{}JI_{az}I_{bz}`$. For the conditional phases of qubit $`a`$, $`\omega _{a0}`$ shifted to to $`\omega _{a0}+JI_{bz}=\omega _{a0}\pm J/2`$, depending on the basis state $`|I_{bz}`$ of $`b`$. One can obtain the non-adiabatic corrections in the two subspaces, with the substitution of $`\omega _{a0}\pm J/2`$ for $`\omega _0`$ in $`R`$ above. For a gate as in lit1 , the gaps are of the order of several hundred Hertz, while $`T`$ is of the order of second, hence the non-adiabatic corrections are of the order of $`10^2`$.

This method was also applied to a Josephson junction circuit falci . The effective Hamiltonian is still as that for NMR, now with $`𝐑=(E_J\mathrm{cos}\alpha ,E_J\mathrm{sin}\alpha ,E_c(12n_{off}))`$, where $`E_J`$ and $`\alpha `$ are decided by the Josephson couplings of two junctions, $`E_c`$ is charging energy, $`2en_{off}`$ is the offset charge. In the charging regime, as used, $`E_cE_J`$. Thus the non-adiabatic correction is of the order of $`\mathrm{}/\sqrt{E_J^2+E_c^2(12n_{off})^2}T`$. Hence if $`12n_{off}`$ is not too small, the adiabatic condition is $`\mathrm{}/E_cT1`$, more relaxed than previously thought falci .

An implementation of WZ phase gate was proposed for trapped ions duan ; unanyan ; recati . The one-bit gates are based on the Hamiltonian $`H=\mathrm{}|e(\omega _00|+\omega _11|+\omega _aa|+h.c.)`$, One can find that the eigenstates are: $`|\varphi ^1=(\omega |e+\omega _0^{}|0+\omega _1^{}|1+\omega _a^{}|a)/\sqrt{2}\omega `$ with eigenvalue $`\mathrm{}\omega `$, where $`\omega =\sqrt{\omega _0^2+\omega _1^2+\omega _a^2}`$, $`|\varphi _\alpha ^0=(\omega _1|0\omega _a|1)/\sqrt{|\omega _0|^2+|\omega _1|^2}`$ and $`|\varphi _\beta ^0=(\omega _a\omega _0^{}|0+\omega _a\omega _1^{}|1(|\omega _0|^2+|\omega _1|^2)|a)/(\omega \sqrt{|\omega _0|^2+|\omega _1|^2})`$ with eigenvalue $`0`$, and $`|\varphi ^1=(\omega |e+\omega _0^{}|0+\omega _1^{}|1+\omega _a^{}|a)/\sqrt{2}\omega `$ with eigenvalue $`\mathrm{}\omega `$. The WZ phase gates are based on $`U^{00}`$, in terms of our notation. Using the instantaneous eigenstates and eigenvalues, the non-adiabatic correction is obtained as $`_{n=1,1}_{x=\alpha ,\beta }a_x(0)U_x^{n0}(T)|\varphi _x^0(0)`$, whose order of magnitude is of $`1/\omega T`$. The two-bit gate proposed there is only a Berry phase gate under the Hamiltonian duan  $`H_{jk}=\frac{\eta ^2}{\delta }[|\mathrm{\Omega }_1|^2\sigma _{j1}^{\varphi _1}\sigma _{k1}^{\varphi _1}+|\mathrm{\Omega }_a|^2\sigma _{ja}^{\varphi _a}\sigma _{ka}^{\varphi _a}]`$, where $`\sigma _{j\mu }^{\varphi _\mu }e^{i\varphi _\mu }|e_{jj}\mu |+h.c.`$, $`\varphi _1\varphi _a=\varphi /2`$, using the notations therein. The eigenstates are $`\varphi ^1=(|\mathrm{\Omega }_1|^2e^{i\varphi }|11+|\mathrm{\Omega }_a|^2|aa+\sqrt{|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4}|ee)/\sqrt{2(|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4)}`$ with eigenvalue $`\frac{\eta ^2}{\delta }\sqrt{|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4}`$, $`|\varphi ^0=(|\mathrm{\Omega }_a|^2|11+|\mathrm{\Omega }_1|^2e^{i\varphi }|aa)/\sqrt{|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4}`$ with eigenvalue $`0`$, and $`\varphi ^1=(|\mathrm{\Omega }_1|^2e^{i\varphi }|11+|\mathrm{\Omega }_a|^2|aa\sqrt{|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4}|ee)/\sqrt{2(|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4)}`$ with eigenvalue $`\frac{\eta ^2}{\delta }\sqrt{|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4}`$. It was proposed to use $`|\varphi ^0`$ to implement the phase gate. The non-adiabatic correction is of the order of $`\mathrm{}/T(\eta ^2/\delta )\sqrt{|\mathrm{\Omega }_1|^4+|\mathrm{\Omega }_a|^4}`$.

Similar proposals were also made in Josephson junction charge qubits choi ; faoro . For the Hamiltonian used in choi , there are an eigenstate with eigenvalue $`\sqrt{h^2+|J_1|^2+|J_2|^2}`$, two degenerate eigenstates with eigenvalue $`h`$, two degenerate eigenstates with eigenvalue $`h`$, which are used to implement the WZ phase gate, and one ground state with eigenvalue $`\sqrt{h^2+|J_1|^2+|J_2|^2}`$, where $`h=E_c(12n_{off})/2`$. Thus the non-adiabatic correction is of the order of $`\mathrm{}/(\sqrt{h^2+|J_1|^2+|J_2|^2}h)T`$. In the two-bit implementation, the eigenvalues are $`\sqrt{|J_b|^2+(2h)^2}`$, $`2h`$, $`0`$, $`2h`$, $`\sqrt{|J_b|^2+(2h)^2}`$. The eigenstates with eigenvalue $`2h`$ are used as the qubit states. The non-adiabatic correction is of the order of $`\mathrm{}/\mathrm{\Delta }T`$, where $`\mathrm{\Delta }`$ is the smaller one of $`\sqrt{|J_b|^2+(2h)^2}2h`$ and $`2h`$. Suppose the order of magnitude of Josephson energy is $`J`$. Then if $`12n_{off}`$ is close to $`1`$, the energy gap is of the order of $`J^2/E_c`$ in the single-bit gate, and is of the order of $`J^2/4E_c`$ in the two-bit gate. Since $`E_cJ`$, the energy gap is smaller than $`J`$. Hence compared with the case of falci , the adiabatic condition is harder to meet, i.e. the non-adiabatic correction is larger. On the other hand, if $`12n_{off}`$ is tuned to be very small, then the energy gap for the cases of both choi  and falci  are of the order of the Josephson energy.

In the one-bit gate in faoro , the energy eigenvalues are $`\delta E_c+\sqrt{(\delta E_c)^2+2J^2}`$, $`0`$, which is with twofold degeneracy and is used to implement the WZ gate, and $`\delta E_c\sqrt{(\delta E_c)^2+J^2}`$, where $`\delta E_c`$ is some charging energy difference, $`J^2=|J_L|^2+|J_M|^2+|J_R|^2`$, using the notations there. Thus the non-adiabatic correction is of the order of $`\mathrm{}/\mathrm{\Delta }T`$, where $`\mathrm{\Delta }`$ is the smaller one of $`|\delta E_c+\sqrt{(\delta E_c)^2+J^2}|`$ and $`|\delta E_c\sqrt{(\delta E_c)^2+J^2}|`$. Hence $`\mathrm{\Delta }`$ is of the order of $`\delta E_c`$ if $`\delta E_c>J`$, and is of the order of $`J`$ if $`\delta E_c<J`$. For the two-bit gate, the three energy eigenvalues are $`\sqrt{|J_X|^2+|J_M^{(2)}|^2}/2`$, $`0`$, and $`\sqrt{|J_X|^2+|J_M^{(2)}|^2}/2`$, where the parameters are as defined there. Hence the non-adiabatic correction is of the order of $`2\mathrm{}/\sqrt{|J_X|^2+|J_M^{(2)}|^2}T`$. For both one-bit and two-bit gates, the energy gap is at most of the order of Josephson energy. Therefore, though WZ gate has more advantages over Berry gate, the non-adiabatic error for choi ; faoro  is larger than that for falci .

A quantum computing process in a gate array consists of many gate operations on a large number of qubits, hence a complete estimation of error must include its scaling with the number of gate operations. Suppose from time $`0`$ to $`T`$, $`M_1`$ adiabatic gates, denoted as $`U^j(t)`$, $`(j=1,\mathrm{},M_1)`$, are in parallel operation, each on a small number of (say, one or two) qubits. For $`0tT`$, the entire quantum computer evolves as $`|\mathrm{\Psi }(t)=U^1(t)\mathrm{}U^{M_1}(t)|\mathrm{\Psi }(0)=_{i_1\mathrm{}i_{M_1}i_r}U^1(t)|\varphi _{i_1}^1(0)\mathrm{}U^{M_1}(t)|\varphi _{i_{M_1}}^{M_1}(0)|\varphi _{i_r}^r`$, where $`|\varphi _{i_j}^j`$ is a basis state of the qubits acted by the $`j`$-th gate, $`r`$ denotes the rest qubits, which are not operated by any gate during this period. We know that $`U^j(t)|\varphi _{i_j}^j(0)=U_0^j(t)|\varphi _{i_j}^j(0)+ϵ^j(t)`$, where $`U_0^j(t)`$ represents the adiabatic limit of $`j`$-th gate, while $`ϵ^j(t)<1`$ is its non-adiabatic correction. The state of the quantum computer at $`T`$ is then $`|\mathrm{\Psi }(T)=|\mathrm{\Psi }_0(T)+\sigma (0,T)`$, where $`|\mathrm{\Psi }_0(T)=U_0^1(T)\mathrm{}U_0^{M_1}(T)|\mathrm{\Psi }(0)`$, $`\sigma (0,T)_{j=1}^{M_1}ϵ^j(T)`$ is the first-order error of the entire quantum computer accumulated from time $`0`$ to $`T`$. Afterwards, during $`T<t<2T`$, the quantum computer is operated in parallel by $`M_2`$ gates, labelled as $`M_1+1,\mathrm{},M_1+M_2`$, to which the qubits are allocated in a way usually different from the period $`0<t<T`$. Using a derivation similar to the above, it can be obtained that $`|\mathrm{\Psi }(2T)=U_0^{M_1+1}(T)\mathrm{}U_0^{M_1+M_2}(T)|\mathrm{\Psi }(T)+\sigma (T,2T)`$, where $`\sigma (T,2T)_{j=M_1+1}^{M_1+M_2}ϵ^j(T)`$ is the error of the entire quantum computer accumulated from $`T`$ to $`2T`$. Therefore $`|\mathrm{\Psi }(2T)=U_0^{M_1+M_2}(T)\mathrm{}U_0^{M_1+1}(T)\mathrm{}U_0^{M_1}(T)\mathrm{}U_0^1|\mathrm{\Psi }(0)+\sigma (0,2T)`$, where $`\sigma (0,2T)\sigma (0,T)+\sigma (T,2T)_{j=1}^{M_1+M_2}ϵ^j(T)`$ is the total error at $`2T`$. Therefore, for a quantum computing process consisting of many gate operations, no matter how they are arranged in space and time, the total non-adiabatic error is, to the first-order, just the sum of the errors of all these gates.

Suppose for each adiabatic gate, the time duration $`T`$, and the minimum energy gap with other eigenstates $`\mathrm{\Delta }`$. Thus the lower bound of the non-adiabatic error for each gate is $`ϵ(T)\mathrm{}/\mathrm{\Delta }T`$. Hence the lower bound of the total error is $`\sigma Mϵ(T)`$, where $`M`$ is is total number of gate operations. For a quantum computing process to make sense, it is constrained that $`\sigma <1`$. Therefore, $`M<1/ϵ(T)=\mathrm{\Delta }T/\mathrm{}`$. In Shor’s algorithm, to factor a number $`N`$, $`M300(\mathrm{log}_{10}N)^3`$ steane . Therefore, adiabatic quantum computing can at most factor $`N10^{(1/300ϵ(T))^{1/3}}`$. For $`ϵ(T)`$ of the order $`10^2`$, $`N10`$.

Let us switch to the adiabatic quantum algorithm farhi , which is based on adiabatically varying the Hamiltonian from a beginning Hamiltonian $`H_b`$ at $`t=0`$ to a final one $`H_p`$ at $`t=T`$. Under the adiabatic limit, if the system starts with the ground state of $`H_b`$, it ends up as the ground state of $`H_p`$, which gives the solution to an optimization problem. For a finite varying rate of the Hamiltonian, to the first order, the more accurate state is as given in Eq. (3). Hence to the first order approximation, the non-adiabatic correction at time $`T`$ is of the order of $`1/\mathrm{\Delta }_pT`$, where $`\mathrm{\Delta }_p`$ is the energy gap of $`H_p`$, which is independent of the specific path in which $`H_b`$ is evolved to $`H_p`$.

According to (3), as far the perturbative approach is valid, i.e. $`|\varphi _0(t)|_t\varphi _n(t)/(E_nE_0)|<1`$, the measurement shows that one of the eigenstates appears with probability clearly the largest. Then one can know that this state corresponds to the ground state and thus the solution to the problem. This is consistent with vandam .

To summarize, we developed a general, perturbative, formulation of the adiabatic quantum computing schemes, which perturbatively describe the accurate evolution of the state. It leads to a deeper understanding of related issues. The formalism is applied to analyze both the the adiabatic geometric quantum computation and the adiabatic quantum algorithm. The order of magnitude of the first-order non-adiabatic error is the inverse of the executing time times the minimum gap with other eigenstates. Several proposed physical implementations of the former are considered from this point of view. Different proposals based on charging Josephson junctions are compared. We also consider an entire quantum computing process consisting of many adiabatic gates, obtaining the lower bound of the non-adiabatic error, as an interesting constraint on the power of the quantum computation based on adiabatic geometric gates. One needs to enlarge the energy gap in order to reduce the non-adiabatic error and thus improve the computational power. For the adiabatic quantum algorithm, it is noted that it can be realized as far as the perturbative approach, rather than the rigorous adiabatic limit, is valid, hence the computational time may be appropriately shortened.

Y.S. was supported by a Humboldt Fellowship when visiting Germany.


text:
# The enhancement of ferromagnetism in uniaxially stressed diluted magnetic semiconductors

## I Introduction

Recently (see Tanaka ; Dietl  and references therein), great advances have been made in the problem of ferromagnetism (FMm) of p-doped diluted magnetic semiconductors (DMS). The question about influence of different physical phenomena in DMS on the critical temperature $`T_c`$ of ferromagnetic (FM) phase transition is of prime interest for searching for the future study trends. After Pashitskii and Ryabchenko prediction of FMm in DMS,PR  the competition between FM correlations mediated by indirect long-range spin-spin interaction and direct short-range antiferromagnetic (AFM) interaction is considered to be decisive in the formation of FM state in DMS. In other words, to obtain FM state with high enough $`T_c`$, we need to inhibit the AFM contribution.

Due to short-range character of AFM interaction, only close pairs of magnetic ions contribute to it. This contribution can be minimized by decreasing the magnetic ion concentration $`n_i`$. On the other hand, when $`n_i`$ is too small, ferromagnetism can be destroyed. To retain FM ordering in this case, we may increase the carriers concentration $`n_c`$ (in p-type DMS this corresponds to holes concentration $`n_h`$ ). The calculations of $`T_c`$ performed in a mean field approximation (MFA) PR  support this standpoint predicting an increase of $`T_c`$ as $`n_c^{1/3}n_i`$.

However, at large $`n_c`$, the Friedel oscillations of carrier spin polarization become significant so that MFA becomes inapplicable. More thorough calculations with respect to Friedel oscillations (i.e. beyond MFA) corroborate above statement and show that in DMS undergoing FM phase transition $`n_c`$ cannot exceed some critical value related to $`n_i`$. The reason for that is the oscillations of the RKKY interaction at the scale of $`1/k_Fn_c^{1/3}`$, which makes impossible long-range FM correlations if $`\overline{r}k_F1`$, ($`\overline{r}`$ is an inter-ion mean distance) SS .

To properly account for above Friedel oscillations, which is indeed a spatial dispersion of inter-ion interaction, we developed so-called random field approximation (RFA) in Ref.SS, . In that work, the ion-ion interaction has been considered in the context of Ising model. This model can be applied for axially symmetric semiconductors with magnetic ions interacting indirectly via holes (i.e. with RKKY interaction). The reason for Ising model usage for RKKY interaction is a complex structure of valence band in DMS that picks out the hole angular momenta projections $`J_Z=\pm 3/2`$ for lowest heavy hole subbands in crystals with distorted cubic or uniaxial lattice.

Additionally to spatial dispersion of inter-ion interaction, there are fluctuations of local magnetic field direction that also cannot be described in terms of MFA. This effect stems from the contribution of transversal spin components in effective Heisenberg - like Hamiltonian of spin-spin interaction. Namely, RKKY-interaction in p-doped undistorted DMS with cubic lattice represents this situation.

In the present paper we pay attention to the fact, that transition from Heisenberg spin-spin interaction to Ising one, i.e. exclusion of transversal spin components from Heisenberg Hamiltonian decreases the system entropy and therefore can enhance $`T_c`$. For quantitative description of this effect we present a comparative analysis of the RFA-theories for critical temperature $`T_c`$ in Heisenberg and Ising models. This analysis has been made to determine the role of directional fluctuations (inherent to Heisenberg model) of localized spins in a random magnetic field. The uniaxial stresses in typical DMS-structures grown on a substrate with some mismatch of lattice constants is shown to be the factor responsible for appearance of Ising-like interaction between magnetic ions spins. Thus, we predict a new effect implying that increasing of $`T_c`$ can be controlled within certain bounds by the uniaxial stresses in DMS. In other words, we expect that strain engineering can efficiently control the value of ferromagnetic phase transition temperature resulting from the hole-mediated exchange interaction between magnetic ions in DMS.

## II Theoretical background

### II.1 Heisenberg model

The Hamiltonian of Heisenberg model for DMS reads

$$=\underset{j<j^{}}{}J(\stackrel{}{r}_{j,j^{}})\stackrel{}{S}_j\stackrel{}{S}_j^{}+\underset{j}{}\stackrel{}{H}_0\stackrel{}{S}_j,$$
(1)
where external magnetic field $`\stackrel{}{H}_0`$ and interaction $`J(\stackrel{}{r}_{j,j^{}})`$ is measured in energy units (i.e. $`g\mu =1`$, $`\mu `$ is the Bohr magneton). The transition to Ising model means, hereafter, keeping only $`S_{jZ}S_{j^{}Z}`$ \- components in scalar product $`\stackrel{}{S}_j\stackrel{}{S}_j^{}`$.

The Hamiltonian (1) incorporates two kinds of randomness. First, (the spatial disorder) is that spin can be randomly present or absent in the specific $`j`$-th cite of a host semiconductor. Second, (the thermal disorder) is a random quantum state of a spin in $`j`$-th cite. These spatial and thermal fluctuations can be taken into consideration by introduction of random field rather than mean field.

In the random field approximation, we consider every spin $`\stackrel{}{S}_j`$ as a source of fluctuating (random) field

$$\stackrel{}{H}_{ri,j}=J(\stackrel{}{r}_i\stackrel{}{r}_j)\stackrel{}{S}_j$$
(2)
affecting other spin at the sites $`\stackrel{}{r}_i`$. In other words, every spin is subjected to some random (rather then mean) field, created by all other spins. So, all thermodynamic properties of the system will be determined by the distribution function $`f(\stackrel{}{H}_r)`$ of the random field $`\stackrel{}{H}_r`$. Namely, any spin dependent macroscopic quantity (like magnetization) $`<<A>>`$ reads

$$<<A>>=<A>_{\stackrel{}{H}_r}f(\stackrel{}{H}_r)𝑑\stackrel{}{H}_r,$$
(3)
where

$$<A>_{\stackrel{}{H}_r}=\frac{\mathrm{Tr}\left\{A\mathrm{exp}\left(_Z/T\right)\right\}}{\mathrm{Tr}\mathrm{exp}\left(_Z/T\right)}$$
(4)
is a single particle thermal average with temperature $`T`$ and effective Zeeman Hamiltonian $`_Z=\stackrel{}{H}_r\stackrel{}{S}`$.

The distribution function $`f(\stackrel{}{H}_r)`$ is defined as

$$f(\stackrel{}{H}_r)=\overline{\delta \left(\stackrel{}{H}_r+\underset{j(i)}{}J(\stackrel{}{r}_i\stackrel{}{r}_j)\stackrel{}{S}_j\stackrel{}{H}_0\right)},$$
(5)
where the bar means averaging over spatial disorder. Our RFA approach is based on micro-canonical statistical theory of magnetic resonance line shape.Stoneham  Latter theory assumes the additivity of local molecular field contributions $`\stackrel{}{H}_r=_j`$ $`\stackrel{}{H}_{ri,j}`$ of each particle $`j`$ (Eq. (2)) as well as the non-correlative spatial distributions of magnetic ions.

Latter assumptions with respect to spectral representation of $`\delta `$ function permit to transform Eq. (5) to the non-linear integral equation for $`f(\stackrel{}{H}_r)f(\stackrel{}{H})`$ in thermodynamic limit. Introducing the probability $`n_i(\stackrel{}{r})d^3\stackrel{}{r}`$ for small volume $`d^3\stackrel{}{r}`$ to be occupied by a particle, we obtain

$`f(\stackrel{}{H})={\displaystyle \mathrm{exp}\left[i\stackrel{}{\tau }(\stackrel{}{H}\stackrel{}{H}_0)𝒢\left(\stackrel{}{\tau }\right)\right]\frac{d^3\stackrel{}{\tau }}{(2\pi )^3}};`$ (6a)
$`𝒢\left(\stackrel{}{\tau }\right)={\displaystyle _V}\mathrm{\Psi }(\stackrel{}{r})n_i(\stackrel{}{r})d^3\stackrel{}{r},`$ (6b)
$`\mathrm{\Psi }(\stackrel{}{r})<<1\mathrm{exp}[iJ(\stackrel{}{r})\stackrel{}{S}\stackrel{}{\tau }]>>=`$
$`={\displaystyle W(\stackrel{}{H})f(\stackrel{}{H})d^3H},`$ (6c)
$`W(\stackrel{}{H})=<1\mathrm{exp}[iJ(\stackrel{}{r})\stackrel{}{S}\stackrel{}{\tau }]>_\stackrel{}{H}.`$ (6d)
Eqs (6) represent the integral equation for distribution function $`f(\stackrel{}{H})`$. In general case this equation can be solved only numerically.

However, in many cases (e.g. for $`T_c`$ or magnetization calculations) it is possible to avoid the solution of the integral equation since in these cases it is exactly reducible to the set of transcendental equations for macroscopic quantities like $`<<\stackrel{}{S}^n>>`$, $`1n2S`$ of the system. Simplest situation corresponds to the case $`S=1/2`$, where only magnetization $`\stackrel{}{}=g\mu <<\stackrel{}{S}>>`$ (or in dimensionless units $`\stackrel{}{m}=2<<\stackrel{}{S}>>`$; $`g`$ is g-factor of a magnetic ion) is a unique order parameter to be found.

Thus, in the case of $`H_0=0`$ and $`\stackrel{}{S}=\frac{1}{2}\stackrel{}{\sigma }`$ ($`\stackrel{}{\sigma }`$ are the Pauli matrices) Eq. (6a) takes the form

$`f(\stackrel{}{H})={\displaystyle }\mathrm{exp}\left[i\stackrel{}{\tau }\stackrel{}{H}\right]\times `$
$`\times \mathrm{exp}\left\{_0\left({\displaystyle \frac{\tau }{2}}\right)i{\displaystyle \frac{\stackrel{}{\tau }}{\tau }}\stackrel{}{m}_1\left({\displaystyle \frac{\tau }{2}}\right)\right\}{\displaystyle \frac{d^3\stackrel{}{\tau }}{(2\pi )^3}},`$ (7)
where

$$\stackrel{}{m}=\frac{\stackrel{}{H}}{H}\mathrm{tanh}\frac{H}{2T}f(\stackrel{}{H})d^3\stackrel{}{H};$$
(8)
$$_0(x)=_Vn(\stackrel{}{r})\left[1\mathrm{cos}\left(J(\stackrel{}{r})x\right)\right]d^3\stackrel{}{r};$$
(9)
$$_1(x)=_Vn(\stackrel{}{r})\mathrm{sin}\left(J(\stackrel{}{r})x\right)d^3\stackrel{}{r}.$$
(10)
Here $`\tau =|\stackrel{}{\tau }|\sqrt{\tau _x^2+\tau _y^2+\tau _z^2}`$. To derive (7) we have used following property of Pauli matrices:

$$\mathrm{exp}(\stackrel{}{\sigma }\stackrel{}{b})=\mathrm{cosh}b+\frac{\stackrel{}{\sigma }\stackrel{}{b}}{b}\mathrm{sinh}b,b|\stackrel{}{b}|.$$
(11)
Substitution of equation (7) into Eq.(8) results in a single closed equation for order parameter $`\stackrel{}{m}`$:

$`\stackrel{}{m}={\displaystyle \frac{1}{\left(2\pi \right)^3}}{\displaystyle }d^3\stackrel{}{H}{\displaystyle }d^3\stackrel{}{\tau }{\displaystyle \frac{\stackrel{}{H}}{H}}\mathrm{tanh}{\displaystyle \frac{H}{2T}}\times `$
$`\times \mathrm{exp}\left\{_0\left({\displaystyle \frac{\tau }{2}}\right)+i\stackrel{}{\tau }\stackrel{}{B}\right\},`$ (12)
where

$$\stackrel{}{B}=\stackrel{}{H}\stackrel{}{m}\frac{_1\left(\frac{\tau }{2}\right)}{\tau }.$$
(13)
To simplify the vector equation (12), we scalarwise multiply its left- and right-hand sides by $`\stackrel{}{m}`$ and integrate the resulting equation for scalar quantity $`m|\stackrel{}{m}|`$ over the angle between vectors $`\stackrel{}{H}`$ and $`\stackrel{}{m}`$. The final result (see Appendix A for details of its derivation) reads

$`m=6{\displaystyle _0^{\mathrm{}}}_{1/2}^H(t)e^{_0(t/2T)}\times `$
$`\times R_2\left(m_1\left({\displaystyle \frac{t}{2T}}\right)\right)dt,`$ (14)
$`_{1/2}^H(t)={\displaystyle \frac{1}{3}}\mathrm{csch}\pi t(1+\pi t\mathrm{coth}\pi t),`$ (15)
$`R_n(x)={\displaystyle \frac{x\mathrm{cos}x\mathrm{sin}x}{x^n}}.`$ (16)
Trivial solution $`m=0`$ of the equation (14) corresponds to paramagnetic phase. Under certain system parameters and temperatures, the equation (14) has nontrivial solution that determines the phase transition to the state with spontaneous magnetization.

To find the critical temperature $`T_c`$, we use the Landau theory with $`m`$ as an order parameter. For this purpose, we may derive (see Appendix B for details of derivation) the free energy of the system in the form

$`F^H=F_0+{\displaystyle \frac{1}{2}}m^2+6{\displaystyle _0^{\mathrm{}}}_{1/2}^H(t)e^{_0(t/2T)}\times `$
$`\times {\displaystyle \frac{\mathrm{sin}\left(m_1\left(\frac{t}{2T}\right)\right)}{m_1^2\left(\frac{t}{2T}\right)}}dt.`$ (17)
In the vicinity of $`T_c`$, the free energy (17) can be substituted by a Landau expansion

$`F_L^H`$ $`=`$ $`F_0+{\displaystyle \frac{1}{2}}m^2\left(12A_1^H\right)+{\displaystyle \frac{1}{20}}m^4A_3^H;`$ (18)
$`A_n^H`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}_{1/2}^H(t)e_1^{_0(t/2T)}_1^n\left({\displaystyle \frac{t}{2T}}\right)𝑑t,`$
where $`F_0`$ is a system free energy in a paramagnetic phase.

It should be noted here that contrary to conventional phenomenological Landau expansions of a free energy, the coefficients $`A_n^H`$ in the function (18) have been derived microscopically within our RFA approach. Free energy functions (17) and (18) give a possibility to describe the experimentally observable equilibrium thermodynamic characteristics (like magnetic susceptibility, specific heat etc) of the DMS both in paramagnetic and in ferromagnetic phases.

According to Landau theory of phase transitions, the phase transition temperature $`T_c`$ is reached, when coefficient $`12A_1^H=0`$ in Eq. (18). This is because $`T_c`$ is defined as a temperature, where nonzero infinitesimal magnetization appears. Obviously, the same equation can be obtained from the Eq. (14) for magnetization in the limit $`m0`$. The explicit form of the equation for $`T_cT_c^H`$ reads

$$1=2_0^{\mathrm{}}_S^H(t)_1\left(\frac{t}{2T_c^H}\right)e^{_0\left(t/2T_c^H\right)}𝑑t,$$
(19)
where $`S=1/2`$ in our case.

Actually, the Eq. (19) determines the $`T_c^H`$ as an implicit function of system parameters (like $`n_i`$, $`n_c`$ etc). This function can be considered as a phase diagram that separates the region of parameters where the ferromagnetic phase with $`m0`$ exists from that where $`m=0`$. Latter phase may be paramagnetic or spin glass phase. In principle, our RFA method permits to investigate this question. This study, however, is beyond the scope of the present paper.

The limit $`T_c^H0`$ in (19) gives the relation between parameters of the system, which determines the condition for ferromagnetic ordering to occur in DMS at $`T=0`$. The explicit form of this condition reads

$$1<\frac{4}{3\pi }_0^{\mathrm{}}\frac{_1\left(x\right)}{x}e^{_0\left(x\right)}𝑑x.$$
(20)
### II.2 Ising model

Let us consider now Ising model. In this case all effective magnetic fields are directed along $`OZ`$ axis, so that scalar product reduces to $`\stackrel{}{S}\stackrel{}{\tau }=S_Z\tau `$, $`S_Z=1/2,1/2`$ and Eq.(6a) becomes

$`f(H)={\displaystyle \frac{1}{2\pi }}{\displaystyle _{\mathrm{}}^{\mathrm{}}}e^{iH\tau 𝒢(\tau )}𝑑\tau ,`$
$`𝒢(\tau )=1{\displaystyle _V}n_i(\stackrel{}{r})\left(e^{i\tau J(\stackrel{}{r})S_Z}\right)d^3\stackrel{}{r}=`$
$`=`$ $`_0\left({\displaystyle \frac{\tau }{2}}\right)+im_1\left({\displaystyle \frac{\tau }{2}}\right),`$ (22)
where definition of $`m`$ is similar to Eq.(8):

$$m=_{\mathrm{}}^{\mathrm{}}\mathrm{tanh}\left(\frac{H}{2T}\right)f(H)𝑑H.$$
(23)
Multiplying Eq.(II.2) by $`\mathrm{tanh}(H/2T)`$ and integrating over $`H`$, we obtain the transcendental equation for order parameter $`m`$. The explicit form of this equation reads

$`m={\displaystyle _{\mathrm{}}^{\mathrm{}}}_{1/2}^I(t)\mathrm{exp}(_0\left({\displaystyle \frac{t}{2T}}\right))\times `$
$`\times \mathrm{sin}\left(m_1\left({\displaystyle \frac{t}{2T}}\right)\right)dt,`$ (24a)
where

$$_{1/2}^I(x)=\frac{1}{2\pi }_{\mathrm{}}^{\mathrm{}}\mathrm{tanh}\left(\frac{h}{2}\right)\mathrm{sin}\left(xh\right)𝑑h=\frac{1}{\mathrm{sinh}\pi x}.$$
(24b)
The equations for free energy and critical temperature $`T_c`$ can be obtained similarly to those for Heisenberg model. They read

$`F^I`$ $`=`$ $`F_0+{\displaystyle \frac{1}{2}}m^2{\displaystyle _{\mathrm{}}^{\mathrm{}}}_{1/2}^I(\pi t)\mathrm{exp}(_0\left({\displaystyle \frac{t}{2T}}\right))\times `$
$`\times `$ $`{\displaystyle \frac{1\mathrm{cos}\left(m_1\left(\frac{t}{2T}\right)\right)}{_1\left(\frac{t}{2T}\right)}}dt,`$
$`F_L^I`$ $`=`$ $`F_0+{\displaystyle \frac{1}{2}}m^2\left(1A_1^I\right)+{\displaystyle \frac{1}{24}}m^4A_3^I,`$
$`A_n^I`$ $`=`$ $`2{\displaystyle _0^{\mathrm{}}}_{1/2}^I(t)\mathrm{exp}\left(_0\left({\displaystyle \frac{t}{2T}}\right)\right)_1^n\left({\displaystyle \frac{t}{2T}}\right).`$
Similar to Eq. (19), the explicit form of equation for $`T_cT_c^I`$ reads

$$1=2_0^{\mathrm{}}_{1/2}^I(t)_1\left(\frac{t}{2T_c^I}\right)\mathrm{exp}\left(_0\left(\frac{t}{2T_c^I}\right)\right)𝑑t,$$
(25)
while the equation for FMm region in the phase diagram at $`T=0`$ has following form

$$1<\frac{2}{\pi }_0^{\mathrm{}}_1\left(x\right)e^{_0\left(x\right)}\frac{dx}{x}.$$
(26)
## III Discussion of RKKY interaction

### III.1 Non Gaussian fluctuations for $`S=1/2`$

Let us analyze the equations for critical temperatures for Heisenberg (19) and Ising (25) models in more details. The difference between them consists only in the form of kernels of integrals for Heisenberg (Eq. (15)) and Ising (Eq. (24b)) cases, so that equation for critical temperature in Heisenberg model can be transformed to that in Ising model (and vice versa) by replacement of $`_{1/2}^H(t)`$ with $`_{1/2}^I(t)`$.

We start the analysis of these equations from their MFA asymptotics. To get this asymptotics, the functions $`_0(\xi )`$ and $`_1\left(\xi \right)`$ in the Eqs (9), (10) should be expanded up to linear terms: $`_0(\xi )0`$, $`_1\left(\xi \right)\xi \overline{J}`$, where $`\overline{J}=_Vn_i(\stackrel{}{r})J(\stackrel{}{r})d^3r`$. After some algebra, the latter approximation allows to reduce the Eqs(19), (25) to the expressions for critical temperatures in Heisenberg $`T_M^{MF}(M=H)`$ or Ising $`(M=I)`$ models:

$$T_M^{MF}=2\overline{J}_0^{\mathrm{}}_{1/2}^M(t)t𝑑t=\frac{1}{4}\overline{J},$$
(27)
One can see that latter expression is identically the same to well-known MFA result for $`S=1/2`$

$$T_c^{MF}=\frac{1}{3}S(S+1)_Vn_i(\stackrel{}{r})J(\stackrel{}{r})d^3r.$$
(28)
The Eqs (27), (28) demonstrate also that MFA is independent of the choice of Heisenberg or Ising model, $`T_H^{MF}=T_I^{MF}=T_c^{MF}`$.

Next terms of expansion of the Eqs (9), (10) correspond to Gaussian asymptotics for distribution function of local fields. The purpose of subsequent analysis is to compare the (actual, i.e. non Gaussian) fluctuations of longitudinal components of random field with those of transversal ones.

Since our theory permits to find the distribution function $`f(H)`$ when a spatial dependence of $`J(\stackrel{}{r})`$ is assigned, we should specify a magnetic interaction in the system. Usually in the problems of carrier-induced ferromagnetism in DMS, the RKKY interactionKittelBook  is considered as an effective spin-spin exchange interaction resulting in FM ordering. To clarify the role of transversal spin fluctuations, here we use the simplest possible form of the interaction and neglect all possible factors that can influence $`J(\stackrel{}{r})`$ (such as nonparabolicity of carrier dispersion law etc, see RefsDietl01, ; McDon01, ; Mattis, ; Mauger, ; UmeKas,  for more details). Also, the stresses may change the form of $`J(\stackrel{}{r}),`$ see below for discussion.

In the case of simple one band carrier structure, the RKKY interaction reads

$$J(\stackrel{}{r})=J_0x_e^{4/3}R_4(2k_Fr),$$
(29)
where $`x_e=n_c/N_0`$, $`J_0`$ $`=\left(\frac{3}{\pi }\right)^{1/3}\frac{3}{2\mathrm{}^2}J_{ci}^2\mathrm{\Omega }^{2/3}m_d`$, $`J_{ci}`$ is a carrier-ion exchange constant, $`N_0=1/\mathrm{\Omega }`$ is a concentration of the cation cites, $`m_d`$ is the density of states effective mass. Note that in our single band approximation, the effects of stress may influence on $`x_e`$ and $`k_F,`$ see Dietl  for details. The threshold temperature of ferromagnetic ordering in MFA now can be found by evaluation of integrals (28) with respect to (29):

$$T_c^{MF}=\frac{1}{24\pi }J_0x_i^{4/3}\nu ^{1/3}.$$
(30)
Here, the factor $`J_0x_i^{4/3}`$ is independent of carrier concentration, $`x_i=n_i\mathrm{\Omega }`$ is a molar fraction of the magnetic ions. The ratio of electron and magnetic ion concentrations $`\nu =n_c/n_i=x_e/x_i`$ plays a crucial role in our theory as a parameter separating the cases of relatively small fluctuations with $`\nu \nu _c`$ and that of large ones with $`\nu \nu _c`$; parameter $`\nu _c`$ is indeed a dimensionless critical concentration (corresponding to equality sign in expressions (20) and (26) for Heisenberg and Ising models respectively).

The functions $`_0(x)`$ and $`_1(x)`$ in Eqs.(9), (10) (with respect to Eq.(29)) assume following form in the case of homogeneous magnetic ions distribution, $`n_i(\stackrel{}{r})=n_i=const`$,

$`_{0,1}\left(\xi \right)`$ $`=`$ $`\phi _{0,1}\left(\xi \right)/6\pi \nu ,`$
$`\phi _0\left(\xi \right)`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}\left\{1\mathrm{cos}\left(\xi R_4(y)\right)\right\}y^2𝑑y,`$ (31)
$`\phi _1\left(\xi \right)`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}\mathrm{sin}\left(\xi R_4(y)\right)y^2𝑑y.`$ (32)
In the case of spin $`S=1/2`$, the result (30) of MFA can be recovered from Eqs (31), (32) by their expansion up to linear terms $`\phi _0\left(\xi \right)0`$; $`\phi _1\left(\xi \right)\xi `$ (of course, with their further substitution into Eq. (14) or Eq. (24a)). Gaussian asymptotics of distribution function corresponds to the next term of expansion of the Eq. (31), $`\phi _0\left(\xi \right)\pi \xi ^2/30,\phi _1\left(\xi \right)\xi `$.

To account for real (non-Gaussian) distribution of fluctuating local field, we do not expand Eqs (31),(32) and calculate them numerically.

In dimensionless variables, the equations for critical temperatures for both above models assume following form

$`{\displaystyle \frac{1}{6\pi \nu }}{\displaystyle _{\mathrm{}}^{\mathrm{}}}_{1/2}^M(t)\phi _1\left({\displaystyle \frac{t}{2\theta _M}}\right)_{1/2}\left({\displaystyle \frac{t}{\theta _M}}\right)𝑑t=1,`$
$`_{1/2}\left({\displaystyle \frac{t}{\theta _M}}\right)=\mathrm{exp}\left({\displaystyle \frac{1}{6\pi \nu }}\phi _0\left({\displaystyle \frac{t}{2\theta _M}}\right)\right),`$ (33)
where $`\theta _M=T_c^M/(J_0x_e^{4/3})`$; $`M`$ stands for $`H`$ (Heisenberg model) or $`I`$ (Ising model). The result of calculation of $`T_c^M/T_c^{MF}`$ with the help of Eq. (33) as a function of $`\nu `$ is reported in the Fig.1a. It is seen, that there are curves that separate the areas of system parameters (including temperature) where FM or non-FM phases occur.

Our results show the limited area of concentrations $`0<\nu <\nu _c`$ which allow FM ordering in both considered models. So, we have found $`\nu _c=0.2473`$ for the Ising model and $`0.0989`$ for the Heisenberg model. Our results also show the critical character of dependence $`T_c=T_c(\nu )`$ that can be well approximated by the function

$$T_cT_c^{MF}\left(1\nu /\nu _c\right)^\lambda =\frac{J_0}{24\pi }x_i^{4/3}\nu ^{1/3}\left(1\nu /\nu _c\right)^\lambda $$
(34)
with $`\lambda 0.47`$ and $`\lambda 0.63`$ for Heisenberg and Ising models respectively.

Thus, the fluctuations of transversal spin components suppress a tendency towards FM ordering in the range of concentration ratios $`0<\nu <0.099`$. Moreover, in the range $`0.099<\nu <0.247`$, our results predict impossibility of FMm in DMS with Heisenberg-like spin-spin interaction while in this interval of $`\nu `$ FMm can still occur in DMS with Ising- like interaction between spins. Since Ising model is inherent to uniaxially stressed semiconductors, latter conclusion means that uniaxial distortion can effectively inhibit transversal spin fluctuations thus enhancing $`T_c`$.

### III.2 Gaussian fluctuations for $`S=5/2`$

The equation (19) for critical temperature $`T_c`$ is exact in the framework of our model, which means that it incorporates non-Gaussian fluctuations. At the same time it is restricted by the case of ion spins $`S=1/2`$. Practically important case of $`Mn^{2+}`$ ions with $`S=5/2`$ needs a special consideration. Mathematically, the problems arise for $`S>1/2`$, when exponential function of spin operator is no more a linear function so that Eq. (11) is no more valid. In this case we should use Sylvester theorem for exponent of a Hermitian matrix $`(2S+1)\times (2S+1)`$(see, e.g.Korn  for more details). Appendix C presents the final result of such calculations for $`W(\stackrel{}{H})`$ and $`S=5/2`$ (Eq.(6d)).

It is apparent from Appendix C that the case $`S=5/2`$ involves much more order parameters (such as $`a_{m,n}=c^m\eta _n=c^m\eta _nf(\stackrel{}{H})d^3H`$, $`m=0,2,4`$; $`n=0,1,2`$, and $`b_{m,n}=c^m\varphi _n`$, $`m=1,3,5`$; $`n=0,1,2`$) than in the case $`S=1/2`$. Not all of them play a crucial role in a formation of FM phase in DMS. For example, if the random field fluctuations are almost Gaussian, we can find the expectation values of spin operators (3) with the aid of only two parameters, representing first and second moments of the distribution function $`f(H)`$.

Our calculations show that Gaussian approximation for $`f(H)`$ adequately describes the actual phase diagram of DMS (except for close vicinity of critical concentration which is not practically important). That is why the analysis of the case of $`S=5/2`$ can be performed with sufficient accuracy in Gaussian approximation. To do this, the Eq. (67) should be substituted by its expansion in $`t=J(\stackrel{}{r})\tau `$ up to the second order. The result reads

$$W(\stackrel{}{H})=ic\left(3\mathrm{coth}\frac{3H}{T}\frac{1}{2}\mathrm{coth}\frac{H}{2T}\right)J(\stackrel{}{r})\tau +\left(\frac{5}{4}\left(1+4c^2\right)\frac{\left(3c^21\right)\left(1+4\mathrm{cosh}\frac{H}{T}\right)}{1+2\mathrm{cosh}\frac{2H}{T}}\right)\frac{\left(J(\stackrel{}{r})\tau \right)^2}{2}.$$
(35)
We apply this result to obtain the equation for critical temperature. The magnetization and mean value of random field $`\stackrel{}{H}`$ are negligibly small at this temperature that suggests that $`f(\stackrel{}{H})`$ is isotropic. Thus, performing angular averaging over the directions of $`\stackrel{}{H}`$ in Eq.(6c), we can put $`c^2=1/3`$ in the Eq. (35).

To relate this result to the case $`S=1/2`$, we note that first term of Eq. (35) includes Brillouin function for spin $`S=5/2`$

$$B_S\left(\frac{SH}{T}\right)=(1+\frac{1}{2S})\mathrm{coth}\frac{(S+1/2)H}{T}\frac{1}{2S}\mathrm{coth}\frac{H}{2T}.$$
(36)
After simple algebra, the Eq.(6c) takes following form for arbitrary spin $`S`$

$$\mathrm{\Psi }(\stackrel{}{r})=i\frac{\stackrel{}{\tau }}{\tau }\stackrel{}{m}SJ(\stackrel{}{r})+\frac{1}{3}S(S+1)\frac{\left(J(\stackrel{}{r})\tau \right)^2}{2},$$
(37)
where

$$\stackrel{}{m}=\frac{\stackrel{}{H}}{H}B_S\left(\frac{SH}{T}\right)f(\stackrel{}{H})d^3H.$$
(38)
One can see that the above equations are formally similar to those for $`S=1/2`$ if we expand Eqs (9) and (10) up to first nonvanishing terms. The Eq. (37) determines the components of Fourier image of distribution function (7) in Gaussian approximation

$`_0\left({\displaystyle \frac{\tau }{2}}\right)={\displaystyle \frac{1}{3}}S(S+1)\overline{J^2}{\displaystyle \frac{\tau ^2}{2}};_1\left({\displaystyle \frac{\tau }{2}}\right)=S\overline{J}\tau ,`$
$`\overline{J^n}={\displaystyle _V}n(\stackrel{}{r})J^n(\stackrel{}{r})d^3r.`$ (39)
Equations (38) and (39) reduce the problem of $`T_c`$ determination in Gaussian approximation for arbitrary ion spin to the case of non-Gaussian fluctuations for spin $`S=1/2`$ considered above. Namely, after substitution of these equations to Eq. (19) along with proper generalization (for the case of arbitrary spin) of the function $`_{1/2}^H(t)`$, we can use this equation to find $`T_c`$ for any $`S>1/2`$. Calculations performed in the same manner as Eqs. (56) and (57) yield

$$_S^H(t)=\frac{1}{6S}\left(\mathrm{coth}\frac{\pi t}{2S+1}\mathrm{coth}\pi t\right)+\frac{\pi t}{6S}\left(\frac{\mathrm{csch}^2\frac{\pi t}{2S+1}}{2S+1}\mathrm{csch}^2\pi t\right).$$
(40)
It can be readily shown that in the case $`S=1/2`$ the Eq. (40) reduces to simpler form (15).

Substitution of (39) to (20) permits to obtain the necessary condition for FMm formation at zeroth temperature in the form (see also Eqs. (20) and (26))

$$\overline{J}>\sqrt{\frac{3\pi }{8}\frac{S+1}{S}\overline{J^2}}.$$
(41)
This inequality has clear physical meaning: small dispersion and positive (FM) inter-ion spin-spin interaction favors creation of FMm in DMS. Also, large spin is also preferable for FMm formation due to relatively small (quantum) fluctuations of its transversal components.

For comparison, we also consider Ising model for $`S1/2`$ in Gaussian approximation (despite the possibility to account for non-Gaussian fluctuations of random field for arbitrary $`S1/2`$).SS  Expansion of the function $`𝒢(\tau )`$ (Eq.(22)) with arbitrary $`S`$ and $`m=_{\mathrm{}}^{\mathrm{}}\mathrm{B}_S\left(\frac{SH}{T}\right)f(H)𝑑H`$ up to second order leads to the distribution function in the form (II.2) with $`𝒢(\tau )`$ in the form (22), where

$$_0\left(\frac{\tau }{2}\right)=M_2\overline{J^2}\frac{\tau ^2}{2};_1\left(\frac{\tau }{2}\right)=S\overline{J}\tau .$$
(42)
Parameters

$$m=_{\mathrm{}}^{\mathrm{}}\mathrm{B}_S\left(\frac{SH}{T}\right)f(H)𝑑H$$
(43)
and

$`M_2=S_Z^2={\displaystyle _{\mathrm{}}^{\mathrm{}}}\mathrm{Q}_S\left({\displaystyle \frac{H}{T}}\right)f(H)𝑑H,`$ (44)
$`\mathrm{Q}_S\left({\displaystyle \frac{H}{T}}\right)=S(S+1)`$
$`S\mathrm{coth}\left({\displaystyle \frac{H}{2T}}\right)\mathrm{B}_S\left({\displaystyle \frac{SH}{T}}\right),`$ (45)
should be found self-consistently with respect to Eqs. (II.2), ( 22) and (42). Thus, in Gaussian approximation, Ising model needs self-consistent determination of two parameters, $`S_Z`$ and $`S_Z^2`$, while in the case of Heisenberg model single parameter $`<<|\stackrel{}{S}|>>`$ needs to be determined. This is because for Heisenberg-like interaction $`<<\stackrel{}{S}^2>>=S(S+1)=const`$ (Eq. (39)) so that this parameter depend neither on temperature nor on field distribution.

The critical temperature can be found from the Eqs (43) and (44) as $`m0`$

$`1={\displaystyle \frac{\overline{J}\sqrt{8}}{\sqrt{\pi M_2\overline{J^2}}}}\times `$
$`\times {\displaystyle _0^{\mathrm{}}}xS\mathrm{B}_S\left({\displaystyle \frac{3\sqrt{2M_2\overline{J^2}}}{S(S+1)\overline{J}}}{\displaystyle \frac{x}{\theta }}\right)e^{x^2}dx;`$ (46)
$`M_2={\displaystyle \frac{2}{\sqrt{\pi }}}{\displaystyle _0^{\mathrm{}}}\mathrm{Q}_S\left({\displaystyle \frac{3\sqrt{2M_2\overline{J^2}}}{S(S+1)\overline{J}}}{\displaystyle \frac{x}{\theta }}\right)e^{x^2}𝑑x,`$ (47)
where $`\theta =T_c/T_c^{MF}`$ is a ratio of actual critical temperature and that obtained in MFA approximation (Eq. (28)). In this approximation, the necessary condition to form FM state at $`T=0`$ in DMS has the form of following inequality

$$\overline{J}>\sqrt{\frac{\pi }{2}\overline{J^2}}.$$
(48)
The independence of Eq. (48) of spin can be thought of (see Eq. (41) for comparison) as a lack of transversal spin components contribution to random field in Ising model.

Now we are ready to compare Heisenberg and Ising models for specific case of RKKY interaction (29). In Gaussian approximation we have to evaluate only two integrals

$$\overline{J}=\frac{1}{6\pi }J_0x_i^{4/3}\nu ^{1/3},\overline{J^2}=\frac{2\pi ^2}{5}\nu \left(\overline{J}\right)^2.$$
(49)
Note that the equations (49) demonstrate the relation between Friedel oscillations and first ($`\overline{J}`$) and second ($`\overline{J^2}`$) moments of distribution function of random magnetic fields. Namely, at frequent Friedel oscillations (i.e. large $`n_c`$ and $`\nu `$) the influence of dispersion $`(\overline{J^2})^{1/2}`$ (which is ”responsible” for disorder in the system) prevail over trend to order the system due to mean value $`\overline{J}`$.SS  This means that frequent Friedel oscillations at the scale of mean inter-ion distance inhibit ferromagnetism in DMS.

Substitution of (49) into (46), (47) permits to find the dependence $`\theta (\nu )`$ for practically important case $`S=5/2`$ (see Fig.1b). One can see that qualitatively situations for $`S=5/2`$ and $`S=1/2`$ are similar: the region of FM state is significantly larger for Ising model than that for Heisenberg model. But there is also a difference. Comparison of Fig.1a and 1b shows that the area between curves $`T_c(n_c/n_i)`$ for Heisenberg (H) and Ising (I) models at low temperatures is smaller for $`S=5/2`$. This means that for $`S=5/2`$ quantum fluctuations (which is the only possible fluctuations at $`T=0`$) are inhibited as compared to the case $`S=1/2`$.

This is actually a reflection of the well known fact that the larger the magnitude of the quantum number of the spin, the ”more classical” it is, i.e. the smaller is the contribution of quantum fluctuations of its transversal components. At $`T0`$ additional thermal fluctuations appear. It is seen from the Figure, that at $`TT_c`$ the extension of FM phase due to enhancement of $`T_c`$ is almost the same for $`S=5/2`$ and $`S=1/2`$. This means that the thermal fluctuations of the spin do not sensitive to its value.

## IV Conclusions

In this paper we have proposed a new mechanism of the enhancement of FM phase transition temperature $`T_c`$ by the uniaxial distortion of DMS. This prediction is based on comparative analysis of Heisenberg (inherent to undistorted DMS) and Ising (inherent to uniaxially distorted DMS) models. The analysis of above models has been carried out in the framework of our recently developed formalism, SS  random field method. This method, which can be regarded as a substitution of conventional MFA for disordered systems with given $`J(\stackrel{}{r})`$, permits to derive self-consistently the equations for order parameter $`m`$ and the free energy functions of DMS.

Now we discuss in more details the influence of stresses on magnetic spin Hamiltonian of DMS (1). We consider such influence in two steps. In first step we consider the influence of the stress on the effective spin-spin interaction potential and in second step we consider the operator part of the Hamiltonian.

Since the effective potential of indirect spin-spin interaction strongly depends on the band structure of specific semiconductor sample, which is affected by the stress, this potential by itself may also depend on stress. This influence Dietl  is manifested via both density of states (in our case we use the effective mass of density of states $`m_d`$ (29)) at the Fermi level and concentration of free holes, related to pinning of the Fermi level by defects and impurities of different nature. As it was noted in Ref. Dietl01 , the influence of pressure on the density of states is small. The influence via pinning centers depends entirely on their nature. We can imagine the situation when the influence of pressure on the concentration is also small. For example, in the cases of the absence (or small number) of pinning centers or ”synchronous movement” of pinning centers with the valence band edge shift due to deformation potential, this effect is negligible and our mechanism of influence of the pressure will be decisive. Here we would like to emphasize that there are no general problems to incorporate the possible dependence of the concentration on the stress into our consideration (this is simply one more modification of $`J(\stackrel{}{r}_{j,j^{}})`$ in (1)). If we do so, the considered effect of elimination of transversal spin components by the stress, which has not been discussed in the literature, is an additional factor enhancing $`T_c.`$

To discuss the stress effect on the operator part, we note that in the typical case of p-type DMS with cubic lattice the transversal spin-spin interaction can be substantially inhibited by the uniaxial stress of a crystal. Really, such stress splits the valence band edge to light and heavy hole subbands. The resulting heavy hole (HH) states are characterized by angular momentum projection $`J_Z=\pm 3/2`$. For such states, the spin-flip scattering of these holes on magnetic ions is forbidden due to angular momentum conservation. Hence only longitudinal (along the distortion axis) $`Z`$-components of the spins will be present in the resulting indirect spin-spin interaction via above heavy holes.

However, aforementioned spin-flip processes are not forbidden both for light hole (LH) states (with spin projection $`\pm 1/2`$) and for transitions between LH and HH states thus contributing to the indirect interaction of transversal spin components of magnetic ions. Hence, the deformational splitting of a valence band edge, leading to preferential occupancy of the heavy hole subbands gives the anisotropy of indirect spin-spin interaction in the form $`=_{i>j}\left[J_{i,j}S_{Zi}S_{Zj}+J_{i,j}\left(S_{Xi}S_{Xj}+S_{Yi}S_{Yj}\right)\right]`$. The ratio $`\gamma =J_{i,j}/J_{i,j}`$ ($`0<\gamma <1`$) should be monotonically decreasing function of the ratio of heavy holes concentration $`n_{HH}`$ to their complete concentration $`n_H`$. Thus, if deformational HH - LH splitting exceeds the Fermi energy of the holes, the ratio $`n_{HH}/n_H`$ reaches its maximal value ($`n_{HH}/n_H1`$), corresponding to Ising Hamiltonian ($`\gamma 0`$) of the indirect spin-spin interaction. In this case the region of ferromagnetic state of DMS expands substantially towards both higher carriers (holes) concentration and higher temperatures.

The detailed theoretical description of all above effects, which is intimately related to the parameters of specific DMS sample and experimental conditions, can be developed within the framework of presented theory for any particular case. However, such calculations can be done only numerically. Note, that the magnitude of effect which we predicted is very sensitive to the holes concentration $`n_h`$. The problem of its correct determination, to the best of our knowledge, is related to the influence of anomalous Hall effect. Latter effect lowers reliability of $`n_h`$ determination from Hall effect data Ohno1 ; Ohno2 .

For better illustration of our effect, we estimate now the magnitude of $`T_c`$ increase for typical ferromagnetic DMS $`Ga_{1x}Mn_xAs`$ deposited on $`GaAs`$ and $`GaP`$ substrates. The mismatch $`\mathrm{\Delta }a`$ of lattice constant $`a`$ leads to biaxial strain that splits the valence band with deformation potential $`b=1.7`$ eV by the value $`\delta E_{1,2}=2\left|bϵ_{zz}^{}\right|`$,BirPikus  where $`zz`$\- component of strain tensor $`ϵ_{zz}^{}=2(\mathrm{\Delta }a/a)c_{12}/c_{11}`$, and the ratio of elastic moduli in GaAs $`c_{12}/c_{11}=0.453`$. For $`x=0.035`$ (or concentration $`n_i=7.7610^{20}`$ cm<sup>-3</sup>) , the relative mismatch $`\left|\mathrm{\Delta }a/a\right|=0.002`$ for $`GaAs`$ substrate and $`\left|\mathrm{\Delta }a/a\right|=0.036`$ for $`GaP`$ substrate.Dietl01  We can see that for $`GaAs`$ substrate small valence band splitting $`\delta E_{1,2}6`$ meV cannot suppress interaction of transversal $`Mn`$-spin components for typical concentration $`n_c=10^{20}`$ cm<sup>-3</sup> which corresponds to Fermi energy $`\epsilon _F80`$ meV, whereas for $`GaP`$ substrate $`\delta E_{1,2}109`$ meV $`>\epsilon _F`$. Thus, our mechanism predicts the enhancement of $`T_c`$ for $`Ga_{1x}Mn_xAs`$ on $`GaP`$ substrate by the factor 1.64 (see Fig.1b for $`\nu =n_c/n_i=0.13`$) as compared to the same DMS but on $`GaAs`$ substrate.

Let us finally note that in the present paper we considered the enhancement of $`T_c`$ due to RKKY interaction only. But there are also other mechanisms, which can lead to appearance of ferromagnetism in DMS, see Dietl  for details. These mechanisms will be eventually reduced to the Hamiltonian (1) with modified potential $`J(\stackrel{}{r}_{j,j^{}})`$. Thus for quantitative discussion of these mechanisms it is sufficient to substitute the corresponding modified potential to our self-consistent equations.

## Appendix A

We begin with equation (12) for magnetization.

$$\stackrel{}{m}=\frac{1}{\left(2\pi \right)^3}d^3\stackrel{}{H}d^3\stackrel{}{\tau }\frac{\stackrel{}{H}}{H}\mathrm{tanh}\frac{H}{2T}\mathrm{exp}\left\{_0\left(\frac{\tau }{2}\right)+i\stackrel{}{\tau }\stackrel{}{B}\right\},$$
(50)
First, we pass to spherical system for $`\stackrel{}{\tau }`$ and integrate the Eq.(50) over $`\stackrel{}{\tau }`$ directions. This yields

$`\stackrel{}{m}`$ $`=`$ $`{\displaystyle \frac{1}{2\pi ^2}}{\displaystyle d^3\stackrel{}{H}_0^{\mathrm{}}\tau ^2𝑑\tau \frac{\stackrel{}{H}}{H}\mathrm{tanh}\frac{H}{2T}\mathrm{exp}\left[_0\left(\frac{\tau }{2}\right)\right]\frac{\mathrm{sin}B\tau }{B\tau }},`$ (51)
$`B`$ $`=`$ $`|\stackrel{}{B}|=\sqrt{H^22\mathrm{cos}\theta mH_1(\tau /2)/\tau +\left(m_1(\tau /2)/\tau \right)^2},`$ (52)
where $`\theta `$ is an angle between vectors $`\stackrel{}{m}`$ and $`\stackrel{}{H}`$.

Next step is a scalar multiplication of its both sides by $`\stackrel{}{m}`$ that yields

$$m^2=\frac{1}{2\pi ^2}d^3\stackrel{}{H}_0^{\mathrm{}}\tau ^2𝑑\tau \frac{\stackrel{}{H}\stackrel{}{m}}{H}\mathrm{tanh}\frac{H}{2T}\mathrm{exp}\left[_0\left(\frac{\tau }{2}\right)\right]\frac{\mathrm{sin}B\tau }{B\tau }$$

or

$$m=\frac{1}{2\pi ^2}d^3\stackrel{}{H}_0^{\mathrm{}}\tau ^2𝑑\tau \mathrm{cos}\theta \mathrm{tanh}\frac{H}{2T}\mathrm{exp}\left[_0\left(\frac{\tau }{2}\right)\right]\frac{\mathrm{sin}B\tau }{B\tau }$$
(53)
The rotational invariance of scalar product permits to point $`\stackrel{}{H}`$ along $`z`$ axis and integrate over angular variables in (53). This yields

$$m=\frac{1}{2\pi ^2}_0^{2\pi }𝑑\phi _0^\pi \mathrm{sin}\theta d\theta _0^{\mathrm{}}H^2𝑑H_0^{\mathrm{}}\tau ^2𝑑\tau \mathrm{cos}\theta \mathrm{tanh}\frac{H}{2T}\mathrm{exp}\left[_0\left(\frac{\tau }{2}\right)\right]\frac{\mathrm{sin}B\tau }{B\tau },$$

where $`B`$ is defined by (52). Consider

$$I=\frac{1}{\tau }_0^{2\pi }𝑑\phi _0^\pi \mathrm{sin}\theta \mathrm{cos}\theta \frac{\mathrm{sin}B\tau }{B}d\theta .$$

Change of variables

$$\mathrm{cos}\theta =z=\frac{H^2+\left(m_1(\tau /2)/\tau \right)^2B^2}{2mH_1(\tau /2)/\tau }$$

with the help of (52) reduces it to the form

$`I`$ $`=`$ $`{\displaystyle \frac{\pi }{\tau \left(mH_1(\tau )/\tau \right)^2}}{\displaystyle _{H+m_1(\tau )/\tau }^{Hm_1(\tau )/\tau }}\left[H^2+\left(m_1(\tau /2)/\tau \right)^2B^2\right]\mathrm{sin}B\tau dB=`$ (54)
$`=`$ $`4\pi R_2\left(m_1(\tau /2)\right)R_2(H\tau ),R_n(x)={\displaystyle \frac{x\mathrm{cos}x\mathrm{sin}x}{x^n}}.`$
Substitution of (54) into (53) gives

$$m=\frac{2}{\pi }_0^{\mathrm{}}𝑑H_0^{\mathrm{}}𝑑\tau \mathrm{tanh}\frac{H}{2T}\mathrm{exp}\left[_0\left(\frac{\tau }{2}\right)\right]R_2\left(mH_1(\tau /2)\right)R_0(H\tau ).$$
(55)
It is also possible to integrate over $`H`$ in (55)

$$_0^{\mathrm{}}\mathrm{tanh}\frac{H}{2T}R_0(H\tau )𝑑H=3\pi T_{1/2}^H(t),t=\tau T.$$
(56)
$$_{1/2}^H(t)=\frac{\mathrm{sinh}\pi t+\pi t\mathrm{cosh}\pi t}{3\mathrm{sinh}^2\pi t}.$$
(57)
With respect to substitution $`\tau T=t`$ this gives final equation for $`m`$ in the form (14) of the text.

## Appendix B

We start the derivation of the free energy from the equation (14) for magnetization (order parameter). We rewrite it in the form

$$m+6_0^{\mathrm{}}_{1/2}^H(\pi t)e^{_0(t/2T)}R_2\left(m_1\left(\frac{t}{2T}\right)\right)𝑑t=0.$$
(58)
Now we recollect that if we have a free energy $`F`$ of a system, then the equation for order parameter (in our case Eq. (58)) should minimize it. In other words, Eq. (58) should be equivalent to condition

$$\frac{F}{m}=0.$$
(59)
Condition (59) is a simple differential equation for $`F`$, its solution yields

$`F`$ $`=`$ $`{\displaystyle 𝑑m\left\{m+6_0^{\mathrm{}}_{1/2}^H(\pi t)e^{_0(t/2T)}R_2\left(m_1\left(\frac{t}{2T}\right)\right)𝑑t\right\}}=`$ (60)
$`=`$ $`F_0+{\displaystyle \frac{1}{2}}m^2+6{\displaystyle _0^{\mathrm{}}}_{1/2}^H(\pi t)e^{_0(t/2T)}{\displaystyle \frac{\mathrm{sin}\left(m_1\left(\frac{t}{2T}\right)\right)}{m_1^2\left(\frac{t}{2T}\right)}}𝑑t.`$
This is indeed equation (17) from the text.

To get Landau expansion of (60), we simply expand $`\mathrm{sin}\left(m_1\left(\frac{t}{2T}\right)\right)`$ in Taylor series at small $`m`$. This gives

$`F`$ $`=`$ $`F_0+6{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{_{1/2}^H(\pi t)e^{_0(t/2T)}}{_1\left(\frac{t}{2T}\right)}}𝑑t+{\displaystyle \frac{1}{2}}m^2m^2A_1^H+{\displaystyle \frac{1}{20}}m^4A_3^H+\mathrm{},`$ (61)
$`A_n^H`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}_{1/2}^H(\pi t)e_1^{_0(t/2T)}_1^n\left({\displaystyle \frac{t}{2T}}\right)𝑑t.`$
Paying attention that second term in (61) does not depend on $`m`$ and hence just renormalizes $`F_0`$, we easily obtain Eq. (18) from the text.

## Appendix C

We are looking for the expression (6d) for the case $`S=5/2`$. We introduce notation $`\stackrel{}{t}=J(\stackrel{}{r})\stackrel{}{\tau }`$ and assume $`T=1`$ (i.e. $`\stackrel{}{H}`$ means $`\stackrel{}{H}/T`$). In these notations Eq.(6d) reads

$$W(\stackrel{}{H})=1\frac{\mathrm{Tr}\left\{\mathrm{exp}\left[i\stackrel{}{S}\stackrel{}{t}\right]\mathrm{exp}\left[\stackrel{}{S}\stackrel{}{H}\right]\right\}}{\mathrm{Tr}\left\{\mathrm{exp}\left[\stackrel{}{S}\stackrel{}{H}\right]\right\}}.$$
(62)
The denominator in Eq.(62) can be immediately evaluated in a reference frame rotating around quantization axis

$$\mathrm{Tr}\left\{\mathrm{exp}\left[\stackrel{}{S}\stackrel{}{H}\right]\right\}=\frac{\mathrm{sinh}3H}{\mathrm{sinh}H/2}.$$
(63)
We introduce the functions

$`\varphi _n={\displaystyle \frac{\mathrm{cosh}\left(nH\right)}{1+2\mathrm{cosh}\left(2H\right)}},`$
$`\eta _n={\displaystyle \frac{\mathrm{cosh}\left(nH\right)\left(\mathrm{cosh}H1\right)}{\mathrm{sinh}3H}}=`$
$`=`$ $`{\displaystyle \frac{\mathrm{cosh}\left(nH\right)\mathrm{sinh}H/2}{\mathrm{cosh}H/2+\mathrm{cosh}3H/2+\mathrm{cosh}5H/2}}`$
and cosine of the angle between vectors $`\stackrel{}{H}`$ and $`\stackrel{}{\tau }\stackrel{}{t}`$,

$$c=\mathrm{cos}\left(\widehat{\stackrel{}{H},\stackrel{}{\tau }}\right)=\frac{\left(\stackrel{}{H}\stackrel{}{\tau }\right)}{H\tau }.$$
(66)
In these notations after lengthy calculations we can get the expression for the trace in numerator of Eq. (62). With respect to Eq. (63), this expression assumes the form

$`W(\stackrel{}{H})`$ $`=`$ $`1{\displaystyle \frac{1}{4}}\mathrm{cos}{\displaystyle \frac{t}{2}}\left\{\left(314c^2+15c^4\right)\varphi _04\left(16c^2+5c^4\right)\varphi _1+5\left(1c^2\right)^2\varphi _2\right\}`$ (67)
$`{\displaystyle \frac{1}{8}}\mathrm{cos}{\displaystyle \frac{3t}{2}}\left\{\left(1+38c^245c^4\right)\varphi _0+4\left(112c^2+15c^4\right)\varphi _1+5\left(1+2c^23c^4\right)\varphi _2\right\}`$
$`{\displaystyle \frac{1}{8}}\mathrm{cos}{\displaystyle \frac{5t}{2}}\left\{\left(310c^2+15c^4\right)\varphi _0+4\left(15c^4\right)\varphi _1+\left(1+10c^2+5c^4\right)\varphi _2\right\}+`$
$`+{\displaystyle \frac{ic}{4}}\mathrm{sin}{\displaystyle \frac{t}{2}}\left\{\left(314c^2+15c^4\right)\eta _04\left(16c^2+5c^4\right)\eta _1+5\left(1c^2\right)^2\eta _2\right\}`$
$`{\displaystyle \frac{ic}{8}}\mathrm{sin}{\displaystyle \frac{3t}{2}}\left\{\left(326c^2+15c^4\right)\eta _04\left(34c^2+5c^4\right)\eta _15\left(32c^2c^4\right)\eta _2\right\}+`$
$`+{\displaystyle \frac{ic}{8}}\mathrm{sin}{\displaystyle \frac{5t}{2}}\left\{\left(1510c^2+3c^4\right)\eta _0+4\left(5c^4\right)\eta _1+\left(5+10c^2+c^4\right)\eta _2\right\}.`$