**Pseudopotential hyperfine calculations through perturbative core-level polarization**

Mohammad Saeed Bahramy,1,*Marcel H. F. Sluiter,2 and Yoshiyuki Kawazoe1
1_{Institute for Materials Research, Tohoku University, Sendai 980-8577, Japan}

2* _{Department of Materials Science and Engineering, Delft University of Technology, Mekelweg 2, 2628CD Delft, The Netherlands}*
共Received 23 November 2006; revised manuscript received 12 February 2007; published 30 July 2007兲

Within density functional theory, an efficient and accurate method for calculating the hyperfine parameters
in the context of pseudopotential formalism is proposed. The spin density at and in the vicinity of the nucleus
is evaluated in two steps. First, a transformation due to Blöchl**关Phys. Rev. B 50, 17953 共1994兲兴 is applied to**
reconstruct the frozen-core all-electron wave functions in the core regions. Second, the contributions of core
orbitals to the charge density at the nucleus are evaluated through first-order perturbation theory in which the
perturbing potential is defined as a functional of charge and spin densities. The current pseudopotential based
method makes it possible to predict hyperfine parameters of complex molecular assemblies and crystal defects
with an accuracy as good as current all-electron method with less computational cost.

DOI:10.1103/PhysRevB.76.035124 PACS number共s兲: 71.15.⫺m, 31.15.Ew, 31.30.Gs, 32.10.Fn

**I. INTRODUCTION**

Much of our knowledge on the geometrical and electronic
structure of molecules and clusters as well as the pointlike
defects in crystals has been obtained through experimental
techniques such as electron paramagnetic resonance and
nuclear magnetic resonance.1–3 _{Common to all these }
tech-niques is the measurement of the so-called hyperfine
struc-ture, which describes the interaction of the electronic wave
functions with the nuclear magnetic moments. In general, the
hyperfine structure is decomposed into the isotropic
共Fermi-contact兲 parameter A*iso* and anisotropic 共dipolar兲 parameter

*Aaniso*. The former is directly proportional to the spin density

*at the nucleus, arising from the spin-polarized s-like orbitals,*
while the latter has contributions from the other
*spin-polarized states, e.g., p- or d-like orbitals, and is proportional*
to the volume integral over the spin density divided by the
distance to the nucleus to the third power. In actuality,
*mo-lecular motion causes Aaniso* to average out to zero so that

*only Aiso*can be observed in the spectra of gases and liquids.

Although the theory underlying the hyperfine structure is
well understood and was already developed in the early days
of quantum mechanics,1,2 _{the first-principles calculations of}
hyperfine parameters共HFPs兲 have proven to be a challenging
task.4,5 The main challenge is the necessity of an accurate
representation of the spin density at and in the vicinity of the
nucleus. Already a half century ago, works based on the
un-restricted Hartree-Fock method6,7 _{revealed that the }
spin-polarized core states contribute substantially to the spin
*den-sity at the nuclei of system with non-s-like singly occupied*
molecular orbitals共SOMOs兲.8–10_{Therefore, contributions of}
both core and valence states in the vicinity of the nucleus
must be accurately accounted for.

Among the computational approaches, all-electron 共AE兲
methods based on density functional theory共DFT兲,11_{such as}
the linearized muffin-tin orbital,12 _{the full-potential }
linear-ized augmented-plane-wave,13 _{and the all-electron }
mixed-basis共AEMB兲14 _{methods have proven to give accurate }
rep-resentations of the spin density near the nucleus. However,
AE methods are computationally less efficient than
plane-wave pseudopotential 共PP兲 methods. Additionally, new

de-velopments in DFT are typically more easily coded within
PP methods than in AE methods. Therefore, already some
time ago it was realized that it is desirable to be able to
compute HFPs with PP methods.15_{However, because of the}
inherent shortcoming of PP methods in approximating the
electronic wave functions within the core region, PP methods
cannot be used in their original form for the calculation of
HFPs. The difficulties originate from共1兲 the nodeless
behav-ior of pseudo-wave-functions in the core region, resulting in
incorrect spin densities in the vicinity of the nucleus, and共2兲
the complete elimination of core orbitals and, consequently,
the complete elimination of the core spin polarization and
corresponding contributions to the HFPs.

So far, there have been several efforts in developing a
pseudopotential treatment for the calculation of HFPs.
How-ever, most of these efforts have been concentrated on the first
challenge, that is, the reconstruction of the full-nodal form of
AE wave functions.16–18 _{In this regard, the so-called}
projector-augmented-wave 共PAW兲 method proposed by
Blöchl19 * _{has shown to yield highly accurate frozen-core}*
wave functions in the core regions and hence is considered as

*a reliable method for calculating HFPs of systems with s-like*SOMO, e.g., hydrogen defects in Si.15

_{Another important}aspect of HFPs computed through PP formalism has only very recently received attention: the effect of core spin

*po-larization was considered by Yazyev et al.*20

_{and Declerck et}*al.*21_{The former}20_{is similar in spirit to the current work, but}
approaches the treatment of the core spin polarization
differ-ently. In the so-called core spin-polarization correction
共CSPC兲 method,20_{a reconstruction of the AE wave functions}
and the frozen-valence spin-density approximation are used
to solve the Kohn-Sham equations for core states only.
De-spite the simplicity of the method, it seems to be quite
*accu-rate for calculating Aiso* of simple molecular radicals,

the remaining atoms. Nontheless, this method appears sensi-tive to the choice of basis set.

Here, we present a general and accurate method for cal-culating HFPs, within the DFT pseudopotential formalism. The approach, which we call perturbative core-level polar-ization共PCLP兲, is based on the evaluation of HFPs by taking into account the contributions of both core and valence spin-polarized orbitals to the spin density at and in the vicinity of the nucleus. In this regard, we first reconstruct the AE wave functions from the corresponding pseudo 共PS兲-wave-functions with the assumption that core electrons are unper-turbed by the presence of other atoms, also known as the frozen-core approximation. Then, the reconstructed AE wave functions are utilized to develop a local spin-density func-tional potential. Finally, using this perturbing potential, the spin-polarized core levels are estimated by means of first-order perturbation theory.

