Magnetic hopfions in solids

Topological solitons^1,2 are localized finite energy solutions to classical non-linear field equations appearing in many fields of science, from nuclear physics³ to cosmic string theory.⁴ Topological solitons are characterized by their ability to move and interact with each other as ordinary particles. They are one answer to Heisenberg’s question how countable particles can appear in continuous fields. Recent celebrated examples in condensed matter are magnetic chiral skyrmions,^5–8 two-dimensional vortex-like field configurations in magnets without inversion symmetry, which have received considerable interest as fundamental objects as well as promising candidates for future spintronic applications.⁹ In contrast, models, which allow three-dimensional (3D) topological solitons,^2,10 to which we shall refer to as hopfions according to a classical topological concept due to Hopf,¹¹ seem very rare in nature. It is worth mentioning that the term “hopfion” is often used for the states characterized by variants of the Hopf map.^11–13 Such textures resembling Hopf fibration were experimentally observed in liquid crystals,^14–19 Bose–Einstein condensate,²⁰ and magnetic multilayers²¹ and theoretically predicted in ferromagnetic chiral magnets,^22–27 antiferromagnetic frustrated chiral magnets,²⁸ and ferroelectrics.²⁹ 

The first statically stable magnetic hopfions have been reported by Bogolubsky in Refs. 30 and 31, where an isotropic micromagnetic model possessing O(3) symmetry with higher-order derivatives of the order parameter has been used. The hopfion solutions found by Bogolubsky are stable at zero applied field when the ground state of the system represents a collinear ferromagnetic state. Recently, Sutcliffe has shown that hopfion solutions for isotropic frustrated magnets described by O(3) symmetric Hamiltonian can also be stabilized in a strong external magnetic field that breaks the symmetry and suppresses the exchange spirals representing the ground state of the system.³² 

Here, we identify the spatially isotropic functional used in Refs. 30–32 as a degenerate case of a more general micromagnetic functional necessary to describe magnetic texture of crystals with competing exchange interactions and cubic symmetry. We present this advanced micromagnetic functional derived from spin-lattice Hamiltonian for simple, face-centered, and body-centered cubic lattices. Following our approach, similar functionals can be derived for materials of arbitrary crystal symmetry.³³ In addition to the presented hopfion solutions, we also provide a criterion for the existence of hopfions in the systems described by such a micromagnetic functional and give an estimate of a typical real space dimensions of such hopfions. This criterion is based on single inequality, where all variables can be deduced from quantum mechanical models of solids typically with the help of density functional theory (DFT).³⁴ Such a computational materials science approach opens a multiscale vista for the discovery of hopfions in real materials by bridging quantum mechanics to higher-level materials specific models.

Topological solitons occur if the order parameter takes values in a curved manifold, e.g., a Lie group or a magnetization sphere $S2$⁠. In this case, the topological concept of homotopy may be used to distinguish certain classes of field configurations from the topologically trivial uniform state. More precisely, for field configurations representing a topological soliton, there is no continuous deformation (homotopy) to any collinear (constant) field configuration.

In case of magnetic hopfions, the relevant order parameter of the system is a unit vector field n(r) = (n_x, n_y, n_z), |n(r)| = 1, defined at any point $r∈R3$⁠. Field configurations n attaining a uniform background state at infinity n(r) → n₀ as |r| → ∞ can be classified according to the linkage of their fibers {n = p}, which, for regular values $p∈S2$⁠, are collections of closed loops in $R3$⁠. Examples of simple and more intricate hopfion configurations are illustrated in Fig. 1.

The Hopf invariant H = H(n) of a field n onto $S2$ is defined as the linking number of two generic fibers, i.e., the oriented number of times fibers wind around each other. The crucial point is that H is a classifying homotopy invariant. The existence of a homotopically non-trivial map with unit H, nowadays known as Hopf fibration, is a classical result due to Hopf¹¹ in 1931 also emerging implicitly in Dirac’s work on quantized singularities of the same year.³⁵ Hopfions are topological solitons of non-zero Hopf invariants, whose existence or non-existence relies on the structure of the governing energy functional.

