text:
# SKYRMIONS, FULLERENES AND RATIONAL MAPS

## 1 Introduction

### 1.1 Overview

The possibility that solitons can be used to represent particles is an attractive one, very much at the heart of current ideas in high energy physics. The first such model, known as the Skyrme model , was proposed in 1961 as a theory for the strong interactions of pions. The resulting non-linear field theory admits topological soliton solutions, which became known as Skyrmions. In the context of nuclear physics, the topological charge which stabilizes these solitons was identified with baryon number and hence the solitons themselves were identified with baryons.

This model was set aside after the advent of gauge theories and Quantum Chromodynamics (QCD) in the late 1960’s, but much later  it was shown to be a low-energy effective action for QCD, in the limit of the number of colours $`(N_\mathrm{c})`$ being large. Subsequent work has shown the model to be capable of describing at least some aspects of the low-energy behaviour of hadrons. In particular, it was shown that the properties of the proton and neutron could be adequately described using the simple quantization of the solution with a single unit of topological charge , and that the other static configurations appeared to be most stable when considering charges corresponding to $`{}_{}{}^{4}\mathrm{He}`$ and $`{}_{}{}^{7}\mathrm{Li}`$ .

Although understanding the properties of light nuclei is the main motivation for this paper, we shall only comment very briefly at various points in our discussion on the implications of our results for this application. Instead, we will concentrate on the Skyrmions themselves, which are of interest in their own right. As examples of three-dimensional topological solitons, they have been seen to be very similar to BPS monopoles . From a mathematical point of view they can be thought of as maps between 3-spheres and, therefore, the minimum energy configurations which we compute are the minimum energy maps relative to the Skyrme energy functional. Although there are no doubt many other possibilities for an energy functional on the space of maps between 3-spheres, the Skyrme functional is the simplest which supports topological structures and, therefore, it is also interesting to speculate as to their generality in this context.

In the following we will present an extensive study of minimum energy solitons using two very different numerical methods111The bare essentials of this work were first reported in ref. . Here, we give a more detailed exposition of our results.. With a small number of caveats, the two approaches come to the same conclusion: that there is a connection with fullerene cages  familiar in carbon chemistry, as conjectured in ref. , and that the solutions can be represented in terms of the rational map ansatz . For each topological charge, we will discuss the structure and symmetries of the solution in detail, and present an explicit analytic formula for the approximate solution. In the small number of cases where there is some deviation from the fullerene hypothesis, we will discuss qualitatively, and sometimes quantitatively, possible reasons. Finally, in the context of the two applications mentioned above, nuclear physics and generalized harmonic maps between 3-spheres, we will attempt to provide some heuristic insight into why our conclusions as to the structure of Skyrmions are consistent with the particular application.

### 1.2 Skyrmions

In terms of the algebra valued currents $`R_\mu =(_\mu U)U^{}`$ of an SU(2) valued field $`U`$, the Skyrme Lagrangian density is,

$$=\frac{1}{24\pi ^2}\left[\mathrm{Tr}\left(R_\mu R^\mu \right)+\frac{1}{8}\mathrm{Tr}\left([R_\mu ,R_\nu ][R^\mu ,R^\nu ]\right)\right].$$
(1.1)
At first sight the domain is $`\text{I}\text{R}^3`$ and the target space is the group manifold of SU(2), $`S^3`$, but the finite energy boundary condition, $`U(\mathrm{})=I`$, means that $`U`$ is in fact a map from compactified $`\text{I}\text{R}^3S^3S^3`$. Such mappings have non-trivial homotopy classes characterized by $`\pi _3(S^3)=\text{Z}\text{Z}`$, which has the explicit representation

$$B=\frac{1}{24\pi ^2}ϵ_{ijk}d^3x\mathrm{Tr}\left(R_iR_jR_k\right).$$
(1.2)
A more geometrical description of the model can be made in terms of the strain tensor defined at each point $`𝐱`$ in the domain by

$$D_{ij}=\frac{1}{2}\mathrm{Tr}\left(R_iR_j\right),$$
(1.3)
which can be thought of as quantifying the deformation induced by the map between 3-spheres. This symmetric, positive definite tensor can be diagonalized with non-negative eigenvalues $`\lambda _1^2`$, $`\lambda _2^2`$, $`\lambda _3^2`$, and the static energy, $`E`$, and baryon number, $`B`$, can be computed as integrals over $`\text{I}\text{R}^3`$ of the corresponding densities $``$ and $``$ given by

$$=\frac{1}{12\pi ^2}(\lambda _1^2+\lambda _2^2+\lambda _3^2+\lambda _1^2\lambda _2^2+\lambda _2^2\lambda _3^2+\lambda _3^2\lambda _1^2),=\frac{1}{2\pi ^2}\lambda _1\lambda _2\lambda _3.$$
(1.4)
A simple manipulation of these expressions allows one to deduce the Faddeev-Bogomolny bound $`E|B|`$, but in contrast to monopoles and vortices this bound cannot be saturated for any non-trivial finite energy configuration. In fact, it can be attained only when all the eigenvalues of the strain tensor are equal to one at all points in space — an isometry — and this is clearly not possible since $`\text{I}\text{R}^3`$ is not isometric to $`S^3.`$ In practice we shall see that minimization of the Skyrme energy functional associated with the Lagrangian density requires that the solution be as close to this bound as possible. The energy minimization can be thought of as finding the map which is as close to an isometry as possible, when averaged over space .

The boundary condition breaks the chiral symmetry ($`\mathrm{SU}(2)\times \mathrm{SU}(2)`$) of the Skyrme model to an SO(3) isospin symmetry $`U𝒪U𝒪^{}`$, where $`𝒪`$ is a constant element of SU(2). When we refer to a spatial symmetry of a Skyrmion, such as spherical symmetry, the fields are not invariant under a spatial rotation, but rather there is an equivariance property that the effect of a spatial rotation can be absorbed into an isospin transformation. This implies that both the energy density, $``$, and baryon density, $``$, are strictly invariant under the symmetry.

The $`B=1`$ Skyrmion is spherically symmetric , the maximally allowed symmetry of the Lagrangian, and its energy is $`E=1.232`$ . In fact, spherically symmetric solutions exist at all charges, but they are not the minimum energy configurations, which are less symmetric. The $`B=2`$ solution is axially symmetric  and the higher charge solutions all have point symmetries  which are subgroups of O(3). For $`B=3,4,7`$ the Skyrmions have the Platonic symmetries of the tetrahedron ($`T_d`$), the cube ($`O_h`$) and the dodecahedron ($`Y_h`$) respectively, while for $`B=5,6,8`$ the Skyrmions have the dihedral symmetries $`D_{2d}`$, $`D_{4d}`$ and $`D_{6d}`$ respectively (see the discussion in section 1.5 if you are unfamiliar with these point groups). For $`B=9`$ a tetrahedrally symmetric Skyrmion has been found, but as we shall discuss later this appears not to be the minimum energy configuration, and hence is probably a low-energy saddle point.

In all the above cases the baryon density (and also the energy density) is localized around the edges of a polyhedron. From these known results we were able, in a previous paper , to formulate some simple geometrical rules222The original Geometric Energy Minimization (GEM) rules were that all the solutions computed for $`B9`$ had the symmetries of almost spherical trivalent polyhedra with $`4(B2)`$ vertices, $`2(B1)`$ faces and $`6(B2)`$ edges. for the structure of Skyrmions which led us to conjecture that higher charge Skyrmions ($`B7`$) would resemble trivalent polyhedra formed from 12 pentagons and $`2B14`$ hexagons. We will refer to such structures as fullerene-like and to the conjecture as the fullerene hypothesis since precisely the same structures arise in carbon chemistry where carbon atoms sit at the vertices of such polyhedra, known as fullerenes . On the basis of this, it was suggested that the minimum energy Skyrmion of charge $`B`$, $`𝒮_B`$, would have the same symmetry as a fullerene from the family $`\mathrm{C}_{4(B2)}`$. For low charges ($`B=7`$, $`B=8`$) this leads to a unique prediction for $`𝒮_B`$, and indeed this was what we found in our original simulations. But as the charge increases the number of possible structures increases, in particular for $`B=9`$ there are 2 possibilities with $`D_2`$ and $`T_d`$ symmetries respectively, for $`B=10`$ there are 6, for $`B=11`$ there are 15, with a rapid increase for $`B>11`$. Using this simple analogy, it was possible to predict that there would be an icosahedral configuration with $`B=17`$ corresponding to the famous Buckminsterfullerene structure of $`\mathrm{C}_{60}`$, and given its highly symmetric structure we suggested that this would be the minimum energy configuration at that charge.

Our original work on computing minimum energy Skyrmions  and also that of ref.  used the full non-linear field equations, or the full non-linear energy functional; a very resource hungry procedure even in the modern era of parallel supercomputers. It is also likely to be somewhat imprecise since (1) it can be very difficult to identify the particular symmetries of the solution that one computes in this way, even using sophisticated visualization packages, and (2) the choice of initial conditions is somewhat arbitrary: even when they are chosen to have no particular symmetry, it is impossible to guarantee that they will relax to the global minimum. Fortunately, help is at hand in the form of the rational map ansatz , which was devised as an approximate representation of the solutions computed in ref. . Although this representation is not exact it reduces the degrees of freedom in the problem to be a finite, manageable number, allowing computations to take place in an acceptable amount of time. Of course, the original approach based on the full non-linear equations still has a role to play in firstly checking and then fully relaxing the approximate solutions computed in the rational map approach. The development of this subject now treads an interesting interface between numerically generated solutions and this analytic approximation.

In this paper we shall first compute what are the minimum energy solutions of the Skyrme model upto $`B=22`$, under the assumption that they can be adequately represented by the rational map ansatz. Then by simulating collisions of well separated Skyrmion clusters using the full non-linear field equations, followed by a numerical relaxation of the resulting coalesced cluster, we shall attempt to verify that these structures are in fact the minimum energy Skyrmions. We find precisely the fullerene structures conjectured in ref.  for all but a small number of cases where there are interesting caveats. We also demonstrate that not only are the minimal energy Skyrmions fullerene-like, but that there are also several other fullerene Skyrmions at a given baryon number whose energies are only very slightly higher than the minimal value. Finally, using the full non-linear equations once again, we relax the rational map generated solutions with specified symmetries to compute accurately the energies of the configurations.

### 1.3 Rational map ansatz

The rational map ansatz was introduced in ref. , and is a way to construct approximate Skyrmions from rational maps between Riemann spheres. Briefly, we use spherical coordinates in $`\text{I}\text{R}^3`$, so that a point $`𝐱\text{I}\text{R}^3`$ is given by a pair $`(r,z)`$, where $`r=|𝐱|`$ is the distance from the origin, and $`z`$ is a Riemann sphere coordinate giving the point on the unit two-sphere which intersects the half-line through the origin and the point $`𝐱`$.

Now, let $`R(z)`$ be a degree $`B`$ rational map between Riemann spheres, that is, $`R=p/q`$ where $`p`$ and $`q`$ are polynomials in $`z`$ such that $`\mathrm{max}[\text{deg}(p),\text{deg}(q)]=B`$, and $`p`$ and $`q`$ have no common factors. Given such a rational map the ansatz for the Skyrme field is

$$U(r,z)=\mathrm{exp}\left[\frac{if(r)}{1+|R|^2}\left(\begin{array}{cc}1|R|^2& 2\overline{R}\\ 2R& |R|^21\end{array}\right)\right],$$
(1.5)
where $`f(r)`$ is a real profile function satisfying the boundary conditions $`f(0)=\pi `$ and $`f(\mathrm{})=0`$, which is determined by minimization of the Skyrme energy of the field (1.5) given a particular rational map $`R`$. It can be shown that this Skyrme field has charge $`B`$, and for $`1B9`$ rational maps were presented in ref.  which reproduced Skyrmions with the same symmetries as those computed in ref. . Furthermore, they were shown to have energies which are only about one or two percent above the numerically calculated values.

Substitution of the rational map ansatz (1.5) into the Skyrme energy functional results in the following expression for the energy

$$E=\frac{1}{3\pi }\left(r^2f^2+2B(f^2+1)\mathrm{sin}^2f+\frac{\mathrm{sin}^4f}{r^2}\right)𝑑r,$$
(1.6)
where $``$ denotes the integral

$$=\frac{1}{4\pi }\left(\frac{1+|z|^2}{1+|R|^2}\left|\frac{dR}{dz}\right|\right)^4\frac{2idzd\overline{z}}{(1+|z|^2)^2}.$$
(1.7)
To minimize the energy (1.6), therefore, one first determines the rational map which minimizes $``$, which may be thought of as an energy functional on the space of rational maps. Then given the minimum value of $``$ it is a simple exercise to find the profile function which minimizes the energy (1.6) using a gradient flow method to solve the appropriate boundary value problem. Thus, within the rational map ansatz, the problem of finding the minimal energy Skyrmion reduces to the simpler problem of calculating the rational map which minimizes the function $``$. Computing the map which minimizes this set up is the essence of our procedure for finding the minimal energy Skyrmion, and in section 2 we shall describe our numerical techniques used to address this problem.

The baryon density is proportional to the derivative of the rational map, and (counting multiplicities) this will have $`2B2`$ zeros, giving the points on the Riemann sphere for which the baryon density vanishes along the corresponding half-lines through the origin. In terms of a baryon density isosurface plot these correspond to holes in a shell-like structure which resembles a polyhedron and the holes correspond to the face centres.

### 1.4 Symmetric maps : general discussion

Since the maps we shall be dealing with describe symmetric Skyrmions, let us recall what it means for a rational map (and hence the associated Skyrmion) to be symmetric under a group $`GSO(3).`$ Consider a spatial rotation $`gSO(3)`$, which acts on the Riemann sphere coordinate $`z`$ as an $`SU(2)`$ Möbius transformation

$$zg(z)=\frac{\gamma z+\delta }{\overline{\delta }z+\overline{\gamma }},\text{where}|\gamma |^2+|\delta |^2=1.$$
(1.8)
Similarly a rotation, $`DSO(3)`$, of the target two-sphere (which corresponds to an isospin transformation) will act in the same way

$$RD(R)=\frac{\mathrm{\Gamma }R+\mathrm{\Delta }}{\overline{\mathrm{\Delta }}R+\overline{\mathrm{\Gamma }}},\text{where}|\mathrm{\Gamma }|^2+|\mathrm{\Delta }|^2=1.$$
(1.9)
A map is $`G`$-symmetric if, for each $`gG`$, there exists a target space rotation, $`D`$, which counteracts the effect of the spatial rotation, that is,

$$R(g(z))=D(R(z)).$$
(1.10)
Note that in general the rotations on the domain and target spheres will not be the same, so that $`(\gamma ,\delta )(\mathrm{\Gamma },\mathrm{\Delta }).`$

Since we are dealing with $`SU(2)`$ transformations the set of target space rotations will form a representation of the double group of $`G`$, which is the group of order $`2|G|`$ obtained from $`G`$ by the addition of an element $`\overline{E}`$ which squares to the identity. The fact that we are dealing with the double group is important since it has representations which are not representations of $`G.`$ From now on it is to be understood that when we refer to a group $`G`$ we shall actually mean its double group.

