The chiral biquadratic pair interaction

TL;DR

This paper introduces and analyzes a new chiral biquadratic pair interaction in magnetic systems, derived from microscopic models and studied through first-principles calculations, with implications for complex magnetic structures.

Contribution

The paper derives the chiral biquadratic interaction from microscopic principles and explores its properties in magnetic dimers and extended systems.

Findings

Identified a new chiral biquadratic pair interaction.

Demonstrated its significance in magnetic dimers and spin structures.

Discussed potential impact on magnetic skyrmion lattices.

Abstract

Magnetic interactions underpin a plethora of magnetic states of matter, hence playing a central role both in fundamental physics and for future spintronic and quantum computation devices. The Dzyaloshinskii-Moriya interaction, being chiral and driven by relativistic effects, leads to the stabilization of highly-noncollinear spin textures such as skyrmions, which thanks to their topological nature are promising building blocks for magnetic data storage and processing elements. Here, we reveal and study a new chiral pair interaction, which is the biquadratic equivalent of the Dzyaloshinskii-Moriya interaction. First, we derive this interaction and its guiding principles from a microscopic model. Second, we study its properties in the simplest prototypical systems, magnetic dimers deposited on various substrates, resorting to systematic first-principles calculations. Lastly, we discuss its importance and implications not only for magnetic dimers but also for extended systems, namely one-dimensional spin spirals and complex two-dimensional magnetic structures, such as a nanoskyrmion lattice.

Tables

Table 1. Table 1: Spin moments, magnetic interaction parameters and opening angles of the magnetic ground state for Cr, Mn, Fe, Co and Ni dimers deposited on Pt(111). M 𝑀 M is the spin magnetic moment of one atom in the dimer. The CBI and the DMI are represented by their dominant vector component C y subscript 𝐶 𝑦 C_{y} and D y subscript 𝐷 𝑦 D_{y} , respectively. The biquadratic and bilinear isotropic interactions are given by B 𝐵 B and J 𝐽 J , respectively. All interaction values are in meV. The canting angles are found by minimizing Eq. ( 4 ) with all the interactions ( Δ ​ α 4 ​ s Δ superscript 𝛼 4 s \Delta\alpha^{\mathrm{4s}} ) or keeping only J 𝐽 J and D y subscript 𝐷 𝑦 D_{y} ( Δ ​ α 2 ​ s Δ superscript 𝛼 2 s \Delta\alpha^{\mathrm{2s}} ). The sign of Δ ​ α Δ 𝛼 \Delta\alpha represents the sign of ( 𝐞 ^ 1 × 𝐞 ^ 2 ) y subscript subscript ^ 𝐞 1 subscript ^ 𝐞 2 𝑦 (\hat{\mathbf{e}}_{1}\times\hat{\mathbf{e}}_{2})_{y} , the vector chirality of the magnetic ground state. For comparison, we also list the effective bilinear interactions defined by D y eff = D y − sgn ​ ( J ) ​ C y superscript subscript 𝐷 𝑦 eff subscript 𝐷 𝑦 sgn 𝐽 subscript 𝐶 𝑦 D_{y}^{\mathrm{eff}}=D_{y}-\mathrm{sgn}(J)C_{y} and J eff = J − 2 ​ sgn ​ ( J ) ​ B superscript 𝐽 eff 𝐽 2 sgn 𝐽 𝐵 J^{\mathrm{eff}}=J-2\,\mathrm{sgn}(J)B .

Dimer		$M$ (${\mu }_B$)	$C_y$	$D_y$	$B$	$J$	$\Delta {\alpha }^{2s}$	$\Delta {\alpha }^{4s}$	$D_y^{\mathrm{eff}}$	$J^{\mathrm{eff}}$
Pt(111)	Cr	$3.26$	$2.5$	$8.5$	$-11.7$	$35.8$	$-13°$	$-6°$	$6.0$	$59.2$
	Mn	$4.05$	$-0.6$	$-3.3$	$0.8$	$58.5$	$3°$	$3°$	$-2.7$	$56.9$
	Fe	$3.32$	$2.6$	$-7.3$	$-2.2$	$-43.0$	$10°$	$6°$	$-4.7$	$-47.4$
	Co	$2.12$	$-1.6$	$7.1$	$0.8$	$-76.8$	$-5°$	$-4°$	$5.5$	$-75.2$
	Ni	$0.62$	$0.5$	$0.8$	$-1.1$	$-5.4$	$-8°$	$-10°$	$1.3$	$-7.6$
Pt(001)	Cr	$2.53$	$2.5$	$11.2$	$-9.7$	$-35.3$	$-18°$	$-14°$	$13.7$	$-54.7$
	Fe	$3.24$	$-0.2$	$-9.5$	$-1.5$	$15.0$	$32°$	$28°$	$-9.3$	$12.0$
Ir(111)	Cr	$3.02$	$3.2$	$10.7$	$-12.1$	$29.5$	$-20°$	$-8°$	$7.5$	$53.7$
	Fe	$3.06$	$1.3$	$-14.6$	$-3.6$	$-16.3$	$42°$	$32°$	$-13.3$	$-23.3$
