Tuning the size of skyrmion by strain at the Co/Pt3 interfaces.

Based on density functional theory calculations, we elucidated the tunability of the atomic structures and magnetic interactions of Co/Pt3 interface (one layer of hcp(0001) Co and three layers of fcc(111) Pt) and thus the skyrmion sizes using strain. The dispersion relations of the spin spiral in the opposite directions, E(q) and E(-q), were evaluated based on generalized Bloch equations. Effective exchange coupling (EC) and Dzyaloshinsky-Moriya interaction (DMI) parameters between different neighbors J i and d i at different lattice constants were derived by fitting the resulting spin spiral dispersion E(q) to EC model with DMI and E(q)-E(-q) formula, respectively. We observed an increase in DMI and a significant decrease in EC with an increase in strain. Hence, the size of Néel-type skyrmions determined by the ratio of EC/DMI can be controlled by applying strain, leading to an effective approach to tailor the formation of skyrmion lattices by inducing slight structural modifications on the magnetic thin films.

Introduction

Magnetic skyrmions, which are topologically protected particle-like nanoscale eddy magnetic structures (Bogdanov and Yablonskii, 1989a, 1989b; Bogdanov et al., 1989), have attracted much attention because of their potential application as skyrmions in high-density, low-power, nonvolatile computing and storage devices in future. So far, they have been observed in many asymmetric chiral magnetic materials or magnetic thin-film materials with broken spatial symmetry (Romming et al., 2013; Gong and Zhang, 2019; Fert et al., 2013). The mechanism for generating magnetic skyrmions is Dzyaloshinsky-Moriya interaction (DMI) (Dzyaloshinsky, 1958; Moriya, 1960; Li et al., 2018; Huang et al., 2018; Emori et al., 2013; Li et al., 2019; Mutter et al., 2019; Mühlbauer et al., 2009; Yu et al., 2010). The key to DMI lies in the spin-orbit coupling and the breaking of symmetry of space inversion, depending on the direction of DMI. The two-dimensional interface with non-centrosymmetry in the vertical direction can generate stable Néel-type skyrmions (Bode et al., 2007). More and more two-dimensional materials have been used to explore the generation of the skyrmions. For example, first-principles calculations and experiments suggest that Co on graphene (Co/Gra) may also have a strong DMI — up to 1.14 meV per atom (Yang et al., 2018a, 2018b) — and it provides a broad direction for the study of spin-orbit motion based on two-dimensional materials. For the well-studied theoretically and experimentally system Pd/Fe/Ir(111), the calculated skyrmion size is 6.5 nm (Zhu et al., 2020a, 2020b), which is very close to the experimental value, 5 nm (Romming et al., 2013). Co layers on Pt(111) surface (Co/Pt) is a kind of classical metal interface with strong DMI, located in the interfacial Co layer, originating from spin-orbit energy provided by the adjacent nonmagnetic Pt layer (Yang et al., 2015; Zimmermann et al., 2019; Deger, 2020; Gusev et al., 2020).

