Phase transition in frustrated thin films -- physics at phase boundaries

TL;DR

This review explores the physics of frustrated systems and thin films, emphasizing exact solutions in low dimensions and phenomena at phase boundaries caused by competing interactions, with applications to spin systems and skyrmions.

Contribution

It provides a comprehensive overview of exact solutions in 2D frustrated spin systems and discusses phenomena at phase boundaries in thin films, including new results on skyrmions.

Findings

Exact solutions reveal high ground-state degeneracy and multiple phase transitions.

Spectacular effects like reentrance and partial disorder are observed.

Surface phenomena and skyrmion behavior in thin films are analyzed.

Abstract

In this review, I outline some principal theoretical knowledge on the properties of frustrated systems and thin films. The two points I would like to emphasize: i) the physics in low dimensions where exact solutions can be obtained, ii) the physics at phase boundaries where spectacular phenomena can occur due to competing interactions of the two phases around the boundary. This competition causes a frustration. I will concentrate my attention to thin films and phenomena occurring near the boundary of two phases of different symmetries. The case of two-dimensional (2D) systems is in fact the limiting case of thin films with a monolayer. Naturally, I will treat this case at the beginning. After a short introduction on frustrated spin systems, I show several 2D frustrated Ising spin systems which can be exactly solved by using vertex models. These systems contain most of the spectacular effects due to the frustration: high ground-state degeneracy, existence of several phases in the ground-state phase diagram, multiple phase transitions with increasing temperature, reentrance, disorder lines, partial disorder at equilibrium. Evidences of such effects in non solvable models are also shown and discussed. Thin films are next presented with different aspects: surface elementary excitations (surface spin-waves), surface phase transition and criticality. Several examples are shown and discussed. New results on skyrmions in thin films and superlattices are also displayed.

Tables

Table 1. Table 1: Critical exponents, effective dimension and critical temperature at infinite x ​ y 𝑥 𝑦 xy limit as obtained in this paper.

$N_z$	$\nu$	$\gamma$	$\alpha$	$\beta$	$d_{\mathrm{eff}}$	$T_c\left(L=\infty ,N_z\right)$
1	$0.9990\pm 0.0028$	$1.7520\pm 0.0062$	$0.00199\pm 0.00279$	$0.1266\pm 0.0049$	$2.0000\pm 0.0028$	$2.2699\pm 0.0005$
3	$0.9922\pm 0.0019$	$1.7377\pm 0.0035$	$0.00222\pm 0.00192$	$0.1452\pm 0.0040$	$2.0135\pm 0.0019$	$3.6365\pm 0.0024$
5	$0.9876\pm 0.0023$	$1.7230\pm 0.0069$	$0.00222\pm 0.00234$	$0.1639\pm 0.0051$	$2.0230\pm 0.0023$	$4.0234\pm 0.0028$
7	$0.9828\pm 0.0024$	$1.7042\pm 0.0087$	$0.00223\pm 0.00238$	$0.1798\pm 0.0069$	$2.0328\pm 0.0024$	$4.1939\pm 0.0032$
9	$0.9780\pm 0.0016$	$1.6736\pm 0.0084$	$0.00224\pm 0.00161$	$0.1904\pm 0.0071$	$2.0426\pm 0.0016$	$4.2859\pm 0.0022$
11	$0.9733\pm 0.0025$	$1.6354\pm 0.0083$	$0.00224\pm 0.00256$	$0.1995\pm 0.0088$	$2.0526\pm 0.0026$	$4.3418\pm 0.0032$
13	$0.9692\pm 0.0026$	$1.6122\pm 0.0102$	$0.00226\pm 0.00268$	$0.2059\pm 0.0092$	$2.0613\pm 0.0027$	$4.3792\pm 0.0034$

Equations

$$P=\prod _{⟨i,j⟩}\mathrm{sign}(J_{i,j}),$$

$$E$$

$$\frac{\partial E}{\partial {\theta }_1}$$

$$\frac{\partial E}{\partial {\theta }_2}$$

$$\frac{\partial E}{\partial {\theta }_3}$$

$$E=J({\mathbf{S}}_1\cdot {\mathbf{S}}_2+{\mathbf{S}}_2\cdot {\mathbf{S}}_3+{\mathbf{S}}_3\cdot {\mathbf{S}}_1)=-\frac{3}{2}J{S}^2+\frac{J}{2}({\mathbf{S}}_1+{\mathbf{S}}_2+{\mathbf{S}}_3{)}^2.$$

$${\theta }_2-{\theta }_1={\theta }_3-{\theta }_2={\theta }_4-{\theta }_3=\frac{\pi }{4}\text{ and }{\theta }_1-{\theta }_4=\frac{3\pi }{4}$$

$$\cos {\theta }_{32}=\cos {\theta }_{43}=\cos {\theta }_{14}\equiv \theta =\frac{1}{2}[\frac{\eta +1}{\eta }{]}^{1/2}$$

$$E$$

$$\frac{\partial E}{\partial \theta }$$

$$\sin \theta =0⟶\theta =0\text{(ferromagnetic solution)}$$

$$\cos \theta =\frac{J_1}{4|J_2|}⟶\theta =\pm \arccos \left(\frac{J_1}{4|J_2|}\right).$$

$$H=-J_1\sum _{(ij)}{\sigma }_i{\sigma }_j-J_2\sum _{(ij)}{\sigma }_i{\sigma }_j-J_3\sum _{(ij)}{\sigma }_i{\sigma }_j$$

$$W=\exp [K_1({\sigma }_1{\sigma }_2+{\sigma }_2{\sigma }_3+{\sigma }_3{\sigma }_4+{\sigma }_4{\sigma }_5+{\sigma }_5{\sigma }_6+{\sigma }_6{\sigma }_1)+$$

$$K_2\sigma ({\sigma }_1+{\sigma }_2+{\sigma }_4+{\sigma }_5)+K_3\sigma ({\sigma }_3+{\sigma }_6)]$$

$$Z=\sum _{\sigma }\prod _{c}W$$

$$W^{\prime }=2\exp [K_1({\sigma }_1{\sigma }_2+{\sigma }_2{\sigma }_3+{\sigma }_3{\sigma }_4+{\sigma }_4{\sigma }_5+{\sigma }_5{\sigma }_6+{\sigma }_6{\sigma }_1)]\times$$

$$\cosh [K_2({\sigma }_1+{\sigma }_2+{\sigma }_4+{\sigma }_5)+K_3({\sigma }_3+{\sigma }_6)]$$

$$\omega$$

$$\overline{\omega }$$

$${\omega }_{56}$$

$${\overline{\omega }}_{56}$$

$${\omega }_{15}$$

$${\overline{\omega }}_{15}$$

$${\omega }_{46}$$

$${\overline{\omega }}_{46}$$

$${\omega }_{13}$$

$${\overline{\omega }}_{13}$$

$${\omega }_{24}$$

$${\overline{\omega }}_{24}$$

$${\omega }_{14}$$

$${\overline{\omega }}_{14}$$

$${\omega }_{23}$$

$${\overline{\omega }}_{23}$$

$${\omega }_{25}$$

$${\overline{\omega }}_{25}$$

$${\omega }_{36}$$

$${\overline{\omega }}_{36}$$

$${\overline{\omega }}_{34}$$

$${\omega }_{34}$$

$${\overline{\omega }}_{35}$$

$${\omega }_{35}$$

$${\overline{\omega }}_{45}$$

$${\omega }_{45}$$

$${\overline{\omega }}_{12}$$

Full text

• •

PACS numbers: 75.10.-b ; 75.10.Hk ; 64.60.Cn

• •

Keywords: Frustration; Phase Transition; Reentrance ; Disorder Lines;

                  Surface Spin-Waves ; Thin Films; Theory; Simulation.

††preprint: APS/123-QED

I Introduction

Materials science has been extensively developed during the last 30 years. This is due to an enormous number of industrial applications which drastically change our life style. The progress in experimental techniques, the advance on theoretical understanding and the development of high-precision simulation methods together with the rapid increase of computer power have made possible the spectacular development in materials science. Today, it is difficult to predict what will be discovered in this research area in ten years.

The purpose of this review is to look back at early and recent results in the physics of frustrated systems at low dimensions: 2D systems and thin films. We would like to connect these results, published over a large period of time, on a line of thoughts: physics at phase boundaries. A boundary between two phases of different orderings is determined as a compromise of competing interactions each of which favors one kind of ordering. The frustration is thus minimum on the boundary (see reviews on many aspects of frustrated systems in Ref. DiepFSS ). When an external parameter varies, this boundary changes and we will see in this review that many interesting phenomena occur in the boundary region. We will concentrate ourselves in the search for interesting physics near the phase boundaries in various frustrated systems in this review.

The study of order-disorder phase transition is a fundamental task of equilibrium statistical mechanics DiepSP ; Zinn . Great efforts have been made to understand the basic mechanisms responsible for spontaneous ordering as well as the nature of the phase transition in many kinds of systems. We will show methods to study the phase transition in thin films where surface effects when combined with frustration effects give rise to many new phenomena. Surface physics has been intensively developed in the last two decades due to its many applications zangwill ; bland-heinrich ; DiepTM ; Baibich ; Grunberg ; Fert ; Fert2013 .