Re(0001)	Cr	$2.18$	$0.4$	$-18.1$	$-3.4$	$-16.4$	$48°$	$40°$	$-17.7$	$-23.2$
	Fe	$2.29$	$0.3$	$0.5$	$0.1$	$-2.3$	$-12°$	$-19°$	$0.8$	$-2.5$

Equations

$$\Omega$$

$$={\Omega }^0+{\Omega }^{\mathrm{soc}}+\sum _p\sum _{k=1}^{p/2}{\Omega }^{p,2k}[\{{\mathbf{S}}_i\}].$$

$$(1d)\to \frac{1}{2}\sum _{i,j}{\mathbf{C}}_{ij}\cdot \left({\mathbf{S}}_i\times {\mathbf{S}}_j\right)\left({\mathbf{S}}_i\cdot {\mathbf{S}}_j\right).$$

$${\mathcal{H}}^{\mathrm{pair}}=\frac{1}{2}\sum _{i,j}\sum _{\alpha ,\beta }{J}_{ij}^{\alpha \beta }{S}_i^{\alpha }{S}_j^{\beta }+\frac{1}{2}\sum _{i,j}\sum _{\alpha ,\beta ,\gamma ,\delta }{B}_{ij}^{\alpha \beta \gamma \delta }{S}_i^{\alpha }{S}_j^{\beta }{S}_i^{\gamma }{S}_j^{\delta }.$$

$$E(\alpha )=J\cos \alpha +D_y\sin \alpha +B{\cos }^2\alpha +C_y\sin \alpha \cos \alpha .$$

$$E({\alpha }_{\min })$$

$$\displaystyle\Rightarrow\Delta\alpha=\frac{$180\text{\,}\mathrm{\SIUnitSymbolDegree}$}{\pi}\,\frac{\mathrm{sgn}(J)D_{y}-C_{y}}{|J|-2B}\quad,$$

$$E_{\mathrm{DFT}}[\{\hat{\mathbf{e}}\}]=\sum _i\sum _{\alpha ,\beta }{K}_i^{\alpha \beta }{\hat{e}}_i^{\alpha }{\hat{e}}_i^{\beta }+\frac{1}{2}\sum _{i,j}\sum _{\alpha ,\beta }{J}_{ij}^{\alpha \beta }{\hat{e}}_i^{\alpha }{\hat{e}}_j^{\beta }+\dots$$

$$E_{\mathrm{cDFT}}[\{\hat{\mathbf{e}}\}]=E_{\mathrm{DFT}}[\{\hat{\mathbf{e}}\}]-\sum _i{\mathbf{b}}_i\cdot \int \mathrm{d}\mathbf{r}{\mathbf{m}}_i(\mathbf{r}).$$

$$\frac{\delta E_{\mathrm{cDFT}}[\{\hat{\mathbf{e}}\}]}{\delta {\hat{\mathbf{e}}}_i}=\frac{\delta E_{\mathrm{DFT}}[\{\hat{\mathbf{e}}\}]}{\delta {\hat{\mathbf{e}}}_i}-{\mathbf{b}}_i|{\mathbf{m}}_i|=0.$$

Full text

{affiliations}

Peter Grünberg Institut and Institute for Advanced Simulation, Forschungszentrum Jülich & JARA, 52425 Jülich, Germany

Department of Physics, RWTH Aachen University, 52056 Aachen, Germany

Introduction

Starting from the seminal work of Heisenberg[1], magnetic materials are often described by bilinear isotropic magnetic interactions, $J_{ij}\,\mathbf{S}_{i}\cdot\mathbf{S}_{j}$. However, a wealth of complex spin-textures were discovered over the last century that called for the enrichment of the original Heisenberg model with various other types of interactions (see e.g. Refs.[2, 3, 4, 5, 6, 7, 8, 9, 10, 11]). The magnetism of 3He is a striking example, being dominated by higher-order isotropic interactions[12] which can be derived from the Hubbard model at half-filling[13, 14, 15] or from Kondo-lattice models[16, 17, 18]. These interactions, such as the biquadratic interaction $B_{ij}\,(\mathbf{S}_{i}\cdot\mathbf{S}_{j})^{2}$ and the related three- and four-site interactions, introduce nonlinear effects into the Heisenberg model. An important consequence is that different spin spirals, characterized by a wavevector $\mathbf{Q}$, can be combined into lower-energy multiple-$Q$-states, as the higher-order interactions invalidate the superposition principle. Prominent examples are the antiferromagnetic $uudd$-state (a $2Q$-state)[19, 20] and the $3Q$-state[21]. Interestingly, this $3Q$-state (also magnetic skyrmions[22, 23] and bobbers[24, 25]) is a noncoplanar magnetic state that hosts interesting Berry-phase physics arising from its non-vanishing scalar spin chirality $\mathbf{S}_{i}\cdot(\mathbf{S}_{j}\times\mathbf{S}_{k})$, such as topological orbital ferromagnetism and Hall effects[26, 27, 28, 29, 30].

