Geometrically Constrained Skyrmions

Abstract

Skyrmions are chiral swirling magnetization structures with nanoscale size. These structures have attracted considerable attention due to their topological stability and promising applicability in nanodevices, since they can be displaced with spin-polarized currents. However, for the comprehensive implementation of skyrmions in devices, it is imperative to also attain control over their geometrical position. Here we show that, through thickness modulations introduced in the host material, it is possible to constrain three-dimensional skyrmions to desired regions. We investigate skyrmion structures in rectangular FeGe platelets with micromagnetic finite element simulations. First, we establish a phase diagram of the minimum-energy magnetic state as a function of the external magnetic field strength and the film thickness. Using this understanding, we generate preferential sites for skyrmions in the material by introducing dot-like “pockets” of reduced film thickness. We show that these pockets can serve as pinning centers for the skyrmions, thus making it possible to obtain a geometric control of the skyrmion position. This control allows for stabilization of skyrmions at positions and in configurations that they would otherwise not attain. Our findings may have implications for technological applications in which skyrmions are used as units of information that are displaced along racetrack-type shift register devices.

1. Introduction

Magnetic skyrmions, predicted by theory almost 30 years ago, have advanced to a central topic of research [1,2,3] in nanoscale magnetism over the last decade following their experimental observation [4,5]. Their particular topological properties [6], which impart them high stability and particle-like behavior [7,8,9], combined with their room-temperature availability [10,11], reduced dimensions [12], and their unique dynamic properties [13], render these magnetic structures promising candidates for future spintronic applications [14]. Skyrmions are formed in non-centrosymmetric magnetic materials exhibiting a sufficiently strong Dzyaloshinsky–Moriya Interaction (DMI) [4,5,15], that is, an antisymmetric energy term that favors the arrangement of the magnetization in helical spin structures with a specific handedness and spiral period. In extended thin films, in addition to the DMI, the formation and stability of skyrmions depends sensitively on various parameters, such as the strength of an externally applied magnetic field, the film thickness, temperature, and the magnetic history of the sample. Phase diagrams have been reported in the literature [10,16,17,18], displaying the parameter ranges within which skyrmions are stable and where they may take different forms. Skyrmions may typically either develop individually or in the form of a hexagonal skyrmion lattice. While the occurrence of individual skyrmions makes them attractive candidates for units of information that can be displaced in a controlled way by spin-polarized electric currents, their spontaneous arrangement in the form of a periodic lattice could be interesting for magnonic applications [19], where these point-like magnetic structures could play the role of scattering centers of planar spin waves. One drawback of possible applications of skyrmions is the difficulty of controlling their position. For instance, if skyrmions are to be used as units of information in race-track-type shift-register devices, it is not only necessary to be able to displace them in a controlled way, but also to make sure that they are shifted between well-defined positions along the track. In domain-wall-based concepts for race-track memory devices [20], which preceded the skyrmion-based variants, this control of the position was typically achieved by inserting indentations (“notches”) into the strips [21,22]. It was shown that such notches represent preferential sites for domain walls, making it possible to trap domain walls at these specific positions, from which they could only be detached after overcoming a certain depinning energy [23]. Although the possibility to capture skyrmions at specific sites has been addressed in the case of ultrathin films, geometric control analogous to the pinning of domain walls at notches does not yet seem to be firmly established for skyrmions. In a recent study, Suess and coworkers reported that skyrmions could be pinned in a racetrack-type geometry by introducing semicircular notches at the lateral boundaries [24]. However, the authors also found that skyrmions have a tendency to lose their topological stability and to dissolve when they interact with these notches. In two-dimensional systems, strategies for the pinning of skyrmions include the insertion of point-like defects [25,26] or atomic-scale vacancies [27]. Remarkably, randomly distributed point defects in ultrathin films have also been reported to have little effect on the current-driven skyrmion dynamics [7]. Motivated to further the discussion on this topic by addressing the three-dimensional case of “bulk” DMI, we use finite-element micromagnetic simulations to study the extent to which geometric control of the skyrmion position in a thin-film element can be achieved by introducing a patterning in the form of local variations of the film thickness. Compared to holes or lateral notches, these geometric features have the advantage of preserving the topological stability of skyrmions interacting with them. Our simulations show that, by locally lowering the film thickness in sub-micron-sized, dot-shaped regions, skyrmions can in fact be “captured” at these geometrically defined sites. We find that, by using geometrical constrains of this type, skyrmions can be stabilized at positions that they would otherwise not adopt. For instance, these geometric manipulations make it possible to generate regular, square lattices of skyrmions which, apart from exceptional situations [28], are not observed naturally in non-centrosymmetric ferromagnets. It is argued that this control of skyrmion positions in magnetic thin films can open up new possibilities for skyrmionic devices, as well as for concepts of magnonic metamaterials.