It is worth mentioning that other perturbation treatments
of core spin polarization for the computation of isotropic
HFPs have been derived within the context of Hartree-Fock
共HF兲 formalism.22,23 _{These studies are mainly concentrated}
on simple atomic and molecular models. This is
understand-able because of the inherent neglect of electron-electron
cor-relation effects and the leading spin contamination errors in
HF approximations make their utilization a nontrivial task
*for larger systems with, e.g., d-type atoms. In contrast, in the*
current DFT-based approach, these problems are essentially
negligible and of little importance.

**II. FORMULATION**

*The hyperfine parameters Aiso* *and Aaniso*are defined as

*Aiso共I兲 =*
2
30*ge**egI**I**s***共R***I*兲 共1兲
and
*Aaniso共I兲 =*
1
40*ge**egI**I*

### 冕

**

**dr**I*s*

**共r**

*I*兲 3 cos2− 1

*2rI*3 , 共2兲

**where R**

*I*

*is the positional*

**is the position of nucleus I, r**Ivector relative to the nucleus**共r***I***= r − R***I*兲, is the angle

**be-tween r***I* and symmetry axes, 0 is the permeability of

vacuum共4⫻10−7_{T}2_{m}3_{J}−1_{兲, g}

*eis the electron g factor,**e*

*is the Bohr magneton, and gI* and *I* are the gyromagnetic

ratio and the magnetic moment of the nucleus. Throughout
*this work, gI*and*I* values are taken from Ref.24. The spin

density*s***共r兲 is the difference between the up and **

spin-down charge densities, *s***共r兲=***↑***共r兲−***↓***共r兲. Using the **

nota-tion for spin sign *共↑ or ↓兲,* **共r兲 can be considered as a**
sum of the core contribution *c*_{共r兲 and the valence }

contri-bution*v*_{共r兲,}

_{共r兲 =}_{}*c*_{共r兲 +}_{}*v*_{共r兲.}_{共3兲}

In the conventional PP formalism, the core contribution
*c*_{共r兲 is not available. Additionally, the valence contribution}

is expressed as
*˜v*_{共r兲 =}

### 兺

*n*具⌿

*˜*

*n*

_{兩r典具r兩⌿}*˜*

*n*

_{典,}

_{共4兲}

where the PS-wave-functions兩⌿*˜ _{n}*典 differ from the exact AE
wave functions 兩⌿

*n*典 within the core regions. Thus, any

at-tempt to implement HFP calculations into PP methods must
consider共1兲 the inclusion of*c*_{in the charge density and}_{共2兲}

the reconstruction of 兩⌿*n*典 from 兩⌿*˜n*典 inside the atomic

spheres.

Below, we first discuss the second task. The兩⌿*n*典 is

re-constructed from兩⌿*˜ _{n}*典 using the PAW19

_{method. In this }con-text, 兩⌿

*n*典 are obtained by applying the following

transfor-mation on兩⌿*˜n*典:
兩⌿*n*典 = 兩⌿*˜n*典 +

### 兺

共兩典 − 兩*˜ 典兲具p*

*˜ 兩⌿*

*˜*

*n*

_{典,}

_{共5兲}

where AE partial waves 兩_{}典 are the solution of the radial
Schrödinger equation and PS-partial-waves兩*˜ 典 are smooth*_{}
functions that coincide with the corresponding AE partial
*waves outside the cutoff radius Rc*. The index is a

*short-hand denoting the particular angular momenta l and m of*
level**, localized at the atomic site R***I*. The functions*兩p˜ 典 are*

projector functions which are well localized inside the core region and defined such that they fulfill the condition

*具p˜ 兩**˜ 典 =* ␦. 共6兲

To determine *兩p˜ 典, we follow the approach proposed by*_{}
*Hetényi et al.*18 _{The procedure is briefly discussed below.}
Initially, we perform a simple AE calculation共e.g., using the
Herman-Skillman method25_{兲 for an isolated atom in the}
ground state in order to compute the partial waves兩_{}典 and
the corresponding eigenvalues_{}as well as the radial
*effec-tive potential Vef f共r兲,*

## 冉

−ⵜ2

2 *+ Vef f*

## 冊

兩典 = 兩典. 共7兲*Next, an arbitrary potential Vloc共r兲 is constructed such that it*

*behaves smoothly inside and coincides with Vef f共r兲 outside*

the core region. A possible choice is to find an exponential
*function for Vloc共r兲, that joins smoothly to Vef f共r兲 at Rc*共see

Ref.19*兲. Using Vloc共r兲, *_{}, and PS-partial-waves兩*˜ 典, a set*_{}

of functions, localized within the core region, is defined as

兩*˜ 典 =*

## 冉

+ⵜ 22 *− Vloc*

## 冊

兩*˜ 典.* 共8兲 Finally, taking a linear combination of 兩

*˜ 典, the projector*

_{}functions

*兩p˜ 典 are constructed as*

_{}

*兩p˜ 典 =*

### 兺

*共B*

−1_{兲}

兩*˜ 典,* 共9兲

*where the matrix elements B*_{} are defined so that the
com-puted*兩p˜ 典 fulfills the orthogonality condition 关Eq. 共*_{} 6兲兴,

Once the projector functions*兩p˜ 典 are determined, transfor-*_{}
mation共5兲 can be applied to obtain 兩⌿*n*典. Similarly, the

va-lence contribution *v*_{共r兲 is reconstructed from}_{˜}_{}*v*_{共r兲 as}

follows:19
*v*_{共r兲 =}_{˜}_{}*v*_{共r兲 +}_{1}_{共r兲 −}_{˜}_{}
1
_{共r兲,}_{共11兲}
where
1_{共r兲 =}