Among these Co/Pt researches, it is proved in both calculations (Deger, 2020) and experiments (Gusev et al., 2020) that strain could enhance the DMI at the interface. A strong variation (from 0.1 to 0.8 mJ/m2) of the DMI constant is demonstrated at 0.1 in in-plane uniaxial deformation of the films and it even changes the sign of the DMI in the direction perpendicular to the strain direction (Gusev et al., 2020). In the experimental results of other systems, anisotropic strain as small as 0.3 modulates the DMI and induces very large deformations of skyrmions in FeGe (Shibata et al., 2015); however, in the case of MnSi, it has also been found that the magnetic modulation period in the skyrmion crystal phase remains almost unchanged on the application of a uniaxial stress and the change is 0.5%, so the effect of strain on the DMI is less relevant in the case of MnSi (Nii et al., 2015); in addition, mechanical strain is an effective means to control the skyrmion in some works (Wu et al., 2015; Liu et al., 2017; Yu et al., 2018). Wu et al. explore physical compression and chemical expansion effects on the skyrmion phase in Cu2OSeO3 and they suggest that tuning of the skyrmion phase under different kinds of pressures might be related to the interatomic bond lengths and bond angles that modulate the competing Heisenberg exchange interactions and DMI and the pressure effect can be used to tune the skyrmion phase (Wu et al., 2015). Liu et al. demonstrate that it is possible to cut a skyrmion from a chiral stripe domain by applying an in-plane uniaxial strain (Liu et al., 2017). Yu et al. find that the hexagonal skyrmion lattice survives at lower temperatures, and no topological phase transition occurs with decreasing temperature, although the shape of the individual skyrmion becomes deformed, possibly owing to the effect of anisotropic strain (Yu et al., 2018). These researches tell us that the strain would provide an efficient way to vary DMI to tune the size of skyrmions. Generally, the size of Néel skyrmions can be effectively reduced by enhancing DMI. Furthermore, it would be more convenient to reduce exchange coupling (EC) to reduce the size of skyrmions (Zhu et al., 2020a, 2020b). However, studies on the overall effect of deformation on skyrmions are still lacking for Co/Pt interface (Deger, 2020; Gusev et al., 2020) and other systems (Shibata et al., 2015; Nii et al., 2015; Wu et al., 2015; Liu et al., 2017; Yu et al., 2018; Koretsune et al., 2015). Therefore, we would investigate how the deformation in the Co/Pt interface simultaneously changes both EC and DMI, and jointly affects the size of the skyrmion. This can lay a solid foundation for future deformation research of skyrmion.

In this work, we calculate the dispersion relation between energy and wave vector (q) of the spin spiral in opposite directions: E(q) and E(-q) of Co/Pt3 interface (one layer of hcp(0001) Co and three layers of fcc(111) Pt) under generalized Bloch condition by first-principles calculations, EC, and DMI at different lattice constants were obtained by fitting. As the lattice constant increases from 2.72 Å to 2.84 Å, d1 increases, d2 decreases, and J decreases, it can be seen that the radius of skyrmion is reduced obviously with strain. The variation of EC and DMI were observed with changes in the lattice constants of Co/Pt3 films, resulting in changes of the size of skyrmions. This implies that the strain-induced modification of both DMI and EC dominates the formation of skyrmions.

In most VASP calculations in the past, the EC J1 between the nearest neighbor atoms is only considered and J between other farther neighbors was not obtained (Yang et al., 2015, 2018a, 2018b; Zhu et al., 2014; Pan et al., 2015). In fact, M. Marsman (Marsman and Hafner, 2002) et al. has calculated the dispersion relationship of γ-Fe spin spiral in VASP for the first time, and the J of each neighbor was obtained and used to explain the magnetic structure of γ-Fe. This paper attempts to calculate the dispersion relation of the spin spiral by nonlinear calculation (Pan et al., 2015), taking spin orbit coupling into account and using the generalized Bloch condition (Bode et al., 2007; Zimmermann et al., 2019; von Malottki et al., 2017; Sandratskii, 1991; Knöpfle et al., 2000; Zhu et al., 2019, 2020a, 2020b). By means of this, the EC parameter J and the DMI parameter d of each neighbor could be obtained.