It has generally been believed that in magnetically ordered crystals, 3D textures may only arise as precessing or moving states,^36–38 while energy dissipation processes eventually lead to instability. In certain cases, the dynamical states under the influence of various external stimuli such as electric current or AC magnetic fields can be stabilized on a quite long time scale.^39,40 Here, we do not consider such dynamic states and study the statically stable solutions only.

In contrast to the classical Skyrme model,³ the energy density of the conventional micromagnetic theory,⁴¹ $E∝|∇n|2≡∑α(∇nα)2$⁠, where α ∈ {x, y, z}, lacks stabilizing terms, either in form of higher powers of derivatives of the magnetization density, n(r), or higher order derivatives. Below, we show that energy functionals incorporating such terms can be derived from the classical isotropic Heisenberg model. This model is microscopically defined through the following Hamiltonian:

which describes an exchange interaction between magnetic spin moments μ located on atomic sites i and j, where n are unit vectors, n = μ/|μ|. ∑_i>j indicates the summation over all pairs of interacting spins. Here, we go beyond typical model approaches where the strength of pair interactions $Jij(rij)$ is set to zero at the distances above a particular value. We take into account both cases when the absolute values of $Jij(rij)$ decay with the distance between spins as well as the case when $Jij(rij)$ do not vanish or decay very slowly with the distance, which is typical for some metals.⁴² Of interest to us are systems with competing ferro- and antiferromagnetic interactions as a function of distance. The various potential energy terms such as magnetocrystalline anisotropy or Zeeman energy can be taken into account. However, as shown in Sec. II, they are not required for the stability of hopfions.

We consider crystals with cubic Bravais lattices: simple (sc), body-centered (bcc), and face-centered (fcc). For a correct continuum representation [n_i → n(r)] of the spin-lattice Hamiltonian (1), which is able to reflect the essence of the competing exchange interactions, we go beyond the conventional micromagnetic approximation and take into account the higher order terms in the series expansion (see Sec. 1 of the supplementary material for details). We obtain the advanced micromagnetic energy functional

where each index runs over x, y, and z. The derivation provides linear relations between micromagnetic and exchange constants J_s, symmetry equivalent representatives of $Jij$⁠,

where a is the lattice constant. The positive coefficients $as$⁠, $bs$⁠, and $cs$ depend on the crystal lattice type (see Tables S1–S3 of the supplementary material and Ref. 33).

For real materials, the micromagnetic constants $A$⁠,$B$⁠, and $C$ can be obtained either experimentally measuring spin-wave spectra along different crystallographic directions or using density functional theory (DFT) calculating the energy density for flat spin spiral textures,⁴³ n_s = (cos(q · r), sin(q · r), 0), with the wave vector q, which reads

Figure 2 illustrates that, in contrast to conventional ferromagnets characterized by a quadratic dispersion relation $Es∼q2$⁠, functional (2) features a quartic behavior of $Es$ near the Γ-point along at least one of the high-symmetry directions.

For negative $B$ or $C$⁠, the energy density (4) is unbounded from below when q → ∞, which indicates the tendency of the system to antiferromagnetic order where functional (2) becomes irrelevant. The occurrence of antiferromagnetic hopfions, which are allowed by the spin-lattice model (1), will be considered elsewhere, while below we assume that $B≥0$⁠, $C≥0$⁠, and $B+C>0$⁠.

In the particular case of $C=6B$⁠, functional (2) describes a spatially isotropic system and can be reduced to a simple form

where Δ = ∇² denotes the Laplacian. We refer to this particular case as the model of Bogolubsky,^30,31 who first reported on the hopfion as a statically stable solution of (5) for $A>0$⁠. Recently, in Ref. 32, for a particular case of model (5) when $A<0$⁠, hopfions stabilized by an external magnetic field were reported. Notably, functional (5) has been also considered in framework of a pure 2D model with $A>0$ in Ref. 44 and with $A<0$ in Ref. 45. Note that functional (5) can be written in the following equivalent form:

