Nematic Bond Theory of Frustrated Heisenberg Models on Triangular and Honeycomb Lattices

(1)

Nematic Bond Theory of

Frustrated Heisenberg Models on Triangular and Honeycomb Lattices

Cecilie Glittum

Thesis submitted for the degree of Master of Science

Department of Physics

Faculty of Mathematics and Natural Sciences UNIVERSITY OF OSLO

May 2020

(2)

(3)

Nematic Bond Theory of Frustrated Heisenberg Models on Triangular and

Honeycomb Lattices

Cecilie Glittum

(4)

Nematic Bond Theory of Frustrated Heisenberg Models on Triangular and Honeycomb Lattices

http://www.duo.uio.no/

Printed: Reprosentralen, University of Oslo

(5)

Abstract

The newly developed Nematic Bond Theory enables an efficient way of studying phases and phase transitions in classical frustrated magnets with Heisenberg interactions. We extend here the Nematic Bond Theory to hold for non-Bravais lattices, and develop in addition an approximate expression for the free energy, which can be calculated at low computational cost. Using these extensions, we study classical frustrated Heisenberg models on triangular and honeycomb lattices. On the triangular lattice, we derive the phase diagram for the J₁-J₂-J₃ Heisenberg model with ferromagnetic nearest neighbour interactions. In addition to detect lattice-nematic phases, we also find a novel symmetry broken state, which shows an extended maximum in the spin correlation function. We conjecture that this state is caused by domain wall excitations. The honeycomb lattice is studied to benchmark the Nematic Bond Theory on a non-Bravais lattice. For the antiferromagneticJ₁-J₂ Heisenberg model, we find that the Nematic Bond Theory shows good agreement with the literature for both symmetry broken lattice-nematic phases and symmetric spin liquid states.

i

(6)

(7)

Takk

Min fantastiske veileder, Olav Fredrik Syljuåsen, fortjener verdens største takk. Jeg setter utrolig stor pris på all hjelpen du har gitt meg og at du alltid har tatt deg tid til å prate med meg, til tross for at det til tider har opplevdes som en lyddiode. Tusen takk for at jeg har fått samarbeide så tett med deg de siste årene, og virkelig fått forske sammen med deg. Jeg kunne ikke ha drømt om en bedre veileder.

Takk til teorigruppa med ansatte og studenter for å inkludere meg, selv under flyt- tingen til Ullevål. Jeg kommer til å savne lunsjene med dere, både i virkelig og virtuell utgave. Takk også til alle ansatte på Fysisk institutt for å bidra til å gjøre fysikkutdan- ningen så god som den er. Tiden min her har vært helt uerstattelig, og det blir med et tungt hjerte at jeg forlater mitt andre hjem.

Takk Anne, Elin, Marianne og resten av realistene på Bamble videregående skole for å gi meg tro på meg selv og selvtillit nok til å være den nerden jeg er. En spesiell takk til Anne for alle gode vibber gjennom studietiden og for å alltid minne meg på at godt nok er godt nok.

Hallooo Wanja. Jeg er så glad for at en underlig fascinasjon resulterte i vårt vennskap.

Takk for at du har vært der for meg de siste tre årene, både gjennom oppturer, nedturer og shoppingturer. Du er min beste venn og jeg kommer til å savne deg.

Jeg hadde ikke hatt mye til sosialt liv, hadde det ikke vært for Metin, Ivar, Jonas, Ylva og Lasse. Takk for alle middager, turer og samtaler vi har hatt de siste årene. En ekstra takk til Ivar for våre gode samarbeid. Jeg gleder meg til jeg kan se dere alle igjen.

I want to thank Taylor Swift for creating marvellous tunes and for inspiring the colour scheme used in this thesis.

Til slutt vil jeg takke Torbjørn. Takk for at du alltid viser interesse for det jeg driver med, og for at du gledelig hører på alt jeg har å si. Du betyr alt.

Oslo, 15.05.20 Cecilie Glittum

iii

(8)

(9)

Chapter 1 Introduction

Magnetic materials are made up of microscopic magnetic moments (spins). In ferro- magnets, the spins align at low temperatures, causing a macroscopic magnetic moment.

In another, more complex type of magnet, the antiferromagnet, spins anti-align at low temperatures with no macroscopic magnetic moment. Despite their different microscopic and macroscopic properties, ferro- and antiferromagnets both have in common that the ground state is uniquely determined by minimizing all pairwise interactions between the spins.

In a frustrated magnet, the spins are not able to minimize all interactions simultaneously due to the geometry of the lattice and/or competing interactions. This frustration gives rise to a large degeneracy of the system’s ground state, making frustrated magnets particularly interesting, as it may lead to emergence of intriguing phases. Examples of such are lattice-nematic phases [1, 2], which break the symmetries of the underlying lattice, and spin liquids [3], being symmetric phases with strong correlations. Frustrated magnetism is also believed to play an important role in the understanding of various unconventional superconductors [4, 5, 6].

Little is known about frustrated magnets, even in the classical limit, and there is thus a fundamental interest in studying such systems. As ordinary expansion techniques are only trustworthy when studying systems with a uniquely defined ground state, frustrated magnets are usually studied using Monte Carlo simulations. This is however time-consuming and demands high computational power.

Recently, Schecter, Syljuåsen and Paaske developed a field-theoretic framework, the Nematic Bond Theory (NBT) [7]. This framework enables a more efficient, but approximate, way of calculating critical temperatures and correlation functions for classical frustrated magnets with Heisenberg interactions.

In the literature the NBT is developed for Bravais lattices. It has been applied to the square and cubic lattices [7, 8] and to systems with complicated power-law interactions [9], showing great success in detecting both lattice-nematic phases and phases breaking the global spin symmetry. However, the NBT is only approximate and has been found to systematically overestimate transition temperatures compared to Monte Carlo

1

(14)

simulations [8]. In addition, there are interpretational issues with NBT as multivaluedness of order parameters close to first order phase transitions and lack of a well-defined way to distinguish between stable and metastable states [8].

