Field and anisotropy driven transformations of spin spirals in cubic skyrmion hosts

Within the phenomenological theory introduced by Dzyaloshinskii Dz64 the magnetic energy density of a non-centrosymmetric ferromagnet with spatially dependent magnetization $\mathbf{M}$ can be written as

where $A$ and $D$ are coefficients of exchange and Dzyaloshinskii-Moriya interactions (DMI); $\mathbf{H}$ is an applied magnetic field; $x_{i}$ are the Cartesian components of the spatial variable. $w_{D}$ is composed of Lifshitz invariants $\mathcal{L}^{(k)}_{i,j}=M_{i}\partial M_{j}/\partial x_{k}-M_{j}\partial M_{i}/\partial x_{k}$ that are energy terms involving first derivatives of the magnetization with respect to the spatial coordinates. In the following, all calculations will be done for cubic helimagnets with $w_{D}=\mathbf{m}\cdot\nabla\times\mathbf{m}$ although the results may be applied for magnets with other symmetry classes JETP89 including different combinations of Lifshitz invariants. $W_{0}(\mathbf{M})$ includes only basic interactions essential to stabilize skyrmion and helical states and specifies their most general features attributed to all chiral ferromagnets. Next we discuss the physics governed by these dominant interactions. The effects of anisotropy $W_{a}(\mathbf{m})$, usually representing a smaller energy scale, will be discussed later.

In chiral magnets the Dzyaloshinskii-Moriya interactions $w_{D}$ play a crucial role in destabilizing the homogeneous ferromagnetic arrangement and twisting it into a helix. At zero magnetic field the helices are single-harmonic modes forming the global minimum of the functional $W_{0}(\mathbf{M})$ Dz64 , $\mathbf{M}=M_{s}\left[\mathbf{n}_{1}\cos\left(\mathbf{q}\cdot\mathbf{r}\right)+\mathbf{n}_{2}\sin\left(\mathbf{q}\cdot\mathbf{r}\right)\right]$. $\mathbf{n}_{1}$, $\mathbf{n}_{2}$ are the unit vectors in the plane of the magnetization rotation orthogonal to the wave vector $\mathbf{q}=1/(2L_{D})$ ($\mathbf{n}_{1}\perp\mathbf{n}_{2};\mathbf{n}_{1}\perp\mathbf{q};\,\mathbf{n}_{2}\perp\mathbf{q}$). $L_{D}$ is proportional to the ratio of the counter-acting exchange and Dzyaloshinskii constants and introduces a fundamental length characterizing the periodicity of chiral modulations in chiral magnets, $L_{D}=A/D$.

The helical modulations have a fixed rotation sense determined by the sign of Dzyaloshinskii-Moriya constant $D$ and are continuously degenerate with respect to propagation directions of the helical modulations in space. An applied magnetic field lifts this degeneracy and stabilizes two types of one-dimensional modulations: cones and helicoids with propagation directions parallel and perpendicular to the applied magnetic field respectively.

For helicoids, analytical solutions for the polar angle $\theta(x)$ of the magnetization written in spherical coordinates, $\mathbf{M}=M_{s}\left(\sin\theta(x)\cos\psi,\sin\theta(x)\sin\psi,\cos\theta(x)\right)$, are derived by solving a pendulum equation $Ad^{2}\theta/dx^{2}-H\cos\theta=0$. Such solutions are expressed as a set of elliptical functions Dz64 and describe a gradual expansion of the helicoid period with increased magnetic field (see the set of angular profiles $\theta(x)$, e.g., in Fig. 4.1. of Ref. thesis, ). In a sufficiently high magnetic field $H_{H}$ the helicoid infinitely expands and transforms into a system of isolated non-interacting 2$\pi$-domain walls (kinks) separating domains with the magnetization along the applied field Dz64 ; JMMM94 . The dimensionless value of this critical field is $h_{H}=H_{H}/H_{D}=\pi^{2}/8=0.30843$ with $H_{D}=D^{2}/AM$. The azimuthal angle $\psi$, on the contrary, is fixed by the different forms of the Lifshitz invariants. In particular $\psi=\pi/2$ for Bloch helicoids in cubic helimagnets and $\psi=0$ for Neel helicoids or cycloids in polar skyrmion hosts with C_nv symmetry such as the lacunar spinels GaV₄S₈ Kezsmarki15 ; Ruff15 , GaV4Se8 Butykai17 , GaMo4S8 Butykai19 or VOSe2O5 Kurumaji17 . In this latter class of materials, the polar axial symmetry gives rise to Lifshitz invariants different from that of the cubic helimagnets (Eq. (1)) and the modulation vectors are restricted to the plane perpendicular to the high-symmetry axis.

