A point particle model of lightly bound skyrmions

A simple model of the dynamics of lightly bound skyrmions is developed in which skyrmions are replaced by point particles, each carrying an internal orientation. The model accounts well for the static energy minimizers of baryon number $1\leq B\leq 8$ obtained by numerical simulation of the full field theory. For $9\leq B\leq 23$, a large number of static solutions of the point particle model are found, all closely resembling size $B$ subsets of a face centred cubic lattice, with the particle orientations dictated by a simple colouring rule. Rigid body quantization of these solutions is performed, and the spin and isospin of the corresponding ground states extracted. As part of the quantization scheme, an algorithm to compute the symmetry group of an oriented point cloud, and to determine its corresponding Finkelstein-Rubinstein constraints, is devised.


Introduction
The Skyrme model is an effective theory of nuclear physics in which nucleons emerge as topological solitons in a field whose small amplitude travelling waves represent pions. It thus provides a unified treatment of both nucleons and the mesons which, in the Yukawa picture, are responsible for the strong nuclear forces between them. While the Skyrme model has been superceded as a fundamental model of strong interactions by QCD, interest in the model revived once it was recognized to be a possible low energy reduction of QCD in the limit of large N c (number of colours) [17,16], and much work has been conducted to extract phenomenological predictions about nuclei from standard versions of the model [11,12,10,3]. Many of these predictions are in good qualitative agreement with experiment, and recent improvements in skyrmion quantization schemes offer hope of significant further improvement to come [6,7].
One area in which standard versions of the model perform poorly, however, is that of nuclear binding energies: typically, classical skyrmions are much more tightly bound than the nuclei they are meant to represent (by a factor of 15 or so). In recent years, no fewer than three variants of the model have been proposed which seek to remedy this problem. In each case, the model is, by design, a small perturbation of a Skyrme model in which the binding energies vanish exactly. Perhaps the most radical proposal, due to Sutcliffe and motivated by holography, couples the Skyrme field to an infinite tower of vector mesons [15]. Small but nonvanishing binding energies are (conjecturally) introduced by truncating this infinite tower at some high but finite level. This proposal, while elegant, has so far not been amenable to detailed analysis. A second proposal, due to Adam, Sanchez-Guillen and Wereszczynzki, starts with a model which is invariant under volume preserving diffeomorphisms of space, then perturbs it by mixing with a small fraction of the conventional Skyrme energy [1]. Skyrmions in this model have the attractive feature of being somewhat akin to liquid drops. However, the large (in fact, infinite dimensional) symmetry group of the unperturbed model is extremely problematic for numerical simulations, and the shapes and symmetries of classical skyrmions, even for rather low baryon number (B ≥ 3) are, so far, not known in this model in the regime of realistically small binding energy [5].
In this paper we will study the third (and arguably least radical) proposal, originally due to one of us [8]. This amounts to making a nonstandard choice of potential term in the standard Skyrme lagrangian and, more importantly, radically shifting the weighting of the derivative terms from the quadratic to the quartic. The resulting model is still amenable to numerical simulation, but its classical solutions are quite different from conventional skyrmions: the lowest energy Skyrme field of baryon number B now resembles a loosely bound collection of B spherically symmetric unit skyrmions, rather than a tightly bound object in which the skyrmions have merged and lost their individual identities. In the terminology of [14], which studied a (2 + 1) dimensional analogue of the model, skyrmions in this lightly bound Skyrme model prefer to hold themselves aloof from one another. Numerical analysis reveals [5] that they also prefer to arrange themselves on the vertices of a face centred cubic spatial lattice, with internal orientations dictated by their lattice position. This suggests that, unlike conventional skyrmions, lightly bound skyrmions can be modelled as point particles, each carrying an internal orientation, interacting with one another through some pairwise interaction potential whose minimum encourages them to sit at a fixed separation with their internal orientations correlated. The aim of this paper is to derive such a simple point particle model, compare its predictions with numerical simulations of the full field theory, and use it to extract, via rigid body quantization, phenomenological predictions about nuclei with baryon number 2 ≤ B ≤ 23. A similar programme (minus quantization) for the (2 + 1) dimensional analogue model was completed in [14].
As we shall see, the point particle model accounts almost flawlessly for static skyrmions with 1 ≤ B ≤ 8, where comparison with simulations of the full field theory is available. For B ≥ 9, it predicts a rapid proliferation of nearly degenerate skyrmions as B grows, all rather close to size B subsets of the face centred cubic lattice. In comparison with conventional skyrmions, these typically have rather little symmetry, and anisotropic mass distribution. Determining the symmetries of these configurations is an interesting and important task, nonetheless, as they determine the Finkelstein-Rubinstein constraints on quantization. Usually, symmetries of skyrmions are determined by ad hoc means: one looks at suitable pictures of the skyrmion, predicts a symmetry by eye, then checks it by operating on the numerical data. By contrast, we will develop an algorithm which automatically computes the symmetry group of any point particle configuration. This allows us to completely automate the rigid body quantization scheme. The result is, as a phenomenological model of nuclei, moderately successful: rigid body ground states plausibly account for the lightest nucleus of baryon number B for 12 of the 23 values considered. Presumably this can be improved by replacing rigid body quantization by something more sophisticated.
The rest of the paper is structured as follows. In section 2 we review the lightly bound Skyrme model, focussing on its spin-isospin symmetry and associated inertia tensors. In section 3 we introduce the point particle model, then in section 4 we describe a numerical scheme to find its energy minimizers, and present the results of this scheme. In section 5 we formulate the rigid body quantization of our classical energy minimizers, focussing particularly on the Finkelstein-Rubinstein constraints. Some concluding remarks and possible future directions of development are presented in section 6.

