Quantum Spin-Wave Theory for Non-Collinear Spin Structures, a Review

Abstract

In this review, we trace the evolution of the quantum spin-wave theory treating non-collinear spin configurations. Non-collinear spin configurations are consequences of the frustration created by competing interactions. They include simple chiral magnets due to competing nearest-neighbor (NN) and next-NN interactions and systems with geometry frustration such as the triangular antiferromagnet and the Kagomé lattice. We review here spin-wave results of such systems and also systems with the Dzyaloshinskii–Moriya interaction. Accent is put on these non-collinear ground states which have to be calculated before applying any spin-wave theory to determine the spectrum of the elementary excitations from the ground states. We mostly show results obtained by the use of a Green’s function method. These results include the spin-wave dispersion relation and the magnetizations, layer by layer, as functions of T in 2D, 3D and thin films. Some new unpublished results are also included. Technical details and discussion on the method are shown and discussed.

1. Introduction

In a solid the interactions between its constituent atoms or molecules gives rise to elementary excitations from its ground state (GS) when the temperature increases from zero. One has examples of elementary excitations due to atom-atom interactions, known as phonons, or due to spin-spin interactions, known as magnons. Note that magnons are spin waves (SW) when they are quantized. Elementary excitations are defined also for interactions between charge densities in plasma, or for electric dipole-dipole interactions in ferroelectrics, among others. Elementary excitations are thus collective motions which dominate the low-temperature behaviors of solids in general.

For a given system, there are several ways to calculate the energy of elementary excitations from classical treatments to quantum ones. Since those collective motions are waves, their energy depends on the wave vector $\mathbf{k}$. The $\mathbf{k}$-dependent energy is often called the SW spectrum for spin systems. Note that though the calculation of the SW spectrum is often for periodic crystalline structures, it can also be performed for symmetry-reduced systems such as in thin films or in semi-infinite solids in which the translation symmetry is broken by the presence of a surface.

In this review we focus on the SW excitations in magnetically ordered systems. The history began with ferromagnets and antiferromagnets with collinear spin GSs, parallel or antiparallel configurations in the early 1950s. Most of the works on the SW used either the classical method or the quantum Holstein–Primakoff transformation. The Green’s function (GF) technique has also been introduced in a pioneering paper of Zubarev [1]. The first application of this method to thin films has been done [2]. Note that unlike the SW theory, the GF can treat the SW up to higher temperatures. We will come back to this point later.

Let us recall some important breakthroughs in the study of non-collinear spin configurations. The first discovery of the helical spin configuration has been published in 1959 [3,4]. Some attempts to treat this non-collinear case have been done in the 1970 and 1980. Let us cite two noticeable works on this subject in Refs. [5,6]. In these works, a local system of spin coordinates have been introduced in the way that each spin lies on its quantization axis. One can therefore use the commutation relations between spin deviation operators. These works took into account magnon-magnon interactions by expanding the Hamiltonian up to three-operator terms at temperature T = 0 [5] or up to four-operator terms at low T [6]. Nevertheless, since these works used the Holstein–Primakoff method, the case of higher T cannot be dealt with. In Ref. [7], the GF method has been employed for the first time to calculate the SW spectrum in a frustrated system where the GS spin configuration is non-collinear. Using the SW spectrum, the local order parameters and the specific heat were calculated. Since this work, we have applied the GF method to a variety of systems where the GS is non-collinear. In this review, we will recall results of some of these published works.

Let us comment on the frustration which is the origin of the non-collinear GS. The frustration is caused by either the competing interactions in the system or a geometry frustration as in the triangular lattice with only the antiferromagnetic interaction between the nearest neighbors (NN) (see Ref. [8]). The frustration causes high GS degeneracy, and for the vector spins (XY and Heisenberg cases) the spin configurations are non-collinear, making the calculation of the SW spectrum harder. A number of examples will be shown in this review paper.

In addition to competing interactions, the Dzyaloshinskii–Moriya (DM) interaction [9,10] is also the origin of non-collinear spin configurations in spin systems. While the Heisenberg model between two spins is written as $-J_{ij}{\mathbf{S}}_i·{\mathbf{S}}_j$ giving rise to two collinear spins in the GS, the DM interaction is written as ${\mathbf{D}}_{ij}·{\mathbf{S}}_i\times {\mathbf{S}}_j$ giving rise to two perpendicular spins. The DM model was historically proposed to explain the phenomenon of weak ferromagnetism observed in Mn compounds [11]. However, the DM interaction is at present known in various materials, in particular at the interface of a multilayer [12,13,14,15,16]. Although in this review we do not show the effect of the DM interaction in a magnetic field which gives rise to topological spin swirls known as skyrmions, we should mention a few of the important works given in Refs. [17,18,19,20,21]. Skyrmions are among the most studied subjects at the time being due to their potential applications in spin electronics [22]. We refer the reader to the rich biography given in our recent papers in Refs. [23,24].

Since this paper is a review on the method and the results of published works on SW in non-collinear GS spin configurations, it is important to recall the method and show main results of some typical cases. We would like to emphasize that, on the GF technique, to our knowledge there are no authors other than us working with this method. Therefore, the works mentioned in the references of this paper are our works published over the last 25 years. The aim of this review is two-fold. First we show technical details of the GF method by selecting a number of subjects which are of current interest in research: helimagnets, systems including a DM interaction, and surface effects in thin films. Second, we show that these systems possess many striking features due to the frustration.

This paper is organized as follows. In Section 2, we express the Hamiltonian in a general non-collinear GS and define the local system of spin coordinates. Here, we also present the calculation of the GS and the foundation of the self-consistent GF technique and the calculation of the SW dispersion relation and layer magnetizations at arbitrary temperature (T). We show in Section 3 the numerical results obtained from the GF. Section 4 shows interesting examples using various kinds of interaction including the DM interaction in a variety of systems from two dimensions, to thin films and superlattices. Section 5 treats a case where the DM interaction competes with the antiferromagnetic interaction in the frustrated antiferromagnetic triangular lattice. Section 6 presents the surface effect in a thin film where its surface is frustrated. Concluding remarks are given in Section 7.

2. Hamiltonian of a Chiral Magnet—Local Coordinates

Chiral order in helimagnets has been subject of recent extensive investigations. In Ref. [25], the surface structure of thin helimagnetic films has been studied. In Ref. [26] exotic spin configurations in ultrathin helimagnetic holmium films have been investigated. In Refs. [27,28] chiral structure and spin reorientations in MnSi thin films have been theoretically studied. In these works, the chiral structures have been considered at $T=0$, but not the SW even at $T=0$. The main difficulty was due to the non-collinear, non-uniform spin configurations. We have shown that this was possible using the GFs generalized for such spin configurations given in Ref. [7].

To demonstrate the method, let us follow Ref. [29]: we consider the body-centered tetragonal (bct) lattice with Heisenberg spins. Each spin interacts with its nearest neighbors (NN) via the exchange constant $J_1$ and with its next NN (NNN) on the c-direction via the exchange $J_2$ (see Figure 1).

We consider the simplest model of a helimagnet, given by the following Hamiltonian

$$\mathcal{H}=-J_1\sum _{i,j}{\mathbf{S}}_i·{\mathbf{S}}_j-J_2\sum _{i,k}{\mathbf{S}}_i·{\mathbf{S}}_k$$

(1)

where ${\mathbf{S}}_i$ is a quantum spin of magnitude 1/2, the first sum is performed over all NN pairs, and the second sum over pairs on the c-axis (cf. Figure 1).

In the case of an infinite crystal, the chiral state occurs when $J_1$ is ferromagnetic and $J_2$ is antiferromagnetic and $|J_2|/J_1$ is larger than a critical value, as will be shown below.

Let us suppose that the energy of a spin $E_C$ in a chiral configuration when the angle between two NN spin in the neighboring planes is $\theta$, one has (omitting the factor $S^2$)

$$E=-8J_{1}cos\theta -2J_{2}cos\left(2\theta \right)$$

(2)

The lowest-energy state corresponds to

$$\begin{array}{ccc}\hfill \frac{dE}{d\theta }& =& 0\hfill \\ \hfill \to 8J_{1}sin\theta +4sin\left(2\theta \right)& =& 0\hfill \\ \hfill 8J_{1}sin\theta (1+\frac{J_2}{J_1}cos\theta )& =& 0\hfill \end{array}$$

(3)

There are two solutions, $sin\theta =0$ and $cos\theta =-\frac{J_1}{J_2}$ The first solution corresponds to the ferromagnetic state, and the second solution exists if $-\frac{J_1}{J_2}\le 1$ which corresponds to the chiral state.

For a thin helimagnetic film, the angle between spins in adjacent layers varies due to the surface. We can use the method of energy minimization for each layer, then we have a set of coupled equations to solve (see Ref. [29]). Figure 2 displays an example of the angle distribution across the film thickness $N_z$.

In order to calculate the SW spectrum for systems of non-collinear spin configurations, let us emphasize that the commutation relations between spin operators are established when the spin lies on its quantization z. In the non-collinear cases, each spin has its own quantization axis. It is therefore important to choose a quantization axis for each spin. We have to use the system of local coordinates defined as follows. In the Hamiltonian, the spins are coupled two by two. Consider a pair ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$. As seen above, in the general case these spins make an angle ${\theta }_{i,j}={\theta }_j-{\theta }_i$ determined by the competing interactions in the systems. For quantum spins, in the course of calculation we need to use the commutation relations between the spin operators $S^z,S^{+},S^{-}$. As said above, these commutation relations are derived from the assumption that the spin lies on its quantization axis z. We show in Figure 3 the local coordinates assigned to spin ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$. We write

$$\begin{array}{ccc}\hfill {\mathbf{S}}_i& =& S_i^x{\widehat{\xi }}_i+S_i^y{\widehat{\eta }}_i+S_i^z{\widehat{\zeta }}_i\hfill \end{array}$$

(4)

$$\begin{array}{ccc}\hfill {\mathbf{S}}_j& =& S_j^x{\widehat{\xi }}_j+S_j^y{\widehat{\eta }}_j+S_j^z{\widehat{\zeta }}_j\hfill \end{array}$$

(5)

Expressing the axes of ${\mathbf{S}}_j$ in the frame of ${\mathbf{S}}_i$ one has

$$\begin{array}{ccc}\hfill {\widehat{\zeta }}_j& =& cos{\theta }_{i,j}{\widehat{\zeta }}_i+sin{\theta }_{i,j}{\widehat{\xi }}_i\hfill \end{array}$$

(6)

$$\begin{array}{ccc}\hfill {\widehat{\xi }}_j& =& -sin{\theta }_{i,j}{\widehat{\zeta }}_i+cos{\theta }_{i,j}{\widehat{\xi }}_i\hfill \end{array}$$

(7)

$$\begin{array}{ccc}\hfill {\widehat{\eta }}_j& =& {\widehat{\eta }}_i\hfill \end{array}$$

(8)

so that

$$\begin{array}{ccc}\hfill {\mathbf{S}}_j& =& S_j^x(-sin{\theta }_{i,j}{\widehat{\zeta }}_i+cos{\theta }_{i,j}{\widehat{\xi }}_i)\hfill \\ & & +S_j^y{\widehat{\eta }}_i+S_j^z(cos{\theta }_{i,j}{\widehat{\zeta }}_i+sin{\theta }_{i,j}{\widehat{\xi }}_i)\hfill \end{array}$$

(9)

Using Equation (9) to express ${\mathbf{S}}_j$ in the $({\widehat{\xi }}_i,{\widehat{\eta }}_i,{\widehat{\zeta }}_i)$ coordinates, we calculate ${\mathbf{S}}_i·{\mathbf{S}}_j$, we get the following Hamiltonian from (28):