The concept of vector spin chirality is embodied by the antisymmetric bilinear Dzyaloshinskii-Moriya interaction (DMI), $\mathbf{D}_{ij}\cdot(\mathbf{S}_{i}\times\mathbf{S}_{j})$[3, 4], which arises due to the combination of spin-orbit coupling and absence of spatial inversion symmetry. The DMI lifts the energy degeneracy of magnetic spirals with opposite vector spin chirality, $\mathbf{S}_{i}\times\mathbf{S}_{j}$, thus stabilizing magnetic structures of well-defined rotational sense, such as chiral spin spirals[31, 32] and magnetic skyrmions[22, 23]. The intricate interplay of higher-order and anisotropic bilinear magnetic interactions generates various magnetic states: conical spin spirals[33] and more complex magnetic structures[30, 34, 35], such as an intricate nanoskyrmion lattice for a monolayer of Fe on the Ir(111) surface[36].

In this work, we utilize a microscopic model combined with first-principles-based simulations to introduce and characterize a new kind of spin-orbit-driven magnetic pair interaction, the chiral biquadratic interaction (CBI). It has the form $\mathbf{C}_{ij}\cdot(\mathbf{S}_{i}\times\mathbf{S}_{j})\,(\mathbf{S}_{i}\cdot\mathbf{S}_{j})$. Like the DMI, this is a unidirectional interaction which is linear in the spin-orbit coupling, and so it is governed by the magnitude and orientation of the CBI vector $\mathbf{C}_{ij}$. We demonstrate that this vector obeys the same symmetry rules as the DMI[4, 37, 38, 39]. Like the isotropic biquadratic interaction, it couples twice a pair of magnetic moments. After systematic investigations on magnetic dimers made of $3d$ elements on various surfaces with strong spin-orbit coupling, namely Pt(111), Pt(001), Ir(111) and Re(0001) surfaces, we find that the CBI can be comparable in magnitude to the DMI. Lastly, we explore the implications of the CBI for magnetic structures in one and two dimensions.

Results

0.1 Systematic microscopic derivation of higher-order interactions

The benefits of studying the properties of the magnetic interactions starting from a microscopic model are well-illustrated by the case of the DMI. Although phenomenological arguments completely determine the form and symmetry properties of the DMI[3], the microscopic analysis of Moriya[4] and later on the intuitive picture proposed by Fert and Lévy[37, 38] have clarified the main ingredients that underpin this interaction. We thus begin by introducing a generic model of the electronic structure of the magnetic material, and then outline how one can systematically extract all kinds of magnetic interactions from the electronic grand potential.

Microscopic model. The microscopic hamiltonian that we consider has three contributions: $\mathcal{H}=\mathcal{H}^{0}+\mathcal{H}^{\mathrm{mag}}+\mathcal{H}^{\mathrm{soc}}.$ Here $\mathcal{H}^{0}$ contains all spin-independent contributions, $\mathcal{H}^{\mathrm{mag}}=\sum_{i}U_{i}\,\mathbf{S}_{i}\cdot\mathbf{\upsigma}$ is the local exchange coupling of strength $U_{i}$ between the magnetic moment $\mathbf{S}_{i}$ on site $i$ and the electronic spin $\mathbf{\upsigma}$, and $\mathcal{H}^{\mathrm{soc}}=\sum_{a}\lambda_{a}\,\mathbf{L}_{a}\cdot\mathbf{\upsigma}$ is the atomic spin-orbit coupling of strength $\lambda_{a}$ on site $a$ between the electron spin and its atomic orbital angular momentum $\mathbf{L}_{a}$. Grouping the spin-dependent terms into $\Delta\mathcal{H}=\mathcal{H}^{\mathrm{mag}}+\mathcal{H}^{\mathrm{soc}}$, it is straightforward to derive a formal power series for the electronic grand potential (see Supplementary Note 1),

[TABLE]

Here $\Omega^{0}$ is the contribution to the grand-canonical potential from the spin-independent $\mathcal{H}^{0}$, and $G^{0}(E)=(E-\mathcal{H}^{0})^{-1}$ is the corresponding retarded Green function. The contributions arising solely from spin-orbit coupling are collected in $\Omega^{\mathrm{soc}}$, and the terms that depend on the magnetic moments are given by $\Omega^{p,2k}[\{\mathbf{S}_{i}\}]$. The Fermi-Dirac distribution for energy $E$ and chemical potential $\mu$ is given by $f(E;\mu)$, and the trace is over all sites, orbitals and spin degrees of freedom.

Diagrammatic rules.

Only a subset of the terms contained in $\Omega^{p,2k}[\{\mathbf{S}_{i}\}]$ are of interest for the purpose of identifying the possible types of magnetic interactions. As detailed in Supplementary Note 1, these can be represented by prototypical diagrams for all kinds of magnetic interactions, from magnetocrystalline anisotropies to pair interactions or many-site interactions. Each diagram contain $p$ vertices connected by $p$ lines. The vertices in a prototypical diagram must correspond to spatially distinct sites, and the lines represent connections between the sites through $G^{0}$. Each diagram contains $2k$ magnetic sites and $p-2k$ spin-orbit sites. A link between two magnetic sites is denoted by a solid line, while a link between a magnetic site and a spin-orbit site is marked by a dashed line. A magnetic site cannot appear consecutively (i.e. a line cannot close on itself), and two consecutive spin-orbit sites (distinct or not) are also excluded. In this work we focus on interactions involving two magnetic sites, and the corresponding prototypical diagrams up to fourth order in the magnetic sites and first order in spin-orbit coupling are shown in Fig. 1a-d.