### 兺

*n*

### 兺

具⌿*˜*

*n*

*典具*

_{兩p˜}**兩r典具r兩**

*典具p˜*兩⌿

*˜n*典 共12兲 and

*˜*

_{1}

**共r兲 =**

### 兺

*n*

### 兺

具⌿*˜*

*n*

*典具*

_{兩p˜}*˜*

**兩r典具r兩**

*˜*

*典具p˜*兩⌿

*˜n*典. 共13兲

Having*v*_{共r兲, the next step is to construct the core }

con-tribution *c*_{共r兲. Recently, we proposed that}_{}*c*_{共r兲 can be}

estimated using first-order perturbation theory.14 _{In this }
con-text, the spin-polarized core orbitals_{}are constructed from
a linear combination of non-spin-polarized eigenstates_{}as

兩_{}_{典 = 兩}_{}

典 +

### 兺

⫽*C*兩典, 共14兲
where we have used the fact that the perturbing potential
*⌬V* differs just in sign for the spin-up and spin-down
chan-nels关see Eq. 共18兲兴 so that

*C*_{}=具*兩⌬V*兩典

− 共15兲

does not depend on. Equations.共14兲 and 共15兲 suggest that
only levelsneed to be taken into account, which are
ener-getically close to the given level and for which
具*兩⌬V*兩典 is considerable. It is to be noted that both
and_{}were calculated already for the reconstruction of the
AE wave functions.

Using Eq.共14兲, the contribution of the core orbitals to the charge density is expressed as

*c*_{共r兲 =}

### 兺

*core*具

_{}

_{兩r典具r兩}_{}

_{典,}

_{共16兲}

where the summation is over the core levels only. Corre-spondingly, one can show that the contribution of the core levels to the spin density becomes

*s*
*c*_{共r兲 =}_{}*c _{↑}*

_{共r兲 −}_{}

*c*

_{↓}

_{共r兲 = 4}### 兺

*core*Re

## 冋

具_{}

**兩r典**

### 兺

⫽*C*

_{}

**具r兩**典

## 册

, 共17兲 where “Re” indicates the real part. If the summation is only*over the core s levels*共e.g., in the case of isotropic HFPs兲, the imaginary part vanishes.

Returning to Eq. 共15兲, using the local spin-density
ap-proximation共LSDA兲,26_{it is possible to define a density }
func-tional form for*⌬V*_{}关see Eq. 共105兲 in Ref.14兴,

*⌬V*_{}⯝ −

*SOMO*_{共r兲}

*ˆ*2/3** _{共r兲}** , 共18兲
where

*SOMO*=

*v↑*

**共r兲−**

*v↓*

**共r兲. The total charge density**

*ˆ*

**共r兲,**

within the frozen-core approximation, is expressed as

*ˆ***共r兲 =***v↑***共r兲 +***v↓***共r兲 + 2**

### 兺

*core*

具**兩r典具r兩**典. 共19兲
In defining a DFT-based perturbing potential in Eq.共18兲,
one may consider the various generalized gradient
approxi-mation共GGA兲 formulations instead of the LSDA. However,
as Asada and Terakura27 _{have shown, GGAs produce a }
sin-gularity in the exchange potential *⌬V*_{} at the nucleus
posi-tion. A possible way to eliminate this singularity is to
repre-sent the nucleus not as a point but as a sphere with a finite
radius within which the nuclear charge is homogeneously
*distributed. However, the work by Battocletti et al.*28 _{has}
shown that the “finite nucleus model+ GGA” does not
pro-duce significant improvement over LSDA within the context
of hyperfine calculations.

It is to be noted that the perturbing potential, defined in
Eq.共18兲, differs for the two spin channels in sign only. Thus,
it is more convenient to replace*⌬V*_{}in Eq.共15兲 by an
effec-tive perturbing potential,

*⌬Vef f*=*⌬V↑*−*⌬V↓*= −

2

*SOMO*_{共r兲}

*ˆ*2/3**共r兲** , 共20兲

so that it is necessary to estimate the core spin polarization in one spin channel only 共here, we arbitrarily selected the spin-up channel兲.

Equation共20*兲 reveals that ⌬Vef f***共r兲 depends sensitively on**

*the type of SOMO. That is, if s-like orbitals contribute *
sub-stantially to the SOMO,*⌬Vef f***共r兲 expands significantly **

out-side the core region共see Figs. 3 and 12 in Ref. 14兲. In that
case, the effect of core spin polarization is expected to be
*small. As a result, s-like SOMOs lead to a dominant isotropic*
*HFP Aiso, whereas Aaniso* is likely to be negligible. On the

*other hand, if the SOMOs are dominated by p- or d-like*
orbitals,*⌬Vef f***共r兲 becomes more localized inside the core **

re-gion resulting in a large spin polarization among core
orbit-als 共see Figs. 3 and 12 in Ref. 14兲. Then, both core and
*valence spin-polarized levels contribute significantly to Aiso*.

*However, as the SOMO is non-s-like, Aaniso*is substantially

dominated by SOMO with a relatively small contribution from spin-polarized core levels. Hence, in both situations,

*s-like and non-s-like SOMOs, the neglect of core spin *

*polar-ization may not affect seriously the Aaniso*values, while the

*neglect of core spin polarization for non-s-like SOMOs *
*re-sults in erroneous Aiso*values.

Figure 1 illustrates the radial part of *⌬Vef f***共r兲 obtained**

directly from an all-electron Kohn-Sham DFT calculation
and as obtained from Eq.共18*兲 for C with a 2p-type SOMO.*
It clearly shows that the perturbing potential is strongly
lo-calized inside the core region*共Rc*for carbon is about 0.6 Å兲.

The figure also shows that*⌬Vef f***共r兲 from Eq. 共**18兲 coincides