$$\begin{array}{ccc}\hfill {\mathcal{H}}_e& =& -\sum _{<i,j>}{J}_{i,j}\{\frac{1}{4}\left(cos{\theta }_{i,j}-1\right)\left(S_i^{+}{S}_j^{+}+S_i^{-}{S}_j^{-}\right)\hfill \\ & +& \frac{1}{4}\left(cos{\theta }_{i,j}+1\right)\left(S_i^{+}{S}_j^{-}+S_i^{-}{S}_j^{+}\right)\hfill \\ & +& \frac{1}{2}sin{\theta }_{i,j}\left(S_i^{+}+S_i^{-}\right)S_j^z-\frac{1}{2}sin{\theta }_{i,j}{S}_i^z\left(S_j^{+}+S_j^{-}\right)\hfill \\ & +& cos{\theta }_{i,j}{S}_i^z{S}_j^z\}\hfill \end{array}$$

(10)

This explicit Hamiltonian in terms of the angle between two NN spins is common for a non-collinear spin configuration due to exchange interactions $J_{i,j}$. For other types of interactions such as the DM interaction, the explicit Hamiltonian in terms of the angle will be different as shown in Section 4.

We define the following GFs for the above Hamiltonian:

$$\begin{array}{ccc}\hfill G_{i,j}(t,t^{\prime })& =& <<S_i^{+}\left(t\right);S_j^{-}\left(t^{\prime }\right)>>\hfill \\ & =& -i\theta (t-t^{\prime })<\left[S_i^{+}\left(t\right),S_j^{-}\left(t^{\prime }\right)\right]>\hfill \end{array}$$

(11)

$$\begin{array}{ccc}\hfill F_{i,j}(t,t^{\prime })& =& <<S_i^{-}\left(t\right);S_j^{-}\left(t^{\prime }\right)>>\hfill \\ & =& -i\theta (t-t^{\prime })<\left[S_i^{-}\left(t\right),S_j^{-}\left(t^{\prime }\right)\right]>\hfill \end{array}$$

(12)

Writing their equations of motion we have

$$\begin{array}{ccc}\hfill i\hslash \frac{d}{dt}{G}_{i,j}\left(t,t^{\prime }\right)& =& 〈\left[S_i^{+}\left(t\right),S_j^{-}\left(t^{\prime }\right)\right]〉\delta \left(t-t^{\prime }\right)\hfill \\ & -& 〈〈\left[\mathcal{H},S_i^{+}\left(t\right)\right];S_j^{-}\left(t^{\prime }\right)〉〉,\hfill \end{array}$$

(13)

$$\begin{array}{ccc}\hfill i\hslash \frac{d}{dt}{F}_{i,j}\left(t,t^{\prime }\right)& =& 〈\left[S_i^{-}\left(t\right),S_j^{-}\left(t^{\prime }\right)\right]〉\delta \left(t-t^{\prime }\right)\hfill \\ & -& 〈〈\left[\mathcal{H},S_i^{-}\left(t\right)\right];S_j^{-}\left(t^{\prime }\right)〉〉,\hfill \end{array}$$

(14)

where

$$\begin{array}{ccc}\hfill S_j^{\pm }& =& S_j^x{\widehat{\xi }}_j\pm i{S}_j^y{\widehat{\eta }}_j\hfill \\ \hfill \left[S_j^{+},S_l^{-}\right]& =& 2S_j^z{\delta }_{j,l}\hfill \\ \hfill \left[S_j^z,S_l^{\pm }\right]& =& \pm S_j^{\pm }{\delta }_{j,l}\hfill \end{array}$$

Note that the equation of motion of the G Green’s function generates the F Green’s functions, and vice-versa. Performing the commutators in Equations (13) and (14), and using the Tyablikov approximation [30] for higher-order GFs, for instance $<<S_{i^{\prime }}^z{S}_i^{+}\left(t\right);S_j^{-}\left(t^{\prime }\right)>>\simeq <S_{i^{\prime }}^z><<S_i^{+}\left(t\right);S_j^{-}\left(t^{\prime }\right)>>$ etc., we obtain

$$\begin{array}{ccc}\hfill i\hslash \frac{d{G}_{i,j}(t,t^{\prime })}{dt}& =& 2<S_i^z>{\delta }_{i,j}\delta (t-t^{\prime })\hfill \\ & -& \sum _{i^{\prime }}{J}_{i,i^{\prime }}[<S_i^z>(cos{\theta }_{i,i^{\prime }}-1)\times \hfill \\ & \times & F_{i^{\prime },j}(t,t^{\prime })\hfill \\ & +& <S_i^z>(cos{\theta }_{i,i^{\prime }}+1)G_{i^{\prime },j}(t,t^{\prime })\hfill \\ & -& 2<S_{i^{\prime }}^z>cos{\theta }_{i,i^{\prime }}{G}_{i,j}(t,t^{\prime })]\hfill \\ & +& 2\sum _{i^{\prime }}{I}_{i,i^{\prime }}<S_{i^{\prime }}^z>cos{\theta }_{i,i^{\prime }}{G}_{i,j}(t,t^{\prime })\hfill \end{array}$$

(15)

$$\begin{array}{ccc}\hfill i\hslash \frac{d{F}_{i,j}(t,t^{\prime })}{dt}& =& \sum _{i^{\prime }}{J}_{i,i^{\prime }}[<S_i^z>(cos{\theta }_{i,i^{\prime }}-1)\times \hfill \\ & \times & G_{i^{\prime },j}(t,t^{\prime })\hfill \\ & +& <S_i^z>(cos{\theta }_{i,i^{\prime }}+1)F_{i^{\prime },j}(t,t^{\prime })\hfill \\ & -& 2<S_{i^{\prime }}^z>cos{\theta }_{i,i^{\prime }}{F}_{i,j}(t,t^{\prime })]\hfill \\ & -& 2\sum _{i^{\prime }}{I}_{i,i^{\prime }}<S_{i^{\prime }}^z>cos{\theta }_{i,i^{\prime }}{F}_{i,j}(t,t^{\prime })\hfill \end{array}$$

(16)

Note that the Tyablikov decoupling scheme is equivalent to the so-called “random-phase-approximation” (RPA).

For the sake of clarity, we write separately the NN and NNN sums, we have

$$\begin{array}{ccc}\hfill i\hslash \frac{d{G}_{i,j}(t,t^{\prime })}{dt}& =& 2<S_i^z>{\delta }_{i,j}\delta (t-t^{\prime })\hfill \\ & -& \sum _{k^{\prime }\in NN}{J}_{i,k^{\prime }}[<S_i^z>(cos{\theta }_{i,k^{\prime }}-1)\times \hfill \\ & \times & F_{k^{\prime },j}(t,t^{\prime })\hfill \\ & +& <S_i^z>(cos{\theta }_{i,k^{\prime }}+1)G_{k^{\prime },j}(t,t^{\prime })\hfill \\ & -& 2<S_{k^{\prime }}^z>cos{\theta }_{i,k^{\prime }}{G}_{i,j}(t,t^{\prime })]\hfill \\ & +& 2\sum _{k^{\prime }\in NN}{I}_{i,k^{\prime }}<S_{k^{\prime }}^z>cos{\theta }_{i,k^{\prime }}{G}_{i,j}(t,t^{\prime })\hfill \\ & -& \sum _{i^{\prime }\in NNN}{J}_{i,i^{\prime }}[<S_i^z>(cos{\theta }_{i,i^{\prime }}-1)\times \hfill \\ & \times & F_{i^{\prime },j}(t,t^{\prime })\hfill \\ & +& <S_i^z>(cos{\theta }_{i,i^{\prime }}+1)G_{i^{\prime },j}(t,t^{\prime })\hfill \\ & -& 2<S_{i^{\prime }}^z>cos{\theta }_{i,i^{\prime }}{G}_{i,j}(t,t^{\prime })]\hfill \end{array}$$

(17)

$$\begin{array}{ccc}\hfill i\hslash \frac{d{F}_{k,j}(t,t^{\prime })}{dt}& =& \sum _{i^{\prime }\in NN}{J}_{k,i^{\prime }}[<S_k^z>(cos{\theta }_{k,i^{\prime }}-1)\times \hfill \\ & \times & G_{i^{\prime },j}(t,t^{\prime })\hfill \\ & +& <S_k^z>(cos{\theta }_{k,i^{\prime }}+1)F_{i^{\prime },j}(t,t^{\prime })\hfill \\ & -& 2<S_{i^{\prime }}^z>cos{\theta }_{k,i^{\prime }}{F}_{k,j}(t,t^{\prime })]\hfill \\ & -& 2\sum _{i^{\prime }\in NN}{I}_{k,i^{\prime }}<S_{i^{\prime }}^z>cos{\theta }_{k,i^{\prime }}{F}_{k,j}(t,t^{\prime })\hfill \\ & +& \sum _{k^{\prime }\in NNN}{J}_{k,k^{\prime }}[<S_k^z>(cos{\theta }_{k,k^{\prime }}-1)\times \hfill \\ & \times & G_{k^{\prime },j}(t,t^{\prime })\hfill \\ & +& <S_k^z>(cos{\theta }_{k,k^{\prime }}+1)F_{k^{\prime },j}(t,t^{\prime })\hfill \\ & -& 2<S_{k^{\prime }}^z>cos{\theta }_{k,k^{\prime }}{F}_{k,j}(t,t^{\prime })]\hfill \end{array}$$

(18)

For simplicity, we suppose in the following $J_{k,k^{\prime }}$ are all equal to $J_1$ for NN interactions and to $J_2$ for NNN interactions. $I_{k,k^{\prime }}$ is taken to be $I_1$ for NN pairs. In addition, in the film coordinates defined above, we denote the Cartesian components of the spin position ${\mathbf{R}}_i$ by three indices $({\ell }_i,m_i,n_i)$ in three directions x, y and z.

Since there is the translation invariance in the $xy$ plane, the in-plane Fourier transforms of the above equations in the $xy$ plane are

$$\begin{array}{ccc}\hfill G_{i,j}\left(t,t^{\prime }\right)& =& \frac{1}{\Delta }\int {\int }_{BZ}d{\mathbf{k}}_{xy}\frac{1}{2\pi }{\int }_{-\infty }^{+\infty }d\omega e^{-i\omega \left(t-t^{\prime }\right)}\hfill \\ & & \times g_{n_i,n_j}\left(\omega ,{\mathbf{k}}_{xy}\right)e^{i{\mathbf{k}}_{xy}·\left({\mathbf{R}}_i-{\mathbf{R}}_j\right)},\hfill \end{array}$$

(19)

$$\begin{array}{ccc}\hfill F_{k,j}\left(t,t^{\prime }\right)& =& \frac{1}{\Delta }\int {\int }_{BZ}d{\mathbf{k}}_{xy}\frac{1}{2\pi }{\int }_{-\infty }^{+\infty }d\omega e^{-i\omega \left(t-t^{\prime }\right)}\hfill \\ & & \times f_{n_k,n_j}\left(\omega ,{\mathbf{k}}_{xy}\right)e^{i{\mathbf{k}}_{xy}·\left({\mathbf{R}}_k-{\mathbf{R}}_j\right)},\hfill \end{array}$$

(20)

where $\omega$ is the SW frequency, ${\mathbf{k}}_{xy}$ the wave-vector parallel to $xy$ planes and ${\mathbf{R}}_i$ the position of ${\mathbf{S}}_i$. $n_i$, $n_j$ and $n_k$ denote the z-components of the sites ${\mathbf{R}}_i$, ${\mathbf{R}}_j$ and ${\mathbf{R}}_k$. The integral over ${\mathbf{k}}_{xy}$ is performed in the first Brillouin zone ($BZ$) whose surface is $\Delta$ in the $xy$ reciprocal plane. $n_i=1$ denotes the surface layer, $n_i=2$ the second layer etc.

In the 3D case, the Fourier transformation of Equations (17) and (18) in the three $(x,y,z)$ directions yields the SW spectrum in the absence of anisotropy:

$$\hslash \omega =\pm \sqrt{A^2-B^2}$$

(21)

where

