PHYSICAL REVIEW B 89, 235114 (2014)

Mean-field model for the quadrupolar phases of UPd3

Manh Duc Le,1,2,3,* Keith A. McEwen,3 Martin Rotter,4,5 Mathias Doerr,6 Alexander Barcza,5,7 Je-Geun Park,1,2
James Brooks,8 Eric Jobiliong,8,9 and David Fort10
1
Center for Correlated Electron Systems, Institute for Basic Science (IBS), Seoul 151-747, Korea
2
Department of Physics and Astronomy, Seoul National University, Seoul 151-747, Korea
3
Department of Physics and Astronomy, and London Centre for Nanotechnology, University College London,
Gower Street, London WC1E 6BT, United Kingdom
4
Max Planck Institute for Chemical Physics of Solids, Nöthnitzer Str. 40, D-01187 Dresden, Germany
5
Institut für Physikalische Chemie, Universität Wien, Währingerstr. 42, 1090 Wien, Austria
6
Institut für Festkörperphysik, Technische Universität Dresden, D-01062 Dresden, Germany
7
Vakuumschmelze GmbH, Grünerweg 37, D-63450 Hanau, Germany
8
National High Magnetic Field Laboratory, 1800 Paul Dirac Dr., Tallahassee, Florida 32306, USA
9
Department of Industrial Engineering, Universitas Pelita Harapan Karawaci, Banten 15811, Indonesia
10
Department of Metallurgy and Materials, University of Birmingham, Birmingham B15 2TT, United Kingdom
(Received 19 March 2014; revised manuscript received 6 May 2014; published 10 June 2014)

UPd3 is known to exhibit four antiferroquadrupolar ordered phases at low temperatures. We report measur-
ements of the magnetization and magnetostriction of single crystal UPd3 , along the principal symmetry directions,
in fields up to 33 T. These results have been combined with recent inelastic neutron and x-ray resonant scattering
measurements to construct a mean-field model of UPd3 including up to fourth nearest-neighbor interactions.
In particular, we find that anisotropic quadrupolar interactions must be included in order to explain the low-
temperature structures derived from the scattering data.

DOI: 10.1103/PhysRevB.89.235114 PACS number(s): 75.80.+q, 75.40.Cx, 75.10.Dg, 71.10.Hf

I. INTRODUCTION magnetic susceptibility [8], and electrical resistivity [9] of

UPd3 at low temperatures, which are indicative of phase tran-
Low-temperature phase transitions in condensed matter
sitions. Polarized neutron diffraction measurements revealed
systems are usually driven by the cooperative actions of
superlattice peaks at Q = ( 21 0 l), which appear only in the
the system’s electronic degrees of freedom. In many cases,
non-spin-flip channel, and thus can be attributed to a structural
the exchange interactions between the spins of electrons on
transition [10], from hexagonal to orthorhombic symmetry.
neighboring ions causes them to order below a characteristic
X-ray resonant scattering (XRS) measurements [11]
transition temperature. However, in lanthanide compounds, showed that this transition arises from the ordering of the
the localized 4f electrons on each ionic site may have electric quadrupole moment of the 5f 2 electrons. In addition
nonspherical charge or current distributions, which may be to this phase transition at T0 = 7.8 K, there are three further
described by a multipole expansion of their electric or mag- transitions to different antiferroquadrupolar (AFQ) ordered
netic fields [1]. These multipoles may interact, and in certain states at T+1 = 6.9 K, T−1 = 6.7 K, and T2 = 4.4 K. Below
cases these interactions may be stronger than the spin exchange T0 a superlattice peak at ( 21 0 l) where l is odd is observed,
interactions, driving a phase transition to a multipolar ordered whilst below T−1 there are additional peaks at ( 12 0 l) where
phase. There have been many examples of electric quadrupolar l is even. The l odd peaks show that there is antiferro-
order observed, such as in CeB6 [2,3], PrPb3 [4], and TmTe [5], quadrupolar ordering along the c direction, also denoted
whilst higher-order electric multipoles were observed in as an antiphase stacking of quadrupoles. The l even peaks
DyB2 C2 [6]. show an additional ordering of quadrupole moments in-phase
In contrast to the lanthanides, there have not been many along c.
cases of multipolar order observed in actinide compounds. Measurements of the order parameter using x-ray resonant
This is due partly to the larger spatial extent of the 5f wave scattering show that the l odd order is associated mainly with
functions, which causes them to become delocalized, and Qzx quadrupoles [12], whilst the l even order is associated
hence invalidates any description of their electric or magnetic with Qxy quadrupoles [13]. The Qzx ordering is accompanied
fields about a particular point. The insulating actinide oxides by a component from the Qx 2 −y 2 quadrupoles, whilst the Qxy
UO2 and NpO2 , however, have been studied extensively, is accompanied by a Qyz component. The directions x and
and exhibit electric quadrupolar and high order magnetic z are equivalent to the a and c crystallographic directions
multipolar order, respectively [1]. In contrast, UPd3 is one of and y is perpendicular to both. In addition, there is also an
the few metallic actinide compounds that has well-localized additional ordering of the Qyz quadrupoles in antiphase below
5f electrons, and it was one of the first compounds to T+1 . This sequence of phase transitions is summarized in
be shown to exhibit quadrupolar ordering. Anomalies were Fig. 1.
observed in the heat capacity [7], thermal expansion [7], UPd3 adopts the double hexagonal close packed TiNi3
4
crystal structure (space group D6h , P 63 /mmc, no. 194) [14]
with lattice parameters a = 5.76 Å and c = 9.62 Å. This
*
Corresponding author: mducle@snu.ac.kr means that the nearest-neighbor U-U distance, 4.1 Å, is larger

1098-0121/2014/89(23)/235114(11) 235114-1 ©2014 American Physical Society

MANH DUC LE et al. PHYSICAL REVIEW B 89, 235114 (2014)

FIG. 1. (Color online) The quadrupolar phases of UPd3 . Labels below the line indicate order parameter as determined from XRS data. “In
phase” indicates the quadrupole moments are aligned in parallel along the c axis. This corresponds to a superlattice reflection at Q = (1 0 4) in
orthorhombic notation, or ( 12 0 4) in hexagonal notation. “Antiphase” indicates an antiparallel alignment along c and corresponds to a reflection
at Q = (1 0 3)orth or ( 21 0 3)hex . (After Walker et al. [12,13]).