2. Results

Before addressing the question of how preferential sites for skyrmions can be generated through nanopatterning, we first investigate, as a preliminary study, the field- and thickness dependence of skyrmionic structures forming in rectangular FeGe platelets.

2.1. Chiral Magnetization States in a Helimagnetic Rectangular Platelet

We consider rectangular thin-film elements with a lateral size of $180 \mathrm{n}\mathrm{m}\times 310 \mathrm{n}\mathrm{m}$ and thicknesses ranging between 5 $\mathrm{n}$$\mathrm{m}$ and 75 $\mathrm{n}$$\mathrm{m}$, and simulate the magnetic structures forming in the presence of a perpendicular external magnetic field with a flux density varying between 0 $\mathrm{m}$$\mathrm{T}$ and 900 $\mathrm{m}$$\mathrm{T}$. The simulations yield a large variety of possible magnetization states in this thickness and field range, which can be classified into six non-trivial types, summarized in Figure 1. The three main types are the helical state shown in panel (a), which is characterized by the presence of regular spin spirals extending over large parts of the sample, the bimeron state (c), which can be interpreted either as a particular type of skyrmion structure that is stretched along one axis or, alternatively, as a helical state in which the extension of the helices is limited, and finally, the skyrmion lattice state (e), which is characterized by a regular, hexagonal arrangement of skyrmions.

In addition to these three fundamental states, there are also hybrid states in which two different types of structures coexist [29], such as the helical-bimeron state shown in Figure 1b and the bimeron-skyrmion state shown in Figure 1d. These mixed states can be considered as transitional configurations between one fundamental state and another, which appear with changes in the external field value or in the film thickness.

At sufficiently large field values, the platelet can also sustain individual skyrmions which are not arranged in the form of a lattice. Figure 2a shows an example of a magnetic structure containing a single skyrmion, located at the center of a 60 $\mathrm{n}$$\mathrm{m}$ thick platelet in a 650 $\mathrm{m}$$\mathrm{T}$ field. At this film thickness, the single-skyrmion state does not represent the minimum energy configuration. Nevertheless, the skyrmion is stable and possesses the usual properties of topological protection. In particular, the skyrmion could be manipulated and displaced, such as by means of spin-polarized electrical currents, as is done in racetrack-type devices. We will show later that it is possible to influence the position of such freely moving skyrmions by introducing geometric variations in the platelet, and thereby generate effective pinning sites for skyrmions. In order to obtain more information on the size and the profile of the skyrmion, we display the magnetization components along a line scan through the center, as shown in Figure 2b. The graphs show that the skyrmion is of Bloch type, with vanishing radial magnetization component along the central line. The diameter of the skyrmion core, defined as the distance between two diametrically opposed zero-crossings of the z-component of the magnetization, is almost exactly 30 $\mathrm{n}$$\mathrm{m}$. We found that the skyrmion diameter does not display noticeable variations if the film thickness changes, such as from 60 $\mathrm{n}$$\mathrm{m}$ to 30 $\mathrm{n}$$\mathrm{m}$.