To determine the existence and compute particular symmetric rational maps is, therefore, a matter of classical group theory. We are concerned with degree $`B`$ polynomials which form the carrier space for $`\underset{¯}{B+1}`$, the $`(B+1)`$-dimensional irreducible representation of $`SU(2).`$ Now, as a representation of $`SU(2)`$ this is irreducible, but if we only consider the restriction to a subgroup $`G`$, $`\underset{¯}{B+1}|_G`$, this will in general be reducible. What we are interested in is the irreducible decomposition of this representation and tables of these subductions can be found, for example, in ref. .

The simplest case in which a $`G`$-symmetric degree $`B`$ rational map exists is if

$$\underset{¯}{B+1}|_G=E+\mathrm{}$$
(1.11)
where $`E`$ denotes a two-dimensional representation. In this case a basis for $`E`$ consists of two degree $`B`$ polynomials which can be taken to be the numerator and denominator of the rational map. A subtle point which needs to be addressed is that the two basis polynomials may have a common root, in which case the resulting rational map is degenerate and does not correspond to a genuine degree $`B`$ map.

More complicated situations can arise, for example, if

$$\underset{¯}{B+1}|_G=A_1+A_2+\mathrm{}$$
(1.12)
where $`A_1`$ and $`A_2`$ denote two one-dimensional representations, then a whole one-parameter family of maps can be obtained by taking a constant multiple of the ratio of the two polynomials which are a basis for $`A_1`$ and $`A_2`$ respectively. An $`m`$-parameter family of $`G`$-symmetric maps can be constructed if the decomposition contains $`(m+1)`$ copies of a two-dimensional representation, that is,

$$\underset{¯}{B+1}|_G=(m+1)E+\mathrm{}$$
(1.13)
where the $`m`$ (complex) parameters correspond to the freedom in the decomposition of $`(m+1)E`$ into $`m+1`$ copies of $`E.`$

Explicit examples corresponding to the above types of decompositions will be constructed in section 1.5. For a detailed explanation of how to calculate these maps by computing appropriate projectors see ref. .

### 1.5 Symmetric maps : specific examples

As we shall discuss in section 2.2 the basic procedure for finding the minimum energy Skyrmion in the rational map ansatz will be to minimize the function $``$ subject to the map having degree $`B`$. However, this will find the map in an arbitrary spatial orientation, preventing identification of the symmetry directly from the map. One can always compute the corresponding baryon density isosurface and identify the symmetry by eye, a procedure which is often helpful, but this is also fraught with difficulties, particularly when, for example, the solution only has a small number of symmetry generators. Therefore, in order to be sure of the symmetry identification we will also search maps which are restricted to have a particular symmetry.

Since all the point groups that we shall consider are subgroups of O(3), they must be either cyclic groups, $`C_n`$, which involve invariance under rotations by $`(360/n)^{}`$ about some axis, dihedral groups, $`D_n`$, which are obtained from the cyclic group by the addition of a $`C_2`$ axis which is perpendicular to the main symmetry axis, tetrahedral groups ($`T`$), which are the symmetries associated with the tetrahedron, octahedral groups ($`O`$), those associated with the octahedron/cube, or icosahedral groups ($`Y`$), those associated with the icosahedron/dodecahedron. Each of these symmetry groups can be extended by the inclusion of reflections. All the icosahedral maps presented in this paper have already been discussed in ref. , while all the octahedral maps used are easily deduced from the tetrahedral maps discussed below, so we shall only concentrate on understanding the details of the dihedral and tetrahedral maps.

In terms of the Riemann sphere coordinate $`z`$ the generators of the dihedral group $`D_n`$ may be taken to be $`ze^{2\pi i/n}z`$ and $`z1/z.`$ This can be extended by the addition of a reflection symmetry in two ways: by including a reflection in the plane perpendicular to the main $`C_n`$ axis, which is represented on the Riemann sphere by invariance under $`z1/\overline{z}`$, and the group $`D_{nh}`$ is obtained. Alternatively, a reflection symmetry may be imposed in a plane which contains the main symmetry axis and bisects the $`C_2`$ axes obtained by applying the $`C_n`$ symmetry to the $`C_2`$ axis. This reflection is represented on the Riemann sphere as invariance under $`ze^{\pi i/n}\overline{z}`$, and the resulting group is $`D_{nd}.`$

To construct $`D_n`$ symmetric maps does not require any group theory formalism discussed in section 1.4 since it is a simple task to explicitly apply the two generators of $`D_n`$ to a general degree $`B`$ rational map to determine a family of symmetric maps. Explicitly, an $`s`$-parameter family is given by333There are other $`D_n`$ symmetric families of maps in addition to those of the form (1.14), but they will not be needed in this paper.

$$R(z)=\left(\underset{j=0}{\overset{s}{}}a_jz^{jn+r}\right)/\left(\underset{j=0}{\overset{s}{}}a_{sj}z^{jn}\right),$$
(1.14)
where $`r=B\text{mod }n`$ and $`s=(Br)/n.`$ Here $`a_s=1`$ and $`a_0,..,a_{s1}`$ are arbitrary complex parameters. Clearly, this map satisfies the conditions for it to be symmetric under $`D_n`$,

$$R(e^{2\pi i/n}z)=e^{2\pi ir/n}R(z),R(1/z)=1/R(z),$$
(1.15)
and imposing a reflection symmetry constrains the otherwise complex coefficients $`a_j`$ to either be real, or pure imaginary. In the case of $`D_{nh}`$ symmetry the condition is that all $`a_j`$ are real, whereas for $`D_{nd}`$ symmetry the coefficient $`a_j`$ is either real or purely imaginary depending on whether $`(sj)\text{mod 2 }`$ is zero or one respectively

The procedure for constructing tetrahedrally symmetric maps is more difficult than for dihedral symmetries and the systematic group theory approach discussed in section 1.4 must be employed. No simple formula such as (1.14) exists, so for later reference we shall, therefore, need to recall the basic facts about the irreducible representations of the tetrahedral group $`T.`$

$`T`$ has three one-dimensional representations, which are the trivial representation, $`A`$, and two conjugate representations $`A_1`$ and $`A_2.`$ There is also a three-dimensional representaion, $`F`$, which is obtained as $`\underset{¯}{3}|_T.`$ In addition to these representations there are three two-dimensional representations of the double group of $`T`$, which we denote by $`E^{},E_1^{},E_2^{}`$, where the prime signifies that these are not representations of $`T`$, but only of the double group of $`T.`$ $`E^{}`$ is obtained as $`\underset{¯}{2}|_T`$ and $`E_1^{}`$ and $`E_2^{}`$ are conjugate representations.

## 2 Numerical minimization algorithms

### 2.1 Overview

Minimization of an energy functional is a classical numerical problem, with no hard and fast optimum method. Methods which are tailored for a particular application can work very badly in others. In this section we will outline the basic features of the numerical methods which we have employed to compute minimum energy Skyrmions with and without using the rational map ansatz. We will attempt to discuss both the advantages and disadvantages of the two methods.

Our original approach to this problem was to use the code first used in ref.  to evolve the full field equations for well-separated Skyrmions, and this is discussed in section 2.3. It worked well for $`B<9`$  and its results were the main motivation for the rational map ansatz. Its disadvantages are that it is very slow, requiring many hours of CPU time on a parallel computer, and in circumstances where there is the possibility of two or more minima separated by a small energy gap, dependence on the choice of initial conditions is also an issue. It is, however, the only way in which the results of approximate methods, such as the rational map ansatz, can be checked. Creating a particular configuration from many different initial conditions using this method can be thought of as strong evidence for it being the minimum energy solution, irrespective of other considerations.

The simulated annealing of rational maps as discussed in section 2.2 is by contrast fast — it only requires a serial processor — by virtue of the fact that the number of degrees of freedom have been substantially reduced, is relatively independent of the initial conditions and much less sensitive to local minima. But, of course, the results are only as good as the rational map approach to describing Skyrmions. Thankfully, as we shall see, the rational map approach in general works very well, allowing us to generate symmetric Skyrmions with large baryon number.

### 2.2 Simulated Annealing of Rational Maps

Simulated annealing is a fairly recent numerical method for obtaining the global minimum of an energy function, and is based on the way that a solid cools to form a lattice . Let $`E(𝐚)`$ be an energy function which depends on a number of parameters $`𝐚.`$ The idea is that at a given temperature, $`T`$, the system is allowed to reach thermal equilibrium, characterized by the probability of being in a state with energy $``$ given by the Boltzmann distribution

$$\text{Pr}(E=)=\frac{1}{Z(T)}\mathrm{exp}(/T),$$
(2.1)
where $`Z(T)`$ is the partition function.

In practice, this is achieved by applying a Metropolis method: starting with a given configuration, $`𝐚`$, a small random perturbation $`\delta 𝐚`$ is made and the energy of the resulting configuration is computed. If the change in energy, $`\delta E=E(𝐚+\delta 𝐚)E(𝐚)`$, is negative then the new configuration is accepted, that is, $`𝐚`$ is replaced by $`𝐚+\delta 𝐚.`$ However, if the change results in an increase in energy then the probability of accepting the new configuration is $`e^{\delta E/T}.`$ By performing a large number of such perturbations thermal equilibrium can be achieved at temperature $`T`$.

The procedure to minimize the energy is to start at a high temperature, bring the system into thermal equilibrium and then lower the temperature before regaining the equilibrium. As the temperature is decreased, the system is more likely to be found in a state with lower energy and in the limit as $`T0`$ the configuration will move toward a minimum of $`E.`$ In the limit of infinitesimally slow variations in the temperature this can be shown to be the global minimum.

From the above description it is immediately clear that the simulated annealing method has a major advantage over other conventional minimization techniques in that changes which increase the energy are allowed, enabling the algorithm to escape from minima that are not the global minimum. Of course, in practice one is not guaranteed to find the global minimum since the number of iteration loops used to bring the system into thermal equilibrium at a fixed temperature and also the number of times the temperature can be decreased are both restricted by computational resources. However, it does provide the most efficient means for searching for a global minimum and with sufficient computational resources, plus sufficient care in applying them, one can be fairly confident of the final result.

For our application to rational maps we obviously take the energy function to be $``$ and the parameters $`𝐚`$ to be the constants in the rational map. To compute $``$ involves a numerical integration over the sphere, which can be performed with standard methods, and in a typical simulation this needs to be calculated approximately a million times for a full simulated annealing run. In each case we take the initial rational map to be the axially symmetric one $`R=z^B,`$ which for large $`B`$ has a very high value of $`.`$

Note that since we are using the Riemann sphere coordinate $`z`$, the two-sphere metric is the Fubini-Study metric in this coordinate system. This means that in terms of moving energy around the sphere there is a bias between points near the south pole as compared to those near the north pole, in terms of small variations of the rational map parameters. To counteract this discrepancy we perform a spatial rotation of the configuration, $`z1/z`$ each time the temperature is decreased.

To identify the rational map, and more importantly its symmetries, produced by the simulated annealing algorithm is a two stage process. First, we compute the minimum energy map assuming no particular symmetry, that is, we allow a general map of degree $`B`$. As already pointed in section 1.5 the end result will be a rational map which is in a random orientation in both the domain and target two-spheres. The target space orientation could easily be fixed (in fact, it is more convenient not to do this, due to the above comments regarding a periodic spatial rotation of the map during the minimization procedure), but there is no simple way to make the spatial orientation such that the symmetry generators of the map are conveniently represented. Once one has the minimizing degree $`B`$ rational map and the corresponding minimum value of $`,`$ the Skyrmion is then constructed and its baryon density is plotted and examined in an attempt to identify its symmetries by eye. This conjectured symmetry is then confirmed by constructing the most general map with this symmetry and minimizing within this constrained symmetric family to check that the same minimum value of $``$ is recovered. As a final check the corresponding Skyrmion is constructed and its baryon density examined to confirm that it is identical to the one obtained previously.

For each charge we have also performed several simulated annealing runs within constrained symmetric families, such as $`D_2`$, $`D_3`$ and $`D_4`$. Since all the symmetry groups must be subgroups of O(3), the number of possibilities is finite and checking just these three possibilities, allows one to rule out a large fraction of them; the rest often being possible by eye. This not only provides an additional check that the minimizing map was found, but also allows us to obtain other low energy maps, which may be either saddle points or local minima.

Finally, it is perhaps worth mentioning that a simulated annealing method has recently been used in a rather different way to study Skyrmions . These authors used a simulated annealing algorithm on a discretized version of the full Skyrme energy, taking the parameters to be the field values at the discretized lattice sites. This is a much more computationally expensive approach than using the rational map ansatz, which is probably the reason these authors only considered Skyrmions upto charge $`B=4`$, where, of course, the results were already known. Nonetheless, it is interesting to know that a simulated annealing approach is viable in this manner and demonstrates yet another application of this versatile technique.

### 2.3 Full Field Dynamics

#### 2.3.1 Sigma Model Formulation

In the $`SU(2)`$ form, the equations of motion are cumbersome to handle numerically, so we convert to the notation of a non-linear sigma model (NLSM), which has Lagrangian,

$$=_\mu \varphi ^\mu \varphi \frac{1}{2}\left(_\mu \varphi ^\mu \varphi \right)^2+\frac{1}{2}\left(_\mu \varphi _\nu \varphi \right)\left(^\mu \varphi ^\nu \varphi \right)+\lambda (\varphi \varphi 1),$$
(2.2)
with the Lagrange multiplier $`\lambda `$ introduced to maintain the constraint $`\varphi \varphi =1`$.

The Euler-Lagrange equations are given by

$$\left(1_\mu \varphi ^\mu \varphi \right)\text{    }\text{ }\varphi \left(^\nu \varphi _\mu _\nu \varphi _\mu \varphi \text{    }\text{ }\varphi \right)^\mu \varphi +\left(^\mu \varphi ^\nu \varphi \right)_\mu _\nu \varphi \lambda \varphi =0,$$
(2.3)
where the Lagrange multiplier can be calculated by contracting (2.3) with $`\varphi `$ and using the second derivative of the constraint,

$`\lambda =`$ $`\left(1_\mu \varphi ^\mu \varphi \right)\varphi \text{   }\text{ }\varphi +\left(^\mu \varphi ^\nu \varphi \right)\left(\varphi _\mu _\nu \varphi \right)`$ (2.4)
$`=`$ $`\left(_\mu \varphi _\nu \varphi \right)\left(^\mu \varphi ^\nu \varphi \right)(1_\mu \varphi ^\mu \varphi )_\nu \varphi ^\nu \varphi .`$ (2.5)
Denoting differentiation with respect to time as a dot, these equations can be recast as

$$M\ddot{\varphi }\alpha (\dot{\varphi },_i\varphi ,_i\dot{\varphi },_i_j\varphi )\lambda \varphi =0,$$
(2.6)
where the symmetric matrix $`M`$ has elements

$$M_{ab}=\left(1+_j\varphi _j\varphi \right)\delta _{ab}_j\varphi _a_j\varphi _b,$$
(2.7)
and $`\alpha `$ is given by