than the Hill limit (≈3.5 Å) [15] and thus the 5f 2 electrons are to the magnetic dipoles in the system, and the energies are far
well localized. The U4+ ions are in the 2a and 2d sites, which too low for inelastic x-ray scattering to resolve. However, the
have, respectively, D3d (3̄m) and D3h (6̄m2) point symmetry, wealth of physical property measurements available on UPd3
which we shall refer to as quasicubic and hexagonal. The XRS has encouraged us to try to synthesize this into a mean-field
data show that the ordering involves primarily the ions on the model which explains, to a large extent, these varied mea-
quasicubic sites [12,13]. surements. We hope that this may encourage the construction
Knowledge of the crystal-field (CF) interactions is essential of models to explain the properties of similar quadrupolar
in determining the more complex two-ion interactions. In (or higher multipolar) ordered compounds thus deepening the
particular, it is crucial to know the CF ground state. Moreover, understanding of what lies behind these phenomena.
the excited states and the matrix elements of the angular
momentum operators, Jˆi , between these and the ground state
II. CRYSTAL-FIELD INTERACTIONS
determine the intensities of excitations observed by inelastic
neutron scattering and also to some extent the magnitude The single-ion properties of uranium are generally found to
of the magnetization in the ordered phase, as explained by be close to the LS-coupling limit [16]. So, in order to simplify
McEwen et al. [8]. Thus we shall first consider in Sec. II the the analysis, we shall ignore any mixing with higher order
CF level scheme deduced from inelastic neutron scattering multiplets in determining the CF parameters. As mentioned
measurements in the paramagnetic phase. in Sec. I, there are two inequivalent sites for the U4+ ions in
We then present high-field magnetization and magnetostric- the crystal structure of UPd3 . The different point symmetry of
tion measurements in Sec. III, from which the critical fields these sites gives rise to different crystal fields, but it happens
and magnetic phase diagrams were determined experimentally. that both split the ninefold degenerate J = 4 ground multiplet
Finally, a mean-field model of the quadrupolar ordering is into three singlets and three doublets. The energies and wave
constructed in Sec. IV. The order parameters of each of the functions of these levels, however, are different for the two
quadrupolar phases (determined by resonant x-ray scattering) sites.
and their transition temperatures, were used to constrain the The nature of the ground state may be deduced from single
two-ion quadrupolar interaction parameters for the model, crystal susceptibility measurements and estimates of the mag-
whilst the dipolar exchange parameters were determined from netic entropy determined from heat capacity measurements.
the measured critical fields. The results of the model are then These results indicate a singlet ground state on the hexagonal
compared to the measured high-field data. sites and a doublet on the quasicubic sites [8].
Whilst the model, which we shall present, is empirically The CF split energy levels were determined from previously
based on the measured physical properties of UPd3 rather than reported inelastic neutron scattering measurements made on
directly on measured exchange constants, we note that it is the time-of-flight spectrometer HET at the ISIS Facility, United
the first attempt to explain comprehensively this fascinating Kingdom [17], and on the triple-axis-spectrometer IN8 at ILL,
compound with its many competing ordered phases from a mi- Grenoble [18]. We identified magnetic excitations at 4.1, 9.7,
croscopic point of view and hope to stimulate further ab initio 12.3, 16.8, 20.4, and 30 meV. The 16 meV peak exhibits
studies of the exchange interactions involved. The quadrupolar considerable dispersion and is assigned to the transition
interactions, which drive the many phase transitions in UPd3 , between the |Jz = 0 singlet ground state and the |Jz = ±1
are difficult to measure directly, because neutrons couple only doublet excited state of the hexagonal site ions. Its dispersion

235114-2
MEAN-FIELD MODEL FOR THE QUADRUPOLAR PHASES . . . PHYSICAL REVIEW B 89, 235114 (2014)

was used to determine the exchange interactions between TABLE I. Crystal-field parameters in Stevens normalization in
hexagonal sites [19]. The remaining peaks are assigned to meV. The parameters were deduced from fitting to inelastic neutron
the quasicubic site ions. spectra and the constraints on the wave functions of the lowest
Henceforth we shall be concerned mainly with the qua- three energy levels of the quasicubic sites as described in the
sicubic sites, as the quadrupolar ordering primarily involves text.
the uranium ions on these sites. These sites have trigonal, 3̄m
(D3d ), point symmetry, so that the crystal-field Hamiltonian is B02 0.035 B06 −0.00012
  B04 −0.012 B36 0.0025
Hcf = Bk0 Ok0 + Bk3 Ok3 + B66 O66 , (1) B34 −0.027 B66 0.0068
k=2,4,6 k=4,6
q q
where Bk are crystal-field parameters and Ok are Stevens
operators. The quantization (z) axis is taken to be the trigonal initial set of parameters to fit the measured energy levels. These
axis, which in this case is parallel to c. refined parameters are subsequently used to calculate the b, b(2)
From the measured transition energies and with the re- and e matrix elements, from which the simulated annealing
striction of a doublet ground state, a crystal-field fitting “energy” is obtained, and hence minimized. Figure 2 shows
program [20] was used to obtain initial estimates of the CF the resulting crystal-field energy splitting and wave functions
parameters for the quasicubic sites. This program relies on for the U4+ ions on the quasicubic sites.
the orthogonality of the spherical harmonic functions from We note that these energy levels, and the crystal-field
which the CF operators are constructed. It allows one to find parameters from which they arise, do not show any cubic
q
a set of CF parameters, Bk , given the energy levels and wave symmetry. Indeed, even a point charge calculation assuming
functions produced by the crystal field. The fitting algorithm Pd+ and U4+ ions yields a crystal-field splitting quite far from
may thus vary either the wave functions to fit a particular set the cubic expectation: in common with Buyers et al. [22],
of energy levels, or vice versa. In this case, however, we also we find a singlet ground state, doublets at ≈3.2, 8.8, and
face constraints on the wave functions. 28.6 meV above this, and two singlets at ≈5.8 and 12.7 meV.
From symmetry considerations, the doublet ground states Were the crystal-field symmetry truly cubic, the two low-lying
have the wave functions doublets and singlets would be degenerate, as deduced by the
|d1  = a|4 + b|1 + c| − 2, LLW [23] parameters W = 0.408 meV and x = 0.44, which
most closely approximate the point charge parameters. Thus
|d2  = a| − 4 − b| − 1 + c|2, (2) the term “quasicubic” is somewhat of a misnomer and has been
retained here only for continuity with the literature, while the
where for brevity the kets denote states |J = 4,Jz . The singlet
true symmetry of the crystal field at this site is trigonal.
wave functions have the forms
|s = d|3 + e|0 − d| − 3,
1
|s   = √ (|3 + | − 3) . (3)
2
In order to ensure that the T0 = 7.8 K transition is accompa-
nied by only a very small entropy change, as deduced from the
heat capacity data, the Landau theory analysis of McEwen
et al. [8] requires that the matrix elements d1 |Qzx |d2  =
d2 |Qzx |d1  ≈ 0, where Qzx = 12 (Jx Jz + Jz Jx ). As shown in
the reference, this implies that bc ≈ 0.
In addition, we note that the basal plane susceptibilities
χx,y increase with decreasing temperature through the T−1 =
6.7 K and T2 = 4.4 K transitions. This may be explained if
the first excited state is a singlet and there is a large Jˆx,y
matrix element between it and some higher-energy state which
increases the x- or y-direction susceptibility in the ordered
phases as progressively more of the singlet state is mixed in
with the doublet ground state [21]. As the Jˆx,y matrix elements
between singlet states are zero, this coupling must be to a
(2) (2)
higher lying doublet, |d1,2 . The condition that s|Jˆx,y |d1,2 
be large, whilst s|Jˆx,y |d1,2  is small is thus satisfied if e ≈ 1,
b(2) ≈ 1 and b ≈ 0. FIG. 2. U4+ quasicubic site crystal field level scheme. The
These requirements are satisfied by the crystal-field pa- crystal-field energy levels and corresponding wave functions in the
rameters in Table I, which yield b = 0.02, b(2) = 0.99 and paramagnetic phase expressed in the |J = 4,Jz  basis. Arrows denote
e = 0.92. The parameters were obtained using a simulated transitions from the ground state with non-negligible dipole matrix
annealing minimization procedure whereby at each iteration, elements whose squared values are shown as numbers near the
the algorithm mentioned above [20] was used to refine an arrow.