At elevated field values and larger film thicknesses, a quasi-homogeneous state is formed (not shown), where the magnetization is largely aligned along the external field direction. At fields below saturation, chiral bobber (ChB) [30,31] structures are also observed. These complex configurations of the magnetization can be considered as variants of skyrmions which do not traverse the entire thickness of the sample. Instead, they have a skyrmion-like structure only on one surface, which evolves into a quasi-saturated configuration on the opposite surface on a path along the film thickness. The apex of the ChB contains a Bloch point at which the magnetic structure changes in a discontinuous way. ChB structures have interesting micromagnetic properties, and have recently been discussed as magnetic structures that could be attractive in the context of spintronic devices [31], but they are not of primary interest for our study. We display an example of a ChB structure only for completeness in the upper right of Figure 1f, where it coexists alongside four ordinary skyrmions. The magnetic structures shown in Figure 1 were obtained by starting from a random initial configuration and by a subsequent energy minimization. For more details, see Section 4.

2.2. Phase Diagram of the Magnetization States

The various magnetic states described in the previous section are possible equilibrium configurations of the magnetization forming in the FeGe platelets at different values of the external field and the film thickness. It is important to note that these magnetic structures are not uniquely determined by the film thickness and the field strength. Because of this, in order to avoid possible misunderstandings, we did not specify the values of the thickness and the field strength at which the states shown in Figure 1 occur. In fact, several metastable states that can be significantly different from each other are often possible under identical conditions, depending only on the magnetic history of the sample or, in a numerical experiment, on the initial conditions of the simulation. While it is generally not possible to identify a unique magnetization state that develops in the thin-film element, micromagnetic simulations can be used to determine the type of magnetic structure that has the lowest energy. The results of these calculations are summarized in the phase diagram shown in Figure 3a.

Although the magnetic structure at a specific thickness and field value is generally not unique, the phase diagram helps identifying the most preferable structure as far as the total energy is concerned. While at lower field values (below about 400 $\mathrm{m}$$\mathrm{T}$) the phase diagram is rather complex, evidencing a multitude of possible magnetic structures showing neither any clearly dominating state nor a significant thickness dependence, the situation becomes simpler at larger field strengths (above about 600 $\mathrm{m}$$\mathrm{T}$). Two main states emerge in these ranges of larger field values—the skyrmion configuration and the quasi-saturated states. Moreover, these states are separated by a clearly defined boundary in the phase diagram, showing a distinct impact of the film thickness. Specifically, if a field of 650 $\mathrm{m}$$\mathrm{T}$ is applied, the formation of skyrmion structures will be energetically favorable if the film thickness is below 50 $\mathrm{n}$$\mathrm{m}$, while a quasi-saturated state will be the lowest-energy configuration at larger thickness values, as shown in Figure 3b. This observation represents the fundamental of the concept of geometrically constrained skyrmions that we present in this study. The idea is the following. If the film thickness is locally modulated within a small dot-shaped region such that, at a given field, the skyrmion structure is favorable in that thinner part while in the rest of the sample, the thickness is large enough to favor a quasi-homogeneous state, these thickness modulations can be designed to capture skyrmions. As we will show, this patterning makes it possible to generate pinning sites for skyrmions and, to some extent, to achieve a geometric control of the skyrmion position within the thin-film element.

2.3. Geometrically Constrained Skyrmions