The Hamiltonian of a Heisenberg model system with DMI can be expressed as:where 0 is the normalized spin vector. For thin films, we apply d0 = d0(z×u0), with u0 being the unit vector between sites 0 and i, z being the unit vector normal to the plane as shown in Figure 1 and d0 is the size of DMI between each neighbor. A schematic diagram of a single Néel skyrmion for 2D hexagonal structures is shown in Figure 1A. The radius (Bogdanov and Hubert, 1994; Han et al., 2010) of the skyrmion defined according to the principle of minimum average energy is expressed as r = Ns∗ × a, where a is the lattice constant and Ns∗ is the number of magnetic moments along the radius to the edge of skyrmion. For a Néel skyrmion, the magnetic moment is gradually rotated from the downward direction of the center (sz = −1) to the upward direction of the edge (sz = +1) in the hexagonal lattice. The actual diameter of a skyrmion is described as D, where D = 2 × a × Ns∗/10 (the unit of lattice constant a is Å and the unit of D is converted to nanometers). The change in EC energy caused by magnetic moment rotation is the relative value compared to the ferromagnetic (FM) state, therefore, for convenience of discussions, the energy of the FM ground state is set to zero (Zhu et al., 2020a, 2020b). If the magnetic moment is set to be changed in a plane perpendicular to the two-dimensional interface to form a spin spiral, and the value of d0 is not zero, then DMI breaks the degeneracy of the spin spiral, resulting in the spatial asymmetry of the spin spiral dispersion relationship. VASP is used to calculate the dispersion relation E(q) and E(-q) of the spin spiral in the opposite direction, where q is the direction vector of the spin spiral. The ΔDMI(q) = E(q)-E(-q) is twice the DMI energy of the calculated system along as q. Generally speaking, under the same wave vector q, the larger the absolute value of DMI is, the larger the d0, the smaller the skyrmions, and then it is more convenient for practical application.

Figure 1

Magnetic structure of Néel skyrmions and atomic structure of Co/Pt3

(A) A schematic representation for the magnetic structure of Néel skyrmions.

(B) Atomic structure of Co/Pt3: top view of the extended 3 × 3 supercell used to describe the EC and DMI for the interfacial Co layer between different neighbors in Co/Pt3, as shown the yellow numbers on the atom; a1 and a2 are basis vectors; b1 and b2 are reciprocal lattice vectors; Γ, M and K are high symmetry k-points in first Brillouin zone (Momma and Izumi, 2011).

(C) Side view of Co/Pt3 atomic structure, different atoms are labeled, blue and gray balls are for Co and Pt atoms respectively; the distance at z-axis between the Co atom and Pt3 atom is labeled as hz (Momma and Izumi, 2011).

Results and discussion

The interface structure for calculation is composed of three layers of fcc(111) Pt and one layer of hcp(0001) Co atoms. The formula of E(q) is closely related to the neighbor number of the chosen Co atom. In the top view, the central Co atom is labeled as 0 (named Co_0) and the other numbers represent the neighbors of this atom. The number of atoms in the first to third nearest neighbor is 6 and the number of atoms in the fourth nearest neighbor is 12.

We set the position of Co_0 as the zero point (whose magnetic moment is S(0)) and the position of the Co atom at position j is represented by Rj = m1 + n2, where 1 and 2 are the basis vectors. The orientation of the spin spiral is described as q = q11 + q22, where 1 and 2 are the reciprocal vectors. In VASP calculations, all magnetic moments are set in the xz plane, with S(0) along the z axis. Under the generalized Bloch condition, in the spin wave with wave vector q, the magnetic moment S(Rj) of the j-th adjacent Co atom is:

Then, according to Equation (1), considering the first to fourth nearest neighbor merely, the sum of EC between S(0) and the magnetic moment S(Rj) of the j-th nearest neighbor is as follows:

The sum of EC is as follows:

Among them, J1 to J4 are EC parameters, including S(0)2.

DMI is considered to the second neighbor and the sum of DMI between S(0) and S(Rj) is:

The sum of DMI is:

Among them, d1 and d2 are DMI parameters, including S(0)2. The total energy of spin spiral E(q) includes EC and DMI:

Obviously, EEC(q) is an even function of q, and EDMI(q) is an odd function of q, so:

VASP is used to calculate the dispersion relation E(q) and E(-q) described by Equation (4) for selected lattice constants, as shown in Figures 2A, 2C and 2E and 2G. The wave vector q of the spin spiral selects the boundary of the simplified Brillouin zone shown in Figure 1A. For EC, KΓ and -K-Γ, ΓM and -Γ-M are equivalent shown in Equation (2), respectively. However, DMI breaks the symmetry, so the energy of E(q) and E(-q) is not exactly the same, as shown in Figures 2A, 2C, 2E, and 2G. For the selected lattice constant, the energy differences of ΔDMI(q) can be regularly displayed in discrete points in Figures 2B, 2D, 2F, and 2H.