A large part of this review, section II, is devoted to the definition of the frustration and to models which are exactly solved. We begin with exactly solved models in order to have all properties defined without approximation. As seen, many spectacular phenomena are exactly uncovered such as partial disorder, reentrance, disorder lines and multiple phase transitions. Only exact mathematical techniques can allow us to reveal such beautiful phenomena which occur around the boundary separating two phase of different ground-state orderings. These exact results permit to understand similar behaviors in systems that cannot be solved such as 3D systems.

In section IV, an introduction on surface effects in thin films is given. In order to avoid a dispersion of techniques, I introduce only the Green’s function method which can be generalized in more complicated cases such as non-collinear spin states. Calculations of the spin-wave spectrum and the surface magnetization are in particular explained.

In section V a number of striking results obtained mainly by my group are shown on a number of frustrated thin films including helimagnetic films. We show in particular the surface phase transition, quantum fluctuations at low temperature, and the existence of partial phase transition. Results obtained by Monte Carlo simulations are also shown in most cases to compare with the Green’s function technique.

The question of the criticality in thin films is considered in section VI. Here, the high-precision multi-histogram techniques are used to show that critical exponents in thin films are effective exponents having values between those of the 2D and 3D universality classes.

Section VII is devoted to skyrmions, a hot subject at the time being due to their numerous possible applications. Here again, we show only results obtained in the author’s group, but we mention a large bibliography. Skyrmions are topological excitations. They are a kind of circular domain walls involving a number of spins. Skyrmions are shown to result from the competition of different antagonist interactions under an applied magnetic field. We find the existence of a skyrmion crystal, namely a network of periodically arranged skyrmions. Results show that such a skyrmion crystal is stable up to a finite temperature. The relaxation time of skyrmions is shown to follow a stretched exponential law.

Concluding remarks are given in section VIII.

II Physics in two dimensions: frustration effects

II.1 Frustration

During the last 30 years, much attention has been paid to frustrated models DiepFSS . The concept of ”frustration” has been introduced Tou ; Villain1 to describe the situation where a spin (or a number of spins) in the system cannot find an orientation to fully satisfy all the interactions with its neighboring spins (see below). This definition can be applied to Ising spins, Potts models and vector spins. In general, the frustration is caused either by competing interactions (such as the Villain model Villain1 ) or by lattice structure as in the triangular, face-centered cubic (fcc) and hexagonal-close-packed (hcp) lattices, with antiferromagnetic nearest-neighbor (nn) interaction. The effects of frustration are rich and often unexpected (see DiepFSS ).

In addition to the fact that real magnetic materials are often frustrated due to several kinds of interactions, frustrated spin systems have their own interest in statistical mechanics. Recent studies show that many established statistical methods and theories have encountered many difficulties in dealing with frustrated systems. In some sense, frustrated systems are excellent candidates to test approximations and improve theories.

Since the mechanisms of many phenomena are not understood in real systems (disordered systems, systems with long-range interaction, three-dimensional systems, etc), it is worth to search for the origins of those phenomena in exactly solved systems. These exact results will help to understand qualitatively the behavior of real systems which are in general much more complicated.

Let us give here some basic definitions to help readers unfamiliar with these subjects.

Consider two spins $\mathbf{S}_{i}$ and $\mathbf{S}_{j}$ with an interaction $J$. The interaction energy is $E=-J\left(\mathbf{S}_{i}\cdot\mathbf{S}_{j}\right)$. If $J$ is positive (ferromagnetic interaction) then the minimum of $E$ is $-J$ corresponding to the configuration in which $\mathbf{S}_{i}$ is parallel to $\mathbf{S}_{j}$. If $J$ is negative (antiferromagnetic interaction), the minimum of $E$ corresponds to the configuration where $\mathbf{S}_{i}$ is antiparallel to $\mathbf{S}_{j}$. It is easy to see that in a spin system with nn ferromagnetic interaction, the ground state (GS) of the system corresponds to the spin configuration where all spins are parallel: the interaction of every pair of spins is fully satisfied. This is true for any lattice structure. If $J$ is antiferromagnetic, the spin configuration of the GS depends on the lattice structure: i) for lattices containing no elementary triangles, i.e. bipartite lattices (such as square lattice, simple cubic lattices, …) the GS is the configuration in which each spin is antiparallel to its neighbors, i.e. every bond is fully satisfied. ii) for lattices containing elementary triangles such as the triangular lattice, the fcc lattice and the hcp lattice, one cannot construct a GS where all bonds are fully satisfied (see Fig. 1). The GS does not correspond to the minimum of the interaction of every spin pair. In this case, one says that the system is frustrated.

The first frustrated system which was studied in 1950 is the triangular lattice with Ising spins interacting with each other via a nn antiferromagnetic interaction Wan . For vector spins, non collinear spin configurations due to competing interactions were first discovered in 1959 independently by Yoshimori Yos , Villain Vill and Kaplan Kapl .

Consider an elementary cell of the lattice. This cell is a polygon formed by faces hereafter called ”plaquettes”. For example, the elementary cell of the simple cubic lattice is a cube with six square plaquettes, the elementary cell of the fcc lattice is a tetrahedron formed by four triangular plaquettes. Let $J_{i,j}$ be the interaction between two nn spins of the plaquette. According to the definition of Toulouse Tou the plaquette is frustrated if the parameter $P$ defined below is negative

[TABLE]

where the product is performed over all $J_{i,j}$ around the plaquette. Two examples of frustrated plaquettes are shown in Fig. 1: a triangle with three antiferromagnetic bonds and a square with three ferromagnetic bonds and one antiferromagnetic bond. $P$ is negative in both cases. One sees that if one tries to put Ising spins on those plaquettes, at least one of the bonds around the plaquette will not be satisfied. For vector spins, we show below that in the lowest energy state, each bond is only partially satisfied.

One sees that for the triangular plaquette, the degeneracy is three, and for the square plaquette it is four, in addition to the degeneracy associated with returning all spins. Therefore, the degeneracy of an infinite lattice composed of such plaquettes is infinite, in contrast to the unfrustrated case.

We emphasize that the frustration can be created with other kinds of interaction such as the Dzyaloshinski-Moriya interaction $E=-\mathbf{D}\cdot\left(\mathbf{S}_{i}\wedge\mathbf{S}_{j}\right)$ Dzyaloshinskii ; Moriya which favors the perpendicular spin configuration in competition with a Heisenberg exchange model which favors a collinear one. We will return to this interaction in the section on skyrmions later in this paper.

The determination of the GS of some frustrated spin systems as well as discussions on their properties are shown in the following.

II.2 Non collinear spin configurations

Consider as examples the plaquettes shown in Fig. 1. In the case of $XY$ spins, the GS configuration is obtained by minimizing the energy of the plaquette $E$. In the case of the triangular plaquette, suppose that spin $\mathbf{S}_{i}$ $(i=1,2,3)$ of amplitude $S$ makes an angle $\theta_{i}$ with the $\mathbf{Ox}$ axis. Writing $E$ and minimizing it with respect to the angles $\theta_{i}$, one has

[TABLE]

A solution of the last three equations is $\theta_{1}-\theta_{2}=\theta_{2}-\theta_{3}=\theta_{3}-\theta_{1}=2\pi/3$. One can also write

[TABLE]

The minimum corresponds to $\mathbf{S}_{1}+\mathbf{S}_{2}+\mathbf{S}_{3}=0$ which yields the $120^{\circ}$ structure. This is true also for the case of Heisenberg spin.

We can do the same calculation for the case of the frustrated square plaquette. Suppose that the antiferromagnetic bond connects the spins $\mathbf{S}_{1}$ and $\mathbf{S}_{2}$. We find

[TABLE]

If the antiferromagnetic bond is equal to $-\eta J$, the solution for the angles is Berge

[TABLE]

and $|\theta_{21}|=3|\theta|$, where $\cos\theta_{ij}\equiv\cos\theta_{i}-\cos\theta_{j}$. This solution exists if $|\cos\theta|\leq 1$, namely $\eta>\eta_{c}=1/3$. One can check that when $\eta=1$, one has $\theta=\pi/4$, $\theta_{21}=3\pi/4$.

We show the GS spin configurations of the frustrated triangular and square lattices in Fig. 2 with $XY$ spins ($N=2$).

At this stage, we note that the two GS found above have a two-fold degeneracy resulting from the equivalence of clockwise or counter-clockwise turning angle (noted by $+$ and $-$ in Fig. 3) between adjacent spins on a plaquette in Fig. 2. Therefore the symmetry of these plaquettes is of Ising type O(1), in addition to the symmetry SO(2) due to the invariance by global rotation of the spins in the plane.

It is expected from the GS symmetry of these systems that the transitions due to the respective breaking of O(1) and SO(2) symmetries, if they occur at different temperatures, belong respectively to the 2D Ising universality class and to the Kosterlitz-Thouless universality class Zinn . The question of whether the two phase transitions would occur at the same temperature and the nature of their universality remains at present an open question.