We now consider magnetic structures forming in an FeGe platelet of 60 $\mathrm{n}$$\mathrm{m}$ thickness containing dot-like cylindrical cavities within which the thickness is locally reduced to 30 $\mathrm{n}$$\mathrm{m}$. The phase diagram displayed in Figure 3 suggests that, at external field values of about 650 $\mathrm{m}$$\mathrm{T}$, the insertion of these cavities results in a geometry with specific regions favoring the stability of skyrmion structures in a thin-film element which would otherwise favor a quasi-homogeneous magnetic configuration. This can lead to the formation, or the trapping, of skyrmions that are geometrically constrained to the regions in which the pockets have been introduced. Figure 4 shows such a geometrically constrained skyrmion in a 60 $\mathrm{n}$$\mathrm{m}$ thick platelet. The skyrmion remains confined to the small region in which the thickness is reduced by 50 through two cylindrical pockets with a depth of 15 $\mathrm{n}$$\mathrm{m}$ and radius of $r=20 \mathrm{n}\mathrm{m}$, inserted symmetrically on both the top and bottom surfaces. In order to accommodate a single skyrmion, the pocket radius is chosen to be slightly larger than the radius of the skyrmion core (15 $\mathrm{n}$$\mathrm{m}$), obtained by a line-scan of the magnetization components of an isolated skyrmion, as shown in Figure 2.

The geometrically constrained skyrmion, shown in Figure 4, is stabilized by the geometry for two reasons. Firstly, as discussed before, in this field range the skyrmion state is generally favored because of the reduced film thickness. Secondly, the vortex-like magnetic configuration forming on the interior cylinder surfaces of the cavity helps in pinning the position of the skyrmion to the center of the pocket. This cylindrical flux-closure structure thereby provides boundary conditions, albeit not in a mathematical sense, which constrain the skyrmion to this dot-like geometry. By forming such a cylindrical vortex structure, the magnetization finds a nearly optimal way to adapt to competing micromagnetic interactions. It thereby satisfies both the tendency of the DMI to introduce chiral, swirling patterns, as well as the tendency imposed by the magnetostatic interaction to form flux-closure structures with the magnetization aligned along the surfaces. The simulations show that a symmetric insertion of these pockets on both the top and bottom surfaces is necessary to obtain the desired stability and localization of skyrmions. If the thickness variation is introduced only on one of the surfaces, the pinning of skyrmions appears to be much less effective. This result suggests that the previously mentioned swirling of the magnetization at the interior of the cylindrical pockets contributes significantly to the stabilization of the skyrmion position, and that sufficiently strong pinning can only be achieved if these swirling boundary conditions are imposed on both sides of the film.

If geometric modifications of the sample surface, as described above, can stabilize a skyrmion, the question arises whether this effect can be used to place skyrmions at specific positions where they might be generated or removed in a controlled way through external manipulation. This could be of interest, such as for device concepts in which skyrmions are utilized as binary units of information, in a context similar to that of dot-patterned magnetic media for high-density data storage [32]. In this case, the skyrmion pockets would take the role of the magnetic nanodots in bit-patterned media. While it is beyond the scope of this study to discuss the technical feasibility of such storage media or to explore the ability to write and delete individual skyrmion patterns into the pockets, we can show that, indeed, it is possible to stabilize skyrmions in various geometrically predefined locations that could be addressed individually.

Figure 5a–e shows several examples of simulations in which the position of skyrmions in a thin-film element can be predetermined by introducing several pockets of the type discussed before. As shown in Figure 5e, our simulations predict the possibility to stabilize six skyrmions at well-defined positions, placed on a regular grid, in our sub-micron FeGe platelet. Although the results shown in Figure 5 may suggest a nearly optimal geometric control of the skyrmion positions, it is important to note that the pockets discussed here merely provide preferential sites for skyrmions. The latter may or may not form or remain pinned at those sites. In particular, it is not sufficient to thin-out a part of the sample to ensure the appearance of geometrically constrained skyrmions. The purpose of such pockets could rather be to capture existing skyrmions and to fix them at well-defined positions, similar to the domain-wall pinning role that is played by notches in conventional racetrack-memory devices [20,33]. It should also be noted that the geometric trapping of skyrmions with such pockets does not always work, in particular when the pockets are too closely packed. As a rule of thumb, the material must observe a characteristic minimal distance between the skyrmions that is given by the material-dependent, long-range helical period $l_D$, which in the case of FeGe, is about 70 $\mathrm{n}$$\mathrm{m}$ (see Section 4).