ad-equate for estimating core spin polarization. At this stage,

*Aiso* *and Aaniso* are calculated through Eqs. 共1兲 and 共2兲 by

substituting the spin density*s***共r兲 obtained from the charge**
densities,*c***共r兲 and***v***共r兲. In the following, the accuracy of**
our method is examined by calculating the HFPs of various
systems and comparing the results with other computational
methods and with experiment.

**III. COMPUTATIONAL DETAILS**

In this work, the DFT calculations have been carried out
using the gradient-corrected PW9129,30_{exchange-correlation}

functional. For sake of comparison, calculations for the free
atoms have been performed also using the LSDA31 _{}
func-tional共see TableI兲. For the reconstruction of AE wave
func-tions 兩⌿*n*典, the PS-orbitals 兩*˜*典 and the PS-wave-functions

兩⌿*˜*

*n*

_{典 were generated using ultrasoft pseudopotentials}
共USPPs兲,32_{as implemented in the}_{VASP}_{code.}33_{Additionally,}
the AE orbitals兩*˜*_{}*典 and corresponding eigenvalues 兩˜*_{}典 were
computed using the ATOMDEF package34 _{with }
gradient-corrected exchange35 _{and Vosko-Wilk-Nusair correlation}36
functionals. For the atomic and molecular systems, we
con-sidered a cubic supercell with a side of 16 Å and integrations
in reciprocal space used the⌫ point only. For the crystalline
point defect calculation, a 6⫻6⫻6 Monkhorst-Pack mesh
has been used to sample the irreducible Brillouin zone.
Stop-ping criterion for structural relaxations was a magnitude of
the force less than 0.001 eV/ Å on each atom.

**IV. RESULTS AND DISCUSSIONS**

As a first step to evaluate the accuracy of our method, the
feasibility of reconstructing AE wave functions via Eq.共5兲 is
examined. Figure2*displays the radial part of the carbon 2s*
orbital as obtained directly from an AEMB calculation and as
reconstructed from the corresponding PS-orbital by means of
Eq. 共5兲. The figure clearly shows the excellent agreement
*between the reconstructed 2s orbital and the corresponding*
AE orbital. In fact, the lines are essentially indistinguishable.
To further test the accuracy of our method, the
contribu-tion of core levels to the spin density at the nucleus,*s*

*c*_{共0兲, as}

well as the contribution of both core and valence levels to the
*isotropic HFP, Aiso*

*c*

*and Aisov* , respectively, for various atoms,

*including a series of first-row elements and 3d transition*
metals, in their ground state electronic configurations have
been calculated using both LSDA31 _{and GGA}29 _{functionals}
共see TableI兲.
-5
-4
-3
-2
-1
0
0 0.3 0.6 0.9 1.2
∆
Veff
(r
)(
V
)
r (A)°
AEMB
analytical

FIG. 1.*共Color online兲 Comparison of perturbing potential ⌬V _{ef f}*
obtained directly from an all-electron mixed-basis calculation
共marked with “AEMB”兲 and as computed from Eq. 共20兲 共marked

with “analytical”兲 for a single isolated carbon atom in its ground
state configuration. The vertical dotted line indicates the cutoff
ra-dius of carbon*共R _{c}*= 0.6 Å兲.

TABLE I. Contribution of core共and valence兲 levels to the spin density at the nucleus, * _{s}c共0兲 共in e/Å*3

_{兲, and, correspondingly, to the}

*isotropic HFP’s, A*

_{iso}c*共and A*兲 共in MHz兲, as computed using the perturbative core-level polarization 共PCLP兲 method and as directly obtained from an all-electron mixed-basis

_{iso}v*共AEMB兲 calculation, for a series of first-row elements and 3d transition metals. The spin configurations all*follow Hund’s rule.

Atom Configuration

*s*
*c*_{共0兲}

*A _{iso}c*

*A*

_{iso}vPCLPa _{PCLP}b _{AEMB} _{PCLP}a _{PCLP}b _{AEMB} _{PCLP}a _{PCLP}b _{AEMB}

C *2s*2* _{2p}*2

_{兵2↑ ,0↓其}_{−1.344}

_{−1.350}

_{−1.351}

_{−111.9}

_{−112.4}

_{−112.4}

_{130.2}

_{132.4}

_{133.8}N

*2s*2

*3*

_{2p}

_{兵3↑ ,0↓其}_{−2.971}

_{−2.985}

_{−3.001}

_{−47.4}

_{−47.6}

_{−47.8}

_{52.7}

_{56.6}

_{57.3}O

*2s*2

*4*

_{2p}

_{兵3↑ ,1↓其}_{−2.830}

_{−2.849}

_{−2.853}

_{127.0}

_{127.8}

_{128.0}

_{−151.3}

_{−158.0}

_{−159.8}F

*2s*2

*2p*5

*兵3↑ ,2↓其*−1.533 −1.546 −1.550 −954.8 −962.9 −965.4 1237.4 1257.7 1261.5 Sc

*3d*1

*2*

_{兵1↑ ,0↓其4s}_{−0.581}

_{−0.588}

_{−0.610}

_{−93.5}

_{−94.7}

_{−98.2}

_{78.7}

_{72.5}

_{75.5}Ti

*3d*2

*兵2↑ ,0↓其4s*2 −1.131 −1.145 −1.168 21.1 21.4 21.8 -12.8 −11.7 −10.7 V

*3d*3

*2*

_{兵3↑ ,0↓其4s}_{−2.505}

_{−2.563}

_{−2.637}

_{−145.6}

_{−149.0}

_{−153.3}

_{108.4}

_{98.5}

_{100.6}Mn

*3d*5

*兵5↑ ,0↓其4s*2 −3.817 −3.906 −4.090 −124.9 −127.8 −133.8 81.3 69.4 70.8 Cu

*3d*10

*1*

_{4s}

_{兵1↑ ,0↓其}_{−0.241}