The next natural lattice to apply the NBT to after the square lattice is the triangular lattice, as it is the quintessential example of a geometrically frustrated antiferromagnet [10]. Furthermore, frustrated magnetism has also shown interesting features on so-called non-Bravais lattices. The simplest non-Bravais lattice based on the triangular lattice is the honeycomb lattice. The honeycomb lattice is in fact a bipartite lattice, and is thus unfrustrated considering nearest neighbour interactions. Interests in frustrated Heisenberg models on the honeycomb lattice have however increased recently [11, 12, 13]

after experiments showed that the S = 3/2 material Bi3Mn4O12(NO3) does not order down to a low temperature [14], pointing in the direction of competing interactions.

Even though the honeycomb lattice is not geometrically frustrated, other non-Bravais lattices have the potential of exerting an even higher degree of geometrical frustration than the triangular lattice. Examples of such are the highly frustrated kagome and pyrochlore lattices [15]. However, to study non-Bravais lattices using the NBT, an extension taking the sublattice nature of the non-Bravais lattices into account would be required.

This thesis has two main goals, the first being to further develop the NBT. The development will include an extension to non-Bravais lattices, as well as attempting to solve the multivaluedness problem by deriving an expression for the free energy. The second goal is to apply the NBT to the triangular and honeycomb lattices. In particular, we will study the J1-J2-J3 Heisenberg model on the triangular lattice with a goal of both demonstrating the strength of the NBT and searching for novel phases by scanning the full range of J2 and J3 couplings. For the honeycomb lattice, we will study the antiferromagnetic J1-J2 Heisenberg model and compare our results with Monte Carlo studies [11, 12] to benchmark the NBT-extension to non-Bravais lattices.

We begin the thesis with an overview of general and necessary background theory in Chapter 2. In Chapter 3 we extend the NBT to non-Bravais lattices and derive an expression for the free energy. The numerical implementation of the NBT is outlined in Chapter 4. Furthermore, in Chapter 5, we outline the relevant details of the lattices and interactions which we study. The numerical results obtained using the NBT are presented in Chapter 6. In Chapter 7 we discuss both the results of the numerical calculations and the developments of the framework. Lastly, we conclude in Chapter 8.

(15)

Chapter 2 Background

The purpose of this chapter is to introduce relevant background concepts. We start by a short description of the origin of magnetic interactions and frustrated magnetism. Then we describe details of the Heisenberg model, before going into details on phase transitions and the low-temperature phases expected when studying frustrated Heisenberg models.

2.1 The Origin of Magnetic Interactions

The classical dipolar interactions between the electronic moments are of order 10⁻⁵ eV, and are thus too weak to explain the observed magnetic transition temperatures which in many materials are of order 10²-10³ K. The much stronger magnetic interactions are caused by a coupling mechanism derived from the following fundamental properties of electrons [16]:

1. The electron’s spin.

2. The electron’s kinetic energy.

3. Pauli’s exclusion principle.

4. Coulomb repulsion between electrons.

Magnetism is thus truly a quantum phenomenon.

The spin-dependent part of the interaction between two electrons can be modelled by a Heisenberg interaction [16]

H =J ~S₁·S~₂. (2.1)

Here, S~₁ and S~₂ are the spin operators for each of the two electrons and J is referred to as the exchange coupling between the two spins. For J <0, alignment of the spins is en- ergetically favoured. This is referred to this as a ferromagnetic interaction. J >0favours anti-alignment of the spins and is analogously referred to as an antiferromagnetic interaction. Whether the interaction is ferromagnetic or antiferromagnetic depends on detailed

3

(16)

1 2

3 ?

Figure 2.1: Three spins interacting antiferromagnetically arranged in a triangle. If spins 1 and 2 are anti-aligned, spin 3 cannot be simultaneously anti-aligned with both spin 1 and spin 2.

atomic mechanisms occurring in the material. For more details on such mechanisms, see for example Ref. [16].

2.2 Frustrated Magnetism

Frustration in physics in general means that a system cannot satisfy all interactions simultaneously. In magnetic systems, frustration can arise as a consequence of the geometry of the lattice and/or competing interactions.

The simplest example of a geometrically frustrated system is three antiferromagnetically interacting spins arranged in a triangle as shown in Fig. 2.1. Once two of the spins are anti-aligned, the third can only satisfy the interaction with one of the other two.

Hence, antiferromagnetic order is incompatible with the triangle.

Other magnets are frustrated due to competing interactions. A spin in a magnetic material may have couplings to its nearest neighbours, next-nearest neighbours and so forth. These interactions can be competing. Consider for instance a spin chain with antiferomagnetic nearest and next nearest neighbour interactions. Both these interactions cannot be satisfied simultaneously, i.e. they are competing interactions leading to frustration.

2.3 The Heisenberg Model on Non-Bravais Lattices

Crystalline structures are divided into two types of lattices: Bravais lattices and non- Bravais lattices. In a Bravais lattice all lattice points are equivalent, and the whole lattice is spanned by a set of lattice vectors. A non-Bravais lattice can be described as a set of interpenetrating Bravais lattices referred to as sublattices. This can also be viewed as a Bravais lattice with an additional basis. In a non-Bravais lattice, all lattice points are not equivalent, and the position of an atom has to be specified both by the unit cell of the underlying Bravais lattice and the position within the unit cell.

(17)

2.4. Phase Transitions 5 Generalizing to a non-Bravais lattice, the Heisenberg Hamiltonian takes the form

H = 1 2

X

R, ~~R⁰

X

ij

JR ~~R⁰,ijS~R,i~ ·S~R~⁰,j, (2.2)

whereS~R,i~ =

SˆR,i,x~ ,SˆR,i,y~ ,SˆR,i,z~

is the quantum mechanical spin operator for a spin on sublattice i in the unit cell positioned at R~. J_{R ~}_~_R0,ij is the exchange coupling between a spin at sublatticei in the unit cell positioned at R~ and a spin at sublattice j in the unit cell positioned at R~⁰.

By the definition of a non-Bravais lattice, the environments of spins on different sublattices are different, and thus the interactions are different for spins on different sublattices. Consequently, the exchange interaction is in general anm×m matrix, where m is the number of sublattices.

2.3.1 Classical spins

Despite magnetism being a quantum mechanical phenomenon, we will in this thesis consider the classical limit. For large spins it is a reasonable approximation to substitute the quantum mechanical spin operators with classical vectors of fixed length. In the case of the Heisenberg model, the spins are then represented as points on the sphere. In other words, the spins are modelled as classical vectors with a unit-length constraint.