Another example is the case of helimagnets. Consider a chain of Heisenberg spins with ferromagnetic interaction $J_{1}(>0)$ between nn and antiferromagnetic interaction $J_{2}(<0)$ between nnn. When $\varepsilon=|J_{2}|/J_{1}$ is larger than a critical value $\varepsilon_{c}$, the spin configuration of the GS becomes non collinear. One shows that the helical configuration displayed in Fig. 4 is obtained by minimizing the interaction energy:

[TABLE]

where one has supposed that the angle between nn spins is $\theta$.

The two solutions are

[TABLE]

and

[TABLE]

The last solution is possible if $-1\leq\cos\theta\leq 1$, i.e. $J_{1}/\left(4|J_{2}|\right)\leq 1$ or $|J_{2}|/J_{1}\geq 1/4\equiv\varepsilon_{c}$.

This is shown in Fig. 4. There are two degenerate configurations corresponding to clockwise and counter-clockwise turning angles as in the previous examples.

Let us enumerate two frequently encountered frustrated spin systems where the nn interaction is antiferromagnetic: the fcc lattice and the hcp lattice. These two lattices are formed by stacking tetrahedra with four triangular faces. The frustration due to the lattice structure such as in these cases is called ”geometry frustration” DiepFSS .

III Exactly solved frustrated models

The 2D Ising model with non-crossing interactions is exactly soluble. Instead of finding the partition function one can map the model on a 16-vertex model or a 32-vertex model. The resulting vertex model will be exactly soluble. We have applied this method for finding the exact solution of several Ising frustrated models in 2D lattices with non-crossing interactions shown in Figs. 5-7.

Details have been given in Ref. Giacomini . I outline here only a simplified formulation of a model for illustration. The aim is to discuss the results. As we will see these models possess spectacular phenomena due to the frustration.

We take the case of the centered honeycomb lattice. The Hamiltonian of this model is as follows :

[TABLE]

where $\sigma_{i}=\pm 1$ is an Ising spin occupying the lattice site i , and the first, second, and third sums run over the spin pairs connected by heavy, light, and doubly light bonds, respectively (see Fig. 7). When $J_{2}=J_{3}=0$, one recovers the honeycomb lattice, and when $J_{1}=J_{2}=J_{3}$, one has the triangular lattice.

Let us denote the central spin in a lattice cell, shown in Fig. 7, by $\sigma$, and number the other spins from $\sigma_{1}$ to $\sigma_{6}$. The Boltzmann weight associated to the elementary cell is given by

[TABLE]

[TABLE]

The partition function is written as

[TABLE]

where the sum is performed over all spin configurations and the product is taken over all elementary cells of the lattice. Periodic boundary conditions are imposed. Since there is no crossing-bond interaction, the model is exactly soluble. To obtain the exact solution, we decimate the central spin of each elementary cell of the lattice. In doing so, we obtain a honeycomb Ising model with multispin interactions.

After decimation of each central spin, the Boltzmann factor associated to an elementary cell is given by

[TABLE]

[TABLE]

This model is equivalent to a special case of the 32-vertex model on the triangular lattice that satisfies the free-fermion condition as seen in the following.

Let us consider the dual lattice of the honeycomb lattice, i.e. the triangular lattice Bax . The sites of the dual lattice are placed at the center of each elementary cell and their bonds are perpendicular to bonds of the honeycomb lattice, as it is shown in Fig. 8.

Each site of the triangular lattice is surrounded by 6 sites of the honeycomb lattice. At each bond of the triangular lattice we associate an arrow. We take the arrow configuration shown in Fig. 9 as the standard one. We can establish a two-to-one correspondence between spin configurations of the honeycomb lattice and arrow configurations in the triangular lattice. This can be done in the following way : if the spins on either side of a bond of the triangular lattice are equal ( different ), place an arrow on the bond pointing in the same ( opposite ) way as the standard. If we do this for all bonds, then at each site of the triangular lattice there must be an even number of non-standard arrows on the six incident bonds, and hence an odd number of incoming ( and outgoing ) arrows. This is the property that characterizes the 32 vertex model on the triangular lattice.

In Fig. 10 we show two cases of the relation between arrow configurations on the triangular lattice and spin configurations on the honeycomb lattice.

In consequence, the Boltzmann weights of the 32-vertex model will be a function of the Boltzmann weights $W^{\prime}(\sigma_{1},\sigma_{2},\sigma_{3},\sigma_{4},\sigma_{5},\sigma_{6})$ , associated to a face of the honeycomb lattice. By using the relation between vertex and spin configurations described above and expression Eq. (9), we find

[TABLE]

Using the above expressions, the critical temperature of the model is determined from the equation (see details in Ref. Diep91a ):

[TABLE]

The solutions of this equation are given below for some special cases.

Following the case studied above, we can study other 2D Ising models without crossing bonds shown in Figs. 5-6: after decimation of the central spin in each square, these models can be mapped into a special case of the 16-vertex model which yields the exact solution for the critical surface (see details in Ref. Giacomini ).

Before showing some results in the space of interaction parameters, let us introduce the definitions of disorder line and reentrant phase.

III.1 Disorder line, reentrance

It is not the purpose of this review to enter technical details. I would rather like to describe the physical meaning of the disorder line and the reentrance. A full technical review has been given in Ref. Giacomini .

Disorder solutions are very useful for clarifying the phase diagrams of anisotropic models and also imply constraints on the analytical behavior of the partition function of these models.

A great variety of anisotropic models (with different coupling constants in the different directions of the lattice) are known to posses remarkable sub-manifolds in the space of parameters where the partition function is computable and takes a very simple form. These are the disorder solutions.

All the methods applied for obtaining these solutions rely on the same mechanism : a certain local decoupling of the degrees of freedom of the model results in an effective reduction of dimensionality for the lattice system. Such a property is provided by a simple local condition imposed on the Boltzmann weights of the elementary cell generating the lattice Ste . On a disorder line, a 2D system can behave as a 1D one: the dimension reduction is due to the decoupling of a degree of freedom in one direction, for instance.

This is very important while interpreting the system behavior: on one side of the disorder line, pre-ordering fluctuations have correlation different from those of the other side. Crossing the line, the system pre-ordering correlation changes. The dimension reduction is often necessary to realize this.

Disorder solutions have recently found interesting applications, for example in the problem of cellular automata (for a review see Rujan Rujan ). Moreover, they also serve to built a new kind of series expansion for lattice spin systems Mail .

Let us give now a definition for the reentrance. A reentrant phase lies between two ordered phases. For example, at low temperature ($T$) the system is in an ordered phase I. Increasing $T$, it undergoes a transition to a paramagnetic phase $R$, but if one increases further $T$, the system enters another ordered phase II before becoming disordered at higher $T$. Phase $R$ is thus between two ordered phases I and II. It is called ”reentrant paramagnetic phase” or ”reentrant phase”.

How physically is it possible? At a first sight, it cannot be possible because the entropy of an ordered phase is smaller than that of an disordered phase so that the disordered phase $R$ cannot exist at lower $T$ than the ordered phase II. In reality, as we will see below, phase II has always a partial disorder which compensates the loss of entropy while going from $R$ to II. The principle that entropy increases with $T$ is thus not violated.

III.2 Phase diagram

III.2.1 Kagomé lattice

A model of great interest is the Kagomé lattice shown in Fig. 5. The Kagomé Ising lattice with nn interaction $J_{\rm 1}$ has been solved a long time agoKa/Na showing no phase transition at finite $T$ when $J_{1}$ is antiferromagnetic. Taking into account the nnn interaction $J_{2}$, we have solved Aza87 this model by transforming it into a 16-vertex model which satisfies the free-fermion condition. The critical surface is given by

[TABLE]

For the whole phase diagram, the reader is referred to Ref. Aza87 . We show in Fig. 11 (bottom) only the small region of $J_{2}/J_{1}$ in the phase diagram which has the reentrant paramagnetic phase and a disorder line. This region lies around the phase boundary between two phases IV (partially disordered) and I (ferromagnetic) in Fig. 11 (top).

The phase X indicates a partially ordered phase where the central spins are free. The nature of ordering was determined by Monte Carlo (MC) simulations Aza87 .

Here again, the reentrant phase takes place between a low-$T$ ordered phase and a partially disordered phase. This suggests that a partial disorder in the high-$T$ phase is necessary to ensure that the entropy is larger than that of the reentrant phase.

When all the interactions are different in the model shown in Fig. 5, i.e. the horizontal bonds $J_{3}$, the vertical bonds $J_{2}$ and the diagonal ones are not equal, the phase diagram becomes complicated with new features Diep91b : in particular, we show that the reentrance can occur in an infinite region of phase space. In addition, there may be several reentrant phases occurring for a given set of interactions when $T$ varies.

The Hamiltonian is written as

[TABLE]

where $\sigma_{i}=\pm 1$ is an Ising spin occupying the lattice site i , and the first, second, and third sums run over the spin pairs connected by diagonal, vertical and horizontal bonds, respectively. When $J_{2}=0$ and $J_{1}=J_{3}$, one recovers the original nn Kagomé lattice Ka/Na . The effect of $J_{2}$ in the case $J_{1}=J_{3}$ has been shown above.