If two skyrmions are constrained in adjacent pockets, the magnetization rotates smoothly by 360 ${}^{\circ }$ along a line connecting the cores of these skyrmions. Such a continuous rotation of the magnetization direction involves the formation of a spin spiral that must be compatible with the properties of the material. In particular, two skyrmions cannot be constrained geometrically at distances that would be too small to accommodate spin spirals with a periodicity that is much smaller than $l_D$. We emphasize that, in this context, $l_D$ only serves as an estimate for the minimal distance between skyrmions, and that inter-skyrmion distances may be larger than $l_D$, in particular in the case of stronger external fields [34]. Beyond the simple analytic estimate of the minimal skyrmion spacing based on $l_D$, detailed studies have established inter-skyrmion repulsion as a general effect [35,36,37], which ensures that a minimum distance between neighboring skyrmions is preserved. As an example of a dense array of constrained skyrmions approaching this limit, the nearest-neighbor spacing in the array (measured as the distance between the centers of the pockets) of constrained skyrmions shown in Figure 5a is 90 $\mathrm{n}$$\mathrm{m}$. This spacing is compatible with the analytically estimated minimum distance of 70 $\mathrm{n}$$\mathrm{m}$ given by the material’s helical long-range period, $l_D$. In this geometry, adding more pockets would reduce the minimal distance below this value, and thereby destabilize the pinning effect. We also found that skyrmions cannot be stabilized at positions too close to the lateral sample boundaries. This observation is consistent with the skyrmion-edge repulsion effect [8,35,37,38,39], which leads to an increase in energy as a skyrmion approaches the edges of a thin-film element. An example of such a failed attempt to stabilize skyrmions close to the lateral boundaries is shown in Figure 5f. In spite of these limitations, the ability to geometrically constrain skyrmions provides an attractive way to obtain control over the skyrmion position in thin-film elements, which could have important technological implications.

3. Discussion

By means of micromagnetic finite-element simulations, we have presented a possibility to control the position of magnetic skyrmions at predefined positions within a thin-film element by introducing cylindrical nano-pockets graved into the surface. Our concept for geometrically constraining skyrmions via such dot-like thickness variations is, in many ways, analogous to the idea of geometrically constrained domain walls [40] in cylindrical nanowires, or to studies in which indentations have been introduced in rectangular strips in order to capture head-to-head domain walls in racetrack-type memory devices [22]. In those cases, too, the desired effect of the geometric constraint is to define preferential sites for specific micromagnetic structures, such that the magnetic structures constrained at those artificial pinning sites require a certain activation energy in order to detach from them. The pockets described in this work could effectively play this role in the case of skyrmions driven along magnetic strips by means of spin-polarized electrical currents. Such geometric control of their position would allow for shifting skyrmions between well-defined points on the track. With the perspective of a possible implementation of this concept in skyrmion-based racetrack devices, future research could address the precise pinning strength that would result from those pockets and the current density that would be required to depin skyrmions captured in such pockets. If the skyrmion pinning effect resulting from pockets, as described in this study, should turn out to be too strong for application purposes, the pinning potential could be reduced in a gradual and straightforward way by modifying the geometric parameters of the pockets. For instance, a weaker pinning effect could be achieved by reducing the depth of the pockets, or by smoothing the edges of the cylinder, such as to yield rounded or cone-shaped pockets.