Prototypical diagrams. It is a simple matter to extract the form of the magnetic interactions from each prototypical diagram, by using the properties of traces of Pauli matrices. The derivations and the forms of the coupling coefficients can be found in Supplementary Note 1. The first and simplest diagram is given in Fig. 1a and translates into the isotropic bilinear exchange interaction $\frac{1}{2}\sum_{i,j}J_{ij}\,\mathbf{S}_{i}\cdot\mathbf{S}_{j}$. Attaching one spin-orbit site to the diagram of Fig. 1a results in the diagram shown in Fig. 1b. This generates the DMI, $\frac{1}{2}\sum_{i,j}\mathbf{D}_{ij}\cdot\left(\mathbf{S}_{i}\times\mathbf{S}_{j}\right)$. The structure of this diagram is identical to the third-order perturbation theory developed by Fert and Lévy[37, 38]. The DMI vector is determined by the properties and geometrical arrangement of the spin-orbit sites, $\mathbf{D}_{ij}=\sum_{a}\mathbf{D}_{ij,a}$. Due to the cross product form, it favors magnetic structures with a definite vector chirality. The next diagram is shown in Fig. 1c and leads to the isotropic biquadratic interaction, $\frac{1}{2}\sum_{i,j}B_{ij}\left(\mathbf{S}_{i}\cdot\mathbf{S}_{j}\right)^{2}$. The diagrams with the same number of lines but connecting either three or four different magnetic sites lead to the isotropic 4-spin 3-site and 4-spin 4-site interactions, respectively (see Supplementary Note 1). Lastly, we find a new kind of magnetic interaction from the prototypical diagram shown in Fig. 1d:

[TABLE]

We name it the chiral biquadratic interaction (CBI), as it is an antisymmetric 4-spin 2-site interaction generated by an additional spin-orbit site. It thus combines the isotropic scalar product $\mathbf{S}_{i}\cdot\mathbf{S}_{j}$ with the chiral coupling $\mathbf{C}_{ij}\cdot\left(\mathbf{S}_{i}\times\mathbf{S}_{j}\right)$ defined by the CBI vector $\mathbf{C}_{ij}=\sum_{a}\mathbf{C}_{ij,a}$, which is generated by the spin-orbit sites. This is our main quantity of interest and its properties will be discussed in detail in this paper. The diagrams with the same number of lines but connecting either three or four different magnetic sites lead to the chiral 4-spin 3-site and 4-spin 4-site interactions, respectively (see Supplementary Note 1).

Symmetry rules. We next study what are the properties of the newly-found CBI vector, $\mathbf{C}_{ij}$, by comparison with those of the DMI vector, $\mathbf{D}_{ij}$. Given a pair of magnetic sites $i$ and $j$ connected with the vector $\mathbf{R}_{ij}$, there are five relevant symmetries, which are illustrated in Fig. 1e-i. Within the picture of the prototypical diagrams, these symmetry operations are a combination of local transformations at each site (e.g. a rotation or a mirroring) and a permutation of the sites. Importantly, symmetry dictates what is the spatial arrangement of the spin-orbit sites around the pair of magnetic sites. One can then relate the diagrams connecting the pair of magnetic sites to each spin-orbit site, noting that the orbital angular momentum operator transforms as a pseudovector, and from this derive the symmetry rules for each magnetic interaction. For the DMI vector, these symmetries lead to the so-called Moriya’s rules[4, 38, 39]. These rules are, for each symmetry operation shown in Fig. 1e-i: (e) $\mathbf{D}_{ij}=0$, (f) $\mathcal{P}_{\hat{\mathbf{n}}}^{\parallel}\,\mathbf{D}_{ij}=0$, (g) $\mathcal{P}_{\hat{\mathbf{n}}}^{\parallel}\,\mathbf{D}_{ij}=0$, (h) $\mathcal{P}_{\hat{\mathbf{n}}}^{\perp}\mathbf{D}_{ij}=0$, and (i) $\mathcal{P}_{\hat{\mathbf{n}}}^{\perp}\mathbf{D}_{ij}=0$ (see Supplementary Note 2). The vanishing components of the DMI vector are those either parallel or perpendicular to $\hat{\mathbf{n}}$, which represents either the rotation axis or the normal to the mirror plane. It follows naturally from comparing the structure of the prototypical diagrams for the DMI and the CBI that precisely the same rules apply to the CBI vector, $\mathbf{C}_{ij}$. The two vectors do not have to be collinear, notably if the only applicable symmetry is of type (f).

Connection to a phenomenological model. Another advantage of of our approach is apparent if we consider the appropriate phenomenological model for the magnetic interactions. To illustrate this point, we consider the most general spin model containing only bilinear and biquadratic pair interactions:

[TABLE]