The phase diagram at temperature $T=0$ is shown in Fig. 12 in the space ($\alpha=J_{2}/J_{1}$, $\beta=J_{3}/J_{1}$) for positive $J_{1}$. The ground- state spin configurations are also displayed. The hatched regions indicate the three partially disordered phases (I, II, and III) where the central spins are free. Note that the phase diagram is mirror-symmetric with respect to the change of the sign of $J_{1}$. With negative $J_{1}$ , it suffices to reverse the central spin in the spin configuration shown in Fig. 12. Furthermore, the interchange of $J_{2}$ and $J_{3}$ leaves the system invariant, since it is equivalent to a $\pi/2$ rotation of the lattice. Let us consider the effect of the temperature on the phase diagram shown in Fig. 12. Partial disorder in the ground state often gives rise to the reentrance phenomenon as in systems shown above. Therefore, similar effects are to be expected in the present system. As it will be shown below, we find a new and richer behavior of the phase diagram: in particular, the reentrance region is found to be extended to infinity, unlike systems previously studied, and for some given set of interactions, there exist two disorder lines which divide the paramagnetic phase into regions of different kinds of fluctuations with a reentrant behavior.

Following the decimation method Giacomini , one obtains a checkerboard Ising model with multispin interactions. This resulting model is equivalent to a symmetric 16-vertex model which satisfies the free-fermion condition Gaff ; Suzu ; Wu72 . The critical temperature of the model is given by

[TABLE]

Note that Eq. (14) is invariant when changing $K_{1}\rightarrow-K_{1}$ and interchanging $K_{2}$ and $K_{3}$ as stated earlier. The phase diagram in the three-dimensional space ($K_{1},K_{2},K_{3}$) is rather complicated to show. Instead, we show in the following the phase diagram in the plane ($\beta=J_{3}/J_{1},T$) for typical values of $\alpha=J_{2}/J_{1}$. We just show now some interesting results in the interval $0>\alpha>-1$. In this range of $\alpha$, there are three critical lines. The critical line separating the F and P phases and the one separating the PD phase I from the P phase have a common horizontal asymptote as $\beta$ tends to infinity . They form a reentrant paramagnetic phase between the F phase and the PD phase I for positive b between a value $\beta_{2}$ and infinite $\beta$ (Fig. 23). Infinite region of reentrance like this has never been found before this model. As $\alpha$ decreases, $\beta_{2}$ tends to zero and the F phase is contracted. For $\alpha<-1$, the F phase disappears together with the reentrance.

In the interval $0>\alpha>-1$, the phase diagram possesses two disorder lines (see equations in Ref. Diep91b ). These two disorder lines are issued from a point near $\beta=-1$ for small negative $\alpha$; this point tends to zero as $\alpha$ tends to -1 (see Fig. 13).

III.2.2 Centered honeycomb lattice

To obtain the exact solution of our model, we decimate the central spin of each elementary cell of the lattice as outlined above. The resulting model is equivalent to a special case of the 32-vertex model Sacco on a triangular lattice that satisfies the free-fermion condition. The general treatment has been given in Ref. Diep91a . Here we take a particular case when $K_{2}=K_{3}$. The critical line obtained from Eq.(LABEL:eq30) is

[TABLE]

In the case $K_{2}=0$, the critical line is given by

[TABLE]

The case $K_{3}=0$ shows on the other hand a reentrant phase. The critical lines are determined from the equations

[TABLE]

[TABLE]

Fig. 14 shows the phase diagram obtained from Eqs. (17) and (18) around the phase boundary $\alpha=-0.5$. The reentrant paramagnetic phase goes down to zero temperature at the boundary $\alpha=-0.5$ separating GS phases II and III (see Fig. 14 right).

Note that the honeycomb model that we have studied here does not present a disorder solution with a dimensional reduction.

III.2.3 Periodically dilute centered square lattices

In this paragraph, we show the exact results on several periodically dilute centered square Ising lattices by transforming them into 8-vertex models of different vertex statistical weights that satisfy the free-fermion condition. The dilution is introduced by taking away a number of centered spins in a periodic manner. For a given set of interactions, there may be five transitions with decreasing temperature with two reentrant paramagnetic phases. These two phases extend to infinity in the space of interaction parameters. Moreover, two additional reentrant phases are found, each in a limited region of phase space Diep92 .

Let us consider several periodically dilute centered square lattices defined from the centered square lattice shown in Fig. 6.

The Hamiltonian of these models is given by

[TABLE]

where $\sigma_{i}=\pm 1$ is an Ising spin occupying the lattice site i , and the first, second and third sums run over the spin pairs connected by diagonal, vertical and horizontal bonds, respectively. All these models have at least one partially disordered phase in the ground state, caused by the competing interactions.

In each dilute square model shown in Fig. 6, the reentrance occurs along most of the critical lines when the temperature is switched on. This is a very special feature of the models which has not been found in other models.

Let us show in Fig. 15 the phase diagrams, at $T=0$, of the models shown in Figs. 6a, 6b and 6d, in the space ( $a,b$ ) where $a=J_{2}/J_{1}$ and $b=J_{3}/J_{1}$. The spin configurations in different phases are also displayed. The three-center case (Fig. 15a), has six phases (numbered from I to VI), five of which (I, II, IV, V and VI) are partially disordered (with, at least, one centered spin being free), while the two-center case (Fig. 15b) has five phases, three of which (I, IV, and V) are partially disordered. Finally, the one-center case has seven phases with three partially disordered ones (I, VI and VII).

Without showing the detailed calculation, let us describe in Fig. 16 the phase diagram of the three-center model (see Fig. 6a) in the space ( $a=J_{2}/J_{1}$, $T$) for typical values of $b=J_{3}/J_{1}$.

For $b<-1$, there are two reentrances. Fig. 16a shows the case of $b=-1.25$ where the nature of the ordering in each phase is indicated using the same numbers of corresponding ground state configurations (see Fig. 15). Note that all phases (I, II and VI) are partially disordered: the centered spins which are disordered at $T=0$ (Fig. 15a) remain so at all $T$. As seen, one paramagnetic reentrance is found in a small region of negative $a$ (schematically enlarged in the inset of Fig. 16a), and the other on the positive $a$ extending to infinity. The two critical lines in this region have a common horizontal asymptote.

For $-1<b<-0.5$, there are three reentrant paramagnetic regions as shown in Fig. 16b: the reentrant region on the negative $a$ is very narrow (inset), and the two on the positive $a$ become so narrower while $a$ goes to infinity that they cannot be seen on the scale of Fig. 16. Note that the critical lines in these regions have horizontal asymptotes. For a large value of $a$, one has five transitions with decreasing $T$: paramagnetic state - partially disordered phase I - reentrant paramagnetic phase - partially disordered phase II - reentrant paramagnetic phase- ferromagnetic phase (see Fig. 16b ). So far, this is the first model that exhibits such successive phase transitions with two reentrances.

For $-0.5<b<0$, there is an additional reentrance for $a<-1$: this is shown in the inset of Fig. 16c. As $b$ increases from negative values, the ferromagnetic region (III) in the phase diagram ”pushes” the two partially disordered phases (I and II) toward higher $T$.

At $b=0$, these two phases disappear at infinite $T$, leaving only the ferromagnetic phase.

For positive b, there are thus only two reentrances remaining on a negative region of $a$, with endpoints at $a=-2$ and $a=-1$, at $T=0$ (see Fig. 16d).

In conclusion, we summarize that in dilute square lattice models shown in Fig. 6a, we have found two reentrant phases occurring on the temperature scale at a given set of interaction parameters. A striking feature is the existence of a reentrant phase between two partially disordered phases which has not been found so far in any other model (we recall that in other models, a reentrant phase is found between an ordered phase and a partially disordered phase).

III.3 Summary and discussion

The present section shows spectacular phenomena due to the frustration. What to be retained is those phenomena occur around the boundary of two phases of different ground states, namely different symmetries. These phenomena include

1. the partial disorder at equilibrium: disorder is not equally shared on all particles as usually the case in unfrustrated systems,

2. the reentrance: this occurs around the phase boundary when $T$ increases $\rightarrow$ the phase with larger entropy will win at finite $T$. In other words, this is a kind of selection by entropy.

3. the disorder line: this line occurs in the paramagnetic phase. It separates the pre-ordering zones between two nearby ordered phases.

The partial disorder and the reentrance which occur in exactly solved Ising systems shown above are expected to occur also in models other than the Ising one as well as in some three-dimensional systems. Unfortunately, these systems cannot be exactly solved. One has to use approximations or numerical simulations to study them. This renders difficult the interpretation of the results. Nevertheless, in the light of what has been found in exactly solved systems, we can introduce the necessary ingredients into the model under study if we expect the same phenomenon to occur.

As seen above, the most important ingredient for a partial disorder and a reentrance to occur at low $T$ in the Ising model is the existence of a number of free spins in the ground state.

