Table Of ContentDipolar-broadening of I = 1-NMR spectra of solids
2
Jens Jensen
Ørsted Laboratory, Niels Bohr Institute, Universitetsparken 5, 2100 Copenhagen, Denmark
5
(February 1, 2008)
9
9
1
n
a
J Abstract
3
1
The equations of motion of dipolar-coupled spins of I = 1/2 placed on a
3
v
9 rigid lattice are solved approximately in the high-temperature and high-field
3
0
limit. The NMR-spectra predicted by this theory are in close agreement
1
0
5 with both the theoretical spectral moments, up to the eighth in the case of a
9
/ simple cubic lattice, and the extremely accurate experimental results which
t
a
m
have been obtained in the case of CaF . The theory is compared with the
2
-
d
recent NMR-experiment of Lefmann et al. on 13C diamond. It predicts a
n
o
c double-peak splitting of the NMR-spectrum when the field is applied along
:
v
i [111], in accordance with the experiment, though the widthsof the calculated
X
r
a resonance lines are smaller than observed.
76.20.+q 76.60.-b
Typeset using REVTEX
1
I. INTRODUCTION
Recently Lefmann et al.1 have measured the NMR-spectra on single crystals of the new
material 13C diamond. The nuclear spin of the 13C isotope is 1/2, and in the case of a
field applied along [111], the square of the magnetic dipole–dipole interaction between two
nearest-neighboring spins in the direction of the field is an order of magnitude stronger than
the interaction between any other neighbors. This leads to a splitting of the resonance into
two peaks, a phenomenon which has only previously been observed in dense systems of
molecules,2–4 The 13C diamond constitutes a system of I = 1/2 nuclei interacting through
the dipole force between their magnetic moments and, because the Debye temperature of
the diamond crystal is high, the lattice may be considered to be rigid. The same model
applies for the NMR in CaF except that the I = 1/2 nuclei, the F-ions, are placed on a
2
simple cubic lattice.
The simplest problem of NMR in solids is the one where identical nuclei of I = 1/2 are
placed on a rigid lattice. In this case no quadrupole moment or quadrupolar broadening
effects have to be considered, and lattice vibrations are neglected. The calculation of the
dipolar broadening of the NMR is, nevertheless, a complicated many-body problem with no
small expansion parameter. Van Vleck calculated the second and fourth moments of the
spectral distribution,5 which indicate that the line shape is neither Lorentzian nor Gaussian.
The real-time response function F(t), the so-called free-induction-decay (fid) curve, derived
from experiments6 on CaF , was found by Abragam7 to be accurately described by the
2
function
F(t) = e−a2t2/2sin(bt)/bt, (1.1)
when the parameters a and b were adjusted so to get the right second and fourth moments
of the spectral function. In the early seventies the calculation of the moments was extended
to include the sixth and eighth,8 and the method of Abragam was generalized in order to
account for the values of these higher moments.4,9 Experimentally, the fid curves of the stan-
dardNMR-systemCaF weredeterminedwithveryhighprecisionbyEngelsbergandLowe.10
2
2
Important progress in determining the response function from the equations of motion was
made by Becker et al.11 and by Lundin.12,13 The method of Lundin, which has recently been
developed further by Shakhmuratov,14 is reminiscent of the simple physical interpretation
made by Lefmann and collaborators of the double peak structure in 13C diamond, where the
splitting is considered to be due to the strongly-coupled pairs of nearest neighbors, and the
width to a Gaussian-like distributed random-field from other neighbors.1,15
In the present work the equations of motion of the Green function are analyzed. One
main result is that a solution of these equations can be found when only the z-components
of the spins are coupled. The Green function in this case is determined by an infinite
continued-fraction, which may be solved accurately by numerical methods. This result is
then generalized to include the effects of the x- and y-components of the dipole–dipole in-
teraction. In this step, we utilize the analytic results obtained previously for the higher
moments. A numerical analysis is then carried through in the cases corresponding to CaF
2
and 13C diamond. The results obtained by the present method agree closely with the ex-
perimental results obtained for CaF . In the case of 13C diamond there are discrepancies,
2
which are probably mostly of experimental origin.
II. EQUATIONS OF MOTION
We want to consider the linear magnetic response of a lattice of nuclear moments sub-
jected to a large stationary field. The normal situation in the case of a NMR experiment
is that kT is much larger than the mean-field splitting between the nuclear levels, which is
again much larger than the energy shifts due to the coupling between the nuclear moments.
These conditions we shall assume to be satisfied throughout. We shall only treat the case
of I = 1/2 nuclei placed in a non-metallic crystal. In this case we need only to consider the
classical coupling of the nuclear moments and, neglecting non-linear effects, the Hamiltonian
may effectively be reduced to7,16
3
∆ 1
H = − σz − J (ij)s+s− − J (ij)σzσz. (2.1)
2 i x i j 4 z i j
X X X
i ij ij
∆ = gµ H is the Zeeman splitting between the I = ±1 states, g is the nuclear g-factor,
N z z 2
and the gyromagnetic ratio γ = gµ /~. The operators s± = I ±iI and σz = 2I . An
N j xj yj j zj
effective dipole-coupling parameter is defined as
3(1−3cos2θ )
D = ij (gµ )2, (2.2)
ij 4R3 N
ij
where R is the distance between the ith and jth nucleus, and θ is the angle between
ij ij
the line connecting the two nuclei and the z-axis. The contribution of the dipole–dipole
interaction is then included in Eq. (2.1) via
J (ij) = 1D ; J (ij) = −2D . (2.3)
x 3 ij z 3 ij
The frequency-dependent Green functions are defined as
∞
G(Aˆ) = lim (−i)θ(t)h[Aˆ(t), s−]iei(ω+iη)tdt, (2.4)
η→0+Z−∞ ~ Pj j
where θ(t) is the step function, and the equation-of-motion17,18 is
~ωG(Aˆ)−G([Aˆ,H]) = h[Aˆ, s−]i. (2.5)
j j
P
Due to our general assumptions, we can neglect the dipole coupling in the determination of
the thermal average on the r.h.s. of this equation, and only h[s+, s−]i = hσzi ≃ ∆/kT
i j j i
P
is non-zero to leading order in 1/kT. The linear response observed in NMR experiments
is proportional to the imaginary part of the Green function G(s+). Utilizing the special
i
properties of the spin-1/2 operators, e.g. (σz)2 = 1, Eqs. (2.1) and (2.5) lead to the following
i
equation:
εG(s+) = −hσzi+ [J (ij)G(σzs+)−J (ij)G(s+σz)] (2.6)
i i x i j z i j
X
j
with ε = ∆−~ω. The next-order Green functions are determined by (i 6= j)
4
εG(s+σz) = J (ij)G(s+)−J (ij)G(s+)
i j x j z i
+ ′ J (ik)G(σzσzs+)−J (ik)G(s+σzσz)
x i j k z i j k
X h
k
+2J (jk){G(s+s+s−)−G(s+s−s+)} , (2.7)
x i j k i j k
i
where the primed sum means that the summation index (k) has to be different from the
other indices (i or j). Next we find, when all indices i, j, and k are different from each other,
εG(s+σzσz) = J (ik)G(σzs+)+J (ij)G(s+σz)
i j k x j k x j k
−J (ik)G(s+σz)−J (ij)G(s+σz)
z i j z i k
+ ′ J (il)G(σzσzσzs+)−J (il)G(s+σzσzσz)
x i j k l z i j k l
X h
l
+2J (jl){G(s+s+σzs−)−G(s+s−σzs+)}
x i j k l i j k l
+2J (kl){G(s+σzs+s−)−G(s+σzs−s+)} (2.8)
x i j k l i j k l
i
and
2εG(s+s+s−) = J (ik){G(σzs+)−G(szσz)}
i j k x i j j k
+J (jk){G(s+σz)−G(s+σz)}
x i j i k
+2 ′ J (il)G(σzs+s−s+)
x i j k l
X h
l
+J (jl)G(s+σzs−s+)−J (kl)G(s+s+σzs−)
x i j k l x i j k l
−{J (il)+J (jl)−J (kl)}G(s+s+s−σz) . (2.9)
z z z i j k l
i
The major complications in the equations of motion above are due to J (ij). The Green
x
function G (s+), determined by these equations when J (ij) ≡ 0, only depends on higher-
z i x
order Green functions G (s+σzσz···) involving one s+ operator, and in the next order:
z i j k
εG (s+σzσzσz) = −J (il)G (s+σzσz)
z i j k l z z i j k
−J (ik)G (s+σzσz)−J (ij)G (s+σzσz)
z z i j l z z i k l
− ′J (im)G (s+σzσzσzσz ), (2.10)
z z i j k l m
X
m
5
where the first three terms all contribute the same to G (s+) and thus may effectively
z i
be replaced by −3J (il)G (s+σzσz). Similarly, the two terms in front of the summation
z z i j k
sign of Eq. (2.8) may be replaced by −2J (ik)G (s+σz). Hence the hierarchy of Green
z z i j
functions determining G (s+) has a transparent structure. Introducing a vector with the
z i
components G (s+), J (ij)G (s+σz), J (ij)J (ik)G (s+σzσz), etc. the equations of motion
z i z z i j z z z i j k
may be written as an infinite-dimensional matrix equation, which has the solution of an
infinite continued-fraction
−hσzi
G (sz) = i (2.11)
z i
J2(ij)
z
ε−
Xj 2J2(ik)
′ z
ε−
Xk 3J2(il)
′ z
ε−
ε−Qz(ijkl)
Xl 4
with
4J2(im)
Qz(ijkl) = ′ z (2.12)
4 ε−Qz(ijklm)
Xm 5
etc. Notice that Qz(ijkl) depends implicitly on the indices j, k, and l, because the summa-
4
tion index m has to be different from all the preceding ones, as indicated by the prime.
The conjecture is now that the J (ij) coupling may be included as a perturbation, so
x
that the result may still be written like the infinite continued-fraction above, except that
the dependence on the coupling parameters is altered, i.e. Qz in (2.11) is replaced by
p
X (ij···)
′ p
Q (ij···) = , (2.13)
p ε−Q (ij···)
X p+1
where the parameters X (ij···) are determined by the equations of motion when the J (ij)
p x
coupling is included. Expanding the infinite continued-fraction in powers of 1/ε, we get
hσzi 1 1
G(s+) = − i 1+ ′X (ij)+ ′X (ij)
i ε h ε2 X 1 ε4 X 1
j j
′ ′ ′
×{ X (ij )+ X (ijk)}+··· . (2.14)
1 2
Xj′ Xk i
6
Assuming all nuclei to be placed in equivalent surroundings, i.e. G(s+) = G(s+), then we
i j
get from Eqs. (2.6)–(2.9)
hσzi 1 1
G(s+) = − i 1+ D2 + D2 (2.15)
i ε h ε2 X ij ε4 X ij
j j
1
×{ D2 + ′(7D2 +5D2 +6D D )}+···
Xj′ ij′ 9 Xk ik jk ik kj i
utilizing that the summation indices may be interchanged (the prime on the j-sums is
unnecessary as D = 0). A comparison of (2.14) and (2.15) then allows the identification of
ii
X and X . The distinction between D2 and D2 in the last sum is arbitrary, and because
1 2 ik jk
G(s+) = (1/N) G(s+) the result should be same if the sum over j, in the first step of the
i i i
P
infinite continued-fraction, is replaced by a sum over i. In order to ensure this, also when
introducing the higher-order terms, it is sufficient to assume the solution to be symmetric
in i and j:
X (ij) = D2 (2.16a)
1 ij
2
X (ijk) = {D2 +D2 +D D }. (2.16b)
2 3 ik jk ik kj
The moments of the spectral energy-distribution, defined with respect to the mean value ∆,
are
∞ ∞
M = hεni = εnG′′(s+;ε)dε G′′(s+;ε)dε, (2.17)
n Z−∞ i .Z−∞ i
where G′′(s+;ε) denotes the imaginary part of the Green function (M is (2π~)n times the
i n
usual frequency moments). The oddmoments vanish by symmetry, andanelementary result
is
hσzi 1 1
G(s+) = − i 1+ M + M +··· , (2.18)
i ε h ε2 2 ε4 4 i
which shows that (2.15) is equivalent to an expansion of G(s+) in terms of the spectral
i
moments. Introducing the notation S = Dn, the moments are
n j ij
P
M = S (2.19a)
2 2
1
M = {7S2 −4S +2S D D } (2.19b)
4 3 2 4 2 ik kj
X
k
7
in accordance with previous results.5,8
The complexity inperforming theexpansion (2.15)ofG(s+) increases dramaticallyinthe
i
next order. Therefore we make use of the additional assumption that terms which involve
products of two different D-couplings can be neglected in X , when p ≥ 3. The coupling
p
D takes on both positive and negative values and D integrated over a sphere vanishes,
ij ij
whereas in contrast D2 is always positive. Advancing one more step in the hierarchy of
ij
Green functions, beyond that of Eqs. (2.8) and (2.9), and keeping track of only the squared
couplingterms, wefindthelinearcombinationX (ijkl) = aD2+bD2+cD2 tobedetermined
3 il jl kl
by a+b+c = 41/18, and further that a solution symmetric in i and j (a = b) implies c = 0,
thus
41
X (ijkl) = {D2 +D2}. (2.20)
3 36 il jl
This result predicts the sixth moment M to be
6
1
{229S3 −373S S −171S +41 D2D2 D2 }
27 2 2 4 6 ij ik kj
X
jk
plusthetermsderiving fromD D inX (ijk), which isconsistent withtheanalyticalresult
ik kj 2
for M derived by Knak Jensen and Kjærsg˚ard Hansen.8 It is possible to include most of the
6
terms in M depending on the linear factors, but the price is that Q (ij) has to be divided
6 2
into at least two continued fractions, which makes the procedure somewhat arbitrary. More
importantly, in principle, is the occurrence8 of the term −(8/27) D D D D D D in
ij ik il lk kj jl
P
M . It falls outside the present scheme and cannot be incorporated into the continued
6
fraction in any simple way, thus indicating a limit to the present procedure. However, as
discussed above, the contributions to M which are neglected in Eq. (2.20) are expected to
6
be small, and this is supported by the numerical analysis discussed below.
In the higher order it is not important to discriminate between the different quadratic
contributions, e.g. between D2 , D2 , ··· in X (ijklm). In any case, such a separation plays
im jm 4
no role in the following numerical analysis, and we may assume,
α
X (ijk···γ) = p{D2 +D2 } ; p ≥ 4 (2.21)
p 2 iγ jγ
8
corresponding to Eq. (2.20). Using the result of Knak Jensen and Kjærsg˚ard Hansen8 that
the leading order term in M is (11031/243)S4, we find α = 2819/738. Introducing an
8 2 4
effective α for the lower-order terms also, we get the following sequence of values for α :
p p
4 41 2819
1, , ,
3 18 738
for p = 1, 2, 3, and 4. The final Green function is not much dependent on how the series
is continued, except that α increases with p (which is guarantied by the J (ij)-coupling).
p z
Althoughthe fourmoments M –M do not determine the response, theyput strong limitson
2 8
the kind of variation which may be achieved by varying the higher-order terms. However,
our ignorance of how the α series continues leaves us a single degree of freedom, which
p
we have utilized for a small over-all adjustment of the calculated response functions, using
for this purpose the accurate experimental results for CaF . We have tried to use a linear
2
extrapolation (with a slope between 1.15 and 1.2), but the α -values above are described
p
closely by α = 1 + 0.033(p − 1) + 0.302(p − 1)2, i.e. α seems to increase quadratically
p p
rather than linearly with p, and we obtained a slightly better result by using the following
quadratic extrapolation:
α = 1+0.327(p−1)2 ; p ≥ 5. (2.22)
p
In the numerical evaluation of the Green function we have to make some additional
approximations. The infinite continued-fraction
−hσzi
G(sz) = i (2.23)
i
D2
ij
ε−
ε−Q (ij)
Xj 2
becomes less and less dependent on the actual values of Q (ij···) the larger p is, and to a
p
first approximation we assume that the summation index only needs to be different from i
and j. In this approximation ′X (ijk···γ) ≃ α V(ij) with V(ij) = S −D2 for p ≥ 3,
p p 2 ij
P
hence
9
4V(ij)+ 2 D D
Q (ij) ≃ 3 3 k ik kj . (2.24)
2 P
α V(ij)
3
ε−
α V(ij)
4
ε−
ε− ...
The infinite continued-fraction has to be determined for each value of the summation index
j. This is done by assuming α to be constant for p ≥ n, i.e. Q (ij) ≃ α V(ij)/{ε−Q (ij)}
p n n n
or
1
Q (ij) = ε− ε2 −α V(ij) .
n 2n p n o
The minus sign in front of the square root is the only choice which leads to the right sign of
G′′(s+). Even for rather large values of n, this termination of the infinite continued-fraction
i
gives rise to spurious oscillations in the calculated response function (as function of n or ε),
which are, however, foundtocanceloutinthemeanvalueoftheresponsefunctioncalculated
for two successive values of n. With the use of this averaging procedure it is found that the
result becomes independent of n, when n is larger than about 20–30 (we have used n = 50
and 51 in the final calculations). The approximation made in Eq. (2.24) for Q (ij) may be
2
improved, and to next order we get
2 D2 +D2 +D D
′ ik jk ik kj
Q (ij) = (2.25)
2 3
X
k α W(ijk)
3
ε−
α W(ijk)
4
ε−
ε− ...
withW(ijk) = S −D2 −(D2 +D2 )/2. Theeffective coupling parameterD2 isproportional
2 ij ik jk ij
to R−6, but although this coupling decreases much faster than the dipole coupling itself,
ij
the j sum has to be extended over 10–20,000 neighbors in order to obtain an acceptable
accuracy. Fortunately, the modifications introduced by the replacement of Eq. (2.24) by Eq.
(2.25) only influence the final response function weakly, and it is only necessary to apply
the improved expression for Q (ij) when j is one of the 20–30 neighbors which are coupled
2
most strongly to the ith site. The result (2.25) may be considered to correspond to the
10