The bilinear interactions are described by a rank-2 cartesian tensor $J_{ij}^{\alpha\beta}$ (9 parameters), which contains the isotropic pair interaction given by $J_{ij}\,\mathbf{S}_{i}\cdot\mathbf{S}_{j}$ (1 parameter), the DMI given by $\mathbf{D}_{ij}\cdot\left(\mathbf{S}_{i}\times\mathbf{S}_{j}\right)$ (3 parameters), and the remaining five parameters describe the symmetric bilinear pair anisotropy. The biquadratic interactions are described by a rank-4 cartesian tensor $B_{ij}^{\alpha\beta\gamma\delta}$ (81 parameters), and are not straightforward to classify. The number of independent elements of the biquadratic tensor is reduced to 25 by noting that $B_{ij}^{\alpha\beta\gamma\delta}=B_{ij}^{\gamma\beta\alpha\delta}=B_{ij}^{\alpha\delta\gamma\beta}=B_{ij}^{\gamma\delta\alpha\beta}$ and that excluding terms which are independent of the spin orientation requires $\sum_{\alpha}B_{ij}^{\alpha\beta\alpha\delta}=\sum_{\beta}B_{ij}^{\alpha\beta\gamma\beta}=0$. The same conclusion as to the number of independent parameters can be arrived at via the spin cluster expansion of the magnetic energy[40, 41] (see Supplementary Note 3). Making use of the prototypical diagrams, we already recovered the isotropic biquadratic interaction $B_{ij}\left(\mathbf{S}_{i}\cdot\mathbf{S}_{j}\right)^{2}$ (1 parameter), and we uncovered the CBI given by $\mathbf{C}_{ij}\cdot\left(\mathbf{S}_{i}\times\mathbf{S}_{j}\right)\left(\mathbf{S}_{i}\cdot\mathbf{S}_{j}\right)$ (3 parameters). Considering prototypical diagrams with more spin-orbit sites is a constructive approach to populate the rest of the $B_{ij}^{\alpha\beta\gamma\delta}$ tensor, from which the form of the magnetic interactions will also follow. If spin-orbit coupling is weak in some sense, we then also obtain a natural classification of the various magnetic interactions in powers of this small parameter. This would justify considering only the isotropic biquadratic and the CBI as the most important interactions among all biquadratic ones. Instead, we shall turn to realistic calculations of the magnitude and properties of the magnetic interactions.

0.2 Magnetic interactions of dimers on various surfaces

In order to quantify the properties and significance of the CBI in relation to the other magnetic interactions, we present a systematic study of a series of prototypical systems: magnetic dimers on several surfaces for which the spin-orbit effects are strong. To do so, we construct a complete magnetic model containing all relevant interactions up to four-spin couplings by defining a mapping from a set of self-consistent constrained DFT calculations, as explained in the Methods and Supplementary Note 4. The parametrizations of the complete magnetic model for all considered systems are given in Supplementary Note 5.

Simplified magnetic model. We focus on the following interactions: the isotropic bilinear interaction $J\,\hat{\mathbf{e}}_{1}\cdot\hat{\mathbf{e}}_{2}$, the DMI $\mathbf{D}\cdot(\hat{\mathbf{e}}_{1}\times\hat{\mathbf{e}}_{2})$, the isotropic biquadratic interaction $B\,(\hat{\mathbf{e}}_{1}\cdot\hat{\mathbf{e}}_{2})^{2}$, and the CBI $\mathbf{C}\cdot(\hat{\mathbf{e}}_{1}\times\hat{\mathbf{e}}_{2})\,(\hat{\mathbf{e}}_{1}\cdot\hat{\mathbf{e}}_{2})$. These are defined in terms of the unit vectors $\hat{\mathbf{e}}_{i}$ representing the orientation of the spin magnetic moment of the $i$-th atom in the dimer. The axis of the dimer is chosen as the $x$-axis, while the normal to the surface is chosen as the $z$-axis, as shown in Fig. 2a. Symmetry then restricts the DMI and CBI vectors to lie in the $yz$-plane for (111) and (0001) surfaces, illustrated in in Fig. 2b, or to lie along the $y$-axis for the (001) surface. In all cases the $y$-component is the dominant one. To understand the interplay between the different interactions, we consider the simplified model obtained from Eq. (3) by confining the magnetic moments to the $xz$-plane and keeping only the mentioned interactions:

[TABLE]

Here $\alpha=\theta_{2}-\theta_{1}$ is the opening angle between the two magnetic moments. The angle that minimizes the energy can be written as $\alpha_{\mathrm{min}}=\alpha_{J}+\Delta\alpha$, where $\alpha_{J}=$0\text{,}\mathrm{\SIUnitSymbolDegree}$$ if $J<0$ (ferromagnetic) or $180\text{\,}\mathrm{\SIUnitSymbolDegree}$ if $J>0$ (antiferromagnetic), and $\Delta\alpha$ is the canting induced by the remaining magnetic interactions. The energy can then be expanded as

[TABLE]

where the last line gives an approximation to the canting angle.