In three dimensions, apart from a particular exactly solved case Horiguchi showing a reentrance, a few Ising systems such as the fully frustrated simple cubic lattice Blan ; Diep85b , a stacked triangular Ising antiferromagnet Blan85 ; Nagai and a body-centered cubic (bcc) crystal Aza89b exhibit a partially disordered phase in the ground state. We believe that reentrance should also exist in the phase space of such systems though evidence is found numerically only for the bcc case Aza89b .

In two dimensions, a few non-Ising models show also evidence of a reentrance. For the $q$-state Potts model, evidence of a reentrance is found in a study of the two-dimensional frustrated Villain lattice (the so-called piled-up domino model) by a numerical transfer matrix calculation foster1 ; foster2 . It is noted that the reentrance occurs near the fully frustrated situation, i.e. $\alpha_{c}=J_{AF}/J_{F}=-1$ (equal antiferromagnetic and ferromagnetic bond strengths), for $q$ between $\simeq 1.0$ and $\simeq 4$. Note that there is no reentrance in the case $q=2$. Below (above) this $q$ value, the reentrance occurs above (below) the fully frustrated point $\alpha_{c}$. For $q$ larger than $\simeq 4$, the reentrance disappears foster2 .

The necessary condition for the occurrence of a partial disorder at finite $T$ is thus the existence of several kinds of site with different energies in the ground state. This has been so far verified in a number of systems as shown above.

IV Physics of thin films: surface magnetism, background

IV.1 Surface parameters

Surface physics has been rapidly developed in the last 30 years thanks to the progress in the fabrication and the characterization of films of very thin thickness down to a single atomic layer. A lot of industrial applications have been made in memory storage, magnetic censors, … using properties of thin films.

Theory and simulation have also been in parallel developed to understand these new properties and to predict further interesting effects. In the following I introduce some useful microscopic mechanisms which help understand macroscopic effects observed in experiments.

The existence of a surface on a crystal causes a lot of modifications at the microscopic levels. First, the lack of neighbors of atoms on the surface causes modifications in their electronic structure giving rise to modifications in electron orbital and atom magnetic moment by for example the spin-orbit coupling and in interaction parameters with neighboring atoms (exchange interaction, for example). In addition, surfaces can have impurities, defects (vacancies, islands, dislocations, …). In short, we expect that the surface parameters are not the same as the bulk ones. As a consequence, we expect physical properties at and near a surface are different from those in the bulk. For the fundamental theory of magnetism and its application to surface physics, the reader is referred to Ref. DiepTM .

In the following we outline some principal microscopic mechanisms which dominate properties of thin films.

IV.2 Surface spin-waves: simple examples

In magnetically ordered systems, spin-wave (SW) excitations dominate thermodynamic properties at low $T$. The presence of a surface modifies the SW spectrum. We show below that it gives rise to SW modes localized near the surface. These modes lie outside the bulk SW spectrum and modify the low-$T$ behavior of thin films.

Let us calculate these modes in some simple cases. We give below for pedagogical purpose some technical details.

We consider a thin film of $N_{T}$ layers stacked in the $z$ direction. The Hamiltonian is written as

[TABLE]

where $J_{ij}$ is the exchange interaction between to nn Heisenberg quantum spins, and $D_{ij}>0$ denotes an exchange anisotropy. $S_{i}^{+}$ and $S_{j}^{-}$ are the standard spin operators $S_{j}^{\pm}=S_{j}^{x}\pm iS_{j}^{y}$.

For simplicity, we suppose no defects and impurities at the surface and all interactions are identical for surface and bulk spins. The microscopic mechanism which governs thermodynamic properties of magnetic materials at low temperatures is the spin waves. The presence of a surface often causes spin-wave modes localized at and near the surface. These modes cause in turn a diminution of the surface magnetization and the magnetic transition temperature. There are several methods to calculate the spin-wave spectrum such as (see examples given in Ref. DiepTM ) the method of equation of motion, the Holstein-Primakoff method and the Green’s function method using a correlation function between two spin operators. Here we use for illustration the Green’s function method which the author has developed for thin films (see details in Ref. Diep1979 ; DiepTF91 ). This method shall be generalized below for helimagnets and other systems with non-collinear spin configurations.

We define the following double-time Green’s function

[TABLE]

The equation of motion of $G_{i,j}(t,t^{\prime})$ is written as

[TABLE]

where $[...]$ is the boson commutator and $\langle...\rangle$ the thermal average in the canonical ensemble defined as

[TABLE]

with $\beta=1/k_{B}T$. The commutator of the right-hand side of Eq. (22) generates functions of higher orders. In the first approximation, these functions can be reduced with the help of the Tyablikov decoupling Tyablikov as follows

[TABLE]

We obtain then the same kind of Green’s function defined in Eq. (21). As the system is translation-invariant in the $xy$ plane, we use the following Fourier transforms

[TABLE]

where $\omega$ is the SW (magnon) pulsation (frequency), $\vec{k}_{xy}$ the wave vector parallel to the surface, $\vec{R}_{i}$ the position of the spin at the site $i$, $n$ and $n^{\prime}$ are respectively the indices of the planes to which $i$ and $j$ belong ($n=1$ is the index of the surface). The integration on $\vec{k}_{xy}$ is performed within the first Brillouin zone in the $xy$ plane. Let $\Delta$ be the surface of that zone. Equation (22) becomes

[TABLE]

where the factors $(1-\delta_{n,1})$ and $(1-\delta_{n,N_{T}})$ are added to remove $C_{n}$ and $B_{n}$ terms for the first and the last layer. The coefficients $A_{n}$, $B_{n}$ and $C_{n}$ depend on the crystalline lattice of the film. We give here some examples:

• •

Film of simple cubic lattice

[TABLE]

where $C=4$ (in-plane coordination number) and $\gamma_{k}=\frac{1}{2}[\cos(k_{x}a)+\cos(k_{y}a)]$.

• •

Film of body-centered cubic lattice

[TABLE]

where $\gamma_{k}=\cos(k_{x}a/2)\cos(k_{y}a/2)$

Writing Eq. (26) for $n=1,2,...,N_{T}$, we obtain a system of $N_{T}$ equations which can be put in a matrix form

[TABLE]

where ${\bf u}$ is a column matrix whose n-th element is $2\delta_{n,n^{\prime}}<S_{n}^{z}>$.

For a given $\vec{k}_{xy}$ the magnon dispersion relation $\hbar\omega(\vec{k}_{xy})$ can be obtained by solving the secular equation $det|{\bf M}|=0$. There are $N_{T}$ eigenvalues $\hbar\omega_{i}$ ($i=1,...,N_{T}$) for each $\vec{k}_{xy}$. It is obvious that $\omega_{i}$ depends on all $\langle S_{n}^{z}\rangle$ contained in the coefficients $A_{n}$, $B_{n}$ and $C_{n}$.

To calculate the thermal average of the magnetization of the layer $n$ in the case where $S=\frac{1}{2}$, we use the following relation (see chapter 6 of Ref. DiepTM ):

[TABLE]

where $\langle S_{n}^{-}S_{n}^{+}\rangle$ is given by the following spectral theorem

[TABLE]

$\epsilon$ being an infinitesimal positive constant. Equation (34) becomes

[TABLE]

where the Green’s function $g_{n,n}$ is obtained by the solution of Eq. (33)

[TABLE]

$|{\bf M}|_{n}$ is the determinant obtained by replacing the n-th column of $|{\bf M}|$ by ${\bf u}$.

To simplify the notations we put $\hbar\omega_{i}=E_{i}$ and $\hbar\omega=E$ in the following. By expressing

[TABLE]

we see that $E_{i}$ ($i=1,...,N_{T}$) are the poles of the Green’s function. We can therefore rewrite $g_{n,n}$ as

[TABLE]

where $f_{n}(E_{i})$ is given by

[TABLE]

Replacing Eq. (39) in Eq. (36) and making use of the following identity

[TABLE]

we obtain

[TABLE]

where $n=1,...,N_{T}$.

As $<S_{n}^{z}>$ depends on the magnetizations of the neighboring layers via $E_{i}(i=1,...,N_{T})$, we should solve by iteration the equations (42) written for all layers, namely for $n=1,...,N_{T}$, to obtain the layer magnetizations at a given temperature $T$.

The critical temperature $T_{c}$ can be calculated in a self-consistent manner by iteration, letting all $<S_{n}^{z}>$ tend to zero.

Let us show in Fig. 17 two examples of SW spectrum, one without surface modes as in a simple cubic film and the other with surface localized modes as in body-centered cubic ferromagnetic case.

Note that a surface mode has a damping SW amplitude when going from the surface to the interior. The SW amplitudes for each mode are in fact their eigenvectors calculated from Eq. (37). It is very important to note that acoustic surface localized spin waves lie below the bulk frequencies so that these low-lying energies will give larger integrands to the integral on the right-hand side of Eq. (42), making $<S_{n}^{z}>$ to be smaller. The same effect explains the diminution of $T_{c}$ in thin films whenever low-lying surface spin waves exist in the spectrum.

Figure 18 shows the results of the layer magnetizations for the first two layers in the films considered above with $N_{T}=4$.