_{−0.246}

_{−0.250}

_{−42.3}

_{−43.2}

_{−43.9}

_{6044.7}

_{5984.3}

_{5979.0}Zn

*3d*10

*2*

_{4s}_{0.0}

_{0.0}

_{0.0}

_{0.0}

_{0.0}

_{0.0}

_{0.0}

_{0.0}

_{0.0}

For both *s*

*c_{共0兲 and A}*

*iso*

*c*

, the table reveals that both
LSDA31_{and GGA}29_{functionals give comparable results: for}
most atoms 1% or less, and about 2% for V and Mn. It is also
apparent that there is excellent agreement between our PCLP
and the AEMB methods: for first-row elements, the
differ-ence is about 1%*共both LSDA兲; for 3d transition metals also,*
in spite of the different character of the SOMO, good
agree-ment between PLCP and AEMB exists with the largest
*dif-ferences occurring when the number of unpaired 3d electrons*
is largest. Even in the extreme case of Mn, with the
maxi-mum number of unpaired electrons, the difference between
PLCP and AEMB is still only 7%共LSDA兲.

On the other hand, the contribution of valence levels to
the isotropic HFP turns out to be more sensitive to the choice
of exchange-correlation functional as also recently reported
*by Declerck et al.*21 _{According to Table}_{I}_{, for first-row }
ele-ments, the difference between the LSDA and GGA results of

*A _{iso}v* varies from 2%共for C兲 to 7% 共for N兲. Such a difference

*becomes even more pronounced for 3d transition metal*

*at-oms. In the extreme case of Mn, the LSDA value of Aisov*

differs about 17% from the corresponding GGA value. The
comparison between PCLP and AEMB results共both LSDA兲
reveals a comparable difference for first-row elements, while
a better agreement between the two sets of results can be
*achieved for 3d transition metal elements. An interesting *
ob-servation is that the agreement between GGA PCLP and
LSDA AEMB results is systematically better than that
be-tween LSDA PCLP and the LSDA AEMB results.

Both PLCP and AEMB show that the value of*s*
*c*_{共0兲 }

in-creases with the number of unpaired electrons. For Zn, with
fully occupied valence orbitals *共3d*10*4s*2兲, *s*

*c*_{共0兲 becomes}

zero, whereas for Mn*共3d*5* _{4s}*2

_{兲, the absolute value of}

_{}

*s*
*c*_{共0兲 is}

the largest among the elements considered here.
Further-more, TableIshows that *s*

*c _{共0兲 for elements with non-s-like}*

SOMO 共e.g., F or Sc兲 is considerably larger than that for

*s-like SOMO*共e.g., Cu兲. These are all in accordance with our

earlier assumptions that the strength of induced core spin

polarization depends sensitively on the number and type of
unpaired electrons in SOMOs关see Eq. 共17兲兴. It is to be noted
that*s*

*c*** _{共0兲 is always negative.}**8–10

_{For the first-row elements,}it can be explained as follows: according to the Pauli exclu-sion principle, the exchange interaction induced by unpaired electrons is attractive but applicable to electrons in the same spin channel only.8–10,37

*As a result, the core 1s electrons of*spin majority type are pulled a little outward, thereby leaving behind a slight depletion of their corresponding charge den-sity at the vicinity of nucleus. Thus, the spin denden-sity

*associ-ated with core electrons becomes negative. For 3d transition*metals, a more complicated mechanism is needed for a proper interpretation of the negative sign of the core spin densities共see our discussion in Sec. IV B of Ref.14兲.

The hyperfine calculations of molecular radicals and
tran-sition metal complexes provide a more practically relevant
test than atomic calculations. Table II shows the isotropic
HFPs as computed with our method, the recently proposed
CSPC method,20 _{and the AEMB method, as well as the }
ex-perimental data, for a series of molecular radicals containing
first-row elements. At all the nuclei in each of the molecules,
our results are similar to the CSPC results and are also in
good agreement with the AEMB and the experimental
re-sults. Obviously, there is no contribution from core levels to

*Aiso*of 1H. For the other atoms, the comparison between the

contributions of core and valence levels to the isotropic HFPs reveals that the former are always lower in magnitude and with an opposing sign as compared to the corresponding valence contributions. A more detailed comparison of the two contributions allows us to divide such molecules, on the basis of their SOMO type, into two main groups.

The first group, including CH3, C3H5, and H2CN, is char-acterized by an unpaired electron occupying a pure

*p*共兲-type orbital. The characteristic feature of this group of

molecules is that the contributions of both core and valence
*levels to Aiso*are relatively moderate and comparable to each

other. In the second group, including C2H3, HCO, FCO, and
NO2*, the unpaired electron occupies an sp-hybridized*
*or-bital. This implies that the s character of SOMO contributes*
directly to the spin density at the nuclei of the radical. As a
*result, Aiso* of these species is substantially dominated by

valence orbitals with large values, as shown in TableII.
To further clarify the above classification, a sample from
each group, C_{3}H_{5} and HCO, is analyzed in more detail. In
C3H5, all the molecular orbitals 共MOs兲 are doubly
occu-pied. Thus, the unpaired electron can only be distributed
among threeMOs where eachMO is localized on one of
the carbon atoms. According to TableII, the distribution of
*SOMO*_{is such that the sign of the spin density at the position}

of the central carbon atom, 13C␣, becomes negative while it
is positive for the two outer carbon atoms, 13C. The
sim-plest configuration that can be used to describe such a
distri-bution for*SOMO*is47

*SOMO*
=2
3兩
*L* _{兩}2_{−}1
3兩␣兩
2_{+}2
3兩
*R* _{兩}2_{,} _{共21兲}
*where the superscripts R and L indicate the left and the right*
13_{C} _{atoms, respectively. This implies that the induced spin}
-0.6
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
1.2
0 0.5 1 1.5 2
R2s
(r)r
r (A)°
AEMB
USPP
reconstructed