$$\begin{array}{ccc}\hfill A& =& J_1〈S^z〉[cos\theta +1]Z\gamma +2Z{J}_1〈S^z〉cos\theta \hfill \\ & & +J_2〈S^z〉[cos\left(2\theta \right)+1]Z_{c}cos\left(k_{z}a\right)\hfill \\ & & +2Z_c{J}_2〈S^z〉cos\left(2\theta \right)\hfill \\ \hfill B& =& J_1〈S^z〉(cos\theta -1)Z\gamma \hfill \\ & & +J_2〈S^z〉[cos\left(2\theta \right)-1]Z_{c}cos\left(k_{z}a\right)\hfill \end{array}$$

where $Z=8$ is the NN coordination number, $Z_c=2$ the NNN number on the c-axis and $\gamma =cos(k_{x}a/2)cos(k_{y}a/2)cos(k_{z}a/2)$ where a is the lattice constant taken the same in three directions. Note that $\hslash \omega$ is zero when $A=\pm B$. This is realized at two points as expected in helimagnets: $k_x=k_y=k_z=0$ ($\gamma =1$) and $k_z=2\theta$ along the helical axis. It is interesting to note that we recover the SW dispersion relation of ferromagnets (antiferromagnets) [2] with NN interaction only by putting $cos\theta =1$$(-1)$ in the above coefficients.

In the case of a thin film, the in-plane Fourier transformation yields the following matrix equation

$$\mathbf{M}\left(\omega \right)\mathbf{h}=\mathbf{u},$$

(22)

where $\mathbf{h}$ and $\mathbf{u}$ are given by

$$\mathbf{h}=\left(\begin{array}{c}{g}_{1,n^{\prime }}\\ f_{1,n^{\prime }}\\ ⋮\\ g_{n,n^{\prime }}\\ f_{n,n^{\prime }}\\ ⋮\\ g_{N_z,n^{\prime }}\\ f_{N_z,n^{\prime }}\end{array}\right),\mathbf{u}=\left(\begin{array}{c}2〈S_1^z〉{\delta }_{1,n^{\prime }}\\ 0\\ ⋮\\ 2〈S_{N_z}^z〉{\delta }_{N_z,n^{\prime }}\\ 0\end{array}\right),$$

(23)

We take $\hslash =1$ hereafter. Note that $\mathbf{M}\left(\omega \right)$ is a $\left(2N_z\times 2N_z\right)$ matrix given by Equation (24) where

$$\mathbf{M}\left(\omega \right)=\left(\begin{array}{cccccccccccc}\omega +A_1& 0& B_1^{+}& C_1^{+}& D_1^{+}& E_1^{+}& 0& 0& 0& 0& 0& 0\\ 0& \omega -A_1& -C_1^{+}& -B_1^{+}& -E_1^{+}& -D_1^{+}& 0& 0& 0& 0& 0& 0\\ \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots \\ \cdots & D_n^{-}& E_n^{-}& B_n^{-}& C_n^{-}& \omega +A_n& 0& B_n^{+}& C_n^{+}& D_n^{+}& E_n^{+}& \cdots \\ \cdots & -E_n^{-}& -D_n^{-}& -C_n^{-}& -B_n^{-}& 0& \omega -A_n& -C_n^{+}& -B_n^{+}& -E_n^{+}& -D_n^{+}& \cdots \\ \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots \\ 0& 0& 0& 0& 0& 0& D_{N_z}^{-}& E_{N_z}^{-}& B_{N_z}^{-}& C_{N_z}^{-}& \omega +A_{N_z}& 0\\ 0& 0& 0& 0& 0& 0& -E_{N_z}^{-}& -D_{N_z}^{-}& -C_{N_z}^{-}& -B_{N_z}^{-}& 0& \omega -A_{N_z}\end{array}\right)$$

(24)

$$\begin{array}{ccc}\hfill A_n& =& -8J_1(1+d)[〈S_{n+1}^z〉cos{\theta }_{n,n+1}\hfill \\ & & +〈S_{n-1}^z〉cos{\theta }_{n,n-1}]\hfill \\ & -& 2J_2[〈S_{n+2}^z〉cos{\theta }_{n,n+2}\hfill \\ & & +〈S_{n-2}^z〉cos{\theta }_{n,n-2}]\hfill \end{array}$$

$$\begin{array}{ccc}\hfill B_n^{\pm }& =& 4J_1〈S_n^z〉(cos{\theta }_{n,n\pm 1}+1)\gamma \hfill \\ \hfill C_n^{\pm }& =& 4J_1〈S_n^z〉(cos{\theta }_{n,n\pm 1}-1)\gamma \hfill \\ \hfill E_n^{\pm }& =& J_2〈S_n^z〉(cos{\theta }_{n,n\pm 2}-1)\hfill \\ \hfill D_n^{\pm }& =& J_2〈S_n^z〉(cos{\theta }_{n,n\pm 2}+1)\hfill \end{array}$$

where we recall that n denotes the layer number, namely $1,2,\cdots ,N_z$ and $d=I_1/J_1$. Note that ${\theta }_{n,n\pm 1}$ denotes the angle between a spin in the layer n and its NN spins in adjacent layers $n\pm 1$ etc., and $\gamma =cos\left(\frac{k_{x}a}{2}\right)cos\left(\frac{k_{y}a}{2}\right).$

In order to obtain the SW frequency $\omega$, we solve the secular equation $det\left|\mathbf{M}\right|=0$ for each given ($k_x,k_y)$. Since the linear dimension of the square matrix is 2$N_z$, we obtain 2$N_z$ eigen-values of $\omega$, half positive and half negative, corresponding to two opposite spin precessions as in antiferromagnets. These values depend on the input values $<S_n^z>$ ($n=1,\cdots ,N_z$). Thus, we have to solve the secular equation by iteration until the convergence of input and output values. Note that, even at $T=0$, $<S_n^z>$ are not equal to $1/2$ due to the zero-point spin contraction [31]. In addition, because of the film surfaces, the spin contractions are not uniform.

The solution for $g_{n,n}$ can be calculated (see Ref. [29]). The spectral theorem [1] can be used to obtain, after a somewhat lengthy algebra (see [29]),:

$$\langle S_n^z\rangle =\frac{1}{2}-\frac{1}{\Delta }\int \int d{k}_{x}d{k}_y\sum _{i=1}^{2N_z}\frac{D_{2n-1}\left({\omega }_i\right)}{e^{\beta {\omega }_i}-1}$$

(25)

where $n=1,\cdots ,N_z$, and

$$D_{2n-1}\left({\omega }_i\left({\mathbf{k}}_{xy}\right)\right)=\frac{{\left|\mathbf{M}\right|}_{2n-1}\left({\omega }_i\left({\mathbf{k}}_{xy}\right)\right)}{{\prod }_{j\ne i}\left[{\omega }_j\left({\mathbf{k}}_{xy}\right)-{\omega }_i\left({\mathbf{k}}_{xy}\right)\right]}.$$

(26)

As $<S_n^z>$ depend each other in ${\omega }_i(i=1,\cdots ,2N_z)$, their solutions should be obtained by iteration at a given temperature T. In the particular case where $T=0$ one has

$$\langle S_n^z\rangle (T=0)=\frac{1}{2}+\frac{1}{\Delta }\int \int d{k}_{x}d{k}_y\sum _{i=1}^{N_z}{D}_{2n-1}\left({\omega }_i\left({\mathbf{k}}_{xy}\right)\right)$$

(27)

Note that the sum is performed over $N_z$ negative ${\omega }_i$ since positive ${\omega }_i$ yields the zero Bose–Einstein factor at $T=0$).

The transition temperature $T_c$ can be calculated self-consistently when all $<S_n^z>$ tend to zero.

We show in the following section, the numerical results using the above formulas.

3. Results for Helimagnets Obtained from the Green’s Function Technique

We use the ferromagnetic interaction between NN as unit, namely $J_1=1$. Take the helimagnetic case where $J_2$ is negative with $|J_2|>J_1$. We have determined above the spin configuration across the film for several values of $p=J_2/J_1$. Replacing the angles ${\theta }_{n,n\pm 1}$ and ${\theta }_{n,n\pm 2}$ in the matrix elements of $\left|\mathbf{M}\right|$, then calculating ${\omega }_i(i=1,\cdots ,2N_z)$ for each ${\mathbf{k}}_{xy}$. For the iterative procedure, the reader is referred to Ref. [29]. The solution $\langle S_n^z\rangle (n=1,\cdots ,N_z)$ is obtained when the input and the output are equal with a desired precision P.

3.1. Spectrum

We calculate the SW spectrum as described above for each a given $J_2/J_1$. The SW spectrum depends on T. We show in Figure 4 the SW spectrum $\omega$ versus $k_x=k_y$ for an eight-layer film with $J_2/J_1=-1.4$ at $T=0.1$ and $T=1.02$ (in units of $J_1/k_B=1$). We observe that

(i)

There are opposite-precession SW modes. Unlike ferromagnets, SW in antiferromagnets and non-collinear spin structures have opposite spin precessions [31]. The negative sign does not mean SW negative energy, but it indicates just the precession contrary to the trigonometric sense,

(ii)

There are two degenerate acoustic “surface” branches one on each side. These degenerate “surface” modes stem from the symmetry of the two surfaces. These surface modes propagate parallel to the film surface but are damped when going to the bulk,

(iii)

With increasing T, layer magnetizations decrease as seen hereafter, this reduces therefore the SW frequency (see Figure 4b),

(iv)

Surface and bulk SW spectra have been observed by inelastic neutron scattering in collinear magnets (ferro- and antiferromagnetic films) [32,33]. However, such experiments have not been reported for helimagnetic thin films.

3.2. Zero-Point Spin Contraction and Transition Temperature

It is known that, in antiferromagnetic materials, quantum fluctuations cause a contraction of the spin length, namely the spin length is shorter than the spin magnitude, at $T=0$ [31]. We demonstrate here that a spin with a stronger antiferromagnetic interaction has a stronger contraction: spins in the first and in the second layers have only one antiferromagnetic NNN on the c-axis while interior spins have two NNN. The contraction at a given $J_2/J_1$ is thus expected to be stronger for interior spins. This is shown in Figure 5: with increasing $|J_2|/J_1$, i.e., the antiferromagnetic interaction becomes stronger, the contraction is stronger. Of course, there is no contraction when the system is ferromagnetic, namely when $J_2\to -1$.

3.3. Layer Magnetizations

We show now the layer ordering in Figure 6 and Figure 7 where $J_2/J_1=-1.4$ and −2, respectively, in the case of $N_z=8$. Consider first the case $J_2/J_1=-1.4$. We note that the surface magnetization, having a large value at $T=0$ as seen in Figure 5, crosses the interior layer magnetizations at $T\simeq 0.42$ to become much smaller than interior magnetizations at higher temperatures. This crossover phenomenon is due to the competition between quantum fluctuations, which dominate low-T behavior, and the low-lying surface SW modes which reduce the surface magnetization at higher T. Note that the second-layer magnetization makes also a crossover at $T\simeq 1.3$ which is more complicated to analyze. Similar crossovers have been observed in other quantum systems such as antiferromagnetic films [34] and superlattices [35].

Similar remarks are also hold for $J_2/J_1=-2$ shown in Figure 7.

Note that the results shown above have been calculated with an in-plane anisotropy interaction $d=0.1$. Larger d yields stronger layer magnetizations and larger $T_c$.

To close this section on SW in helimagnetic bct thin films, we mention that the same investigation was done in the case of simple-cubic helimagnetic films where the surface spin reconstruction and the surface SW have been shown [36]. We have also studied the frustrated bct Heisenberg helimagnet in which the SW spectrum of the non-collinear spin configuration has been calculated [37].

4. Dzyaloshinskii–Moriya Interaction in Thin Films

Let us consider a thin film made of N square lattices stacked in the y direction perpendicular to the film surface. The results for this system have been published in Ref. [38]. Hereafter, we review some of these important results. The Hamiltonian is given by

$$\begin{array}{ccc}\hfill \mathcal{H}& =& {\mathcal{H}}_e+{\mathcal{H}}_{DM}\hfill \end{array}$$