Calculations for antiferromagnetic thin films and other cases with non-collinear spin configurations can be performed using generalized Green’s functions Diep1979 ; DiepTF91 ; Quartu1998 with the general Hamiltonian defined for two spins $\mathbf{S}_{i}$ and $\mathbf{S}_{j}$ forming an angle $\cos\theta_{ij}$: one can express the Hamiltonian in the local coordinates as follows Diep2015

[TABLE]

The last term is an anisotropy added to facilitate a numerical convergence for ultra thin films at long-wave lengths since it is known that in 2D there is no ordering for isotropic Heisenberg spins at finite temperatures Mermin .

The determination of the angles in the ground state can be done either by minimizing the interaction energy with respect to interaction parameters NgoSurface ; NgoSurface2 . Using their values, one can follow the different steps presented above for the collinear magnetic films, one then obtains a matrix which can be numerically diagonalized to get the spin-wave spectrum which is used in turn to calculate physical properties in the same manner as for the collinear case presented above.

V Frustrated thin films: surface phase transition

Having given the background in the previous section, we can give some results here. The reader is referred to the original papers for lengthy technical details. Our aim here is to discuss physical effects due to the conditions of the surface.

As said earlier, the combination of the frustration and the surface effect gives rise to drastic effects as seen in the examples shown in the following.

V.1 Frustrated surfaces

We consider an example in this section: a ferromagnetic film with frustrated surfaces Ngo2007 . We study, by the analytical Green’s function method and extensive Monte Carlo simulations, effects of frustrated surfaces on the properties of thin films made of stacked triangular layers of atoms bearing Heisenberg spins with an Ising-like interaction anisotropy. We suppose that the in-plane surface interaction $J_{s}$ can be antiferromagnetic or ferromagnetic while all other interactions are ferromagnetic. We show that the ground-state spin configuration is non linear 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 related to disordering of surface and interior layers. There is a good agreement between Monte Carlo and Green’s function results.

V.2 Model

We consider a thin film made up by stacking $N_{z}$ planes of triangular lattice of $L\times L$ lattice sites.

The Hamiltonian is given by

[TABLE]

where $\mathbf{S}_{i}$ is the Heisenberg spin at the lattice site $i$, $\sum_{\left<i,j\right>}$ indicates the sum over the nearest neighbor spin pairs $\mathbf{S}_{i}$ and $\mathbf{S}_{j}$. The last term, which will be taken to be very small, is needed to ensure that there is a phase transition at a finite temperature for the film with a finite thickness when all exchange interactions $J_{i,j}$ are ferromagnetic. Otherwise, it is known that a strictly two-dimensional system with an isotropic non-Ising spin model (XY or Heisenberg model) does not have a long-range ordering at finite temperatures Mermin .

We suppose that the interaction between two nearest neighbors on the surface is equal to $J_{s}$, and all other interactions are ferromagnetic and equal to $J=1$ for simplicity. The two surfaces of the film are frustrated if $J_{s}$ is antiferromagnetic ($J_{s}<0$), due to the triangular lattice structure.

V.3 Ground state

In this paragraph, we suppose that the spins are classical. The classical ground state can be easily determined as shown below. Note that for antiferromagnetic systems, even for bulk materials, the quantum ground state though not far from the classical one, cannot be exactly determined because of the quantum fluctuations DiepTM .

For $J_{s}>0$ (ferromagnetic interaction), the magnetic ground state is ferromagnetic. However, when $J_{s}$ is negative there is a competition between the non collinear surface ordering and the ferromagnetic ordering of the spins of the beneath layer.

We first determine the ground state configuration for $I=I_{s}=0.1$ by using the steepest descent method : starting from a random spin configuration, we calculate the magnetic local field at each site and align the spin of the site in its local field. In doing so for all spins and repeating until the convergence is reached, we obtain in general the ground state configuration, without metastable states in the present model. The result shows that when $J_{s}$ is smaller than a critical value $J_{s}^{c}$ the magnetic ground state is obtained from the planar $120^{\circ}$ spin structure in the $XY$ plane, by pulling them out of the $xy$ plane by an angle $\beta$. The three spins on a triangle on the surface form thus an “umbrella” with an angle $\alpha$ between them and an angle $\beta$ between a surface spin and its beneath neighbor (see Fig. 19). This non planar structure is due to the interaction of the spins on the beneath layer, just like an external applied field in the $z$ direction. Of course, when $|J_{s}|$ is smaller than $|J_{s}^{c}|$ one has the collinear ferromagnetic ground state as expected: the frustration is not strong enough to resist the ferromagnetic interaction from the beneath layer.

We show in Fig. 20 $\cos(\alpha)$ and $\cos(\beta)$ as functions of $J_{s}$. The critical value $J_{s}^{c}$ is found between -0.18 and -0.19. This value can be calculated analytically by assuming the “umbrella structure”. For ground state analysis, it suffices to consider just a cell shown in Fig. 19. This is justified by the numerical determination discussed above. Furthermore, we consider as a single solution all configurations obtained from each other by any global spin rotation.

Let us consider the full Hamiltonian (44). For simplicity, the interaction inside the surface layer is set equal $J_{s}$ $(-1\leq J_{s}\leq 1)$ and all others are set equal to $J>0$. Also, we suppose that $I_{i,j}=I_{s}$ for spins on the surfaces with the same sign as $J_{s}$ and all other $I_{i,j}$ are equal to $I>0$ for the inside spins including interaction between a surface spin and a nn spin on the beneath layer.

The spins are numbered as in Fig. 19: $S_{1}$, $S_{2}$ and $S_{3}$ are the spins in the surface layer (first layer), $S^{\prime}_{1}$, $S^{\prime}_{2}$ and $S^{\prime}_{3}$ are the spins in the internal layer (second layer). The Hamiltonian for the cell is written as

[TABLE]

Let us decompose each spin into two components: an $xy$ component, which is a vector, and a $z$ component $\mathbf{S}_{i}=(\mathbf{S}_{i}^{\parallel},S_{i}^{z})$. Only surface spins have $xy$ vector components. The angle between these $xy$ components of nearest neighbor surface spins is $\gamma_{i,j}$ which is chosen by ($\gamma_{i,j}$ is in fact the projection of $\alpha$ defined above on the $xy$ plane)

[TABLE]

The angles $\beta_{i}$ and $\beta^{\prime}_{i}$ of the spin $\mathbf{S}_{i}$ and $\mathbf{S}^{\prime}_{i}$ with the $z$ axis are by symmetry

[TABLE]

The total energy of the cell (45), with $S_{i}=S^{\prime}_{i}=\frac{1}{2}$, can be rewritten as

[TABLE]

By a variational method, the minimum of the cell energy corresponds to

[TABLE]

We have

[TABLE]

For given values of $I_{s}$ and $I$, we see that the solution (49) exists for $J_{s}\leq J_{s}^{c}$ where the critical value $J_{s}^{c}$ is determined by $-1\leq\cos\beta\leq 1$. For $I=-I_{s}=0.1$, one obtains $J_{s}^{c}\approx-0.1889J$ in excellent agreement with the numerical result.

The classical ground state determined here will be used as input ground state configuration for quantum spins in the calculation by the Green’s method.

V.4 Results from the Green’s function method

Let us consider the quantum spin case. The details of the method in the case of non collinear spin configuration have been given in Ref. Ngo2007 . We just show the results on the surface phase transition and compare with the Monte Carlo results performed on the equivalent classical model.

V.4.1 Phase transition and phase diagram of the quantum case

We first show an example where $J_{s}=-0.5$ in Fig. 21. As seen, the surface-layer magnetization is much smaller than the second-layer one. In addition there is a strong spin contraction at $T=0$ for the surface layer. This is due to the antiferromagnetic nature of the in-plane surface interaction $J_{s}$ DiepTM . One sees that the surface becomes disordered at a temperature $T_{1}\simeq 0.2557$ while the second layer remains ordered up to $T_{2}\simeq 1.522$. Therefore, the system is partially disordered for temperatures between $T_{1}$ and $T_{2}$. This result is very interesting because it confirms again the existence of the partial disorder in quantum spin systems observed earlier in the bulk Rocco ; santa2 . Note that between $T_{1}$ and $T_{2}$, the ordering of the second layer acts as an external field on the first layer, inducing therefore a small value of its magnetization. A further evidence of the existence of the surface transition will be provided with the surface susceptibility in the Monte Carlo results shown below.

Figure 22 shows the non frustrated case where $J_{s}=0.5$, with $I=I_{s}=0.1$. As seen, the first-layer magnetization is smaller than the second-layer one. There is only one transition temperature. Note the difficulty for numerical convergency when the magnetizations come close to zero.

We show in Fig. 23 the phase diagram in the space $(J_{s},T)$. Phase I denotes the ordered phase with surface non collinear spin configuration, phase II indicates the collinear ordered state, and phase III is the paramagnetic phase. Note that the surface transition does not exist for $J_{s}\geq J_{s}^{c}$.

V.4.2 Monte Carlo results

The Green’s function method can go up to $T_{c}$ but due to the decoupling scheme, it cannot give a correct critical behavior at $T_{c}$. An alternative method for high temperatures is to consider the counterpart classical spins and to use Monte Carlo simulations to obtain the phase diagram for comparison. This is somewhat justified because the quantum nature of spins is no more important at high $T$.

