Quantum Spin-Wave Theory for non-collinear Spin Structures, a Review

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\'e 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.


I. INTRODUCTION
In a solid the interaction 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 interaction known as phonons or due to spin-spin interaction known as magnons. Note that magnons are spin waves (SW) when they are quantized. Elementary excitations are defined also for interaction between charge densities in plasma, or for electric dipole-dipole interaction 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, its energy depends on the wave vector k. The 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 50's. 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 pioneer 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 70's and 80's. 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 parameter, 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 noncollinear spin configurations in spin systems. While the Heisenberg model between two spins is written as −J ij S i · S j giving rise to two collinear spins in the GS, the DM interaction is written as D ij ·S i ×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 discoververed 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]. Skyrmions are among the most studied subjects at the time being due to their potental applications in spinelectronics. [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, 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 II, 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 III the numerical results obtained from the GF. Section IV 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 V treats a case where the DM interaction competes with the antiferromagnetic interaction in the frustrated antiferromagnetic triangular lattice. Section VI presents the surface effect in a thin film where its surface is frustrated. Concluding remarks are given in section VII.

II. 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 Fig.1).
We consider the simplest model of helimagnet given by the following Hamiltonian where 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. Fig. 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 θ, one has (omitting the factor S 2 ) The lowest-energy state corresponds to There are two solutions, sin θ = 0 and cos θ = − J1 J2 The first solution corresponds to the ferromagnetic state, and the second solution exists if − J1 J2 ≤ 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]).  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 noncollinear 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 S i and S j . As seen above, in the general case these spins make an angle θ i,j = θ j − θ 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 Fig. 3 the local coordinates assigned to spin S i and S j . We write Expressing the axes of S j in the frame of S i one haŝ so that FIG. 3: Spin Si lies along the ζ axis (its quantization axis), while spin Sj lies along its quantization axis ζ which makes an angle θ with the ζ axis. The axes ξ and ξ are perpendicular respectively to ζ and ζ . T he perpendicular axesηi andηj coincide with the c axis, perpendicular to the basal plane of the bct lattice.
Using Eq. (9) to express S j in the (ξ i ,η i ,ζ i ) coordinates, we calculate S i · S j , we get the following Hamiltonian from (28): 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 IV. We define the following GFs for the above Hamiltonian: Writing their equations of motion we have where For simplicity, we suppose in the following J k,k are all equal to J 1 for NN interactions and to J 2 for NNN interactions. I k,k 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 R i by three indices ( 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 where ω is the SW frequency, k xy the wave-vector parallel to xy planes and R i the position of S i . n i , n j and n k denote the z-components of the sites R i , R j and R k . The integral over k xy is performed in the first Brillouin zone (BZ) whose surface is ∆ 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 Eqs. (17)- (18) in the three (x, y, z) directions yields the SW spectrum in the absence of anisotropy: where where Z = 8 is the NN coordination number, Z c = 2 the NNN number on the c-axis and γ = 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 ω is zero when A = ±B. This is realized at two points as expected in helimagnets: k x = k y = k z = 0 (γ = 1) and k z = 2θ 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 θ = 1 (−1) in the above coefficients.
In order to obtain the SW fequency ω, we solve the secular equation det |M| = 0 for each given (k x , k y ). Since the linear dimension of the square marix is 2N z , we obtain 2N z eigen-values of ω, half positive and half negative, corresponding to two opposite spin precessions as in antiferromagnets. These values depend on the input values < S z n > (n = 1, ..., 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 z n > 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]),: where n = 1, ..., N z , and As < S z n > depends each other in ω i (i = 1, ..., 2N z ), their solutions should be obtained by iteration at a given temperature T . In the particular case where T = 0 one has Note that the sum is performed over N z negative ω i since for positive ω i yield the zero Bose-Einstein factor at T = 0). The transition temperature T c can be calculated self-consistently when all < S z n > tend to zero. We show in the following section, the numerical results using the above formulas.

III. 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 θ n,n±1 and θ n,n±2 in the matrix elements of |M|, then calculating ω i (i = 1, ..., 2N z ) for each k xy . For the iterative procedure, the reader is referred to Re. [29]. The solution S z n (n = 1, ..., N z ) is obtained when the input and the output are equal with a desired precision P .

A. 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 Fig. 4 the SW spectrum ω versus k x = k y for an 8-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 Fig. 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.