Magnetic dimers on Pt(111). We first compare the magnetic properties of five different homoatomic dimers on the Pt(111) surface, with the corresponding data collected in Table 1. All dimers except Ni possess large spin magnetic moments, which depend very weakly on the various imposed magnetic structures. Comparing the CBI to the DMI, we see that the magnitude of $C_{y}$ is around $20$–$30\%$ of the one of $D_{y}$, even reaching $60\%$ for Ni. For most dimers, $B$ is similar in magnitude to the CBI, and is even stronger than the DMI for Cr and Ni. According to $J$, which is the dominant interaction, Cr and Mn are antiferromagnetic, while Fe, Co and Ni are ferromagnetic. Considering only $J$ and $D_{y}$ leads to a canting of the magnetic structure given by $\Delta\alpha^{\mathrm{2s}}$ in Table 1, while considering also $B$ and $C_{y}$ we obtain $\Delta\alpha^{\mathrm{4s}}$. The difference between these values is the largest for Cr and Fe, so these are the dimers for which the biquadratic interactions are most important. Lastly, we also include the values of the effective bilinear interactions defined by the coefficients of $\Delta\alpha$ and $(\Delta\alpha)^{2}/2$ in Eq. (5). These correspond to $J^{\mathrm{eff}}=J-2\,\mathrm{sgn}(J)B$ and $D_{y}^{\mathrm{eff}}=D_{y}-\mathrm{sgn}(J)C_{y}$. The vector chirality of the magnetic ground state is set by the combination of the DMI and CBI vectors. These can be parallel, antiparallel, or substantially noncollinear (shown in Fig. 2c), in particular for the Co dimer. This shows that the CBI has not only the potential to impose the opposite vector spin chirality to the one favoured by the DMI ($\Delta\alpha$ changing sign in Eq. (5)), but also to tilt in away from the direction defined by the DMI vector.

Cr and Fe dimers on other surfaces. The Cr and Fe dimers on Pt(111) were found to have the most important contributions from the CBI. To ascertain whether this is particular to the Pt(111) surface, we placed these dimers on other surfaces with strong spin-orbit coupling, namely Pt(001), Ir(111) and Re(0001). We see from Table 1 that the CBI is generally a sizeable fraction of the DMI. On the Pt(001) and Ir(111) surfaces, the two dimers display a very large DMI, even in relation to its isotropic bilinear interaction $J$, leading to a strong canting of the magnetic structure. This canting is substantially modified when the biquadratic interactions are accounted for. The same behavior is found for the Cr dimer on Re(0001), while for the Fe dimer on this surface the interactions are found to be surprisingly weak, but still support a strongly noncollinear magnetic structure.

Electronic origin of the magnetic interactions. The origin of the different magnetic interactions can be further understood by comparing their dependence on the filling of the electronic states with the corresponding density of states of each dimer. This is shown in Fig. 3 for the Cr and Fe dimers on Pt(111). $D_{y}$ is largest in the energy range of the Pt $d$-states, which shows that the DMI is strongly enhanced by hybridization of the magnetic $d$-states of the dimer with the $d$-states of the Pt surface. $B$ and $C_{y}$ (and also $J$) have their largest values in the small energy range of the minority $d$-states of the Fe dimer, which suggests that the biquadratic interactions require less hybridization with Pt and so more localized magnetic $d$-states in the dimer. This can be understood from the microscopic theory. According to the prototypical diagrams in Fig. 1a,c, the isotropic interactions are direct interactions between the dimer atoms, resulting in a strong dependence on the Fe $d$-states. The prototypical diagram of Fig. 1b shows that the DMI is an interaction mediated by a spin-orbit site, which are supplied by the Pt surface atoms, and so this interaction is strongly dependent on the Pt $d$-states. The CBI involves both a direct exchange between the magnetic sites and an excursion through a spin-orbit site, Fig. 1d, so it can be amplified in those two ways, leading to a more complicated dependence.

0.3 Implications of the chiral biquadratic interaction

The CBI has different important implications for a broad class of noncollinear magnetic nanostructures. For a magnetic dimer, we already found that the CBI influences the opening angle and the vector spin chirality of the magnetic structure. To gain further understanding, we return to the previous example of a CBI vector in the $y$-direction, for which the interaction energy has the form $E_{\mathrm{CBI}}(\alpha)=C_{y}\sin\alpha\cos\alpha$ (see Eq. (4)). The $\cos\alpha$ term comes from the $\hat{\mathbf{e}}_{1}\cdot\hat{\mathbf{e}}_{2}$ part of the interaction, while the $\sin\alpha$ terms comes from $(\hat{\mathbf{e}}_{1}\times\hat{\mathbf{e}}_{2})_{y}$, with $\alpha$ the opening angle. As the dot product is isotropic, the opening is favored in the plane perpendicular to the CBI vector. Fixing $C_{y}>0$ for definiteness, there are two energy minima for $\alpha_{\mathrm{min}}\in\{$-45\text{,}\mathrm{\SIUnitSymbolDegree}$,$135\text{,}\mathrm{\SIUnitSymbolDegree}$\}$, and two maxima for $\alpha_{\mathrm{max}}\in\{$45\text{,}\mathrm{\SIUnitSymbolDegree}$,$-135\text{,}\mathrm{\SIUnitSymbolDegree}$\}$. Strikingly, the two values of $\alpha_{\mathrm{min}}$ have opposite signs, which means that the sign of their vector spin chirality (projected on the $y$-axis) is also opposite. Thus, and in contrast to the DMI, the CBI favors both possible rotational senses at once (although with different opening angles). Starting from a ferromagnetic or antiferromagnetic structure (set by $J$), the DMI will induce a canting of the same rotational sense for both cases, while the CBI will favor cantings for each structure which have opposite rotational senses.