For Monte Carlo simulations (see methods in Refs. Metropolis ; Binder ; Ferrenberg ; Ferrenberg2 ; DiepSP ), we use the same Hamiltonian (44) but the spins are the classical Heisenberg model of magnitude $S=1$. The film sizes are $L\times L\times N_{z}$ where $N_{z}=4$ is the number of layers (film thickness) taken as in the quantum case presented above. We use here $L=24,36,48,60$ to study finite-size effects. Periodic boundary conditions are used in the $xy$ planes. The equilibrating time is about $10^{6}$ Monte Carlo steps per spin and the averaging time is $2\times 10^{6}$ Monte Carlo steps per spin. $J=1$ is taken as unit of energy in the following.

Figure 24 shows the layer magnetizations of the first two layers as a function of $T$ , in the case $J_{s}=0.5$ (no frustration) with $N_{z}=4$ (the third and fourth layers are symmetric). In this case, there is clearly no surface transition just as in the quantum case.

Figure 25 shows a frustrated case where $J_{s}=-0.5$. The surface layer in this case becomes disordered at a temperature much lower than that for the second layer. Note that the surface magnetization is slightly smaller than 1 at $T=0$ (not seen with the scale of the figure). This is because the surface spins make an angle with the $z$ axis so their $z$ component is less than 1 in the ground state.

The phase diagram is shown in Fig. 26 in the space $(J_{s},T)$. This phase diagram resembles remarkably to that obtained for the quantum counterpart model shown in Fig. 23.

V.5 Frustrated thin films

We have also studied frustration effects in antiferromagnetic fcc Heisenberg films NgoSurface2 . In this case, the whole film is frustrated due to the geometry of the lattice.

Let us consider a film of fcc lattice structure with (001) surfaces. To avoid the absence of long-range order of isotropic non Ising spin model at finite $T$ when the film thickness is very small, i.e. quasi 2D system Mermin , we add in the Hamiltonian an Ising-like uniaxial anisotropy term. The Hamiltonian is given by

[TABLE]

where $\mathbf{S}_{i}$ is the Heisenberg spin at the lattice site $i$, $\sum_{\left<i,j\right>}$ indicates the sum over the nn spin pairs $\mathbf{S}_{i}$ and $\mathbf{S}_{j}$.

In the following, the interaction between two nn surface spins is denoted by $J_{s}$, while all other interactions are supposed to be antiferromagnetic and all equal to $J=-1$ for simplicity.

The ground state is shown (see demonstration in Ref. NgoSurface2 ) to depend on the surface in-plane interaction $J_{s}$ with a critical value $J^{c}_{s}=-0.5$ at which ordering of type I coexists with ordering of type II (see Fig. 27).

For $J_{s}<J^{c}_{s}$, the spins in each $yz$ plane are parallel while spins in adjacent $yz$ planes are antiparallel (Fig. 27a). This ordering will be called hereafter ”ordering of type I”: in the $x$ direction the ferromagnetic planes are antiferromagnetically coupled as shown in this figure. Of course, there is a degenerate configuration where the ferromagnetic planes are antiferromagnetically ordered in the $y$ direction. Note that the surface layer has an antiferromagnetic ordering for both configurations. The degeneracy of type I is therefore 4 including the reversal of all spins.

For $J_{s}>J^{c}_{s}$, the spins in each $xy$ plane is ferromagnetic. The adjacent $xy$ planes have an antiferromagnetic ordering in the $z$ direction perpendicular to the film surface. This will be called hereafter ”ordering of type II”. Note that the surface layer is then ferromagnetic (Fig. 27b). The degeneracy of type II is 2 due to the reversal of all spins.

Monte Carlo simulations have been used to study the phase transition in this frustrated film. We just show below three typical cases, at and far from $J^{c}_{s}$. Figure 28 shows the sublattice layer magnetizations at $J_{s}^{c}=-0.5$ where one sees that the surface layer undergoes a transition at a temperature lower than the interior ones. Far from this value there is a single phase transition as seen in Fig. 29. However, when $J_{s}$ is negatively stronger, we have a hard surface, namely the surface undergoes a phase transition at a $T$ higher than the interior layer, as seen in Fig. 30

The phase diagram is shown in Fig. 31

Note that near the phase boundary $J_{s}^{c}$ ($-0.5\leq J_{s}\leq-0.43$) a reentrant phase is found between phases I and II (not seen with the figure scale). As said in the 2D exactly solved models above, one has to be careful to examine a very small region near the phase boundary where unexpected phenomena can occur. This is the case here.

The nature of the phase transition is also studied by a histogram technique Ferrenberg ; Ferrenberg2 . Critical exponents are found to have values between 2D and 3D universality classes. The reader is referred to Ref. NgoSurface2 for details. We will return to this point in section VI below.

V.6 Helimagnetic films

Bulk helimagnets have been studied a long time ago Harada ; Diep1989a ; Diep1989b . A simple helimagnetic order resulting from the competition between the nn and nnn interactions is shown in section II.2. Helimagnetic films are seen therefore as frustrated films.

We have recently used the Green’s function method and Monte Carlo simulations to study helimagnetic films in zero field Diep2015 ; Diep2016 and in a perpendicular field DiepField2017 . We summarize here some results and emphasize their importance.

Consider the following helimagnetic Hamiltonian

[TABLE]

where $J_{i,j}$ is the interaction between two spins $\mathbf{S}_{i}$ and $\mathbf{S}_{j}$ occupying the lattice sites $i$ and $j$ and $\mathbf{H}$ denotes an external magnetic field applied along the $c$ axis. To generate helical angles in the $c$ direction, we suppose an antiferromagnetic interaction $J_{2}$ between nnn in the $c$ direction in addition to the ferromagnetic interaction $J_{1}$ between nn in all directions. For simplicity, we suppose that $J_{1}$ is the same everywhere. For this section we shall suppose $J_{2}$ is the same everywhere for the presentation clarity. Note that in the bulk in zero field, the helical angle along the $c$ axis is given by $\cos\alpha=-\frac{J_{1}}{4J_{2}}$ for a simple cubic lattice DiepTM with $|J_{2}|>0.25J_{1}$. Below this value, the ferromagnetic ordering is stable.

In zero field the helical angle has been shown to be strongly modified near the surface as shown in Fig. 32

Some results from the laborious Green’s function are shown in Fig. 33. Note the crossover of the layer magnetizations at low $T$. This is due to quantum fluctuations which are different for each layer, depending on the antiferromagnetic interaction strength (namely the so-called zero-point spin contractions, see Ref. DiepTM ). Without such a theoretical insight, it would be difficult to analyze experimental data when this happens.

In an applied field DiepField2017 , we have observed a new phenomenon, namely a partial phase transition in the helimagnetic film: contrary to what has been shown above (surface phase transition below or above the bulk one), here we have each single interior layer undergoes a separate transition. Synthetically, we can summarize that under the applied magnetic field, the spins in the GS make different angles between them across the film. When the temperature increases, the layers with large $xy$ spin-components undergo a phase transition where the transverse (in-plane) $xy$ ordering is destroyed. This ”transverse” transition is possible because the $xy$ spin-components do not depend on the field. Other layers with small $xy$ spin-components, namely large $z$ components, do not make a transition because the ordering in $S^{z}$ is maintained by the applied field. The transition of a number of layers with large $xy$ spin-components, not all layers, is a new phenomenon discovered here with our present model. Experiments have been performed on materials with helical structures often more complicated than the model considered in this paper. However, the clear physical pictures given in our present analysis are believed to be useful in the search for the interpretation of experimental data.

VI Criticality of thin films

One of the important fundamental questions in surface physics is the criticality of the phase transition in thin films.

To answer this question, we studied the critical behavior of magnetic thin films as a function of the film thickness DiepCritical . We used the ferromagnetic Ising model with the high-resolution multiple histogram Monte Carlo simulation Ferrenberg ; Ferrenberg2 . We showed that though the 2D behavior remains dominant at small thicknesses, there is a systematic continuous deviation of the critical exponents from their 2D values. We explain these deviations using the concept of ”effective” exponents suggested by Capehart and Fisher Fisher in a finite-size analysis. The shift of the critical temperature with the film thickness obtained here by Monte Carlo simulation is in an excellent agreement with their prediction.

We summarize here this work.

Let us consider the Ising spin model on a film made from a ferromagnetic simple cubic lattice. The size of the film is $L\times L\times N_{z}$. We apply the periodic boundary conditions (PBC) in the $xy$ planes to simulate an infinite $xy$ dimension. The $z$ direction is limited by the film thickness $N_{z}$. If $N_{z}=1$ then one has a 2D square lattice.

The Hamiltonian is given by

[TABLE]

where $\sigma_{i}$ is the Ising spin of magnitude 1 occupying the lattice site $i$, $\sum_{\left<i,j\right>}$ indicates the sum over the nn spin pairs $\sigma_{i}$ and $\sigma_{j}$.