B. Zero-point spin contraction and transition temperature
Ii isknown 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 Fig. 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 → −1.

C. Layer magnetizations
We show now the layer ordering in Figs. 6 and 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 Fig. 5, crosses the interior layer magnetizations at T 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 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 Fig. 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 a 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] IV. 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 where J i,j and D i,j are the exchange and DM interactions, respectively, between two quantum Heisenberg spins S i and S j of magnitude S = 1/2. We supppose 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 S i makes an angle θ i,j with its neighbor S j . It is obvious that the quantization axes of S i and S j are different. Therefore, the transformation using the local coordinates, Eqs. (4)-(9), is necessary. Let us suppose that the vector D i,j is along the y axis, namely theη i axis. We write where e i,j =+1 (-1) if j > i (j < i) for NN j on thex orẑ 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 Fig. 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 θ on J 1 : the energy of the spin S i is written as where θ = |θ i,j | minimizing E i with respect to θ one obtains The result is in agreement with that obtained by the steepest descent method. An example has been shown in Fig. 8. We rewrite the DM term of Eq. (30) as From Eq. (31), we obtain where we have replaced S x by (S + + S − )/2. Note that e i,j sin θ i,j is always positive since for a NN on the positive axis direction, e i,j = 1 and sin θ i,j = sin θ where θ is positively defined, while for a NN on the negative axis direction, e i,j = −1 and sin θ i,j = sin(−θ) = − sin θ.