A softer constraint would be to require the sum of all the spins on a sublattice to sum up to the volume V, i.e.

1 V

X

R~

S~R,i~ ·S~R,i~ = 1, (2.3) In reciprocal space it can be written

1 V

X

~ q

S~_~_q,i·S~_~_q,i = 1. (2.4)

2.4 Phase Transitions

A phase transition is a point in parameter space at which a system goes from one phase to another, i.e. several physical parameters change in a non-analytic way. Such a transition can be thermal, meaning that the transition occurs upon varying the temperature, or it can be quantum, meaning that the transition occurs at zero temperature upon varying a non-thermal parameter. We will here consider thermal phase transitions.

The temperature at which the phase transition occurs is called the transition temperature, which we denoteT_c. The phase transition is characterized by an order parameter, which is non-zero in the ordered phase (T < T_c) and zero in the disordered phase (T > T_c). It is usually taken to be the thermal average of a first derivative of the free

(18)

O

T /T_c 1

1

F

T /T_c 1

Figure 2.2: Illustration of the order parameter (left) and the free energy (right) near a first order phase transition. The order parameter makes a discontinuous jump at the transition temperature T_c. The free energy is continuous, but has a discontinuous slope at T_c. The dashed line in the free energy shows an overshoot of metastable states, characteristic for first order transitions.

energy. An example of an order parameter is the net magnetization, characterizing the ferromagnetic-paramagnetic transition.

There are two classes of phase transitions: first order and continuous. A phase transition is said to be first order if of the free energy has discontinuities in its first derivatives at the transition temperature T_c. Consequently, the order parameter changes discontinuously at the phase transition. The typical shape of an order parameter for a first order transition is illustrated in Fig. 2.2 left panel.

The entropy

S =− ∂F

∂T

V

(2.5) is also a first derivative of the free energy and changes discontinuously at a first order phase transition. This gives rise to a latent heat and also a discontinuity in the slope of the free energy at the transition temperature. An illustration of the typical shape of the free energy as a function of temperature for a first order transition is shown in Fig. 2.2 right panel. The dashed line shows metastable states, which may occur close to a first order phase transition. These are metastable as there exist other states at the same temperature with a lower free energy.

A continuous phase transition is characterized by continuous first derivatives and singular second derivatives of the free energy. Consequently, the order parameter and the entropy are continuous at all temperatures, including the transition temperature.

Examples of the order parameter and the free energy of a continuous phase transition are illustrated in Fig. 2.3. The order parameter approaches zero continuously as the temperature approaches the transition temperature. The free energy has a continuous slope, and does not show signs of metastable states.

(19)

2.4. Phase Transitions 7

O

T /T_c 1

1

F

T /T_c 1

Figure 2.3: Illustration of the order parameter (left) and the free energy (right) near a continuous phase transition. The order parameter goes continuously to zero as the temperature increases towards the transition temperature T_c. Both the free energy and its slope are continuous at T_c.

2.4.1 Spontaneous Symmetry Breaking and the Mermin-Wagner The- orem

If the Hamiltonian is symmetric under a symmetry operation, while the state of the system is not, the considered symmetry is said to be spontaneously broken. Such spontaneous symmetry breaking implies a phase transition.

Conventional phase transitions in magnetic systems break the global spin symmetry.

One example of such is the ferromagnetic-paramagnetic phase transition. In the paramagnetic state, the spin configuration is symmetric under global rotations of the spins. In the ferromagnetic phase the spins on the other hand align in one specific direction. Thus, the ferromagnetic phase has broken the global spin symmetry.

As the magnetization can point anywhere on the sphere for the Heisenberg model, the global spin symmetry is a continuous symmetry. It was proven by Mermin and Wagner in 1966 that the spontaneous breaking of continuous symmetries is precluded for short-ranged interactions on lattices of dimension less than three [17]. In the cases of one- and two-dimensional lattices there are thus no finite-temperature phase transitions breaking continuous symmetries if the systems have short-ranged interactions. This is a consequence of fluctuations being enhanced in low dimensions, where the number of spin waves generated at finite temperatures are diverging [16, 18].

However, discrete symmetries of the Hamiltonian, such as lattice symmetries, are not affected by the Mermin-Wagner theorem. The breaking of lattice symmetries may thus be possible at finite temperature also in low dimensions, leading to lattice-nematic phases.

We will discuss lattice-nematic phases further in Section 2.5.1.

(20)

2.5 Spiral States

For Bravais lattices with a Heisenberg Hamiltonian, the states of lowest energy are in general known to be helimagnetic [19, 20], i.e. spiral states given by¹

S~_~_r(Q) =~ ~ucos Q~ ·~r

+~vsin Q~ ·~r

, (2.6)

where~uand ~v are orthogonal vectors of unit length in R³. These vectors can be rotated continuously, and preserve the global spin symmetry. Q~ is one of the absolute minima of the exchange coupling in reciprocal space J_~_q. As the spin configuration in Eq. (2.6) is specified by one wave vector, we refer to this as a single-~q spiral state.

J_~_q has in general several absolute minima. Two minima are said to be equivalent if they are related by a reciprocal lattice vector. The set of non-equivalent minima of J_~_q constitute what is referred to as the star ofQ~ [1]. The star ofQ~ can be either a discrete or continuous set.

If the star ofQ~ has a discrete set of elements and the corresponding spin configurations S~_~_r are different under global spin rotations, the single-~q ground state is said to have a discrete degeneracy. We then have the possibility of breaking a discrete symmetry, allowing for a finite-temperature phase transition also in low-dimensional systems. There are two cases in which the single-~q ground state does not possess a discrete degeneracy, despite the star of Q~ being a discrete set:

• If the star of Q~ has only one element, all single-~q ground states can be reached by global spin rotations. Thus, the ground state does not have a discrete degeneracy and we do not expect a phase transition². Such cases include the ferromagnet.

• If the star of Q~ is reduced to Q~ and −Q~, all single-~q ground states can be reached by global spin rotations. Inverting Q~ is equivalent to rotating~v an angle π about

~

u. Thus, the spin configurations S~_~_r(Q)~ and S~_~_r(−Q)~ are equal up to a global spin rotation, and we do not expect a phase transition.