$`\alpha =\left(\dot{\varphi }_i_i\varphi _i\varphi _i\dot{\varphi }\right)\dot{\varphi }+2\left(\dot{\varphi }_i\varphi \right)_i\dot{\varphi }\left(\dot{\varphi }_i\dot{\varphi }\right)_i\varphi \dot{\varphi }^2_i_i\varphi `$ (2.8)
$`+\left(_i\varphi _i_j\varphi _j\varphi _i_i\varphi \right)_j\varphi +\left(1+_j\varphi _j\varphi \right)_i_i\varphi \left(_i\varphi _j\varphi \right)_i_j\varphi .`$ (2.9)
Quite clearly these equations of motion are not analytically tractable. In subsequent sections we will discuss our numerical methods for evolving the equations of motion for spatially discretized initial conditions and for obtaining minimal energy static Skyrmions. As we shall see this is still a highly non-trivial task and can only be done for a specialized, but ill-defined, set of initial conditions.

#### 2.3.2 Discretization and boundary conditions

There are three aspects common to almost all numerical approaches to solving non-linear PDE’s. The first is a spatial discretization and an approximation for spatial derivatives. We discretized on a regular, cubic grid with $`N`$ points in each of the Cartesian directions and the array $`\varphi _{i,j,k}\varphi (i\mathrm{\Delta }x,j\mathrm{\Delta }x,k\mathrm{\Delta }x)`$. The choice of $`N`$ and $`\mathrm{\Delta }x`$ is of critical importance, since the soliton configurations we wish to represent are localized. We found that grids with $`N=100`$ and $`\mathrm{\Delta }x=0.1`$ were convenient for representing all the configurations which we study in this paper, although larger grid spacings $`(\mathrm{\Delta }x=0.2)`$ also give sensible results, and larger grids ($`N=200`$) were used to obtain more accurate calculations of energies.

The spatial derivatives used were fourth order, so as to accurately represent the large spatial gradients of the solitonic configurations. Since the reader may not be totally familiar with this procedure, various expressions for derivatives are presented below: for first order derivatives,

$$\frac{\varphi }{x}=\frac{\varphi _{i+2,j,k}+8\varphi _{i+1,j,k}8\varphi _{i1,j,k}+\varphi _{i2,j,k}}{12\mathrm{\Delta }x}+𝒪(\mathrm{\Delta }x^4),$$
(2.11)
for second order derivatives,

$$\frac{^2\varphi }{x^2}=\frac{\varphi _{i+2,j,k}+16\varphi _{i+1,j,k}30\varphi _{i,j,k}+16\varphi _{i1,j,k}\varphi _{i2,j,k}}{12\mathrm{\Delta }x^2}+𝒪(\mathrm{\Delta }x^4),$$
(2.12)
and for mixed second order derivatives

$`{\displaystyle \frac{^2\varphi }{x^2}}+2{\displaystyle \frac{^2\varphi }{xy}}+{\displaystyle \frac{^2\varphi }{y^2}}=𝒪(\mathrm{\Delta }x^4)+`$ (2.13)
$`(\varphi _{i+2,j+2,k}+16\varphi _{i+1,j+1,k}30\varphi _{i,j,k}+16\varphi _{i1,j1,k}\varphi _{i2,j2,k})/(12\mathrm{\Delta }x^2).`$
The next part of the procedure is a method for time evolution. The equations of motion can be transformed into first order form,

$$M\dot{\psi }\alpha (\psi ,_i\varphi ,_i\psi ,_i_j\varphi )\lambda \varphi =0,$$
(2.14)
by defining $`\psi =\dot{\varphi }`$, and this can be solved using a leapfrog method. This involves replacing

$$\dot{\varphi }=\frac{\varphi ^+\varphi ^{}}{2\mathrm{\Delta }t}+𝒪(\mathrm{\Delta }t^2),\dot{\psi }=\frac{\psi ^+\psi ^{}}{2\mathrm{\Delta }t}+𝒪(\mathrm{\Delta }t^2),$$
(2.15)
where $`+`$ and $``$ correspond to the values of $`\psi `$ and $`\varphi `$ at one step after and one before respectively. In a much simpler case such as the wave equation, this creates a decoupling of the arrays containing the discretized versions $`\varphi `$ and $`\psi `$. However, the non-linear dependence of the function $`\alpha `$ on $`\psi `$ requires the storing of two copies of each array, and hence four in total, requiring 64Mb of core memory for $`N=100`$. The choice of $`\mathrm{\Delta }t`$ is also crucial, since it can create numerical instability. The standard Courant condition for a linear equation in three dimensions states that $`\sqrt{3}\mathrm{\Delta }t<\mathrm{\Delta }x`$. The non-linear nature of the Skyrme equations leads to a non-trivial modification to this relation, which still stands as the best possible due to reasons of causality. We find, essentially by trial and error, that $`\mathrm{\Delta }t\mathrm{\Delta }x/10`$ leads to a stable algorithm. We should note at this stage that there is another, potentially more pathological, instability of the Skyrme model which is discussed in a subsequent section.

Finally, one is required to specify boundary conditions for the finite grid employed. Since the fourth order spatial derivatives require a five point wide stencil, one point away from the boundary it is necessary to use second order spatial approximations. But on the boundary itself the spatial derivatives cannot be evaluated, since the second order spatial approximation requires a three point wide stencil. We experimented with various different types of boundary conditions, such as Neumann (zero normal derivative), Dirichlet (fixed) and periodic. The results presented here are for Dirichlet boundary conditions, although we believe that the use of Neumann boundaries would have little effect on the results.

#### 2.3.3 Imposing the constraint

In the previous section we have discussed all aspects of the numerical solution of the Skyrme equations of motion which are common to the numerical solution of most non-linear PDE’s. There is, however, an added extra which makes life much more difficult in the case of a NLSM. In the previous section we included the Lagrange multiplier $`\lambda `$, assuming that it can be calculated from $`\psi `$, $`\varphi `$ and their spatial derivatives. In fact, this leads to a numerical scheme which becomes unstable for any choice of $`\mathrm{\Delta }t`$ in under ten timesteps. The problem is that we have ignored the reason for its introduction; to maintain the constraint $`\varphi \varphi =1`$, which is manifest in the NLSM. A number of approaches have been developed to deal with such constraints.

Firstly, one could modify the numerical scheme to calculate $`\lambda `$ so that it explicitly maintains the constraint for the discretized equations of motion. At each step, this is seen to be almost equal to calculating $`\lambda `$ from the formula (2.5), but the two differ by numerical discretization effects at a level well below 1%, which nonetheless cause the solution to slip off the unit sphere, if allowed to accumulate.

An alternative approach, which was found to work well in simulations of Baby Skyrmions , is to simply rescale the field to have unit modulus, that is, to continually make the replacement

$$\varphi \frac{\varphi }{\sqrt{\varphi \varphi }},$$
(2.16)
at each point on the discretized grid after each timestep. While appearing ugly from a purist numerical analysis point of view, this technique is effective and does not become unstable except in the most extreme circumstances. One is effectively projecting the field back onto the sphere along the field itself, which has no particular physical motivation, but if the modification is small, which can be arranged by choosing a sufficiently small timestep, then this should be as good as any other arbitrary choice.

Another possibility is to require that the derivative of the constraint is zero: the relation $`\varphi \varphi =1`$ not only implies that the solution lies on the unit sphere, but that it cannot come off, that is, all the derivatives of the constraint are also satisfied. This leads to an infinite hierarchy of relations which must hold. Obviously, with a second order time evolution one cannot hope to maintain them all explicitly in the discretized system. However, if one manages to satisfy the first one, it may be possible to construct an effective code. It is possible to impose that $`\varphi \psi =0`$ by computing $`\lambda `$ in order to satisfy $`\varphi ^+\psi ^+=0`$ which requires that,

$$\lambda =\frac{\psi ^{}\psi +\varphi M^1\alpha +2\psi M^1\alpha \mathrm{\Delta }t}{\varphi M^1\varphi +2\psi M^1\varphi \mathrm{\Delta }t}.$$
(2.17)
Although we should note that this does not implicitly imply the imposition of the constraint on the discretized equations.

Largely by trial and error, we find that the best method for our numerical scheme is a hybrid of the last two methods. It just so happens that this is possible within our scheme, since there are two discretized grids with essentially their own time evolution. Firstly, we calculate $`\lambda `$ to satisfy $`\varphi ^+\psi ^+=0`$, followed by the rescaling transformation (2.16) on the field $`\varphi ^+`$. We find that this hybrid methods maintains the constraint for many thousands of timesteps.

#### 2.3.4 Non-Hyperbolic regions

We have already discussed the potential for our numerical scheme to become unstable because of Courant-type instability and also due to the imperfect imposition of the constraint $`\varphi \varphi =1`$. However, there is a much more pathological instability which comes about since the equations of motion are not manifestly hyperbolic and their precise nature, hyperbolic, parabolic or elliptic, depends on the particular configuration being evolved .

The problem is that the specific numerical scheme we have designed will work only for the hyperbolic case and we know of no way of treating all configurations within a single numerical scheme. Fortunately, the hyperbolic regime is the only physically meaningful one. One can understand this by thinking of the Skyrme model as a low energy effective action, with higher order terms ignored. At higher energies, where the terms which are ignored would be large, the Skyrme equations of motion become non-hyperbolic.

It was suggested in ref. that the equations of motion became non-hyperbolic whenever the kinetic energy is greater than the potential energy. Empirically, we find that this is at least partially true, with an instability associated with a very large kinetic energy density relative to that of potential energy. However, this statement is a little imprecise since we found that sometimes the kinetic energy density rose above the potential locally on the discretized grid without creating instability. Unfortunately, the complicated nature of the equations of motion makes it almost impossible to say for certain, which configurations will eventually lead to an instability and which will not, although experience has taught us that most of the configurations we wish to evolve for physical applications, such as low energy nuclear physics, are possible.

#### 2.3.5 Locating minima

In the preceding sections we have discussed how to construct a numerical scheme which evolves the full non-linear equations of motion for the Skyrme model. This allowed us to simulate the dynamics of Skyrmion collisions in ref. . However, for the purposes of this paper one would also like to create static multi-soliton configurations. We are assisted in this by the observation that when well separated Skyrmions coalesce they seem to create low energy symmetric multi-soliton states.

The procedure that we use is to set up initial conditions with the required topological charge involving a collision between one or more Skyrmions of some particular charge, which can be done using the rational map ansatz, described in section 1.3, using charges for which the rational map is already known. We then evolve the configuration until they visibly coalesce or until the potential energy begins to increase. At this point all the kinetic energy is removed, that is, $`\dot{\varphi }=0`$ at all points on the grid, and the evolution is continued until once again the potential energy rises. This procedure is repeated many times until the energy is no longer decreasing. Once the solution is sufficiently relaxed it often pays to just evolve the solution under the full equations of motion without removing kinetic energy since this can prevent the minute oscillations required to achieve the minimum from gaining momentum.

Even this procedure can be very slow, but can be speeded up by adding a dissipative term to the equations of motion,

$$M\ddot{\varphi }\alpha (\dot{\varphi },_i\varphi ,_i\dot{\varphi },_i_j\varphi )\lambda \varphi =ϵ\dot{\varphi },$$
(2.18)
for $`ϵ>0`$. Clearly, the static solutions are still the same. However, the dissipation causes the solution to roll down the potential well to the minimum much quicker in certain circumstances, with $`ϵ=0.5`$ seeming to work well. We should note that the addition of this dissipation can also effect the Courant instability of the algorithm, and, in particular, very near to the minimum one is plagued by instability. Experience, has shown us that a combination of running with dissipation and then without helps speed up the process.

Obviously, one should be concerned that this process might not necessarily lead one to the global minima. One might, for example, relax down to a metastable local minima or the initial conditions may have some symmetry which is maintained by the equations of motion and hence the final solution. We attempt to ensure that we do not encounter the latter possibility by creating initial conditions using the product ansatz, which is manifestly asymmetric $`(U_1U_2U_2U_1)`$. The possibility of local minima, however, can never be totally excluded, but one can build up confidence in the minima by using different initial configurations.

For low charges ($`B4`$), the attractive channel configurations discussed in ref. are particularly good initial conditions. But for higher charge no such maximally attractive channels exist for $`B`$ well separated Skyrmions and only a small number of attractive configurations are known. We find that sensible initial conditions can be produced for any charge $`B>10`$ by using two clusters, one with charge $`Bn`$ and the other with charge $`n`$, such that $`n=1,2,3,4,5`$ but no larger. These are Lorentz boosted together with a velocity of $`v=0.3`$ in a collision which has a small but non-zero impact parameter.

## 3 Skyrmion identification

In this section we will present the results of an extensive set of simulations performed with the intention of identifying the symmetry and associated polyhedron of the minimum energy configurations for $`B22`$. We should note that making definitive statements as to an identification for a particular charge has, historically, been fraught with difficulties. In particular, the identifications of the $`B=5`$ and $`B=6`$ solutions in ref. and the $`B=9`$ configuration in ref.  were incorrect to varying degrees. Suffice to say when the two numerical methods agree for a wide range of initial conditions and simulation parameters, there is strong grounds to believe that we have created the correct configuration; identifying the symmetry is then the only complication. Conversely, when they disagree, or we know more than one very low energy solution for a particular charge, it is a matter of some debate which is the true minimum energy solution, or whether there is some more complicated symmetry allowing local, or even degenerate minima. We will engage in this debate at various stages, but the reader should make up their own mind as to the strength of our arguments.

We have argued strongly in the previous sections on numerical methods that the simulated annealing of rational maps is the simplest and probably cleanest approach to compute low energy Skyrmion configurations. Clearly, its veracity for computing the true minima depends on the viability of the ansatz for describing Skyrmions, and that the energy functional based on $``$ is a good approximation to the true energy, or at least the relative energies of particular configurations. With these caveats in mind, we present our first attempt at identification of the minima as the rational maps which minimize $``$, before discussing other possibilities.

The results of the simulated annealing algorithm applied to a general rational map of degree $`B(22)`$ and the symmetry identification procedure discussed in section 1.5 are presented in table 1. In each case, we tabulate the identified symmetry group $`G`$, the minimum value of $``$, the quantity $`/B^2`$ — which is strikingly uniform at around 1.2-1.3 — and the value of $`E/B`$ for a profile function which minimizes the energy functional (1.6) for the particular map.

We should first comment that for $`B8`$ the rational maps which minimize $``$ are exactly those presented in ref.  to approximate the results of the full non-linear simulations ; thus the simulated annealing algorithm provides a nice numerical check that the same maps are reproduced by searching the full parameter space of rational maps. Also for $`B7`$ all the symmetry groups with the exception of $`B=9`$, $`B=10`$ and $`B=13`$ are compatible with the fullerene hypothesis: that $`𝒮_B`$ has $`4(B2)`$ trivalent vertices and is constructed from $`2B14`$ hexagons and 12 pentagons. The baryon density isosurfaces for each of the solutions are displayed in fig. 1 along with a model of the associated polyhedron in fig. 2, which confirm that for the most part they are indeed of the fullerene type. The symmetry groups of $`B=9`$, $`B=10`$ and $`B=13`$ all contain the cyclic subgroup $`C_4`$, which is not compatible with them being of the fullerene type, since the associated polyhedron of such a solution must contain either a four-valent bond, or a square. These are the first Skyrme solutions found which do not comply with the Geometric Energy Minimization (GEM) rules suggested in ref. . We shall discuss these solutions in more detail in the subsequent section, but it is gratifying to note that all the other solutions appear qualitatively to comply with our expectations based on the GEM rules.