(28)

$$\begin{array}{ccc}\hfill {\mathcal{H}}_e& =& -\sum _{〈i,j〉}{J}_{i,j}{\mathbf{S}}_i·{\mathbf{S}}_j\hfill \end{array}$$

(29)

$$\begin{array}{ccc}\hfill {\mathcal{H}}_{DM}& =& \sum _{〈i,j〉}{\mathbf{D}}_{i,j}·{\mathbf{S}}_i\times {\mathbf{S}}_j\hfill \end{array}$$

(30)

where $J_{i,j}$ and ${\mathbf{D}}_{i,j}$ are the exchange and DM interactions, respectively, between two quantum Heisenberg spins ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$ of magnitude $S=1/2$.

We suppose in this section the in-plane and inter-plane exchange interactions between NN are both ferromagnetic and denoted by $J_1$ and $J_2$, respectively. The DM interaction is defined only between NN in the plane for simplicity. The J term favors the collinear spin configuration while the DM term favors the perpendicular one. This will lead to a compromise where ${\mathbf{S}}_i$ makes an angle ${\theta }_{i,j}$ with its neighbor ${\mathbf{S}}_j$. It is obvious that the quantization axes of ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$ are different. Therefore, the transformation using the local coordinates, Equations (4)–(9), is necessary. Let us suppose that the vector ${\mathbf{D}}_{i,j}$ is along the y axis, namely the ${\widehat{\eta }}_i$ axis. We write

$${\mathbf{D}}_{i,j}=D{e}_{i,j}{\widehat{y}}_i$$

(31)

where $e_{i,j}$ = +1 (−1) if $j>i$ ($j<i)$ for NN j on the $\widehat{x}$ or $\widehat{z}$ axis. One has by definition $e_{j,i}=-e_{i,j}$.

The easiest way to determine the GS is to minimize the local energy at each spin: taking a spin and calculating the local field acting on it from its neighbors. Then, we align the spin in its local-field direction to minimize its energy. Repeating this procedure for all spins, we say we realize one sweep. We have to make a sufficient number of sweeps to obtain the convergence with a desired precision (see details in Ref. [39]). This local energy minimization is called “the steepest descent method”. We show in Figure 8 the configuration obtained for $D=-0.5$ using $J_1=J_2=1$.

We see that each spin has the same angle with its four NN in the plane (angle between NN in adjacent planes is zero). We demonstrate now the dependence of $\theta$ on $J_1$: the energy of the spin ${\mathbf{S}}_i$ is written as

$$E_i=-4J_1{S}^{2}cos\theta -2J_2{S}^2+4D{S}^{2}sin\theta$$

(32)

where $\theta =|{\theta }_{i,j}|$ minimizing $E_i$ with respect to $\theta$ one obtains

$$\frac{d{E}_i}{d\theta }=0\Rightarrow -\frac{D}{J_1}=tan\theta \Rightarrow \theta =arctan(-\frac{D}{J_1})$$

(33)

The result is in agreement with that obtained by the steepest descent method. An example has been shown in Figure 8.

We rewrite the DM term of Equation (30) as

$$\begin{array}{ccc}\hfill {\mathbf{S}}_i\times {\mathbf{S}}_j& =& (-S_i^z{S}_j^y-S_i^y{S}_j^{x}sin{\theta }_{i,j}+S_i^y{S}_j^{z}cos{\theta }_{i,j}){\widehat{\xi }}_i\hfill \\ & & +(S_i^x{S}_j^{x}sin{\theta }_{i,j}+S_i^z{S}_j^{z}sin{\theta }_{i,j}){\widehat{\eta }}_i\hfill \\ & & +(S_i^x{S}_j^y-S_i^y{S}_j^{z}sin{\theta }_{i,j}-S_i^y{S}_j^{x}cos{\theta }_{i,j}){\widehat{\zeta }}_i\hfill \end{array}$$

(34)

From Equation (31), we obtain

$$\begin{array}{ccc}\hfill {\mathcal{H}}_{DM}& =& \sum _{〈i,j〉}{\mathbf{D}}_{i,j}·{\mathbf{S}}_i\times {\mathbf{S}}_j\hfill \\ & =& D\sum _{〈i,j〉}(S_i^x{S}_j^x{e}_{i,j}sin{\theta }_{i,j}+S_i^z{S}_j^z{e}_{i,j}sin{\theta }_{i,j})\hfill \\ & =& \frac{D}{4}\sum _{〈i,j〉}[(S_i^{+}+S_i^{-})(S_j^{+}+S_j^{-})e_{i,j}sin{\theta }_{i,j}\hfill \\ & & +4S_i^z{S}_j^z{e}_{i,j}sin{\theta }_{i,j}]\hfill \end{array}$$

(35)

where we have replaced $S^x$ by $(S^{+}+S^{-})/2$. Note that $e_{i,j}sin{\theta }_{i,j}$ is always positive since for a NN on the positive axis direction, $e_{i,j}=1$ and $sin{\theta }_{i,j}=sin\theta$ where $\theta$ is positively defined, while for a NN on the negative axis direction, $e_{i,j}=-1$ and $sin{\theta }_{i,j}=sin(-\theta )=-sin\theta$.

4.1. Formulation of the Green’s Function Technique for the Dzyaloshinskii–Moriya System

Using the transformation into the local coordinates, Equations (4)–(9), one has

$$\begin{array}{ccc}\hfill \mathcal{H}& =& -\sum _{<i,j>}{J}_{i,j}\{\frac{1}{4}\left(cos{\theta }_{i,j}-1\right)\left(S_i^{+}{S}_j^{+}+S_i^{-}{S}_j^{-}\right)\hfill \\ & +& \frac{1}{4}\left(cos{\theta }_{i,j}+1\right)\left(S_i^{+}{S}_j^{-}+S_i^{-}{S}_j^{+}\right)\hfill \\ & +& \frac{1}{2}sin{\theta }_{i,j}\left(S_i^{+}+S_i^{-}\right)S_j^z-\frac{1}{2}sin{\theta }_{i,j}{S}_i^z\left(S_j^{+}+S_j^{-}\right)\hfill \\ & +& cos{\theta }_{i,j}{S}_i^z{S}_j^z\}\hfill \\ & +& \frac{D}{4}\sum _{〈i,j〉}[(S_i^{+}+S_i^{-})(S_j^{+}+S_j^{-})e_{i,j}sin{\theta }_{i,j}\hfill \\ & & +4S_i^z{S}_j^z{e}_{i,j}sin{\theta }_{i,j}]\hfill \end{array}$$

(36)

Note that the quantization axes of the spins are in the $xz$ planes as shown in Figure 3.

We emphasize that, while the sine terms of the DM Hamiltonian, Equation (35), remain after summing over the NN, the sine terms of ${\mathcal{H}}_e$, the 3rd line of Equation (36), are zero after summing over opposite NN because there is no $e_{i,j}$ term.

It is very important to emphasize again that the commutation relations between spin operators $S^z$ and $S^{\pm }$ are valid when the spin lies on its local quantization axis. Therefore, it is necessary ro use the local coordinates for each spin.

In two dimensions (2D) there is no long-range order at non-zero T for isotropic spin models with short-range interaction [40]. Thin films have very small thickness, not far from 2D systems. Thus, in order to stabilize the ordering at very low T, we use a very small anisotropy interaction between between ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$ as follows

$${\mathcal{H}}_a=-\sum _{<i,j>}{I}_{i,j}{S}_i^z{S}_j^{z}cos{\theta }_{i,j}$$

(37)

where $I_{i,j}(>0)$ is positive, small compared to $J_1$, and limited to NN in the $xz$ plane. For simplicity, we suppose $I_{i,j}=I_1$ for all such NN pairs. As we will see below, the small value of $I_1$ does stabilize the SW spectrum when D becomes large. The Hamiltonian is finally given by

$$\mathcal{H}={\mathcal{H}}_e+{\mathcal{H}}_{DM}+{\mathcal{H}}_a$$

(38)

Using the two GF’s in the real space given by Equations (11) and (12) and using the same method, we study the effect of the DM interaction. For the DM term, the commutation relations $[\mathcal{H},S_i^{\pm }]$ lead to:

$$D\sum _{l}sin\theta [\mp S_i^z(S_l^{+}+S_l^{-})+\pm 2S_i^{\pm }{S}_l^z]$$

(39)

which gives rise, using the Tyablikov decoupling, to the following GF’s:

$$<<S_i^z{S}_l^{\pm };S_j^{-}>>\simeq <S_i^z><<S_l^{\pm };S_j^{-}>>$$

(40)

These functions are in fact the G and F functions. There are thus no new GF’s generated by the equations of motion.

As in Section 2, the Fourier transforms in the $xz$ plane $g_{n,n^{\prime }}$ and $f_{n,n^{\prime }}$ of the G and F lead to the matrix equation

$$\mathbf{M}\left(E\right)\mathbf{h}=\mathbf{u},$$

(41)

$\mathbf{M}\left(E\right)$ being given by Equation (42) below

$$\left(\begin{array}{ccccccccc}E+A_1& B_1& C_1& 0& 0& 0& 0& 0& 0\\ -B_1& E-A_1& 0& -C_1& 0& 0& 0& 0& 0\\ \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots \\ \cdots & 0& C_n& 0& E+A_n& B_n& C_n& 0& 0\\ \cdots & 0& 0& -C_n& -B_n& E-A_n& 0& -C_n& 0\\ \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots & \cdots \\ 0& 0& 0& 0& 0& C_N& 0& E+A_N& B_N\\ 0& 0& 0& 0& 0& 0& -C_N& -B_N& E-A_N\end{array}\right)$$

(42)

where $E=\hslash \omega$ is the SW energy and the matrix elements are given by

$$\begin{array}{ccc}\hfill A_n& =& -J_1[8<S_n^z>cos\theta (1+d_n)\hfill \\ & & -4<S_n^z>\gamma (cos\theta +1)]\hfill \\ & & -2J_2(<S_{n-1}^z>+<S_{n+1}^z>)\hfill \\ & & -8Dsin\theta <S_n^z>\gamma \hfill \\ & & +8Dsin\theta <S_n^z>\hfill \end{array}$$

(43)

$$\begin{array}{ccc}\hfill B_n& =& 4J_1<S_n^z>\gamma (cos\theta -1)\hfill \\ & & -8Dsin\theta <S_n^z>\gamma \hfill \end{array}$$

(44)

$$\begin{array}{ccc}\hfill C_n& =& 2J_2<S_n^z>\hfill \end{array}$$

(45)

where $n=1,2,\cdots ,N$ denotes the layer numbers, $d_n=I_1/J_1$, $\gamma =(cos{k}_{x}a+cos{k}_{z}a)/2$, $k_x$ and $k_z$ are the wave-vector components in the $xz$ planes, a being the lattice constant. Remarks: (i) if $n=1$ (surface layer) then there are no $n-1$ terms in the $A_n$, (ii) if $n=N$ then there are no $n+1$ terms in $A_n$.

For a thin film, the SW frequencies at a given wave vector $\overrightarrow{k}=(k_x,k_z)$ are obtained by diagonalizing (42).

The magnetization of the layer n at finite T is calculated as in the helimagnetic case shown in the previous section. The formula of the zero-point spin contraction is also presented there. The transition temperature $T_c$ can be also calculated by the same method. Let us show in the following the results.

4.2. Results for 2D and 3D Cases

In the 2D case, one has only one layer. The matrix (42) is

$$\begin{array}{ccc}\hfill (E+A_n)g_{n,n^{\prime }}+B_n{f}_{n,n^{\prime }}& =& 2<S_n^z>\delta (n,n^{\prime })\hfill \\ \hfill -B_n{g}_{n,n^{\prime }}+(E-A_n)f_{n,n^{\prime }}& =& 0\hfill \end{array}$$

(46)

where $A_n$ is given by (43) but without $J_2$ term for the 2D case. Coefficient $B_n$ is given by (44) and $C_n=0$. The SW frequencies are determined by the following secular equation