235114-3
MANH DUC LE et al. PHYSICAL REVIEW B 89, 235114 (2014)

III. HIGH-MAGNETIC FIELD MEASUREMENTS 30

Single crystals were grown by the Czochralski technique at
the University of Birmingham and cut with faces perpendicular
to the orthogonal axes x, y, and z, where xa and zc. 20

H||x (T)
HF
These were used in magnetization and magnetostriction
measurements in fields up to 14 T at Birkbeck College [24]
by two of us (KAM and JGP). Subsequently, high-field 10
AFQ1
magnetization measurements were carried out (by JGP) at
the Grenoble High Magnetic Field Laboratory, and magne- AFQ3 AFQ2
PM
tostriction measurements at the National High Magnetic Field 0
Laboratory, Tallahassee. The magnetostriction was measured 30
using a miniature capacitance dilatometer [25] in which the
single-crystal samples were mounted with either the x, y,
or z faces parallel to the capacitor plates. The dilatometer 20

H||y (T)
could be rotated so that the magnetic field is perpendicular
to the capacitor plate allowing the transverse components of
magnetostriction to be measured. 10
AFQ1
Figure 3 shows the magnetization at several different PM
applied fields, whilst Fig. 4 collates this and other data [26] AFQ3 AFQ2
to construct the magnetic phase diagrams of UPd3 . The data 0
show the T0 transition increasing in temperature with increas- 30
ing field, which is characteristic of an antiferroquadrupolar
transition. In general, three ordered phases can be identified
20 PM
from the data, as the phase between T−1 and T+1 cannot be
H||z (T)

distinguished from the magnetization data.

Figure 5 shows the forced magnetostriction data at 4.2 K,
plotted as l/ l = [l(H,T = 4.2K) − l(H = 0,T = 4.2 K)]/ 10
l(H = 0,T = 4.2 K). We measured the longitudinal compo-
nents of magnetostriction with the field in the x, y, and z AFQ3 AFQ2 AFQ1
directions, and also the transverse components y/y and 0
2 4 6 8 10
z/z with applied field parallel to x and x/x with field Temperature (K)

FIG. 4. (Color online) Magnetic phase diagrams of UPd3 for

0.8 fields parallel to a, b, and c. Filled squares are from previous
H||a
Moment (μ /ion)

23T measurements by McEwen et al. [24]. Filled triangles are from

22T
high-field magnetization measurements described in this work, whilst
B

21T
0.6 19T filled diamonds are from high-field magnetostriction. The shaded
17T background shows the MCPHASE calculated phase diagrams, with
15.5T
14T
zero-field phases AFQ1-AFQ3 as described in the text. PM is
0.4 paramagnetic (or ferromagnetically polarized high-field) phase, and
0.8 HF is the high-field phase for H ||a, which has the ordering wave
H||b
Moment (μ /ion)

23T vector ( 21 21 0) with all quadrupolar moments ordered.

20T
B

0.6 17T
14T parallel to y. The measurements were repeated to confirm
the reproducibility of the data, and subsequently, binned and
0.4 averaged. In addition, the signal was corrected for artifacts
0.3 due to eddy currents. The mechanical noise from the magnet
H||c
Moment (μ /ion)

13T cooling system and electrical noise in the leads meant that we
10T obtained a resolution of 10−6 in l/ l.
B

0.25
For H ||x, the magnetostriction parallel (perpendicular) to
0.2 the applied field first decreases (increases) until approximately
3 T, then increases (decreases) to about 15 T before decreasing
0.15 (increasing) slightly. Similar, but less pronounced, behavior is
0 5 10 15 20 25
Temperature (K)
also observed for H ||y.
Transitions at high field were observed when the field was
FIG. 3. (Color online) The high-field magnetization of UPd3 as a applied in the basal plane, with a slight anisotropy between
function of temperature. Solid lines are calculated from the mean-field the x and y directions. With applied field parallel to the x
model. The calculated values in the bottom panel (for field parallel direction we see a steplike change in both the longitudinal and
to c) have been divided by 3. transverse magnetostriction at 28 T, whereas for field parallel

235114-4
MEAN-FIELD MODEL FOR THE QUADRUPOLAR PHASES . . . PHYSICAL REVIEW B 89, 235114 (2014)

20 to the y direction a change of slope is observed at 27 T. For

H||x both field directions, the magnetostriction parallel to the field
−5

10 decreases (the crystal contracts) at the phase transition, whilst

Δ l/l × 10

the magnetostriction perpendicular to the field increases (the

0 crystal expands).
When the field is applied parallel to the z direction, we
observed no high-field transitions above 20 T, but instead see
−10
40 anomalies (indicated by arrows in Fig. 5) in the longitudinal
H||y magnetostriction, at 7, 11, and 17 T, in agreement with our
−5

20 magnetization data and the phase diagram of Tokiwa et al. [9].

Δ l/l × 10

0
IV. QUADRUPOLAR TWO-ION INTERACTIONS
−20 AND MEAN-FIELD MODEL FOR UPD3
−40 The ordered quadrupolar structures were calculated from a
10 mean-field model, with quadrupolar interactions between the
H||z Δ x/x
quasicubic site ions, using the package MCPHASE [27,28]. In
−5

0 Δ y/y
Δ l/l × 10

order to determine the stable ordered structure, a set of super-