In the case where the star ofQ~ is a continuous set, the single-~qground state has a continuous degeneracy, and it should in principle not be possible to break any symmetries.

However, the continuous degeneracy of the ground state may be lifted due to thermal fluctuations, a principle known as order by disorder [2]. Low-temperature thermal excitations in magnetic systems are spin waves. Different spiral states may have different spin wave dispersions, and consequently different spin wave entropies.

As the free energy is

F =E−T S, (2.7)

1There are some exceptions, which we will not consider in this thesis. Some of the cases where there exist additional ground states are presented in Ref. [1].

2An exception is the celebrated Kostelitz-Thouless transition occuring in the XY-model [21].

(21)

2.5. Spiral States 9

Figure 2.4: An example of three different spiral states on the triangular lattice where the wave vectors are related by lattice symmetries. The nearest neighbouring bonds with strongest correlations are drawn as blue lines.

the continuous degeneracy of the ground state may be lifted at finite temperature by differences in entropy. It is then the ground states maximizing the entropy that will have the lowest free energy. If thermal fluctuations lift the continuous degeneracy of the ground state such that the degeneracy becomes discrete, it should theoretically be possible to have a phase transition breaking a discrete symmetry also when the star ofQ~ is a continuous set.

For the rest of this thesis, we will refer toQ~ and −Q~ as the same element of the star of Q~. Thus, we will also refer to peaks in the correlation function at Q~ and −Q~ as one peak.

2.5.1 Lattice-Nematic Phases

Discrete symmetries of the ground state are in general point group symmetries of the lattice, i.e. rotational symmetries and/or mirror symmetries. If the rotational symmetry of the lattice is broken spontaneously, we get so-called lattice-nematic phases. In the same way nematics in liquid crystals break rotational symmetry by aligning in one direction, lattice-nematic phases break rotational symmetry by having a direction of strongest correlations.

Three different spin configurations on the triangular lattice are shown in Fig. 2.4.

These configurations are related by lattice symmetries. For each spin configuration, blue lines are drawn to mark the directions of strongest nearest neighbour correlations.

For high temperatures, we have a thermal mixture of all spirals, and the system has all symmetries preserved. Upon lowering the temperature, we may encounter a phase transition spontaneously breaking the lattice symmetries. In the low temperature phase, the system may be governed by one of the states shown in Fig. 2.4. In such a case, the correlations are clearly stronger in one specific direction, i.e. the rotational symmetry of the lattice is broken. It may look like also the global spin symmetry is broken, but indeed the vectors ~u and~v are free, preserving the global spin symmetry.

For such transitions, breaking the lattice symmetries, the transition temperature is typically of the order of the free energy barriers between the different ground states. If the temperature is low, the system gets locked into one of the ground state configurations,

(22)

as thermal fluctuations are not large enough to cross the energy barriers. For higher temperatures, the energy barriers can be crossed, and the system fluctuates between the distinct ground state configurations.

2.5.2 Entropy

The spin wave entropy per spin is given by s=−1

N X

~k

lnω(~k), (2.8)

whereω(~k)is the spin wave dispersion. In Ref. [22], Seabra et al. derive an expression for the entropy of a spiral state. This expression is obtained by expanding around a planar spiral state characterized by a momentum vector Q, given by Eq. (2.6). They find that~ the spin wave dispersion for a spiral state is given by

ω(~k) = r

β²1 2

h

J(Q~ +~k) +J(Q~ −~k)−2J(Q)~ i h

J(~k)−J(Q)~ i

. (2.9)

(23)

Chapter 3 The Nematic Bond Theory

The NBT [7] is a field-theoretic framework to study phases and phase transitions of classical Heisenberg models with N spin components. The partition function consists of integrals over all spins, where the spins are constrained to have unit length. These integrals are impossible to treat exactly. In the NBT, the unit length constraint on the classical spins is imposed in form of delta-functions with a constraint field. When calculating the partition function, the integrals over the spins can then be treated as individual Gaussian integrals. The partition function is then left as integrals over the constraint field.

The Self-Consistent Gaussian Approximation (SCGA) is a method similar to the NBT.

In the SCGA, the constraint field is set to be spatially homogeneous, and the unit length constraint is then ensured by adjusting the value of the constraint field [23]. The SCGA gives good results for the spin correlation function in the symmetric phases. However, it is unable to predict phase transitions with broken lattice symmetries as the spin correlation function inherits all the symmetries of the exchange couplings.

The NBT differs from the SCGA in the treatment of the constraint field by allowing for a spatial dependence. In the NBT, the constraint field is separated into two parts:

the spatially homogeneous part and the spatial fluctuations from the former. The unit length constraint is enforced, similarly to the SCGA, by adjusting the homogeneous part of the constraint field. However, by a careful treatment of the spatial fluctuations in the constraint field, the NBT is able to capture also symmetry breaking phases.

The spatial fluctuations are treated by a diagrammatic perturbation theory in an1/N expansion. This leads to a momentum-dependent self-energy, which enables symmetry breaking. The aim of the NBT is to compute this self-energy and the spin correlation function through a set of coupled self-consistent equations. These equations can be solved numerically by iteration to give information about phases and phase transitions.

We will in the following extend the NBT to hold also for non-Bravais lattices, following closely the approach of Schecter et al. [7]. In addition, we derive an approximate expression for the free energy.

11

(24)

3.1 Notation

In extending the NBT to non-Bravais lattices it will be convenient to write observables as matrices in both momentum space and sublattice space. These matrices have four indices: Two momentum indices and two sublattice indices. We denote such matrices as bold symbols. An example we will encounter is the exchange coupling matrix J. An element of this matrix is specified as

J~q~q⁰,ij =J~q,ijδ~q~q⁰. (3.1) We will also refer to matrices that are matrices only in the sublattice space. These matrices willnot be in bold. An example of such a matrix is shown in Eq. (3.1), J_~_q. An element of this matrix is specified by two sublattice indices, i.e. J_~_q,ij.

For dot products, the convention we use is to complex conjugate the first vector, i.e.

A~·B~ =X

α

A^α∗B^α. (3.2)

We will in general setk_B =~=a= 1, wherea is the lattice spacing. Thus entropies, spins, lengths and momenta are unit-less. Also, β = 1/T.

3.1.1 Fourier Conventions

We specify here the Fourier conventions used throughout this thesis. For the spins, we will use