While racetrack-type skyrmion devices are probably the most straightforward application in which our idea to capture skyrmions at geometrically defined sites could be implemented, the concept of constrained skyrmions could also be useful for other purposes. As mentioned before, the trapping of skyrmions at dot-like sites could serve as a principle for skyrmion-based data storage devices, without necessarily involving any displacement or depinning processes. If skyrmions can be selectively generated and dissolved at such preferential sites, such as by means of the field of a magnetized nano-tip or through a localized spin-polarized current traversing the film thickness, the geometrically constrained skyrmions could represent units of information that could be written and erased. Perhaps such skyrmionic dot material could even be stacked in three dimensions for ultra high-density storage purposes. To address such possibilities, future research directions could explore ways to reversibly insert skyrmions in these geometrically defined regions. Another potentially interesting use of our concept concerns magnonic applications [19,41]. Since skyrmions can act as point-like scattering centers for spin waves, the ability to arrange them at regular lattice sites, as described in this study, could open up new perspectives, as this could result in a new type of magnonic metamaterial in the form of artificial magnon Bragg lattices consisting of skyrmions. Such artificial structures could be tailored to yield specific scattering and interference properties for spin waves that could not be obtained otherwise.

4. Materials and Methods

The material modelled in this study is FeGe. Due to its well-known helimagnetic properties, this B20-type, non-centrosymmetric material serves as a prototype for materials hosting chiral magnetic structures that develop due to the “bulk” DMI effect, as opposed to certain systems of ultra-thin magnetic films and substrates that can generate an “interfacial” DMI [1]. The micromagnetic parameters of FeGe are [18,42] $A=8.78\times {10}^{-12}\mathrm{J}{\mathrm{m}}^{-1}$, $M_{\mathrm{s}}=384\mathrm{k}\mathrm{A}{\mathrm{m}}^{-1}$, and $D=1.58\times {10}^{-3}\mathrm{J}{\mathrm{m}}^{-2}$, where A is the ferromagnetic exchange constant, $M_{\mathrm{s}}$ the saturation magnetization, and D the DMI constant. We neglect any magnetocrystalline anisotropy of the material, setting the uniaxial anisotropy to zero, $K_{\mathrm{u}}=0\mathrm{J}{\mathrm{m}}^{-3}$. A characteristic length scale of this material is the long-range helical period $l_d=4\pi A/\left|D\right|\simeq 70 \mathrm{n}\mathrm{m}$ [43]. This length scale describes the typical period length of magnetic spirals forming as a result of the competing interactions of the ferromagnetic exchange on one hand, and the DMI on the other. The ordering temperature of FeGe is $278.7$ $\mathrm{K}$ [43]. Therefore, skyrmion devices based on this material are not expected to be operational at room temperature. Chiral structures in FeGe have been observed experimentally, such as at 100 $\mathrm{K}$ [37], but recent studies have evidenced that in this material system, skyrmions can be stable also near room-temperature [10,44]. These developments suggest that, although thermal stability is a topic of central importance in skyrmion research [3], temperature limits do not represent a significant obstacle in the case of FeGe.

With the material parameters defined above, the total energy $E_{\mathrm{tot}}$ of the system is given by the sum of the Zeeman term, the ferromagnetic exchange, the DMI interaction, and the magnetostatic energy:

$$E_{\mathrm{tot}}=\underset{V}{\int }\left({\mu }_0{\mathit{H}}_{\mathrm{ext}}·\mathit{M}+A·\sum _{i=x,y,z}{\left(\nabla {\mathit{m}}_i\right)}^2+D\mathit{m}·\left(\nabla \times \mathit{m}\right)-\frac{{\mu }_0}{2}\mathit{M}·\nabla u\right)\mathrm{d}V$$

(1)

Here, V is the sample volume, ${\mathit{H}}_{\mathrm{e}xt}$ is the externally applied magnetic field, ${\mu }_0=4\pi \times {10}^{-7}\mathrm{V}\mathrm{s}{\mathrm{A}}^{-1}{\mathrm{m}}^{-1}$ is the vacuum permeability, $\mathit{m}=\mathit{M}/M_{\mathrm{s}}$ is the reduced (normalized) magnetization, and u is the magnetostatic scalar potential. The magnetostatic (demagnetizing) field ${\mathit{H}}_{\mathrm{d}}=-\nabla u$ is the gradient field of the magnetostatic potential. We calculate the magnetostatic potential u, which accounts for the long-range dipolar interaction, by using the hybrid finite-element method/boundary element method (FEM/BEM) introduced by Fredkin and Koehler [45,46]. The dense matrix occurring in the boundary integral part of this formalism is represented using ${\mathcal{H}}^2$-type hierarchical matrices [47]. This data-sparse representation effectively overcomes size limitations arising from the boundary element method, as it yields a linear scaling of the computational resources required for the calculation of the magnetostatic term, which would otherwise grow quadratically with the number of degrees of freedom on the surface.