The conical state combines the properties of the ferromagnetic and the helical states as a compromise between the Zeeman and DM energies. The conical spirals retain their single-harmonic character with $\psi=z/2L_{D}$ and $\quad\cos\theta=|\mathbf{H}|/2H_{D}$. The magnetization component along the applied field has a fixed value $M_{\bot}=M\cos\theta=MH/2H_{D}$ and the magnetization vector $\mathbf{M}$ rotates within the surface of a cone. The critical value $h_{d}=2H_{D}$ marks the saturation field of the conical state. For the functional $W_{0}(\mathbf{M})$ in Eq. (1) the conical phase is the global minimum over the whole range of applied fields where modulated states exist ($0<h<h_{d}$), whereas helicoids and skyrmion lattices can exist only as metastable states.

In this section we present solutions for spiral states, which we obtained by exact numerical minimization of the energy functional, given by Eq. (1), including the cubic anisotropy $K_{c}$ (2). For this purpose, we consider the dimensionless anisotropy constant, $k_{c}=K_{c}A/D^{2}$. The magnetic field $\mathbf{H}$ is applied along the high-symmetry directions [110], [111], and [001]. The energy density is measured in non-dimensional units, $w(\mathbf{m})=W(\mathbf{M})/H_{D}M$.

a. First- or second-order phase transitions between the conical and the ferromagnetic state.

In the present section, we consider the solutions for the conical state in the field applied along [001] or [111] for both signs of $k_{c}$. As it will become evident later, first-order phase transitions occur when the magnetic field is applied along the easy anisotropy axes. In this case, on its way to saturation, the magnetization is forced to embrace hard anisotropy axes, as shown in Fig. 1 (a) and (c). This is an energetically unfavourable configuration, that forces a discontinuous transition of the magnetisation towards the magnetic field. On the contrary, if the magnetic field points towards a hard anisotropy axis, the minima of the cubic anisotropy gradually align along the field and underlie the continuous second-order phase transition between the conical and the ferromagnetic state, as illustrated by Fig. 1 (b) and (c).

In the following, we carry out a comprehensive analysis of the magnetization rotation in the conical phase by considering the surface energy plots shown in first row of Fig. 2, and the three-dimensional polar plots of the cubic anisotropy energy depicted in the second row of Figs. 2. The different topologies of these surface plots underly the particular types of the phase transitions.

For $k_{c}>0$ and $h=0$ (Fig. 2 (a)) the equilibrium states of the homogeneous state correspond to the easy axes of cubic anisotropy oriented along the $\langle 001\rangle$ crystallographic directions. Maxima of the cubic anisotropy energy correspond to the hard $\langle 111\rangle$ directions. Let us consider the helical mono-domain state where the magnetization for $h=0$ rotates in the plane (001). While rotating, the magnetization leaves one energy minimum corresponding to $\langle 001\rangle$ directions and, rotating through the saddle point between hard axes $\langle 111\rangle$, gets into another energy minimum along a $\langle 001\rangle$ direction. The trace of the magnetization in the conical phase is shown by the yellow line in Fig. 2 . When applying a magnetic field along the [001] axis, the magnetization in the conical phase is forced to rotate in the vicinity of hard $\langle 111\rangle$ axes, as shown in Fig. 2 (b) for $h=0.2$. This underlies the subsequent first-order phase into the ferromagnetic state, which becomes evident, by e.g. the magnetization jumps in the magnetization curves shown in Fig. 1 (e). In Fig. 1 (d) we plot the difference between the energy densities of the conical spiral and the ferromagnetic states. In this plot the point $B$ for $k_{c}=0.1$ stands for the first-order phase transition that takes place when the energy difference between the two phases vanishes. Logically the same situation should also occur for $h||[111]$ and $k_{c}<0$, where three hard anisotropy axes $\langle 001\rangle$ impose distortions to the conical phase and underlie the magnetization jump into the ferromagnetic state.