Δ z/z cells and corresponding wave vectors is generated. From the
−10 wave vector, a configuration of moments (dipole, quadrupole,
etc.) is generated and used as an initial configuration for a
−20 self-consistent mean-field calculation. For each solution of
0 5 10 15 20 25 30 35
Applied Field (Tesla) the mean-field iteration, the free energy is calculated. The
self-consistent ordered structure with the lowest free energy is
FIG. 5. (Color online) The high-field longitudinal and transverse taken to be stable and used for the computation of the physical
magnetostriction of UPd3 at 4.2 K. Solid vertical lines indicate the properties. In this way, the phase boundaries between the
observed high-field transitions. Blue circles, green squares and red different quadrupolar ordered structures were determined in
triangles denote the magnetostriction parallel to the x, y, and z order to construct the magnetic phase diagrams of UPd3 .
directions, respectively, for the indicated field directions. Solid lines The Hamiltonian
are the results of the mean-field calculations.
⎧
 1⎨  
ij ij
 ij  i j ij  
H= Hcf
i
+ HZi − J11 Jx Jx + Jy Jy + K11
ˆ ˆ ˆ i ˆj
cos(2φij ) Jˆxi Jˆxj − Jˆyi Jˆyj − sin(2φij ) Jˆxi Jˆyj + Jˆyi Jˆxj
i
2⎩ ij ij

ij

ij
ij

ij
ij

ij

+ J10 Jˆzi Jˆzj + J21 Q̂zx Q̂zx + Q̂yz Q̂yz + K21
i j i j
cos(2φij ) Q̂izx Q̂jzx − Q̂iyz Q̂jyz
ij ij ij

 ij

ij
j
+ sin(2φij ) Q̂izx Q̂jyz + Q̂iyz Q̂jxz + J22 Q̂ix 2 −y 2 Q̂x 2 −y 2 + Q̂ixy Q̂jxy
ij
⎫

ij
 i  ⎬
ij j j
+ K22 cos(4φij ) Q̂x 2 −y 2 Q̂x 2 −y 2 − Q̂ixy Q̂jxy − sin(4φij ) Q̂ixy Q̂x 2 −y 2 + Q̂ix 2 −y 2 Q̂jxy (4)
⎭
ij

was employed, where the site indices i and j run over symmetry of the interactions is satisfied. Such expressions
nearest- (ij = nn) and next-nearest neighbors (ij = nnn) within were used to explain the properties of elemental Pr [31], and
an ab plane and nearest- (ij = nnc) and next-nearest neighbors the dispersion of crystal-field excitations due to interactions
(ij = n3c) between planes, and φij is the angle in the ab plane between the hexagonal sites in UPd3 [19]. We note that the
between the projection of the vector between i and j in this form of the local anisotropic interaction chosen is higher than
plane and the x (a) axis. Hcf
i
is the crystal-field Hamiltonian of that demanded by the symmetry of the midpoint between two
the ith quasicubic ion given in Eq. (1) and HZi is the Zeeman quasicubic sites (C2h ), but is required in order to ensure that
Hamiltonian, −gJ μB Ji · H. The form of the two-ion exchange the dispersive crystal-field modes remain degenerate at certain
Hamiltonian was derived by considering a high-symmetry high-symmetry q vectors as observed in the neutron scattering
anisotropic interaction between each pair of neighbors in data [19].
a local coordinate system defined by an x  axis along the The Hamiltonian (4) suffices to describe the ordered phases
projection of the bond in the ab plane, and z axis along c, and of UPd3 . In particular, the anisotropic exchange terms (∝ K)
then rotating them into a global coordinate system [29,30] are required since the first (second) order isotropic quadrupolar
as described in Appendix. In this way, the full hexagonal interactions couple the Q̂zx (Q̂xy ) and Q̂yz (Q̂x 2 −y 2 ) operators

235114-5
MANH DUC LE et al. PHYSICAL REVIEW B 89, 235114 (2014)

equally, but the principal order parameters measured by reso-

nant x-ray diffraction are Q̂zx and Q̂xy . Thus the anisotropic
terms are required to favour these (Q̂zx , Q̂xy ) interactions over
the others (Q̂yz , Q̂x 2 −y 2 ).
In order to stabilize a structure with antiphase (in-phase)
stacking along c, we require J nnc < 0 (J nnc > 0), whilst the 1
4

AFQ order in the ab plane requires J nn,nnn < 0. Thus

Jzxnnc = J21
nnc
+ K21
nnc
cos 2φnnc < 0,
Jxy
nnc
= J22nnc
− K22
nnc
cos 4φnnc > 0, (5)
 nnc   nnc 
J  > J ,
zx xy
 nn  nn 1
J21 − K21 cos 4φnn < J22 + K22
nn nn
cos 2φnn < 0 (6) 4

