Transition between distinct hybrid skyrmion textures through their hexagonal-to-square crystal transformation in a polar magnet

Magnetic skyrmions, topological vortex-like spin textures, garner significant interest due to their unique properties and potential applications in nanotechnology. While they typically form a hexagonal crystal with distinct internal magnetisation textures known as Bloch- or Néel-type, recent theories suggest the possibility for direct transitions between skyrmion crystals of different lattice structures and internal textures. To date however, experimental evidence for these potentially useful phenomena have remained scarce. Here, we discover the polar tetragonal magnet EuNiGe3 to host two hybrid skyrmion phases, each with distinct internal textures characterised by anisotropic combinations of Bloch- and Néel-type windings. Variation of the magnetic field drives a direct transition between the two phases, with the modification of the hybrid texture concomitant with a hexagonal-to-square skyrmion crystal transformation. We explain these observations with a theory that includes the key ingredients of momentum-resolved Ruderman–Kittel–Kasuya–Yosida and Dzyaloshinskii-Moriya interactions that compete at the observed low symmetry magnetic skyrmion crystal wavevectors. Our findings underscore the potential of polar magnets with rich interaction schemes as promising for discovering new topological magnetic phases.

There could be different reasons for the unexpected finite   in phases I and IV.For example, Equation 1, which is used to calculate    , includes the assumption that the normal Hall coefficient R0, which is also included in the topological Hall coefficient, is unchanged across different magnetic phases.However, this hypothesis is only valid if the helical structure has a long period.In other words, in the presence of weak DMI, the magnetic order cannot modify the electronic structure dramatically, and the normal Hall coefficient, R0, is not affected by variation in the magnetic order.Since the itinerant electrons are undoubtedly coupled to the magnetism in this system, we thus conclude pronounced changes likely take place in the electronic structure between the different phases.As a result, we can expect R0 to have a Hdependency between the ground state phase I and the saturated regime.Therefore, the standard analysis of the topological Hall signal based on a H-independent R0 cannot apply here.This likely leads to the formation of an uncompensated topological signal in phase I and IV.
On the other hand, if the conventional analysis approach for   data would in fact be appropriate, the additional Hall signal observed in phases I and IV may then have its origin in an extrinsic mechanism.A possible explanation is skew scattering from topologically nontrivial features of domain walls separating different magnetic domains S1-S3 .Since phases I and IV are determined to naturally form multi-domain states according to our SANS measurements, topological defects may be expected to exist.Thus, the magnitude of skew scattering-related Hall resistivity may exceed the intrinsic topological Hall effect.Furthermore, such scattering events are less likely to occur in the single-domain multi-Q phases II and III

Supplementary Note III : Spin Polarisation Factor
Here we attempt a quantitative analysis of the THE, whereby we assume the magnetic order in phases II and III correspond to skyrmion phases.According to the continuum approximation, the emergent magnetic field due to a skyrmion is given by Bem = -Φ0/λSk 2 = -(h/e) nSk, where λSk represents the magnetic period of the skyrmion lattice (SkL), and nSk is the skyrmion density.
From Fig. 3d, we obtain nSk ~ 0.22 nm -2 in Phase II and nSk ~ 0.184 nm -2 in Phase III, which yields Bem ~ -900 T and -750 T, respectively.Such enormous effective fields can be related empirically to the topological Hall resistivity through the normal Hall coefficient R0, and the effective spin polarisation P of conduction electrons given by S4,S5 : After roughly extracting topological Hall resistivity    ~ 0.02 μΩ.cm and R0 ~ 0.03 μΩ.cm/T from Fig. 2c for both Phases II and III, we estimate P = 7 × 10 -4 and 8 × 10 -4 , respectively.
In general, P is relatively strong for d-d coupling in transition metal compounds (e.g. the B20 family).For example P ≃ 0.1 for MnSi S4 , P ~ 0.25-0.38 for MnSi under pressure S5 , and for slightly-doped Mn1-xFexSi, P ~ 0.3-0.45S6 .However, it is expected to be small for rare earth systems with moderate f-d coupling, for example P = 0.07 for Gd2PdSi3 29 .We performed a spin density functional theory (SDFT) calculation of the density of states (DOS) near the Fermi level of EuNiGe3, and estimate a small value of P ~ 0.015.The small estimated size of P can be attributed to the dominant contribution of Ni and Ge to the DOS compared with Eu.While the estimated value for P from SDFT is larger than those estimated above from the experimental data, the small sizes of both theoretical and experimental P estimates prevent a meaningful quantitative comparison.Nonetheless, the SDFT result suggests the small contribution of Eu electron bands at the Fermi level, and is thus generally consistent with a small value of P, and consequently the small topological Hall resistivity observed in this compound.