For $k_{c}<0$ and $h=0$ (Fig. 2 (c), (d)), the equilibrium states of the cubic anisotropy energy correspond to the $\langle 111\rangle$-type crystallographic directions. Thus, the maxima of the cubic anisotropy are the hard $\langle 001\rangle$-type directions, and if a magnetic field is applied along a $\langle 001\rangle$ axis, the angle between the easy directions and the field is $70.6^{\circ}$. In this case, the magnetization in the conical phase rotates so that it sweeps the easy directions $\langle 111\rangle$ (Fig. 1 (b)). Even at zero field the magnetic structure is conical, thus, it has a non-zero uniform magnetization (see magnetization curves in Fig. 1 (g)). However, this low-field conical phase is in reality a metastable solution overshadowed by the oblique spiral states with $\mathbf{q}$-vectors along the easy $\langle 111\rangle$-type directions. In a finite magnetic field the global minima of the cubic anisotropy gradually approach the field direction and thus display a second-order phase transition to the ferromagnetic state (Fig. 2 (d)). The same principle is valid also for $\mathbf{h}||[111]$ and $k_{c}>0$, where the angle between the field direction and the easy $\langle 001\rangle$ axes is $54.7^{\circ}$, whereas for the hard $\langle 111\rangle$ axes the angle is $70.6^{\circ}$. In this case the topology of the corresponding polar plots would allow the gradual alignment of the magnetization along the hard $\langle 111\rangle$ direction with increasing magnetic field.

b. Conical domains in the ferromagnetic state.

As a consequence of the first-order phase transition that occurs for either $k_{c}>0$ and $h||[001]$ or $h||[111]$ and $k_{c}<0$, domains of the conical phase may persist in the field polarised state adopting the configuration illustrated in Fig. 3. Such domains exhibit both types of the rotational sense at their boundaries. At the left side of the depicted conical droplet (Fig. 3 (a)-(c)) the magnetization in the conical phase adjusts to the homogeneous state with the ”correct” twist supported by DMI. This leads to energy gain as shown by the negative energy density difference in the $xy$ plane, i.e., for a cross-section perpendicular to the magnetic field $\mathbf{h}$ and the conical $\mathbf{q}$-vector, shown in Fig. 3 (d)). At the right side of the conical droplet, the magnetization acquires the ”wrong” in-plane twist, which results in an energetic ”penalty”, as shown in Fig. 3 (d).

A similar behaviour is also found along the $z$ direction, which is parallel to $\mathbf{h}$ and $\mathbf{q}$. As illustrated in Fig. 3 (e), (f) the energy density of the conical droplet exhibits alternating negative and positive contributions corresponding to ”correct” and ”wrong” senses of the magnetization rotation at the boundary to the field polarised state. When the magnetic field is increased from the conical to the ferromagnetic phase, the conical droplet may either continuously expand for $h<h_{B}$, as seen in Fig. 1 (d) or gradually shrink for $h>h_{B}$, a feature that underlies the kinetic character of such first-order phase transition. Remarkably, when the magnetic field is decreased from the field polarized to the conical state, the energy density shown in Fig. 3 (f) leads to an ”in-phase” merge of neighboring conical droplets that initially may be formed with some phase shift with respect to each other. Then, the region with the positive energy density at the boundary of one conical droplet is annihilated by the corresponding part with the negative energy of another droplet, a feature that underlies in-phase ”zipping” and ”unzipping” of conical droplets.