There are two straightforward approaches that allow one to find hopfion solutions for particular material parameters. The first approach is based on direct energy minimization of the corresponding spin-lattice Hamiltonian (1). This approach requires the exact values of the exchange coupling constants $Jij(rij)$⁠, which can generally be found by DFT methods. An obvious obstacle in this approach is the case when a set of $Jij(rij)$ represents an infinite series. The second approach assumes an energy minimization of the micromagnetic functional (2), for instance, employing various finite-difference schemes. In this case, the material parameters $A$⁠, $B$⁠, and $C$ can be extracted from dispersion curves (Fig. 2). The main disadvantage of this approach is that the chosen type of finite-difference scheme and the density of the discrete mesh affect the precision and performance of the numerical calculations. Implementating a finite-difference scheme for the functional containing forth order terms may prove computationally less efficient than the approach discussed below.

For an efficient numerical approximation of hopfion solutions for the advanced micromagnetic functional, we propose a minimal effective model based on a spin-lattice Hamiltonian for a simple cubic lattice, where each spin interacts with four types of neighbors—the first four shells. The constants of these exchange interactions are denoted by $J̃1,2,3,4$ and may be adapted to the effective micromagnetic constants $A$⁠, $B$⁠, and $C$⁠.

Let us consider the following example. Suppose, a magnetic crystal has fcc cubic symmetry with the lattice constant a_fcc = 0.4 · 10⁻⁹ m. Let us assume that the material parameters are $A=1.5⋅10−14$ J m⁻¹, $B=0.8⋅10−32$ J m, and $C=5.2⋅10−32$ J m. In this case, those exact values of $Jij(rij)$ of the corresponding Heisenberg model (1) are not essential for finding the solutions. The constants $A$⁠, $B$⁠, and $C$ can be found either from their relations to $Jij(rij)$ provided in Tables S1–S3 of the supplementary material or can be directly extracted from the experimental or theoretical dispersion curves (Fig. 2).

Let us now consider a simple cubic lattice with the lattice constant $ã=afcc$ and coupling constants $J̃1,2,3,4$⁠, which give the same material parameters $A$⁠, $B$⁠, and $C$⁠. According to (3) and Table S1 of the supplementary material, these constants must satisfy the following system of equations:

It is also necessary to take into account an additional constraint that guarantees that the ferromagnetic state is lower in energy than the antiferromagnetic state; therefore,

For the considered example here, one can choose $J̃1=10−21$ J (⁠$≈6.2$ meV), $J̃2=0.188J̃1$⁠, $J̃3=−0.274J̃1$⁠, and $J̃4=−0.161J̃1$⁠. It is easy to show that this set of $J̃1,2,3,4$ satisfies (7) and (8). The results of direct energy minimization of such an effective model with those parameters $J̃1,2,3,4$ will correspond to the solution of both micromagnetic functional (2) and the corresponding Heisenberg model (1) of the fcc crystal.

In the case of a simple cubic crystal with nonzero exchange interactions only in the first four shells, the minimal effective model is exactly the Heisenberg model (1). Below, whenever we refer to the results obtained with the minimal effective model, we use a tilde symbol for the exchange interaction constants.

For positive $A$⁠, the energy density in (2) is non-negative and provides a true extension of a conventional ferromagnet governed by $E∝|∇n|2$⁠. In this case, functional (2) satisfies the Derrick–Hobart criterion¹ necessary for the existence and stability of solutions localized in 3D space. From a mathematical perspective, the advanced micromagnetic model features striking analogies with the Skyrme–Faddeev model.⁴⁶ The latter one is accompanied by a well-established theory on the occurrence of string-like solitons.