ao
ah
should be satisfied in order to yield two phases with Q̂zx
antiphase (Q̂xy in-phase) order along c at higher (lower) tem-
bo bh
peratures, as measured. Furthermore, because the Q̂zx (Q̂xy )
and Q̂x 2 −y 2 (Q̂yz ) operators share the same symmetry [8], an
FIG. 6. (Color online) The double-hexagonal close packed struc-
ordering of one of these pairs will induce a secondary ordering ture. Only U atoms are shown, each of which has sixfold coordination
of the other quadrupole of the pair on the same site. That is, with neighboring Pd atoms. Circles indicate the quasicubic, whilst
a nonzero expectation value Q̂zx  implies Q̂x 2 −y 2  = 0 also hexagons the hexagonal sites. Labels in brackets apply to atoms at
(angled
 brackets denote the thermal expectation values Ô = z = 12 . Dashed lines show the exchange interaction pathways. Figures
n n| Ô|n exp( −E
kB T
n
)/Z, where the states |n are eigenstates on the right illustrate the combination of order parameters, which
of the Hamiltonian (4) and Z is the partition function). lead to a different quadrupolar moments on each quasicubic site in
This contributes to the effective field acting on an ion and the unit cell, and thus to four transitions observed by inelastic neutron
combines with the exchange interaction to reinforce (or scattering at 2 K, as explained in the text.
suppress) the ordering of some particular quadrupole. The
strength of this contribution is dependent on the crystal-field
wave functions, and we found that for UPd3 , the Q̂x 2 −y 2  sites C1-C4, and thus the splitting of the ground state doublet,
moment induced by a Q̂zx order is of the same order of be unique. This is illustrated schematically on the right side
magnitude as the primary Q̂zx  moment and acts to reinforce of Fig. 6.
the Q̂x 2 −y 2 ordering (induced Q̂x 2 −y 2 >0). This means that, Finally, we calculate that for the CF scheme outlined in
unfortunately, in our model, there is always a large Q̂x 2 −y 2 Sec. II, d1,2 |Q̂zx |s and d1,2 |Q̂yz |s are an order of magnitude
moment in disagreement with the measured XRS azimuthal lower than d1,2 |Q̂xy |s and d1,2 |Q̂x 2 −y 2 |s, which means that
dependence, which indicates a 10%–15% contribution. the exchange parameters J21 should be an order of magnitude
The other induced moments are generally an order of larger than J22 to give similar ordering temperatures for the
magnitude weaker than their primary order parameter. They first and second order quadrupoles as observed [32]. From
are needed, though, to account for the observation by inelastic the above considerations, we first determined the order of
neutron scattering of four almost dispersionless excitations magnitude of exchange parameters that result in ordering
below 4 meV at 1.8 K, which arise from transitions between temperatures below 10 K. Subsequently, a simulated annealing
the levels of the ground state doublet on the quasicubic and particle swarm optimization [33] search was carried out
sites, whose degeneracy is lifted by the quadrupolar order. to find sets of parameters that yield at least two transitions
As there are four quasicubic sites in the ordered unit cell, at ≈4 and ≈8 K and a splitting of the ground-state doublet
this implies that the splitting on each site is different, at 2 K close to the measured values 1.28, 1.68, 2.20, and
which may only occur if the magnitude of the quadrupolar 2.60 meV. Additional criteria, including the requirement that
moments on each site is different. In-phase ordering of the the calculated magnetization with an in-plane applied field
quadrupoles along c (denoted CS for c same in the following should increase with decreasing temperature, and that there
for brevity) means that the moments on sites C1-C4 of should be nonzero moments of all quadrupoles in the lowest
Fig. 6 are Q̂C1  = Q̂C2  = −Q̂C3  = −Q̂C4 , whilst an- temperature phases were then used to sort the candidate sets
tiphase ordering (henceforth denoted CD, c-different) implies of parameters found by the search. We then calculated the
that Q̂C1  = −Q̂C2  = −Q̂C3  = Q̂C4 . However, since magnetic phase diagram for the four best sets of parameters,
|Q̂xy | = |Q̂x 2 −y 2 | and |Q̂zx | = |Q̂yz | [8], the combina- and selected the set which have phase boundaries most similar
tion of Q̂zx CD, Q̂yz CD, Q̂xy CS, Q̂x 2 −y 2 CS ordering imposed to those measured. Finally, the parameters were refined by
by the Hamiltonian (4) and conditions (5) and (6) will result hand to better match the transition temperatures and fields.
in Q̂C1  = −Q̂C3  = Q̂xy  and Q̂C2  = −Q̂C4  = Q̂zx  so We found during this procedure that unless the dipolar in-
that only two excitations may be expected. Only by including teractions are included, structures with ordering wave vectors
the induced moments, which yields Q̂x 2 −y 2 CD and Q̂yz CS ( 12 12 0) are favoured over the ( 12 00) observed (indexed with
ordering (amongst others), will the moments on each of the respect to the dhcp cell). Thus small values of J1m and K1m

235114-6
MEAN-FIELD MODEL FOR THE QUADRUPOLAR PHASES . . . PHYSICAL REVIEW B 89, 235114 (2014)

TABLE II. Deduced exchange parameters in meV.

J10nn −0.017 J10nnc −0.04

J11nn −0.01 K11
nn
−0.02 J11nnc −0.01 K11
nnc
0
J21nn 0.01827 K21
nn
0.00107 J21nnc −0.06088 K21
nnc
−0.00583
J21nnn −1.16×10−5 K21
nnn
0.00778 J21n3c −0.00973 K21
n3c
−2.25×10−4
J22nn −8.42×10−4 K22
nn
−3.65×10−4 J22nnc −2.55×10−5 K22
nnc
0.52×10−5
J22nnn −0.00348 K22
nnn
0.00173 J22n3c 0.00318 K22
n3c
2.12×10−4

were included, but not varied in the search procedure. They not obtain good fits with simple analytical forms of the RKKY
were subsequently refined by hand along with the quadrupolar exchange. A more sophisticated approach, using the measured
interaction parameters to better fit the measured critical fields. band structure of the Pd-U conduction band, may give better
The final parameters are shown in Table II. At 2 K, the full results. This consideration may also apply to the quadrupolar
mean-field model with these parameters yields a splitting of the interaction parameters, and may account for the discrepancies
doublet ground state on ions C1-C4 of 0.69, 1.65, 1.82, and between the measured and calculated phase boundaries.
2.08 meV, respectively, [at Q = ( 32 23 0)], which is somewhat Finally, we note that ultrasound measurements [35] and a
lower than the experimentally measured values (1.28, 1.68, symmetry analysis of the XRS data [36] showed that there
2.20, and 2.60 meV [34]). The dispersion of the levels along is a sequence of structural transitions from hexagonal to
[00l], ≈0.5 meV, is close to the measured value but the in-plane orthorhombic symmetry at T0 and from orthorhombic to mon-
dispersion of ≈1 meV is in stark contrast to the measurements, oclinic at T+1 . As the parameters in Table II have the symmetry
which showed the modes to be almost dispersionless. of the high-temperature hexagonal structure, it is possible
Furthermore, the calculated order parameters differ from that using temperature dependent exchange parameters, which
the XRS measurements: in the model, the dominant order incorporate deviations from hexagonal symmetry proportional
parameters are Q̂x 2 −y 2 and Q̂xy rather than Q̂zx and Q̂xy . The to the order parameter may yield a better fit.
phase denoted AFQ1 in Fig. 4 has a large Q̂x 2 −y 2 moment The magnetoelastic strain is proportional to the strain
(ordered in antiphase along c) inducing a small (≈2%) Q̂zx derivative of the free energy [37],
moment. In the AFQ2 phase, Q̂xy moments become ordered,
 s αβ ∂F
inducing some Q̂yz quadrupoles; both these quadrupolar α = − , (7)
moments double below the transition to the AFQ3 phase. β
V ∂ β
Thus although the calculated sequence of ordering agrees
with experimental data, the type of quadrupolar order does where s αβ is the elastic compliance and the indices α and
 that F = −kB T ln Z,