The lightly bound Skyrme model
The field theory of interest is defined as follows. There is a single Skyrme field U : R 3,1 → SU(2), required to satisfy the boundary condition U (t, x) = 1 as |x| → ∞ for all t. Such a field, if smooth, has at each t, a well-defined integer valued topological charge the topological degree of the map U (t, ·) : R 3 ∪ {∞} → SU(2) ∼ = S 3 . Since the field is smooth, B(t) is smooth and integer valued, hence automatically conserved. Physically it is interpreted as the baryon number of the field U . The right invariant current associated with U is R µ = (∂ µ U )U † , in terms of which the lagrangian density is Here F π is the pion decay constant, m π the pion mass, and e > 0, 0 ≤ α < 1 are dimensionless parameters. In [5] the following values were chosen for these parameters so that classical binding energies in the model are comparable with experimentally-measured nuclear binding energies: 1 F π = 36.1 MeV, m π = 303 MeV, e = 3.76, α = 0.95.
There is certainly room for improvement in this calibration: for example, obtaining the correct pion mass was not a priority in [5], and we expect that a more thorough analysis could result in a parameter set for which m π is closer to its experimental value of 137MeV. However, the aim in the present paper is not to fine-tune the parameters, but rather to study qualitative properties of static solutions, which we expect to be insensitive to details of the calibration. 1 The value for F π recorded here corrects a typographical error in [5] It will be convenient to use F π /4e √ 1 − α as a unit of energy and 2 √ 1 − α/F π e as a unit of length; in these units the lagrangian takes the form and m := (2m π √ 1 − α/F π e). In the parameter set given above, m = 1.00. Note that when α = 0, L is the lagrangian of the conventional Skyrme model with pion mass, while for α = 1 this is a completely unbound model [8]: there is a topological energy bound of the form V ≥ const|B|, but this is attained only when |B| ≤ 1.
The first approximation to a nucleus containing B nucleons is a static Skyrme field U : R 3 → SU(2) of degree B which minimizes the potential energy V . Thus it is important to identify static classical energy minimizers. These are referred to as skyrmions. A better approximation to a nucleus is obtained by allowing solitons to carry spin and isospin. The lagrangian is invariant under a left action of the group G := SU(2) I × SU(2) J , defined by where we have identified physical space R 3 with the Lie algebra su(2) via x ∼ = ix j σ j , σ 1 , σ 2 , σ 3 being the Pauli matrices, to define the action of h on x. Equivalently, The conserved quantities associated with these symmetries are isospin and spin. We refer to transformations g ∈ SU(2) I as isorotations, in analogy with rotations h ∈ SU(2) J .
Every ω ∈ g := su(2) I ⊕ su(2) J defines a one-parameter subgroup {exp(tω) : t ∈ R} of G isomorphic to S 1 , whose action on a static skyrmion U generates a rigidly isorotating and rotating skyrmion, U ω = exp(tω) · U , of constant kinetic energy T [U ω ]. The mapping ω → T [U ω ] is a quadratic form on g, and hence defines a unique symmetric bilinear form Λ : g × g → R called the inertia tensor of the skyrmion U . By its definition, Λ vanishes on the subspace of g tangent to the isotropy group G U of U (that is, the subgroup G U := {(g, h) ∈ G : (g, h) · U = U } < G which leaves U unchanged). If G U is discrete, as is the case for all the skyrmions studied in this paper except when B = 1, then Λ is a positive bilinear form, and thus defines a left invariant Riemannian metric on G. In order to identify spin and isospin quantum numbers of skyrmions corresponding to those of nuclei, isorotations and rotations needed to be treated quantum mechanically rather classically. The inertia tensor plays an important role in the simplest quantization scheme, known as rigid body quantization, which will be reviewed in section 5, and amounts to quantizing geodesic motion on (G, Λ), subject to certain symmetry constraints required to give skyrmions fermionic exchange statistics. Clearly, by choosing a basis for su(2), we obtain a basis for g which can be used to represent Λ as a real symmetric 6 × 6 matrix. We shall consistently represent inertia tensors in this way, having chosen the basis [− i 2 σ 1 , − i 2 σ 2 , − i 2 σ 3 ] for su(2).

The point particle model
Extensive numerical simulations reported in [5] showed that skyrmions in the lightly bound Skyrme model with B > 0 invariably resemble collections of B particles. Encouraged by this observation, we have developed a point particle model in which a Skyrme field U with baryon number B is replaced by B oriented point particles in R 3 .
To explain how the model is derived, we begin by recalling the structure of the simplest skyrmion, which has B = 1, and is of "hedgehog" form with f (r) a real function satisfying f (0) = π, f (r) → 0 as r → ∞, and r = |x|. The profile function is determined by solving (numerically) the Euler-Lagrange equation for V restricted to fields of hedgehog form, a certain nonlinear second order ODE for f . One finds that U H has total energy M H := V [U H ] ≈ 87.49, and its energy density is monotonically decreasing with r and concentrated around the origin. The 1-skyrmion has a high degree of symmetry: In other words, G U H is the diagonal subgroup of SU(2) I × SU(2) J . This basic skyrmion can be moved and rotated using symmetries of the model. A 1skyrmion with position x 0 ∈ R 3 and orientation q 0 ∈ SU(2) is given by The energy-minimizers with 2 ≤ B ≤ 8 resemble superpositions of fields of this type [5]. More precisely, their energy densities are concentrated at B well-separated points x 1 , . . . , x B , and near each such point x a the field U is approximately of the above form for some q a . These positions and orientations are the basic degrees of freedom in our point particle model, and will be allowed to depend on time t. The lagrangian for this point particle model takes the form where |q| 2 := 1 2 Tr(qq † ) and is an interaction potential.
The terms involving time derivatives of x a and q a represent the kinetic energy of a moving skyrmion. Their coefficients could be deduced from the Skyrme model. It is known that the 1-skyrmion has inertia tensor From this it follows that the kinetic energy of a rigidly rotating skyrmion should take the form