By means of the effective model introduced in Sec. II C, we have found a wide spectrum of solutions with different Hopf indices. Figure 3 illustrates hopfions for the limiting cases $C=0$⁠, $C=6B$⁠, and $B=0$⁠. The axial symmetry of the toroidal hopfions with Hopf indices H = 1, 2 reflects the isotropic structure of functional (2) in case of $C=6B$⁠. Note that the hopfion solution with H = 1 for the particular case of $C=6B$ was first suggested by Bogolubsky.³⁰ The numerical solution in Ref. 30 was found by energy minimization of the spin-lattice Hamiltonian, which is conceptually similar to our minimal effective model (see Sec. II C). In particular, the calculations were performed for a simple cubic lattice with four shells, $J̃1=61$⁠, $J̃2=−10$⁠, $J̃3=0$⁠, and $J̃4=−5$⁠, which corresponds $C≈6B$⁠. This set of parameters has also been used in the simulations of dynamics of magnetic hopfions.⁴⁷ 

By applying mathematical transformations, we show that (2) is bounded from below by a variant of the Skyrme–Faddeev functional (the details are provided in Sec. 5 of the supplementary material). The latter one has a topological lower bound, also known as Vakulenko–Kapitanski inequality.⁴⁶ After adapting constants,^48,49 we find

with an explicit constant $η=(32π2/31/8)2Amin(C,6B)$⁠, which is not expected to be optimal. This lower energy bound can be matched qualitatively by suitable trial configurations of the prescribed Hopf invariant as in Ref. 50, proving the previously mentioned fractional energy law. During completion of this work, we became aware of the work⁵¹ reporting similar results related to energy bounds but for the case $A<0$⁠, which is not considered here.

The behavior of the numerical solutions for hopfions with different Hopf indices confirms the fractional power law for minimal energies E ∼ |H|^3/4 (see Fig. 4). By virtue of the arguments in Ref. 52, the sublinear energy growth essentially implies the attainment of minimal energies in infinitely many homotopy classes. Such behavior of the energy of hopfions with increasing topological index is in contrast to the 2D systems. The linear energy growth for positive topological charge Q is a crucial element in the rigorous proof of minimal energy attainment for chiral skyrmion with Q = −1.⁷ In contrast to Belavin–Polyakov solutions,⁵³ however, the energy bounds are only attained approximately by collapsing Belavin–Polyakov lumps. It is also interesting to note that the physics of two-dimensional skyrmions with diverse topological charges can be quite different.⁵⁴ 

The characteristic scale of inhomogeneity, l₀, which can be seen as a measure of smoothness of the hopfion’s vector field, represents an important parameter related to the stability of such states in the spin-lattice model. We define l₀ as

where n₂ and n₁ are vectors at arbitrary positions r₂ and r₁, respectively, and ∠(n₂, n₁) is angle between these vectors. The characteristic scale, l₀, can be assumed to be nearly independent on Hopf index, H, while all hopfions can be decomposed on similar structural elements (for details, see Sec. 6 of the supplementary material). The characteristic scale of inhomogeneity was estimated analytically by analyzing the structure of an infinite hopfion. This technique is based on the approach proposed by Faddeev² in which the hopfion is considered as a twisted and looped skyrmion tube (see also Refs. 55 and 56). For some σ-models,^48,49 this approach allows us to perform an exact analytical analysis assuming that the length of such a closed tube tends to infinity. We adapted this technique for the case of the spatially anisotropic Hamiltonian (2) and derived an analytical expression for the estimate of the characteristic scale of inhomogeneity (see Sec. 6 of the supplementary material),

If we denote the distance between nearest neighbors of the lattice as d (sc: d = a, bcc: $d=a3/2$⁠, and fcc: $d=a2/2$⁠), then δ_m = d/l₀ defines the maximal angle between neighboring spins. Equation (11) together with (7) allows us to estimate the dependence of δ_m on the exchange constants. Figure 5 illustrates such a dependence for the minimal effective model. According to our estimate, δ_m should not exceed ∼45°, half of the critical angle for spin-lattice models,⁵⁷ because the continuum functional (2) may become an insufficient approximation of the lattice Hamiltonian (1). An approximate inequality d/l₀ ⪅ π/4 together with estimation (11) provide the following criterion for the existence of hopfions:

For the case of nonzero external magnetic field, B_ext, one has to add an additional term to functional (2),

where E is the exchange part of the energy defined by (2). For positive $A$⁠, the ground state of the system is a spin-polarized state with n(r) ⇈ B_ext for each r. On the other hand, for negative $A$ and 0 ≤ |B_ext| < |B_s|, the ground state of the system is a spin spiral. Here, B_s is the critical field value above which the system undergoes a transition into the spin-polarized state. Let us consider the case of $Bext⇈êz$⁠. In the most general case, functional (13) admits the solution corresponding to equilibrium spin spiral, $n(r)=cos(q⋅r)sin(θ),sin(q⋅r)sin(θ),cos(θ)$⁠, of two types,

and

The saturation field, B_s, is defined by the condition θ(B_ext) = 0. Taking into account relation (7), one can extrapolate solutions (14) and (15) to the case of the spin-lattice model. Indeed, Eqs. (15) and (7) correlate well with some well-known classical results for lattice models. For instance, as follows from (15) and relations (7) in the one-dimensional case, the condition for appearance of the spin spiral in the ground state is $J̃4<−J̃1/4$⁠, which agrees with earlier results.⁵⁸ 

The degenerate case when $C=6B$ and the first term E in (13) becomes spatially isotropic has been considered earlier in Ref. 32. It is easy to see that in this case, solutions (14) and (15) coincide, and the orientation of the q-vector of spiral is not restricted to the directions of 111-type and 100-type—all directions for spin spiral propagation are energetically equivalent. The saturation field in the degenerate case is $Bs=A2/(8MsB)$⁠.

Hopfion solutions for the degenerate case have been reported only for the external magnetic field, which is twice as large as the saturation field B_ext = 2B_s.³² The study of hopfion solutions for the general case (13) at moderate fields B_ext ∼ B_s goes beyond the scope of this work and will be presented elsewhere.

In order to prove that hopfions remain stable in realistic systems, one has to take into account the effect of dipole–dipole interactions (DDI), which may become significant in the case of finite size samples as films, plates, and discs. We consider a 30 nm thick film of a simple cubic crystal with a lattice constant $ã=0.5$ nm. The normal of the film is parallel to the [001] direction. We assume that the magnetic moment at each atom is equal to one Bohr magneton, |μ| = 1μ_B. The nearest-neighbor ferromagnetic exchange interaction is $J̃1=10$ meV, and the antiferromagnetic exchange for the neighbors in the fourth shell is $J̃4=−0.245J̃1$⁠. According to (7), the Heisenberg coupling constants correspond to the micromagnetic constants as $A≈3.2⋅10−14$ J m⁻¹, $B≈2.4⋅10−32$ J m, and $C≈4.9⋅10−32$ J m. Direct numerical energy minimization was performed on a domain of 256 × 256 × 60 lattice sites and with periodic boundary conditions along the x- and y-axes in the plane of the film. Figure 6(a) illustrates a stable configuration of the hopfion with Hopf number H = 1 in the center of the film. In order to estimate the contribution of the DDI to the total energy of the hopfion, we performed independent energy minimizations for hopfions with and without DDI and compared the energies of these states, E_h, with the energies of the ground state E_g with and without DDI, respectively. The corresponding differences, ΔE_h = E_h − E_g, are ΔE_h = 30.6J with DDI and ΔE_h = 29.4J without DDI. Thereby, one may conclude that in our model with quite realistic parameters, the contribution of the DDI is very small.

Another important aspect of hopfion stability in geometrically confined systems is related to the interaction of the hopfion with the edges of the sample. When the initial hopfion configuration is shifted closer to one of the free surfaces, the system first convergences quickly to a near-equilibrium configuration as found in the bulk. After a surprisingly large number of iterations, the particle starts to move gradually toward the nearest free surface and eventually escapes from the film. Such a behavior represents one of the manifestations of the particle-like properties of hopfions. By centering the initial hopfion guess to the middle plane of the film, one can prevent the escape of the hopfion through the boundary.