In table 2, fig. 3 and fig. 4, we present the results of our extensive search for the minimum energy maps with particular symmetries, usually dihedral groups or chosen from the extensive tables of fullerenes presented in ref. , which lend further weight to our conclusions that those presented in table 1 are in fact the minima relative to the energy functional $``$. They do, however, turn up the possibility that in certain cases the minima of $``$ may not necessarily be the minimum of the Skyrme energy, since some of them have values of $``$ very close to the values presented in table 1. When the values are so close it is difficult to make a guess as to how the relaxation to the true solution might effect their relative positions; an issue to which we shall return in subsequent sections. For the moment we will denote them by $``$, and conclude at least that they are not a global minima of $``$, but are believed to represent other critical points.

For $`9B22`$ we shall now describe in detail the rational maps we have obtained, the structure of the associated Skyrmions and make a comparison with the results from full field simulations. Our study of the rational maps has already turned up a few oddities and we shall attempt to interpret these at the relevant charge. Charges where the fullerene hypothesis appears to break down are $`B=9`$ and $`B=13`$, while the rational map approach to representing the minimum energy Skyrmion appears to need careful consideration for $`B=10`$, $`B=14`$, $`B=16`$ and $`B=22.`$

### 3.1 $`B=9`$

In ref.  it was suggested that the $`B=9`$ minimum energy configuration had $`T_d`$ symmetry, a symmetric configuration of $`\mathrm{C}_{28}`$ (it corresponds to configuration 28:2 in ref. ). The polyhedron to which it corresponds comprizes of 12 pentagons, fuzed into 4 triplets, placed at the vertices of a tetrahedron, with 4 hexagons placed at the vertices of a dual tetrahedron. This, plus the solutions for $`B=7`$ and $`B=8`$, was one of the main motivations of the fullerene hypothesis. Unfortunately, it appears that the original identification of the symmetry was incorrect and further relaxation of this configuration using the full non-linear field equations lead to a somewhat different solution.

The minimizing map in this case has $`D_{4d}`$ symmetry and $`=109.3`$, with the functional form of the map being given by

$$R=\frac{z(a+ibz^4+z^8)}{1+ibz^4+az^8},$$
(3.1)
where $`a=3.38,b=11.19.`$ It should be noted that this is slightly lower than the value $`=112.8`$ of the tetrahedral map

$$R=\frac{5i\sqrt{3}z^69z^4+3i\sqrt{3}z^2+1+az^2(z^6i\sqrt{3}z^4z^2+i\sqrt{3})}{z^3(z^63i\sqrt{3}z^4+9z^25i\sqrt{3})+az(i\sqrt{3}z^6+z^4+i\sqrt{3}z^21)},$$
(3.2)
where $`a=1.98.`$

Amazingly, the solution which was created by the relaxation of initially well-separated $`B=8`$ and $`B=1`$ configurations (henceforth such initial conditions will denoted $`8+1`$) using the full non-linear field equations, and confirmed using a number of different initial conditions (for example, $`7+2`$ and $`6+3`$), is precisely that which corresponds to the rational map (3.1). The associated polyhedron is not a fullerene, since the symmetry group $`D_{4d}`$ is incompatible with pentagons and hexagons forming a trivalent polyhedron. In fact, it has two four-valent links which occur between four pentagons forming the top and bottom pseudo-faces444We shall use the term pseudo-face to refer, rather loosely, to a set of connected polygons, which act from the point of view of symmetry of the associated polygon as a single face. of a rather flat polyhedron, linked by a belt of eight alternately up and down pointing pentagons; the top and bottom being rotated relative to each other by $`45^{}`$ to give the $`D_{4d}`$ symmetry. In section 5.2 we will discuss how this solution can be formed by the symmetry enhancement of the other known fullerene corresponding to $`\mathrm{C}_{28}`$ which has $`D_2`$ symmetry (this is labelled 28:1 in ref. ). On the basis of this we conclude that $`𝒮_9`$ has $`D_{4d}`$ symmetry and not $`T_d`$ as previously suggested, but that the known $`T_d`$ symmetric solution is a saddle point.

### 3.2 $`B=10`$

A glance at table 1 and table 2 shows that there are (at least) four maps whose $``$ values are very close together. Minimizing over all maps yields the value $`=132.6`$ and this can be obtained from a $`D_{4d}`$ symmetric map of the form

$$R=\frac{z^2(a+ibz^4+z^8)}{1+ibz^4+az^8},$$
(3.3)
where $`a=8.67`$ and $`b=14.75.`$ This is not a fullerene Skyrmion, the four-fold symmetry this time manifesting itself in the existence of a square on the top and bottom of the associated polyhedron.

There is, however, a fullerene Skyrmion with $`=132.8`$, very close to that of the minimum, which corresponds to a $`D_3`$ symmetric map of the form

$$R=\frac{z(1+az^3+bz^6+cz^9)}{c+bz^3+az^6+z^9},$$
(3.4)
where $`a,b,c`$ are complex parameters. The minimum is obtained for the values $`a=4.401.72i,b=2.38+3.10i,c=0.12+0.19i.`$ The symmetry can be increased to $`D_{3d}`$ by choosing $`b`$ to be real and $`a`$ and $`c`$ to both be purely imaginary, and within this class the minimum is very slightly higher at $`=133.5,`$ which is attained when $`a=20.40i,b=30.22,c=4.69i.`$ Finally, if $`a,b,c`$ are all real then the symmetry is $`D_{3h}`$ and the minimum in this class has $`=143.2`$ when $`a=5.14,b=2.20,c=0.36.`$

The baryon density of the $`D_{4d}`$ symmetric map is presented in fig. 1 and for the other three maps, $`D_3,D_{3d},D_{3h}`$, in fig. 3. The polyhedron associated with the $`D_{4d}`$ solution can be constructed by taking two squares each surrounded by 4 hexagons, connected via a band of 8 pentagons alternately pointing up and down; the two squares being rotated relative to each other by $`45^{}`$. Each of the links is trivalent, but instead of comprizing of 12 pentagons and 6 hexagons as would have been suggested by the fullerene hypothesis, it contains 8 pentagons, 8 hexagons and 2 squares, although the number of vertices $`324(B2)`$ and faces $`182(B1)`$ are still compatible with the GEM rules The other three maps give Skyrmions of fullerene type, with the baryon density isosurface comprising the requisite number of pentagons and hexagons arranged in a trivalent polyhedron. The associated polyhedron for the $`D_{3h}`$ solution (which corresponds to 32:5 in ref. ) comprizes of two copies of a hexagonal triple linked by a belt of 12 pentagons which can be thought of as being made of 3 sets of four fuzed pentagons in a $`C_3`$ arrangement. The $`D_3`$ and $`D_{3d}`$ solutions are very similar: to make each of the associated polyhedra (which correspond to configurations 32:6 and 32:4 in ref.  respectively) first start with two pentagon triples. There are six places on each triple to which one can add another polygon, which fall into two types — one can connect to a single pentagon edge, or between two pentagons connecting to an edge of each. To make the two different configurations, one must add three hexagons and three pentagons alternately around each of the triples; the difference being which polygon (pentagon or hexagon) connects to the two different sites. In particular, the $`D_{3d}`$ configuration has hexagons connected to the single pentagon, and pentagons between the two pentagons; vice versa for the $`D_3`$ configuration. Once one has added these 6 polygons, the two identical copies are then connected, the two being rotated relative to each other by $`60^{}`$ in the $`D_{3d}`$ configuration, and at no particular fixed angle in the case of $`D_3`$. As fullerenes these three structures are tabulated in ref. , along with three other possibilities which have less symmetry ($`2\times C_2`$ and $`D_2`$). Unfortunately, these lower symmetry solutions are impossible to find using the simulated annealing algorithm since their symmetry groups are subgroups of $`D_{4d}`$ and the minimum energy rational map with this symmetry is already known.

In order to try and understand which of the four configurations is the true minimum we have relaxed all the different initial configurations made from two individual Skyrmion clusters whose baryon numbers sum to 10, that is, $`9+1`$, $`8+2`$, $`7+3`$, $`6+4`$ and $`5+5`$. None of these form the $`D_{4d}`$ configuration, nor that with $`D_{3h}`$, suggesting — but not proving — that neither of them are the minimum energy solution, and are hence likely to be saddle points. However, it appears that it is possible to produce both the $`D_{3d}`$ and $`D_3`$ configurations from collisions. In particular, the $`7+3`$ relaxation appears to give the more symmetric $`D_{3d}`$ configuration, while all the others give one which only has $`D_3`$ symmetry. We have already found solutions which we believe to be saddle points, for example, the $`B=9`$ configuration with $`T_d`$ symmetry, but here we have evidence for a new phenomena — local minima. Given that there is no symmetry which can be invoked to explain why there might be degenerate minima, it seems likely that the energies of the two configurations must differ by a minute amount. As we shall discuss in section 4, given the uncertainties we are unable to ascertain which is the global minimum.

The reason that the rational map ansatz is so successful in describing Skyrmions is that they appear to prefer to be as spherical555A sensible quantification of how spherical a solution is might be to consider the eigenvalues of the moments of the baryon density. A distribution which is isotropic, and hence almost spherical, would have all the eigenvalues approximately equal, whereas a more elongated solution would have one which is substantially different to the other two. as possible. Examination of the models of the associated polyhedra sheds some light on the preference of the rational map ansatz for the $`D_{4d}`$ and $`D_3`$ configurations, which are very spherical, as opposed to the more elongated $`D_{3d}`$ configuration. It might be that this oddity is not reproduced in the full non-linear energy functional, and the configuration with $`D_{4d}`$ symmetry is impossible to reproduce.

The phenomenon of many different configurations at a given baryon number, often with energies very close to the minimal value, is a feature of the fullerene hypothesis which one might have been able to predict since there are many possible polyhedra which contain 12 pentagons and $`2B14`$ hexagons for $`B9`$. The possibility of four-valent vertices and also trivalent configurations containing squares only make things worse. This was the main motivation for the choice of a simulated annealing algorithm as our minimization scheme for rational maps; a choice which appears to have been vindicated. It is also one reason why we have used two very different numerical techniques, the rational map approach and full field simulations, to try and confirm the results we obtain, thereby increasing the confidence that the solutions we construct are the global minima. Unfortunately, in this case we were unable to make a definitive identification of the minimum energy Skyrmion, $`𝒮_{10}`$, but have presented evidence that it is one of two configurations which are almost indistinguishable.

### 3.3 $`B=11`$

The minimum value at this charge is $`=161.1`$ and this is obtained from the $`D_{3h}`$ symmetric map

$$R=\frac{z^2(1+az^3+bz^6+cz^9)}{c+bz^3+az^6+z^9},$$
(3.5)
where $`a,b,c`$ are real parameters, taking the values $`a=2.47,b=0.84,c=0.13`$ at the minimum. The associated polyhedron (which corresponds to 36:13 in ref. ) can be constructed by considering a hexagon to which 3 pentagons and 3 hexagons are connected alternately with a $`C_3`$ symmetry. Each of the pentagons is part of a set of four fuzed together, each of which is placed in the $`C_3`$ arrangement. Each of these is then connected to another hexagon, which is directly below the first one and can be thought of as being equivalent to the original one from the point of view of symmetry. The spaces in between are filled up with hexagons, the whole structure comprizing of 12 pentagons and 8 hexagons.

The exact same configuration was produced by the collision and then relaxation of two Lorentz-boosted $`B=3`$ Skyrmions and a stationary $`B=5`$ solution in a linear arrangement ($`3+5+3`$). Given that the fullerene hypothesis is clearly not the whole story for $`B=9`$ and $`B=10`$, it is reassuring that things get back on track at $`B=11`$ with what appears to be the unique global minimum, $`𝒮_{11}`$, being a fullerene type solution describable by the rational map ansatz.

### 3.4 $`B=12`$

Considering all degree 12 maps the minimum is found to be $`=186.6`$ and this can be reproduced from a $`T_d`$ symmetric map constructed as follows. Decomposing $`\underset{¯}{13}`$ as a representation of $`T`$ gives

$$\underset{¯}{13}|_T=2A+A_1+A_2+3F.$$
(3.6)
Now let $`p_\pm `$ be the Klein polynomials

$$p_\pm =z^4\pm 2\sqrt{3}iz^2+1,$$
(3.7)
associated with the vertices and faces of a tetrahedron. On applying the $`C_3`$ generator contained in the tetrahedral group to these polynomials they acquire the multiplying factors $`p_\pm e^{\pm 2\pi i/3}p_\pm .`$ Thus, the degree 12 polynomials $`p_\pm ^3`$ are strictly invariant, forming a basis for the representation $`2A`$ in the above decomposition and the polynomials $`p_+^2p_{}`$ and $`p_+^2p_{}`$ are bases for the representations $`A_1`$ and $`A_2`$ respectively. Explicitly, the rational map

$$R=\frac{ap_+^3+bp_{}^3}{p_+^2p_{}},$$
(3.8)
is $`T_d`$ symmetric for all real $`a`$ and $`b`$, with the minimal value $`=186.6`$ obtained for $`a=0.53,b=0.78.`$ We should note that there are other maps with $`T_d`$ symmetry (the denominator in (3.8) can be replaced by $`p_+^3`$, for example), but it appears that all these have a larger value for $`.`$

As for $`B=11`$, this fullerene-like configuration was reproduced in non-linear field theory simulations, this time from initially well-separated $`B=7`$ and $`B=5`$ solutions ($`7+5`$), allowing us to conclude that it is the unique $`𝒮_{12}`$. The associated polyhedron (which corresponds to configuration 40:40 in ref. ) is in some ways similar to the $`T_d`$ solution at $`B=9`$: there being four pentagon triplets positioned on the vertices of a tetrahedron. Each of these triplets is completely surrounded by hexagons forming a polyhedron well-known in fullerene chemistry , where it is one of 40 configurations with 12 pentagons and 10 hexagons which are candidates for a $`\mathrm{C}_{40}`$ cage.

### 3.5 $`B=13`$

The minimal map of degree 13, deduced from simulated annealing of general maps, has cubic symmetry, another with four-fold symmetry which is incompatible with the fullerene hypothesis. It is interesting to note that the fullerene hypothesis would have predicted a trivalent polyhedron made from 12 pentagons and 12 hexagons. We have not discussed representations of the cubic group666The cubic symmetry group $`O`$ is that of the octahedron/cube, without all the reflection symmetries contained in the full symmetry group $`O_h`$., $`O`$, so we shall describe the group theory of this example by embedding it into the tetrahedral group, whose representations were reviewed earlier.

One finds that

$$\underset{¯}{14}|_T=3E^{}+2E_1^{}+2E_2^{},$$
(3.9)
so there is a two parameter family of $`T`$ maps associated with the first component in (3.9). Setting one of these parameters to zero extends the symmetry to $`O`$ and results in the one parameter family of maps

