Size and profile of skyrmions in skyrmion crystals

Size is a fundamental quantity of magnetic skyrmions. A magnetic skyrmion can be a local circular object and in an isolated form. A skyrmion can also coexist with a group of its siblings in a condensed phase. Each skyrmion in a condensed phase takes a stripe shape at low skyrmion density and a circular shape at high skyrmion density. Skyrmions at high density form a skyrmion crystal (SkX). So far, skyrmion size in an SkX has not been seriously studied. Here, by using a generic chiral magnetic film, it is found that skyrmion size in an SkX has a different parameter dependence as those for isolated skyrmions and stripes. A size formula and a good spin profile for skyrmions in SkXs are proposed. These findings have important implications in searching for stable smaller skyrmions at the room temperature. Understanding the parameters which govern the size of a skyrmion are important to gain greater control over their stability and dynamics but the underlying mechanisms differ from individual skyrmions to skyrmion crystals. To this end, using micromagnetic simulations, the authors theoretically investigate the relationship between the underlying material parameters and skyrmion radius within a skyrmion crystal

In this paper, a generic chiral magnetic film with exchange stiffness constant A, the DMI strength D, and perpendicular magnetic anisotropy K is used to show that κ ≡ π 2 D 2 / (16AK) = 1 separates isolated skyrmions (κ < 1) from condensed skyrmion states (κ > 1). In contrast to isolated skyrmions whose size increases with D/K and is insensitive to κ ≪ 1 and stripes whose width increases with A/D and is insensitive to κ ≫ 1, the skyrmion size in SkXs is inversely proportional to the square root of skyrmion number density and decreases with A/D. This finding provides general guidance for skyrmion size manipulation.