Supplementary Note IV : Experimental data for H parallel [100]
In Supplementary Fig. 3 we show experimental data obtained for progressively stronger magnetic fields applied along the [100] direction, i.e. perpendicular to the polar axis.
Supplementary Fig. 3a shows the H-T phase diagram constructed from temperature and H- Turning the H-dependence of the observable magnetic scattering at 1.9 K, as seen in Supplementary Fig. 3b we do not observe the magnetic Q-vectors to display a drastic rearrangement up to saturation.Instead, a smooth variation of the overall scattering intensity is observed as the field is increased (Supplementary Fig. 3c), with a small increase seen in the low field region.Notably there is no obvious feature in the high field portion related to the phase line determined by susceptibility shown in Supplementary Fig. 3a.Therefore, further work is needed to clarify the physical origin of this feature.The variation in SANS intensity in the low field region may have its origin in a differing behaviour of the unobserved single Qdomains which have Q-vectors nearly perpendicular to H, or the deformation of the observed ground state into a conical-like structure.Again, further experiments are needed to elucidate the fate of all single-Q domains under in-plane fields.Thus, here we focus on the behaviour of the two observable single-Q domains.Since they both survive until saturation, this provides evidence that the modulation of the ground state phase I has a mainly helical character, since the stability of helical structures benefit from the generally enhanced susceptibility for fields applied at angles to the helical plane.This contrasts with the expectation for a cycloidal modulation with moments rotating in the plane containing Q, and for which the susceptibility is enhanced normal to the cycloidal plane.In the cycloidal modulation scenario, the observed single-Q domains could be expected to be destabilised (or polarised) well before saturation, and instead it would be the unobserved single-Q domains that would be more likely to survive until saturation.Therefore, the data at hand provide support that the ground state phase I in EuNiGe3 has a significant helical character.