$$R=\frac{z(a+(6a39)z^4(7a+26)z^8+z^{12})}{1(7a+26)z^4+(6a39)z^8+az^{12}},$$
(3.10)
whose minimum occurs at $`a=0.40+5.18i`$ when $`=216.7.`$ The associated polyhedron is in many ways similar to a cube comprizing of six pseudo-faces, each of which are made of four pentagons with a four-valent bond, very similar to those in the $`B=9`$ configuration. Clearly, in order for the them to fit together, with all the other bonds being trivalent, each of these pseudo-faces must be rotated slightly relative to the one diametrically opposite, which removes the possibility of the reflection symmetries of the cube and, hence, the symmetry group $`O_h`$. The polyhedron comprizes of a total of 24 pentagons, as opposed to the 12 pentagons and 12 hexagons that would have been expected had the fullerene hypothesis been correct for this charge. As we shall discuss in section 5.2, this solution, which is reproduced in relaxation of a range of initially well separated clusters ($`3+7+3`$ and $`12+1`$, for example), is rather special; it being obtainable via a multiple symmetry enhancement of a $`D_2`$ fullerene polyhedron (probably from either configurations 44:75 or 44:89 in ref. ).

If $`a`$ is real then the symmetry can be extended to $`O_h`$, but the minimum in this class is quite a bit higher at $`=265.1`$ for $`a=7.2.`$ This Skyrmion, which is probably a saddle point, has recently been computed in ref.  from the relaxation of a single Skyrmion surrounded by 12 others in an initially face centred cubic array. The polyhedron associated with this configuration is more akin to the octahedron than the cube, comprizing of eight triangular pseudo-faces. It contains mainly four-valent bonds; the only trivalent ones being placed in the centre of the pseudo-faces.

There is a degree 13 map with $`D_{4d}`$ symmetry whose $``$ value is extremely close to the minimal one, in fact $`=216.8.`$ This map is

$$R=\frac{z(ia+bz^4+icz^8+z^{12})}{1+icz^4+bz^8+iaz^{12}}$$
(3.11)
where $`a,b,c`$ are real and take the values $`a=5.15,b=50.46,c=46.31`$ at the minimum. This Skyrmion looks similar to the $`O`$ symmetric minimum, but it only has two four-valent vertices as opposed to the six in the cubic configuration, and can be thought of as being an extension of the $`D_{4d}`$ configuration with $`B=9`$. The two four valent bonds are part of two pseudo-faces forming the top and bottom, which are linked by eight copies of a single hexagon connected to a pentagon, alternately arranged pointing up and down, so that four hexagons and four pentagons connect to both the top and bottom pseudo-faces. This configuration contains 8 hexagons and 16 pentagons, breaks the GEM rules since the number of vertices is 42 rather than the predicted number 44 (the number of faces is still $`242(B1)`$) and can be created by a single symmetry enhancement of the same $`D_2`$ fullerene as the $`O`$ symmetric configuration (see section 5.2 for details). Given the similarity of this configuration to the minimum, it is no surprize that the values of $``$ are very close. We conclude, therefore, that $`𝒮_{13}`$ has $`O`$ symmetry, can be approximated by the rational map (3.10), and that the GEM rules and the fullerene hypothesis breakdown at this charge as they did for $`B=9`$.

### 3.6 $`B=14`$

The minimizing map of degree 14 has only a relatively small symmetry, that of $`D_2.`$ The map can be written in the form

$$R=\left(\underset{j=0}{\overset{7}{}}a_jz^{2j}\right)/\left(\underset{j=0}{\overset{7}{}}a_{7j}z^{2j}\right),$$
(3.12)
where $`a_7=1`$ and $`a_0,..,a_6`$ are complex parameters. The minimum is $`=258.5`$ which occurs when the parameters are those given in table 3. This configuration has much less symmetry than any of the others previously described and the associated polyhedron is difficult to visualize in detail. It can be constructed by arranging the 12 pentagons in 4 sets of 3. The 4 sets should be split up into two pairs, each of which is connected by three hexagons, one in the gap between two pentagons on each side, and the other two either side of the first. The two pairs should then be connected by a band of a further 8 hexagons, making 14 in total. The configuration is of the fullerene type (it corresponds to configuration 48:144 in ref. ) and is one of 192 possibilities containing 12 pentagons and 14 hexagons.

Attempts to reproduce this solution by the relaxation of well-separated clusters have proved unsuccessful. We have tried initial conditions which comprize of $`7+7`$, $`12+2`$ and $`13+1`$ and in each case the same configuration, shown in fig. 5 was the end-product. This configuration has even less symmetry than the minimum energy rational map; it having just $`C_2`$ symmetry. The associated polyhedron is almost impossible to describe due to the lack of symmetry suffice to say that it contains 12 pentagons and 14 hexagons, and corresponds to configuration 48:83777This was identified by computing the pentagon index for the associated polyhedron and checking it against the table in ref. . in ref. .

Note that this Skyrmion is very elongated and so it is not surprising that the rational map approximation does not describe this configuration very well, since it assumes that the baryon density has the same angular distribution on concentric spherical shells. Presumably there is a rational map which describes a distorted, more spherical, version of this Skyrmion, but its $``$ value will be larger than the minimal one. Unfortunately we are unable to find this rational map using simulated annealing since its symmetry group is contained within that of the minimizing map. We believe that $`𝒮_{14}`$ is the elongated configuration shown in fig. 5, based on the fact that it was created easily in the collision of well-separated clusters. It is of the fullerene type, but at this stage we do not have a good description of it in terms of a rational map. In subsequent sections, when we use the rational map ansatz as the starting point for a quantitative investigation of Skyrmion properties we will be forced to use the $`D_2`$ configuration instead of what we believe to be the minimum. However, this should not lead to substantial errors.

### 3.7 $`B=15`$

Considering all degree 15 maps the minimum is found to be $`=296.3`$ which has tetrahedral symmetry888As for the cubic group $`O`$, the group $`T`$ is that of the tetrahedron, but without the reflection symmetries of $`T_d`$., $`T`$. To construct the map the relevant decomposition is

$$\underset{¯}{16}|_T=2E^{}+3E_1^{}+3E_2^{}.$$
(3.13)
At first sight it may appear that there is a one (complex) parameter family of tetrahedral maps corresponding to the $`2E^{}`$ component. However, this is not the case since this family of maps is degenerate, having common factors. From

$$\underset{¯}{8}|_T=2E^{}+E_1^{}+E_2^{},$$
(3.14)
it follows that there is a one parameter family of degree seven tetrahedral maps (this family is constructed explicitly in ). Furthermore,

$$\underset{¯}{9}|_T=A+A_1+A_2+2F,$$
(3.15)
so there is a strictly invariant degree eight tetrahedral polynomial, which is given by $`p_+p_{}=1+14z^4+z^8`$ and is the vertex polynomial of a cube. A basis for the $`2E^{}`$ component of (3.13) is obtained by multiplying each basis polynomial for the $`2E^{}`$ component in (3.14) by $`p_+p_{}`$, and hence the corresponding map is degenerate, being only degree seven rather than 15.

The $`3E_1^{}`$ component in (3.13) does correspond to a genuine two (complex) parameter family of degree 15 tetrahedral maps. Using the methods described in ref.  we find that this family of maps is given by $`R=p/q`$ where

$`p=`$ $`i\sqrt{3}\left(1+ab\right)z^{15}+\left(7799a5b\right)z^{13}`$ (3.16)
$`+`$ $`i\sqrt{3}\left(637+21a+35b\right)z^{11}+\left(1001+561a65b\right)z^9`$
$`+`$ $`i\sqrt{3}\left(429+99a+45b\right)z^7+\left(1001297a127b\right)z^5`$
$``$ $`i\sqrt{3}\left(273+185a+15b\right)z^3+\left(115+27a+5b\right)z,`$
and $`q(z)=z^{15}p(1/z)`$. The value $`=296.3`$ is obtained when $`a=0.16+2.06i`$, $`b=4.478.57i.`$ If $`a`$ and $`b`$ are both real then the symmetry extends to $`T_d`$, but the minimum in this class is higher at $`=313.7,`$ when $`a=4.64,b=20.45.`$

The polyhedron associated with the $`T`$ symmetric solution is of the fullerene type. It contains 12 pentagons and 16 hexagons which can be thought of as being arranged in 8 pseudo faces: 4 of these comprize of hexagon triples, whereas the other 4 can be made from a hexagon connected to 3 pentagons in a $`C_3`$ arrangement. The 4 hexagon triples can be thought of as being placed on the vertices of a tetrahedron, and the other 4 pseudo-faces, which are connected to the others, can be thought of as being on the vertices of another tetrahedron which is not dual to the first, removing the possibility of reflection symmetries.

There is no fullerene polyhedron with $`T_d`$ symmetry , and hence this configuration must contain four-valent bonds. In fact it contains more four valent bonds than trivalent ones, in an essentially similar way to the $`O_h`$ configuration for $`B=13`$. The $`T`$ symmetric solution was reproduced in the relaxation of clusters containing $`12+3`$ and $`13+2`$, initially well separated, and hence we conclude that $`𝒮_{15}`$ has $`T`$ symmetry, is of the fullerene type and can be reproduced by a rational map.

### 3.8 $`B=16`$

$`B=16`$ is another interesting situation where it appears that the minimum energy rational map may not in fact be the minimum energy Skyrmion. The minimizing map of degree 16 has $`D_3`$ symmetry which takes the form

$$R=\left(\underset{j=0}{\overset{5}{}}a_jz^{3j+1}\right)/\left(\underset{j=0}{\overset{5}{}}a_{5j}z^{3j}\right),$$
(3.17)
where $`a_5=1`$ and $`a_0,..,a_4`$ are complex parameters. The minimum is $`=332.9`$ which is attained when the parameters take the values presented in table 4. The associated polyhedron can be constructed by first taking two sets of 3 hexagons, each of which is almost flat. Now connect a total of 6 pentagons and 3 hexagons, in 3-fold cyclic order around the flat structure, such that each of the gaps between the original hexagons is filled by a pentagon flanked on one side by another pentagon and on the other side by a hexagon. A further 6 pentagons, split into 3 pairs are used to fuze the two halves; each of the pairs connecting sets of 2 pentagons. Within the structure, there are a number of symmetrically placed groupings of polygons which are exactly those which were found to lead to symmetry enhancement as observed for $`B=9`$ and $`B=13`$. However, in this case a close examination of the baryon density isosurface shows that all the bonds remain trivalent and the associated polyhedron is of the fullerene type, containing 12 pentagons and 18 hexagons.

There exists a family of $`D_2`$ symmetric rational maps with $`B=16`$ of the form

$$R=\left(\underset{j=0}{\overset{8}{}}a_jz^{2j}\right)/\left(\underset{j=0}{\overset{8}{}}a_{8j}z^{2j}\right),$$
(3.18)
where $`a_8=1`$ and $`a_0,..,a_7`$ are complex parameters. The minimum in this class takes place when $`=333.4`$, very close to that with $`D_3`$, and the parameters take the values presented in table 5. Since the solution has very little symmetry it is difficult to describe the associated polyhedron as for that with $`B=14`$. It is of the fullerene type, comprizing of 12 pentagons and 18 hexagons, and can be thought of being formed from two identical half shells, connected together. To construct each of the two shells, start with a hexagon and attach to it, in cyclic order, two pentagons, a hexagon, another two pentagons and another hexagon. Then connect a pentagon, pointing downwards to the edge of the two hexagons, and add three hexagons to each side connecting the two pentagons.

Given that there are at least two rational maps with very similar values of $``$, it is interesting to see if we can create them both via the relaxation of initial well-separated clusters. To this end we have tried a number of different initial conditions ($`7+9`$, $`12+4`$ and $`13+3`$). In contrast to the $`B=10`$ case where we were able to create both the $`D_3`$ and $`D_{3d}`$ configurations, in each case the end-product of the relaxation process had $`D_2`$ symmetry. This strongly suggests, but does not prove, that this is the global minimum energy solution and that the $`D_3`$ configuration may be a saddle point solution. This is interesting since it suggests that the relative ordering of the $`D_3`$ and $`D_2`$ solutions is probably different when considering the full non-linear energy functional and that for the rational map ansatz.

### 3.9 $`B=17`$

The case of $`B=17`$ is interesting since it was conjectured in ref.  on the basis of the fullerene hypothesis that there might exist a Skyrmion configuration of this charge with the same structure as that of Buckminsterfullerene, $`\mathrm{C}_{60}`$, which comprizes of 20 hexagons and 12 pentagons in an icosahedral configuration. This is well known as the standard design for a football since it is almost spherical, and also in civil engineering where it was championed by Buckminster Fuller as a candidate for a geodesic dome. It is the isolated pentagon structure (each pentagon is isolated by connecting it to 5 hexagons) with the lowest number of vertices, and appears to be the most stable carbon cage; it being the subject of much interest in chemistry in the recent past.

An icosahedrally symmetric rational map was found in ref.  and this exact same map is reproduced by minimizing over all degree 17 maps. The value $`=363.4`$ is given by the $`Y_h`$ symmetric buckyball map999In ref. the value of $``$ quoted for the $`B=17`$ buckyball map was the result of a typographical error and should read $`363.41`$, not $`367.41.`$

$$R=\frac{17z^{15}187z^{10}+119z^51}{z^2(z^{15}+119z^{10}+187z^5+17)}.$$
(3.19)
Further confirmation that this is indeed the minimum energy configuration at this charge comes from non-linear field theory simulations. We have performed relaxations of initial configurations of $`12+5`$, $`13+4`$ and $`5+7+5`$ all of which relax very quickly to the buckyball structure.

In addition to the buckyball map there is also a non-fullerene map which has a low value of $`.`$ This map has $`O_h`$ symmetry, and again we shall describe its group theory construction by embedding it into the tetrahedral group. The tetrahedral decomposition in this case reads

$$\underset{¯}{18}|_T=3E^{}+3E_1^{}+3E_2^{},$$
(3.20)
with the two-parameter family of $`T`$ maps corresponding to the $`3E^{}`$ component, reducing to a one-parameter family of $`O`$ maps by setting one of the two parameters to zero. If the remaining parameter is real then the symmetry extends to $`O_h`$ and the map is given by $`R=p/q`$ where

$$p=(129+a)+(2380+116a)z^4+(24310+286a)z^8+(6188156a)z^{12}+(17+9a)z^{16},$$
(3.21)
and $`q(z)=z^{17}p(1/z).`$ The minimal $`O_h`$ map in this class has $`=367.2`$ when $`a=280.9.`$

Given that Buckyball map is the minimum map with $`B=17`$ and it is reproduced in the numerical field theory relaxations, we conclude that the fullerene hypothesis and the conjecture of ref.  are spectacularly confirmed at this charge; $`𝒮_{17}`$ has $`Y_h`$ symmetry and the associated polyhedron is the buckyball.

### 3.10 $`B=18`$

After the particularly high symmetry of the $`B=17`$ solution, the minimizing map of degree 18 is relatively unremarkable having only $`D_2`$ symmetry, and takes the form