S~R,i~ = 1

√V X

~ q

S~~q,ie^i~^q·^R^~, S~~q,i= 1

√V X

R~

S~R,i~ e^−i~^q·^R^~. (3.3) Analogously for the constraint field, we will use

~λR,i~ = 1

√V X

~ q

~λ_~_q,ie^i~^q·^R^~, ~λ_~_q,i = 1

√V X

R~

~λR,i~ e^−i~^q·^R^~. (3.4) When it comes to the exchange coupling, we have chosen a rather different convention in order to get the Hamiltonian in momentum space on a simpler form. The chosen convention is

J_R,ij_~ = 2 V

X

~ q

J_~_q,ije^−i~^q·^R^~, J_~_q,ij = 1 2

X

R~

J_R,ij_~ e^i~^q·^R^~. (3.5)

3.2 The Heisenberg Model and the Unit Length Constraint

As described in Section 2.3, the classical Heisenberg Hamiltonian for non-Bravais lattices reads

H = 1 2

X

R ~~R⁰

X

ij

J_R_~0−R,ij~ S~_R,i_~ ·S~_R_~0,j, (3.6)

(25)

3.2. The Heisenberg Model and the Unit Length Constraint 13 where we have taken the exchange coupling to depend on the relative separation of the spins. The corresponding partition function reads

Z = Z

D ~S e^−βH. (3.7)

As the spins are classical, the integral over each spin is constrained to the unit sphere.

These constraints can be enforced throughδ-functions δ

|S~_R,i_~ | −1

=

Z βdλR,i~

π e^−iβλ^R,i^~

S~R,i~ ·S~R,i~ −1

, (3.8)

where the integration variable λ_R,i_~ enforces unit length of the spin S~_R,i_~ . We refer to λ as the constraint field. The constraint field has a spatial dependency, as opposed to the SCGA, and is also dependent on the sublattice index. When enforcing the unit length constraints by δ-functions, the partition function can be written

Z = Z

D ~SDλe^−βHe^−iβ

P

R,i~ λ_R,i_~

S~_R,i_~ ·S~_R,i_~ −1

, (3.9)

where the constant factors in front of the partition function have been absorbed into the integration measure. These constants will not be of importance as long as we only consider expectation values. However, when studying the free energy, they will play an important role. The constant factors are thus retrieved in Section 3.6.1.

As the constraint field enforces unit length of the spins, we can now integrate over each spin component separately over the range (−∞,∞). The action is quadratic in the spin components, which implies that we can carry out the integrals over the spin components as independent Gaussian integrals. This is easiest seen by rewriting the partition function as

Z = Z

D ~SDλe^−β

P

~ q

J~q,ijδ_~_q,~_q0+i^√¹

Vλ_~_q−~_q0,iδij

S~~q,i·S~_~_q0,j+iβ√ VP

iλ~q=0,i

, (3.10)

where we in addition have transformed to reciprocal space. By scaling the spin components by 1/√

β and performing the Gaussian integrals over the spin components, we obtain

Z = Z

Dλe⁻^N² Tr ln[(J+Λ)]+iβ√ VP

iλ_~_q=0,i

, (3.11)

where we have defined the matrices

J_~_q~_q⁰_,ij =J_~_q,ijδ_~_q~_q⁰; Λ_~_q~_q⁰_,ij = 1

√V iλ_~q−~q⁰,iδ_ij. (3.12)

3.2.1 Separating the Constraint Field

To evaluate the integrals over the constraint field, we expand around the spatially homogeneous value of the constraint field as

λ_~_q,i =−i∆_i√

V δ_~_q,0+√

Vλ˜_~_q,i, λ˜_~_q=0,i = 0, (3.13)

(26)

where we have one homogeneous value ∆_i for each sublattice, and λ˜_~_q,i are the spatial fluctuations from ∆_i. ˜λ_~_q,i does therefore not have any ~q = 0 component. We will later determine the value of ∆_i by a saddle point approximation to ensure unit length of the spins on sublattice i.

Inserting the expansion of the constraint field into the partition function, we obtain the effective constraint field action

S[∆, λ] = N

2 Tr lnh β

J+∆−Λ˜i

−βV X

i

∆_i, (3.14)

with

∆_~_q~_q⁰_,ij = ∆_iδ_~_q~_q⁰δ_ij; Λ˜_~_q~_q⁰_,ij = (−i)˜λ_~q−~q⁰,iδ_ij. (3.15)

3.3 Diagrammatic Treatment of the Fluctuations

We treat the integrals over the spatial fluctuations of the constraint field by diagrammatic perturbation theory. When doing perturbation theory, one does an expansion of observables in powers of a small parameter, which in our case will be1/N, whereN is the number of spin components. The number of spin components is usually1(Ising),2(XY) or 3(Heisenberg), but we will anyhow develop a self-consistent theory in the limit of an asymptotically large N. We do so, hoping that we capture the important features of the systems also for the relatively small values of N.

We have the following partition function and action Z =

Z

D∆Dλ e˜ ^−S

h

∆,˜λ i

, Sh

∆,λ˜i

= N

2 Tr lnh β

K−Λ˜i

−βV X

i

∆_i, (3.16) where we have defined the matrix

K=J+∆⇒K_~_q~_q⁰_,ij ≡K_~_q,ijδ_~_q~_q⁰ ≡J_~_q,ijδ_~_q~_q⁰+ ∆_iδ_~_q~_q⁰δ_ij. (3.17) From this, we observe that K is diagonal in the momentum indices, as opposed to Λ˜, which is diagonal in the sublattice indices. We aim to rewrite the action to find the Gaussian action and the interaction terms, i.e. the quadratic terms in ˜λ and the higher order terms in λ˜, respectively. By expanding the action in powers of Λ˜, we get

Sh

∆,λ˜i

=−βV X

i

∆_i+ N

2 Tr lnK+S₂h

∆,λ˜i

+S_int.h

∆,λ˜i

, (3.18)

where the Gaussian action S₂ is given by

S2

h

∆,λ˜ i

= 1 2

X

~q6=0

X

ij

˜λ−~q,i



 N

2 X

~ p

K_~_q+~⁻¹_p,ijK_~_p,ji⁻¹



λ˜~q,j, (3.19)

(27)

