Table Of ContentSoliton Systems at Finite Temperatures and Finite Densities
Oliver Schwindt and Niels R. Walet
Dept. of Physics, UMIST, Manchester M60 1QD, UK
(Dated: February 1, 2008)
Thefinite-densityandfinite-temperaturephaseportraitsofthebaby-SkyrmeandSkyrmemodels
areinvestigated. Bothgrand-canonicalandcanonicalapproachesareemployed. Thegrand-canonical
approach can be used to find the “natural” crystal structure of the baby-Skyrme model and it is
showntohavetriangularsymmetry. Thephaseportraitsincludesolid,liquidandphase-coexistence
between solids and vacuum states. Furthermore, a chiral phase transition can be observed for
both models. The phase portrait of the Skyrmemodel is compared with that expected of strongly
interacting matter, and we also contrast ourresults to othermodels.
PACSnumbers: 12.39.Dc,21.65.+f,05.10.-a,73.43.-f
2
0
0 I. INTRODUCTION of the phase diagram is the presence of a critical point.
2 This separates an area where a first-order phase tran-
n Oneofthegreatquestionsinthetheoryofstronginter- sition occurs, from one where a continuous path can be
a takenfromthehadronicphasetothequark-gluonplasma
actions is the behavior of baryonic matter at finite den-
J phase. Atthecriticalpointitself,itisbelievedthatthere
sity and finite temperatures. Even though it is widely
2 is a second-order phase transition. The position of the
believed that quantum chromo-dynamics (QCD) is the
2 criticalpointisextremelydifficulttomodelandhasbeen
correct theory to describe the strong interactions, its
1 non-perturbative nature makes it extremely difficult to estimatedtohaveatemperatureTE ≈140−190MeVand
a chemical potential of µ ≈ 200−800MeV [3, 7]. The
v describe the low energyconsequencesof the theory. Lat- E
3 ticeQCD[1]isthepreferrednumericaltechniquetostudy color superconducting phase cannot be modeled by the
0 the QCD Lagrangian directly. Here, QCD is modeled Skyrme model because color-flavor locking is observed,
2 and in order to model it, three colors need to be consid-
non-perturbatively by discretising the theory on a lat-
1 ered rather than the N = ∞ approximation underlying
tice, and a numericalapproximationis made to the path c
0
the Skyrme model, which will be discussed later. The
integral describing either the Euclidean time evolution
2
color superconducting phase will therefore not be dis-
0 or the thermodynamics of the theory. The masses of
cussed further.
/ mesons, baryons, and glueballs have already been cal-
h
culated to reasonable accuracy [2]. The limitations of The Skyrme model has a venerable history [8] in the
p
- the approach are both numerical (momentum cut-off, fi- non-perturbativedescriptionofnucleonstructureandthe
p nite lattice size), and more theoretical. The most im- low-energybehaviorofbaryonicmatter,sinceitcontains
e
portant difficulty is the description of fermions on the a good description of the long-wave length behavior of
h
lattice. Since lattice QCD relies on using Monte Carlo thedynamicsofhadrons. Ashasbeenarguedby’tHooft
:
v methods, computing power is another of the key limi- andWitten[9,10,11],thisiscloselyrelatedtothelarge-
i
X tations. If sufficient powerful computers were available, number-of-colors limit of QCD, in which baryons must
the inaccuracies due to discretising QCD could be re- emerge as solitons, in much the same way as happens in
r
a duced by adopting larger gridsizes. Therefore, as com- the Skyrme model. Alternatively, we can interpret the
puting power increases in future, more accurate results model in terms of chiral perturbation theory, in terms
will be achieved. Even then, dealing with a finite den- what is much like a gradient expansion in terms of the
sity offermions is anespecially tricky problem,since the pion field. The model can be used to describe, with due
fermionicdeterminanttobecalculatedbecomescomplex, care[12,13],systemsofafewnucleons,andhasalsobeen
leading to various instabilities. applied to nuclear and quark matter. Within the stan-
Due to these difficulties, it is quite common to use dard zero-temperature Skyrme model description there
models to study the phase diagram of strongly interact- are signatures of chiral symmetry restoration at finite
ing matter. Recent examples of such approaches can for density, but in a rather special way, where a crystal of
instancebefoundinstudiesbyRajagopaletal[3,4],who nucleons turns into a crystal of half nucleons at finite
describe the phase diagram using a Nambu-Jona-Lasino density,whichischirallysymmetriconlyonaverage[14].
model [5] and a random matrix model [6]. The QCD Thequestionofthefinite-temperaturehasneverbeenad-
phase diagram, obtained from a model with two mass- dressedandwouldbeofsomeinterest. Weshallalsolook
less quarks in Ref. [3]. A typical phase diagram shows atthe two-dimensionalSkyrmemodel,whichhasthead-
a hadronic-matter phase, the quark-gluonplasma, and a vantageofbeingeasiertointerpret,butalsohasphysical
colorsuperconductingphase. Earth’ssurface,ingeneral, relevance; a special form of the two-dimensional Skyrme
has a low temperature (0.025eV) and a chemical poten- model has recently been developed for use in quantum
tialmuchlessthan1GeV,andthereforeweliveina(liq- Hall systems [15]. This model is obtained as an effective
uid) hadronic-matter phase. The most striking feature theory when the excitations relative to the ν = 1 ferro-
2
magnetic quantum Hall state are described in terms of tuations (with care, because nonremormalisability plays
(a gradient expansion in) the spin density, a field with a role in the thermodynamics as well!). Furthermore, we
properties analogous to the pion field in nuclear physics might well wish to study the physics of these systems at
[16]. Apart from obvious changes due to the number of temperatures where thermal fluctuations dominate the
dimensions, the new approach differs from the histori- physics. Therefore, we shall concentrate on the finite
cal Skyrme model by having a different time-dependent temperature phases of these classical field theories.
term in the Lagrangian, and the appearance of a non- The most direct approach to the problem is to per-
localinteraction,wherethe topologicalchargedensity at form a Monte-Carlo (Metropolis algorithm based) study
different points interacts through the Coulomb force. In of the partition function of each of these models. As we
the limit of large Skyrmions this last term can be ap- shallarguebelow the leastobviousaspectin suchanap-
proximated by a more traditional local “Skyrme term”, proach is how to deal with the topological conservation
which is quartic in the fields, leading to the standard laws. Rather than immediately tackle the 3D model, we
baby-Skyrme model with local interactions. Even with shall first concentrate on the (local) 2D model. This has
the non-local complication the model is on the whole the advantage that the visualizationof, and also the un-
remarkably similar to the nuclear Skyrme model. The derstandinggainedby,theresultsismuchmorestraight-
Skyrme-field effective degrees of freedom describe the forward. The extensions to quantum-HallSkyrmionsare
ground state of such systems, and probably also the low under way and will be presented elsewhere [19]. We ex-
energy dynamics and thermodynamics, so that we can pect the current for the QHE solitons to be very similar
ask similar questions about the finite-density physics as to the results reported here.
for baryonic matter. Thispaperisorganizedasfollows. InSec.IIwediscuss
In this paper we shall consider the local baby Skyrme the formalismsuccinctly. In the next section,Sec. III we
model, which can be regarded as the 2D version of the discuss in detail the results obtained in our simulations
traditional 3D nuclear one. We can also extrapolate to (a smallexcerpt has alreadybeen published in Refs. [20,
one dimension. The model obtained in that way is the 21]). We draw some conclusions,andgive anoutlook for
sine Gordon model expressed in terms of an unit com- possible future work in Sec. IV.
plex field. All of these models contain a topologically
conserved charge, which in the three-dimensional case is
identified with baryon number. This gives rise to topo- II. FORMALISM
logicalsolitons,whichinnuclearphysicsareidentifiedas
baryons. Since the behavior of nuclear matter at high A. Lagrange density
densities is expected to reveal, probably by one or two
phasetransitions,thesubstructureofbaryonicmatter,it The Lagrangedensity forSkyrme models canbe given
is interesting to study this phenomenon in Skyrme mod- as [8, 22]
els. Atthesametimeitallowsustolookatthephasedi-
agram of quantum Hall Skyrmionic systems, which may 1 1
L = (∂ φ )2− (∂ φ ∂µφ ∂ φ ∂νφ )
well be the easiest way to calculate parts of the phase 2 µ k 4 µ k k ν l l
diagram of the underlying electronic model. 1 m2
+ (∂ φ ∂µφ ∂ φ ∂νφ )+ π (φ −1) , (1)
At finite density but zero temperature the Skyrme 4 µ k l ν k l f2e2 0
π
models and the sine Gordon model have a crystalline
structure, which consist of regularly spaced solitons. In which as long at we take the space index µ=0,1,... ,d
the one-dimensional case, where exact solutions exist todescribed+1dimensionalspace,andthefieldcompo-
[17], thisis knownto be azerotemperature artifact,and nents k,l =1,... ,d+1 to describe a d+1 dimensional
wehavealiquidatanyfinitetemperature,nomatterhow unitvector,candescribetheSine-Gordontheory(d=1),
small. InthenuclearSkyrmemodelwewouldliketofind the baby Skyrme model (d=2) and the nuclear Skyrme
afluid,tomimicthequantumliquidbehaviorexpectedin model (d = 3). The Lagrange density can be separated
nuclearmatter [18]. Inorderto appreciatethe subtleties into a kinetic and a potential term
involved, one must understand that the Skyrme model,
as a classical field theory, can be understood as a semi- L=T[φ]−V[φ] , (2)
classical(largeaction)limitofaquantumtheory. Clearly
the quantum fluctuations in the underlying theory could where
be large enough to washout the crystalline structure, as
1
happens for the sine Gordon model. In the special case T[φ] = ∂ φ ∂ φ M [φ] , (3)
t a t b ab
2
ofonespacedimensionitisalsowellknownthatthermal
M [φ] = δ +∂ φ ∂ φ δ −∂ φ ∂ φ , (4)
andquantumfluctuationsplayexactlythesamerole,and ab ab i c i c ab i a i b
thermalfluctuations also breakthe crystallinestate. For V[φ] = 1∂ φ ∂ φ + m2π (1−φ )
theSkyrmeandbabySkyrmemodelsitisnoteasytoac- 2 i a i a f2e2 0
π
cessthequantumfluctuations,sincethefieldtheoriesare 1
non-renormalisable, but we can access the thermal fluc- + ([∂iφa∂iφa]2−[∂iφa∂jφa][∂iφb∂jφb])(5.)
4
3
For our later calculation we need we need the expres- which we shall refer to as the scaled energy, is
sionfor the energydensity, T +V, which for the Skyrme
dMd
model becomes E =E −E =E − . (9)
scaled total vib total
β
1
E = ∂tφa∂tφbMab(φ(x))+V(φ(x)) . (6) It is clear that these vibrational contributions are un-
2
physical and need to be removed, since in the contin-
The two sets of constraints (implemented by the delta- uum limit their contributionwouldbecome infinite. The
functions in the integrals) are φ2−1 = 0 and φ˙ ·φ = 0 sine-Gordon model is renormalisable, and by subtract-
ateachpointinspace,wherethefirstconstraintforcesφ ing the energy due to the lattice vibrations, we obtain
to have unit length and the second constraint follows by the contribution due entirely to the solitons in a similar
taking the time derivative of the first constraint. way as was done analytically in Ref. [17]. The baby-
In3Dthetopologicalchargedensity,whichwasidenti- Skyrme and Skyrme models are non-renormalisable,and
fiedbySkyrmewiththebaryon-numberdensityoperator, we find that it is therefore not possible to remove all
is given by thelattice-dependencyfromtheenergyexpression. How-
ever,byremovingtheharmoniccontribution,thelattice-
1 dependency is reduced. We will not use the energy of a
B(x) = ε0σρνεαβγδφ ∂ φ ∂ φ ∂ φ
12π2 α ν β σ γ ρ δ system to examine phase transitions, because of these
1 lattice-dependencies. Instead, we use other, less lattice-
= det(φ,∂ φ,∂ φ,∂ φ) , (7)
2π2 x y z dependent methods, which we will discuss in Secs. IIIB
and IIIC.
with simpler expressions in 1D and 2D. The energy den- Itisusefultousedensity-likequantitiesratherthanto-
sity(6)andbaryon-numberdensity(7)aretheninserted tal quantities, and therefore we often use quantities like
into the grand-canonicalpartition function (8). the average potential-energy density V = V and the
hdMd
average baryon density B = B , where hdMd is the
hdMd
volume of the system. Results are even more compara-
B. Partition Function ble when measuring V (or V), which gives the potential
B B
energy per soliton. By measuring density-like quanti-
In order to study the thermodynamics of a theory, we ties,wearemeasuringquantitieswhichremainfiniteand
must construct a partition function. We shall discuss meaningful in the limit of an infinite simulation volume.
both grand-canonical and canonical partition functions, The chemical potential µ determines the particle den-
whicharealsodiscussedinRefs.[23,24]. Theexplicitex- sity in a particular system. We initially expected that a
pressions for the thermodynamic partition functions will simulation of the grand-canonical partition function au-
be given for the baby-Skyrme and the Skyrme models. tomatically includes the probability of adding or remov-
The Metropolis principle is used to evaluate the parti- ing a soliton from the system. Unfortunately, because
tion functions, see Refs. [23, 24] for more details. we are working with extended particles, a soliton can
only be added or removed from the system if a number
of field vectors on our discretised lattice change simulta-
neously. The probability of this happening is extremely
1. The Grand-Canonical Partition Function
small, and has never been observed in our simulations.
Theonlyexamplesofthe totalparticlenumberchanging
The general thermodynamic partition function for a
is if the numerics become unstable. In this case, neigh-
general unitary field theory on the lattice is given by
boring vectors are pointing in quite different directions
such that the derivatives, and therefore the energy and
Md
baryon densities, are not calculated correctly. The nu-
Z = dd+1φ dd+1φ˙ δ(φ˙ ·φ )δ(φ2−1)
p p p p p merical breakdown, i.e. the unphysical loss of winding
Z pY=1 number, is dependent on the lattice spacing, because it
×exp(−β(E−µB)) , (8) occurs at lower densities and temperatures if the lattice
spacing is increased. Such numerical breakdown must
where an integration is to be performed over the field
occur when we have less than the minimum number of
and field derivative at each lattice point. The quantities
lattice points that can describe a Skyrmion.
E and B are the total energy and baryon number of the
To solve the problem that the number of particles in
system.
a given simulation is fixed by the topology, we use an
Ifthefieldmakessmallvibrations(whichmustbehar- opensystem,whereparticlescanflowthroughthebound-
monic) at eachlattice point then eachdegree of freedom aries freely and therefore we have a heatbath that acts
contributes 12kBT to the free energy. There are Md lat- as a source of particles. The method we use to imple-
tice sites, each with d potential and d kinetic degrees ment our open system will be discussed in Secs. IIIB,
of freedom. The total energy due to lattice vibrations and IIIC, when we discuss the thermodynamics of our
is thus Evib = 2d2Mβd. The energy due only to solitons, soliton models.
4
2. The Canonical Partition Function withinthe canonicalapproach. Forsuchstatic solutions,
a renormalisation is irrelevant. Furthermore, the kinetic
The canonical approach is usually applied to a closed energy is irrelevantsince the total energy is equal to the
system, i.e. where particles are trapped in a box, and potential energy.
the number of particles in this box remain constant. As We now use the fact that φ˙ is an eigenvector of the
alreadymentionedinthelastsection,thenumberofsoli- mass matrix (4) to perform the integral over φ˙ analyti-
tons is conserved by the topology of the system. There- cally. The partition function now simplifies to
fore,the canonicalapproachseemsto be verysuitable to
investigate the thermodynamics of soliton systems. The Md
average density of particles is therefore defined by the Z =N−1β−dM2d dd+1φp 1+∂iφc,p∂iφc,p
initial choice of the number of particles B in the volume Z p=1 s det(Mp) !
Y
V and it never changes throughout the simulation. We
×exp(−β(V −µB )) ,
choosenottoimplementfixedboundaries,whichhavethe p p
vacuumvalueonthesimulationedges,butinsteadweuse
(12)
periodic boundaries. A system with periodic boundaries
canbeinterpretedasanapproximationofaninfinitesys-
where V and B are the potential energy and baryon
tem because the simulation box is effectively duplicated p p
number at the lattice site p and represent the potential
an infinite number of times. Periodic boundary condi-
energy and baryon number in the surrounding unit cell.
tions imply thatif a Skyrmionflows throughone bound-
The overall constant factor N−1 is irrelevant when ap-
aryofthesimulation,thenitreappearsthroughtheoppo-
plyingtheMetropolisprinciple,andwillbe ignoredfrom
site side. Periodic boundary conditions have two advan-
now on. Eq. (12) will therefore be used as the grand-
tages over fixed boundaries; firstly that the simulation
canonical partition function. The lattice does not need
can be interpreted as an infinite system, and secondly
to be cubic, but we haveadoptedthe conventionthatall
thatweavoidunphysicalboundaryeffectsoccurringfrom
the sides have equal length in order to avoid the addi-
Skyrmions interacting with the fixed boundaries.
tional notation required when using unequal lengths.
Thegeneralexpressionforthecanonicalpartitionfunc-
tion is The procedure discussed above can be applied to the
canonical formalism, and as expected we obtain the ex-
Mn pression (12) with µ=0.
Z = dn+1φ dn+1φ˙ δ(φ˙ ·φ )δ(φ2−1)
B p p p p p The partition function (12) cannot be evaluated ana-
Z pY=1 lytically,andthereforeaMetropolisalgorithmisapplied,
×exp(−βE) , (10) see Refs. [23, 24] for more details.
The only exception is the sine-Gordon model. The
where all the symbols are again defined as in Sec. IIB1.
grand-canonicalpartitionfunctionisderivedinRef.[17].
The harmonic contributions to the potential energy are
Itcanveryeasilybecomparedtooursimulations,sinceit
removed in a similar manner as for the grand-canonical
is obtained by discretising the field and then taking the
approach,
continuum limit. Like many other field-theory models,
nMn the observables tend to infinity as the continuum limit
Vscaled =Vtotal− . (11) is being reached. It thus illustrates clearly the role of
2β
renormalisation.
In Ref. [24] the zero-temperature minimal-energy solu- The canonical partition function of the sine-Gordon
tions solutions are calculated using simulated annealing model in the angle representation is given by
M
1 1
Z = d2α d2α˙ exp −β( ∂ α ∂ α + ∂ α ∂ α +1−cosα )h
B p p t p t p x p x p p
2 2
Z p=1 (cid:18) (cid:19)
Y
M
1 1
= const d2α exp −β( ∂ α ∂ α +1−cosα )h , (13)
p x p x p p
β 2
r Z p (cid:18) (cid:19)
Y
where α is the field, and 1M originates from evaluat- latticespacinghisrequiredtoconverttheenergydensity
β tothe energyinacellfromthe discretelattice. Thecon-
ing the time derivative anqalytically, see Sec. IIB1. The stant factor in the expression for Z is irrelevant when
B
5
measuringthermodynamicquantities,andisthereforeset
equalto one. The reasonevery lattice point behaves like
a classical harmonic oscillator is that the kinetic energy
and potential energy are quadratic in time derivatives
0.4
ρ
andspatialderivatives,respectively. Thecontributionto 10
0.2
the fluctuations at individual lattice points from the po- 7.5 0
tential 1−cos(α) is negligible. The energy at a single β 5
1155
lattice point due to fluctuations in the field is given by 2.5 1100µ
(αi+21h−2αi)2 ×h, or (αi+12−hαi)2, where i is the lattice index 00 55
and h is the lattice spacing. Therefore,
2π M/2 M (α −α )2
i+1 i
ZB = β exp −β 2h 0.4 ρ
(cid:18) (cid:19) Z i=1 i ! 10 0.2
Y×dα ...dXα 7.5 0
1 M β 5
1155
2π M/2 const×h M/2 2.5 1100µ
≈ . (14) 55
β β 0 00
(cid:18) (cid:19) (cid:18) (cid:19)
FIG. 1: The exact (top figure) an numerical (lower figure)
Using P = 1 ∂ZB, the harmonic oscillator contribution
βL ∂β results for the grand-canonical equation of state of the sine-
to the pressure becomes P = M ln(const×h). Since the Gordon model.
Lβ β
harmonic oscillator at every lattice point gives a finite
contributionto the pressureP,it becomes infinite inthe
III. RESULTS
continuumlimit. Theconstantisnotevaluatedexplicitly
here because its contribution to the energy which we are
interested in is irrelevant. Furthermore, the lattice spac- A. Sine-Gordon model
ingh withinthe logarithmalsovanishes. The energybe-
coming infinite in the continuum limit can be overcome Acheckofthe calculationscanbe made by comparing
by removing the contribution from the harmonic oscil- theexactresultsasdiscussedinthe previoussectiontoa
lators. This way, all quantities are calculated without calculationusingopenboundaryconditions. This means
contributions from the harmonic oscillators,and contain thatweimposenoboundaryconditionsatallonthefinite
only contributions that are not divergent in the thermo- simulationvolume,butallowsolitonstoenterthesystem
dynamic limit. Therefore, the pressure due entirely to throughthoseboundaries. Eventhough,inprinciple,one
the solitons is given by cancreate topologicalchargeon a grid,in practice topo-
logical charge is very well conserved. The only way to
perform a grand-canonicalsimulation in practice is thus
A −iµβ,4β2
′ π using open boundary conditions. We then compare the
P =− +1 . (15)
h 8β2 i simulatedresultforthedensityofsolitonstotheanalytic
result, Eq. (16), in Fig. 1 as a function of β and µ. As
wecanseetheresultsareverysimilar,withsomestatisti-
cal noise in the simulated results. Since it is much more
We find that the density and internal-energy density
straightforwardtoperformlargesimulationsforthesine-
foragivenβ andµ,whichareobtainedbyusingthestan-
Gordon model than for the 2D and 3D Skyrme models,
dard thermodynamic relationships [25, 26] ρ = ∂P
∂µ thisshowsthebestqualityofresultsthatcanbeexpected
β
and u=− ∂(Pβ) respectively, are (cid:16) (cid:17) (apartfromalimitationinsimulationtime,sincewehave
∂β µβ looked at quite a dense set of grid for the results).
(cid:16) (cid:17)
∂A −iµβ,4β2
1 π B. Thermodynamics of the Baby-Skyrme Model
ρ=− (16)
8β2 h ∂µ i
1. The Grand-Canonical Approach
and
The specialfeature of the grand-canonicalapproachis
that we use open-boundary conditions. This means that
∂ 1 A −iµβ,4β2 +β
8β π thenumberofparticlespresentinaparticularsimulation
u= . (17)
(cid:16) h ∂β i (cid:17) canchangewhen particlesenter or leavethe systemover
µβ the boundary. The lattice points at a boundary behave
6
0.175
0.15
0.125
ρ 0.1
0.15 0.075
0.03 0.1 0.05
0.05ρ 0.025
0.02 0
T0.01 4400 5 10 15 20µ25 30 35 40
3300
0 2200 µµ 0.175
0.15
0.125
ρ 0.1
0.075
0.05
0.15 0.025
0.1
0.03 0.05ρ 5 10 15 20µ25 30 35 40
0.02 0
T 4400
0.01
3300 FIG. 3: The hysteresis effect when increasing and then de-
0 2200 µ creasingthechemicalpotentialµ. Top:Thesolidlines,where
thearrowsmarkthedirectionofthesimulation,showthehys-
teresis effect for a low temperature T =0.0001. The dashed
FIG. 2: The µ, ρ, T phase portrait for the baby-Skyrme lines are the results imported from Fig. 2, where the tem-
model. Top: The order in which the points in the plot are perature was changed for constant µ. Bottom: For higher
calculated is in the increasing T direction with constant µ. temperatures, the hysteresis effect is less visible. The axes
Bottom:ThecalculationorderisinthedecreasingT direction are labeled in Skyrmeunits.
with constant µ. The plots differ because the initial field
configuration at each point is the final configuration of the
previously calculated point. The axes are labeled in Skyrme hysteresis effect is more visible than for higher tempera-
units. tures,becausesolitonsdonotmoveasquicklyandhence
it takes longer to reach equilibrium. The magnitude of
thehysteresiseffectdependsonthecomputingtimeused,
slightly differently than the internal ones. If a field vec- asequilibriumisreachedextremelyslowlyandwedidnot
tor is sampled at a boundary, it is accepted or rejected have the time to let the system strictly reach it.
depending only on the change in the integrand on the Each baby-Skyrmion has to have an energy greater
neighboringinternalplaquettes(weusethewordplaque- than the topological lower bound for it to exist in an
tte for a square enclosedby four nearestneighbor lattice open system simulation. Fig. 3 shows that the chemical
points). Thefieldvectorsattheboundariesaretherefore potential,whichdefineshowmuchenergyisgiventoeach
less restricted than those at other lattice points. soliton, must be greater than approximately 18 Skyrme-
Theµ, ρ, T phasediagramforthebaby-Skyrmemodel energyunitsperbaby-Skyrmion,µ&18. Oncethechem-
shown in Fig. 2. Each point is a separate simulation ical potential is greater than the threshold, the density
which has reached equilibrium for a given chemical po- increases approximately linear with the chemical poten-
tential µ anda temperature T. To savecomputing time, tial, ρ∝µ−µ .
threshold
the initial field configuration used to calculate the den- As long as µ is large enough, which is the case for al-
sityρateachpointintheplotwasthefinalconfiguration most all of the finite-ρ region, the state of the system is
forthepreviouslycalculatedtemperature. Thearrowsin crystalline. Theopenboundariesallowcrystalstructures
Fig.2showtheorderthesimulationswherecalculatedin. to form without defects, and therefore we can determine
Two plots of the same phase portrait have been shown the “natural” crystal structure which will be described
where the order of calculation has been reversed. The in the next section. The liquid state and the phase co-
reasonforshowingthisisbecausethecomputingtimere- existencebetweensolidandvacuumstatesexistnearthe
quiredtoreachequilibriumbecomeslargenearthephase phase transition and are too difficult to model using the
transition, and the simulation time was not long enough grand-canonicalapproach. Thesestateswillbeexamined
to reach equilibrium near the phase transition. By com- using the canonical approach in Sec. IIIB2.
paring the two graphs, one can see that it is difficult to For high densities and low temperatures, the
predictthe chemicalpotentialofthe transitionto within Skyrmions merge into a lattice where it is not possible
5Skyrme-energyunitspersoliton. InFig.3,thenetden- to identify individual Skyrmions, see Fig. 4. Since there
sity of baby-Skyrmions ρ is shown against the chemical are no fixed boundaries in this simulation, there are no
potential µ at a fixed temperature T, for two different defects in the crystalline solutions either. The lattice is
temperatures. Ineachplot,theorderinwhichµchanges determinedmosteasilyfromthe structureofthe “holes”
hasbeenshown;oneintheincreasingµdirectionandone with zero baryon-number density. There are two such
in the decreasing µ direction. For low temperatures, the minima present per baryon. In this crystal phase, the
7
12
20.5
2
10 20
E/B 19.5
8 19
18.5 1
6
0.02 0.03 0.04 0.05 0.06 0.07
ρ
4
1
2
0
0 2 4 6 8 10 12
FIG. 4: Baryon-number density plots of the natural crystal
structureofthebaby-Skyrmemodelwithopen-boundarycon-
ditions. The baryon-number density is averaged over 20000
steps. The axes are labeled in Skyrme-length units. Lighter
2
shadesrepresenthigherbaryon-numberdensities. Theparam-
etersforthissimulationareµ=30,T =0.05,andthenumber
ofsolitonsinthesimulation volumeisabout13(rho=xxx).
FIG.5: Thehigh-densityphaseatzerotemperature. Top:En-
average field hφi is zero, and therefore a mean-field chi-
ergy per baryon against density. Middle: A baryon-number
ral symmetry exists.
density plot of the lowest energy per baryon number attain-
The triangular nature of the crystal is evidence that able, i.e. the“natural” crystal structure. Lighter shades rep-
the type of crystal favoredby the baby-Skyrme model is resent higher baryon-number density. Bottom: A baryon-
independent of the lattice used to discretise the model. number density plot at high density. The baryon number is
Eventhoughthelatticeissquare,thecrystalwhichforms more evenly distributed throughout the structure. The mid-
doesnotalignitselftoit. Sincethissimulationallowsthe dle and right plot are drawn to the same scale.
crystal structure having the lowest energy per baryon to
form without being effected by boundary conditions, we
call it the “natural” crystal structure. The field config- 2. The Canonical Approach
uration within a unit cell was copied (since it satisfies
periodic boundary conditions) and simulated annealing Aswehavearguedbefore,withperiodicboundarycon-
withperiodicboundaryconditionswasusedtoaccurately ditions it seems impossible to generate additional topo-
investigate the zero-temperature properties. logical charge. Nevertheless these calculations are much
The correct high-density structure [20] where preferable over the open boundary calculations, since
Skyrmions merge into a lattice is shown in Fig. 5. The the results suffer much less from artifacts due to those
lattice is determined most easily from the structure of boundaries (but the results will now depend much more
the “holes” with zero baryon-number density, since it is stronglyonthe symmetriesofthe initialseedtothe sim-
not possible to identify individual Skyrmions. There are ulation). This is still the preferred way to calculate the
two such minima present per baryon. In the low-density phase diagram. In order to create the ρ, T, u phase
ρ
phase, the averagefield hφki is polarizedin the direction diagram, the lattice-dependent terms must be removed
ofthesigmafield,butforhighdensities,theaveragefield from the internal energy. Unfortunately, the lattice-
is zero. Therefore, a mean-field chiral symmetry exists dependency cannot be removed entirely, because the
in the high-density phase and not in the low-density baby-Skyrme model is non-renormalisable.
phase. The energy per baryon ratio when minimizing The internal energy per unit area u is found by using
the energyofafieldconfigurationhavingthe formofthe the thermodynamic relationship
naturalcrystalis 18.27Skyrme-energyunits ata density
of 0.024 particles per cubic Skyrme-length unit. If the 1 ∂lnZ
u=U/L2 =− , (18)
density is lowered further, then the minimal-energy L2 ∂β
(cid:18) (cid:19)B
solution consist of small numbers of Skyrmions bound
together, i.e. the multi-Skyrmions. In fact, this can be whereL2isthesimulationvolumeandBdenotesthecon-
defined as the point where the phase transition between stant particle number defined by the initial conditions.
the high-density and low-density matter occurs. Therefore,
8
u = M2 + 1 pd3φp 1+∂dieφtc(,Mp∂pi)φc,p Vpexp(−β(Vp)) = 1 + 1 hVi , (19)
L2β L2R Q d3φ(cid:16)q 1+∂iφc,p∂iφc,(cid:17)p exp(−β(V )) h2β L2
p p det(Mp) p
R Q (cid:16)q (cid:17)
where h is the lattice spacing and L12hVi is the potential T 1111..552
energy per unit area. 00..55
00
One can attempt to remove some of the lattice-
dependency by assuming that the fluctuations at lattice
pointsbehavelikeharmonicoscillators,andthereforeas- 50
suming that the contribution from the potential energy u
of the harmonic oscillator, 1 , can be removed, so ρ
h2β 40
1 1
u = hVi− . 30
scaled L2 h2β
Unfortunately uscaled still depends strongly on the lat- 00 00..002200..0044ρ00..006600..0088 0.120
tice spacing. Furthermore,forhightemperatures,u
scaled
can become negative, making its interpretation difficult.
FIG. 6: The ρ, T, u phase diagram. Some of the contribu-
If the number of solitons in the system is zero, uscaled ρ
tion from finite-lattice effects is removed by using Eq. (20)
still has a non-linear behavior that can become nega-
to find the internal-energy density. The density ρ is altered
tive, when we would expect it to be zero. This is clearly
bychangingthelatticespacingandkeepingthesamenumber
an indication of the non-renormalisabilityof the Skyrme
of Skyrmions in the system. The axes are labeled in Skyrme
model, and we can no longer use the simple technique
units. (16 solitons on a 90×90 lattice.)
(countthelatticedegreesoffreedom)usedtorenormalise
the sine-Gordon model
Insteadofsearchingfora detailedtheoreticalexplana-
Atfinite temperatures,vibrationsateachlatticepoint
tion of the lattice-dependency, which would require the
are present, but the contribution to the energy depends
study of counterterms, we choseto comparesimulations
on the number of lattice points per unit volume. Since
withafixednumberofsolitonswithsimulationsthatuse
we vary the lattice spacing to changethe solitondensity,
the samelatticeparametersandareatthe sametemper-
theinternalenergypersolitonasshowninFig.6maynot
ature, but contain zero solitons. Thus, we use
be properly compared for different densities. This phase
diagram will not be discussed further, because it cannot
1 1
uscaled = L2hVi− L2hV0i , (20) be used to interpret the states of matter throughout the
ρ-T plane. In the next two sections, the methods that
where hVi is the potential energy of a simulation with successfully determine the states of matter for the baby-
solitons, and hV i is the potential energy of the same Skyrmemodelthroughouttheρ−T plane arediscussed.
0
simulation without any solitons. Although there is still Along the lines first proposed by Klebanov [27], we
somelattice-dependency,itisnotverystrong. Itislikely analyse the grid-averaged fields hφi. The average of the
that it is not possible to completely remove the depen- pion fields hπi is always zero, but the average sigma
dency on the lattice spacing because the model is non- field hσi changes. The average sigma field hσi is plotted
renormalisable. Nonetheless,weusethismethodbecause against the temperature T and the density ρ in Fig. 7.
the internal energies for different lattice spacing become When few Skyrmions exist in a large volume, hσi>0.
quite comparable, and phenomena such as phase transi- This is because the vacuum, where σ = 1, contributes
tion points do not seem to depend on the lattice scale a significant amount to the average sigma field. For
chosen. high densities, however, the baby-Skyrmions combine in
Wefirstuseacanonicalapproachtostudytheρ, T, u such a way that the individual baby-Skyrmions cannot
ρ
phasediagram. Theinternalenergypersolitonisplotted be identified. Furthermore, an approximately chirally-
against the density ρ and the temperature T in Fig. 6. symmetric configuration is formed where hσi ≈ 0 and
The internal energy is measured using Eq. (20). The the sigma field and the pion fields are interchangeable
density ρ is altered by changing the lattice spacing h without change in energy. For the baby-Skyrme model,
while keeping the number of Skyrmions in the system at chiral symmetry is never exactly satisfied (in a time av-
a constant value. Although Eq. (20) is used to remove erage), because the pion-mass term, which is required
some of the lattice-dependency, there is still a remaining to stabilize the solitons, violates this symmetry. As the
contribution. density of a system is increased, hσi = 0 is approached
9
2
T 11..55
11
00..55
00
1
0.8
<σ>
0.6
0.4
0.2
00 00..002200..0044ρ00..0066 00..0088 0.10
FIG. 7: The chiral symmetry (magnetization) hσi as a func- FIG. 8: A simulation in the solid phase. Pseudo-time-
tion of ρ and T. For fluids and solids, we find hσi ≈ 0. For averaged baryon-number density plot. Lighter shades rep-
the phase coexistence between solid and vacuum,hσi>0. T resenthigherbaryon-numberdensityandtheunitsshownare
andρaregiveninSkyrmeunitsandhσiisdimensionless. The in Skyrme-length units. The density ρ=0.04, which implies
parameterhσihasbeenmeasuredfromasinglefieldconfigura- thateachlight-coloredcirclecontainstwounitsoftopological
tionandhasnotbeenaveragedthroughtime,andtheresults charge.
show large fluctuations. All simulations were performed for
16solitonsona90×90lattice,andthedensitywasvariedby
changing thelattice spacing.
32
asymptotically. InFig.7,aphasetransitioninhσicanbe
24
observed to occur between the densities 0.02<ρ < 0.04
for various temperatures.
Thefluctuationsinhσiarelargestintheregion0.02< 16
ρ < 0.08, T > 0.5. In the next section, we identify a
liquid and a solid phase in the chirally symmetric phase,
8
andthe liquidregionis locatedwherethe fluctuations in
hσi are largest. The liquid and solid states possibly ex-
tend into the broken chiral symmetry phase, depending 8 16 24 32
on where one defines the chiral phase transition to be.
These fluctuations are present because in Fig. 7 hσi was FIG. 9: A simulation in a liquid phase, as a pseudo-time-
calculatedfroma single field configurationand not aver- averagedbaryon-densityplotwithanaveragetakenwith2000
agedintime. Also,thehσiphaseportraitseemstobein- pseudo-time steps. The shades from dark to light represents
lowtohighbaryon-numberdensitiesandtheunitsshownare
dependent of the lattice parameters, unlike the internal-
in Skyrme-length units. As in Fig. 8 ρ = 0.04, but the tem-
energy density. When doubling the number of lattice
peratureis higher.
point and keeping the same number of baby-Skyrmions
in the same volume, the hσi phase diagram is not signif-
icantly different than the one shown in Fig. 7.
can be observed although defects are present, which are
By examining the pseudo-time-averaged baryon-
a consequence of using periodic boundary conditions, in
density plots at various temperatures and densities, one
contrast to the absence of defects seen when using the
can easily identify three different states that a particu-
grand-canonicalapproach.
lar simulation may be in, namely in a solid, in a liquid,
orinaphase-coexistencestate. Furthermore,correlation Apseudo-time-averageofaliquidshowsanalmostcon-
functions may be used to draw the same conclusions. In stant baryon-number density, because structures move,
this section, the number of solitons in a system is varied see Fig. 9. The visible structure in a pseudo-time-
tocreatedifferentdensities,whilethelatticesizeandthe averagedbaryon-densityplotisdependentonthenumber
lattice spacing remainfixed. Although the results shown of iterations used to create it, in this case 2000 pseudo-
have the same lattice parameters, the same conclusions time steps.
canbedrawninsimulationswithdifferentlatticeparam- In the phase-coexistence region, one observes series of
eters, as has been checked for a number of examples. multi-solitons that do not change with pseudo-time, i.e.
An example of the baryon-number density distribu- a solid structure surroundedby large regions of vacuum.
tion for a solid is given in Fig. 8. On the left we show These can either be interpreted as a crystal of nuclei, or
a pseudo-time-averagedbaryon-numberdensity plot and asapercolatingnetwork. AnexampleisshowninFig.10.
on the right a snapshot where thermal fluctuations are The solid, liquidand phase-coexistencestates canalso
evident. The approximate triangular crystal structure be identified by examining correlation graphs. The cor-
10
(a) Solid
C
32 0.04
0.03
0.02
24 0.01
r
2.5 5 7.5 10 12.5 15 17.5
-0.01
16 -0.02
(b) Phase Coexistence
8 C
0.04
0.03
0.02
8 16 24 32
0.01
r
2.5 5 7.5 10 12.5 15 17.5
FIG. 10: Pseudo-time averaged baryon-number density plot -0.01
in the phase-coexistence region. The shades from dark to -0.02
light represents low to high baryon-number densities. As in
Figs. 8,9 ρ = 0.04, but the temperature is intermediate be-
(c) Liquid
tween thetwo. C
0.04
0.03
relation graphs are created by plotting the correlation 0.02
between the baryon number at the lattice points, 0.01
r
2.5 5 7.5 10 12.5 15 17.5
allploaitnttisce alploliantttsice -0.01
C(r)= hB2B2i−hB2ihB2i -0.02
a b a b
a=1 b=1
X X (cid:0) (cid:1) FIG.11: Thecorrelationfunctionsforasolid,astateinphase
×δ(r−|r −r |), (21)
a b
coexistence,andaliquid,respectively. (a):Forthesolid,there
arepeaksadistanceofseveralSkyrme-lengthunitsawayfrom
against r, which is the separation between lattice points
theoriginbecausethecrystalisrepetitive. (b):Forphaseco-
a and b. The correlation function C at r is calculated
existence, there are also peaks a distance of several Skyrme-
by taking the average of the correlation with respect to
lengthunitsawayfromtheoriginbecausethebaby-Skyrmions
all a and b separated by a distance r. To keep track of
areboundtogether,butthepeaksareweakerthanforasimu-
the value Ba2Bb2 for each separation distance, we bin the lationinthesolidphase,andlackthefinestructure. (c):Fora
values, and therefore the distances are rounded to the liquid,thecorrelation functiondecreasesrapidlytozeroafter
nearest half-lattice spacing, i.e. to the nearest h/2. a few Skyrme-length units. The non-zero component exists
Thecorrelationgraphscanbeusedtoidentifythestate because of the correlation of Skyrmions with themselves. (r
ofasystem. Forsolidsandphasecoexistence(whichis a is in Skyrme-lengthunitsand C is dimensionless.)
solidsurroundedby vacuum), correlationfunctions show
peaks a distance of several Skyrme-length units away
from the origin, see Fig. 11. For liquids, however, the through either method. The simulations whose baryon
correlation function levels to zero quickly. The non-zero densitiesareplottedonthelefthaveallbeenidentifiedas
component is the correlation of Skyrmions with them- solidsandthoseplottedontherighthavebeenidentified
selves, and is therefore non-zero for approximately the as liquids. The pseudo-time average for these plots was
length of the radius of a baby-Skyrmion. Some struc- set to 2000 iterations, and therefore it can be seen that
ture is visible to a distance of 10 Skyrme-length units the baryons move around more for simulations at higher
becauseinteractionswithneighboringSkyrmionsarestill temperature. When using more iterations, the baryon
present. Itisbecauseofthevisiblelong-rangeinteraction density approaches the average everywhere, but because
that we interpret the state to be a liquid rather than a there is less structure, it defeats the purpose of these il-
gaseous state. The advantage of using correlation func- lustrations. It is not possible to identify the structure at
tions rather than pseudo-time-averaged baryon-density T = 0.2 directly from the baryon-density plot, but the
plots to identify liquids is that they are not dependent correlation function clearly shows that it is a solid. The
on the number of iterations used to create them. apparent motion of the solitons are lattice vibrations.
The transition from solid to liquid as the temperature Fig. 13 shows the final phase portrait for the baby-
is raised at a constant density is shown in Fig. 12. Both Skyrme model which has been created by examining the
the pseudo-time-averaged baryon-density plots and the pseudo-time-averaged baryon-density plots and correla-
corresponding correlation functions are shown to high- tionfunctionsatvarioustemperaturesanddensities. We
lightthefactthatthestateofthesystemcanbeidentified identify solid, liquid, and phase-coexistence regions, by