In this section we present qualitative results on the spiral reorientation processes that occur for $k_{c}<0$ and for different directions of the applied magnetic field (Fig. 2 (e)). Numerically exact solutions for the spiral states are obtained by direct minimization of the energy functional, given by Eq. (1), including $K_{c}$ in Eq. (2), and are discussed elsewhere chapter . Since $k_{c}>0$ favours easy axes of $\langle 111\rangle$ type, the ground state at zero field consists of domains of spiral states with the corresponding directions of the $\mathbf{q}$-vectors. Fig. 2 (d) shows a schematic layout of the $\mathbf{q}$-vectors in stable and metastable spiral states with respect to the field and the crystallographic directions.

For $\mathbf{H}\|[001]$ (first panel in Fig. 2 (e)), all four equivalent spiral states with $\mathbf{q}\|\langle 111\rangle$ are at an angle of $35.2^{\circ}$ with respect to the field. As also discussed in Ref. chapter, , with increasing magnetic field, these wave vectors first rotate gradually and then jump by a first-order phase transition towards the field direction.

For $\mathbf{H}\|[111]$ (second panel in Fig. 2 (e)), only the spiral that is co-aligned with the magnetic field is stable, whereas the spirals along the other equivalent wave vectors $\mathbf{q}\|\langle 111\rangle$ become metastable. With decreasing magnetic field only one spiral (or conical) domain with the $\mathbf{q}$-vector along the field persists and at zero field the other spiral domains are not recovered.

For $\mathbf{H}\|[110]$ (third panel in Fig. 2 (e)), there are two sets of $\langle 111\rangle$ directions composing smaller and larger angles with respect to the field. Thus the field-driven transition into the conical phase is a two-step process: the energetically unfavorable spiral states with the angle $90^{\circ}$ reorient first, whereas the stable spirals with the angle $27.35^{\circ}$ gradually rotate and eventually flip along the field. With decreasing magnetic field the transitions at both critical fields should occur in a reproducible manner, possibly with some hysteresis.

For $\mathbf{H}||[110]$, independently on the sign of the cubic anisotropy constant $k_{c}$, the $\mathbf{q}$-vector of the conical state is surrounded by two hard and two easy anisotropy axes $\langle 111\rangle$ and $\langle 100\rangle$. This rosette of anisotropy axes leads to the characteristic elliptical distortion of the conical state, which moreover changes when the magnetization rotation occurs in the exterior or interior of the rosette. In particular for $k_{c}<0$ and low field values, the trajectory of the rotating magnetization is elongated in-plane along the hard $\langle 100\rangle$ axes (Fig. 4 (a), (b)). With increasing magnetic field the magnetization approaches the rosette, penetrates the rosette’s interior and changes its elliptical distortion being now elongated along the easy $\langle 111\rangle$ axes (Fig. 4 (c), (d)).

For $k_{c}>0$, the same process takes place, but with inverse switch of the elliptical distortion (Fig. 5). The ratios of the ellipse major and minor axes for both signs of the anisotropy constant are plotted as a function of the magnetic field for different values of $k_{c}$ in Fig. 4 (e) and Fig. 5 (e). From these curves as well as from the magnetization curves (Fig. 4 (f)), it becomes clear that the phase transition between these two conical phases is more abrupt for the negative values of $k_{c}$. This might be explained by the different angles the easy and hard anisotropy axes compose with respect to the field – $35.3^{\circ}$ for $\langle 111\rangle$ and $45^{\circ}$ for $\langle 100\rangle$. Indeed, when $k_{c}<0$, the rotating magnetization tries to avoid two hard axes [100] and [010] and this results in the jump of the magnetization from the rosette’s exterior to its interior. For $k_{c}>0$ [100] and [010] axes are easy anisotropy directions. Thus, although the sense of ellipticity changes, no jump of the magnetization is observed. An important point is that the phase transition between the two different conical states leaves the periodicity of the spiral state almost unaffected. In fact, the period changes just by a few percent for low values of the cubic anisotropy (less than $4\%$ in Fig. 4 (g) for $k_{c}=-0.2$), and this change is hardly detectable by e.g. neutron scattering experiments. The field dependence of $dm_{h}/dh$ (Figs. 4 inset of (f) and Fig. 5 (f)), on the contrary, provides a clear hint for this transition as it shows a $\lambda$-shaped spike only for $k_{c}<0$.