3.3. Diagrammatic Treatment of the Fluctuations 15 where the~q = 0component is excluded from the sum, asλ˜has no~q= 0component. The interaction S_int. is given by

S_int.h

∆,λ˜i

=−N 2

∞

X

j=3

1 j Tr

K⁻¹Λ˜j

. (3.20)

This partitioning of the action enables us to define two types of expectation values:

The expectation value with respect to the full action S

h. . .i=

R Dλ˜(. . .)e^−S

h

∆,λ˜i

R Dλ e˜ ^−S

h

∆,λ˜i (3.21)

and the expectation value with respect to the Gaussian action S₂

h. . .i₂ =

R D˜λ(. . .)e^−S²

h

∆,˜λi

R Dλ e˜ ^−S²

h

∆,λ˜i . (3.22)

From the full action in Eq. (3.18), we see that the two expectation values are related by h. . .i= h. . . e^−S^int.i₂

he^−S^int.i₂ . (3.23)

3.3.1 Free Propagators

Having separated the action into a Gaussian term and an interaction term, we are in position to introduce the diagrammatic language. We begin by defining the free propagators of the theory, also referred to as the bare propagators.

The bare constraint field propagator is defined as the two-point correlation function with respect to the Gaussian action

D_0~_q,ij =D

λ˜_~_q,iλ˜−~q,j

E

2, (3.24)

which is a matrix in sublattice space of the same dimensions as J_~_q. From Eq. (3.19), we get

D_0~_q,ij = 2 N



 X

~ p

K_~_q+~⁻¹_p,ijK_~_p,ji⁻¹





−1

. (3.25)

The mathematical derivation is shown in Appendix A.1. From the expression, we see that the bare constraint field propagator is proportional to the perturbation parameter 1/N.

We represent the bare constraint field propagatorD_0~_q,ij by a bare wavy line with momentum~q flowing from sublattice index ito sublattice index j, as illustrated in Fig. 3.1.

In Appendix A.2.1, we show that the spin correlation function to zeroth order inλ˜ is 2

N Tχ_~_q,ij ≡ 2 N T

DS~_~_q,j·S~_~_q,iE

=K_~_q,ij⁻¹ +O λ˜

, (3.26)

(28)

D_0~_q,ij = i j

~ q

K_~_q,ij⁻¹ = i j

~ q

Figure 3.1: The bare propagators of the theory. The wavy line represents the bare constraint field propagator, and the solid line represents the bare spin propagator.

where we also have introduced the susceptibility, χ_~_q. From Eq. (3.26), it is clear that K_~_qij⁻¹ should be defined as the bare spin propagator. We represent it as a solid line with momentum ~q flowing fromi toj, as illustrated in Fig. 3.1.

The propagators have here been defined as matrices in sublattice space. We can also define them as matrices in both momentum space and sublattice space

D_0~_q~_q⁰_,ij =D_0~_q,ijδ_~_q,~_q⁰; K⁻¹_~_q~_q0,ij =K_~_q,ij⁻¹δ_~_q,~_q⁰. (3.27) As momentum has to be conserved, the propagator matrices are diagonal in momentum space.

3.3.2 Interaction Terms

Having defined the free propagators of the theory, we can move on to study the interaction terms S_int., defined in Eq. (3.20). The j-th term in the sum is

S_j =− N 2·j Tr

K⁻¹Λ˜j

, j ≥3. (3.28)

All the terms have a similar form. We will first study the term S3, which is cubic in λ˜ S₃ =− N

2·3TrK⁻¹ΛK˜ ⁻¹ΛK˜ ⁻¹Λ.˜ (3.29) Remembering the expressions for the matrix elements of K⁻¹ andΛ, which can be found˜ in Eq. (3.17) and Eq. (3.15), respectively, we end up with

S₃ =− N

2·3(−i)³ X

~ q1,~q2,~q3

X

i,j,k

K_~_q⁻¹

1,ijλ˜_~_q₁−~q2,jK_~_q⁻¹

2,jkλ˜_~_q₂−~q3,kK_~_q⁻¹

3,ki˜λ_~_q₃−~q1,i. (3.30) This term can be visualized diagrammatically as a spin propagator loop (spin loop) with three vertices, as shown in Fig. 3.2. Each vertex is labelled by a sublattice index and each line carries a momentum. In this diagram, each vertex is connecting two bare spin propagators K⁻¹ and a factor ˜λ. ˜λ is represented by a dangling wavy line. As already described, the bare spin propagatorsK_~_q,ij⁻¹ carries a momentum ~q from sublattice index i to sublattice index j. The convention used for the dangling wavy lines is that λ˜_~_q,i has momentum ~q going out from a vertex assigned a sublattice index i. From the mathematical expression for the diagram, Eq. (3.30), we see that each vertex is associated with a factor(−i)and that we have momentum conservation at each vertex. This can be

(29)

3.3. Diagrammatic Treatment of the Fluctuations 17 seen by for example considering all factors conneting to the vertex j, K_~_q⁻¹

1,ijλ˜_~_q₁_−~_q₂_,jK_~_q⁻¹

2,jk. Here, we have momentum~q₁ coming in through one bare spin propagator and momentum

~

q₂ going out through the other bare spin propagator. As the dangling wavy line carries momentum ~q₁−~q₂ out from the vertex, we have momentum conservation. This holds for all vertices. The diagram also has a numerical factor _2·3^N, which we associate to the spin loop.

j i

~ k

q₂ −~q₃ ~q₁−~q₂

~ q₃−~q₁

~

q₃ ~q₁

~q₂

Figure 3.2: The cubic term in the action, (−S3).

We generalize this to higher order terms. Due to the trace, combined withK⁻¹ being diagonal in momentum space and Λ˜ being diagonal in sublattice space, S_j forms a spin loop with j vertices, i.e. a spin loop withj dangling wavy lines. This spin loop is associated with a factor _2·j^N. Again, each vertex is associated with a factor(−i)and a sublattice index, and there is momentum conservation at each vertex. All undetermined sublattice indices and momenta should be summed over when writing out the mathematical expressions for the diagrams. From here on, we will omit labelling of such undetermined sublattice indices and momenta in the diagrams.

As each dangling wavy line represents a factorλ˜, it is clear that the expectation value of two dangling wavy lines with respect to the Gaussian action connects the two lines to give a bare constraint field propagator, as illustrated in Fig. 3.3.