$$R=\left(\underset{j=0}{\overset{9}{}}a_jz^{2j}\right)/\left(\underset{j=0}{\overset{9}{}}a_{9j}z^{2j}\right),$$
(3.22)
where $`a_9=1`$ and $`a_0,..,a_8`$ are complex parameters. The minimum value is $`=418.7`$ which is attained when the parameters take the values given in table 6. The associated polyhedron, which is of fullerene type, but is not an isolated pentagon structure101010In fact, no isolated pentagon structures exist with 12 pentagons and 22 hexagons ., is difficult to describe to similar reasons to the $`B=14`$ and $`B=16`$ solutions with $`D_2`$ symmetry. It contains 12 pentagons and 22 hexagons, and can best be described, as for $`B=16`$, in terms of two half shells which fit together to create the whole polyhedron. Each shell can be created by first taking a hexagon and connecting to it, in cyclic order, two hexagons, a pentagon, two hexagons and another pentagon. Now connect a pentagon, pointing downward to each of the hexagons, then fill in the gaps, of which there are 6, with hexagons.

Given the rather low symmetry of the minimum energy rational map, one might think that as, for example, with the cases of $`B=14`$ and $`B=16`$ that there might be some confusion. However, by relaxing initial conditions comprizing of $`9+9`$ and $`17+1`$ we have reproduced the $`D_2`$ configuration. Hence, we conclude that $`𝒮_{18}`$ is of the fullerene type and can be well approximated by the $`D_2`$ rational map above.

### 3.11 $`B=19`$

The minimum value at degree 19 is $`=467.9`$ and is attained by a $`D_3`$ symmetric map of the form

$$R=\left(\underset{j=0}{\overset{6}{}}a_jz^{3j+1}\right)/\left(\underset{j=0}{\overset{6}{}}a_{6j}z^{3j}\right),$$
(3.23)
where $`a_6=1`$ and $`a_0,..,a_5`$ are complex parameters given in table 7. The associated polyhedron can be constructed by taking two sets of three hexagons which are almost flat. To each, connect a pentagon in the gap connecting two hexagons and another pentagon next to it connected to only one of the hexagons in the triple. These two structures are $`C_3`$ symmetric and can be thought of as forming the top and bottom of the polyhedron. They are connected together by 9 sets of two hexagons, one which connects to the top and the other which connects to the bottom. This is a fullerene type polyhedron containing 12 pentagons and 24 hexagons.

There is a more symmetric map, with $`T_h`$ symmetry, whose $``$ value is only slightly higher than the minimum at $`=469.8.`$ Computing the relevant decomposition

$$\underset{¯}{20}|_T=4E^{}+3E_1^{}+3E_2^{},$$
(3.24)
shows that there is a three parameter family of tetrahedral maps of degree 19 corresponding to the first component in (3.24). This family of maps is given by $`R=p/q`$ where

$`p=`$ $`(2399b)z+(503a25c)z^3+(5508+460b)z^5+(1300a+284c)z^7`$ (3.25)
$`+(4862286b)z^9+(1794a+210c)z^{11}+(9996196b)z^{13}`$
$`+(484a+44c)z^{15}+(135+31b)z^{17}+(ac)z^{19},`$
and $`q(z)=z^{19}p(1/z).`$ The minimum in this family is $`T_h`$ symmetric with $`a=5.5`$, $`b=6.3`$, $`c=37.3`$ and produces the $``$ value given above. The associated polyhedron contains many four valent bonds and is most definitely not of the fullerene type. In fact a cursory glance at the configuration might convince one that the polyhedron has cubic symmetry with there being eight sets of three fuzed pentagons effectively situated at the corners of a cube. However, the two pentagons which are situated at the centre of each face break the cubic symmetry since they point alternately in different directions. In total the configuration comprizes of a total of 36 pentagons.

The $`D_3`$ configuration was reproduced in the collision and subsequent relaxation of $`17+2`$ and therefore we conclude that it corresponds to $`𝒮_{19}`$. It is of the fullerene type and can be approximated using the rational map (3.23).

### 3.12 $`B=20`$

For $`B=20`$ the minimum value is $`=519.6`$ and is reproduced by the $`D_{6d}`$ symmetric map

$$R=\frac{z^2(ia+bz^6+icz^{12}+z^{18})}{1+icz^6+bz^{12}+iaz^{18}},$$
(3.26)
with $`a=16.8`$,$`b=288.3`$ and $`c=215.8.`$ The associated polyhedron can be thought of as being created from two half shells connected together and rotated relative to each other by $`30^{}`$ to give $`D_{6d}`$ symmetry. Each half can be constructed from a single hexagon, surrounded by another six forming an almost flat hexagonal structure, which are then surrounded by 6 hexagons and 6 pentagons. The flat structure has 12 positions for attaching another polygon, 6 places to connect to one hexagon and 6 places to connect to two. The pentagons connect to two and the hexagons connect to one forming a structure which is of the isolated pentagon fullerene type (corresponding to configuration 72:1 in ref. ). It contains 12 pentagons and 26 hexagons, was reproduced in the relaxation of $`17+3`$ and, hence, we conclude that it is $`𝒮_{20}`$.

### 3.13 $`B=21`$

The $`B=21`$ minimizing map is $`T_d`$ symmetric. From the decomposition

$$\underset{¯}{22}|_T=3E^{}+4E_1^{}+4E_2^{}$$
(3.27)
there is a three parameter family of $`T`$ symmetric maps corresponding to the $`4E_1^{}`$ component and this family of maps is given by $`R=p/q`$ where

$`p=`$ $`1025+3a+b+c+i\left(210+890a+74b10c\right)\sqrt{3}z^2`$ (3.28)
$`+\left(5985+6327a+1433b75c\right)z^4`$
$`+i\left(542647752a680b+392c\right)\sqrt{3}z^6`$
$`+\left(203490+5814a1598b+690c\right)z^8`$
$`+i\left(352716+16796a+2652b+260c\right)\sqrt{3}z^{10}`$
$`+\left(29393025194a+442b+130c\right)z^{12}`$
$`+i\left(1162807752a1768b120c\right)\sqrt{3}z^{14}`$
$`+\left(20349+14535a+221b243c\right)z^{16}`$
$`+i\left(1330646a+234b10c\right)\sqrt{3}z^{18}`$
$`+\left(21+51a+13b+9c\right)z^{20},`$
and $`q(z)=z^{21}p(1/z).`$ The minimum of $`=569.9`$ is obtained when $`a=20.8`$, $`b=102.0`$, and $`c=570.1`$, for which the symmetry extends to $`T_d`$ since $`a,b,c`$ are all real. The associated polyhedron can be thought of in terms of four copies of two different pseudo faces, one set is placed on the vertices of a tetrahedron and the other on the vertices of a tetrahedron dual to the first. In this respect it is very similar to the $`T_d`$ configuration with $`B=9`$, and different to the $`T`$ configuration with $`B=15`$. One set of pseudo faces comprize of a hexagon triple, whereas the others consists of a hexagon surrounded alternately by hexagons and pentagons.

Note that this map is the latest in an infinite family of tetrahedral maps, corresponding to charges $`B=6n+3`$, where $`n=0,1,2,3,\mathrm{}`$. This is because

$$\underset{¯}{6n+4}|_T=nE^{}+(n+1)E_1^{}+(n+1)E_2^{}$$
(3.29)
so there is an $`n`$ parameter family of tetrahedral maps corresponding to the middle component in the above. For $`n=0,2,3`$ $`(B=3,15,21)`$ we have seen that this family includes the minimal map, and for $`n=1`$ $`(B=9)`$ this family includes a map which is very close to the minimal value. Thus it seems possible that other members of this family will be minimal maps, for example, for $`B=27`$, although this configuration must have only $`T`$ symmetry if it is to be of the fullerene type and will therefore be similar to that with $`B=15`$.

This solution was reproduced in the relaxation of $`17+4`$ and $`20+1`$, and hence we conclude that it corresponds to $`𝒮_{21}`$. It is of the isolated pentagon fullerene type (corresponding to configuration 76:2 in ref. ), comprizing of 12 pentagons and 28 hexagons and can be reproduced by the rational map (3.28).

### 3.14 $`B=22`$

The minimum value at degree 22 is $`=621.6`$, obtained from a $`D_{5d}`$ symmetric map

$$R=\frac{az^2+ibz^7+cz^{12}+idz^{17}+z^{22}}{1+idz^5+cz^{10}+ibz^{15}+az^{20}},$$
(3.30)
where $`a=24.8,b=814.6,c=2000.3,d=320.3.`$

The polyhedron associated with this configuration can be constructed in two halves, which fit together as with many of the solutions already described. To construct each half, take a pentagon and surround it by 5 hexagons. There are 10 places to position another polygon, 5 of which connect to two hexagons and 5 which connect to just one. Place 5 pentagons in the gaps connecting to two hexagons and 5 hexagons just connecting to one. Then place a further 5 hexagons, connecting to the pentagons. This configuration, which comprizes of 12 pentagons and 30 hexagons, is of the isolated pentagon fullerene type (it corresponds to configuration 80:1 in ref. ).

For $`B=22`$ we find that there is an interesting phenomenon in that an icosahedrally symmetric fullerene polyhedron exists , but no corresponding rational map generated Skyrmion. The $`Y_h`$ symmetric $`C_{80}`$ fullerene is constructed in a similar manner to the $`D_{5d}`$ fullerene described above, except that one interchanges the pentagons and hexagons at the point at which there was a choice in inserting polygons into the 10 positions. It is easy to check that there are no $`Y`$ symmetric rational maps of degree 22, since $`\underset{¯}{23}|_Y`$ contains only representations of dimension three and higher. Thus there are symmetric fullerene polyhedra which do not correspond to symmetric rational maps.

This apparent puzzle can be understood111111We thank Conor Houghton for pointing this out to us. by realizing that there are rational map generated Skyrmions whose baryon (and energy) density has more symmetry than the Skyrme field itself. As mentioned earlier, the baryon density of a Skyrmion is localized around the edges of a polyhedron with the face centres of the polyhedron given by the vanishing of the derivative of the rational map, or more accurately by the zeros of the Wronskian of the numerator and the denominator

$$w(z)=p^{}(z)q(z)q^{}(z)p(z),$$
(3.31)
which is in general a degree $`2B2`$ polynomial in $`z.`$ For $`B=22`$ the Wronskian is, therefore, a degree 42 polynomial, and although there are no $`Y`$ symmetric degree 22 rational maps there is a $`Y`$ symmetric degree 42 polynomial, given by the product of the Klein polynomials corresponding to the edges and vertices of an icosahedron . Therefore, it appears that the existence of a symmetric fullerene polyhedron coincides with the existence of a rational map whose Wronskian has this symmetry, but that the existence of such a Wronskian does not imply the existence of a symmetric rational map itself. Fortunately, we have not encountered this situation in our study of minimal energy rational maps and Skyrmions for the other charges we have studied, since it would make the problem of identifying and constructing a particular rational map a much more difficult exercise.

There exists a family of $`D_3`$ symmetric rational maps with $`B=22`$ of the form

$$R=\left(\underset{j=0}{\overset{7}{}}a_jz^{3j+1}\right)/\left(\underset{j=0}{\overset{7}{}}a_{7j}z^{3j}\right),$$
(3.32)
where $`a_7=1`$ and $`a_0,..,a_6`$ are complex parameters. The minimum value of $``$ in this class of maps is $`=623.4`$, which is very close to that for the $`D_{5d}`$. The coefficients at this minimum are presented in table 8. This configuration, which is also of the isolated pentagon fullerene type, corresponds to configuration 80:4 in ref. . It can be constructed by first taking two hexagon triples each of which are surrounded by 3 pentagons and 6 hexagons in a $`C_3`$ arrangement, hexagons filling the gaps which connect to 2 of the hexagons in the original triples. Then connect 3 more pentagons to each ‘half’ in between each of the $`C_3`$ symmetric hexagon triples. The two ‘halves’ should then be connected by a band of 12 hexagons around the centre, which is split up into 3 lots of 4 by the $`C_3`$ symmetry. The whole polyhedron contains 12 pentagons and 30 hexagons, as for the $`D_{5d}`$ configuration.

Relaxation of clusters containing $`17+5`$, $`20+2`$ and $`21+1`$ all lead to the same structure, that with $`D_3`$ symmetry. Therefore, we conclude that $`𝒮_{22}`$ is that approximated by the rational map (3.32). This removes, from the point of view of this paper at least, further motivation for attempting to create the Skyrmion with the $`Y`$ symmetric baryon density isosurface.

### 3.15 Summary

The main conclusion of this section on Skyrmion identification is that the fullerene hypothesis appears to apply for a wide range of charges, and that the rational map ansatz can be used to make a good approximation to the solutions. In particular, we have concluded121212We should note that we have also confirmed the results of refs.  for $`1B8`$. that $`𝒮_B`$ is of the fullerene type for $`7B22`$, except when $`B=9`$ and $`B=13`$. For these charges the associated polyhedron contains four-valent bonds, but as we shall discuss in section 5.2, even these solutions can be related to fullerene polyhedra via symmetry enhancement. Clearly, there is a strong correlation between the structure of multi-Skyrmions and that of fullerene polyhedra.

For $`B=17`$, $`B=20`$, $`B=21`$ and $`B=22`$ where there are fullerene polyhedra in which all the pentagons are surrounded by just hexagons, this type of configuration is picked out as that with minimum energy. In fullerene chemistry the isolated pentagon isomers are thought to minimize energy by placing the pentagonal defects as far as possible from each other, and it likely that this also taking place here. It is interesting to speculate that these highly spherical, isolated pentagon configurations will be the minima at higher charge in a similar way to our suggestion of the fullerene hypothesis just on the basis of the $`B=7`$, $`B=8`$ minima and the $`B=9`$ saddle point. Although there will no doubt be caveats, as we have reported in this paper for the fullerene hypothesis for $`B=9`$ and $`13`$, we believe they are likely to be the exception rather than the rule.

Notwithstanding these successes we have turned up a few oddities. Firstly, we have seen that in the cases of $`B=10`$, $`B=14`$, $`B=16`$ and $`B=22`$, that the minimum energy rational maps might not be $`𝒮_B`$. For $`B=10`$ the situation is particularly interesting since the minimum energy rational map is not of the fullerene type, yet only fullerenes are produced in the relaxation of initially well-separated clusters. While for $`B=14`$ we have been unable to produce a rational map which accurately reproduces the configuration which is readily produced in the full relaxation, probably since it is elongated and has little symmetry. The cases of $`B=16`$ and $`B=22`$ are probably more mundane since the differences in the $``$ values for the two fullerene-like configurations are very small, and it is not difficult to imagine that using $``$ as the energy functional rather than the true energy manages to swap around the two configurations. These special cases vindicate our approach of using two different methods to make our identifications.

## 4 Skyrmion energies

It is usual when presenting numerically generated minimum energy Skyrmion configurations to discuss their energy. Before we present what we believe are good representations of the energies, we should just make a point that identifying the symmetry is probably a better way of judging success, rather than just this single number. In particular, making comparisons between different approaches for computing the energy presented in the literature is extremely hazardous, whereas the symmetry identification should be universal.