Supplementary Note V : Field tilted SANS measurements
In the field-tilted SANS measurements performed in Phase II we studied the effect on the hexagonally coordinated multi-Q domains when the magnetic field of 2.6 T was applied at increasing angle to the polar axis within the (110) plane.Experimentally, the sample was always zero-field cooled (ZFC), before the field applied at progressively larger tilt angles up to 15° from the c-axis.Due to the mutually orthogonal alignment of the two domains for H || [001], the tilting of the field away from c in the (110) plane naturally breaks the symmetry for the multi-Q domain formation, and is thus is expected to affect the stability of the two domains differently.
The experimental configuration for these unpolarised SANS measurements is illustrated in Supplementary Fig. 4a, with the magnetic field tilt angle away from the [001] axis labelled as ø.For ø = 0°, Supplementary Fig. 4b shows we observe the 12-spot pattern due to the twodomain triple-Q state, where all Q-vectors are the same as described in Fig. 3j from the main text.For ø = 15° however, Supplementary Fig. 4b shows that indeed one of the hexagonal multi-Q domains is indeed preferentially stabilised compared with the intensity due to other domain severely suppressed.The ø -dependence of the overall SANS intensity due to each of the two domains is shown in Supplementary Fig. 4d.The inverse correlation between the ø -where     and     stand for the polar-and radial-type DMIs;     ⊥   and     ||   .We set 3  = 0 (11) 5  = 0 (13) with  1  = 0.9 and  1  = √ 1 − ( 1  ) 2 .Finally, we set (   ,    ,    ) as follows: where α = 1.05, α' = 1.35, β = 0.12 / √5, and γ = 0.06 / √5.In the calculations, we consider all the ordering wavevectors that are symmetry-equivalent to Q1-Q6.
The role of the above model parameters is as follows.κ1 and κ2 represent the interaction ratios to the interaction channel at Q1; the interaction in the Q1 channel is dominant, while those in Q2-Q6 are subdominant.This choice of the parameters leads to the spiral state with Q1 at zero field.In Eqs. ( 8)-( 15), we set the model parameters in terms of the DMI so that the helicity of the spiral and skyrmion phases observed in experiments is reproduced.In particular, the radial component in Eq. ( 7) plays an important role in determining the helicity in Phase I and Phase III.Meanwhile, the DMI plays a lesser role in the stabilization and the choice of the helicity detail how the helicity of each skyrmion type is quantified and, in addition, how the spatial distribution of helicity varies amongst the different skyrmion types.
To begin, we elaborate on how the real-space spin textures of Néel, Bloch, and hybrid skyrmions are generated in the absence of genuine atomistic models for magnetic textures in EuNiGe3.The real-space schematics of skyrmions are created by employing a Fourier-space description, which involves a summation over different magnetic components: Here, m(r) represents the magnetization vector as a function of spatial coordinates (r), where   is a parameter describing the magnetization contributed by the i th Fourier component.
These components have their orientations and magnitudes determined using a basis scheme indexed by j.Each Fourier component resides at a wavevector   and is assigned a phase parameter   within the basis.Magnetic skyrmion lattices are generated by constructing a two-dimensional grid of points with varying r.Equation 25can then be evaluated once the relevant magnetic wavevectors and their corresponding magnetic Fourier components are known.
In our approach, we adopt a scheme in which we set the in-plane magnetic components to have a phase of   =0, while the out-of-plane directions are given a 90° phase shift between each wavevector.Through such a scheme, spatial variations in the magnitude of the magnetic moments are minimized, which is in alignment with the expectation of our numerical simulations.
Additionally, provided that a 90° phase difference is maintained between the in-plane and outof-plane   values, a global phase shift, denoted as   , results in a translation of the skyrmion lattice when placed on a rectangular grid.Consequently, the positions of skyrmion cores can be arbitrarily adjusted by varying the value of   .Here, we choose   = 0 such that a skyrmioncore is located at  = 0.
Armed with models for the various real-space spin textures, we now move on to quantifying the overall helicity distribution for each type of skyrmion by evaluating χ at over 10000 evenly spaced locations within the magnetic unit cell.The calculated helicity distributions are shown in Supplementary Fig. 7b, d, f, and h, along with the schematic depictions of the various skyrmion types investigated in Supplementary Fig. 7a, c, e, and g.In the insets for Supplementary Figs.7b, d, f, and h, we report estimates for the mean helicity ̅ = ∑   χ  /  , and standard deviation  = √∑   (  − ̅ ) 2 /( − 1)