Figure 2

Energy dispersion relation of spin spiral

(A, C, E, and G) The scattered symbols are calculated energy dispersion E(q) and E(-q) of spin spirals at lattice constants of 2.72 Å, 2.76 Å, 2.80 Å and 2.84 Å respectively in Co/Pt3 as a function of the spiral wave vector q.

(B, D, F, and H) The scattered symbols are ΔDMI(q), obtained from the energy differences between E(q) and E(-q). In (A) to (H), lines are the fitted ones. K, Γ, and M are special K-points in the first Brillouin zone shown in Figure 1.

The dispersion relation E(q) and E(-q) of Co/Pt3 spin-spiral calculated by VASP should satisfy Equation (4) and in fact, in Figures 2A, 2C, 2E, and 2G above, the solid line fitted by the method of least square agrees well with the discrete point of the calculated value. As shown in Figures 2B, 2D, 2F, and 2H, the difference ΔDMI(q) between E (q) and E (-q) is twice of the sum of the DMI between S(0) and other magnetic moments, which should satisfy Equation (5) as well. The fitted solid line is basically consistent with the trend of calculated discrete points. It can be seen from Figure 2 that E(q) and E(-q) gradually decrease with increasing lattice constant. In contrast, ΔDMI(q) increases significantly with the increasing lattice constant. Both EC and DMI vary with the lattice constants. This is the basis of the Co/Pt3 interface to tune skyrmions. The variations of specific magnetic parameters with lattice constants are described in Figure 3.

Figure 3

Trend of changes about the fitted DMI d and EC J

(A) and (B) are the trend of changes about the fitted DMI d and EC J. J is considered to be the fourth neighbor and d is considered to be the second neighbor. Zero is shown as a horizontal dotted line.

The calculated lattice constant of Pt is 2.72 Å, and the experimental value is 2.81 Å. We carefully select the lattice constants close to our calculated and experimental values to discuss the effect of deformation. The EC parameters J1 to J4 and the DMI parameters d1 and d2 at different lattice constants were fitted by the method of least square shown in Figure 3. The sum of J1 and J2 in our work is similar to the sum of J1 and J2 in others' work (Zimmermann et al., 2019). Obviously, Co/Pt3 as classical skyrmion interfaces (Zhu et al., 2020a, 2020b; Belmeguenai et al., 2015; Boulle et al., 2016; Juge et al., 2018) has very large d1 parameters. As the lattice constant increases, d1 increases, d2 decreases and J decreases, as shown in Figure 3. For d1, if we increase the lattice constant of Co/Pt3 from 2.72 Å to 2.84 Å, d1 and d2 increases and decreases by 2.94 and 1.49 meV, respectively, as shown in Figure 3. The decrease in d2 is less than the increase in d1. Hence, it could be seen that the magnitude of DMI for the first neighbor increases with increasing lattice constant, while J decreases all the time. When the lattice constant is 2.84 Å, Co/Pt3 has the largest d1, 6.89 meV, which is about 3 meV larger than Co/Pt3 at 2.72 Å.

For a single Co atom, the d orbital has five degenerate orbits, and a very localized d-states peak is formed in the DOS diagram. When atoms interact to form crystals, the peaks tend to flatten as the distance between them decreases. The DOS is uniformly distributed at a lattice constant of 2.72 Å, but relatively locally distributed at a lattice constant of 2.84 Å as shown in Figure 4. When the distance between Co atoms in Co/Pt3 is small, the magnetic interaction between Co/Pt3 is stronger, and the EC is also stronger, resulting in a larger J value for EC parameter. Conversely, when Co/Pt3 expands, the distance between atoms increases and the interaction would decrease, resulting in a smaller J value.

Figure 4

Partial DOS of the Co atom