For each energy term, an effective field ${\mathit{H}}_{\mathrm{eff}}$ is defined as the variational derivative of the corresponding partial energy E,

$${\mathit{H}}_{\mathrm{eff}}(\mathit{r},t)=-\frac{\delta E\left[\mathit{M}\left(\mathit{r},t\right)\right]}{{\mu }_0\delta \mathit{M}}.$$

(2)

Specifically, the effective field of the ferromagnetic exchange is

$${\mathit{H}}_{\mathrm{eff}}^{\left(\mathrm{xc}\right)}(\mathit{r},t)=\frac{2A}{{\mu }_0{M}_{\mathrm{s}}}\Delta \mathit{m},$$

(3)

and the effective field of the DMI is

$${\mathit{H}}_{\mathrm{eff}}^{\left(\mathrm{DMI}\right)}(\mathit{r},t)=-\frac{2D}{{\mu }_0{M}_{\mathrm{s}}}\left(\nabla \times \mathit{m}\right).$$

(4)

Together with the magnetostatic field and the external (Zeeman) field, these effective fields enter the Landau-Lifshitz-Gilbert (LLG) equation [48], which describes the evolution of the magnetization field $\mathit{M}(\mathit{r},t)$ in time,

$$\frac{\mathrm{d}\mathit{M}}{\mathrm{d}t}=-\gamma \left(\mathit{M}\times {\mathit{H}}_{\mathrm{eff}}\right)+\frac{\alpha }{M_s}\left(\mathit{M}\times \frac{\mathrm{d}\mathit{M}}{\mathrm{d}t}\right),$$

(5)

where $\gamma$ is the gyromagnetic ratio and $\alpha$ is a phenomenological, dimensionless damping constant. We use the LLG equation to calculate equilibrium structures $\mathit{M}\left(\mathit{r}\right)$ of the magnetization, by integrating in time until convergence is reached. In the numerical simulations, convergence is achieved when either the total energy ceases to change over a long period, or when the torque (magnitude of the right hand side of the LLG equation) drops below a user-defined threshold. The skyrmions stabilized in the pockets, as shown in Figure 5, result from a relaxation process starting from an out-of-plane saturated state of the film. As is common practice in micromagnetic simulations, we consider thermal effects only implicitly via the material parameters. Hence, temperature does not enter directly as a parameter in the simulations.

The geometry of the samples is designed with FreeCAD [49] and the discretization into linear tetrahedral elements is performed with Netgen [50]. The visualization of the FEM data was done with ParaView [51]. The finite-element cell size does not exceed $2.5$ $\mathrm{n}$$\mathrm{m}$, which is well below the exchange length $l_{\mathrm{ex}}=\sqrt{2A/{\mu }_0{M}_s^2}\simeq 9.7 \mathrm{n}\mathrm{m}$, in order to avoid discretization errors. The finite element meshes in this study contain typically about one million finite elements. The micromagnetic simulations are done with our proprietary GPU-accelerated finite-element software [47].

The discretized representation of the vector field of the magnetization is given by a value of $M_i$ defined at each node (vertex) i of the finite-element mesh. The magnetostatic field, as well as the effective fields of the ferromagnetic exchange and the DMI, are calculated within each tetrahedral element. The element-based data of these fields are then mapped onto the nodes of the mesh, in order to calculate the effective field acting on the magnetization, and thus to calculate the evolution of the magnetization in time at each node according to the LLG equation. More details on the calculation of the converged magnetization states are given in Appendix A.