Symmetries of the interaction potential
The point particle model inherits an action of G = SU(2) I × SU(2) J from the Skyrme model. The action of (g, h) ∈ G on the field U (x; Therefore the action of (g, h) on a point particle configuration is (x a , q a ) → (R(h)x a , gq a h −1 ), a = 1, . . . , B.
The point particle lagrangian should be invariant under these transformations, and under translations x a → x a + c for c ∈ R 3 . It should be invariant under changes of the signs of any of the q a , because U (x; x 0 , −q 0 ) = U (x; x 0 , q 0 ). It should also be invariant under permutations of the particles, because configurations of particles that are the same up to a re-ordering describe the same Skyrme field. Finally, the Skyrme model is invariant under the inversion which is equivalent, for a field of the form (3.2), to (x 0 , q 0 ) → (−x, q 0 ). Hence, our point particle lagrangian should be invariant under (x a , q a ) → (−x a , q a ). (3.5) The kinetic terms in (3.3) obviously have these symmetries. Demanding that the potential (3.4) is also invariant imposes constraints on the function V int (x 1 , q 1 , x 2 , q 2 ) which we now describe.
Translation symmetry implies that V int (x 1 , q 1 , x 2 , q 2 ) depends on the positions of the skyrmions only through their relative position X := x 1 −x 2 . Isorotation symmetry implies that it depends on q 1 , q 2 only through the isorotation-invariant combination Q = q −1 1 q 2 . Thus while rotational symmetry demands that A permutation (x 1 , q 1 , x 2 , q 2 ) → (x 2 , q 2 , x 1 , q 1 ) changes the sign of X and inverts Q, so permutation invariance implies that Finally, symmetry under inversion (3.5), implies To proceed further, it is helpful to think of V red as a one-parameter family of real functions V ρ on S 2 × SU(2), parametrized by ρ := |X| ∈ (0, ∞). We may expand each such function in a convenient basis for L 2 (S 2 × SU(2)), for example, the basis of eigenfunctions of the Laplacian. A natural truncation to finite dimensions is obtained by keeping only eigenfunctions up to a fixed finite eigenvalue. The effect of this truncation is to exclude from V red terms with fast orientation dependence. This motivates the following definition: for each λ in the spectrum of ∆ S 2 ×S 3 , let E λ denote the corresponding eigenspace, and for any µ ≥ 0, (3.10) Let C ∞ µ denote the space of smooth functions on V : Proof. Recall that the eigenvalues of the Laplacian on S n are λ (n) d = d(d+n−1), d = 0, 1, 2, . . ., and the corresponding eigenspaces, E (n) d , are spanned by (the restrictions to S n ⊂ R n+1 of) harmonic homogeneous polynomials in R n+1 of degree d [2]. It follows that the eigenvalues of ∆ S 2 ×S 3 are λ d . By (3.6), (3.9), V is invariant under both X → −X and Q → −Q, so we may restrict d and d to only even values (homogeneous polynomials of odd degree are parity odd). Further, since V ∈ C ∞ µ with µ < 20, each restriction V ρ lies in 2 ). (3.12) Now SU(2) acts on both E carries the irreducible spin representation of SU (2). In particular, E (3.13) Now the tensor product R 2d+1 ⊗ R 2 +1 contains no trivial subrepresentation if d = , and exactly one if d = .
Hence, of the summands in (3.12), E 0 , E 8 and E 14 each contain a onedimensional subspace on which SU(2) acts trivially (while E 6 does not) and, by (3.7), V ρ lies in the three-dimensional space spanned by these. Clearly E triv 0 = E 0 which is spanned by the constant function (X, Q) → 1. Consider the functions (3.14) These are manifestly SU(2) invariant and extend to homogeneous polynomials on R 3 × R 4 of bidegree (0, 2) and (2, 2) respectively. Furthermore, one may readily check that these polynomials are harmonic (separately with respect to X and Q). Hence, they span E triv 8 and E triv 14 respectively. Noting that |X| 2 ≡ 1 on S 2 , the claim follows. From now on, we assume that V red lies in the truncated function space C ∞ 14 , so that it has the structure prescribed by Proposition 1.
Recall that, in the standard Skyrme model, the interaction potential for well separated skyrmion pairs can be modelled using the dipole formalism [13]: far from its centre, a unit skyrmion looks like the field induced in the linearization of the Skyrme model about the vacuum, U = 1, by an orthogonal triplet of scalar dipoles placed at the skyrmion's centre. The interaction potential for a skyrmion pair with relative position X and orientation Q can then be approximated by the interaction energy of a pair of triplets of dipoles held at relative displacement X and orientation Q, interacting via the linear theory. This approximation introduces another useful constant associated with the unit skyrmion, namely the strength of the (necessarily equal) dipoles. In practice this is determined numerically by reading off a coefficient C in the large r asymptotics of the skyrmion profile function. This formalism is readily adapted to the lightly bound Skyrme model, producing an interaction potential of the form (3.11) with V 0 (r) = 0 The dipole strength (for α = 0.95 and m = 1) is found numerically to be C ≈ 14.58. These formulae reproduce the usual prediction of attractive and repulsive channels for well-separated skyrmions. That is, V red is maximally attractive (increases fastest with |X|) if the orientations of the skyrmions differ by a rotation by π about any direction orthogonal to X, is maximally repulsive if the orientations differ by a rotation by π about X, and is nonmaximally repulsive if their orientations are equal. We refer to these three situations as the attractive, repulsive and product channels respectively. The existence of these three channels allows us to fix the functions V 0 , V 1 , V 2 numerically by conducting scattering simulations of skyrmion pairs in the full field theory, in similar fashion to Salmi and Sutcliffe's work on the (2 + 1) dimensional model [14]. We begin with a Skyrme field of the form where s > 0 is large and U H is a unit hedgehog skyrmion defined (numerically) in a ball of radius less than s/2 (so U H (x) = 1 for all |x| ≥ s/2, and the product above commutes). Such a field represents a pair of skyrmions located at x = (±s/2, 0, 0), that is, with separation s, in the attractive channel. Here, and henceforth, we define the skyrmion positions of a Skyrme field U : R 3 → SU(2) to be those points where U = −1. We now allow U to evolve with time according to the dynamics defined by the lagrangian (2.2), using the fourth order spatial discretization employed by the energy minimization scheme of [5], and a fourth order Runge-Kutta scheme with fixed time step for the time evolution. This numerical scheme conserved total energy E = T + V to extremely high accuracy, for all the dynamical processes presented here. As the dipole model predicts, the skyrmions with these initial data slowly move towards one another, attain a minimum separation, then recede again. By recording their separation s(t) and potential energy V (t) at each time step, we recover a numerical approximation to the attractive channel interaction potential which, according to (3.11 We then repeat the process with intial data which are in the repulsive and product channels respectively. To make the skyrmions approach one another and interact, we now Galilean boost them towards one another at low speed (v=0.1). Note that the reflexion symmetries of the initial data trap these fields in their respective channels for all time. From these numerical solutions we obtain numerical approximations to the repulsive and product channel interaction potentials, which are related It is clear that V a , V r , V p uniquely determine V 0 , V 1 , V 2 and hence, within the ansatz (3.11), V int . Figure 1: Interaction energies of skyrmions pairs with separation s in the attractive (blue), product (red) and repulsive (green) channels. In each case the thick curve represents numerical data extracted from a scattering process, the thin curve is a fit to this, and the dashed curve is the interaction energy predicted by the dipole model.
Graphs of V a , V r , V p , determined numerically as described above, are presented in figure  1. These curves also show the potentials predicted by the dipole model (with dipole strength C = 14.58). Clearly, the dipole formulae (3.15) do not provide an accurate quantitative picture of skyrmion interactions in the lightly bound model at any separation where the interactions are not negligible. This is, perhaps, not surprising, since the dipole formalism replaces the full field theory by terms originating only in the quadratic and pion mass potential terms of the lagrangian, and these are precisely the terms which are given very low weighting, 1 − α, in the lightly bound regime. The qualitative predictions of the dipole picture are reliable however: the interaction potentials appear to decay exponentially fast, and the three channels identified have the behaviour predicted (attractive, repulsive, more weakly repulsive). For later use, it is convenient to have explicit functions which approximate the numerical data for V a , V r , V p . For our purposes, it is important that these functions decay exponentially with s and accurately fit the numerical data for s ≥ s 0 , where s 0 is somewhat smaller then the equilibrium separation defined by V a (that is, the separation at which V a is minimal). The behaviour for s < s 0 is not so important, provided the formulae introduce a repulsive core interaction, and is, in any case, inaccessible to our numerical scheme (since close approach of lightly bound skyrmions is forbidden in low energy scattering processes). Figure 1 also depicts the following fit functions (3.23) Of these, the most elaborate is V a , a Padé approximant on [0, 7.096] spliced to an exponentially decaying tail, the splice being chosen so that V a is continuously differentiable. Unlike V r and V p , V a is well defined at s = 0, where it is chosen to equal the static energy of the axially symmetric B = 2 solution (a saddle point of the Skyrme energy), obtained numerically by a different scheme, a choice made mainly for aesthetic reasons. From now on, we choose V red to be the function defined by (3.11), where and V a , V p , V r are the functions defined in (3.23). It is straightforward to show that this function V red is bounded below as, on physical grounds, it should be.

The FCC lattice
We have seen that the interaction potential V int prefers particles to be in the attractive channel, i.e. such that their relative orientation corresponds to a rotation about an axis perpendicular to their line of separation through angle π. It is therefore desirable to find a way to pack them together such that all neighbouring pairs of particles are in the attractive channel. The face-centred-cubic (FCC) lattice provides a solution to this problem. The face-centred cubic lattice may be defined to be with λ > 0 defining a lattice scale. The underlying cubic lattice is given by points (n 1 λ, n 2 λ, n 3 λ) for which n 1 , n 2 , n 3 are all even. Those points for which some of the coordinates n i are odd lie on faces of the underlying cubic cells.
We assign orientations to these points as follows: those points on the vertices have orientation 1 ∈ SU(2), those on faces perpendicular to the x-axis have orientation i, those on faces perpendicular to the y-axis have orientation j, and those perpendicular to the z-axis have orientation k. Here we have implicitly identified elements q ∈ SU(2) with unit quaternions q ∈ H, such that i = −iσ 1 , j = −iσ 2 , k = −iσ 3 and 1 is the identity matrix. Put differently, the orientation q of a particle at lattice site (n 1 , n 2 , n 3 )λ is such that The reader may verify that any pair of nearest neighbours, separated by a distance λ √ 2, is in the attractive channel.
One might expect that minimizers of the potential energy derived from (3.3) resemble subsets of the FCC lattice. This was certainly true of all global minima of the Skyrme energy identified in [5], and all but one of the local minima.

Inertia tensors
The point particle model (3.3) makes simple predictions for the inertia tensors of lightly bound skyrmions. These are obtained by calculating the kinetic energy of a rotating and isorotating oriented point cloud.
Let {(x a , q a )} be a minimizer of the potential energy derived from (3.3). Choose any pair of angular velocities (ω I , ω J ) ∈ su(2) ⊕ su (2). It is useful to identify each ω ∈ su(2) with a vector ω ∈ R 3 by choosing − i 2 σ j , j = 1, 2, 3, as a basis for su(2) (so ω = − i 2 ω · σ). Consider the following configuration, which is isorotating and rotating at constant angular velocity ω = (ω I , ω J ): We findẋ Therefore the kinetic energy is where the inertia tensor is The point particle model predicts that this is a good approximation to the inertia tensor of a lightly bound degree B skyrmion. We will test this prediction in the next section.
4 Energy minimizers in the point particle model 4

.1 Light nuclei
Having introduced the point particle model for lightly bound skyrmions, in this section we present our results for energy-minimizing configurations of point particles. We begin by discussing our results for eight particles or fewer, where comparison can be made with energy minima in the lightly bound Skyrme model found in [5].
We have developed an iterative zero-temperature annealing algorithm to minimize the energy of a configuration of particles. We applied this algorithm both to randomly-chosen initial ensembles of particles and to initial ensembles that are subsets of the FCC lattice. We ran a large number of simulations for each value of B, typically obtaining several local energy minima, and record here only the lowest local minimum and up to two closest competitors. Energies of these local minima with 2 ≤ B ≤ 8 are presented in table 1. The particle ensembles themselves are depicted in figure 2. The corresponding binding energies in the lightly bound Skyrme model are also recorded in the table. These are defined to be the energy of the B-skyrmion minus B times the energy of the 1-skyrmion.  Table 1: Energies and numbers of bonds of the lowest-energy local minima in the point particle model, and energies of the corresponding lightly bound skyrmions (taken from [5], except those marked * , which result from new simulations conducted with the same numerical scheme).
Our results are almost entirely consistent with the results obtained for the lightly bound Skyrme model in [5]. For 1 ≤ B ≤ 5 we obtained the same global minima as in the lightly bound Skyrme model. For B = 6, 7, 8 multiple local minima were previously obtained in the lightly bound Skyrme model. All of these occured as local minima in the point particle model. For B = 7, 8 the ordering of energies in the point particle also agreed with the ordering of energies in the lightly bound Skyrme model. The only failure of the point particle model is for 6 particles: here the energies of the two lowest-energy local minima appear in the wrong Figure 2: Local energies minimizers in the point particle model. Each ball is centred on a point skyrmion position x a and its colour represents the internal orientation q a . Each picture also depicts the FCC lattice configuration of size B to which the minimizer best fits. Thick grey line segments indicate interskyrmion bonds no more than 10% longer than the FCC bond length, while thin magenta line segments show nearest neighbour bonds in the best fit lattice configuration. In most cases the fit is so good that the thin bonds are not visible. They show quite clearly on 5(a), 5(b), 6(b), 6(c), 7(b) and 8(c) however. order.
In addition to reproducing previously-known minimizers from the lightly bound Skyrme model our point particle model also predicted some new local minima. Most interestingly, the global energy-minimizer in the point particle model for B = 7, labelled 7a in figure 2, did not correspond to any solution of the lightly bound Skyrme model found in [5]. Based on this discovery, we constructed an approximate Skyrme field with a similar shape to the point particle energy-minimizer, and minimized its energy using the same numerical scheme that was used in [5]. After relaxation this Skyrme field had a lower energy than any of the configurations discovered in [5], as predicted by the point particle model. Thus we have a new candidate global energy minimizer at charge seven. Similarly, new simulations find local energy minimizers in the lightly bound model of similar shape to 5b, 6c and 8c, and these have energies ordered exactly as the point particle model predicts (so E 8c > E 8b > E 8a , for example).
In every case, the minimizers found look, to the naked eye, like subsets of the FCC lattice. 2 It is an interesting problem to measure this property quantitatively. Given an oriented point cloud (X, Q) = (x 1 , . . . , x B , q 1 , . . . , q B ), we wish to identify the FCC subset of size B which best approximates it. To do this, we consider the orbit of (X, Q) under the group S of similitudes of R 3 , For each s ∈ S, we define d 2 to be the squared distance from s · X to the FCC lattice, i.e.
Now, given a neighbouring triple of particles in X (a particle x, its nearest neighbour x and next-nearest neighbour x ), we construct a similitude s 0 which maps x to 0, x to (1, 1, 0) and x to the plane spanned by (1, 1, 0), (0, 1, 1). We then solve the gradient flow equation of d 2 : S → R, with s(0) = s 0 , to find a local minimum of d 2 close to s 0 . Repeating over all neighbouring triples, we keep the lowest local minimum s min of d 2 found (note that d 2 never has a global minimum since s · X can be made arbitrarily close to (0, 0, 0) by taking λ sufficiently large). In this way we identify the closest FCC subset to X and its root mean square distance from X, namely d RM S = d(s min ) 2 /B. Having found s min · X, the FCC colouring rule predicts the internal orientations (q 1 , . . . , q B ) the particles should have. These should be compared with (q 1 h −1 min , . . . , q B h −1 min ), bearing in mind that orientations are defined only up to sign, and that the system is isospin invariant. Thus we minimize over g ∈ SU(2) I , again by gradient flow. This gives us a measure of the root mean squared distance of the internal orientations of the configuration (X, Q) from those imposed by the colouring rule applied to its closest FCC approximant, namely d iso RM S = d 2 iso (g min )/B. It also allows us to "coarse grain" the internal orientations, that is, map each q a to the element of {±1, ±i, ±j, ±k} to which gq a h −1 min is closest. We used this method to determine the particle colours and FCC bonds in figure 2. We will present graphs of d RM S and d iso RM S in the next section.
In addition to comparing energies we have also compared inertia tensors in the point particle and lightly bound Skyrme models. Under isorotations and rotations (g, h) ∈ SU(2) I × SU(2) J inertia tensors transform as In comparing the inertia tensors of a charge B skyrmion, obtained by solving the field theory, and a charge B point particle energy minimizer, we must account for the fact that the orientations of these two objects are completely unrelated. We do this by introducing a standard form for inertia tensors which fixes these symmetries. We say that an inertia tensor Λ is in • if λ 1 = λ 2 then µ 1 , µ 2 , µ 3 are either all non-negative or all non-positive and ν 1 , ν 2 , ν 3 are either all non-negative or all non-positive; and • if λ 1 = λ 2 then ν 3 = 0, |µ 1 | > |µ 2 |, and µ 1 , µ 2 , µ 3 , ν 1 , ν 2 are either all non-negative or all non-positive.
Any inertia tensor has a matrix of standard form in its SU(2) I × SU(2) J orbit, and, in generic cases this matrix is unique. Note that we have chosen not to define standard form as being a form in which both the upper-left and lower-right blocks of Λ are diagonal, even though such a form is arguably simpler than the one described above. The reason is that the upper-left block of any inertia tensor obtained in the point particle model is proportional to the identity, so diagonalising the upper left block does not fix the isorotation symmetry. We shall measure the distance between inertia tensors by the distance between their standard forms, using the usual Euclidean norm on the space of real matrices, that is Λ 2 := Tr(Λ T Λ).

Heavier nuclei
When searching for local energy minima with large numbers of particles, one faces the problem that the number of connected subsets of the FCC lattice grows rapidly with the number of name 1 2 3 4 5 6b 7a 8a error 1.96% 5.96% 1.61% 1.15% 4.65% 9.91% 1.97% 3.20% Table 2: Percentage error in inertia tensors calculated in the point particle model, as compared with the lightly bound Skyrme model particles, and hence the number of candidate local minima of the energy grows rapidly. We addressed this problem by seeking only local minima corresponding to FCC lattice subsets with a large number of bonds. More precisely, we used as initial conditions in our relaxation algorithm only lattice subsets whose number of bonds is at most two less than the maximum possible for the given number of particles. In the end we found that global energy minima always had at most one less than the maximum number of bonds, so this restriction seems reasonable.
Even with this simplification, the number of initial conditions to consider is large and it is difficult to be sure that enough simulations have been run to find the global energy minimizer.
To solve this problem we separated our minimization algorithm into two stages: in the first stage, a list of distinct lattice subsets is generated, and in the second stage these subsets are relaxed as before. Our method for telling whether two lattice subsets are distinct is to compute their energy: if two lattice subsets have the same energy to high precision we assume that they are identical and discard one. In doing so we run the risk that a lattice subset whose energy happens to coincide with another is wrongfully discarded. For example, the initial FCC subsets used to generate solutions 6b and 6c have exactly the same spectrum of bonds, and hence exactly the same energy. Only after relaxation away from the FCC lattice do their energies separate. To mitigate against this danger we ran extensive simulations up to 16 particles starting from randomly chosen lattice subsets satisfying the bond number constraint; in all cases we obtained the same minimum energy as when we started with a list of distinct lattice subsets.
One distinct advantage of our method is that it makes it easy to identify not just the global energy-minimizer but also local energy minima. Another is that it allows one to tell with reasonable confidence when sufficiently many lattice subsets have been sampled. Throughout the procedure the number of occurences of each subset is recorded, and when all of these numbers are above a fixed minimum one may assume that all distinct lattice subsets have been found and terminate the algorithm.
In order to generate lattice subsets to use as initial conditions we developed a crystalgrowing algorithm. Again, this algorithm proceeds in two stages. In the first stage a connected subset of the FCC lattice is generated iteratively. This scheme starts with a lattice subset consisting of a single point. At each step of the iteration a member of the lattice subset is chosen at random and one of its twelve nearest neighbours is chosen at random. If the neighbour is not already a member of the lattice subset, it is appended, otherwise it is discarded. This continues until the lattice subset has the required number of particles. In the second stage of the algorithm the subset is modified so as to increase the number of bonds while maintaining a fixed number of particles. At each step the algorithm chooses at random one member of the subset and a neighbour of another member. If the neighbour is not already a member of the subset, and replacing the original member with this neighbour increases the number of bonds, the algorithm makes this replacement; otherwise, nothing happens. This continues for a fixed number of steps. At each step of the second stage the lattice subset is recorded, so running the algorithm once generates a large number of lattice subsets.
The crystal-growing algorithm was run repeatedly and distinct subsets with sufficiently many bonds saved until it was deemed that enough lattice subsets had been sampled, according to the above-defined criteria. The maximum number of bonds and the number of crystals identified satisfying our criteria are recorded in table 3.
The output of our algorithm is recorded in table 4. The total number of local energy minima found was huge; in this table we list all local minima whose energy is within 0.1 of the lowest energy found, together with configurations 5b, 6b and 6c. For the most part, energy minimizers have the maximum number of bonds possible (exceptions in the table are marked by asterisks), and have the most even distribution of particle "colours" (after coarse graining) possible. A notable exception to both these rules is 23a, which has one less bond than maximal, and a rather uneven colour distribution (8,5,5,5) but is, nonetheless, the lowest energy B = 23 configuration found. This minimizer also has unusually high symmetry, as can be seen from figure 3, which also depicts the highly symmetric minimizers 10b and 19b. One should note, however, that the point particle model does not always favour highly symmetric configurations. The B = 13 configuration, let us call it 13sym, obtained by augmenting a single point by all its nearest neighbours, for example, has the maximal number of bonds, but has energy −8.556, which is much higher than the 13a. It also has a very uneven colour distribution: 4,4,4,1. So for B = 13, unlike B = 23, the model prefers to sacrifice symmetry  Table 4: The lowest energy local energy minima in the point particle model. Asterisks in column 2 indicate that the configuration has one or two fewer bonds than the maximum bond number found for that particle number. Column 3 indicates the number of particles of each internal orientation, after coarse-graining. The classical energy is the potential V of eq. (3.4), and the quantum energy is V + ∆, where ∆ is defined in section 5. The experiment column identifies the lightest nucleus for given baryon number B if this nucleus has the spin and isospin predicted.
in favour of uniform colour distribution. These two charge 13 configurations are also depicted in figure 3. Note that all particles in 13a are contained in just two planes of the FCC lattice, a feature it has in common with all global minimizers for 4 ≤ B ≤ 15. The corresponding predictions for nuclear binding energies per nucleon, defined to be −V int /B, are plotted in figure 4. Here, as in [5], energies in table 4 have been converted to MeV by multiplying with 10.72. The curve shows that ensembles of 4 and 16 particles have unusually high binding energies, in agreement with nuclear experiment, although these effects are less pronounced in the point particle model than in experiment. The energy minimizers corresponding to these two peaks are particularly special: they both have tetrahedral symmetry. Note that our binding energy curve lacks the peak seen at baryon number 12 in the nuclear binding energy curve.
In order to analyze the overall shape of the energy minimizers we have calculated for each 10b 13a 13sym 19b 23a Figure 3: Selected energy minimizers for 10 ≤ B ≤ 23. 23a is a large octohedron with a tetrahedron glued to one face. This unusually symmetric configuration is the global minimizer for B = 23, despite having less than maximal bond number and rather uneven colour distribution. By contrast, the exceptionally symmetric configuration 13sym has much higher energy than 13a. Also depicted are the local minimizers 10b and 19b, a large tetrahedron and octohedron respectively.
its second moment matrix M ij , defined by since M is positive definite. In figure 5 we have plotted M 0 against M 1 for the minimizers listed in table 4. Overall there seems to be a downward trend in M 0 / M 1 as M 1 increases, indicating that larger minimizers are closer to being round than small minimizers. However, even for large nuclei the level of anisotropy is substantial.  The determinant det(M 0 ) measures the qualitative nature of the anisotropy. Let us order the eigenvalues of M so that λ 1 ≤ λ 2 ≤ λ 3 . If the point cloud is long and thin then λ 1 ≤ λ 2 < λ < λ 3 , so M 0 has two negative eigenvalues and one positive, whence det(M 0 ) > 0. By contrast, if the point cloud is flat and round then λ 1 < λ < λ 2 ≤ λ 3 , so det(M 0 ) < 0. It is useful to define µ i = λ i − λ, the eigenvalues of M 0 . By extremizing the function det M 0 = µ 1 µ 2 µ 3 on the circle obtained by intersecting the sphere of radius M 0 with the plane µ 1 + µ 2 + µ 3 = 0, one finds that with equality precisely when two of the eigenvalues (of M 0 or, equivalently, M ) coincide. Just as for 2 ≤ B ≤ 8, we can measure the distance of each local minimum found from its closest FCC lattice approximant, both in space and in internal (orientation) space, as defined in (4.2) and (4.3). The results of this analysis are presented in figure 6. With very few exceptions the minimizers match up very closely with FCC subsets, and their internal orientations are very close to the FCC prediction. Note, however, that the optimal lattice scale varies quite significantly with B (rightmost graph), so it is not a good approximation to fix this at the start (to match the FCC bond length to the optimal separation of a single skyrmion pair, for example) and minimize energy only over FCC subsets of that fixed scale. Any attempt to proceed in this way always gets the relative energy ordering of local minima wrong for several values of B. Figure 6: Comparison of energy minimizers of the point particle model with subsets of the FCC lattice for baryon number 3 ≤ B ≤ 23. In each case we shift, scale and rotate the configuration until it matches, as closely as possible, a connected subset of the standard FCC lattice (with integer coordinates). Plot (a) shows the root mean square distance of the particles in each transformed configuration from the FCC lattice, while plot (b) shows the root mean square distance of their internal orientations from those predicted by the FCC colouring rule. Plot (c) shows the scale factor used. In each case, data corresponding to global energy minima (i.e. configurations labelled "a") are circled in red. The dashed line in (c) marks the optimal lattice scale for B = 2.

Rigid body quantization
Nuclei are inherently quantum-mechanical, so to make a direct comparison between skyrmions and nuclei it is necessary to include quantum effects in the Skyrme model. Traditionally this is done semiclassically, treating the classical skyrmions as rigidly rotating and isorotating bodies. The quantized wavefunction is required to satisfy the Finkelstein-Rubinstein constraints. In practice these constraints restrict the spin and isospin quantum numbers of a quantized skyrmion. For example, they guarantee that quantized skyrmions have either half-integer or integer spin and isospin according to whether the baryon number is odd or even. In cases where solitons have symmetry they yield more nontrivial information. In the present section we describe how to apply rigid body quantization and the Finkelstein-Rubinstein constraints in the point particle model. The procedure reduces to a numerical algorithm which we have implemented and applied to the minima presented in the previous section.

Finkelstein-Rubinstein constraints
We begin by recalling the definition of the Finkelstein-Rubinstein constraints (readers not interested in topological details could skip this subsection and continue reading at the start of the next subsection). The classical configuration space of solitons with baryon number B is the space S B of continuous maps U : R 3 → SU(2) of topological degree B, satisfying the boundary condition U (x) → 1 as |x| → ∞. This space is topologically nontrivial: it contains non-contractible loops, so has a nontrivial fundamental group. In fact, π 1 (S B ) = Z 2 for all B ∈ Z. S B has a universal covering spaceS B together with a two-to-one map π S :S B → S B , such that all loops inS B are contractible and a loop in S B is contractible if and only if it can be lifted to a closed loop inS B . The soliton wavefunction is a function Ψ :S B → C. The Finkelstein-Rubinstein constraint on Ψ states that for every pair y and y of distinct points iñ S B such that π S (y) = π S (y ), Ψ(y) = −Ψ(y ). . . x B in R 3 and B elements q 1 , . . . , q B ∈ SU(2). The energy function disfavours vectors x a from being too close, so for practical purposes we may demand that their separations are greater than some fixed minimum δ > 0. Therefore the naive configuration space for the point particle model is The map to Skyrme configuration space is given by a so-called "relativised product ansatz" U (x; x σ(a) , q σ(a) ) . Applying P to a sum of products of matrices in SU(2) yields an SU(2) matrix, as long as the sum of products is everywhere nonvanishing. The argument of P in eq. (5.4) is nonvanishing if the positions x a are sufficiently well-separated, so the right hand side is a map R 3 → SU(2).
This action leavesF invariant, that is,F • Φ (σ,s) =F , soF descends to a continuous map is the true point particle configuration space. SinceC B is simply connected and the action is free,C B is the universal cover of C B and π 1 ( Choose any pair of points x 0 ∈C B , y 0 ∈S B such thatF (x 0 ) = π S (y 0 ). By a standard theorem of topology (see [9, pp 61-2] for example),F has a unique continuous liftF :C B →S B withF (x 0 ) = y 0 . The situation is summarized in the following commutative diagram Note that F is a lift of F . Any wavefunction Ψ :S B → C defines a wavefunction ψ = Ψ •F onC B , which must satisfy some nontrivial constraints derived from the Finkelstein-Rubinstein constraints: Proposition 2. If x, x are two points inC B such that π C (x) = π C (x ), then where (σ, s) ∈ Σ B (Z 2 ) B is the unique group element that maps x to x .
Proof. This result follows almost directly from two important results of Finkelstein and Rubinstein [4]. First, if x ∈C B , and t a ∈ (Z 2 ) B is the transformation that changes the sign of q a (only) and α is a path inC B from x to Φ (Id,t a ) (x) then F • π C • α is non-contractible. Second, if σ ∈ Σ B is a transposition and β is a path inC B from x to Φ (σ,1) (x) then F • π C • α is also non-contractible. Thus the constraint (5.1) implies that and ψ(Φ (σ,1) (x)) = −ψ(x).
Now any element of Σ B (Z 2 ) B can be written as a product of sign flips and transpositions, so the claim follows.

Rigid body quantization
In rigid body quantization, motion is restricted to the rotation-isorotation orbit of a fixed minimum x = (x 1 , . . . , x B , q 1 , . . . , q B ) of the classical energy. Thus the classical configuration space is taken to be G = SU(2) I × SU(2) J with each (g, h) ∈ G identified with The wavefunction ψ : G → C is required to solve a Schrödinger equationĤψ = Eψ, wherê H is (up to a constant factor) the Laplacian operator on G associated with the left invariant metric Λ, the inertia tensor of x. In order to model a nucleus of definite spin and isospin one assumes that ψ is an eigenstate of the total isospin and spin operators with isospin I and spin J. This is consistent with the Schrödinger equation because the hamiltonian commutes with these operators. By the Peter-Weyl theorem, any such ψ is a finite sum of functions of the form where for ∈ 1 2 N, ρ : SU(2) → SU(2 + 1) denotes the spin-representation of SU (2). For each w ∈ V I,J denote by V (w) the subspace of functions with w fixed. Clearly V (w) ∼ = V I,J for all w = 0. Furthermore,Ĥ preserves V (w) , and its action on every V (w =0) is unitarily equivalent. Hence we may, without loss of generality, fix w = 0, and representĤ by a linear operator H I,J on V (w) ∼ = V I,J . To write this operator down explicitly, it is useful to introduce the usual basis for g = su(2) I ⊕ su(2) J , namely with respect to which Λ is a symmetric 6 × 6 real matrix. Denote its entries Λ ab and those of its inverse Λ ab . ThenĤ acts on V (w) ∼ = V I,J as where ρ * I,J : su(2) I ⊕ su(2) J → su((2I + 1)(2J + 1)) is the Lie algebra representation associated to ρ I ⊗ ρ J : σ(a) ).
This in turn implies that I,J = {0} if I is half integer and B is even, or if I is integer and B is odd. Hence, there are no half integer isospin energy eigenstates when B is even, and no integer isospin energy eigenstates when B is odd. Similar comments apply to spin.

Two particles
We now illustrate the quantization procedure for the simple example of two particles. After rotation and centreing, the energy minimizer 2a is where the lattice scale parameter λ takes the value 2.9 to minimize energy and e 1 = (1, 0, 0). This configuration has D 2 dihedral symmetry, and the nontrivial elements of the symmetry group are (i, j), (i, i), and (1, k). The actions of these transformations, and the corresponding signs χ(g, h), are as follows: We will briefly explain how these signs χ(g 0 , h 0 ) have been determined. The first transformation (i, j) permutes the two particles but does not change any signs. It therefore has sgn(σ) = −1, s 1 = s 2 = 1, and consequently χ(i, j) = −1. The second does not permute the particles but does change the sign of exactly one orientation, so has χ(i, i) = −1 . The third permutes the particles and changes one sign, so has χ(1, k) = (−1) 2 = 1. Now we determine the subspaces V is nontrivial, states with higher spin and isospin certainly have higher energy. We need to use the spin 1 representation ρ 1 of SU(2), for which In the case (I, J) = (1, 0) the constraints (5.8) reduce to In the case (I, J) = (0, 1) the constraints (5.8) say that The inertia tensor of the configuration (5.10) is , whose inverse is The representation ρ * 1 of the Lie algebra su(2) is In the case (I, J) = (1, 0) the hamiltonian defined in (5.7) is In the case (I, J) = (0, 1) it is Therefore, the lowest eigenvalue of the restriction of H to V

Automation of the quantization procedure
Most of the quantization procedure described above is linear algebraic and so easily automated, but determining the symmetries of a configuration and the associated signs arising in the Finkelstein-Rubinstein constraints can be tricky. We have developed an algorithm that finds symmetries of configurations of particles, and hence have been able to automate the entire quantization procedure.
Our symmetry-finding algorithm first identifies rotational symmetries of the set of particle positions, and then determines which of these can be lifted to symmetries in SU(2) I × SU(2) S . To find spatial symmetries it first translates the configuration so that its centre of mass is at the origin and rotates it so that its second moment matrix (4.2) is diagonal. It then treats separately three cases corresponding to different degeneracies of the diagonal entries, i.e. eigenvalues, of the second moment matrix.
If the eigenvalues are all different then the symmetry group is a subgroup of the group D 2 ∈ SO(3) consisting of rotations about the three coordinate axes through π. Each element of this group is applied to the particle positions and a distance between the resulting configuration and the original configuration is measured (taking account of permutations). If the distance is sufficiently small then the group element is accepted as a symmetry.
If two eigenvalues are equal and the third distinct then the spatial symmetry group is a subgroup of O(2). Since, by observation, none of the minimizers with three or more particles have continuous symmetry, the symmetry group is assumed to be dihedral or cyclic. To identify cyclic symmetries, rotations through angle θ are applied to the configuration and the distance from the resulting configuration to the original configuration measured as a function of θ. Minima of this function close to zero are interpreted as symmetries. Dihedral symmetries are identified in a similar way.
If all three eigenvalues are equal then the symmetry group is likely to be a discrete subgroup of SO(3) which is neither dihedral nor cyclic. Since icosahedral symmetry is not compatible with the FCC lattice our algorithm works on the assumption that the symmetry group is either the octahedral group O or the tetrahedral group T < O. It computes the fourth moment tensor x a i x a j x a k x a l and finds minima or maxima of the function S 2 → R defined by ijkl n i n j n k n l , n · n = 1. This is a polynomial function on the sphere of degree less than or equal to 4 containing terms of even degree only. It is known that there are only two linearly independent functions of this type with tetrahedral symmetry, namely the constant function and the function n → n 4 1 + n 4 2 + n 4 3 .
The latter furthermore has octahedral symmetry and its maxima are at the points where the coordinate axes intersect the sphere. Therefore, if the configuration has either octahedral or tetrahedral symmetry and the function constructed from the fourth moment tensor is nonconstant then either its maxima are at mutually-orthogonal points on the sphere or its minima are. The algorithm seeks either a pair of orthogonal maxima or a pair of orthogonal minima and rotates these to lie on two of the coordinate axes. It then tests whether each element of the octahedral group is a symmetry of the configuration. Once the configuration's rotational symmetries are known, the algorithm determines whether each lifts to a full rotation-isorotation symmetry of (x 1 , . . . , x B , q 1 , . . . , q B ). Given a rotational symmetry R, we choose h ∈ SU(2) such that R = R(h). Being a rotational symmetry means precisely that R(h)x a = x σ(a) for some permutation σ. Note that spatial rotations change the orientations also, q a → q a h −1 . The spatial symmetry h lifts to a full symmetry if there exists g ∈ SU(2) such that gq a h −1 = ±q σ(a) for all a. If such an isorotation g exists, it is unique up to sign. In fact it must be g = q σ(1) hq −1 1 (or minus this). Hence, our algorithm computes q a := q σ(1) hq −1 1 q a h −1 for each a = 2, . . . , B and tests whether q a = ±q σ(a) (to some numerical tolerance) for all a. If so, (g, h) is accepted as a full symmetry, and its FR factor χ(g, h) is readily computed from the sign of the permutation σ and the signs occuring in q a = ±q a . If not, the spatial symmetry R(h) is discarded.
The output of our algorithm is recorded in table 4. For each classical energy minimizer we record the spin and isospin quantum numbers corresponding to the ground state, and the quantum mechanical energy (defined to be the sum of the classical energy V and the O( 2 ) correction ∆ calculated by our algorithm). The numerical value of is fixed by the calibration proposed in [5]. We have been using energy and length units F π /4g √ 1 − α and 2 √ 1 − α/F π g. In natural units Planck's constant is 1, so in our units it equals In [5] the dimensionless parameter g was determined to be 3.96 by comparing the charge radii of the one-skyrmion and the proton, so the numerical value for 2 is 4 × 3.96 4 ≈ 799.5.  In most cases the configuration with the lowest quantum energy is the same as the configuration with the lowest classical energy. There are two exceptions to this trend. The quantum corrections to the 6-particle minimizers are relatively large and their order is reversed, so that 6c is the lightest and 6a the heaviest. The two configurations 10a and 10b have almost identical classical energies, and after quantization the order of their energies is reversed.
The table also lists spin and isospin quantum numbers of the lightest nucleus for each mass number. In 12 out of 22 cases there is a quantized point particle configuration with the same quantum numbers. Sometimes the configuration with the correct quantum numbers is not that with lowest energy: for example, 5b has the same spin and isospin as 5 He 2 , but its energy exceeds that of 5a.
Including the quantum corrections gives new predictions for nuclear masses and binding energies, which are plotted in figure 7. The mass of a quantized configuration of B particles is BM + V + ∆, where V + ∆ is the quantum energy recorded in 4. The binding energy is the difference between this quantity and B times the quantized mass of 1 particle. The calculation of the quantum mechanical correction to the mass of one particle is a standard calculation similar to those described above; the end result is The binding energy per nucleon is therefore Note that 3 2 /8L ≈ 5.521. Rigid body quantization has the effect of increasing the binding energy of nuclei, so that they are roughly 5% of the 1-skyrmion mass rather than 1%. It is easy to see why: the quantum correction to the 1-skyrmion mass represents about 5% of its total mass, whereas the quantum corrections to the masses of larger nuclei represent a much smaller percentage of their total mass. This means that the quantum corrections to binding energies are also around 5% of the 1-skyrmion mass, and much larger than the classical binding energies.
The fact that binding energies calculated by rigid body quantization are too large does not represent a failure of the lightly bound Skyrme model, but rather illustrates the pitfalls of rigid body quantization itself. A collection of B point particles has 6B −3 degrees of freedom, but in rigid body quantization at most 6 of these are quantized. Only in the case B = 1 are all degrees of freedom quantized, so rigid body quantization systematically underestimates the mass of configurations with a large number of particles. From the point of view of the lightly bound Skyrme model, the degrees of freedom corresponding to moving particles are almost massless, because 1-skyrmions interact only weakly, so arguably these are of comparable importance to the massless degrees of freedom studied in rigid body quantization.

Concluding remarks
We have constructed a simple point particle model of lightly bound skyrmions which almost flawlessly reproduces the results of numerical field theoretic energy minimization for charges 1 to 8. The only exception is charge 6. Here, the point particle model predicts minimizers with shapes, in order of ascending energy, octahedron, bowtie and pyramid-plus-one, whereas full field simulations find that the correct order is bowtie, octohedron, pyramid-plus-one, albeit with the first two of these very close to degenerate. Alongside this minor blemish one should set some unexpected successes: the point particle model predicted previously unknown energy minimizers at charges 5, 7 and 8, all of which corresponded to local energy minimizers of the field theory with correct energy ordering. This includes the (so far) lowest energy skyrmion at charge 7. The point particle model makes a simple prediction for the inertia tensors of lightly bound skyrmions which, with only two free parameters, fits the field theoretic data for the global minimizers with 1 ≤ B ≤ 8 to within 10%. In judging this, one should bear in mind that an inertia tensor is not a single number, but rather (after accounting for symmetries) 15 independent numbers, so we are actually fitting 120 independent quantities here.
Having checked consistency with field theory simulations for 1 ≤ B ≤ 8, we then proceeded to generate local energy minimizers of the point particle model for 9 ≤ B ≤ 23, where full field simulations are, so far, unavailable. We found that the number of nearly degenerate local energy minimizers grows rapidly with B, that minimizers consistently resemble subsets of a face centred cubic lattice, with internal orientations correlated with lattice position, and that minimizers often have one fewer than the maximum possible number of nearest-neighbour bonds. We have, furthermore, implemented a simple rigid body quantization scheme for all the local minima we found (1 ≤ B ≤ 23). As part of this, we devised an automated algorithm to compute the spin-isospin symmetry group of an oriented point cloud, which simultaneously computes the Finkelstein-Rubinstein constraint associated with each symmetry. This allowed us to compute the spin and isospin of the quantum ground states associated (in rigid body quantization) with each local energy minimizer. Since classical binding energies are so small in the lightly bound model, quantization occasionally altered the energy ordering of local minima with a given charge B. For 12 baryon numbers (out of 23), this simple quantization procedure produced states corresponding to the spin-isospin data of the lightest nucleus of that baryon number.
Our numerical scheme to find energy minimizers of the point particle model had two steps: first a crystal-building algorithm was run to generate a subset of the FCC lattice with sufficiently many (nearest neighbour) bonds. Given such a subset, an initial point particle configuration was constructed with particles at the occupied vertices, their internal orientations being fixed by a simple colouring rule, the lattice length scale being chosen to minimize total interaction energy. The second step was to relax this initial FCC subset using a simulated annealing algorithm which allowed the particle positions, and internal orientations, to vary continuously. The results suggest that, in retrospect, this second step is actually superfluous, since the relaxed configuration always stays very close to some FCC lattice (see figure 6). If one merely wishes to find good approximations to classical energy minimizers of the lightly bound Skyrme model, it would seem that considering only FCC lattice subsets, with the lattice scale left as a free parameter, is a fast and effective strategy. It is possible that low energy FCC subsets may also provide useful sets of initial data for energy minimization in more standard variants of the Skyrme model. Certainly this is a quick and convenient means to generate rather uniform initial data of a qualitatively new kind, not obtainable from rational map or alpha particle clustering methods.
A more interesting problem is to find a better quantization scheme than rigid body quantization. In principle, one could attempt to solve the full Schödinger equation onC B , subject to the FR constraints. For B = 2, there is sufficient symmetry that this may well be tractable. For larger B, however, it is clearly hopeless. Instead, one should attempt to implement some form of "vibrational" quantization scheme, as used for the conventional Skyrme model in [6,7]. This requires one to find a low-dimensional moduli space of configurations, including all relevant local energy minima, but also configurations interpolating between them, which captures the most important vibrational processes of the classical skyrmion. In this regard, the full point particle model, with positions and orientations allowed to leave the set of FCC configurations, will be essential. Indeed the model of point skyrmions introduced here may well prove to be an ideal testing ground for vibrational quantization techniques.