*polarization of both core 1s and valence 2s orbitals of*13C␣
must be lower in magnitude with opposite sign compared to
that of 13C atoms. Interestingly, the calculated isotropic
HFPs confirm the accuracy of this simple model.

For HCO, the unpaired electron, which occupies aMO,
is substantially localized on the HuC bond. Our
*calcula-tions indicate that the s character of**SOMO*_{on}1_{H and}13_{C is}
*26% and 6.9%, respectively. As a result, the Aiso* values of

*the corresponding atoms are significantly dominated by the s*
character of the SOMO共see TableII兲. On the other hand, the
*p character of the SOMO leads to a slight spin polarization*

*of core 1s orbitals, thereby leaving behind a small negative*
contribution to the isotropic HFP of13C.

Having confirmed the accuracy of our method for
first-row elements and the related molecular radicals, we apply it
to the calculation of both isotropic and anisotropic HFPs of a
*series of 3d transition metal oxides. The results along with*
those obtained from AEMB calculations and all-electron

Gaussian-type orbital共AE-GTO兲 calculations by Munzarova
and Kaupp48 _{as well as the experimental data are shown in}
Table III. The overall agreement between our method, the
AEMB and AE-GTO results, and experiment appears
reason-able. The only exception is the ligand component of CuO for
*which both Aiso* *and Aaniso* obtained from our PCLP and

AEMB calculations are in better agreement with experiment than the AE-GTO results. A possible explanation for the anomalous AE-GTO results might be a poor basis set selec-tion.

As for the molecular radicals, TableIIIshows that HFPs
depend sensitively on the type of SOMO. In ScO, TiO, and
VO, the*-type SOMOs have predominantly metal 4s *
*char-acter with some mixing from 3dz*2 *and 4p. This explains the*

large isotropic HFPs on the metal atoms, while the
*aniso-tropy is relatively less pronounced. It is to be noted that gI*of

47_{Ti is negative so that the related HFPs become negative. In}
going to higher spin multiplicities, MnO exhibits another

TABLE II. Comparison of isotropic HFP for a number of molecular radicals共in MHz兲. The values labeled “PCLP” were computed with Eq. 共1兲 using spin densities obtained from Eqs. 共17兲 and 共11兲 for core and valence contributions, respectively, “CSPC” refers to core

spin-polarization correction method results from Ref.20, “AEMB” refers to all-electron mixed-basis calculations, and “Expt.” refers to experimental results reported in Refs.37–46.

Molecule Lewis structure Nucleus

PCLP

CSPC AEMB Expt.

Core Valence Total

CH_{3} _{C}* _{˙ H}*
3
1

_{H}

_{−65.6}

_{−65.6}

_{−65.3}

_{−69.8}

_{−70.3}a

_{;}

_{−64.5}b 13

_{C}

_{−110.1}

_{186.4}

_{76.3}

_{74.3}

_{82.7}

_{79.6}a

_{;}

_{107.3}b C

_{2}H

_{3}

*␣*

_{共HC˙兲}_{v共CH}2兲 1

_{H}␣

_{40.4}

_{40.4}

_{43.7}

_{39.8}

_{35.9}1

_{H}␣ 

_{192.5}

_{192.5}

_{181.9}

_{193.6}

_{184.7}1

_{H}

*s*

_{112.1}

_{112.1}

_{113.8}

_{115.2}

_{111.0}13

_{C}␣

_{−104.8}

_{404.4}

_{299.6}

_{297.1}

_{303.2}

_{301.5}13

_{C}

_{10.1}

_{−26.3}

_{−16.2}

_{−11.5}

_{−22.4}

_{−24.1}C3H5 共H2C兲 .... u共CH兲␣ .... u共CH2兲 1 H␣ 10.1 10.1 10.6 12.0 11.5 1

_{H}1

_{−41.5}

_{−41.5}

_{−42.0}

_{−41.8}

_{−41.5}1

_{H}2

_{−38.4}

_{−38.4}

_{−38.9}

_{−37.8}

_{−38.9}13

_{C}␣

_{18.2}

_{−63.1}

_{−44.9}

_{−40.9}

_{−44.8}

_{−48.2}13

_{C}

_{−68.1}

_{113.5}

_{45.4}

_{47.6}

_{50.2}

_{61.4}H2CN H2C

*vN˙*1 H 236.8 236.8 238.8 237.1 233.2 13

_{C}

_{13.7}

_{−79.9}

_{−66.2}

_{−61.9}

_{−70.6}

_{−81.0}14 N −52.7 67.0 14.3 11.8 21.3 26.1 HCO

_{H}

*1H 375.8 375.8 375.0 378.3 379.5 13 C −63.0 453.7 390.7 411.1 379.7 375.2 17*

_{uC˙vO}_{O}

_{35.0}

_{−66.4}

_{−31.4}

_{−26.6}

_{−32.8}

_{−42.3}FCO

_{F}

*19F −184.4 1109.5 925.1 972.2 882.5 906.0 13*

_{uC˙vO}_{C}

_{−42.9}

_{852.5}

_{809.6}

_{822.8}

_{793.4}

_{803.2}17 O 33.9 −76.3 −42.4 −46.5 −44.8 NO

_{2}

_{N}

*2 14*

_{˙ O}_{N}

_{−8.1}

_{170.1}

_{162.0}

_{160.6}

_{146.0}

_{153.6}17 O 37.5 −85.2 −47.7 −51.8 −54.4 −45.7⬃ −56.9