not. Unfortunately, we found that it is impossible to stabilize β are the Cartesian directions. Noting
the Q̂zx order parameter over the Q̂x 2 −y 2 order parameter where the partition function is Z = n exp(−En /kB T ), we
at higher temperatures whilst maintaining a ground state
β =  ∂ β  where the angled brackets indicate
∂F ∂En
thus find that ∂
with all quadrupolar moments ordered. This may be due to the thermal expectation value.
limitations in the mean-field approximation or may indicate In the context of the Hamiltonian developed above, there
that cooperative effects between the considered quadrupoles are two main contributions to the magnetoelastic strain. These
are important, which might be modeled, for example, by arise from the single-ion crystal field and the two-ion exchange
exact diagonalization or Monte Carlo methods. Alternatively, interactions, both of which depend on the position of the ions,
couplings not considered here, for example, between first and are hence coupled to any change in the lattice. These
(q = |1|) and second (q = |2|) order quadrupoles such as terms two interactions give rise to crystal-field striction [38] and
of the form Q̂zx Q̂x 2 −y 2 , may be required. exchange striction [37], respectively. In order to calculate the
Another discrepancy between the model and experiments magnetostriction, we thus have to find expressions for the
is the very low critical fields when H ||c (seen at the bottom of energy levels En as a function of these two interactions. This
Fig. 4) and the magnitude of the c axis magnetization, which is is done by expanding the Hamiltonian (4) in a Taylor series
some three times smaller than measured (Fig. 3). In principle, and making the harmonic approximation by keeping only the
this can be altered by increasing the J10 exchange parameters, first-order term, which yields
however, we found that raising these from the values in
Table II suppresses the quadrupolar ordering completely,  s αβ ∂B q  q 
cfα = − k
Ok , (8)
in favor of a dipolar order. Alternatively, the quadrupolar V ∂
kq,i
exchange parameters J2m may be altered, but an increase in
the critical field necessitates also an increase in the transition 1  s αβ ∂Jij  ˆβ ˆβ 
temperatures.
α
ex = J J . (9)
2 β,ij V ∂ β i j
Since the exchange interactions are likely to arise from the
RKKY mechanism, one expects that it should be long ranged.
αβ ∂B k ij αβ ∂J
Thus including interactions further than nearest neighbor may The prefactors Aα = sV ∂q and Kαβ = sV ∂ βij may be
give better agreement with the data. However, this vastly taken to be independent of field and temperature and can then
increases the parameter space, and unfortunately we could be fitted to experimental data given the thermal expectation

235114-7
MANH DUC LE et al. PHYSICAL REVIEW B 89, 235114 (2014)

TABLE III. Fitted magnetoelastic parameters. model was then used to calculate the high-field magnetization
and magnetostriction, and the magnetic phase diagram up to
α 30 T. With an applied magnetic field in the basal plane, the
a b c calculations for the x and y directions agree with the data. The
calculated magnetization with the field parallel to z is a factor
Aα,xy × 10−5 26 −100 of 3 too large, however, and the calculated critical field in this
Aα,yz × 10−5 −35 85 direction is much smaller than that measured. Despite this, the
Aα,z2 × 10−5 −1 5 −30 z-direction calculated magnetic phase diagram qualitatively re-
Aα,zx × 10−5 80 −50 −30 produces the measured one. In conclusion, we have developed
Aα,x 2 −y 2 × 10−5 5 −60 −30 a nearest-neighbor mean-field model that reproduces many of
the main features observed in UPd3 .
ab
Kα,xy × 10−5 48 20
Kα,yz × 10−5
ab
−10 30
ab
Kα,zx × 10−5 −150 −45 25 ACKNOWLEDGMENTS
−5
2 −y 2 × 10
ab
Kα,x 1 10 10 The authors thank Jens Jensen, Gillian Gehring, Amanda
c
Kα,xy × 10−5 −48 −14 Gipson, and Helen Walker for fruitful discussions. M.D.L.
c
Kα,yz × 10−5 −3 −40 5
thanks the UK Engineering and Physical Sciences Research
Council for a research studentship, and UCL Graduate School
Kα,zx × 10−5
c
200 35
−5
for a project grant. J.G.P. acknowledges assistance by G.
2 −y 2 × 10 −4
c
Kα,x 8 Chouteau and the support provided by the Grenoble High
Magnetic Field Laboratory. M.D.L. and M.R. acknowledge
q the support of the European Science Foundation’s Cooper-
values Ok  and Jˆi Jˆj  obtained from the mean-field model at
ation in Science and Technology (COST) program, contract
different applied magnetic fields.
number COST-STSM-P16-02275. Experiments performed at
Considering only the nonzero terms up to rank 2 in the
the NHMFL were supported by NSF Cooperative Agreement
Hamiltonian (4), we thus have
⎡ DMR-0084173, by the state of Florida, and the US Department
1      ab  i i+(a+b)/2  of Energy. Support by the Deustche Forschungsgemeinschaft
α = ⎣ Aαβ Q̂iβ + Kαβ Q̂β Q̂β (SFB463), and the Austrian FWF (P17226) is gratefully
N i β β acknowledged.
⎤
  
+ c
Kαβ Q̂iβ Q̂i+c
β
⎦, (10) APPENDIX
β We derive in this appendix the form of the anisotropic
exchange interactions used in Eq. (4). In addition to the
where N = 4 is the number of U ions in the magnetic unit
4+
isotropic terms of the form,
cell, and the indices α = x,y,z, β = xy,yz,z2 ,zx,x 2 − y 2 .
In order to reduce the number of parameters in fitting   ij 
J k
q q
Hiso = Ô(i)k · Ô(j )k , (A1)
equation (10) to the data, we considered only the quadrupolar q
ij k
exchange interactions between nearest-neighbor ions in the c
direction and the basal plane, since the calculated expectation q
where the Ô(i)k denote the Stevens operators of the ith ion,
values of Q̂i Q̂j  Jˆi Jˆj . Table III shows the fitted param- we introduce local anisotropic couplings of the form
eters, whilst the calculated magnetostriction is shown as solid
lines in Fig. 5. Haniso
local

The calculated magnetostriction fits the data well in the ⎧

high-field regions, but does not reproduce the low-field   ⎨k
−q −q
K kq Ĉ(i)k · Ĉ(j )k + Ĉ(i)k · Ĉ(j )k
ij q q
=
behavior, particular for the case where the field is applied ⎩
ij k q=1
along x. It may be that at low fields, some of the measured ⎫
magnetostriction is due to domain rotation, which is not ⎬
+ K k0 Ĉ(i)0k · Ĉ(j )0k
ij
considered in the mean-field model. , (A2)
⎭
V. CONCLUSIONS q
where the spherical tensor operators Ĉk (which are equivalent
We have deduced a mean-field model including up to to spherical harmonics in a similar way to the equivalence of
four nearest-neighbor dipolar exchange and quadrupolar in- Stevens operators to tesseral harmonics) are defined as
teractions between the 5f 2 electrons of UPd3 , which is in
generally good qualitative agreement with a broad range −|q| 1  |q| −|q|
Ĉk = √ Ôk − i Ôk , (A3a)
of experimental results. The interactions between electrons 2
on the quasicubic sites were deduced from resonant x-ray Ĉk0 = Ôk0 , (A3b)
scattering measurements of the quadrupolar order parameters
of each of the four low-temperature ordered phases, and from +|q| (−1)|q|  |q| −|q|
Ĉk = √ Ôk + i Ôk . (A3c)
the measured transition temperatures and critical fields. The 2

235114-8
MEAN-FIELD MODEL FOR THE QUADRUPOLAR PHASES . . . PHYSICAL REVIEW B 89, 235114 (2014)