$$\begin{array}{ccc}& & (E+A_n)(E-A_n)+B_n^2=0\hfill \\ & \to & E^2-A_n^2+B_n^2=0\hfill \\ & \to & E=\pm \sqrt{(A_n+B_n)(A_n-B_n)}\hfill \end{array}$$

(47)

Several remarks are in order:

(i)

when $\theta =0$, the last three terms of $A_n$ and $B_n$ are zero: one recovers the ferromagnetic SW dispersion relation

$$E=2Z{J}_1<S_n^z>(1-\gamma )$$

(48)

where $Z=4$ is the coordination number of the square lattice (taking $d_n=0$),

(ii)

when $\theta =\pi$, one has $A_n=8J_1<S_n^z>$, $B_n=-8J_1<S_n^z>\gamma$. One recovers then the antiferromagnetic SW dispersion relation

$$E=2Z{J}_1<S_n^z>\sqrt{1-{\gamma }^2}$$

(49)

(iii)

when there is a DM interaction, one has $0<cos\theta <1$ ($0<\theta <\pi /2$). If $d_n=0$, the quantity in the square root of Equation (47) becomes negative at $\gamma =1$ when $\theta$ is not zero. The SW spectrum is not stable at $k_x=k_y=0$ because the energy is not real. The anisotropy $d_n$ can remove this instability if it is larger than a threshold value $d_c$. We solve the equation $(A_n+B_n)(A_n-B_n)=0$ to find $d_c$. In Figure 9 we show $d_c$ versus $\theta$. As seen, $d_c$ increases from zero with increasing $\theta$.

As we have anticipated, we need to include an anisotropy in order to allow for SW to be excited even at $T=0$ and for a long-range ordering at non-zero T in 2D as seen below.

We show in Figure 10 the SW dispersion relation calculated from Equation (47) for $\theta =0.2$ and 0.6 (radian). As seen, the spectrum is symmetric for positive and negative wave vectors. It is also symmetric for left and right precessions. One observes that for small $\theta$, namely small D, $E\left(k\right)$ is proportional to $k^2$ at low k (see Figure 10a). This behavior is that in ferromagnets. For large $\theta$, one observes that $E\left(k\right)$ becomes linear in k as seen in Figure 10b. This behavior is similar to that of antiferromagnets. Note that the change of behavior is progressive with increasing $\theta$, we do not observe a sudden transition from $k^2$ to k behavior. This behavior is also observed in 3D and in thin films as well.

As said earlier, the inclusion of an anisotropy d permits a long-range ordering at $T\ne 0$ in 2D: Figure 11 displays the magnetization M ($\equiv <S^z>$) calculated by Equation (2) where in each case the limit value $d_c$ has been used. We note that M depends strongly on $\theta$: at high T the larger $\theta$ the stronger M. However, at $T=0$ the spin length is smaller for larger $\theta$ due to the zero-point spin contraction [31] calculated by Equation (27). As a consequence there is a cross-over of layer magnetizations at low T as shown in Figure 11b. The spin length at $T=0$ is shown in Figure 12 for several $\theta$.

We now consider the 3D case. The crystal is infinite in three direction. The Fourier transform in the y direction, namely $g_{n\pm 1}=g_n{e}^{\pm i{k}_{y}a}$ and $f_{n\pm 1}=f_n{e}^{\pm i{k}_{y}a}$ reduces the matrix (23) to two coupled equations of g and f functions. One has

$$\begin{array}{ccc}\hfill (E+A^{\prime })g+Bf& =& 2<S^z>\hfill \\ \hfill -Bg+(E-A^{\prime })f& =& 0\hfill \end{array}$$

(50)

where

$$\begin{array}{ccc}\hfill A^{\prime }& =& -J_1[8<S^z>cos\theta (1+d)\hfill \\ & & -4<S^z>\gamma (cos\theta +1)]\hfill \\ & & -4J_2<S^z>\hfill \\ & & +4J_2<S^z>cos\left(k_{y}a\right)\hfill \\ & & -8Dsin\theta <S^z>\gamma \hfill \\ & & +8Dsin\theta <S^z>\hfill \end{array}$$

(51)

$$\begin{array}{ccc}\hfill B& =& 4J_1<S^z>\gamma (cos\theta -1)\hfill \\ & & -8Dsin\theta <S^z>\gamma \hfill \end{array}$$

(52)

The spectrum is given by

$$E=\pm \sqrt{(A^{\prime }+B)(A^{\prime }-B)}$$

(53)

In the ferromagnetic case, $cos\theta =1$, thus $B=0$. Arranging the Fourier transforms in three directions, one gets the 3D ferromagnetic dispersion relation $E=2Z<S^z>(1-{\gamma }^2)$ where $\gamma =[cos\left(k_{x}a\right)+cos\left(k_{y}a\right)+cos\left(k_{z}a\right)]/3$ and $Z=6$, coordination number of the simple cubic lattice.

As in the 2D case, we find a threshold value $d_c$ for which is the same for a given $\theta$. This is rather obvious because the DM interaction operates in the plane making an angle $\theta$ between spins in the plane, therefore its effects act on SW in each plane, not in the y direction perpendicular to the “DM planes”. Using Equation (53), we calculate the 3D spectrum displayed in Figure 13 for a small and a large value of $\theta$. As in the 2D case, we observe $E\propto k$ when $k\to 0$ for large $\theta$. The main properties of the system are thus governed by the in-plane DM interaction.

Figure 14 displays the magnetization M versus T for several values of $\theta$. As in the 2D case, when the DM interaction is included, the spins undergo a zero-point contraction which increases with increasing $\theta$. The competition between quantum fluctuations at $T=0$ and thermal effects at high T gives rise to magnetization cross-over shown in Figure 14b. The spin length at $T=0$ vs. $\theta$ is shown in the inset of Figure 14b. Comparing these results to those of the 2D case, we see that the spin contraction in 2D is stronger than in 3D. This is physically expected because quantum fluctuations are stronger at lower dimensions.

4.3. The Case of a Thin Film

As in the 2D and 3D cases, in the case of a thin film it is necessary to use a value for $d_n$ larger or equal to $d_c$ given in Figure 9 to stabilize the SW at long wave-length. Note that for thin films with more than one layer, the value of $d_c$ calculated for the 2D case remains valid.

Figure 15 displays the SW spectrum of a film of eight layers with $J_1=J_2=1$ for a small and a large $\theta$. As in the previous cases, E is proportional to k for large $\theta$ (cf. Figure 15b) but only for the first mode. The higher modes are proportional to $k^2$.

Figure 16 shows the layer magnetizations of the first four layers in a 8-layer film (the other half is symmetric) for two values of $\theta$. One observes that the surface magnetization is smaller than the magnetizations of other interior layers. This is due to the lack of neighbors for surface spins [2].

The spin contraction at $T=0$ is displayed Figure 16c.

The effects of the surface exchange and the film thickness have been shown in Ref. [38].

To close this section, let us mention our work [41] on the DM interaction in magneto-ferroelectric superlattices where the SW in the magnetic layer have been calculated. We have also studied the stability of skyrmions at finite T in that work and in Refs. [42,43].

5. Effect of Dzyaloshinskii–Moriya Interaction in a Frustrated Antiferromagnetic Triangular Lattice

The results of this section are not yet published [44]. We will not present this model in details. We show the Hamiltonian, the GS and the SW spectrum.

5.1. Model—Ground State

We consider a triangular lattice occupied by Heisenberg spins of magnitude 1/2. The DM interaction was introduced historically to explain the weak ferromagnetism in compounds MnO. The superexchange between two Mn atoms is modified with the displacement of the oxygen atom between them. If the displacement of the oxygen is in the $xy$ plane (see Figure 17a), then the DM vector ${\mathbf{D}}_{i,j}$ is perpendicular to the $xy$ plane and is given by [45,46]

$${\mathbf{D}}_{i,j}\propto {\mathbf{r}}_{iO}\times {\mathbf{r}}_{Oj}\propto -{\mathbf{r}}_{ij}\times \mathbf{R}$$

(54)

where ${\mathbf{r}}_{iO}={\mathbf{r}}_O-{\mathbf{r}}_i$ and ${\mathbf{r}}_{Oj}={\mathbf{r}}_j-{\mathbf{r}}_O$, ${\mathbf{r}}_{ij}={\mathbf{r}}_j-{\mathbf{r}}_i$. ${\mathbf{r}}_O$ is the position of non-magnetic ion (oxygen) and ${\mathbf{r}}_i$ the position of the spin ${\mathbf{S}}_i$ etc. These vectors are defined in Figure 17a in the particular case where the displacements are in the $xy$ plane. We have therefore ${\mathbf{D}}_{i,j}$ perpendicular to the $xy$ plane in this case.

Note, however, that if the atom displacements are in 3D space, ${\mathbf{D}}_{i,j}$ can be in any direction. In this paper, we consider also the case where ${\mathbf{D}}_{i,j}$ lies in the $xy$ plane as shown in Figure 17b where ${\mathbf{D}}_{i,j}$ is taken along the vector connecting spin ${\mathbf{S}}_i$ to spin ${\mathbf{S}}_j$.

Note that from Equation (54) one has

$${\mathbf{D}}_{j,i}=-{\mathbf{D}}_{i,j}$$

(55)

In the case of perpendicular ${\mathbf{D}}_{i,j}$, let us define ${\mathbf{u}}_{i,j}$ as the unit vector on the z axis. From Equations (54) and (55) one writes

$$\begin{array}{ccc}\hfill {\mathbf{D}}_{i,j}& =& D{\mathbf{u}}_{i,j}\hfill \end{array}$$

(56)

$$\begin{array}{ccc}\hfill {\mathbf{D}}_{j,i}& =& D{\mathbf{u}}_{j,i}=-D{\mathbf{u}}_{i,j}\hfill \end{array}$$

(57)

where D represents the DM interaction strength. Note however that the DM interaction goes beyond the weak ferromagnetism and may find its origin in various physical mechanisms. So, the form given in (56) is a model, a hypothesis.

In the case of in-plane ${\mathbf{D}}_{i,j}$, we suppose that ${\mathbf{D}}_{i,j}$ is given as

$${\mathbf{D}}_{i,j}=D({\mathbf{r}}_j-{\mathbf{r}}_i)/|{\mathbf{r}}_j-{\mathbf{r}}_i|=D{\mathbf{r}}_{ij}$$

(58)

where D is a constant and ${\mathbf{r}}_{ij}$ denotes the unit vector along ${\mathbf{r}}_j-{\mathbf{r}}_i$. The case of in-plane ${\mathbf{D}}_{i,j}$ on the frustrated triangular lattice (see Figure 17b) has been recently studied since this case gives rise to a beautiful skyrmion crystal composed of three interpenetrating sublattice skyrmions in a perpendicular applied magnetic field [44,47,48]. A description of this case is however out of the purpose of this review.

5.2. Ground State with a Perpendicular $\mathbf{D}$ in Zero Field

The Hamiltonian is given by

$$\begin{array}{ccc}\hfill \mathcal{H}& =& -J\sum _{\langle ij\rangle }{\mathbf{S}}_i·{\mathbf{S}}_j-D\sum _{\langle ij\rangle }{\mathbf{u}}_{i,j}·{\mathbf{S}}_i\times {\mathbf{S}}_j\hfill \\ & & -H\sum _i{S}_i^z\hfill \end{array}$$

(59)

where ${\mathbf{S}}_i$ is a classical Heisenberg spin of magnitude 1 occupying the lattice site i. The first sum runs over all spin nearest-neighbor (NN) pairs with an antiferromagnetic exchange interaction J ($J<0$), while the second sum is performed over all DM interactions between NN. H is the magnitude of a magnetic field applied along the z direction perpendicular to the lattice $xy$ plane.