In the following we will show that spiral reorientation processes akin to those at $H_{c1}$ may intervene also near the saturation field $H_{c2}$. Jumps of the spiral $\mathbf{q}$-vectors are less expected although they may be realized for relatively moderate values of the cubic anisotropy constant (Fig. 6). Indeed, the conical phase co-aligned with the field is considered to be the most favorable state that benefits from the Zeeman interaction. In the following, we show theoretically that a tilted spiral state may originate from the interplay between Zeeman and cubic anisotropy interactions, which is generic to chiral magnets, and thus may be important for explaining the experimental data of specific cubic helimagnets. The spiral jumps are also sensitive to the direction of the applied magnetic field and the value of $k_{c}$, since relatively large values of $k_{c}$ are needed to implement particular reorientations.

As an illustrative example, we consider the case with $\mathbf{H}||[111]$ and easy $\langle 111\rangle$ axes of the cubic anisotropy (Fig. 6). The $\mathbf{q}$-vector undergoes a first-order transition from the conical state (Fig. 6 (a),(b)) into a tilted spiral state with almost constant tilt angle for a range of magnetic field values (Fig. 6 (c),(d)), with both the tilt angle and the field range dependent on $k_{c}$. The sketches in Fig. 6 (c)-(d) are instructive to deduce the nature of this tilted spiral state. Indeed, in an increasing magnetic field, the magnetization in the conical state sweeps only two $\langle 111\rangle$ anisotropy axes (Fig. 6 (a), (b)). In tilted state, however, by sacrificing the Zeeman energy, the magnetization may try to embrace all four $\langle 111\rangle$ axes (Fig. 6 (c), (d)). This leads to the maximum tilt angle $70.6^{\circ}$, reached for the higher $k_{c}$ values (Fig. 6 (e)). This phase transition is accompanied by characteristic jumps in the magnetization curve. An interesting point is that such a reorientation phase transition requires a threshold value of $k_{c}$. For $k_{c}=-0.1$ not even a hint of the metastable tilted spiral becomes apparent in the angle-dependent energy density (like that in Fig. 6 (f)).

Conceptually the same reorientation process may occur also for $\mathbf{H}||[110]$.

Remarkably, the internal properties of spiral states are directly imprinted on the functionalities of chiral skyrmions. Indeed, isolated skyrmions acquire a magnetization distribution according to a surrounding ”parental” or host state. In a state that is homogeneously magnetized along the field, axially-symmetric skyrmion particles may form that underlie a repulsive inter-skyrmion potential LeonovNJP16 . On the other hand, when surrounded by a conical phase with the wave vector along the skyrmion axes, skyrmions develop a three-dimensional structure: the central core region nearly preserves the axial symmetry, whereas the domain-wall region, that connects the core with the embedding conical state acquires a crescent-like shape. This asymmetric profile of the cross-section forms a screw-like modulation along the skyrmion axis, matching the rotating magnetization of the conical phase at all cross-sections LeonovJPCM16 . To reduce the energy excess necessarily related to such a rotation, isolated skyrmions tend to form compact clusters with the unimpeded magnetization rotation in their interior.