Partial DOS of the Co atom at lattice constant of 2.72 Å and 2.84 Å. Zero energy indicates the position of the Fermi level shown as a vertical black dotted line. The horizontal black dotted line is used to distinguish the spin up and down states.

In the COHP curve of Figure 5, it is notable that the Fermi level passes through the antibonding area of all the bonds (Song and Zhao, 2014), and there are some antibonding states present below the Fermi level. This is because the atomic radius of Pt is much larger than that of Co. However, the bonding area below the Fermi level is still far larger than the antibonding area, therefore, the Co can grow on the Pt layer. In fact, there has been a Co/Pt interface (Ryu et al., 2013) in experiment. Moreover, with the increase of the lattice constant, the antibonding area below the Fermi level increases slightly but it is still substantially smaller than the bonding area, and this trend can also be seen from the DOS in Figure 4. The results show that by increasing the lattice constant, the distance between Co and Pt layer can be reduced.

Figure 5

-COHP of Co-Pt bond in the Co/Pt3 system

-COHP of Co-Pt bond in the Co/Pt3 system at lattice constants of 2.72, 2.76, 2.80, 2.84Å respectively considering the spin-polarization.

These fitted J and d can still fit well the simple extended NdJ relationship to form skyrmions (Zhu et al., 2020a, 2020b).

This formula indicates that there is a simple relationship between the skyrmion size and fundamental magnetic parameters. From this equation, we can see that the size of skyrmion can be reduced by increasing d in the denominator or decreasing the effective J in the numerator, which may be achieved via mutual cancellation among J (Zhu et al., 2020a, 2020b).

The contribution of the decrease in J to the reduction of the size of the skyrmion is substantially close to the increase in d as shown in Figures 6A and 6B. We increase the lattice constant of Co/Pt3, which can effectively increase J and decrease d at the same time. According to the Equation (6), as shown in Figure 6C, Ns∗ decreases obviously with strain, and as the lattice constant of Co/Pt3 increases from 2.72 Å to 2.84 Å, the actual size of skyrmions at different lattice constants varies from 22 nm to 17 nm as shown in Figure 6D. This indicates that the skyrmion size decreases as the lattice constant increases. Deformation can indeed be used to tune the size of the skyrmion.

Figure 6

Effect of DMI and EC on the formation of skyrmions

(A and B) are effect of DMI and EC on the formation of skyrmions with strain. (C) and (D) are calculated Ns∗ and the actual diameter of Skyrmion described as D at different lattice constants.

What is causing the skyrmion to shrink? When the lattice constant is 2.84 Å, Co is close to the Pt layer, and the height between the Co and Pt layer increases as the lattice constant decreases, as shown in Figure 7. This indicates that the distance between magnet and the substrate is one of the main factors determining EC and DMI. Why can we judge that the magnitude of DMI is determined by distance? It is found that DMI is predominantly located in the interfacial Co layer, originating from spin-orbit coupling provided by the adjacent nonmagnetic layer (Yang et al., 2015, 2018a, 2018b). With the increase of the lattice constant, the spin-orbit coupling enhances as shown in Figure 7, and the magneto crystalline anisotropy energy also increases as shown in Figure S3. The spin-orbit coupling provided by the heavy metal Pt to the adjacent Co will be strong if the Co atom is close to the Pt layer, and weaker otherwise. When we increase the lattice constant, the gap between Pt atoms becomes larger, which makes it for the Co atom layer to fall into the gap, closer to the Pt layer. More interestingly, as the lattice constant increases, the distance between Co atoms also becomes larger, which also reduces the EC between Co atoms. Therefore, increasing the lattice constant of the Co/Pt interface has a dual effect on the reducing the size of skyrmions, and thus is an extremely effective way to tune the size of skyrmions.

Figure 7

Distance and energy difference

Distance at z-axis between Co atom and adjacent Pt atom (black square) and energy difference between without and with spin-orbit coupling (red circles) at different lattice constants.

Conclusions