In the absence of J, unlike the bipartite square lattice where one can arrange the NN spins to be perpendicular with each order in the $xy$ plane, the triangular lattice cannot fully satisfy the DM interaction for each bond, namely with the perpendicular spins at the ends. For this particular case of interest, we can analytically calculate the GS spin configuration as shown in the following. One considers a triangular plaquette with three spins numbered as 1, 2, and 3 embedded in the lattice. For convenience, in a hexagonal (or triangular) lattice, we define the three sublattices as follows: consider the up-pointing triangles (there are three in a hexagon, see the blue triangles in Figure 18). For the first triangle one numbers in the counter-clockwise sense 1, 2, 3 then one does it for the other two up-pointing triangles of the hexagon, one sees that each lattice site belongs to a sublattice. The DM energy of a plaquette is written as

$$\begin{array}{ccc}\hfill H_p& =& -2D[{\mathbf{u}}_{1,2}·{\mathbf{S}}_1\times {\mathbf{S}}_2+{\mathbf{u}}_{2,3}·{\mathbf{S}}_2\times {\mathbf{S}}_3+{\mathbf{u}}_{3,1}·{\mathbf{S}}_3\times {\mathbf{S}}_1]\hfill \\ & =& -2D[sin{\theta }_{1,2}+sin{\theta }_{2,3}+sin{\theta }_{3,1}]\hfill \end{array}$$

(60)

where the factor 2 of the D term takes into account the opposite neighbors outside the plaquette, and where ${\theta }_{1,2}={\theta }_2-{\theta }_1$ is the oriented angle between ${\mathbf{S}}_1$ and ${\mathbf{S}}_2$, etc. Note that the u vectors are in the same direction because we follow the counter-clockwise tour on the plaquette.

The minimization of $H_p$ yields

$$\begin{array}{ccc}\hfill \frac{d{H}_p}{d{\theta }_1}& =& 0=-2D[-cos({\theta }_2-{\theta }_1)+cos({\theta }_1-{\theta }_3)]\hfill \end{array}$$

(61)

$$\begin{array}{ccc}\hfill \frac{d{H}_p}{d{\theta }_2}& =& 0=-2D[cos({\theta }_2-{\theta }_1)-cos({\theta }_3-{\theta }_2)]\hfill \end{array}$$

(62)

$$\begin{array}{ccc}\hfill \frac{d{H}_p}{d{\theta }_3}& =& 0=-2D[cos({\theta }_3-{\theta }_2)-cos({\theta }_1-{\theta }_3)]\hfill \end{array}$$

(63)

The solutions for the above equations are

$$\begin{array}{ccc}\hfill {\theta }_{1,2}& =& {\theta }_{3,1}\mathrm{so} \mathrm{that}{\theta }_{3,2}={\theta }_{3,1}+{\theta }_{1,2}=2{\theta }_{1,2}\hfill \end{array}$$

(64)

$$\begin{array}{ccc}\hfill {\theta }_{2,3}& =& {\theta }_{1,2}\mathrm{so} \mathrm{that}{\theta }_{1,3}={\theta }_{1,2}+{\theta }_{2,3}=2{\theta }_{2,3}\hfill \end{array}$$

(65)

$$\begin{array}{ccc}\hfill {\theta }_{3,1}& =& {\theta }_{2,3}\mathrm{so} \mathrm{that}{\theta }_{2,1}={\theta }_{2,3}+{\theta }_{3,1}=2{\theta }_{3,1}\hfill \end{array}$$

(66)

Note that the second and third lines can be obtained by the circular permutation of the indices 1, 2, and 3 using the first line. These three equations, Equations (64)–(66), should be solved. There is more than one solution. We have from Equation (61) $cos\left({\theta }_{1,2}\right)=cos\left({\theta }_{3,1}\right)$. Using Equation (66) one obtains

$$cos\left(2{\theta }_{3,1}\right)=cos\left({\theta }_{3,1}\right)\to 2{cos}^2\left({\theta }_{3,1}\right)-cos\left({\theta }_{3,1}\right)-1=0$$

(67)

This second-degree equation gives $cos\left({\theta }_{3,1}\right)=\frac{1\pm \sqrt{1+8}}{4}$. Only the negative solution is acceptable so that ${\theta }_{3,1}={\theta }_{2,3}=\pi /6$. From Equation (66), one has ${\theta }_{2,1}=\pi /3$. This is one solution given by Equation (68) below. Note that we have taken one of them, Equation (66), to obtain explicit solutions for the three angles given in Equation (68). We can do the same calculation starting with Equations (64) and (65) to get explicit solutions given in Equations (69) and (70). We note that when we make a circular permutation of the indices of Equation (68) we get Equation (69), and a circular permutation of Equation (69) gives Equation (70). One summarizes the three degenerate solutions below

$$\begin{array}{c}\hfill {\theta }_{3,1}={\theta }_{2,3}=\pi /6,{\theta }_{2,1}=\pi /3\end{array}$$

(68)

$$\begin{array}{c}\hfill {\theta }_{1,2}={\theta }_{3,1}=\pi /6,{\theta }_{3,2}=\pi /3\end{array}$$

(69)

$$\begin{array}{c}\hfill {\theta }_{2,3}={\theta }_{1,2}=\pi /6,{\theta }_{1,3}=\pi /3\end{array}$$

(70)

We show in Figure 18 the spin orientations of the solution (68). The GS energy is obtained by replacing the angles into Equation (60). For the three solutions, one gets the energy of the plaquette

$$H_p=-3D\sqrt{3}$$

(71)

We have three degenerate GSs.

Note that this solution can be numerically obtained by the steepest descent method described above. The result is shown in Figure 19 for the full lattice. We see in the zoom that the spin configuration on a plaquette is what is obtained analytically, with a global spin rotation as explained in the caption of Figure 18.

As said above, to use the steepest descent method, we consider a triangular lattice of lateral dimension L. The total number of sites N is given by $N=L\times L$. To avoid the finite size effect, we have to find the size limit beyond which the GS does not depend on the lattice size. This is found for $L\ge 100$. Most of calculations have been performed for $L=100$.

5.3. Ground State with Both Perpendicular $\mathbf{D}$ and J in Zero Field-Spin Waves

When both J and perpendicular $\mathbf{D}$ are present, a compromise is established between these competing interactions. In zero field, the GS shows non-collinear but periodic in-plane spin configurations. The planar spin configuration is easily understood: when $\mathbf{D}$ is perpendicular and without J, the spins are in the plane. When J is antiferromagnetic without $\mathbf{D}$, the spins are also in the plane and form a 120-degree structure. When $\mathbf{D}$ and J exist together the angles between NN’s change but they remain in the plane in order to keep both D and J interactions as low as possible. An example is shown in Figure 20 where one sees that the GS is planar and characterized by two angles $\theta =102$ degrees and one angle $\beta =156$ degrees formed by three spins on a triangle plaquette. Note that there are three degenerate states where $\beta$ is chosen for the pair (1,2) (Figure 20a) or the pair (2,3) or the pair (3,1). Changing the value of D will change the angle values. Changing the sign of D results in a change of the sense of the chirality, but not the angle values.

In the case of perpendicular ${\mathbf{D}}_{i,j}$ in zero-field, as shown above we find the GS on a hexagon of the lattice is defined by four identical angles $\beta$ and two angles $\theta$ as shown in Figure 20. The values of $\beta$ and $\theta$ depend on the value of D. We take $J=-1$ (antiferromagnetic) hereafter. For $D=0.5$ we have $\beta =156$ degrees and $\theta =102$ degrees. For $D=0.4$ we obtain $\beta =108$ degrees and $\theta =144$ degrees, using $N=60\times 60$.

The periodicity of the GS allows us to calculate the SW spectrum in the following.

The model for the calculation of the SW spectrum uses quantum Heisenberg spins of magnitude $1/2$, it is given by

$$\begin{array}{c}\hfill \mathcal{H}=-J\sum _{〈i,j〉}{\mathbf{S}}_i·{\mathbf{S}}_j-D\sum _{〈i,j〉}{\mathbf{u}}_{i,j}·{\mathbf{S}}_i\times {\mathbf{S}}_j-I\sum _{〈i,j〉}{S}_i^z{S}_j^{z}cos{\theta }_{ij}\end{array}$$

(72)

where ${\theta }_{ij}$ is the angle between ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$ and the last term is an extremely small anisotropy added to stabilize the SW when the wavelength k tends to zero [31,40]. Note that ${\mathbf{u}}_{i,j}$ points up and down along the z axis for respective two opposite neighbors.

As before, in order to calculate the SW spectrum for systems of non-collinear spin configurations, we have to use the system of local coordinates. The Hamiltonian becomes

$$\begin{array}{cc}\hfill \mathcal{H}=& -J\sum _{〈i,j〉}\frac{1}{4}(S_i^{+}{S}_j^{+}+S_i^{-}{S}_j^{-})(cos{\theta }_{ij}-1)+\frac{1}{4}(S_i^{+}{S}_j^{-}+S_i^{-}{S}_j^{+})(cos{\theta }_{ij}+1)\hfill \\ \hfill & +\frac{1}{2}{S}_j^{z}sin{\theta }_{ij}(S_i^{+}+S_i^{-})-\frac{1}{2}sin{\theta }_{ij}{S}_i^z(S_J^{+}+S_j^{-})+S_i^z{S}_j^{z}cos{\theta }_{ij}\hfill \\ \hfill & -D\sum _{〈i,j〉}{S}_i^z{S}_j^{z}sin{\theta }_{i,j}+\frac{1}{4}sin{\theta }_{i,j}(S_i^{+}{S}_j^{+}+S_i^{+}{S}_j^{-}+S_i^{-}{S}_j^{+})+\frac{1}{2}cos{\theta }_{i,j}(S_i^z(S_j^{+}+S_j^{-})-S_j^z(S_i^{+}+S_i^{-}))\hfill \\ \hfill & -I\sum _{〈i,j〉}{S}_i^z{S}_j^{z}cos{\theta }_{i,j}\hfill \end{array}$$

We define the two GFs by Equations (11) and (12) and use the equations of motion of these functions (13) and (14), we obtain

$$\begin{array}{cc}\hfill i\hslash \frac{d{G}_{i,j}(t-t^{\prime })}{dt}& =2<S_i^z>{\delta }_{i,j}\delta (t-t^{\prime })-J\sum _{〈l〉}<S_i^z>F_{l,j}(t-t^{\prime })(cos{\theta }_{i,l}-1)\hfill \\ \hfill & +<S_i^z>G_{l,j}(t-t^{\prime })(cos{\theta }_{i,l}+1)-2cos{\theta }_{i,l}<S_l^z>G_{i,j}(t-t^{\prime })\hfill \\ \hfill & +D\sum _{〈l〉}2sin{\theta }_{i,l}<S_i^z>F_{l,j}(t-t^{\prime })-sin{\theta }_{i,l}<S_i^z>(G_{l,j}(t-t^{\prime })+F_{l,j}(t-t^{\prime }))\hfill \\ \hfill & -2I\sum _{〈l〉}cos{\theta }_{i,l}<S_i^z>F_{l,j}(t-t^{\prime })\hfill \end{array}$$

$$\begin{array}{cc}\hfill i\hslash \frac{d{F}_{i,j}(t-t^{\prime })}{dt}& =J\sum _{〈l〉}<S_i^z>G_{l,j}(t-t^{\prime })(cos{\theta }_{i,l}-1)\hfill \\ \hfill & +<S_i^z>F_{l,j}(t-t^{\prime })(cos{\theta }_{i,l}+1)-2cos{\theta }_{i,l}<S_l^z>F_{i,j}(t-t^{\prime })\hfill \\ \hfill & -D\sum _{〈l〉}2sin{\theta }_{i,l}<S_i^z>G_{l,j}(t-t^{\prime })-sin{\theta }_{i,l}<S_i^z>(G_{l,j}(t-t^{\prime })+F_{l,j}(t-t^{\prime }))\hfill \\ \hfill & +2I\sum _{〈l〉}cos{\theta }_{i,l}<S_i^z>G_{l,j}(t-t^{\prime })\hfill \end{array}$$