For the rank k = 1 case, these are the familiar angular mo- We note that this form of the anisotropy explicitly assumes a
mentum ladder operators Ĉ11 = −Jˆ+ , Ĉ1−1 = Jˆ− , and Ĉ10 = Jˆz . planar character. A more general expression may be obtained
In the same way, the rank k = 1 Stevens operators are the by considering anisotropy along the c direction, by using the
same as the Cartesian angular momentum operators Ô11 = Jˆx , local coordinate system defined by z ||rij [39], y  in the ab
Ô1−1 = Jˆy , and Ô10 = Jˆz . Thus the Cartesian form of Eq. (A2) plane and perpendicular to z and x  perpendicular to both y 
for rank k = 1 is and z . In this case, two rotations, by θij about y  to rotate x 
into the ab plane (and also making z ||c) and φij to rotate x  to
local  ij 
K 11 Jˆxi Jˆxj − Jˆyi Jˆyj + K 10 Jˆzi Jˆzj . become parallel to a are required to obtain the expressions in
ij
H1 aniso = (A4)
ij
the global coordinate system.
The rank k = 1 rotation matrix to accomplish this is [40]
The above anisotropic Hamiltonian (A2) is expressed in the ⎛ 1 −iφ −iφ
1 −iφ
⎞
e (1 + cos θ ) e √sin θ
e (1 − cos θ )
local frame defined by x  , y  , and z , where x  is perpendicular ⎜
2 2 2
⎟
R1 = ⎝ − sin
√θ cos θ sin
√θ ⎠,
to the hexagonal c axis, in the plane formed by the c axis and 2 2
eiφ√sin θ
the vector rij joining ions i and j ; y  is perpendicular to c and
1 iφ
2
e (1 − cos θ ) 2
1 iφ
2
e (1 + cos θ )
x  and z ||c. Thus to obtain the anisotropic couplings in the
(A5)
global coordinate system, one needs only to rotate Eq. (A2) by
the angle φij between the projection of rij into the ab plane where the subscripts ij has been omitted for clarity, and the
and the global x axis (defined as parallel to hexagonal a). rank k = 2 matrix is

⎛ 1 −2iφ ⎞
4
e (1 + cos θ)2 − 12 e−2iφ (1 + cos θ) sin θ 3 −2iφ
2
e sin θ 2 1 −2iφ
2
e (1 − cos θ) sin θ 1 −2iφ
4
e (1 − cos θ)2

⎜− 1 e−iφ (1 + cos θ) sin θ 1 −iφ 3 −iφ 1 −iφ 1 −iφ ⎟
⎜ 2 e (1 + cos θ)(−1 + 2 cos θ) e cos θ sin θ e (1 − cos θ)(1 + 2 cos θ) e (1 − cos θ) sin θ ⎟
⎜ 
2
√
2 2
√
2
 ⎟
R2 = ⎜
⎜
1 3
sin θ 2 − 6 cos θ sin θ 1
(−1 + 3 cos θ 2 ) 6 cos θ sin θ 1 3
sin θ 2 ⎟.
⎟
⎜ ⎟
2 2 2 2 2

⎝ − 1 eiφ (1 − cos θ) sin θ e (1 − cos θ)(1 + 2 cos θ)
1 iφ
− 23 eiφ cos θ sin θ e (1 + cos θ)(−1 + 2 cos θ)
1 iφ
e (1 + cos θ) sin θ ⎠
1 iφ
2 2 2 2

4
e (1 − cos θ)2
1 2iφ
− 12 e2iφ (1 − cos θ) sin θ 3 2iφ
2
e sin θ 2
2
e (1 + cos θ) sin θ
1 2iφ
e (1 + cos θ)2
1 2iφ
4

(A6)

Multiplying the tensor operators vector by this matrix then transforms the local anisotropic Hamiltonian into the global frame,
⎡ ⎛ ⎞⎤ ⎡ ⎛ −k ⎞⎤
K kk Ĉk−k (i)
ij
Ĉk (j )
⎢ ⎜ ⎟ ⎥ ⎢ ⎜ ⎟⎥
Haniso
global
= ⎢Rk ⎜ .. ⎟⎥ · ⎢Rk ⎜ .. ⎟⎥ . (A7)
⎣ ⎝ . ⎠⎦ ⎣ ⎝ . ⎠⎦
ij k
K kk Ĉkk (i)
ij
Ĉkk (j )

In the current case, we chose the local coordinates such that θij = 0 so the matrices Rk simplify to
⎛ ⎞
e−kiφ 0
⎜ .. ⎟
Rk = ⎝ . ⎠ (A8)
kiφ
0 e

and Eq. (A7) to

⎧ ⎫
  ⎨k
−q −q
⎬
K kq e2qiφij Ĉ(i)k · Ĉ(j )k + Ĉ(i)k · Ĉ(j )k + K k0 Ĉ(i)0k · Ĉ(j )0k
ij q q ij
Haniso
global
= . (A9)
⎩ ⎭
ij k q=1

Using the definitions (A3), one may obtain the Hamiltonian (4), which in matrix notation is
 ij
Hk = Jk O(i) ⊗ O(j ), (A10a)
ij
⎛ ⎞
Ôk−k (i)
⎜ . ⎟
O(i) = ⎜
⎝ .. ⎠ ,
⎟ (A10b)
Ôkk (i)

235114-9
MANH DUC LE et al. PHYSICAL REVIEW B 89, 235114 (2014)

where the exchange coupling matrices are, for rank k = 1,

⎛ ij ⎞
J11 − K+ 11 −K− 11
ij ij
0
⎜ ⎟
J1 = ⎜ ⎟,
ij ij
⎝ 0 J10 0 ⎠ (A11)
−K− 11 J11 + K+ 11
ij ij ij
0

where K± kq denote Kkq cos qφij and Kkq sin qφij , respectively, and the rank k = 2 matrices are
ij ij ij

⎛ ij ⎞ ⎛ ⎞
−K− 22
ij ij
J22 − K+ 22 −K− 22 Jxy
ij ij
0 0 0 0 0 0
⎜ ⎟ ⎜ ⎟
⎜ 0 J21 − K+ 21 0
ij ij
−K− 21
ij
0 ⎟ ⎜ 0 ij
Jyz 0 −K− 21
ij
0 ⎟
⎜ ⎟ ⎜ ⎟
ij ⎜ ⎟ ⎜ 0 ⎟
J2 = ⎜ 0 0 0 0 0 ⎟=⎜ 0 0 0 0 ⎟.
⎜ ⎟ ⎜ ⎟
⎜ −K− 21
ij
0 J21 + K+ 21
ij ij ⎟ ⎜ 0 −K− 21
ij
0
ij
Jzx 0 ⎟
⎝ 0 0 ⎠ ⎝ ⎠
−K− 22
ij ij
−K− 22 J22 + K+ 22 Jx 2 −y 2
ij ij ij
0 0 0 0 0 0
(A12)
The restrictions (5) and (6) imply that Jzxnnc < 0, Jxy
nnc
> 0, |Jzxnnc | > |Jxy
nnc
|, and Jzxnn < Jxy
nnc
< 0.