In the following, we shortly discuss some important aspects related to the observation of hopfions in realistic systems with modern experimental techniques.

The inhomogeneity of the magnetization and the corresponding inhomogeneity of the magnetic induction field allow the detection of the hopfions by various electron microscopy techniques such as Lorentz Transmission Electron Microscopy (TEM)⁵⁹ or off-axis electron holography.⁶⁰ The latter allows quantitative measurements of the phase of the electron wave passing through the sample in the projection mode. Figure 6(b) shows the calculated magnetic phase image, which is expected to observe in off-axis electron holography experiment.⁶⁰ Modern instruments allow to achieve spatial resolution up to a few nanometers and a phase resolution of δφ ≈ 2 · 10⁻³ rad.⁶⁰ As seen from Fig. 6(b), that should be enough for hopfion observation.

Feasible alternatives for hopfion observation are also Lorentz electron microscopy with exit wave reconstruction,^61,62 electron tomography,^63,64 and advanced x-ray microscopy and spectroscopy techniques.^21,65,66

In conclusion of this section, we want to discuss in short some effects we observe near the transition between the ferromagnetic ground state when $A>0$ and spin spiral state when $A<0$⁠. The observations listed below have an aim to illustrate rich physics standing behind the presented model here and are expected to motivate further research in this field.

In order to observe stable hopfion solutions for both $A>0$ and $A<0$⁠, we apply an auxiliary external magnetic field B_ext = 50 mT along the plane of the film $(Bext⇈êy)$⁠. For such relatively weak magnetic field, the isolated hopfion with H = 1 remains stable with and without DDI in a finite range of $J̃4$ values, which in reduced units corresponds at least to the range $−0.254≤J̃4/J̃1≤−0.24$⁠. By varying value of $J̃4$⁠, one can effectively change the sign and magnitude of micromagnetic constant $A$⁠, which leads to the listed below effects.

• For $−0.239≤J̃4/J̃1$⁠, the constant $A>0$⁠, and with DDI, the hopfion becomes unstable and collapses. However, without DDI and with $J̃4/J̃1=−0.239$⁠, the hopfion is stable.

• For $J̃4/J̃1=−0.25$⁠, according to (7), the constant $A$ vanishes, i.e., $A=0$⁠. This case may turn to be quite interesting for the following study because according to Ref. 67, in a similar case, the hopfions in the Skyrme–Faddeev model may resemble compactons.⁶⁸ 

• For the range of $−0.252⪅J̃4/J̃1<−0.25$⁠, the constant $A<0$⁠, but its absolute value is small and the applied magnetic field exceeds the saturation field, B_ext > B_s. The lowest energy state in this case is the spin-polarized in-plane state. This case, therefore, is similar to the suppressed spiral mode that was studied in Ref. 32.

• The case $−0.254≤J̃4/J̃1⪅−0.252$ corresponds to $A<0$ in the regime when B_ext < B_s. Despite the fact that the system transforms into a spiral state of type (14), the hopfion remains stable being embedded (“interposed”) in that modulated state. Similar solutions have been recently reported⁶⁹ for skyrmions in the 2D isotropic $(C=6B)$ system, which represents 2D version of functional (13). Those skyrmions can also be considered as embedded into the modulated phase. Remarkably, the energy of these skyrmions does not depend on the displacement. Similar behavior is expected for hopfions presented here in bulk crystals.

• For $J̃4/J̃1=−0.255$⁠, the hopfion becomes elliptically unstable.