type of SOMOs. Comparing with the previous cases, MnO
has two additional SOMOs, antibonding orbitals with metal
*3d*_{} *and 4p*_{} *as well as ligands with 2p*_{} character. Due to
*the large number of d-type SOMOs, spin-polarization effects*
via the core shells are more pronounced. Suffice it to say that
the contribution of the core levels共−230 MHz兲 to the total
value of the isotropic HFP of 55Mn 共+452.4 MHz兲 is very
significant. Thus, the neglect of core spin polarization
defi-nitely leads to a large error in the isotropic HFP. Finally, in
CuO, the situation becomes even more critical. Here, there is
*no contribution from the metal 4s states to the SOMO. *
*In-stead, the SOMO is due to 3d-type orbitals and the negative*
core contribution is actually larger than that of the valence
states.

Since for most of the present transition metal oxides,
*SOMO* _{is significantly localized on the metal atoms,}49 _{both}

*Aiso* *and Aaniso*of the ligand part turn out to be small. The

only exception is CuO, whose spin density is slightly delo-calized on the17O atom.48

As a final verification of the PLCP method, we have
com-puted the HFPs pertaining to a point defect in a II-VI crystal
i.e., a positively charged Zn interstitial共Zn*i*+兲 in ZnSe. A 32

+ 1 atom cell was constructed with 16 Zn and 16 Se on the
cubic zinc-blende ZnSe sites where the cation, Zn+_{, occupies}

a tetrahedral interstitial site surrounded by Se atoms*共T _{d}*Se兲.51
The Zn interstitial caused relaxations: in our calculations the
first 共second兲 neighbor shells Se1NN共Zn2NN兲 move outward
by 0.12 Å共0.07 Å兲 along 具111典 共具100典兲 directions, similar to
relaxations found by van de Walle and Blöchl.15

_{The }calcu-lated HFPs for the corresponding Zn

*i*+ and its first 共Se1NN兲, second共Zn

_{2NN}兲, and third 共Se

_{3NN}兲 are listed in TableIV, and compared with theoretical15

_{and experimental}50

_{data. The}table indicates that the HFPs obtained by PCLP are in good agreement with the experimental data; in fact, the agreement is a little better than the previous calculation.15

_{Our }calcula-tions show that the spin density is substantially localized on the Zn

*i*

+

and its first neighbors, Se1NN, and is significantly
*4s-like. In the case of Se*1NN, it has an additional minor
*con-tribution from its valence 4p orbitals. Therefore, the core*
contribution is rather significant on Se_{1NN}, and it explains
why the earlier calculation without core contribution15
over-estimated the HFP on this site.

**V. CONCLUSIONS**

A computationally expedient, density-functional-theory-based method for calculating hyperfine parameters within pseudopotential formalism has been presented. The accurate

TABLE III. Hyperfine parameters 共in MHz兲 for transition metal complexes as computed using PCLP, AEMB, and AE-GTO 共from Ref. 48兲 methods, and as experimentally measured 共Expt.兲 共as reported

in Ref.48兲.

Complex Nucleus HFP

PCLP

AE-GTOa _{AEMB} _{Expt.}

Core Valence Total

representation of the spin density at and in the vicinity of the
nucleus was achieved by reconstructing the full-nodal
behav-ior of the wave functions with inclusion of core spin
polar-ization via first-order perturbation theory. The method was
successfully applied to calculate the hyperfine parameters of
various molecular radicals, complexes, and a crystalline
*point defect containing first-row of elements and 3d *
transi-tion metals. The current method is highly versatile because it
does not impose restrictions on the exchange-correlation
functional in the pseudopotential calculations. The additional

information regarding the code with which these calculations were carried out can be obtained upon request from Bahr-amy.

**ACKNOWLEDGMENT**

The authors gratefully acknowledge the Center for Com-putational Materials Science at the Institute for Materials Re-search for allocations on the Hitachi SR8000 supercomputer system.

*saeed@imr.edu

1_{A. Abragam and B. Bleaney, Electronic Paramagnetic Resonance}

*of Transition Ions*共Clarendon, Oxford, 1970兲.

2_{J. A. Weil, J. R. Bolton, and J. E. Wertz, Electron Paramagnetic}

*Resonance: Elementary Theory and Practical Applications*

共Wiley, New York, 1994兲.

3_{I. Bertini, C. Luchinat, and G. Parigi, Solution NMR of }

*Paramag-netic Molecules*共Elsevier, Amsterdam, 2001兲.

4_{D. Feller and E. R. Davidson, Theor. Chim. Acta 68, 57}_{共1985兲.}
5_{S. A. Perera, J. D. Watts, and R. J. Bartlett, J. Chem. Phys. 100,}

1425共1994兲.

6_{R. K. Nesbet, Proc. R. Soc. London, Ser. A 230, 312}_{共1955兲.}
7_{G. W. Pratt, Jr., Phys. Rev. 102, 1303}_{共1956兲.}

8_{A. J. Freeman and R. E. Watson, Phys. Rev. Lett. 5, 498}_{共1960兲.}
9_{R. E. Watson and A. J. Freeman, Phys. Rev. 123, 2027}_{共1961兲.}
10_{R. E. Watson and A. J. Freeman, Phys. Rev. Lett. 6, 277}_{共1961兲.}
11_{P. Hohenberg and W. Kohn, Phys. Rev. 136, B864}_{共1964兲.}
12_{H. Overhof, M. Scheffler, and C. M. Weinert, Phys. Rev. B 43,}

12494共1991兲.

13** _{H. Katayama-Yoshida and N. Hamada, Phys. Rev. B 35, 407}**
共1987兲.

14_{M. S. Bahramy, M. H. F. Sluiter, and Y. Kawazoe, Phys. Rev. B}
**73, 045111**共2006兲.

15** _{C. G. Van de Walle and P. E. Blöchl, Phys. Rev. B 47, 4244}**
共1993兲.

16** _{J. Vackář, M. Hyt’ha, and A. Šimůnek, Phys. Rev. B 58, 12712}**
共1998兲.

17_{B. Meyer, K. Hummler, C. Elsässer, and M. Fähnle, J. Phys.:}
**Condens. Matter 7, 9201**共1995兲.