A. Formulation of the Green's function technique for the Dzyaloshinskii-Moriya system
Using the transformation into the local coordinates, Eqs. (4)-(9), one has Note that the quantization axes of the spins are in the xz planes as shown in Fig. 3. We emphasize that while the sinus term of the DM Hamiltonian, Eq. (35), remain after summing over the NN, the sinus terms of H e , the 3rd line of Eq. (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 ± 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 S i and S j as follows 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 Using the two GF's in the real space given by Eqs. (11)-(12) and using the same method, we study the effect of the DM interaction. For the DM term, the commutation relations [H, S ± i ] lead to: which gives rise, using the Tyablikov decoupling, to the following GF's: 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 II, the Fourier transforms in the xz plane g n,n and f n,n of the G and F lead to the matrix equation M (E) being given by Eq. (42) below · · · · · · · · · · · · · · · · · · · · · · · · · · · 0 · · · · · · · · · · · · · · · · · · · · · · · · 0 0 0 0 0 where E = ω is the SW energy and the matrix elements are given by where n = 1, 2, ..., N denoting the layer numbers, d n = I 1 /J 1 , γ = (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 frequecies at a given wave vector k = (k x , kz) 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.

B. Results for 2D and 3D cases
In the 2D case, one has only one layer. The matrix (42) is (E + A n )g n,n + B n f n,n = 2 < S z n > δ(n, n ) −B n g n,n + (E − A n )f n,n = 0 (46) where A n is given by (43) but without J 2 term for the 2D case. Coefficients B n is given by (44) and C n = 0. The SW frequencies are determined by the following secular equation Several remarks are in order: (i) when θ = 0, the last three terms of A n and B n are zero: one recovers the ferromagnetic SW dispersion relation where Z = 4 is the coordination number of the square lattice (taking d n = 0), (ii) when θ = π, one has A n = 8J 1 < S z n >, B n = −8J 1 < S z n > γ. One recovers then the antiferromagnetic SW dispersion relation (iii) when there is a DM interaction, one has 0 < cos θ < 1 (0 < θ < π/2). If d n = 0, the quantity in the square root of Eq. (47) becomes negative at γ = 1 when θ 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 Fig. 9 we show d c versus θ. As seen, d c increases from zero with increasing θ.
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 Fig. 10 the SW dispersion relation calculated from Eq. (47) for θ = 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 θ, namely small D, E(k) is proportional to k 2 at low k (see Fig. 10a). This behavior is that in ferromagnets. For large θ, one observes that E(k) becomes linear in k as seen in Fig. 10b. This behavior is similar to that of antiferromagnets. Note that the change of behavior is progressive with increasing θ, 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 to have a llong-range ordering at T = 0 in 2D: Fig. 11 displays the magnetization M (≡< S z >) calculated by Eq. (II) where in each case the limit value d c has been used. We note that M depends strongly on θ: at high T the larger θ the stronger M . However, at T = 0 the spin length is smaller for larger θ due to the zero-point spin contraction [31] calculated by Eq. (27). As a consequence there is a cross-over of layer magnetizations at low T as shown in Fig. 11b. The spin length at T = 0 is shown in Fig. 12 for several θ.
We now consider the 3D case. The crystal is infinite in three direction. The Fourier transform in the y direction, namely g n±1 = g n e ±ikya and f n±1 = f n e ±ikya reduces the matrix (23) to two coupled equations of g and f functions. One has The spectrum is given by In the ferromagnetic case, cos θ = 1, thus B = 0. Arranging the Fourier transforms in three directions, one gets the 3D ferromagnetic dispersion relation E = 2Z < S z > (1 − γ 2 ) where γ = [cos(k x a) + cos(k y a) + cos(k z a)]/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 θ. This is rather obvious because the DM interaction operates in the plane making an angle θ 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 Eq. (53), we calculate the 3D spectrum displayed in Fig. 13 for a small and a large value of θ. As in the 2D case, we observe E ∝ k when k → 0 for large θ. The main properties of the system are thus governed by the in-plane DM interaction.
FIG. 13: Spin-wave spectrum E(k) versus k ≡ kx = kz for θ = 0.1 (red crosses) and θ = 0.6 (blue circles) in three dimensions. Note the linear-k behavior at low k for the large value of θ. See text for comments. Figure 14 displays the magnetization M versus T for several values of θ. As in the 2D case, when the DM interaction is included, the spins undergo a zero-point contraction which increases with increasing θ. The competition between quantum fluctuations at T = 0 and thermal effects at high T gives rise to magnetization cross-over shown in Fig. 14b. The spin length at T = 0 vs θ is shown in the inset of Fig. 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.

C. The case of a thin film
As in 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 Fig. 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 8 layers with J 1 = J 2 = 1 for a small and a large θ. As in the previous cases, E is proportional to k for large θ (cf. Fig. 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 θ. 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 Fig. 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].

V. 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.

A. 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 Fig. 17a), then the DM vector D i,j is perpendicular to the xy plane and is given by [45,46] where r iO = r O − r i and r Oj = r j − r O , r ij = r j − r i . r O is the position of non-magnetic ion (oxygen) and r i the position of the spin S i etc. These vectors are defined in Fig. 17a in the particular case where the displacements are in the xy plane. We have therefore D i,j perpendicular to the xy plane in this case. Note however that if the atom displacements are in 3D space, D i,j can be in any direction. In this paper, we consider also the case where D i,j lies in the xy plane as shown in Fig. 17b where D i,j is taken along the vector connecting spin S i to spin S j .
Note that from Eq. (54) one has In the case of perpendicular D i,j , let us define u i,j as the unit vector on the z axis. From Eqs. (54)-(55) one writes 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 D i,j , we suppose that D i,j is given as where D is a constant and r ij denotes the unit vector along r j − r i . The case of in-plane D i,j on the frustrated triangular lattice (see Fig. 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.

B. Ground State with a Perpendicular D in Zero Field
The Hamiltonian is given by where 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.
FIG. 18: Perpendicular Di,j: Ground-state spin configuration with only Dzyaloshinskii-Moriya interaction on the triangular lattice (J = 0) is analytically determined. One angle is 120 degrees, the other two are 60 degrees. Note that the choice of the 120-degree angle in this figure is along the horizontal spin pair. This configuration is one GS, the other two GSs have the 120-degree angles on respectively the two diagonal spin pairs. Note also that the spin configuration is invariant under the global spin rotation in the xy plane. For convenience, the spins are decomposed into three sublattices numbered 1, 2 and 3. See text for explanation.
Note that this solution can be numerically obtained by the steepest descent method described above. The result is shown in Fig. 19 for the full lattice. We see in the zoom that the spin configuration on a plaquette is what obtained analytically, with a global spin rotation as explained in the caption of Fig. 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 × 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 ≥ 100. Most of calculations have been performed for L = 100.

C. Ground State with both perpendicular D and J in Zero Field-Spin Waves
When both J and perpendicular 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 D is perpendicular and without J, the spins are in the plane. When J is antiferromagnetic without D, the spins are also in the plane and form a 120-degree structure. When D and J exist together the angles between NN's change but they still in the plane in order to keep both D and J interactions as low as possible. An example is shown in Fig. 20 where one sees that the GS is planar and characterized by two angles θ = 102 degrees and one angle β = 156 degrees formed by three spins on a triangle plaquette. Note that there are three degenerate states where β is chosen for the pair (1,2) (Fig. 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 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 β and two angles θ as shown in Fig. 20. The values of β and θ depend on the value of D. We 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 where θ ij is the angle between S i and 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 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 We define the two GFs by Eqs. (11)- (12) and use the equations of motion of these functions (13)- (14), we obtain Note that < S z i > 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 and where µ i ≡< S z i >, k is the wave vector in the reciprocal lattice of the triangular lattice, and ω the SW frequency. Note that the index z in S z i is not referred to the real space direction z, but to the quantization axis of the spin S i . At this stage, we have to replace θ i,j by either β or θ according on the GS spin configuration given above (see Fig.  20).
As in the previous sections, writing the above equations under a matrix form, we have where M ( ω) is a square matrix of dimension 2 × 2, h and C are given by and the matrix M ( ω) is given by The nontrivial solution of g and f imposes the following secular equation: where A = −J(8µ i cos β(1 + I) + 4µ i cos θ(1 + I) − 4µ i γ(cos β + 1) − 2µ i α(cos θ + 1)) − D(4µ i sin βγ + 2µ i sin θα) + D(8µ i sin β + 4µ i sin θ) (81) where the sum on the two NN on the x axis (see Fig. 20b) is and the sum on the four NN on the oblique directions of the hexagon (see Fig. 20b) is Solving Eq. (80) for each given (k x , k y ) one obtains the SW frequency ω(k x , ky): Plotting ω(k x , ky) in the space (k x , ky) one obtains the full SW spectrum.
The spin length S z i (for all i, by symmetry) is given by (see technical details in Ref. [31]): where E i (i = 1, 2) = ± √ A 2 − B 2 are the two solutions given above, and Q(E i ) is the determinant (cofactor) obtained by replacing the first column of M by C at E i .
The spin length S z at a given T is calculated self-consistently by following the method given in Refs. [31,38]. Let us show the SW spectrum ω (taking = 1) for the case of J = −1 and D = 0.5 in Fig. 21 versus k y with k x = 0 (Fig. 21a) and versus k x for k y = 0 (Fig. 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 Eq. (85) giving rise to imaginary ω. Note that the anisotropy I is contained in A so that increasing I for small k will result in A > B making ω real. We show now in Fig. 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 Fig. 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 0.40, quite small; ii) the magnetic ordering is destroyed at T 1.2.
To close the present section, we note that in the case of perpendicular 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.

VI. 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 c s . The film surfaces are then frustrated. In the frustrated case, there are two phase transitions, one correponds 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 where the first sum is performed over the NN spin pairs S i and 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 unit of energy in the following.

A. 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 configuratrion. [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 Eq. (31). Let us describe qualitatively the GS configuration: when J s is negative and J s < J c s where J c s (< 0) is a critical value, the GS is formed by pulling out the planar 120 • spin structure along the z axis by an angle β. This is shown in Fig. 23). Figure 24 shows cos α and cos β versus J s obtained by the steepest descent method. As seen for J s > J c s , the angles are zero, namely the GS is ferromagnetic. The critical value J c s 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 Fig.  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 , S 2 and S 3 are the We project the spins on the xy plane and on the z axis. One writes S i = (S i , S z i ). One observes that only surface spins have non-zero xy vector components. Let the angle between these xy components of NN surface spins be γ i,j which is in fact the projection of the angle α on the xy plane. By symmetry, we have The angles β i and β i of S i and S i formed with the z axis are by symmetry The total energy of the cell (88), with S i = S i = 1 2 , is thus The minimum of the cell energy verifies this condition: One deduces This solution exists under the condition −1 ≤ cos β ≤ 1. The critical values is determined from this condition. For I = −I s = 0.1, J c s ≈ −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 (87) in the local framework: where cos (θ ij ) 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 the 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 M is defined as . . . · · · · · · · · · . . . 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 × 10 6 MC steps per spin. We show in Fig. 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 0.25 while the second-layer magnetization stays ordered up to T 2 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 Fig. 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.

VII. 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 worth to use it to get results with clear physical mechanisms lying behind observed phenomena in frustrated spin systems.