From a dimer to an infinite chain. Next we relate the magnetic ground state of a dimer to that of an infinite chain, assuming that the interactions present are $J$, $D_{y}$ and $C_{y}$ (we take $B=0$ for simplicity), being nearest-neighbor interactions for the chain. The energy as a function of the opening angle for the dimer is given by Eq. (4), and the same form applies for the energy of a spiral magnetic structure on an infinite chain, with $\mathbf{\hat{e}}_{i}=\sin(Qx_{i})\,\mathbf{\hat{e}}_{x}+\cos(Qx_{i})\,\mathbf{\hat{e}}_{z}$, by setting $Q=\alpha/a$ ($a$ is the nearest-neighbor distance). The canting angle/spiral wavevector that minimizes the energy is shown for different values of $D_{y}/C_{y}$ and $J/C_{y}$ in Fig. 4. The sharp transition in the sign of $\alpha$ found along $J=D_{y}$ is driven by $C_{y}$, and signals a change in the sign of the vector spin chirality. At the edges of the diagram we recover well-known magnetic structures ($\alpha=$0\text{,}\mathrm{\SIUnitSymbolDegree}$,$180\text{,}\mathrm{\SIUnitSymbolDegree} are the ferromagnetic and antiferromagnetic structures for $J$-only, respectively, and $\alpha=\pm$90\text{\,}\mathrm{\SIUnitSymbolDegree} pertain to $D_{y}$-only).

Impact of the CBI on complex 2D magnetic structures. As a final example, we consider two-dimensional magnetic systems. Higher-order isotropic interactions can help stabilizing complex magnetic structures called multiple-$Q$-states[21, 19, 36, 20]. We thus address the potential role of that the CBI might play for such complex magnetic structures, choosing the Fe monolayer on Ir(111) by way of example[36]. The ground state is a nanoskyrmion lattice, which is a type of $2Q$-state made of two symmetry-related wavevectors $\mathbf{Q}_{1}$ and $\mathbf{Q}_{2}$. Other combinations lead to further noncollinear states which were calculated to have a similar energy: single-$Q$ spin spirals, the $\mathbf{Q}_{m}$-star and $\mathbf{Q}_{m}$-vortex states, and the nanovortex lattice. These magnetic structures are visualized in Fig. 5a-f, with the triangular Néel state included for comparison. To show how different magnetic interactions affect these states, we write each energy contribution as $E_{I}=N\,I\,\Gamma_{I}$, where $N$ is the number of atoms in the magnetic unit cell, $I=\{D,B,C\}$ is the interaction strength, and $\Gamma_{I}$ is a factor determined by the lattice and the magnetic structure. Fig. 5g plots the $\Gamma_{I}$ factors computed assuming nearest-neighbor pair interactions following $C_{\mathrm{3v}}$ symmetry. We do not discuss the contribution of $J$, as this interaction is complex and long-ranged for this system[36]. Interestingly, the contribution of the isotropic biquadratic interaction $B$ is similar for every structure, none being particularly favoured. On the other hand, both the DMI as well as the CBI do favour certain magnetic structures. The DMI favours almost equally the single-$Q$ spiral and the nanoskyrmion lattice, while the CBI clearly favours the nanoskyrmion lattice ($2Q$ state) over the spin spirals. This shows that the CBI can play a deciding role in stabilizing noncoplanar multiple-$Q$-states, which have a non-vanishing scalar spin chirality.

Discussion

We presented a comprehensive analysis of a new chiral higher-order magnetic pair interaction, the chiral biquadratic interaction (CBI). Using a microscopic model and a systematic expansion of the electronic grand-potential, we identified the prototypical diagrams behind all kinds of magnetic interactions. This led us to uncover a new chiral interaction, the CBI, which is linear in the spin-orbit coupling and is the biquadratic equivalent of the DMI, following the same symmetry rules. In its most general form, this interaction couples four distinct magnetic sites, and consists of terms of the form $\mathbf{C}_{ijkl}\cdot(\mathbf{S}_{i}\times\mathbf{S}_{j})\,(\mathbf{S}_{k}\cdot\mathbf{S}_{l})$. We note that a recent study has found signatures of higher-order interactions in magnetic chains on the Re(0001) surface[42]. We proved the existence and quantified the importance of this new interaction in relation to other well-known ones by performing systematic first-principles calculations for prototypical systems, magnetic dimers on surfaces with strong spin-orbit coupling. For these dimers, the CBI is typically around $20-30\%$ of the DMI, reaching $60\text{\,}\mathrm{\char 37\relax}$ for a Ni dimer deposited on Pt(111). Furthermore, we contrasted the properties of the CBI with those of the DMI, revealing that the CBI can simultaneously favor structures with opposite vector spin chiralities (canted ferromagnetic vs. canted antiferromagnetic structure), and that the CBI vector is not necessarily parallel to the DMI one. Insight into the electronic origin of the CBI and the other magnetic interactions was obtained via their dependence on the filling of the electronic states, supporting the proposed mechanism derived from the microscopic model. We related the magnetic phase diagram of a dimer to that of an infinite chain, and surveyed the implications of the CBI for magnetic monolayers. Considering several magnetic structures relevant for an Fe monolayer on Ir(111), we showed that the CBI can play a deciding role in favoring a noncoplanar nanoskyrmion lattice over other possible structures. We thus conclude that the CBI can influence simultaneously the vector and the scalar spin chirality, opening new avenues for stabilizing and engineering complex magnetic textures with both fundamental and technological interest. We also expect that there should be a non-trivial impact on the related dynamics, transport and topological properties, which should motivate future investigations.