Note that $<S_i^z>$ is the average of the spin i on its local quantization axis in the local-coordinates system (see Ref. [38]). We use now the time Fourier transforms of the G and F, we get

$$\begin{array}{cc}\hfill \hslash \omega g_{i,j}& =2{\mu }_i{\delta }_{i,j}-J\sum _{〈l〉}{\mu }_i{f}_{lj}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}(cos{\theta }_{i,l}-1)\hfill \\ \hfill & +{\mu }_i{g}_{lj}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}(cos{\theta }_{i,l}+1)-2{\mu }_{l}cos{\theta }_{i,l}{g}_{i,j}\hfill \\ \hfill & -D\sum _{〈l〉}2sin{\theta }_{i,l}{\mu }_l{g}_{i,j}-sin{\theta }_{i,l}{\mu }_i(g_{l,j}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}+f_{l,j}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)})\hfill \\ \hfill & +2I\sum _{〈l〉}{\mu }_{l}cos{\theta }_{i,l}{g}_{i,j}\hfill \end{array}$$

(73)

and

$$\begin{array}{cc}\hfill \hslash \omega f_{i,j}& =J\sum _{〈l〉}{\mu }_i{g}_{lj}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}(cos{\theta }_{i,l}-1)\hfill \\ \hfill & +{\mu }_i{f}_{lj}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}(cos{\theta }_{i,l}+1)-2{\mu }_{l}cos{\theta }_{i,l}{f}_{i,j}\hfill \\ \hfill & +D\sum _{〈l〉}2sin{\theta }_{i,l}{\mu }_l{f}_{i,j}-sin{\theta }_{i,l}{\mu }_i(g_{l,j}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}+f_{l,j}{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)})\hfill \\ \hfill & -2I\sum _{〈l〉}{\mu }_{l}cos{\theta }_{i,l}{f}_{i,j}\hfill \end{array}$$

(74)

where ${\mu }_i\equiv <S_i^z>$, $\mathbf{k}$ is the wave vector in the reciprocal lattice of the triangular lattice, and $\omega$ the SW frequency. Note that the index z in $S_i^z$ is not referring to the real space direction z, but to the quantization axis of the spin ${\mathbf{S}}_i$. At this stage, we have to replace ${\theta }_{i,j}$ by either $\beta$ or $\theta$ according on the GS spin configuration given above (see Figure 20).

As in the previous sections, writing the above equations under a matrix form, we have

$$\mathbf{M}\left(\hslash \omega \right)\mathbf{h}=\mathbf{C},$$

(75)

where $\mathbf{M}\left(\hslash \omega \right)$ is a square matrix of dimension $2\times 2$, $\mathbf{h}$ and $\mathbf{C}$ are given by

$$\mathbf{h}=\left(\begin{array}{c}{g}_{i,j}\\ f_{i,j}\end{array}\right),\mathbf{C}=\left(\begin{array}{c}2〈S_i^z〉{\delta }_{i,j}\\ 0\end{array}\right),$$

(76)

and the matrix $\mathbf{M}\left(\hslash \omega \right)$ is given by

$$\mathbf{M}\left(\hslash \omega \right)=\left(\begin{array}{cc}\hslash \omega +A& B\\ -B& \hslash \omega -A\end{array}\right)$$

The nontrivial solution of g and f imposes the following secular equation:

$$0=\left(\begin{array}{cc}\hslash \omega +A& B\\ -B& \hslash \omega -A\end{array}\right)$$

(77)

where

$$\begin{array}{cc}\hfill A& =-J(8{\mu }_{i}cos\beta (1+I)+4{\mu }_{i}cos\theta (1+I)-4{\mu }_i\gamma (cos\beta +1)-2{\mu }_i\alpha (cos\theta +1))\hfill \\ \hfill & -D(4{\mu }_{i}sin\beta \gamma +2{\mu }_{i}sin\theta \alpha )+D(8{\mu }_{i}sin\beta +4{\mu }_{i}sin\theta )\hfill \end{array}$$

(78)

$$\begin{array}{cc}\hfill B& =J(4{\mu }_i\gamma (cos\beta -1)+2{\mu }_i\alpha (cos\theta -1))-D(4\gamma {\mu }_{i}sin\beta +2{\mu }_i\alpha sin\theta )\hfill \end{array}$$

(79)

where the sum on the two NN on the x axis (see Figure 20b) is

$$\sum _l{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}=2cos\left(k_x\right)\equiv 2\alpha$$

(80)

and the sum on the four NN on the oblique directions of the hexagon (see Figure 20b) is

$$\sum _l{e}^{-i\mathbf{k}·({\mathbf{R}}_i-{\mathbf{R}}_l)}=4cos(k_x/2)cos(\sqrt{3}{k}_y/2)\equiv 4\gamma$$

(81)

Solving Equation (77) for each given ($k_x,k_y$) one obtains the SW frequency $\omega (k_x,k_y)$:

$${(\hslash \omega )}^2=A^2-B^2\to \hslash \omega =\pm \sqrt{A^2-B^2}$$

(82)

Plotting $\omega (k_x,k_y)$ in the space $(k_x,k_y)$ one obtains the full SW spectrum.

The spin length $\langle S_i^z\rangle$ (for all i, by symmetry) is given by (see technical details in Ref. [31]):

$$\langle S^z\rangle \equiv \langle S_i^z\rangle =\frac{1}{2}-\frac{1}{\Delta }\int \int d{k}_{x}d{k}_z\sum _{i=1}^2\frac{Q\left(E_i\right)}{e^{E_i/k_{B}T}-1}$$

(83)

where $E_i(i=1,2)=\pm \sqrt{A^2-B^2}$ are the two solutions given above, and $Q\left(E_i\right)$ is the determinant (cofactor) obtained by replacing the first column of $\mathbf{M}$ by $\mathbf{C}$ at $E_i$.

The spin length $\langle S^z\rangle$ at a given T is calculated self-consistently by following the method given in Refs. [31,38].

Let us show the SW spectrum $\omega$ (taking $\hslash =1$) for the case of $J=-1$ and $D=0.5$ in Figure 21 versus $k_y$ with $k_x=0$ (Figure 21a) and versus $k_x$ for $k_y=0$ (Figure 21b). In order to see the effect of the DM interaction alone we take the anisotropy $I=0$. One observes here that for a range of small wave-vectors the SW frequency is imaginary. The SW corresponding to these modes do not propagate in the system. Why do we have this case here? The answer is that when the NN make a large angle (perpendicular NN, for example), one cannot define a wave vector in that direction. Physically, when k is small the B coefficient is larger than A in Equation (82) giving rise to imaginary $\omega$. Note that the anisotropy I is contained in A so that increasing I for small k will result in $A>B$ making $\omega$ real.

We show now in Figure 22a the spectrum along the axis $k_x=k_y$ at $T=0$ for $I=0$. Again here the frequency is imaginary for small k, as in the previous figure. The spin length $<S^z>$ along the local quantization axis is shown in Figure 22b. Several remarks are in order: (i) At $T=0$, the spin length is not equal to $1/2$ as in ferromagnets because of the zero-point spin contraction due to antiferromagnetic interactions (see Ref. [31]), its length is $\simeq 0.40$, quite small; (ii) the magnetic ordering is destroyed at $T\simeq 1.2$.

To close the present section, we note that in the case of perpendicular $\mathbf{D}$ considered above, we did not observe skyrmion textures when applying a perpendicular magnetic field: all spin configurations are no more planar, making the calculation of the SW spectrum more difficult. This problem is left for a future investigation.

6. Other Systems of Non-Collinear Ground-State Spin Configurations: Frustrated Surface in Stacked Triangular Thin Films

In this section, we study by the GF technique the effect of a frustrated surface on the magnetic properties of a film composed triangular layers stacked in the z direction. Each lattice site is occupied by a quantum Heisenberg spin of magnitude 1/2. Let the in-plane surface interaction be $J_s$ which can be antiferromagnetic or ferromagnetic. The other interactions in the film are ferromagnetic. We show in the following that the GS spin configuration is non-collinear when $J_s$ is lower than a critical value $J_s^c$. The film surfaces are then frustrated. In the frustrated case, there are two phase transitions, one corresponds to the disordering of the two surfaces and the other to the disordering of the interior layers. The GF results agree qualitatively with Monte Carlo simulation using the classical spins (see the original paper in Ref. [39]).

In this section we review some ot the results given in the original paper Ref. [39], emphasizing the SW calculation and the important results. The Hamiltonian is written as

$$\mathcal{H}=-\sum _{〈i,j〉}{J}_{i,j}{\mathbf{S}}_i·{\mathbf{S}}_j-\sum _{<i,j>}{I}_{i,j}{S}_i^z{S}_j^z$$

(84)

where the first sum is performed over the NN spin pairs ${\mathbf{S}}_i$ and ${\mathbf{S}}_j$, the second sum over their z components. $J_{i,j}$ and $I_{i,j}$ are respectively their exchange interaction and their anisotropic one. The latter is small, taken to ensure the ordering at finite T when the film thickness goes down to a few layers, without this we know that a monolayer with vector spin models does not have a long-range ordering at finite T [40].

Let $J_s$ be the exchange between two NN surface spins. We suppose that all other interactions are ferromagnetic and equal to J. We shall use $J=1$ as the unit of energy in the following.

6.1. Ground State

In the case where $J_s$ is ferromagnetic, the GS of the film is ferromagnetic. When $J_s$ is antiferromagnetic, the situation becomes complicated. We recall that for a single triangular lattice with antiferromagnetic interaction, the spins are frustrated and arranged in a 120-degree configuration [8]. This structure is modified when we turn on the ferromagnetic interaction J with the beneath layer. The competition between the non collinear surface ordering and the ferromagnetic ordering of the bulk leads to an intermediate structure which is determined in the following.

The GS configuration can be determined by using the steepest descent method described below Equation (31). Let us describe qualitatively the GS configuration: when $J_s$ is negative and $J_s<J_s^c$ where $J_s^c(<0)$ is a critical value, the GS is formed by pulling out the planar ${120}^{\circ }$ spin structure along the z axis by an angle $\beta$. This is shown in Figure 23.

Figure 24 shows $cos\alpha$ and $cos\beta$ versus $J_s$ obtained by the steepest descent method. As seen for $J_s>J_s^c$, the angles are zero, namely the GS is ferromagnetic. The critical value $J_s^c$ is numerically found between −0.18 and −0.19.

We show in the following that this value can be analytically calculated by assuming the structure shown in Figure 23). We number the spins as in that figure: $S_1$, $S_2$ and $S_3$ are the spins in the surface layer, $S_1^{\prime }$, $S_2^{\prime }$ and $S_3^{\prime }$ are the spins in the second layer. The energy of the cell is

$$\begin{array}{ccc}\hfill H_p& =& -6\left[J_s\left({\mathbf{S}}_1·{\mathbf{S}}_2+{\mathbf{S}}_2·{\mathbf{S}}_3+{\mathbf{S}}_3·{\mathbf{S}}_1\right)\right.\hfill \\ & & +I_s\left(S_1^z{S}_2^z+S_2^z{S}_3^z+S_3^z{S}_1^z\right)\hfill \\ & +& J\left({\mathbf{S}}_1^{\prime }·{\mathbf{S}}_2^{\prime }+{\mathbf{S}}_2^{\prime }·{\mathbf{S}}_3^{\prime }+{\mathbf{S}}_3^{\prime }·{\mathbf{S}}_1^{\prime }\right)\hfill \\ & & +I\left.\left(S_1^{\prime z}{S}_2^{\prime z}+S_2^{\prime z}{S}_3^{\prime z}+S_3^{\prime z}{S}_1^{\prime z}\right)\right]\hfill \\ & -& 2J\left({\mathbf{S}}_1·{\mathbf{S}}_1^{\prime }+{\mathbf{S}}_2·{\mathbf{S}}_2^{\prime }+{\mathbf{S}}_3·{\mathbf{S}}_3^{\prime }\right)\hfill \\ & & -2I\left(S_1^z{S}_1^{\prime z}+S_2^{\prime z}{S}_2^{\prime z}+S_3^z{S}_3^{\prime z}\right),\hfill \end{array}$$