Model.
A two-dimensional (2D) chiral magnetic film of thickness d in the xy-plane is modelled by the following energy where K u , M s , μ 0 , and H d are perpendicular uniaxial crystalline anisotropy constant, the saturation magnetization, the vacuum permeability and demagnetization field, respectively. The energy is set to zero, E = 0, for ferromagnetic state of |m z | = 1. In a theoretical analysis, the demagnetization effect can be included through the effective perpendicular anisotropy constant . This is a good approximation when the film thickness d is much smaller than the exchange length 34 . Films with κ ≡ π 2 D 2 /(16AK) < 1 support isolated skyrmions that are metastable spin textures of energy 8πAd ffiffiffiffiffiffiffiffiffiffiffi 1 À κ p 34 .
Demagnetization field as an effective anisotropy. In order to prove that demagnetization field can accurately be included in the effective anisotropy K ¼ K u À μ 0 M 2 s =2 even for highly inhomogeneous spin structures, we consider a 300 nm × 300 nm × 0.4 nm film with A = 10 pJ m −1 , D = 6 mJ m −2 , Starting from a nucleation domain of 10 nm in diameter and under the dynamics of the LLG equation, a stable ramified stripe has the same skyrmion number as an isolated skyrmion or a skyrmion in an SkX, (see "Methods" section). The stable ramified stripe from MuMax3 is presented in Fig. 1a. In MuMax3 38 simulations, the unit-cell size is 1 nm × 1 nm × 0.4 nm. The black and red dots (90,000 in total) in Fig. 1b denote the crystalline anisotropy energy E a and the sum of the anisotropy and demagnetization energies E a + E d for each cell from an MuMax3 simulation, respectively. Data are sorted by m z . The blue curve is the effective anisotropy energy where V is the volume of each cell. The red dots are almost on the blue curve, demonstrating the excellence of the approximation.
Phase diagram of an isolated skyrmion and a ramified stripe. Figure 2 a plots the energy of one skyrmion in chiral magnets as a function of κ. The squares and diamonds, circles and stars, uptriangles and pentagons, are respectively from three groups of model parameters listed in "Methods" section. E is plotted in arcsinhðEÞ scale because E covers a huge range of values on both sides of zero, and κ is in log ðκÞ scale. Clearly, κ = 1 separates metastable isolated skyrmions from ramified stripes 36 .
For a better comparison, earlier results on spin profile of isolated skyrmions and stripes are summarized first. The spin profile of isolated skyrmions can be approximated by ΘðrÞ ¼ 2 arctan sinhðr=wÞ sinhðR=wÞ h i 34 . Θ is the polar angle of the magnetization at position r measured from the center of a skyrmion where m z = 1. R and w measure respectively the skyrmion size and skyrmion wall thickness. Figure 2b shows the excellent agreement between the theoretical spin profile (the solid curve) with R = 4.513 nm and w = 1.356 nm with numerical simulations from MuMax3 (the symbols) for an isolated skyrmion. Model parameters are The inset shows the skyrmion structure. This profile leads to the skyrmion size formula of R ¼ w= ffiffiffiffiffiffiffiffiffiffiffi 1 À κ p and w = πD/(4K) 34 .
The spin profile of stripes can be approximated by ΘðxÞ ¼ 2 arctan sinhðL=2wÞ sinhðjxj=wÞ h i for m z < 0 and ΘðxÞ ¼ 2 arctan sinhðjxj=wÞ sinhðL=2wÞ h i for m z > 0 36 . x is measured across a stripe and x = 0 is at the stripe centre (m z = 1). L is the stripe width that is the distance between two parallel contour lines of m z = 0. w measures skyrmion wall thickness. Figure 2c is spin profile of stripes with model parameters of The symbols are the numerical data from various locations of stripes shown in the inset and the solid curve is the theoretical fit with L = 8.95 nm and w = 1.82 nm. The excellent agreements between the numerical data and the approximate profiles demonstrate that one can use this theoretical profile to extract the skyrmion size accurately. Data from different stripes, falling on the same curve, demonstrate that stripes, building blocks of irregular spin textures of skyrmion number 1, are identical. The stripe profile leads to the stripe width formula of L = g(κ)A/D, where g(κ) is about 2π for κ ≫ 1 36 , opposite dependence on D as that for an isolated skyrmion.
From an earlier study 36 , it is known that skyrmion-skyrmion interaction is important in SkX formation unlike isolated skyrmions and stripes that do not rely on the interaction. Thus, one should expect a different spin profile and parameter dependence of skyrmion size for an SkX. This is the main issue in this study.
SkX and spin profile. In order to study skyrmions in SkXs, we consider model parameters satisfying κ > 1 below. All initial configurations have sufficient number of nucleation domains in triangular lattices. Each domain of m z = 1 is a disk of 10 nm in diameter. Each initial configuration leads to an SkX in a triangular lattice as shown in Fig Skyrmion size in SkXs. If energy per skyrmion in an SkX is from a skyrmion whose spin texture is described by Θ(r) for 0 ≤ r ≤ L/2, then energy per skyrmion can be expressed as a function of variable R, 3) are defined in the "Methods" section. The skyrmion size is then the solution of dE dR ¼ 0. Under certain reasonable assumptions (see Supplementary Note 2), one has Unlike stripe width that does not depend on n and increases with A/D, Eq. (3) says that R is almost inversely proportional to ffiffiffi n p and decreases with A/D. The interesting dependences agree with the underlying physics of SkXs that result from the balance of skyrmion-skyrmion repulsion and intrinsic stripe nature of skyrmions. R should then be proportional to the lattice constant of an SkX that is the skyrmion-skyrmion distance and, in turn, is inversely proportional to ffiffiffi n p . Since skyrmions in an SkX come from the squeeze of stripe skyrmions, the amount of width deduction from the squeeze should be proportional to the original stripe width proportional to A/D. Thus, R should not only be smaller than the natural stripe width, but also decreases with A/D as what Eq. (3) says. Of course, one should not expect coefficients of 1= ffiffiffi n p and A/D to be very accurate because the obvious approximations involved. For example, it is impossible to fill a  We carried out MuMax3 simulations with very different n, A, D, K in three different parameter regions. For each fixed parameters, skyrmion size R in the corresponding SkX can be obtained either from the size of m z = 0 contour or through fitting of numerical spin profile to the new spin profile. The difference between two extracted values is negligible (see Supplementary  Fig. 1), and all data shown below is from contour measurement. One can test our theoretical predictions of Eq. (3) against numerical simulations. As mentioned above, simulations show that R depends only on n, A/D, and κ. Figure 4a plots R vs. 1= ffiffiffi n p for various fixed A/D and κ from three groups of model parameters. Simulation results (the symbols) are well described by our formula (the solid curves) of Eq. (3) without any fitting parameter. R is almost linear in 1= ffiffiffi n p . The A/D dependence and κ dependence of R can also be well captured by Eq. (3) as shown by good agreement between simulation data (the symbols) and theoretical formula of Eq. (3) (the solid curve) in Fig. 4b (A/D) and Fig. 4c (κ), respectively. Equation (3) says that R does not depend on anisotropy K for κ ≫ 1. Our new spin profile captures successfully the skyrmion-skyrmion repulsion effect that requires R to be proportional to the lattice constant, the inverse of the skyrmion number density n.