Contrary to chiral magnetic skyrmions requiring materials with strong spin–orbit interaction in combination with a lattice lacking inversion symmetry,^5–7,43 the criteria for the existence of magnetic hopfions looks more accessible. The magnetic hopfions discussed in this work can exist in both centrosymmetric and noncentrosymmetric crystals and do not require the presence of spin–orbit coupling at all. The key requirement for the formation/stability of single/isolated hopfions is according to our advanced micromagnetic energy functional a result of the proper competition between the conventional ferromagnetic exchange expressed by the typical spin-stiffness term proportional to square of the wavevector and higher-order terms proportional to a quartic relation of the wavevector. The latter can be a result of competing exchange interaction as put forward in this paper or due to higher order or multi-spin interactions such as the bi-quadratic, three-spin four-site, four-spin,⁷⁰ or chiral–chiral interactions.⁷¹ The proper competition manifests itself in the material parameters of the crystal by a flat magnonic dispersion curve near the gamma point (see Fig. 4). Examples of materials possessing such property are known.³³ The approach for comparing experimental data from magnets with complex interactions and the results of DFT calculations is also well established (see, for instance, Ref. 72).

Knowing criterion (12) motivates a DFT driven virtual design of materials for stable hopfions á la carte⁷³ with traditional approaches or the coevolutionary search method.^74,75 At the moment, the most promising techniques for detection and direct observation of hopfions seem to be x-ray microscopy/nanotomography and off-axis electron holography/tomography. Thus, we expect experimental confirmation of their existence in the near future.

Hopfions offer exciting perspectives in the future information technology in which hopfions appear as the smallest 3D information carrying particles in spintronics. Electrons passing hopfions will experience an emergent magnetic field created by the non-coplanar nature of the hopfion’s magnetization texture, giving rise to a hopfion Hall effect,⁷⁶ enabling the individual readout of Hopf invariants. Spin-currents can be used to move hopfions in all three directions leading eventually to the development of truly 3D memories. The non-linear dynamical response of hopfions in a box offers fascinating perspectives for neuromorphic computing.

See the supplementary material for the series expansion techniques, calculation of the topological charge, numerical energy minimization, continuous formula for the energy functional, bounding from below by a variant of the Skyrme–Faddeev functional, evaluation of characteristic scale of inhomogeneity for hopfion textures.

We thank O. Chugreeva, B. Dupé, P. F. Bessarab, J. Garaud, A. Samoilenka, and E. Babaev for discussions. C.M. acknowledges funding from Deutsche Forschungsgemeinschaft (DFG Grant No. ME 2273/3-1). S.B. acknowledges funding from DFG through SPP 2137 “Skyrmionics” (Grant No. BL 444/16-2) and through the Collaborative Research Center SFB 1238 (Project C01) as well as the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (Grant No. 856538, project “3D MAGiC”). L.D., C.M., and S.B. acknowledge seed-fund support from JARA-FIT. N.S.K. acknowledges Deutsche Forschungsgemeinschaft (DFG) for financial support through SPP 2137 “Skyrmionics” (Grant No. KI 2078/1-1). The work of F.N.R. was supported by the Swedish Research Council Grant Nos. 642-2013-7837, 2016-06122, and 2018-03659, by the Göran Gustafsson Foundation for Research in Natural Sciences and Medicine, by Olle Engkvists Stiftelse, and by the “Roland Gustafssons Stiftelse för teoretisk fysik.”

The authors have no conflicts to disclose.

F.N.R. and N.S.K. derived the extended micromagnetic model, performed the calculations and related analysis, and contributed equally to this work. A.B.B. and L.D. contributed the mathematical background for the governing functional. F.N.R., N.S.K., C.M., and S.B. contributed to writing the manuscript. F.N.R. and S.B. initialized the work independently.

Filipp N. Rybakov: Conceptualization (equal); Investigation (equal); Writing – original draft (equal). Nikolai S. Kiselev: Investigation (equal); Visualization (equal); Writing – original draft (equal). Aleksandr B. Borisov: Investigation (supporting); Writing – review & editing (supporting). Lukas Döring: Investigation (supporting). Christof Melcher: Investigation (equal); Writing – review & editing (equal). Stefan Blügel: Conceptualization (lead); Writing – review & editing (equal).

The data that support the findings of this study are available within the article.