[1] P. Santini, S. Carretta, G. Amoretti, R. Caciuffo, N. Magnani, [18] A. Martin-Martin, Ph.D. thesis, University College London,
and G. H. Lander, Rev. Mod. Phys. 81, 807 (2009). 2000.
[2] F. Yakhou, V. Plakhty, H. Suzuki, S. Gavrilov, P. Burlet, L. [19] M. D. Le, K. A. McEwen, M. Rotter, J. Jensen, R. I. Bewley,
Paolasini, C. Vettier, and S. Kunii, Phys. Lett. A 285, 191 (2001). T. Guidi, and D. Fort, J. Phys.: Condens. Matter 24, 036002
[3] H. Nakao, K. Magishi, Y. Wakabayashi, Y. Murakami, K. (2012).
Koyama, K. Hirota, Y. Endoh, and S. Kunii, J. Phys. Soc. Jpn. [20] D. J. Newman and B. K. C. Ng, Crystal Field Handbook
70, 1857 (2001). (Cambridge University Press, Cambridge, UK, 2000).
[4] T. Onimaru, T. Sakakibara, T. Tayama, N. Aso, H. Y. amd Dai [21] G. A. Gehring (private communication).
Aoki, Y. Ōnuki, T. Kawae, T. Kitai, and T. Takeuchi, Physica B [22] W. J. L. Buyers, A. F. Murray, T. M. Holden, E. C. Svensson, P.
359–361, 935 (2005). de V. DuPlessis, G. H. Lander, and O. Vogt, Physica B+C 102,
[5] T. Nagao and R. Shiina, J. Phys. Soc. Jpn. 79, 034703 (2010). 291 (1980).
[6] Y. Tanaka, T. Inami, S. W. Lovesey,, K. S. Knight, F. Yakhou, [23] K. R. Lea, M. J. M. Leask, and W. P. Wolf, J. Phys. Chem. Solids
D. Mannix, J. Kokubun, M. Kanazawa, K. Ishida, S. Nanao, 23, 1381 (1962).
T. Nakamura, H. Yamauchi, H. Onodera, K. Ohoyama, and Y. [24] K. A. McEwen, U. Steigenberger, K. N. Clausen, J. Kulda, J.-G.
Yamaguchi, Phys. Rev. B 69, 024417 (2004). Park, and M. B. Walker, J. Magn. Magn. Mater. 177–181, 37
[7] S. W. Zochowski, M. de Podesta, C. Lester, and K. A. McEwen, (1998).
Physica B 206&207, 489 (1995). [25] M. Rotter, H. Müller, E. Gratz, M. Doerr, and M. Loewenhaupt,
[8] K. A. McEwen, J.-G. Park, A. J. Gipson, and G. A. Gehring, Rev. Sci. Instr. 69, 2742 (1998).
J. Phys.: Condens. Matter 15, S1923 (2003). [26] S. W. Zochowski and K. A. McEwen, Physica B 199&200, 416
[9] Y. Tokiwa, K. Sugiyama, T. Takeuchi, M. Nakashima, R. Settai, (1994).
Y. Inada, Y. Haga, E. Yamamoto, K. Kindo, H. Harima, and [27] M. Rotter, J. Magn. Magn. Mater. 272–276, E481 (2004).
Y. Onuki, J. Phys. Soc. Jpn. 70, 1731 (2001). [28] M. Rotter et al., http://www.mcphase.de
[10] U. Steigenberger, K. A. McEwen, J. L. Martinez, and D. Fort, [29] R. A. Cowley and J. Jensen, J. Phys.: Condens. Matter 4, 9673
J. Magn. Magn. Mater. 108, 163 (1992). (1992).
[11] D. F. McMorrow, K. A. McEwen, U. Steigenberger, H. M. [30] J. G. Houmann, B. D. Rainford, J. Jensen, and A. R. Mackintosh,
Rønnow, and F. Yakhou, Phys. Rev. Lett. 87, 057201 (2001). Phys. Rev. B 20, 1105 (1979).
[12] H. C. Walker, K. A. McEwen, D. F. McMorrow, S. B. Wilkins, [31] J. Jensen and A. R. Mackintosh, Rare Earth Magnetism
F. Wastin, E. Colineau, and D. Fort, Phys. Rev. Lett. 97, 137203 (Clarendon Press, Oxford, 1991).
(2006). [32] For a two-level system, the mean-field susceptibility is inversely
[13] H. C. Walker, K. A. McEwen, M. D. Le, L. Paolasini, and proportional to ( − 2nsd |d|Ô|s|2 J (q)), which may be
D. Fort, J. Phys.: Condens. Matter 20, 395221 (2008). equated with the Curie-Weiss susceptibility, χ ∝ 1/(T − TN ).
[14] T. J. Heal and G. I. Williams, Acta. Cryst. 8, 494 (1955). At low temperatures, T , the thermal population factor
[15] H. H. Hill, Nucl. Metall. 17, 2 (1970). nsd ≈ /kB T , so that TN ∼ ∝ |d|Ô|s|2 J (q).

[16] K. T. Moore and G. van der Laan, Rev. Mod. Phys. 81, 235 [33] J. Kennedy and R. Eberhart, in Proceedings of IEEE In-
(2009). ternational Conference on Neural Networks, Vol. 4 (IEEE,
[17] K. A. McEwen, H. C. Walker, M. D. Le, D. F. McMorrow, E. Piscataway, NJ, 1995), pp. 1942–1948.
Colineau, F. Wastin, S. B. Wilkins, J.-G. Park, R. I. Bewley, and [34] K. A. McEwen, U. Steigenberger, and J. L. Martinez,
D. Fort, J. Magn. Magn. Mater. 310, 718 (2007). Physica B 186–188, 670 (1993).

235114-10
MEAN-FIELD MODEL FOR THE QUADRUPOLAR PHASES . . . PHYSICAL REVIEW B 89, 235114 (2014)

[35] N. Lingg, D. Maurer, V. Müller, and K. A. McEwen, Phys. Rev. [38] P. Morin and D. Schmitt, in Ferromagnetic Materials, edited
B 60, R8430 (1999). by K. H. J. Buschow and E. P. Wohlfarth, Vol. 5 (Elsevier,
[36] J. Fernández-Rodrı́guez, S. W. Lovesey, and J. A. Blanco, Amsterdam, 1990), pp. 1–132.
J. Phys.: Condens. Matter 22, 022202 (2010). [39] J. Jensen (private communication).
[37] M. Doerr, M. Rotter, and A. Lindbaum, Adv. Phys. 54, 1 [40] H. A. Buckmaster, R. Chatterjee, and Y. H. Shing, Phys. Stat.
(2005). Sol. (a) 13, 9 (1972).

235114-11