Using the high-precision multi-histogram Monte Carlo technique Ferrenberg ; Ferrenberg2 we have calculated various critical exponents as functions of the film thickness using the finite-size scaling Barber described as follows. In Monte Carlo simulations, one calculates the averaged order parameter $\langle M\rangle$ ($M$: magnetization of the system), averaged total energy $\langle E\rangle$, specific heat $C_{v}$, susceptibility $\chi$, first order cumulant of the energy $C_{U}$, and $n^{th}$ order cumulant of the order parameter $V_{n}$ for $n=1$ and 2. These quantities are defined as

[TABLE]

Let us discuss the case where all dimensions can go to infinity. For example, consider a system of size $L^{d}$ where $d$ is the space dimension. For a finite $L$, the pseudo ”transition” temperatures can be identified by the maxima of $C_{v}$ and $\chi$, …. These maxima do not in general take place at the same temperature. Only at infinite $L$ that the pseudo ”transition” temperatures of these respective quantities coincide at the real transition temperature $T_{c}(\infty)$. So when we work at the maxima of $V_{n}$, $C_{v}$ and $\chi$, we are in fact working at temperatures away from $T_{c}(\infty)$. This is an important point to bear in mind for the discussion given below. Let us define the reduced temperature which measures the ”distance” from $T_{c}(\infty)$ by

[TABLE]

This distance tends to zero when all dimensions go to infinity. For large values of $L$, the following scaling relations are expected (see details in Ref. Bunker ):

[TABLE]

[TABLE]

and

[TABLE]

at their respective ’transition’ temperatures $T_{c}(L)$, and

[TABLE]

[TABLE]

and

[TABLE]

where $A$, $C_{0}$, $C_{1}$ and $C_{A}$ are constants. We estimate $\nu$ independently from $V_{1}^{\max}$ and $V_{2}^{\max}$. With this value we calculate $\gamma$ from $\chi^{\max}$ and $\alpha$ from $C_{v}^{\max}$. Note that we can estimate $T_{c}(\infty)$ using the last expression. Then, using $T_{c}(\infty)$, we can calculate $\beta$ from $M_{T_{c}(\infty)}$. The Rushbrooke scaling law $\alpha+2\beta+\gamma=2$ is then in principle verified Barber .

The results are shown in Table 1 where we observe a systematic deviation of the 2D critical exponents with increasing thickness.

An example of the calculation of $\nu$ is shown in Fig. 34 for $N_{z}=11$ to illustrate the precision of the method: the slope of the ”perfect” straight line of our data points gives $1\nu$.

Note that the PBC in the $z$ direction does not change the result if we do not apply the finite-size scaling in that direction DiepCritical .

We have also shown that by decreasing the film thickness, a first-order transition in a frustrated fcc Ising thin film can become a second-order transition DiepCrossover .

VII Skyrmions in thin films and superlattices

Skyrmions are topological excitations in a spin system. They result from the competition between different interactions in an applied magnetic field.

Skyrmions have been discovered by Skyrme Skyrme in the context of nuclear physics. Skyrmions have been shown to exist in condensed matter Bogdanov2003 ; Leonov ; Ackerman ; Ezawa ; Yu2 ; Yu1 ; Muhlbauer ; Bauer .

We consider in this section the case of a sheet of square lattice of size $N\times N$, occupied by Heisenberg spins interacting via a nn ferromagnetic exchange interaction $J$ and a nn Dzyaloshinskii-Moriya (DM) interaction Dzyaloshinskii ; Moriya . The Hamiltonian is given by

[TABLE]

where the $D$ vector of the DM interaction is chosen along the $\hat{z}$ direction perpendicular to the plane.

In zero field we have studied the spin waves and layer magnetizations at $T=0$ and at finite $T$ SahbiS . The results show that the DM interaction strongly affects the first mode of the spin-wave spectrum. Skyrmions appear only when an external field is applied perpendicular to the film, as seen in the following.

With $H\neq 0$, we minimize numerically the above Hamiltonian for a given pair ($H,D$), taking $J=1$ as unit, we obtain the GS configuration of the system. The phase diagram is shown in Fig. 35. Above the blue line is the field-induced ferromagnetic phase. Below the red line is the labyrinth phase with a mixing of skyrmions and rectangular domains. The skyrmion crystal phase is found in a narrow region between these two lines, down to infinitesimal $D$.

Let us show an example of the skyrmion crystal observed at ($D=1,H=0.5)$ and ($D=0.5,H=0.15)$ (Fig. 36 left). We see that the skyrmions form a crystal of triangular lattice. The size of each skyrmion depends on the ratio $H/D$.

Below the red line of Fig. 35 is the labyrinth phase (Fig. 36 right).

We wish to study the stability of the skyrmion crystal phase at finite $T$. To that end, we define an order parameter of the crystal as the projection of the actual spin configuration at the time $t$ at temperature $T$ on the GS configuration. We should average this projection over a large number of Monte Carlo steps per spin. The order parameter $M$ is

[TABLE]

where $\mathbf{S}_{i}(T,t)$ is the $i$-th spin at the time $t$, at temperature $T$, and $\mathbf{S}_{i}(T=0)$ is its state in the GS. The order parameter $M(T)$ is close to 1 at very low $T$ where each spin is only weakly deviated from its state in the GS. $M(T)$ is zero when every spin strongly fluctuates in the paramagnetic state.

We show in Fig. 37 the dependence of $M$ and $M_{z}$ on $T$ which indicate that for the skyrmion crystal remains ordered up to a finite $T$. This stability at finite $T$ may be important for transport applications.

We have studied finite-size effects on the phase transition at $T_{c}$ and we have seen that from $N=800$, all curves coincide: there is no observable finite size effects for $N\geq 800$.

We have studied the relaxation time of skyrmions and found that it follows a stretched exponential law SahbiS2 .

The DM interaction has been shown above to generate a skyrmion crystal in a 2D lattice. But skyrmions have been shown to exist in various kinds of lattices Rosales ; Buhrandt ; Kwon1 ; Kwon2 and in crystal liquids Bogdanov2003 ; Leonov ; Ackerman . Experimental observations of skyrmion lattices have been realized in MnSi in 2009 Muhlbauer ; Bauer and in doped semiconductors in 2010 Yu1 . Also, the existence of skyrmion crystals have been found in thin films Ezawa ; Yu2 and direct observation of the skyrmion Hall effect has been realized Jiang . In addition, artificial skyrmion lattices have been devised for room temperatures Gilbert .

It is noted that applications of skyrmions in spintronics have been largely discussed and their advantages compared to early magnetic devices such as magnetic bubbles have been pointed out in a recent review by W. Kang et al.WangKang . Among the most important applications of skyrmions, let us mention skyrmion-based racetrack memory Parkin , skyrmion-based logic gates Zhang ; ZhouEzawa , skyrmion-based transistor XZhang2015b ; JKim ; YShiota and skyrmion-based artificial synapse and neuron devices YHuang ; SLi .

Finally, we mention that we have found skyrmions confined at the interface of a superlattice composed alternately of a ferromagnetic film and an ferroelectric film Sharafullin ; Sharafullin2 . These results may have important applications.

VIII Conclusion

The aim of this review was to show a number of studied cases on the frustration and the surface effects in two dimensions and in thin films. Interesting phenomena have been found to occur near the frontier of two competing phases of different ground-state orderings. Without frustration, such frontiers do not exist. Exact solutions obtained for 2D Ising frustrated models show many striking features such as ”partial disorder”, namely a number of spins stay disordered in coexistence with ordered spins at equilibrium, ”reentrance”, namely a paramagnetic phase exists between two ordered phases in a small region of temperature, and ”disorder lines”, namely lines on which the system looses one dimension to allow for a symmetry change from one side to the other. Such beautiful phenomena can only be uncovered and understood by means of exact mathematical solutions. They allow us to understand qualitatively systems possessing similar microscopic ingredients but impossible to solve.

The surface effects have been studied by means of the Green’s function method for frustrated non-collinear spin systems. Monte Carlo simulations have also used to elucidate many physical phenomena where analytical methods cannot be used. Surface spin-waves, surface magnetization and surface phase transition have been analyzed as functions of interactions, temperature and applied field.

We have also treated the question of surface criticality. Results of our works show that critical exponents in thin films depend on the film thickness, their values lie between the values of 2D and 3D universality classes.

Recent results on skyrmions have also been reviewed in this paper. One of the striking findings is the discovery of a skyrmion crystal in a spin system with Dzyaloshinskii-Moriya interaction in competition with an exchange interaction, in a field. This skyrmion crystal is shown to be stable at finite temperature.

To conclude, we would like to say that investigations on the subjects discussed above continue intensively today. Note that there is an enormous number of investigations of other researchers on the above subjects and on other subjects concerning frustrated thin films. We mentioned these works in our original papers, but to keep the paper length reasonable we did not present them here.

Acknowledgments

The author wishes to thank his former doctorate students and collaborators for their close collaborations in the works presented in this review. In particular, he is grateful to Hector Giacomini, Patrick Azaria, Ngo Van Thanh, Sahbi El Hog, Aurélien Bailly-Reyre and Ildus F. Sharafullin who have greatly contributed by their works to the understanding of frustrated thin films.