Table Of ContentSELF-CONSISTENT METHODS IN NUCLEAR STRUCTURE
PHYSICS
8 JACEK DOBACZEWSKI
9 Institute of Theoretical Physics, Warsaw University
9 ul. Hoz˙a 69, PL-00681 Warsaw, Poland
1 E-mail: [email protected]
n
a Joint Institute for Heavy Ion Research, Oak Ridge National Laboratory,
J P.O. Box 2008, Oak Ridge, TN 37831, U.S.A.
7
2 Department of Physics, University of Tennessee,
Knoxville, TN 37996, U.S.A.
1
v WepresentaverybriefdescriptionoftheHartree-Fockmethodinnuclearstructure
6 physics,discussthenumericalmethodsusedtosolvetheself-consistentequations,
5 and analyze the precision and convergence properties of solutions. As an appli-
0 cation we present results pertaining to quadrupole moments and single-particle
quadrupolepolarizationsinsuperdeformednucleiwithA∼60.
1
0
8
1 Introduction
9
/
h Self-consistent methods have been used in the low-energy nuclear structure
t studies over many years, and represent a mature field with numerous suc-
-
l cessful applications.1,2,3,4,5 A number of computer codes solving the nuclear
c
u Hartree-Fock (HF) problem have already been developed. Two types of effec-
n tive nucleon-nucleon interactions have been mainly employed. Starting with
: the work of Vautherin and Brink6 many authors have applied the nuclear
v
i HF theory with the Skyrme effective interaction, while the work of Gogny7,8
X
initiated numerous studies with the force which carries his name.
r The methods employed to solve the HF equations depend mainly on the
a
effective force used and on the assumed symmetries of the many-body wave
functions. For the solutions which allow at least triaxial deformations, two
different methods have been applied for the two above mentioned effective in-
teractions. The first one, used in conjuncture with the Skyrme interaction, is
formulated in the spatial coordinates and makes use of the finite-difference,9
orFourier,10 or spline-collocation11,12 methods to approximatedifferentialop-
erators. The solution is then obtained by using the imaginary time evolution
operator.13
The second one, used for the finite-range Gogny interaction, employs a
truncatedharmonicoscillator(HO)basis14,15 andsolvestheproblemeitherby
aniterativediagonalizationofthemean-fieldHamiltonian,orbythegradient16
1
ortheconjugategradient17methods. Recently,themethodwhichincorporates
the advantages of both existing approaches, and combines the robustness of
theCartesianHObasiswiththesimplicityoftheSkyrmeinteraction,hasbeen
implemented.18
The methods using spatial coordinates have several advantages. First of
all,variousnuclearshapes canbe easilytreatedonthe samefooting; the same
cubic lattice of points in three spatial dimensions is suitable to accommodate
wave functions with, in principle, arbitrary deformations restricted only by
a specific symmetrization of the lattice. This allows easy studies of systems
for which the deformation is not a priori known, or is ill defined because
of deformation instabilities or a shape coexistence. Secondly, using spatial
coordinates allows studies which address the asymptotic form of nucleonic
wave functions at large distances. This is particularly important for a precise
description of weakly bound nuclei, where the use of spatial coordinates is a
necessity.19 Third, for the Skyrme zero-range interaction, the mean fields are
local(apartfroma velocitydependence) andcanbe easilyprogrammedin the
spatialcoordinates. Lastbutnotleast,thetreatmentofwavefunctionsonlarge
lattices(12 12 12isatypicalexample)is easilyamenabletovectorizationor
× ×
parallelization of the algorithm.
Methods using the HO basis have other advantages. Firstly, the basis
provides a natural cut-off for many operators which otherwise are unbound
and require particularly delicate treatment in the spatial coordinates. This
concerns in particular the multipole moment and the angular momentum op-
erators which are often used as constraining operators. For the corresponding
constraints the solutions can become unstable when the non-zero probability
amplitudes (wave functions) move towards large distances as it is the case
for e.g. weakly bound nucleons. Secondly, much smaller spaces are usually
requiredtodescribethenuclearwavefunctionswithinagivenprecision. Typi-
callyabasisofabout300HOwavefunctionsissufficientformostapplications.
Third,the iterativediagonalizationofthe meanfield Hamiltoniancanbe used
to find the self-consistent solutions, which provides a rapidly converging al-
gorithm, and, last but not least, scalar or superscalar computers can also be
used, because the typical sizes of the information handled is smaller and the
performance is less dependent on the use of a vector processor.
Detaileddiscussionofthe stability,convergence,andefficiencyofmethods
using the expansion on the Cartesian HO basis has recently been presented
together with the description of the hfodd code.18 In the present communi-
cation we give a brief r´esum´e of the methods employed (Sec. 2), present tests
of the precision and discuss the CPU times required (Sec. 3), and illustrate
possible applications by presenting a few recentresults obtained for the A 60
∼
2
nuclei (Sec. 4).
2 Hartree-Fock method
′
The eigenequations for the HF single-particle Routhians h are called the HF
α
equations,
h′ ψ (rσ)=e′ ψ (rσ), (1)
α i,α i,α i,α
where the neutron (α=n) and proton (α=p) Routhians read
¯h2
h′ = ∆+Γ +δ UCoul+Umult ω Jˆ. (2)
α −2m α αp − y y
Detailed expressions for the nuclear mean field operators Γ , Coulomb poten-
α
tial UCoul, and multipole constraint potential Umult are given in Ref.18 After
solving Eq. (1) one calculates the total energy of the system as a sum of the
kinetic, Skyrme, Coulomb, and pairing terms
= kin+ Skyrme+ Coul+ pair. (3)
E E E E E
′
Thesameenergycanbe calculatedbyusingtheRouthianeigenvaluese and
i,α
occupation probabilities v2 as
i,α
˜= 1 v2 e′ + 1 kin+ pair rear+ 1 Coul mult cran. (4)
E 2X i,α i,α 2E E −E 3Eexch −Ecorr −Ecorr
i,α
Here, the rearrangement energy rear results from the density dependence of
E
the Skyrme interaction. Similarly, the Coulombexchange energy Coul can be
Eexch
considered as resulting from a zero-order interaction term depending on the
density as ρ−2/3. The remaining terms correct for the fact that Routhians
p
contain constraint terms which are not present in the total energy.
The difference between the energies ˜and is exactly equal to zero when
E E
thedensitiesandfieldsdonotchangefromoneiterationtothenextone. Hence
their difference
δ = ˜ (5)
E E −E
provides a useful measure of the quality of convergence, and is called the sta-
bility of the HF energy.
The HF equations(1)canbe solvedbyexpanding the single-particlewave
functions ψ (rσ) onto the deformed HO wave functions ψ (rσ) in the
i nxnynz,sz
Cartesian coordinates, i.e.,
Nx Ny Nz
ψ (rσ)= Anxnynz,szψ (rσ). (6)
i X X X X i nxnynz,sz
nx=0ny=0nz=0sz=−21,12
3
Here N , N , and N are the maximum numbers of the HO quanta corre-
x y z
sponding to the three Cartesian directions. However, the sums over n , n ,
x y
andn areperformedoverthegridofpointswhichformapyramidratherthan
z
a cube.
The HO wave functions have the standard form
ψ (rσ)=ψ (x)ψ (y)ψ (z)δ , (7)
nxnynz,sz nx ny nz szσ
where
ψnµ(xµ)=bµ12Hn(0µ)(ξµ)e−21ξµ2, (8)
and ξ =b x are dimensionless variables scaled by the oscillator constants
µ µ µ
b = mω /¯h. (9)
µ q µ
PolynomialsH(0)(ξ)areproportionaltothestandardHermiteorthogonalpoly-
n
nomials H (ξ),20
n
−1
H(0)(ξ)= √π2nn! 2 H (ξ). (10)
n (cid:0) (cid:1) n
Accuracy of the solution of the HF equations with the wave functions ex-
pandedontotheCartesianHObasis,Eq.(6),dependsonthethreeparameters
¯hω , h¯ω , h¯ω defining the HO frequencies in three Cartesian directions, and
x y z
on the number M of the HO states included in the basis. In the code hfodd
we use the standard prescription21,22 to chose the HO states included in the
basis, namely, the M states with the lowest HO single-particle energies,
ǫ =h¯ω (n + 1)+h¯ω (n + 1)+h¯ω (n + 1), (11)
nxnynz x x 2 y y 2 z z 2
areselectedamongthosewhichhaven N ,n N ,andn N ,whereN is
x 0 y 0 z 0 0
≤ ≤ ≤
the fixed maximum number of HO quanta. It should be noted that in general
both M and N have to be specified to define the basis. Only for large N ,
0 0
the basis is defined solely by M and does not depend on N . In this case,
0
the grid of points (n ,n ,n ) defining the states included in the basis forms a
x y z
pyramidinthreedimensions,withthe inclinedfacedelimitedbythe condition
ǫ const. On the other hand, only for small values of N the basis is
nxnynz≤ 0
defined solely by N and does not depend on the energy cut-off. In this case
0
the corresponding grid of points n n n forms a cube of the size N . In all
x y z 0
intermediate cases the shape of the basis corresponds to a pyramid with the
corners cut off, or to a cube with the corners cut off. Usually N is chosen
0
large enough so that all the states allowed by the energy cut-off are included
in the basis.
4
3 Precision and convergence
TheHFequations(1)canbesolvedbyusingstandarditerativemethods1which
arebasedontheself-consistencyprinciple. Namely,theself-consistentsolution
′
is attained when the mean fields Γ appearing in the mean-field Routhian h
α α
equal to the mean fields calculated from the occupied single-particle states
ψ (rσ). Therefore,by iteratively calculating the meanfields and using them
i,α
toobtainnewsingle-particlestatesonemayeventuallyreachtheself-consistent
solution. A convergence of such a procedure is not guaranteed and, in fact,
in nuclear applications the method such as formulated above almost always
diverges. Thereasonforthatisthefactthattheconsecutivecorrectionsofthe
mean fields turn out to be too violent and lead to too strong modifications of
the states. Since early years of practical implementations of the HF method
in nuclear structure physics the simple method to remedy this situation has
been devised. Namely, only the fraction f of the mean fields calculated from
the single-particlestatesisusedineveryiteration,anditiscombinedwiththe
fraction 1 f of the old mean fields.
−
In practicalcalculations,the optimal value off is determinedby checking
the convergence rate. An example of such an analysis is presented in Fig. 1,
where the stability energy (5) is shown for the values of f ranging from 30%
to 90%. It can be seen that the rate of convergence increases very fast when
f increases from 30% to about 80%, but then it decreases again dramatically.
Unfortunately, the optimal value of f depends critically on the shell structure
of the nucleus; it is rather high for magic nuclei (as superdeformed 152Dy )
66 86
which have large shell gaps, and has to be lowered when there are several
single-particlestatesneartheFermisurface. Inpracticalcalculationsitismore
efficienttouseratherlowvaluesoff,andensureconvergenceindependentlyof
the shell structure, than to risk a divergence. Usually, f=50% can be a safely
recommended value.
The total CPU times required to solve the HF equations depend not only
on the convergence rate governed by f, but also (of course) on the overall
speedofthe computer and(veryimportantly)onthe performancesofthe spe-
cific compiler optimization procedures. This is illustrated in Fig. 2, where the
relative times required to perform specific tasks are shown in percent of the
totaltime. ThetotalCPUtimesrequiredtocalculatethecompletesuperdefor-
medbandin152Dy rangefromone houronthe CrayC-90supercomputerto
66 86
several hours on fast workstations (SG Power Challenge L or IBM RS/6000),
and to about 36 hours ona typicalworkstation(SUN-Hyper). However,these
times result from averaging very different performances pertaining to specific
tasks. Typically, three tasks take most of the time, namely, the calculation of
5
102 30% ---> 25 iter/decade
40% 19
)
V 50% 15
Me 100 60% 12
y ( 70% 10
g 80% 9
er 10-2 90% 23
n
E
f
o
y 10-4
bilit
a
St 10-6
10-8
0 20 40 60 80 100 120 140 160
Number of iterations
Figure1: Stabilityenergyasfunctionofnumberofiterationsforthecalculationsperformed
for the superdeformed state of 16562Dy86. Different curves correspond to different values of
thefractionf (seetext), andleadtodifferentnumbers ofiterations whicharenecessaryto
improvethestabilitybyanorderofmagnitude.
the HO matrix elements of the mean fields, the diagonalizationof the Routhi-
ans, and the calculation of the Coulomb field. On a vector computer, the first
of these tasks takes relatively long time (40%) because it is composed of mul-
tiple short loops,18 whereas the last one contributes very little (5%) due to
rather long loops it requires. On superscalarcomputers these three tasks take
almostthe sametimesofabout20%each. The exampleoftypicalworkstation
shown in Fig. 2 illustrates the importance of the optimization options; here
the compiler failed to optimize the subroutine calculating the mean fields and
therefore this task takes now almost 30% of the CPU time.
4 Superdeformed bands in the A 60 nuclei
∼
Very recently there have been several reports of superdeformed rotational
bandsdiscoveredinnucleiaroundA 60.23,24,25,26 Specific nucleiinthisregion
∼
havealreadybeenaddressedfromthepointofviewoftheNilsson-Strutinsky27
and Relativistic Mean Field28 descriptions. The single-particle structure of
6
Execution times in %
Matrix elements
Diagonalisation
Coulomb
Densities HFODD
S
K
S.p. averages
S
A
T Mean fields
CRAY-C90 (x1)
I/O
SG-PCL (x3)
Start & End RISC/6000 (x4)
SUN (x36)
Other
0 10 20 30 40
Figure2: Percentages oftheCPUtimespentonperformingseparatetasks requiredbythe
hfoddcodetosolvetheHFequations.
the light superdeformed nuclei is very simple and thus amenable to a system-
atic description in terms of a microscopic self-consistent mean-field approach.
Using the code hfodd (v1.75),29 hundreds of bands have recently been calcu-
lated in light Ni-Ga nuclei.30 Based on these results one can find and analyze
the generic features exhibited by rapidly rotating nuclei. In the present com-
munication we report on quadrupole polarizabilities characterizing individual
single-particle orbitals in this region.
Following the method introduced in the A 150 nuclei31 we aim at a de-
∼
scription of the proton quadrupole moments Q in terms of contributions
20
givenby the particle (p) and hole (h) states defined with respect to the magic
superdeformed nucleus 60Zn , i.e.,
30 30
Q =Q 60Zn + δQp n + δQh n (12)
20 20(cid:2)30 30(cid:3) X 20 p X 20 h
p h
7
Deviations from additivity
100
hw w =1MeV
80
s 283 bands
d
n
a 60
b
f
o
o. 40
N
20
0
-0.1 -0.05 0 0.05 0.1
d d Q (b)
20
Figure3: DeviationsofcalculatedHFprotonquadrupolemomentsinA∼60nucleifromthe
sumofindividualsingle-particlecontributions (12).
Thecore60Zn correspondstoalln =0andalln =0,whiletheothersuperde-
30 30 p h
formedstates in this regioncanbe obtained by creatingup to 6 particles with
n =1 and/or up to 6 holes with n =1, both for protons and for neutrons.30
p h
Hence, the sum in Eq. (12) contains altogether 24 unknown coefficients δQ
20
whicharefitted by the least-squaremethod to reproduceallcalculatedproton
quadrupole moments Q , separately at every rotational frequency h¯ω.
20
The quality of describing the proton quadrupole moments by the simple
additive contributions is illustrated in Fig. 3 where we show deviations of the
calculatedHFvaluesfromthesumgiveninEq.(12). Onecanseethatformost
of the bands the additivity hypothesis is correctup to 0.05b. This constitutes
aremarkableagreementwiththe valuesofQ whichareoftheorderof2–3b.
20
InFig.4thevaluesofindividualcontributionsδQ areshownasfunctions
20
of the rotational frequency for the 24 orbitals considered. It seen that both
proton and neutron orbitals significantly contribute to the proton quadrupole
moment. This is so because the direct contributionsof protons are ofthe sim-
ilarorderas the polarizationcontributionsexisting forboth kinds ofparticles.
The positive-parity orbitals [440]1/2 and [431]3/2, originating from the 1g
9/2
intruder orbital, have a larger influence on the quadrupole moment than the
negative-parityorbitals(notethatthe scaleistwiceexpandedinthetoppanel
of Fig. 4). Deformations gradually increase when consecutive 1g orbitals
9/2
are occupied. The influence of proton positive-parity orbitals is about twice
larger than that of their neutron counterparts. The negative-parity [303]7/2
8
0.6
n[431]3/2+
(p)
0.4 n[431]3/2-
0.2 (n) n[440]1/2+
n[440]1/2-
0
p[431]3/2+
(n)
-0.2
p[431]3/2-
-0.4 (p) p[440]1/2+
-0.6 p[440]1/2-
0.1 n[310]1/2+
n[303]7/2+
0
n[321]1/2+
)
b
( n[312]5/2+
0 -0.1
2 n[310]1/2-
Q
d -0.2 n[303]7/2-
d
n[321]1/2-
-0.3
n[312]5/2-
0.1 p[310]1/2+
p[303]7/2+
0
p[321]1/2+
p[312]5/2+
-0.1
p[310]1/2-
-0.2 p[303]7/2-
p[321]1/2-
-0.3
p[312]5/2-
0 0.4 0.8 1.2 1.6
hw w (MeV)
Figure4: ContributionsδQ20 ofindividualsingle-particleorbitalstotheprotonquadrupole
moments of nuclei in the A∼60 region. The orbitals are denoted by the dominant Nilsson
labels. Solid (dashed) lines correspond to particle (hole) states, while the closed and open
symbols correspond to positive (r=i) and negative (r=−i) signatures, respectively. Error
barsaredeterminedbytheleast-squareprocedure.
9
orbitals,originatingfromthe1f extruderorbital,significantlycontributeto
7/2
thequadrupolemomenttoo. Inthiscasetheneutronandprotoncontributions
arealmostequal,i.e.,thedirectprotoncontributionisverysmall. Ontheother
hand, the proton contribution from the [312]5/2 orbital is much larger than
thatofthe neutroncounterpart;hence this orbitalcontributesmostly through
the direct quadrupole moment.
Contributionscomingfromthesignature-partnerorbitalsarealmostequal,
except from those related to the [431]3/2 and [310]1/2 orbitals where the
positive-signature partners contribute significantly more. These two orbitals,
and also [303]7/2,exhibit contributions slightly dependent on the angular fre-
quency. However, these variations are relatively much smaller (do not exceed
0.2b) than the variation of the core quadrupole moment Q 60Zn . In-
20(cid:2)30 30(cid:3)
deed, the latter decreases quite substantially with increasing h¯ω, from 3.28b
at h¯ω=0.2MeV to 2.56b at h¯ω=1.6MeV. This is at variance with the results
calculated for 152Dy , where the core quadrupole moment is fairly indepen-
66 86
dent of spin.
Acknowledgments
ThisworkissupportedinpartbythePolishCommitteeforScientificResearch
(KBN) and by the computational grants from the Interdisciplinary Center for
Mathematical and Computational Modeling (ICM) of the Warsaw University
andfromtheInstitut du D´eveloppement et de Ressources en Informatique Sci-
entifique (IDRIS) of CNRS, France (Project No. 960333).
References
1. P. Ring and P. Schuck, The Nuclear Many-Body Problem (Springer-
Verlag, Berlin, 1980).
2. P.QuentinandH.Flocard,Annu. Rev. Nucl. Part. Sci. 28,523(1978).
3. S. ˚Aberg, H. Flocard and W. Nazarewicz, Ann. Rev. Nucl. Part. Sci.
40, 439 (1990).
4. C. Baktash, B. Haas, and W. Nazarewicz, Annu. Rev. Nucl. Part. Sci.
45, 485 (1995).
5. P. Butler and W. Nazarewicz, Rev. Mod. Phys. 68, 349 (1996).
6. D. Vautherin and D.M. Brink, Phys. Rev. C5, 626 (1972).
7. D. Gogny, in: Proc. Int. Conf. on Nuclear Physics, eds. J. De Boer and
H. Mang (North-Holland, Amsterdam, 1973).
8. D. Gogny, in Nuclear Self-Consistent Fields, ed. G. Ripka and M.
Porneuf (North-Holland, Amsterdam, 1975) p. 333.
10