18_{B. Hetényi, F. D. Angelis, P. Giannozzi, and R. Car, J. Chem.}
**Phys. 115, 5791**共2001兲.

19_{P. E. Blöchl, Phys. Rev. B 50, 17953}_{共1994兲.}

20_{O. V. Yazyev, I. Tavernelli, L. Helm, and U. Röthlisberger, Phys.}
**Rev. B 71, 115110**共2005兲.

21_{R. Declerck, E. Pauwels, V. Van Speybroeck, and M. Waroquier,}
**Phys. Rev. B 74, 245103**共2006兲.

22_{D. A. Goodings, Phys. Rev. 123, 1706}_{共1961兲.}

23** _{G. D. Gaspari, W.-M. Shyu, and T. P. Das, Phys. Rev. 134, A852}**
共1964兲.

24_{D. R. Lide, Handbook of Chemistry and Physics, 85th ed.}_{共CRC,}
Boca Raton, 2004兲, pp. 9–93.

25* _{F. Herman and S. Skillman, Atomic Structure Calculations}*
共Prentice-Hall, Englewood Cliffs, NJ, 1963兲.

26_{G. L. Oliver and J. P. Perdew, Phys. Rev. A 20, 397}_{共1979兲.}
27_{T. Asada and K. Terakura, Phys. Rev. B 46, 13599}_{共1992兲.}
28** _{M. Battocletti, H. Ebert, and H. Akai, Phys. Rev. B 53, 9776}**
TABLE IV. Hyperfine parameters共in MHz兲 for Zn

*+*

_{i}*at the T*Sesite in ZnSe. Results of this work关previous

_{d}pseudopotential calculation 共Ref. 15兲兴 are listed in columns under PCLP 共PP兲. Values are given for the

interstitial at the center of the defect as well as for the first共77Se兲, second 共67Zn兲, and third 共77Se兲 neighbors.

Nucleus HFP

PCLP

PPa _{Expt.}b

Core Valence Total

67
Zn* _{i}*+

*Aiso*−11.1 1093.1 1082.0 1078 1089

*A*−0.1 5.2 5.1 77 Se

_{aniso}_{1NN}

*Aiso*−84.3 693.5 609.2 736 481

*Aaniso*−0.9 15.2 14.3 11 16.8 67

_{Zn}2NN

*Aiso*0.0 10.7 10.7 10

*Aaniso*0.0 0.0 0.0 77

_{Se}3NN

*Aiso*−0.1 41.9 41.8 55 37.5

*Aaniso*0.0 0.6 0.6 a

_{Reference}

_{15}

_{.}

共1996兲.

29* _{J. P. Perdew, Electronic Structure of Solids ‘91, edited by P. }*
Zi-esche and H. Eschrig共Akademie Verlag, Berlin, 1991兲, p. 11; J.

**P. Perdew and Y. Wang, Phys. Rev. B 45, 13244**共1992兲. 30

_{J. P. Perdew, J. A. Chevary, S. H. Vosko, K. A. Jackson, M. R.}

**Pederson, D. J. Singh, and C. Fiolhais, Phys. Rev. B 46, 6671**
**共1992兲; 48, 4978共E兲 共1993兲.**

31_{J. P. Perdew and A. Zunger, Phys. Rev. B 23, 5048}_{共1981兲.}
32_{D. Vanderbilt, Phys. Rev. B 41, 7892}_{共1990兲.}

33_{Vienna ab initio software package}_{共}_{VASP}_{兲, Version 4.6.12, http://}
cms.mpi.univie.ac.at/vasp/

34_{R. Pis Diez, computer code}_{ATOMDEF}_{, Version 1.4, 2003, can be}
downloaded from http://www.quimica.unlp.edu.ar/cequinor/
rpd_en.htm

35** _{J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77,}**
3865共1996兲.

36** _{S. H. Vosko, L. Wilk, and M. Nusair, Can. J. Phys. 58, 1200}**
共1980兲.

37_{D. M. Chipman, Theor. Chim. Acta 82, 93}_{共1992兲.}

38_{E. Hirota and C. Yamada, J. Mol. Spectrosc. 96, 175}_{共1985兲.}
39_{R. W. Fessenden and R. H. Schuler, J. Chem. Phys. 39, 2147}

共1963兲.

40_{R. W. Fessenden, J. Phys. Chem. 71, 74}_{共1967兲.}

41_{H. J. McManus, R. W. Fessenden, and D. M. Chipman, J. Phys.}
**Chem. 92, 3778**共1988兲.

42_{H. J. McManus, R. W. Fessenden, and D. M. Chipman, J. Phys.}
**Chem. 92, 3781**共1988兲.

43** _{F. J. Adrian, B. F. Kim, and J. Bohandy, J. Chem. Phys. 82, 1804}**
共1985兲.

44** _{F. J. Adrian, E. L. Cochran, and V. A. Bowers, J. Chem. Phys. 43,}**
462共1965兲.

45** _{E. L. Cochran, F. J. Adrian, and V. A. Bowers, J. Chem. Phys. 44,}**
4626共1966兲.

46_{Z. Luz, A. Reuveni, R. W. Holmberg, and B. L. Silver, J. Chem.}
**Phys. 51, 4017**共1969兲.

47_{H. M. McConnell, J. Chem. Phys. 28, 1188}_{共1958兲.}

48** _{M. L. Munzarová and M. Kaupp, J. Phys. Chem. A 103, 9966}**
共1999兲.

49_{M. L. Munzarová, P. Kubáček, and M. Kaupp, J. Am. Chem. Soc.}
**122, 11900**共2000兲.

50_{F. Rong and G. D. Watkins, Phys. Rev. Lett. 58, 1486}_{共1987兲.}
51_{The supercell without interstitial has as translation vectors}

*a*_{ZnSe}*共−1,1,1兲, a*_{ZnSe}*共1,−1,1兲, and a*_{ZnSe}*共1,1,−1兲, where a*_{ZnSe}