In conclusion, general magnetic exchange models including both EC and DMI can be well fitted using the calculated spin spiral E(q) spectra with EC parameters J1 to J4 and DMI parameters d1 and d2 for each neighbor at different lattice constants of Pt atom are obtained by fitting with the least square method. As the lattice constant of Co/Pt3 films increases, it is clear that DMI increases and EC decreases. Based on the extended NdJ relationship, the size of skyrmions (Behera et al., 2019) can be tuned by either enhancing d or more conveniently reducing J (Yang et al., 2018a, 2018b; Chen et al., 2013, 2015; Johansen et al., 2019). Furthermore, enlarging the lattice constant of the Co/Pt3 cell can satisfy both of these conditions at the same time, and sharply reduce the diameter of the skyrmion. Our results provide simple guidelines for the design of skyrmion materials, which are still challenging in this field (Schott et al., 2017; Pizzini et al., 2014; Kézsmárki et al., 2015; Ruff et al., 2015; Kurumaji et al., 2017).

Limitations of the study

There are two major limitations in this study that could be addressed in future research. First, we only consider one layer of Co atoms, and there are relatively thick layers of Co atoms in the experiments. But DMI mainly exists at the interface, and HBI will increase with the increase of Co layer, which is not beneficial to skyrmion generation according to our extended NdJ relationship. The experiment can also try the case with a few Co layers, or the case of one layer. Second, Co atoms are relatively small. In the case of tensile deformation, Co atoms may penetrate into the Pt layer.

STAR★Methods

Key resources table

Resource availability

Lead contact

Further information and requests should be directed to and will be fulfilled by the lead contact, Zhu Yan (yzhu@nuaa.edu.cn).

Materials availability

This study did not generate any unique reagents.

Experimental model and subject details

Our study does not use experimental models typical in the life sciences.

Method details

Computational details

We adopt Vienna ab initio simulation package (VASP), a software of first-principles pseudopotential plane wave method (Kresse and Joubert, 1999) based on Density Functional Theory (DFT), for calculation. The Perdew-Burke-Ernzerhof (PBE) potential (Perdew et al., 1996) of VASP is used for the calculation. VASP can employ plane-wave basis functions to solve the Kohn-Sham equation through self-consistent iterative method, and calculate the force and tensor by wave function. Co/Pt3 interface is simulated by a slab model with periodic boundary conditions and along with a vacuum region of 15 Å between adjacent slabs to avoid a spurious dipole moment from image supercells. Atoms in the calculated system are relaxed to the ground state until the forces are less than 10−2 eV/ Å. Energy cutoffs and k-points are tested in the supplementary materials, and 225 eV of energy cutoff and 9 × 9×1 of k-points is sufficient. However, they are chosen as 375 eV and 41 × 41×1 k-point grids respectively for accurate calculations, which are much larger than that typically recommend. The convergence criterion for self-consistent total-energy calculations is 1.0 × 10−6 eV. Methfessel-Paxton smearing with a half-width of 0.05 eV is used to accelerate the convergence for relax calculations and static calculations. We perform fully noncollinear magnetic calculations within the projector-augmented wave (PAW) formalism, as implemented in the VASP code by Hobbs et al. (2000), and spin-orbit coupling (SOC) is also included in the present calculations. In order to express the bonding character and perform a comparative analysis with the density of states (DOS), the crystal orbital Hamilton population (Dronskowski and Bloechl, 1993) (COHP) was calculated (Soler et al., 2002) with numerical atomic orbital basis sets (Junquera et al., 2001) and Troullier-Martins norm-conserving pseudopotentials (Troullier and Martins, 1991).

Quantification and statistical analysis

Our study does not include statistical analysis or quantification.

Additional resources

Our study has not generated or contributed to a new website/forum or if it is not part of a clinical trial.

REAGENT or RESOURCE	SOURCE	IDENTIFIER
Software and algorithms
VESTA	K. Momma and F. Izumi, 2011	https://doi.org/10.1107/S0021889811038970
Online Least Squares Fitting	Yunsuan Network	http://www.yunsuan.info/index.html