Experimentally, clusters of skyrmions have been observed in thin (70 nm) single-crystal samples of Cu₂OSeO₃ using transmission electron microscopy Loudon18 . In a bulk Cu₂OSeO₃ for $\mathbf{H}||[001]$, skyrmion clusters surrounded by the conical state exist up to the field $H_{c2}$ of the first-order phase transition between the cone and the FM state (Fig. 2 (a), (b)). With the onset of the homogeneous state, the inter-skyrmion attraction turns into repulsion and the constituent skyrmions disperse. If the number of isolated skyrmions is sufficiently high, they may form a SkL with an equilibrium period that exists above the $H_{c2}$ line, as observed experimentally Bannenberg2019 ; chacon2018 . However, if the number of skyrmions is low, they keep moving away from each other and the skyrmion signal in the experiment should fade away with time. Thus, the following questions related to the skyrmion structure and skyrmion-skyrmion interaction arise: (i) What would be the impact of the elliptically distorted conical states for $\mathbf{H}||[110]$ (Fig. 4, 5) or the oblique spiral states (Fig. 7)? Does the skyrmion attraction become anisotropic in these cases? (ii) How does the structure of skyrmion clusters distort during the spiral flips in Fig. 6? Do skyrmions avert this sort of spiral jumps or are they instead annihilated via Bloch points?

Recently it was found that skyrmion tubes within the conical phase may orient either along or perpendicular to a wave vector Leonov2018a ; Sohn19 . The first type of skyrmions perfectly blends into the homogeneously saturated state after the cone closes in a strong magnetic field either by the first-order phase transition (Fig. 2 (b)) or by the second-order one (Fig. 2 (d)). The second type of skyrmions perfectly blends into the spiral state and represents two merons with equally distributed topological charge $Q=1/2$. Then, a crossover between the two types of isolated skyrmions occurs for an intermediate value of an applied magnetic field. Thus, the found oblique spiral states (Fig. 6, 7) may underlie new skyrmion states oblique with respect to the field and their controllable switch by the field.

Skyrmion nucleation may occur via different mechanisms specified by the stable and metastable spiral states Bannenberg2019 ; droplets . In particular, the low-temperature stable SkL in the bulk-insulating cubic helimagnet Cu2OSeO3 for $\mathbf{H}||[001]$ may appear via rupture formation of spiral states with wave vectors perpendicular to the field Bannenberg2019 , along [100] and [010]. In this rather conventional mechanism also observed in the first papers on skyrmion visualization yuFeCoSi ; YuFeGe , the pairs of merons behave as if they were free particles. The properties of such skyrmions and merons within the spiral background were recently investigated in Refs. Muller, ; Ezawa, . The essential difference in this case is that the spiral states are metastable solutions, but give rise to a stable SkL.

The experimental findings in Cu₂OSeO₃Bannenberg2019 also hint that skyrmions may appear from the conical and the homogeneous states near the critical field $H_{cr2}$. This may occur via formation of magnetic torons of finite length accompanied by their elongation due to the cubic or exchange anisotropy, which leads to the negative energy density of skyrmion cores and their mutual attraction. Such a homogeneous nucleation mechanism via torons – spatially localized three-dimensional particles composed of a skyrmion filament of finite length cupped by two Bloch points terminating its prolongation – was recently discussed in Refs. leonov2018x, ; Mueller2020, . In this sense, all sorts of the magnetization inhomogeneities provided, e.g., by the domain walls between spiral domains near the flip transition (Fig. 6) or by the superstructure formed from the rotating oblique spiral (Fig. 7), represent a favorable medium for the nucleation and stabbilisation of skyrmions.

Acknowledgements. The authors are grateful to Christos Panagopoulos, Sandor Bordacs, Ulrich Rößler, and Maxim Mostovoy for useful discussions. This work was funded by DFG Priority Program SPP2137, Skyrmionics, under Grant No. KE 2370/1-1. AOL thanks Ulrike Nitzsche for technical assistance.