𝑖
, where in each equation   is the frequency of the  th helicity bin,   is the mid-point of the helicity bin, and  is the number of locations in the unit cell over which the helicity was calculated.
For the archetypal Néel and Bloch skyrmions, the histograms shown in Supplementary Figs.7b and d are dominated by a single peak near either χ = 0 or π/2, as expected for skyrmions with isotropic distributions of helicity.Small contributions to χ slightly away from these peaks are observed, which we attribute to the breakdown of radial symmetry towards the edges of the unit cell and the geometry of the skyrmion lattice (Supplementary Figs.7a and c).This difference in spatial profile between an azimuthally symmetric skyrmion profile (denoted by either a circle or ellipse for Phase II) and the magnetic unit cell plays an observable role in the χ distribution at the extremities of the latter.The result is minor contributions to the histogram at values deviating slightly beyond the range expected for a single isolated skyrmion in a ferromagnetic background.For the Néel and Bloch type skyrmion histograms, the result is a deviation from a truly singular helicity value, a value of ̅ slightly different to either 0 or π/2, and a small, yet finite value for .
Nonetheless, when comparing between the calculated χ distributions of the different skyrmion types in Supplementary Figs.7b, d, f, and h, the histograms for Bloch and Néel skyrmions display a clearly much narrower distribution compared with those of hybrid skyrmion Phases II and III in EuNiGe3.Instead, Phases II and III each display a significantly extended χ distribution due to their hybrid nature.The distribution for Phase II is described by ̅ = 1.428 and  = 0.105, which is skewed towards χ values approaching the pure Bloch value of π/2, but the existence of a competing Néel-type winding leads to a distribution of χ values clearly both less and broader than the pure Bloch case.This leads us to describe the hybrid skyrmions in Phase II as having a 'weak' Bloch character.On the other hand, the χ distribution for Phase III is described by ̅ = 0.581 and  = 0.223, which describes a broader distribution that is skewed towards values below π/4, which is the mid-point between pure Néel and Bloch helicities.Therefore, we describe the hybrid skyrmions in Phase III as displaying a 'weak' Néeltype character.
Next, we next investigate the spatial distribution of helicity within the Néel, Bloch, Phase II and Phase III skyrmions.In the lower parts of main text Fig. 1a, b, e, and f, the variation of the colour in the schematics for Phases II and III already show the calculated helicity may vary with both radial distance from the skyrmion centre, and azimuthal angle within the plane of the two-dimensional skyrmion lattice.In Supplementary Fig. 8, we highlight this finding further by presenting representative one-dimensional (1D) cuts of the magnetisation texture that are dependent magnetic susceptibility measurements.Below saturation, the constructed phase diagram displays two broad regions, with the low field, low temperature portion demarcated by a kink in the susceptibility curves (data not shown).Unpolarised SANS data obtained for H || [100] were measured at the SANS-I beamline, PSI, and are shown in Supplementary Figs.3b and c.In the chosen experimental setup with H perpendicular to the neutron beam, we could observe only two of the expected single-Q domains described by Q1 in phase I, namely those with Q-vectors aligned closest to both H and the horizontal axis.The other two domains that would be expected to show scattering in the top and bottom portions of the detector were shadowed by the more restricted neutron beam access of the horizontal field magnet compared with the H parallel to beam geometry used for the H || [001] study discussed in the main text.Therefore, as seen in Supplementary Fig.3b, just two single-Q domains are observed after zero-field cooling, and in addition with different relative intensities presumable related to the different populations of single-Q domains in the crystal.
hexagonal Bloch-type skyrmion lattice, c Phase II of EuNiGe3, and d Phase III of EuNiGe3 are shown in the upper parts of each panel.The colour map encodes the   component, while arrows indicate the direction of in-plane magnetisation.The white bars with roman numerals indicate cuts along high symmetry directions of the skyrmion lattice, and over which the spatial dependence of the helicity is calculated as shown in panels in e to h.The lower parts of panels a to d show the relevant 2D skyrmion projections originally presented in Figs.1a, b, e and f, also overlaid with the same cut bars as shown in the upper panels.The colour of the discs encodes the spatial distribution of the calculated helicity , and in accord with the relevant colourbar at the bottom of the figure.e The spatial dependence of the calculated helicity along the nearest-neighbour (cut i) and next-nearest neighbour direction (cut ii) of the Néel-type skyrmion lattice.The arrows indicate the local direction of the magnetisation, and they are coloured according to their   component (the colourbar above panels a to d applies).The colour of the strip underneath the arrows encodes the spatial dependence of calculated helicity , according to the colourbar at the bottom of the figure.f The same as for e but for a model Bloch type skyrmion lattice.g and h respectively show the cuts along the high symmetry skyrmion lattice directions for Phases II and III of EuNiGe3.