i j

~

q −~q

h...i₂

−−→ i j

~ q

Figure 3.3: Diagrammatical representation of the Gaussian expectation value of two dangling wavy lines. Mathematically it is D

E

2

=D_0~_q,ij.

When computing expectation values of observables with respect to the full action, i.e.

Eq. (3.21), we see from Eq. (3.23) that this corresponds to taking the Gaussian expectation value of the observable multiplied by exp (−S_int.) and dividing by the Gaussian expectation value of exp (−S_int.). This denominator, R

D˜λ e^−S²

h

∆,˜λ i

, effectively corresponds to omitting all bubble diagrams of the numerator. In other words, computing an expectation value with respect to the full action is reduced to computing the fully connected diagrams of the numerator in Eq. (3.23), i.e. the fully connected diagrams

(30)

of the Gaussian expectation value of the observable multiplied by exp (−S_int.). This is known as thelinked-cluster theorem [24, 25].

Taylor expanding exp (−S_int.) gives terms consisting of a number of loops. The expansion also contributes with an extra numerical factor to the diagrams. This factor is (k₃!k₄!k₅!. . .)⁻¹, where k_j is the number of loops with j wavy lines in the diagram. As an example, the term illustrated in Fig. 3.4 has two loops with three wavy lines and one loop with four wavy lines, thus k₃ = 2 and k₄ = 1 and the diagram gets an additional factor (2!1!)⁻¹.

Figure 3.4: Example of a term in the expansion of exp (−S_int.).

The averaging over a diagram with respect to the quadratic action is performed by connecting all possible pairs of dangling wavy lines, analogous toWick’s theorem[24, 25].

If the number of vertices is odd, it is impossible to pair all the wavy lines, resulting in the average over the diagram being zero. An example of this is shown in Fig. 3.5a. The simplest non-zero examples of the Gaussian expectation value of a spin loop diagram is shown in Fig. 3.5b and c. The combinatorial factor in front of each diagram accounts for all the possible ways to pair up the wavy lines. In Fig. 3.5 the combinatorial factors are explicitly shown, but from here on, we will absorb these factors into the diagrams.

As the spatially homogeneous part of the constraint field is separated out, there is no ~q = 0 component of the constraint field propagator. Therefore all tadpole diagrams, i.e. diagrams where a propagator carries zero momentum, vanish. An example of such a diagram is the first diagram on the right hand side in Fig. 3.5c.

Before continuing, we should summarize how the diagrams depend on the perturbation

a)

* +

2

= 0

b)

* +

2

= 2 +

c)

* +

2

= 9 + 6

Figure 3.5: Examples of expectation values of different diagrams with respect to the Gaussian action. a) The average over a diagram with an odd total number of vertices is zero. b) Example of Wick’s theorem for a single loop with four vertices. c) Example of Wick’s theorem for two loops with three vertices each.

(31)

3.4. Effect of the Interaction Term 19 parameter: Each bare constraint field propagator goes as 1/N and each spin loop goes asN. This will be important to keep in mind in the coming sections, as we need to keep track of the order of the diagrams.

3.4 Effect of the Interaction Term

We have now introduced all formal elements of the diagrammatic perturbation theory, i.e. the Gaussian action, the free propagators, the interaction terms and how these are represented diagrammatically. We have also discussed how expectation values are calculated diagrammatically. The next step is to study how the interaction terms modify the free propagators of the theory.

3.4.1 Dressed Propagators

When computing propagators with respect to the full action, we get so-called dressed propagators. For the dressed constraint field propagator, we have

D_~_q,ij =D

E

=

D˜λ_~_q,i˜λ−~q,je^−S^int.E

2

he^−S^int.i₂ . (3.31)

We see that this gives all fully connected diagrams with two external wavy legs contracted with the expansion of exp (−S_int.). We represent D_~_q,ij diagrammatically as a bold wavy line. The lowest order diagrams contributing to the dressed constraint field propagator are shown in Fig. 3.6.

D_~_q,ij = i j

~ q

= i j

~q

+ i j

~

q ~q

+ i j

~

q ~q

+ i j

~

q ~q

+. . .

Figure 3.6: Lowest order expansion of the dressed constraint field propagator. The bare constraint field propagator is O(1/N) while the lowest order corrections (the three last diagrams shown on the right hand side) are O(1/N²).

When evaluating the dressed spin propagator, we have to take the full expression derived in Appendix A.2 into account. This gives

K_eff⁻¹_~_q,ij ≡ 2 N T

DS~_~_q,j ·S~_~_q,iE

=

*

1−K⁻¹Λ˜−1

K⁻¹

~ q~q,ij

+

. (3.32)

We have assumed inversion symmetry, which we do throughout the whole thesis. The dressed spin propagator is proportional to the spin correlation function, which is related to the structure factor S(~q)and the susceptibility. As the structure factor can be measured

(32)

experimentally, we will strive to calculate the dressed spin propagator as a function of temperature.

Taylor expanding and rewriting as averages with respect to the Gaussian action, we have

K_eff⁻¹_~_q,ij =K_~_q,ij⁻¹ +

P∞ n=1

h

K⁻¹Λ˜n

K⁻¹i

~ q~q,ij

e^−S^int.

2

he^−S^int.i₂ . (3.33)

The first term is the bare spin propagator and the second term contains the corrections.

Diagrammatically, the factorh

K⁻¹Λ˜n

K⁻¹i

~

q~q,ij in the second term is represented by a straight spin line with momentum~q coming in with a sublattice index iand momentum

~

q going out with a sublattice index j. The line has n vertices. An example of this factor is shown diagrammatically in Fig. 3.7 forn = 3.

q q

i j

Figure 3.7: Diagrammatic representation of the factor h

K⁻¹ΛK˜ ⁻¹ΛK˜ ⁻¹ΛK˜ ⁻¹ i

~ q~q,ij. We represent the dressed spin propagator diagrammatically as a bold straight line.

The diagrams that contribute up to second order in the perturbation parameter are shown in Fig. 3.8.

K_eff⁻¹_~_q,ij = i j

~ q

= i j

~ q

a) + i j

~

q ~q

b) + i j

~

q ~q c)

+ i j

~

q ~q

d) + i j