Discussion
Although most previous studies [39][40][41] investigate the lattice constant of SkXs rather than skyrmion size in an SkX, ref. 42 42 ]. Thus, according to our theory, the skyrmion size should be around 22.7 nm that agrees well with the measured value of~21 nm [estimated from Fig. 2h in ref. 42 ]. One should also note that the presence of magnetic field can also affect the size of skyrmions that explains the small discrepancy.
According to the current results, magnetic skyrmions can be isolated, or condensed regular or irregular stripes distributed randomly or arranged as helical states of one-dimensional stripe    lattices, or condensed circular skyrmions in crystal forms. The irregular stripes can even be in ramified, or in dendrite or maze forms 36 . The underlying physics of skyrmion size are different for isolated skyrmions, stripes in different morphologies and structures, and skyrmions in an SkX. It is known that skyrmion formation energies are positive for isolated skyrmions and negative (with respect to the energy of a ferromagentic state) for stripes and SkXs. Positive formation energy results in an extra energy for skyrmion surfaces, similar to the surface tension of a liquid droplet, that favours a circular shape. Negative skyrmion formation energy is similar to a negative surface tension of a liquid droplet that prefers to wet its contacts and to spread out. From this viewpoint, it may not be surprising to see different parameter dependences of skyrmion size for isolated skyrmions and skyrmions in SkXs or stripes, each with skyrmion number 1. It is also not surprising that spin profile of skyrmions in SkXs has little similarity to those of isolated skyrmions and stripes because their shape comes from the balance of the skyrmion squeeze and their stripe nature. Two nearby skyrmions start to repel each other when their distance is order of L ≈ 13A/D. Thus, the critical skyrmion number density for SkX formation is about n c ≈ D 2 /(169A 2 ). For n > n c , skyrmions tend to form a closely packed triangular lattice of lattice constant 1:07= ffiffiffi n p . The skyrmion squeeze leads to Eq. (3) in which skyrmion size decreases with A/D, opposite to the parameter dependence of the width of stripe skyrmions. There is a maximal skyrmion number density beyond which an SkX is not stable any more. The physics around the maximal skyrmion number density surely deserves a study, but it is not an issue in this work.
In applications, creating stable skyrmions of smaller size at the room temperature is a goal. According to current results, one needs to distinguish an isolated skyrmion (κ < 1) from a skyrmion in a condensed skyrmion phase (κ > 1). In the case of κ < 1, materials with large A and smaller D are preferred since larger A and smaller D imply a higher T c and smaller skyrmion size. On the other hand, in the case of κ > 1, materials with larger A and larger D are required so that one may have both higher T c and smaller skyrmion size. Interestingly, the same critical value of κ = 1 was also obtained by many people from comparing energy difference between a one-dimensional magnetic domain wall and the ferromagnetic state 8,28,44,45 .
In conclusion, material parameter dependences of skyrmion size are very different for circular isolated skyrmions, stripes, and circular skyrmions in SkXs. Unlike an isolated skyrmion whose size increases with D/K and κ < 1, or a stripe whose size increases with A/D and decreases with κ > 1, the size of a circular skyrmion in SkXs is inversely proportional to the square root of skyrmion number density and decreases with A/D. Skyrmion size is a very important quantity, and high-density data storage requires, in general, smaller skyrmions. The findings here provide a practical guidance for searching proper materials in different applications.

Methods
Numerical simulations. The dynamics of m is governed by the Landau-Lifshitz-Gilbert (LLG) equation, where γ and α are respectively gyromagnetic ratio and the Gilbert damping constant. The effective magnetic field H eff ¼ 2A the exchange field, the crystalline magnetic anisotropy field, the demagnetization field H d , and the DMI field H DM respectively.
Previous studies 36,37 demonstrate that a small nucleation domain of m z = 1 in the background of ferromagnetic phase of m z = − 1 can develop into a skyrmion in a chiral magnetic film. In the absence of energy source such as an electric current and at zero temperature, the LLG equation describes a dissipative system whose energy can only decrease 46,47 . Thus, the long-time solution of the LLG equation with any initial configuration must be a stable steady spin texture of Hamiltonian Eq. (1). α does not change the stable steady spin structures of the LLG equation.
However, α can change spin dynamics and which stable spin texture to pick because the evolution path and energy dissipation rate depend on damping. This is similar to the sensitiveness of attractor basins to the damping in a macro spin system 48 . In this work, we use a large α = 0.25 to speed up our simulations.
In this study, we choose three groups of model parameters to simulate those chiral magnets of MnSi with A = 0. 27 To obtain the energy of one static isolated skyrmion or one static ramified stripe, we start with one small domain of m z = 1 of 10 nm in diameter. After a few nanoseconds, one can obtain a stable steady spin structure as well as quantities like Q and E from MuMax3 simulator.
In order to obtain an SkX, many small domains of m z = 1 of 10 nm in diameter each are initially arranged in a triangular lattice in the background of m z = −1. Under the dynamics of the LLG equation, each small domain of zero skyrmion number becomes a skyrmion of Q = 1. For each giving n > n c , an SkX in triangular lattice is obtained. We can obtain SkXs with different density n by varying the number of domains in a film of a fixed-size.  16A p 2 ðεÞ 1 κ in the denominator is very small for ε > 0.5. Replacing p 2 (ε) by p 2 (0.5), skyrmion size R, is then Eq.(3).

Data availability
The data that support the plots within this paper and other findings of this study are available from the corresponding author on reasonable request.