{methods}

0.4 Density functional theory

We performed systematic density functional theory calculations with the full-potential Korringa-Kohn-Rostoker (KKR) Green function method[43]. Exchange and correlation effects are treated in the local spin density approximation (LSDA) as parametrized by Vosko, Wilk and Nusair[44], and SOC is added to the scalar-relativistic approximation[45]. The pristine surfaces are modelled by a slab of 22 layers (except for Pt(111) for which 40 layers were used) with the experimental lattice constants with open boundary conditions in the stacking direction, and surrounded by two vacuum regions. No relaxation of the surface layer of the pristine surface is considered, as it was shown to be negligible[46]. We use $150\times 150$ $k$-points in the two-dimensional Brillouin zone, and the angular momentum expansions for the scattering problem are carried out up to $\ell_{\mathrm{max}}=3$. In the next step, we utilize an embedding method to place each dimer on the fcc-like threefold hollow position of Pt(111) and Ir(111), on the fourfold hollow sites on the Pt(001) surface, and on the hcp-like threefold hollow site of Re(0001). The embedding region consists of a spherical cluster around each magnetic adatom including nearest-neighbor Pt atoms. The appropriate structural relaxations of each dimer towards the surface were obtained with the plane-wave code Quantum Espresso[47], using ultrasoft scalar relativistic pseudopotentials[48] with the PBEsol functional[49], and considering a $4\times 4$ supercell with $5$ substrate layers (more details and structural information are given in Supplementary Note 4). All calculations were performed on the supercomputer JURECA at the Forschungszentrum Jülich[50].

0.5 Mapping DFT to an atomistic spin model — The torque method

To map the first-principles calculations to an atomistic spin model we adopted the method of constraining fields[51]. The net magnetic moment for each atom in the dimer is defined by averaging the spin magnetization density over its atomic cell, $\mathbf{m}_{i}=\int\!\mathrm{d}\mathbf{r}\;\mathbf{m}_{i}(\mathbf{r})=|\mathbf{m}_{i}|\,\hat{\mathbf{e}}_{i}$. The orientations of the magnetic moments $\hat{\mathbf{e}}_{i}$ play the role of emergent degrees of freedom that can be used to define an atomistic spin model,

[TABLE]

To determine the energy of a target magnetic configuration, the total energy functional is augmented by a Zeeman term enforcing the constraint,

[TABLE]

The constraining magnetic field is transverse to the orientation of the local magnetic moment, $\mathbf{b}_{i}\cdot\hat{\mathbf{e}}_{i}=0$, and it opposes the magnetic force that acts on it if the magnetic structure is not a stationary point of the total energy functional,

[TABLE]

The induced moments in the surface atoms are allowed to relax without any constraint. The model parameters are then determined by linear least-squares fitting the constraining fields obtained for a set of self-consistent magnetic configurations to the form of the magnetic force supplied by the atomistic spin model. The magnetic configurations for the dimers have been chosen using a Lebedev grid[52] containing $14$ directions for each atom, which is well-suited to describe spherical harmonics up to $\ell=2$, resulting in a total of $14^{2}=196$ configurations, which using symmetry arguments (time-reversal invariance of the magnetic energy plus the spatial symmetries that apply on different surfaces) the number of configurations can be further reduced to $56$ for (111) and (0001) surfaces and to $36$ for (001) surfaces.

{addendum}

This work was supported by the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (ERC-consolidator grant 681405 – DYNASORE). The authors gratefully acknowledge the computing time granted through JARA-HPC on the supercomputer JURECA at the Forschungszentrum Jülich[50].

S.B. and M.d.S.D. derived the magnetic exchange interactions from the microscopic model and defined the mapping procedure from the first-principles calculations. S.B. performed the DFT calculations for all considered systems, the parametrization of the extended magnetic hamiltonian and the model study of several extended magnetic structures. All authors analyzed the results and discussed their significance, and collaborated in writing the manuscript.

The authors declare that the main data supporting the findings of this study are available within the article and its Supplementary Information files. Further data is available from the corresponding author upon reasonable request.

The authors declare that they have no competing financial interests.

Correspondence and requests for materials should be addressed to S.B. (email: s.brinker@fz-juelich.de) or to S.L. (email: s.lounis@fz-juelich.de).