(85)

We project the spins on the $xy$ plane and on the z axis. One writes ${\mathbf{S}}_i=({\mathbf{S}}_i^{\Vert },S_i^z)$. One observes that only surface spins have non-zero $xy$ vector components. Let the angle between these $xy$ components of NN surface spins be ${\gamma }_{i,j}$ which is in fact the projection of the angle $\alpha$ on the $xy$ plane. By symmetry, we have

$${\gamma }_{1,2}=0,{\gamma }_{2,3}=\frac{2\pi }{3},{\gamma }_{3,1}=\frac{4\pi }{3}.$$

(86)

The angles ${\beta }_i$ and ${\beta }_i^{\prime }$ of ${\mathbf{S}}_i$ and ${\mathbf{S}}_i^{\prime }$ formed with the z axis are by symmetry

$$\left\{\begin{array}{c}{\beta }_1={\beta }_2={\beta }_3=\beta ,\\ {\beta }_1^{\prime }={\beta }_2^{\prime }={\beta }_3^{\prime }=0,\end{array}\right.$$

The total energy of the cell (86), with $S_i=S_i^{\prime }=\frac{1}{2}$, is thus

$$\begin{array}{ccc}\hfill H_p& =& -\frac{9(J+I)}{2}-\frac{3(J+I)}{2}cos\beta -\frac{9(J_s+I_s)}{2}{cos}^2\beta \hfill \\ & +& \frac{9J_s}{4}{sin}^2\beta .\hfill \end{array}$$

(87)

The minimum of the cell energy verifies this condition:

$$\frac{\partial H_p}{\partial \beta }=\left(\frac{27}{2}{J}_s+9I_s\right)cos\beta sin\beta +\frac{3}{2}(J+I)sin\beta =0$$

(88)

One deduces

$$cos\beta =-\frac{J+I}{9J_s+6I_s}.$$

(89)

This solution exists under the condition $-1\le cos\beta \le 1$. The critical values are determined from this condition. For $I=-I_s=0.1$, $J_s^c\approx -0.1889J$ which is in excellent agreement with the results obtained from the steepest descent method.

Now, using the GF method for such a film in the way described in the previous sections, we obtain the full Hamiltonian (84) in the local framework:

$$\begin{array}{ccc}\hfill \mathcal{H}& =& -\sum _{<i,j>}{J}_{i,j}\{\frac{1}{4}\left(cos{\theta }_{ij}-1\right)\left(S_i^{+}{S}_j^{+}+S_i^{-}{S}_j^{-}\right)\hfill \\ & +& \frac{1}{4}\left(cos{\theta }_{ij}+1\right)\left(S_i^{+}{S}_j^{-}+S_i^{-}{S}_j^{+}\right)\hfill \\ & +& \frac{1}{2}sin{\theta }_{ij}\left(S_i^{+}+S_i^{-}\right)S_j^z-\frac{1}{2}sin{\theta }_{ij}{S}_i^z\left(S_j^{+}+S_j^{-}\right)\hfill \\ & +& cos{\theta }_{ij}{S}_i^z{S}_j^z\}-\sum _{<i,j>}{I}_{i,j}{S}_i^z{S}_j^z\hfill \end{array}$$

(90)

where $cos\left({\theta }_{ij}\right)$ is the angle between two NN spins. We define the two coupled GF, and we write their equations of motions in the real space. Taking Tyablikov’s decoupling scheme to reduce higher-order GFs, and then using the Fourier transform in the $xy$ plane we arrive at a matrix equation as in the previous section with the matrix $\mathbf{M}$ is defined as

$$\mathbf{M}\left(\omega \right)=\left(\begin{array}{ccccc}{A}_1^{+}& B_1& D_1^{+}& D_1^{-}& \cdots \\ -B_1& A_1^{-}& -D_1^{-}& -D_1^{+}& ⋮\\ ⋮& \cdots & \cdots & \cdots & ⋮\\ ⋮& C_{N_z}^{+}& C_{N_z}^{-}& A_{N_z}^{+}& B_{N_z}\\ \cdots & -C_{N_z}^{-}& -C_{N_z}^{+}& -B_{N_z}& A_{N_z}^{-}\end{array}\right),$$

(91)

where

$$\begin{array}{ccc}\hfill A_n^{\pm }& =& \omega \pm [\frac{1}{2}{J}_n〈S_n^z〉\left(Z\gamma \right)\left(cos{\theta }_n+1\right)\hfill \\ & -& J_n〈S_n^z〉Zcos{\theta }_n-J_{n,n+1}〈S_{n+1}^z〉cos{\theta }_{n,n+1}\hfill \\ & -& J_{n,n-1}〈S_{n-1}^z〉cos{\theta }_{n,n-1}-Z{I}_n〈S_n^z〉\hfill \\ & & -I_{n,n+1}〈S_{n+1}^z〉-I_{n,n-1}〈S_{n-1}^z〉],\hfill \end{array}$$

(92)

$$\begin{array}{ccc}\hfill B_n& =& \frac{1}{2}{J}_n〈S_n^z〉\left(cos{\theta }_n-1\right)\left(Z\gamma \right),\hfill \end{array}$$

(93)

$$\begin{array}{ccc}\hfill C_n^{\pm }& =& \frac{1}{2}{J}_{n,n-1}〈S_n^z〉\left(cos{\theta }_{n,n-1}\pm 1\right),\hfill \end{array}$$

(94)

$$\begin{array}{ccc}\hfill D_n^{\pm }& =& \frac{1}{2}{J}_{n,n+1}〈S_n^z〉\left(cos{\theta }_{n,n+1}\pm 1\right),\hfill \end{array}$$

(95)

where $Z=6$ is the in-plane coordination number, ${\theta }_{n,n\pm 1}$ denotes the angle between two NN spins belonging to the adjacent layers n and $n\pm 1$, while ${\theta }_n$ is the angle between two NN spins of the layer n, and

$$\gamma =\left[2cos\left(k_{x}a\right)+4cos\left(k_{y}a/2\right)cos\left(k_{y}a\sqrt{3}/2\right)\right]/Z.$$

Note that in the above coefficients, we have used the following notations:

(i)

$J_n$ and $I_n$ are the in-plane interactions. $J_n$ is equal to $J_s$ for the two surface layers and equal to J for the interior layers. All $I_n$ are taken equal to I.

(ii)

The interlayer interactions are denoted by $J_{n,n\pm 1}$ and $I_{n,n\pm 1}$. Note that $J_{n,n-1}=I_{n,n-1}$ = 0 if $n=1$ and $J_{n,n+1}=I_{n,n+1}$=0 if $n=N_z$.

As described in the previous sections, the SW spectrum $\omega$ is obtained by solving det$\left|\mathbf{M}\right|=0$. Using $\omega$ we calculate the magnetizations layer by layer for typical values of parameters. The results are shown in the following.

6.2. Quantum Surface Phase Transition

Let us show a typical case in the region of frustrated surface where $J_s=-0.5$ in Figure 25. Several comments are in order:

(i)

The surface magnetization is very small with respect to the magnetization of the second layer,

(ii)

At $T=0$, the length of the surface spin is about 0.425, much shorter than the spin magnitude 1/2. This is due to the antiferromagnetic interaction at the surface which causes a strong spin contraction. For the second layer, the spins are aligned ferromagnetically, their length is fully 0.5,

(iii)

The surface undergoes a phase transition at $T_1\simeq 0.2557$ while the second layer remains ordered up to $T_2\simeq 1.522$. The system is thus disordered at the surface and ordered in the bulk, for temperatures between $T_1$ and $T_2$. This partial disorder is very interesting. It gives another example of the partial disorder observed earlier in bulk frustrated quantum spin systems.

(iv)

One observes that between $T_1$ and $T_2$, the first layer has a small magnetization. This is understood by the fact that the strong magnetization of the second layer acts as an external field on the first layer, inducing therefore a small value of its magnetization.

We plot the phase diagram in the space $(J_s,T)$ in Figure 26. Phase I denotes the surface canted-spin state, phase IIA denotes the partially ordered phase: the surface is disordered while the bulk is ordered. Phase IIB separated from phase IIA by a vertical line issued from $J_s^c\simeq -0.19$ indicates the ferromagnetic state, and phase III is the paramagnetic phase.

6.3. Classical Phase Transition: Monte Carlo Results

In order to compare with the quantum model shown in the previous subsection, we consider here the classical counterpart model, namely we use the same Hamiltonian (28) but with the classical Heisenberg spin of magnitude $S=1$. The aim is to compare their qualitative features, in particular the question of the partial disordering at finite T.

We use Monte Carlo simulations for the classical model where the film dimensions are $N\times N\times N_z$, $N_z$ being the film thickness which is taken to be $N_z=4$ as in the quantum case shown above. We use here $N=24,36,48,60$ to see the lateral finite-size effect. Periodic boundary conditions are used in the $xy$ planes. We discard ${10}^6$ MC steps per spin to equilibrate the system and average physical quantities over the next $2\times {10}^6$ MC steps per spin.

We show in Figure 27 the result obtained in the same frustrated case as in the quantum case shown above, namely $J_s=-0.5$. we see that the surface magnetization falls at $T_1\simeq 0.25$ while the second-layer magnetization stays ordered up to $T_2\simeq 1.8$. This surface disordering at low T is similar to the quantum case. Between $T_1$ and $T_2$ the system is partially disordered.

Figure 28 shows the phase diagram obtained in the space $(J_s,T)$. It is interesting to note that the classical phase diagram shown here has the same feature as the quantum phase diagram displayed in Figure 26. The difference in the values of the transition temperatures is due to the difference of quantum and classical spins.

To close this review, we should mention a few works works where SW in the regime of non-collinear spin configurations have been studied: the frustration effects in antiferromagnetic face-centered cubic Heisenberg films have been studied in Ref. [49], a frustrated ferrimagnet in Ref. [50] and a quantum frustrated spin system in Ref. [51]. These results are not reviewed here to limit the paper’s length. The reader is referred to those works for details.

7. Concluding Remarks

As said in the Introduction, the self-consistent Green’s function theory is the only one which allows to calculate the SW dispersion relation in the case of non-collinear spin configurations, in two and three dimensions, as well as in thin films. The non-collinear spin configurations are due to the existence of competing interactions in the system, to the geometry frustration such as in the antiferromagnetic triangular lattice, or to the competition between ferromagnetic and/or antiferromagnetic interactions with the Dzyaloshinskii–Moriya interaction. We have shown that without an applied magnetic field, the GS spin configuration is non-collinear but periodic in space. We have, in most cases, analytically calculated them. We have checked them by using the iterative numerical minimization of the local energy (the so-called steepest-descent method). The agreement between the analytical method and the numerical energy minimization is excellent. The determination of the GS is necessary because we need them to calculate the SW spectrum: SW are elementary excitations of the GS when T increases.

The double-fold purpose of this review is to show the method and the interest of its results. We have reviewed a selected number of works according to their interest of the community: helimagnets, materials with the Dzyaloshinskii–Moriya interaction, and the surface effects in thin magnetic films. The Dzyaloshinskii–Moriya interaction gives rise not only a chiral order but also the formation of skyrmions in an applied magnetic field. The surface effects in helimagnets and in films with a frustrated surface give rise to the reconstruction of surface spin structure and many striking features due to quantum fluctuations at low T such as the zero-point spin contraction and the magnetization crossover). We have also seen above the surface becomes disordered at a low T while the bulk remains ordered up to a high T. This coexistence of bulk order and surface disorder in a temperature region is also found in several frustrated systems [8].

To conclude, we say that the Green’s function theory for non-collinear spin systems is laborious, but it is worthwhile to use it to get results with clear physical mechanisms lying behind observed phenomena in frustrated spin systems.