The approach that we shall discuss here is based on using the full non-linear dynamics to relax a solution generated using the rational map ansatz for a given charge and symmetry. We have already seen that the rational map ansatz systematically over estimates the ratio of the energy to baryon number ($`E/B)`$ for a given configuration and it is clear that the numerical relaxation is likely to reduce this somewhat.

Specifically, we have performed relaxations on numerical grids for all the solutions listed in table 1 and table 2 with (a) $`\mathrm{\Delta }x=0.02`$ and $`N=100`$ and (b) $`\mathrm{\Delta }x=0.01`$ and $`N=200`$, both using Dirichlet (fixed) boundary conditions at the edge of the grid. These two discretized grids have exactly the same spatial extent and so when computing the initial profile function for the rational map ansatz we have set the profile function to zero at $`r_{\mathrm{}}=10`$. The results of this extensive set of simulations are presented in table 9 for the solutions which are the minimum energy rational maps with respect to the energy functional $``$ and in table 10 for the other critical points of $``$, some of which we have concluded in the previous section are, in fact, the minimum energy Skyrmions. The total simulation time — which is the number of timesteps multiplied by $`\mathrm{\Delta }t`$ — for $`N=100`$ is twice that for $`N=200`$, although making an exact comparison can be misleading since the rate at which the relaxation takes place is somewhat arbitrary due to the periodic removal of kinetic energy. Suffice to say, in both cases we believe that we have run the code for long enough for it to settle down to the minimum.

The first thing to notice is that the computed values of the baryon number for the discrete grid $`B_{\mathrm{dis}}`$ are less than the relevant integer; the value for $`N=200`$ being closer than that for $`N=100`$. There are two possible sources for this error: the first is numerical discretization error from the computation of the spatial derivatives (and the resulting numerical integration) and the other is that the physical size of the box $`r_{\mathrm{}}`$ is not large enough to encompass the whole solution and we have underestimated the gradient energy.

In order to understand which is the dominant effect we have repeated the same relaxations with $`r_{\mathrm{}}=20`$ ($`\mathrm{\Delta }x=0.02`$ and $`N=200`$) for $`B=14,9,13,17`$ and the results are presented in table 11 for a total simulation time which is exactly the same as that for $`\mathrm{\Delta }x=0.02`$, $`N=100`$ and $`r_{\mathrm{}}=10`$. The values of $`B_{\mathrm{dis}}`$ agree to the third decimal place for $`B=14`$ with the agreement being less good for larger values of $`B`$. Clearly discretization error is the dominant effect for the small values of $`B`$, and the truncation error due to $`r_{\mathrm{}}`$ being finite becomes important as it increases.

Remarkably, the value of $`E_{\mathrm{dis}}/B_{\mathrm{dis}}`$ appears to remain constant at the level of around $`\pm 0.001`$ for all the relaxations, irrespective of what is the main source of the error in computing the individual values for the sensible parameters $`\mathrm{\Delta }x`$ and $`N`$ that we have used. Therefore, the main conclusion we will draw from these simulations is the value of $`E_{\mathrm{dis}}/B_{\mathrm{dis}}`$ which we will equate with the true value of $`E/B`$. Clearly, knowledge of $`E/B`$ to a certain level of accuracy allows one to compute the energy of the solution, $`E_B`$, to the same level of accuracy, and the values of the energies based upon this hypothesis are presented in table 12 and table 13 for the minima of $``$ and the other critical points respectively. For the subsequent discussions we will take the value of $`E/B`$ to be that computed when $`\mathrm{\Delta }x=0.01`$, $`N=200`$ and $`r_{\mathrm{}}=10`$ subject to an assumed error of $`\pm 0.001`$, with the relative difference between adjacent values of $`B`$ probably being even more accurate. This error budget is used to include the many other possible systematic uncertainties in computing $`E/B`$ which have not already been discussed and the spread of the values computed.

The precise values of $`E/B`$ we have computed here differ from the those quoted for $`B=19`$ in our earlier work 131313Note that the $`B=9`$ solution of ref.  had $`T_d`$ symmetry and so one should be careful to make the correct comparison, and for $`B=14`$ in recent work which used a simulated annealing algorithm on the full field dynamics . By making comparison with ref.  at the level of accuracy suggested there ($`1\%`$), we see that only the $`B=2`$ and $`B=3`$ are discrepant and then the difference is only of the order of an extra $`1\%`$, but that the actual values quoted are very different at the level of 3 decimal places. We believe the earlier method has a tendency to under estimate the true value of $`E/B`$ since the value of the field at the boundary, which is kept fixed during the relaxation runs, is sensitive to the initial conditions, which are well-separated in contrast to the situation here. This is borne out on inspection of the $`E/B`$ values computed for the relaxed solutions from well-separated clusters for large $`B`$ used to identify the minima. The comparison with the results of ref.  is less good, their results being systematically about $`1\%`$ higher than those presented here and $`2\%`$ higher than those of ref. . Although it is difficult to make strong conclusions, we believe that these higher values could well be a function of using periodic boundary conditions rather than fixed. This assumption was used to make their computed values of $`B`$ almost exactly an integer, but it is clear that such an assumption will modify the scale of the solution. In particular, their quoted value for the energy of a single Skyrmion is larger than the known value $`(E=1.232)`$, and the solution assuming periodic boundary conditions with a domain of finite extent will be somewhat different.

We should comment on the computed values of $`E/B`$ for values of $`B`$ where we have more than one rational map which is of low energy with respect to $``$. Previously in a number of cases ($`B=10,16`$ and $`22`$) we had concluded that the minimum energy rational map is not necessarily the minimum energy Skyrmion and one might hope that the relaxation process might confirm this with their energies being separated by a significant amount.

For $`B=9,15,17`$ and $`19`$ its clear that the computed values of $`E/B`$ are systematically much higher for the solutions which are not minima of the energy functional $``$ and therefore we conclude that these solutions are definitely not the minima of the full Skyrme energy functional (remember that on the basis of relaxation of well-separated clusters we have concluded that the minima with respect to $``$ are in fact the minima). For $`B=13`$ the solution with $`D_{4d}`$ symmetry has exactly the same value as that with $`O`$, while that with $`O_h`$ is much higher. It might seem remarkable that the values for the $`D_{4d}`$ and $`O`$ symmetric solutions are exactly the same, but one should remember that the two solutions are very similar and can be related by symmetry enhancement (see section 5.2). Since the $`O`$ solution was produced from the relaxation of clusters and the $``$ values for the two solutions are very close anyway, we conclude that the $`O`$ solution is probably the minimum. For $`B=10,16`$ and $`22`$, we are unable to tell the different candidate minima apart based on the computed energies since our quoted error of $`\pm 0.001`$ for each of the values of $`E/B`$ encompasses the different solutions under consideration. For $`B=10`$ the $`D_{3h}`$ solution is of higher energy and is clearly not the minima, but the other 3 are well within the range of uncertainty. This is also the case for $`B=16`$ and $`B=22`$, where each of the candidate minima are within the quoted range for $`E/B`$.

In table 14 we have summarized the computed values of $`E/B`$ and $`E_B`$ for the solutions which we have identified as the minima in section 3. Included also are the ionization energy $`I_B=E_{B1}+E_1E_B`$, the energy required to remove a single Skyrmion, and the binding energy per baryon given by $`\mathrm{\Delta }E/B=E_1(E/B)`$, which is the energy required to separate the solution up into single Skyrmions divided by the total baryon number. The accuracy of $`\mathrm{\Delta }E/B`$ will be exactly that of the computed value of $`E/B`$ since the value of $`E_1`$ we compute by this method appears to be exact within the quoted limits, but the errors in computing $`I_B`$ could theoretically be larger since it is the difference of two energies. For $`B>10`$ the worst case errors in computing $`I_B`$ could be as much as $`\pm 0.02`$ (a significant amount on inspection of the quoted values), but since we have already commented that we believe the difference in energies for adjacent values of $`B`$ will be even more accurate than the absolute errors in the energy we suspect that things will be much better. We will comment on this in subsequent sections.

The values of $`E/B`$ computed for these relaxed solutions and also for the original rational map are plotted against $`B`$ in fig. 6. Both start at approximately the same value for $`B=1`$, and both appear to asymptote for large $`B`$, albeit at different values. For the relaxed solutions this appears to be about $`67\%`$ above the Faddeev-Bogomolny bound, which is compatible with that computed for the hexagonal Skyrme lattice , which can be thought of as the infinite limit of a shell-like Skyrmion, while the asymptote for the rational map ansatz is higher at around $`9\%`$. The curve for the relaxed solutions is smoother than that for the rational maps, which has notable dips associated with the highly symmetric solutions with $`B=4,7,13`$ and $`17`$. Although these deviations from what appears to be approximately a smooth curve do not totally disappear after relaxation, one can deduce that the other solutions, not being particularly spherical, do not fit the rational map ansatz as well, but that the relaxation using the full non-linear field dynamics softens these effects.

The binding energy per baryon, $`\mathrm{\Delta }E/B`$, is plotted against $`B`$ in fig. 7, its shape just being the inversion of fig. 6. Interestingly, it increases to an asymptote just as one might expect in a simple model of nuclei which excludes the Coulomb interaction within the nucleus. We shall return to this issue in section 5.6.

The ionization energy, $`I_B`$, is plotted against $`B`$ in fig. 8. We have already commented that our quoted errors in $`E/B`$ might lead to substantial errors in computing $`I_B`$, but that systematic errors in computing the energy of a particular solution are likely to be similar and therefore our computed values for $`I_B`$ could probably be more accurate than one might naively expect. This is borne out by a cursory inspection of fig. 8: we see that the most stable solutions with respect to the removal of a single Skyrmion are those with the most symmetry $`B=4,7`$ and $`17`$, while those with the least symmetry $`B=5,8,14`$ and $`18`$ are much less stable. This is very much as one might expect.

## 5 Discussion

There are some remarkable aspects of the solutions which we have created. In this section we point out, discuss and attempt to explain some of them.

### 5.1 Platonic symmetries

It had been known for sometime that Skyrmion solutions existed whose associated polyhedra are platonic solids ($`B=3,4`$ and 7) and it had been conjectured that the minimum energy solution with $`B=17`$ had $`Y_h`$ symmetry. These platonic symmetry groups are the groups with the highest symmetry (most generators) which are discrete subgroups of $`O(3)`$, the dihedral groups having considerably fewer generators. Here, we have shown that platonic symmetries are even more prevalent in Skyrmion solutions, be they minimum energy solutions $`B=12,13,15`$ and $`21`$, or low energy saddle points $`B=9,13,15,17`$ and $`19`$. The existence of these solutions is truly remarkable.

This could have been expected within the fullerene hypothesis since the platonic groups $`T`$ and $`Y`$ are compatible with the associated polyhedron comprizing of pentagons and hexagons. But we have also seen that the existence of a highly symmetric fullerene polyhedron compatible with the fullerene hypothesis at a particular charge does not necessarily imply that it is a minimum energy configuration. In particular, tetrahedral fullerene structures are compatible with $`B=16`$ and $`19`$, and the group theory decomposition is consistent with the existence of appropriate rational maps. However, in each of these cases the minimum energy Skyrmion has much less symmetry. Clearly, the solution being highly symmetric is not the sole criterion in minimizing the Skyrme energy functional. We have also encountered one case, icosahedral symmetry for $`B=22`$, in which a platonic fullerene polyhedron exists, but no corresponding platonic rational map (and probably no Skyrmion either). However, it appears that a rational map exists (and hence a corresponding Skyrmion) for which the baryon density surface has more symmetry than the Skyrme field and is icosahedrally symmetric. But again this very symmetric structure appears not to be the minimum energy Skyrmion.

For $`B>7`$ the octahedral group is incompatible with the fullerene hypothesis since it requires four-fold symmetry which is impossible in a polyhedron comprizing only of pentagons and hexagons. However, the minimum energy Skyrmion with $`B=13`$ has $`O`$ symmetry and we have also been able to find a Skyrmion with $`O_h`$ symmetry with $`B=17`$ which is a low-energy saddle point. Many such solutions in which there is a four-valent bond connecting 4 pentagons can be related to a fullerene polyhedron by symmetry enhancement as we shall discuss in the next section. This produces an extra twist to the fullerene hypothesis which allows octahedral Skyrmion solutions.

### 5.2 Symmetry enhancement

We have seen that, in keeping with the fullerene hypothesis and the expectations of ref. , most links in the polyhedra associated with Skyrmion solutions are trivalent. In particular, the baryon density isosurface of what we have identified as the minimum energy solutions consist of a trivalent polyhedron for all cases except $`B=9`$ and $`B=13.`$ In these two cases the polyhedra contain four-valent vertices which means that they are not fullerene Skyrmions, since by definition all links are trivalent. However, it turns out that these two exceptional cases can be obtained from fullerenes by the application of a simple rule, which we refer to as symmetry enhancement, and shall now explain.

Consider part of a fullerene which has the form shown in the first figure in fig. 9, consisting of two pentagons and two hexagons with a $`C_2`$ symmetry. The symmetry enhancement operation is to shrink the edge which is common to the two hexagons (the thick line) until it has zero length, which results in the coalescence of two vertices. The final object formed is shown in the second figure in fig. 9. It has a four-valent vertex connecting four pentagons and the symmetry is enhanced from $`C_2`$ to $`C_4.`$ We find, empirically, that symmetry enhancement operations appear to take place in pairs, with a particular operation always being accompanied by the same operation on an opposite face of the associated polyhedron.

There is a $`\mathrm{C}_{28}`$ fullerene with $`D_2`$ symmetry (denoted 28:1 in ref.) that contains two of the structures shown in the first figure in fig. 9. If the above symmetry enhancement operation is performed on both these structures then the resulting object is precisely the $`D_{4d}`$ configuration of the $`B=9`$ Skyrmion described earlier. There are also $`D_2`$ symmetric $`C_{44}`$ (denoted 44:75 and 44:89 in ref. ) fullerenes to which similar statements apply. In this case, which is $`B=13`$, there are an equal number of pentagons and hexagons (12 of each) and so a very symmetric configuration can be obtained by applying the symmetry enhancement operation at all possible vertices (in this case 6) which results in the cubic Skyrmion. Selective application of the symmetry enhancement rule to these fullerenes allows one to create the associated polyhedron for the $`D_{4d}`$ symmetric Skyrmion at this charge, and also that for a $`D_4`$ symmetric Skyrmion, whose rational map we have been unable to deduce since the minimum energy rational map in this class of maps has $`D_{4d}`$ symmetry.

In the context of fullerenes it is, of course, impossible for vertices to coalesce since they correspond to the positions of the carbon atoms, but for Skyrmions the vertices represent local maxima of the baryon density and so there is no restriction that they be distinct; it just appears that in most cases it is energetically favourable to have distinct vertices. Note that, by an examination of the baryon density isosurface by eye, it can often be difficult to identify whether a given vertex is tri-valent or four-valent, since the edge length required to be zero for symmetry enhancement could be small, but non-zero.

### 5.3 Vertices, faces and rational maps

We will now attempt to explain the various features of the Skyrmions we have created by considering the basic properties of the rational map ansatz. Recall that the baryon density of a Skyrmion is localized around the edges of a polyhedron. The face centres of this polyhedron are given by the zeros of the Wronskian of the numerator and the denominator

$$w(z)=p^{}(z)q(z)q^{}(z)p(z),$$
(5.1)
which is in general a degree $`2B2`$ polynomial in $`z.`$ All the solutions which we have created have the property that all the roots of (5.1) are distinct and hence within the rational map ansatz this is a vindication of one of the GEM rules, in that it explains why the number of faces of the polyhedron is $`F=2B2.`$

Often (though not always, as we shall discuss further below) the position of the vertices correspond to local maxima of the density which occur in the integrand defining $``$ in equation (1.7). This density depends on the modulus of the rational map and its derivative so in general it is not possible to obtain such a simple characterization of the location of its maxima. However, in particularly symmetric cases they can be identified with the zeros of a polynomial constructed from the Hessian of the Wronskian . Explicitly, the polynomial is

$$H=(2B2)w(z)w^{\prime \prime }(z)(2B3)w^{}(z)^2,$$
(5.2)
and has degree $`4(B2).`$ As an example in which the above formula does work, consider the degree 7 rational map describing the icosahedrally symmetric minimal energy charge 7 Skyrmion

$$R=\frac{z^77z^57z^21}{z^7+7z^57z^2+1}.$$
(5.3)
The Wronskian and the Hessian in this case are given by

$$w=28z(z^{10}+11z^51),H=8624(z^{20}228z^{15}+494z^{10}+228z^5+1),$$
(5.4)
which are proportional to the Klein polynomials associated with the faces and vertices of a dodecahedron respectively .

Note that if the zeros of the Hessian (5.2) could always be identified with the vertices of the polyhedron then this would explain another of the GEM rules, that is, the fact that the number of vertices is equal to $`V=4(B2).`$ In this case the number of edges $`E=6(B2)`$ by the Euler formula, and hence $`E=3V/2`$, that is, the polyhedron is trivalent. Unfortunately, this is not true in general, though for the minimal energy maps it may be the case that the zeros of (5.2) give good approximations to the positions of the maxima.

We have already discussed explicit cases where the number of vertices and edges is affected by symmetry enhancement, and therefore it comes as no surprize that we cannot make general statements about the $`E,F`$ and $`V`$ based solely on the rational map ansatz. Although we do not have a general global characterization of the vertices it is possible, by a local analysis of the rational map, to check whether a given point is a local maximum and to obtain its valency.141414We thank Nick Manton for suggesting this possibility. By using the freedom to perform rotations of both the domain and target spheres it is always possible to choose the point we are considering to be given by $`z=0`$ and the rational map to have a local expansion about this point of the form

$$R=\alpha (z+\beta z^{p+1}+O(z^{p+2}))$$
(5.5)
where $`\alpha `$ and $`\beta `$ are real positive constants. A possible exception to this is if the derivative of the map is zero at the point we are considering, but since such points correspond to face centres they are clearly not maxima and may be ignored for our purposes here. Substituting the expansion (5.5) into the expression for the baryon density one can obtain the following result.

If $`p>2`$ then $`z=0`$ is a $`p`$-valent vertex if $`\alpha >1.`$ If $`p=1`$ then $`z=0`$ is not a vertex. The remaining case of $`p=2`$ is a little more subtle. In many cases there is a one-to-one correspondence between the vertices and the local maxima of the baryon density polyhedron. However, this is not always true and in some cases (the lowest charge example being $`B=5`$) only some of the local maxima are vertices, whilst others correspond to midpoints of an edge. In this situation some of the edges may appear thicker than others, reflecting their local maxima nature. The rational map description of such a bivalent maximum is the final $`p=2`$ case, where a local maximum requires the more restrictive condition that $`\alpha >\sqrt{1+3\beta }.`$

As an example of this analysis, consider the $`B=9`$ map with $`D_{4d}`$ symmetry given by (3.1). Expanding this map about the point $`z=0`$ gives

$$R=az+ib(1a)z^5+\mathrm{}$$
(5.6)
which can be rotated into the form (5.5) with $`p=4`$, $`\alpha =a=3.38>1`$ and $`\beta =b(1a).`$ Thus the point $`z=0`$ is a four-valent vertex, as we have observed from the baryon density plot. The other minimizing map with four-valent vertices (this time six of them) is the $`B=13`$ $`O`$ map (3.10), which can be checked in a similar way.

### 5.4 Isomerism — local minima and saddles

We have argued strongly that there is an analogy between Skyrmion solutions and polyhedra found in carbon chemistry. Moreover, at some charges we have found more than one solution which has very low energy, and therefore it might seem sensible within the analogy to chemistry to describe these solutions as isomers, whether they be saddle point solutions or local minima 151515The existence of degenerate minima, would probably require some kind of symmetry between the solutions. Although we cannot rule it out, it appears to be very unlikely..

In most cases, for example, $`B=9,16,17,19`$ and $`22`$ the symmetries and structures of the known isomers are unrelated to those of the minimal energy Skyrmions and in these cases it has been easy to identify the minimum energy configurations using the relaxation of initially well-separated clusters. The cases of $`B=10`$ and $`13`$ are interesting since there are known configurations whose associated polyhedra are related by symmetry enhancement. We have already commented that the polyhedron associated with $`O`$ symmetry at $`B=13`$ can be created by 6 symmetry enhancement operations from a $`D_2`$ fullerene polyhedron, and that with $`D_{4d}`$ symmetry requires just 2.

However, we have not discussed the $`B=10`$ solutions in this context. The $`D_{4d}`$ and $`D_{3h}`$ solutions do not appear to be related to this concept, but one can understand the $`D_3`$ and $`D_{3d}`$ solutions in terms of a more symmetric polyhedron which can be created from either by 3 symmetry enhancement operations. Clearly, this highly symmetric configuration, which is likely to be of higher energy, can be thought of as a saddle point in configuration space, with the two minima on either side. It is interesting to speculate that the true minimum energy Skyrmion is more closely associated with a polyhedron in which partial symmetry enhancement has taken place, that is, the bond lengths have shrunk, but not totally to a four valent bond. This might explain our difficulty in identifying the symmetry of the true minima using our methods.

### 5.5 Skyrmion architecture

One of the most important reasons for performing full field simulations in our study of Skyrmions is to verify that the minimal energy Skyrmion (at least for $`B22`$) consists of a single shell structure, which is the main assumption in the rational map ansatz. In this section we speculate on the kinds of structures which may form for Skyrmions of higher charge.

The lowest known value for the energy per baryon in the Skyrme model arises from an infinite three-dimensional cubic crystal , with an energy  only $`3.6\%`$ above the Faddeev-Bogomolny bound. In considering a large single-shell fullerene Skyrmion, where hexagons are dominant, the twelve pentagons may be viewed as defects, inserted into a flat hexagonal structure, in order to generate the required curvature necessary to close the shell. Energetically the optimum structure of this form is an infinite hexagonal lattice and this was constructed in one of our earlier papers , and found to have an energy per baryon which is $`6.1\%`$ above the Faddeev-Bogomolny bound. This is therefore the value to which a bigger and bigger single-shell structure will asymptote. Since this value is higher than for the Skyrme crystal it is reasonable to expect that above some critical charge $`B^{}`$, the minimal energy Skyrmion will resemble a portion cut from the crystal rather than a single shell. However, it is very difficult to estimate the value of $`B^{}`$ (note that we have seen that at least $`B^{}>22`$) since it relies on a delicate comparison of the surface to volume energy of a finite portion of the crystal and this is very sensitive to the way in which the portion of the crystal is smoothed off at the edges. As the crystal is basically composed of stacking $`B=4`$ cubes together then $`B=32`$ is the first charge at which any sizeable chunk of the crystal can emerge, and even then it has a large area to volume ratio, so perhaps the charge will need to be even larger than this before a crystal regime takes over from the shell regime.

An intermediate between a single-shell and a crystal is a multi-shell structure and this has recently been studied by Manton and Piette . For the charges they consider ($`B=12,13,14`$) the relaxation of an initial multi-shell structure produces a single-shell configuration which has relatively high energy in comparison with the minimal energy single-shell Skyrmions we have found. As a shell may be thought of as a spherical domain wall, connecting the two vacua $`U=\pm 1`$, then all configurations with an odd number of shells have $`U=1`$ at the centre, whereas if there are an even number of shells then $`U=1`$ at the centre. Thus a single shell structure obtained from an initial odd number of shells can relax to one of the configurations we have found. In fact, as we have already mentioned, the $`B=13`$ Skyrmion obtained by Manton and Piette is the $`O_h`$ symmetric saddle point solution whereas the minimal energy Skyrmion is only $`O`$ symmetric.

In summary, there are a number of alternatives to a single-shell structure for higher charge Skyrmions and what is remarkable is that none of these alternatives appear to arise at least for $`B22.`$ It seems reasonably clear that single-shells can not be the whole story for large enough charge, but whether this charge is so large as to be irrelevant in applications to nuclear physics has yet to be determined.

### 5.6 Relation to applications

In the introduction we commented on two diverse motivations for creating Skyrmion solutions, namely from a purely mathematical point of view, to study an interesting class of maps between 3-spheres which generalize the harmonic map equations, and from a physical perspective to investigate a phenomenological model of nuclei. Here, we comment briefly on the relevance of our results to these two applications and suggest interesting avenues for future research which we have opened up with this work.

We have already noted that the Skyrme model is the simplest model in which one finds stable solitonic solutions which correspond to maps from $`S^3`$ to $`S^3`$, and so our solutions may have some generality to other extensions. An interesting feature of the solutions which we have found is that, in some sense, they can be thought of as being close to a conformal map between the two 3-spheres, for which the three eigenvalues of the strain tensor, $`\lambda _1^2,\lambda _2^2,\lambda _3^2`$, would all be equal. The rational map ansatz, which we have seen provides a good approximation to the true solutions, has two of these eigenvalues equal and it has been observed  that the shape of the profile function appears to be such that the deviation from a conformal map is minimized when averaged over space.

For a conformal map which is locally an isometry the values of $``$ and $``$ would be exactly equal and so any deviations from the map being locally isometric can be visualized by plotting $``$. When one does this the relevant isosurface is highly localized around the edges of the associated polyhedron and also in the centre of each face where $``$ is close to zero and $``$ is large in comparison. Such an isosurface is a plot of second order effects due to curvature. The fact that the associated polyhedra are generally of the fullerene type is also interesting because in chemistry such structures arise since they minimize what is called steric strain, the overall strain of the delocalized electron distribution. Using this analogy we suggest that the effects of the strain tensor for maps between 3-spheres can be thought of as being analogous to steric strain in fullerene molecules.

Finally, we should note that the existence of a Skyrmion with a particular symmetry, which can be described by a rational map, implies that there also exists an SU(2) BPS monopole with the same symmetry, although, of course, all BPS monopoles of a given charge have the same energy. The fields and Lagrangians of monopoles and Skyrmions are very different but the structures which arise in each case are remarkably similar. This suggests that these types of configurations may be generic as low energy states in a variety of 3-dimensional soliton models and elsewhere.

Although the original motivation for studying the Skyrme model was to make quantitative predictions for the properties of nuclei based on a model which is derived in some limit from QCD, this has historically been a tricky business. Our hope, based on the extensive results that we have presented here, is that some progress can be made in this direction. Part of the problem in achieving such a goal is how one should understand the model in the context of nuclei. Based on the idea that it is a low energy effective action for QCD several studies have attempted to quantize Skyrmion solutions as rigidly rotating spinning tops for $`B=1,2`$ and $`3`$ (see, for example, refs.) quantifying corrections in terms of their position in the $`1/N_\mathrm{c}`$ expansion. This is not only complicated but probably is also too simplified an approach, as demonstrated by the more sophisticated quantization of the $`B=2`$ Skyrmion in ref.. Here we would like to make a few interesting points in terms of thinking of the solitons as a phenomenological model for nuclei which incorporates classical isospin.

The first thing that we would like to discuss is the shell structure of the soliton solutions which we have created. Naively, one might think that such a structure is incompatible with the solutions modelling nuclei which are usually assumed to be solid lumps. What one has to appreciate to reconcile this is, that until the solutions are quantized, the charge can be thought of as a fluid. If a nuclei were assumed to be comprized of a simple positively charge fluid then a shell structure would be expected under the action of the strong force since there is a long range, but weak, attraction at long distances, and a short distance repulsion due to nucleon-nucleon repulsion. Therefore, the hollow shell structure of the multi-Skyrmions which we have observed is largely due to the continuum version of nucleon-nucleon repulsion.

Following on from this point we should note that some of the features of the classical values of $`I_B`$, the ionization energy, are very much in line with expectations based on nuclei. Let us focus in detail on the solutions for $`B=4`$ and $`B=5`$. The $`B=4`$ solution is highly symmetric and $`I_4`$ is relatively large, whereas the $`B=5`$ solution has little symmetry and $`I_5`$ is much smaller. This is exactly as one might have expected since $`{}_{}{}^{4}\mathrm{He}`$ is the most stable nucleus whereas there is no naturally occurring stable nucleus with $`A=5`$. Since the packing structure of the solutions is a feature of the symmetry of the solution this suggests that there may be something even more than just a good model of the strong force potential within the Skyrme model.

There is also an interesting trend in $`\mathrm{\Delta }E/B`$, the classical binding energy per baryon, which appears to asymptote to a value defined by that of an infinite Skyrmion lattice. We have plotted this compared to the experimentally determined values for nuclei with $`A=1`$ to $`22`$  in fig. 10 with an arbitrary normalization factor (which amounts to multiplying the curve in fig. 7 by about 50) accounting for our ability to define the Skyrmion energy units. Although this is crude, it makes the point that the curve has the correct shape. This is very encouraging and is the subject of on-going research.

Since the fullerene polyhedra are clearly very important for our understanding of Skyrmions, and as we have just argued there are a number of appealing features of the model for explaining the properties of nuclei, it is tempting to make an analogy between the delocalized electron distributions in fullerenes and nuclear charge distributions. Although the analogy is not exact, it might be possible to relate the Skyrme model to density functional theory methods (see, for example, ref. ) used in the study of electron distributions.

## 6 Conclusion

We have performed an exhaustive study of minimal energy Skyrmions for all charges upto $`B=22`$, using a variety of methods and involving a substantial amount of CPU time on a parallel machine. At each charge we have discussed in detail the symmetry, structure and energy of the minimal energy Skyrmion (and often several others) in addition to providing an approximate description by presenting its associated rational map. Supplementary to the detailed investigation at specific charges we have found a number of interesting general phenomena. These include the verification of the fullerene hypothesis, which applies to all except two cases (which can be understood in terms of symmetry enhancement), the discovery that there are often several Skyrmions with very different symmetries from the minimal one but nonetheless have energies which are remarkably close to the minimal value, and finally the confirmation that the shell-like structure of Skyrmions continues to large charges (at least $`B=22`$) with the rational map ansatz providing an effective approximation to the true solution. Hopefully this comprehensive piece of work will provide a useful foundation for further studies on Skyrmions, both mathematical and physical, with the ultimate aim being a comparison with experimental data on nuclei.

## Acknowledgements

Many thanks to Conor Houghton, Nick Manton and Tom Weidig for useful discussions. We thank the EPSRC (PMS) and PPARC (RAB) for Advanced Research Fellowships. PMS acknowledges the EPSRC for the grant GR/M57521. The parallel computations were performed on the COSMOS at the National Cosmology Supercomputing Centre in Cambridge.