~

q ~q

h)

+ i j

~

q ~q

e) + i j

~

q ~q

f)

+ i j

~

q ~q

g) + i j

~

q ~q

i) +. . .

Figure 3.8: All terms up to O(1/N²) contributing to the dressed spin propagator. Di- agram a) is O(1), diagram b) is O(1/N) and the rest of the diagrams are O(1/N²).

(33)

3.4. Effect of the Interaction Term 21 The dressed propagators consist of infinitely many diagrams. It would be useful to divide these diagrams into two different sub-classes: diagrams that are so-called one- particle irreducible and diagrams that are not. We will put the one-particle irreducible diagrams of the two dressed propagators into objects known as the proper polarization Π and the proper self-energyΣ.

3.4.2 Polarization

The proper polarization Πconsists of all one-particle irreducible diagrams of the dressed constraint field propagator. By that we mean all the diagrams in the dressed constraint field propagator that cannot be divided into two disconnected pieces by cutting one internal bare constraint field propagator, excluding the external wavy legs. Having a diagrammatic expression for the dressed constraint field propagator up to second order in the perturbation parameter enables us to define the lowest order diagrams contributing to the proper polarization. From Fig. 3.6, we find that the lowest order diagrams in the proper polarization are the ones shown in Fig. 3.9. All the diagrams are O(1) as the external legs are not part of the proper polarization.

ΠΠΠΠΠΠΠΠ Π Π Π Π ΠΠ Π

ΠΠ = + +

Figure 3.9: The lowest order diagrams contributing to the proper polarization. All of these diagrams are of order unity as we have deleted the external wavy lines.

In Fig 3.9, we have not written out the external momenta and sublattice indices explicitly. The objectΠis indeed of the same dimensions as the constraint field propagator D₀, i.e. a matrix in both momentum and sublattice space. We specify an element of the polarization as

Π_~_q~_q⁰_,ij = Π_~_q,ijδ_~_q,~_q⁰, (3.34) where Π_~_q,ij carries momentum ~q from sublattice index i to sublattice index j. As a consequence of momentum conservation, Π is diagonal in momentum space.

3.4.3 Self-Energy

In the same way as the proper polarization is related to the dressed constraint field propagator, the proper self-energyΣis related to the dressed spin propagator. The proper self-energy is thus defined as all the one-particle irreducible diagrams in the dressed spin propagator. That is all diagrams that cannot be divided into two separate pieces by cutting one internal bare spin propagator. Again, we exclude the external legs. Using the diagrammatic expression for the dressed spin propagator found in Fig. 3.8 we obtain the diagrams contributing to the proper self-energy up to second order. These diagrams are shown in Fig. 3.10. The first diagram on the right hand side is O(1/N), while the last six areO(1/N²).

(34)

Σ Σ Σ Σ ΣΣ Σ Σ Σ Σ Σ ΣΣ Σ Σ

ΣΣ = + + + + + +

Figure 3.10: The diagrams contributing to the proper self-energy up to O(1/N²). The first diagram on the right hand side is O(1/N) while the rest are O(1/N²).

Again, we have not included momentum and sublattice indices in the diagrams. The object Σis of the same dimensions as the spin propagator and we specify an element of the proper self-energy as

Σ_~_q~_q⁰_,ij = Σ_~_q,ijδ_~_q,~_q⁰, (3.35) whereΣ_~_q,ij carries momentum~q from sublattice indexito sublattice indexj. Once again δ_~_q,~_q⁰ ensures momentum conservation.

3.4.4 Dyson Equations

For simplicity, we will from now on omit the “proper” in front of the polarization and self-energy. Using expressions for the polarization and the self-energy correct to all orders in the perturbation parameter, the exact expressions for the dressed propagators could be reproduced recursively by writing down the Dyson equations shown in Fig. 3.11. These equations show how the polarization and self-energy renormalize the bare propagators.

a) = + Π

b) = + Σ

Figure 3.11: Dyson equations for a) the dressed constraint field propagatorD and b) the dressed spin propagator K⁻¹_eff.

The Dyson equations can be expressed mathematically as D_~_q,ij =D_0~_q,ij+X

kl

D_0~_q,ikΠ_~_q,klD_~_q,lj (3.36) for the dressed constraint field propagator, and

K_eff⁻¹_~_q,ij =K_0~⁻¹_q,ij +X

kl

K_0~⁻¹_q,ikΣ_~_q,klK_eff⁻¹_~_q,lj (3.37) for the dressed spin propagator. In both cases, we have to sum over the sublattice indices of the polarization and self-energy, meaning that the Dyson equations can be viewed matrix products in sublattice space. Solving these matrix equations, we obtain

D_~_q⁻¹ =D⁻¹_0~_q −Π~q (3.38) and

K_eff_~_q =K_0~_q−Σ_~_q. (3.39)

Nematic Bond Theory of Frustrated Heisenberg Models on Triangular and Honeycomb Lattices

Nematic Bond Theory of

Frustrated Heisenberg Models on Triangular and Honeycomb Lattices

Cecilie Glittum

Thesis submitted for the degree of Master of Science

Department of Physics

Faculty of Mathematics and Natural Sciences UNIVERSITY OF OSLO

Nematic Bond Theory of Frustrated Heisenberg Models on Triangular and

Honeycomb Lattices

Cecilie Glittum

Contents

Chapter 1 Introduction

Chapter 2 Background

2.1 The Origin of Magnetic Interactions

2.2 Frustrated Magnetism

2.3 The Heisenberg Model on Non-Bravais Lattices

2.3.1 Classical spins

2.4 Phase Transitions

2.4.1 Spontaneous Symmetry Breaking and the Mermin-Wagner The- orem

2.5 Spiral States

2.5.1 Lattice-Nematic Phases

2.5.2 Entropy

Chapter 3

The Nematic Bond Theory

3.1 Notation

3.1.1 Fourier Conventions

3.2 The Heisenberg Model and the Unit Length Constraint

3.2.1 Separating the Constraint Field

3.3 Diagrammatic Treatment of the Fluctuations

3.3.1 Free Propagators

3.3.2 Interaction Terms

3.4 Effect of the Interaction Term

3.4.1 Dressed Propagators

3.4.2 Polarization

3.4.3 Self-Energy

3.4.4 Dyson Equations