Comparison of first-principles methods to extract magnetic parameters in ultra-thin films: Co/Pt(111)
Bernd Zimmermann, Gustav Bihlmayer, Marie B\"ottcher and, Mohammed Bouhassoune, Samir Lounis, Jairo Sinova, Stefan Heinze and, Stefan Bl\"ugel, Bertrand Dup\'e

TL;DR
This study compares three first-principles computational methods to accurately determine magnetic interactions in a Co monolayer on Pt(111), revealing differences in parameters depending on the magnetic structure and providing insights into magnetic properties at the nanoscale.
Contribution
It systematically compares different first-principles approaches for magnetic parameter extraction in ultra-thin films, highlighting their consistency and differences in various magnetic regimes.
Findings
Methods (i) and (ii) yield consistent micromagnetic parameters.
Calculated spin stiffness and Curie temperature are significantly higher than bulk Co.
Nearest-neighbor exchange interaction in the monolayer is increased by 50%.
Abstract
We compare three distinct computational approaches based on first-principles calculations within density functional theory to explore the magnetic exchange and the Dzyaloshinskii-Moriya interactions (DMI) of a Co monolayer on Pt(111), namely (i) the method of infinitesimal rotations of magnetic moments based on the Korringa-Kohn-Rostoker (KKR) Green function method, (ii) the generalized Bloch theorem applied to spiraling magnetic structures and (iii) supercell calculations with non-collinear magnetic moments, the latter two being based on the full-potential linearized augmented plane wave (FLAPW) method. In particular, we show that the magnetic interaction parameters entering micromagnetic models describing the long-wavelength deviations from the ferromagnetic state might be different from those calculated for fast rotating magnetic structures, as they are obtained by using (necessarily…
Click any figure to enlarge with its caption.
Figure 1
Figure 2
Figure 1
Figure 4
Figure 5
Figure 6| Co/Pt(111) | ||
|---|---|---|
| fcc | hcp | |
| 2.76 | 2.76 | |
| 2.02 | 2.03 | |
| 2.37 | 2.38 | |
| Total energy | 0 | 134 |
| self-consistent | perturbative | |
|---|---|---|
| KKR | FM state | non-collinearity |
| FLAPW–gBT | any spin-spiral | spin-orbit coupling |
| FLAPW–supercell | large spin-spiral | — |
| stacking | |||||||||
|---|---|---|---|---|---|---|---|---|---|
| position | (meV) | (meV) | (meV) | (meV) | (pJ/m) | (mJ/m2) | (meV) | (meV) | |
| KKR | Co(fcc) | 19.9 | 1.7 | 0.4 | 1.1 | 40.9 (38.8) | 17.6 (14.76) | 29.5 (27.9) | 1.75 (1.47) |
| Co(hcp) | 20.8 | 1.5 | 0.2 | 1.0 | |||||
| FLAPW-gBT (coned) | Co(fcc) | 26.0 | 1.2 | 0.7 | 1.2 | 44.0 | 14.4 | 31.7 | 1.43 |
| Co(hcp) | 27.4 | 0.6 | 0.4 | 1.0 | |||||
| FLAPW-gBT (flat) | Co(fcc) | 27.8 | 2.5 | -0.2 | 1.8 | 44.4 | 14.8 | 32.0 | 1.47 |
| Ref. Simon et al.,2018 | Co(fcc) | (39.86) | (15.11) | (27.2) | (1.43) |
| (mJ/m2) | (meV) | ||
|---|---|---|---|
| coned | scSOC | 20.3 | 2.02 |
| flat | scSOC | 16.2 | 1.60 |
| 1-shot SOC | 21.9 | 2.18 | |
| Ref. Yang et al.,2015 | 19.0 | 2.17 |
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Comparison of first-principles methods to extract magnetic parameters in ultra-thin films: Co/Pt(111)
Bernd Zimmermann
Peter Grünberg Institut and Institute for Advanced Simulation, Forschungszentrum Jülich and JARA, 52425 Jülich, Germany
Gustav Bihlmayer
Peter Grünberg Institut and Institute for Advanced Simulation, Forschungszentrum Jülich and JARA, 52425 Jülich, Germany
Marie Böttcher
Institute of Physics, Johannes Gutenberg-Universität, 55128 Mainz, Germany
Graduate School Materials Science in Mainz, 55128 Mainz, Germany
Mohammed Bouhassoune
Peter Grünberg Institut and Institute for Advanced Simulation, Forschungszentrum Jülich and JARA, 52425 Jülich, Germany
Samir Lounis
Peter Grünberg Institut and Institute for Advanced Simulation, Forschungszentrum Jülich and JARA, 52425 Jülich, Germany
Jairo Sinova
Institute of Physics, Johannes Gutenberg-Universität, 55128 Mainz, Germany
Stefan Heinze
Institute of Theoretical Physics and Astrophysics, Christian-Albrechts-Universität, 24098 Kiel, Germany
Stefan Blügel
Peter Grünberg Institut and Institute for Advanced Simulation, Forschungszentrum Jülich and JARA, 52425 Jülich, Germany
Bertrand Dupé
Institute of Physics, Johannes Gutenberg-Universität, 55128 Mainz, Germany
Abstract
We compare three distinct computational approaches based on first-principles calculations within density functional theory to explore the magnetic exchange and the Dzyaloshinskii-Moriya interactions (DMI) of a Co monolayer on Pt(111), namely (i) the method of infinitesimal rotations of magnetic moments based on the Korringa-Kohn-Rostoker (KKR) Green function method, (ii) the generalized Bloch theorem applied to spiraling magnetic structures and (iii) supercell calculations with non-collinear magnetic moments, the latter two being based on the full-potential linearized augmented plane wave (FLAPW) method. In particular, we show that the magnetic interaction parameters entering micromagnetic models describing the long-wavelength deviations from the ferromagnetic state might be different from those calculated for fast rotating magnetic structures, as they are obtained by using (necessarily rather small) supercell or large spin-spiral wave-vectors. In the micromagnetic limit, which we motivate to use by an analysis of the Fourier components of the domain-wall profile, we obtain consistent results for the spin stiffness and DMI spiralization using methods (i) and (ii). The calculated spin stiffness and Curie temperature determined by subsequent Monte Carlo simulations are considerably higher than estimated from the bulk properties of Co, a consequence of a significantly increased nearest-neighbor exchange interaction in the Co-monolayer (+50%). The calculated results are carefully compared with the literature.
Density functional theory, spintronics
pacs:
73.20.-r 71.15.Mb
I Introduction
In recent years, non-collinear magnetic structures and in particular skyrmions, have attracted a lot of interest due to their peculiar properties and their technological perspective in the field of information technology Fert et al. (2013). Typically these magnetic structures are stabilized by the competition between the Heisenberg exchange, magnetic anisotropy and dipolar interaction. Recently, the Dzyaloshinskii-Moriya interaction (DMI) Dzialoshinskii (1957); Moriya (1960) has emerged as a new key stabilization mechanism. The DMI arises due to spin-orbit coupling (SOC) and is present in every system, which lacks structural inversion symmetry. The presence of DMI explains the stabilization of skyrmions in bulk B20 alloys such as MnSi Mühlbauer et al. (2009) or in thin films of semiconductor Fe*(1-x)Cox*Si Münzer et al. (2010) and in ultra-thin films at low temperature such as in Fe/Ir(111) Heinze et al. (2011) or Pd/Fe/Ir(111) Romming et al. (2013).
As skyrmions are becoming relevant for technological applications Zhang et al. (2015), additional design goals for skyrmions have been formulated Fert et al. (2017). For example, (i) they should be stable above room temperature, (ii) skyrmions should not be too small ( nm diameter), (iii) skyrmions in ultra-thin films and heterostructures thereof are preferred over skyrmions in bulk samples, and (iv) preferably use materials that are simple to integrate into current manufacturing processes. The latter brings Co/Pt(111) ultra-thin film into play, a material that is well known and is used for perpendicular magnetic recording Yang et al. (2015); Moreau-Luchaire et al. (2016); Boulle et al. (2016). Several recent studies focus on Co/Pt based systemsMaccariello et al. (2018); Chauleau et al. (2018) and explore the possibility to tune the material parameters, e.g. through additional buffer layers Vida et al. (2016); Simon et al. (2018); Jia et al. (2018), alloying Shahbazi et al. (2018); Hanke et al. (2018) or dusting Zimmermann et al. (2018) with a third chemical element.
Obviously, it becomes crucial to understand the stabilization mechanism of skyrmions, predict and design their properties in setups for technological use by theoretical models. Ab initio spin-lattice models proved to be a very powerful approach to realistically describe non-collinear magnets, single skyrmions and skyrmions lattices in experimentally realized systems Dupé et al. (2014); Nandy et al. (2016); Dupé et al. (2016). In such a model, magnetic moments are localized at atomic sites and their interactions are described by parameters, typically the Heisenberg exchange constants, , and the DMI vector, , for two-site interactions ( and label magnetic sites), or the magnetic on-site anisotropy . The parameters for an ab initio spin-lattice model are obtained directly from the total energy of the electronic structure by density functional theory, and the magnetic ground state is found for example by spin-dynamics or Monte-Carlo methods.
While this multiscale approach provides a very efficient description of the energy landscape, the quality of the description depends crucially on the parameters obtained by the mapping of the density functional description of the magnetic states onto the model. In this respect it is important to notice that almost all bulk and interface stabilized skyrmion systems are itinerant magnets, i.e. those electrons responsible for the formation of magnetism also participate in the formation of a complex Fermi surface and hop across the lattice. As a consequence the magnetic interactions are typically not short ranged, as it is often assumed in spin-lattice models. Additionally, the size of the magnetic moments, , is not an integer multiple of the Bohr magneton, but depends on the details of the electronic structure. Most importantly, and in difference to the basic assumptions of many spin models, the size of the moments and the interaction parameters between them depend on their relative orientation. This effect increases with the number of magnetic neighbors and is thus stronger for itinerant bulk magnets than for ultra-thin films. The effect of the change of the magnetic moment for a spin-spiral state with wave vector imposed onto the bulk magnets Cr, Mn and Fe was shown Bihlmayer (2007). In the case of Ni, the magnetic moment can even be completely quenched in the antiferromagnetic state which leads to an overestimation of the total energy Singer et al. (2005). Similarly, the magnitude of induced magnetic moments strongly depends on the spin-configuration of neighboring strong moments Polesya et al. (2010). As a consequence the range of validity of the ab initio spin-lattice model depends crucially on the initial spin configuration for which the parameters entering the spin models have been calculated. The choice of the spin-configuration depends then on the purpose of the application of the spin-lattice model for example to explore the ground states at low or high temperatures or the excited states of small or large spin structures.
In case of the skyrmions, and in particular for skyrmions in the Co thin films on 5d substrates Boulle et al. (2016); Corredor et al. (2017); Perini et al. (2018) or Co/Ru(0001) Hervé et al. (2018), one deals with relatively large non-collinear magnetic textures, with sizes in the order of 100 nm. In this case, the variation of the angle between magnetic moments is small from one atom to the next and the ferromagnetic state (FM) is a good initial state to determine the model parameters. Nevertheless, it is interesting to know how sensitive the model parameters are with respect to the initial state.
In this work, we compare three different approaches to extract the model parameters, namely (i) by the method of infinitesimal rotations based on the Korringa-Kohn-Rostoker (KKR) Green function method, or by mapping the energies of spin-spiral states which are calculated (ii) by using the generalized Bloch theorem (gBT) or (iii) in a supercell geometry. The latter two approaches employ the full potential linearized augmented plane wave (FLAPW) method. We obtain consistent results for the micromagnetic spin stiffness and DMI spiralization in the long-wavelength limit (i.e. around the ferromagnetic state) between KKR and gBT calculations. For larger spin-spiral vectors, the details of the computational procedure (e.g. whether flat or coned spin-spirals are used) become relevant. As a consequence, relying only on one data-point is to be taken with caution, as it is typically done in supercell calculations. Our calculated spin-stiffness is considerably higher (more than 50%) as compared to experimentally determined values for Co-thicknesses below 1 nm. We also find a considerably higher Curie temperature as compared to an estimate based on the Curie temperature of bulk Co., which we trace back to an increased nearest-neighbor exchange interaction in the monolayer as compared to bulk Co.
The paper is structured as follows: in Section II we describe the structure of the Co/Pt(111) system that we study, as well as the magnetic models for which we extract the parameters and detail the computational approaches used. In Section III we present the results and discuss them in comparison to the existing literature, before we conclude in Section IV. The Appendices investigate the dependence of the predicted Curie temperature on the computational procedure (Appendix A), give arguments why the micromagnetic limit is the suitable one for Co/Pt based systems (Appendix B) and investigate the importance of the induced Pt moments for the spin-spiral energy dispersion (Appendix C).
II Methods and Computational Details
We have studied a Co monolayer on Pt(111) in both fcc and hcp stacking positions by means of density functional theory (DFT) calculations. For all calculations, the calculated in-plane lattice parameter of bulk fcc-Pt in local-density approximation (LDA) Vosko et al. (1980) was used ( nm). We use a two-dimensional setup, i.e. embedding a finite number of layers between two semi-infinite vacuum regions.
II.1 Structural Relaxations
The structural relaxations were performed using the DFT package FLEUR FLE employing the full potential linearized augmented plane wave (FLAPW) method. We used a mixed density functional, which was introduced in Ref. De Santis et al., 2007 to treat combined systems of 3- and 5-transition metals. It combines the LDA in the muffin tin (MT) spheres of the 5 atoms and the generalized gradient approximation (GGA) Perdew et al. (1992) everywhere else. It was already used to obtain accurate description of surfaces and interfaces containing 3 and 5 materials Dupé et al. (2014, 2015). We used a cutoff parameter for the basis functions of and 72 -points in one twelfth of the two-dimensional Brillouin zone (BZ), where is the Bohr radius. The symmetric slab was composed of five Pt layers and one Co layer on each side. The positions of Co and the top Pt layers were relaxed for both Co stacking positions. As shown in Table 1, the relaxed structural parameters are basically independent of the Co stacking. Our calculations suggest that the ground state of a Co monolayer on Pt(111) is obtained for the fcc stacking position. Hence, in the analysis presented in Sec. III, we focus mostly on the fcc stacking position.
II.2 Magnetic models
II.2.1 Extended Heisenberg model
Magnetic ultrathin films are well described by the general atomistic extended Heisenberg Hamiltonian,
[TABLE]
where and are the magnetic moments of unit length at position and respectively, are the magnetic exchange parameters, are the Dzyaloshinskii-Moriya vectors and is the onsite uniaxial magnetocrystalline anisotropy.
The extended Heisenberg Hamiltonian associates magnetic moments to atomic sites. It was a priori not designed to study itinerant magnets. However, the inclusion of exchange and DM energy parameters beyond the first nearest neighbor approximation allows for an accurate description of the energy landscape of even frustrated 2-dimensional itinerant magnets on Ir(111) Heinze et al. (2011); Simon et al. (2014); Dupé et al. (2014) and Ir(001) Hoffmann et al. (2015) and on W(110) substrates Hoffmann et al. (2017).
A particularly important subset of non-collinear magnetic states are spin spirals,
[TABLE]
with being the spin-spiral vector, is the position of site , the rotation matrix brings the local axis to the rotation axis of the spin spiral, and is the cone angle. For the special value , we obtain flat spin spirals. Within the spin-lattice model, Eq. (II.2.1), the Heisenberg and DM energy contributions of flat spin-spirals are given by
[TABLE]
with .
One may calculate the and parameters in different ways (see Sec. II.3). In any case, Eq. (II.2.1) is only valid under the condition that the size of all magnetic moments is constant. In DFT, this constraint does not exist and particular care must be taken in the extraction of magnetic exchange interactions if this condition is not respected Ležaić et al. (2013); Polesya et al. (2010).
II.2.2 Micromagnetic and effective model
The micromagnetic model employs a continuous vector field, the magnetization . In case of a thin-film geometry, the micromagnetic energy is given by the functional
[TABLE]
with the exchange stiffness (also termed spin stiffness) , uniaxial anisotropy constant and interfacial DMI constant (or spiralization) , and is the unit vector along the direction perpendicular to the film. The material parameters and are related to the parameters of the atomistic model by Schweflinghaus et al. (2016)
[TABLE]
Here, is the volume of the magnetic part of the unit cell, i.e. in our case the Co monolayer. It may become difficult to define the thickness of a layer in the presence of relaxations (see e.g. the discussion in Ref. Simon et al., 2018). Here, we take as thickness of the monolayer (see Table 1).
Finally, we can express the micromagnetic parameters as effective parameters of a nearest-neighbor Heisenberg model (see Eq. (II.2.1)),
[TABLE]
which is designed to reproduce the energy in the long-wavelength limit, but deviations for other magnetic states are expected.
II.3 Extraction of magnetic interaction parameters from DFT
In the next step, we describe how to extract the parameters for a spin-lattice model from DFT calculations employing three different approaches, which are briefly described here and more detailed information is given in the rest of this subsection. A summary of main differences is presented in Table 2
(KKR) The first approach relies on the KKR method. We perform self-consistent calculation for the ferromagnetic state, possibly including SOC. The change in total energy due to infinitesimal rotations of the magnetic moments is related to the Heisenberg exchange parameter Liechtenstein et al. (1987) and DMI vectors Udvardi et al. (2003); Ebert and Mankovsky (2009). This approach is considered to give very accurate model parameters for large non-collinear magnetic textures. 2. 2.
(FLAPW–gBT) Secondly, we work in reciprocal space by means of spin-spiral states, Eq. (2), employing the FLEUR code FLE using the generalized Bloch theorem Kurz et al. (2004) (gBT). Self-consistent calculations (without SOC) are performed for any arbitrary spin-spiral wave-vector . We determine the parameters from a fit to the spin-spiral energies for various spin-spiral vectors in the whole Brillouin zone (or high-symmetry lines thereof). The perturbative inclusion of SOC can be analogously related to the atomistic parameters. This method is flexible in the sense that it allows to access both regimes of slowly and fast rotating non-collinear magnetic structures more realistically through self-consistent calculations. 3. 3.
(FLAPW–supercell) The third approach is conceptually similar to the FLAPW-gBT approach as it compares the energy of different spiraling magnetic states, but now in a supercell geometry. The direction of magnetic moments is fixed by constraining fields. In principle, no perturbative treatment of either SOC or the non-collinearity is needed, but due to the fast increase of computational costs with system size, this approach is practically restricted to small supercell sizes and hence large spin-spiral wave vectors . We use two different magnetic configurations in a supercell, one reproduces the setup of Yang et al. Yang et al. (2015), and the other employs small cone-angles to keep a nearly ferromagnetic alignment of Co-moments.
The parameters and for a micromagnetic model are either evaluated using Eqs. (6) and (7) in case of the KKR method, or can be directly obtained from ab-initio calculations as quadratic and linear terms of the spin-spiral dispersion curve around the ferromagnetic state from FLAPW-gBT calculations Bode et al. (2007); Heide et al. (2008, 2009); Zimmermann et al. (2014).
In all approaches, we used the LDA exchange-correlation functional in the parameterization of Vosko et al. Vosko et al. (1980). For the Co/Pt(111) system, we modeled the Pt substrate by 5 layers taking the structural relaxations as described in Sec. II.1.
II.3.1 The KKR method employing infinitesimal rotations
We use the full-potential Korringa-Kohn-Rostoker Green function (FP-KKR-GF) method juK ; Bauer (2013) to converge the charge and spin densities in scalar-relativistic approximation for the ferromagnetic state using -points in the full Brillouin zone. In a next step, we obtain Heisenberg parameters and DM vectors in real-space by relating the change in energy of infinitesimal rotations of the magnetic moments at lattice sites and Liechtenstein et al. (1987); Ebert and Mankovsky (2009), including spin-orbit coupling only in this step (one-shot SOC). Due to the infinitesimal rotations, this approach should be best to obtain parameters for large non-collinear structures, in particular skyrmions. We obtain interactions for Co pairs up to a distance of seven in-plane lattice constants and sample the full Brillouin zone by -points for this step. We truncate the expansion of the scattering wave-functions into spherical harmonics at . The energy contour integration includes a Fermi-function with an electronic temperature of 473 K, five Matsubara poles and is sampled by 39 points.
This methods also yields interaction parameters between strong Co moments and induced Pt moments. We can simply include these contributions in Eqs. (3), (4), (6) and (7) which implicitly assumes that the size of the induced moments is rigid. A more sophisticated treatment is described in Appendix C. As it turns out, the simple approach is sufficient for our considerations.
Another advantage of this method is the possibility to include SOC self-consistently in all steps and obtain an even more realistic set of parameters. However, the modification of the DM vectors is small (less than 0.1 meV for the nearest-neighbor interaction, ) and therefore a perturbative treatment of SOC is justified to obtain the DMI. We refer to Appendix C for details.
II.3.2 The FLAPW method employing spin-spiral states
The second approach employs the FLAPW method as implemented in the DFT code FLEUR FLE utilizing the generalized Bloch theorem (gBT) Kurz et al. (2004). The main goal is to obtain the total energy of spin-spiral states, Eq. (2), and extract parameters for a magnetic models thereof. By virtue of the gBT, even long wave-length spin spirals meaning small values for (the spin-spiral period length is given by ), can be treated very efficiently in the chemical unit cell, i.e. without the need of large supercells.
Depending on the fitting procedure, we can either extract parameters for long wave-length magnon excitations (low- region) or for spin spiral ground states (large- region). To extract the Heisenberg parameters, we converge the total energy of spin-spiral states in scalar-relativistic approximation using a -point mesh to an accuracy of 0.01 meV. The size of the basis set was determined by , with the Bohr radius .
In a second step, we include spin-orbit coupling (SOC) via first order perturbation theory Heide et al. (2009). The energy correction is used to extract either the micromagnetic DMI spiralization Freimuth et al. (2014) or the two-site DM vectors for a spin-lattice model.
It is important to notice that the dispersion curves scale with the cone angle as by . Hence, the energies for coned spin spirals as we use them later () are a factor of approx. 40 smaller as compared to flat spin-spirals, which requires quite high computational cutoffs.
II.3.3 The supercell method
Additionally, we perform supercell calculations using the FLEUR code FLE to extract the DMI. The supercell contains four magnetic atoms and is spanned by the basis vectors and (see Fig. 1). We then impose magnetic states onto the supercell:
We mimic a flat spin-spiral (see Eq. (2) with ), i.e. the direction of magnetic moments is given by , where labels the position of the magnetic atom in the supercell. This choice reproduces the setup of Yang *et al. * Yang et al. (2015) and is shown in Fig. 1(a). 2. 2.
We extend these states to coned spin-spirals [see Fig. 1(b)], where the direction of moments is given by
[TABLE]
where is the so-called cone-angle. The flat spin-spiral is reproduced for , and all moments point along the direction for .
In both cases, the lines of ferromagnetically aligned moments are parallel to the direction. This corresponds to spin-spiral states described by [see Fig. 1(a)].
A nearest neighbor DMI, , is derived by comparing the energies of right rotating [defined by Eq. (13)] and left rotating magnetic structure [which is obtained by replacing in Eq. (13)],
[TABLE]
We note that the parameter implicitly depends on the choice of through the self-consistent charge and spin densities which enter the Kohn-Sham Hamiltonian of DFT. On the one hand, we expect little modifications of the densities as compared to the FM state if . This configuration is very similar to a perturbative calculation in the vicinity of the FM state. On the other hand, we expect larger changes for .
Concerning the computational details, a basis-set cutoff of was chosen, and the Brillouin-zone was sampled by and -points for coned and flat spin-spirals, respectively. The spin-orbit interaction was included self-consistently in these calculations.
III Results
III.1 Spin stiffness and Heisenberg exchange
Our main results are displayed in Fig. 2 and summarized in Table 3. Let us first compare the non-relativistic spin-spiral dispersion curves in Fig. 2(a-b): KKR and FLAPW data agree very well in the parabolic region around the point, i.e. in the region of interest when large non-collinear magnetic textures are studied. We zoom into this regime and plot the data as a function of in Fig. 2(b) to obtain the spin stiffness (or effective Heisenberg exchange ). Fits to the data of flat or coned spin-spirals yield the same spin stiffness of about 44 pJ/m, which converts to meV (see Table 3). The KKR-derived values are 10% smaller as compared to FLAPW results. To estimate the role of the induced Pt moments for the energies from KKR, we switch their contributions manually off and see basically no change in the spin stiffness as compared to the full calculation. However, taking only the nearest-neighbor interactions of Co into account underestimates the spin stiffness by about 30% [see Fig. 2(b)].
As the -vector approaches the Brillouin zone boundaries at the and points, larger deviations of about 70 meV in energy between flat spin-spirals from FLAPW-gBT calculations and KKR-deduced values occur [see Fig. 2(a)]. The reason lies partly in a different electronic structure: a frozen DFT effective potential from the ferromagnetic state (i.e. =0) is used in the KKR approach, whereas the potential for the spin-spiral is subject to change during the self-consistency as witnessed by a reduction of the Co moment (up to ) and quenching of induced Pt-moments for large (see Fig. 3). In order to perturb the ferromagnetic state less, FLAPW-gBT calculations with coned spin-spirals are preferable. Indeed, the magnetic moments stay practically constant as function of spin-spiral vector , and the spin-spiral energies get corrected to lower values towards the KKR result [see Fig. 2(a)], but a rather large discrepancy (of about e.g. 30 meV at the -point) remains. In order to test whether technical differences between the KKR and FLAPW methods, such as the division of space into Voronoi cells for KKR as opposed to muffin-tin spheres and interstitial region for FLAPW, induce such a discrepancy, we performed self-consistent KKR calculations for the -point (i.e. the row-wise anti-ferromagnetic state) in a simple -supercell geometry. We obtain a perfect agreement to the FLAPW result [compare the green square on top of the blue star at the point in Fig. 2(a)], highlighting the consistency of the two methods if similar approximations are employed. We can finally speculate that the energy difference of about 30 meV at the -point between coned spin-spiral and KKR-deduced energies might be caused by higher order magnetic exchange interactions, such as the four-spin-three-site-interactionKrönlein et al. (2018), which are not captured by the KKR approach based on infinitesimal rotations, but naturally included in the self-consistent spin-spiral calculations. If present, it leads to an effective renormalization of the ’s (concretely, would be increased by 30%), which also reflects in an increased Curie temperature as compared to KKR (see Appendix A).
III.2 Dzyaloshinskii-Moriya interaction
Let us now turn to the DMI [see Figs. 2(c-d)]. We obtain a very good agreement throughout the whole Brillouin zone between spin-spiral energies from KKR and FLAPW-gBT calculations with coned spirals (compare blue line and purple dots). Flat spin-spirals yield considerably higher DM-energies for -vectors in the middle of and , respectively. This discrepancy is again a reflection of the different electronic structure for large -vectors and links directly to the changed hybridization between Co and Pt Sandratskii (2017). We can conclude that it might fail to take a single-point calculation from a rather large value and infer information about the low- regime from this point. In the micromagnetic limit, however, coned and flat FLAPW-gBT spin-spiral calculations yield a similar DMI of (which converts to ; see dashed line). The KKR-obtained data is also well approximated by this line [see Fig. 2(d)]. If we carefully converge the DMI spiralization from KKR according to Eq. (7) with respect to the number of neighbors, we obtain a 20% (10%) higher DMI as compared to FLAPW-gBT values including (neglecting) the interactions between Co and induced Pt moments. This enhancement is due to the fact that by evaluating Eq. (7), we calculate the exact ()-limit. Indeed, if we perform a linear fit to the KKR-derived DMI energies in an interval of , similarly to the procedure used for FLAPW-gBT calculations, we obtain an identical spiralization (up to two digits) as compared to the FLAPW-gBT calculations. Simon et al. Simon et al. (2018) use an analogous fitting-procedure with pair-interaction parameters as calculated by the KKR method, employing a different code, and obtain a value in perfect agreement to our KKR-derived values (see Table 3).
The fact that our DMI data exhibits the same slope from towards and , respectively, shows that the DMI spiralization is isotropic, which is expected for a system with symmetry Bogdanov and Yablonskii (1989). The positive sign of the slope corresponds to a lowering of energies of magnetic states with left-handed (anti-clockwise) chirality.
The effective model, which is obtained from the fits in the low- region, reproduces the energies from coned spin-spirals and via KKR rather satisfactorily even for large (see solid black line in Fig. 2c). We observe pronounced deviations from the simple sine-behavior in the middle of the Brillouin zone for the KKR-obtained curves, which stem from contributions beyond the nearest neighbors. Near the -point (which represents the non-collinear Néel state for a flat spin-spiral), a qualitative change of the energy dispersion is caused by interactions beyond nearest neighbors and cannot be captured by the effective model. Similarly, a different slope between effective model and KKR-derived energy curves appears near the -point, highlighting the limitations of the effective model when extrapolating from the low- region to arbitrary -vectors.
Next, we discuss the Dzyaloshinskii-Moriya interaction as determined from supercell calculations and take the effective model with values as fitted from FLAPW-gBT calculations as reference. A coned supercell-spiral yields a DM energy very close to the corresponding FLAPW-gBT calculation at [see Fig. 2(d)]. Despite this good agreement for this single point (5.7 vs. 6.1 meV), the inferred is 40% higher as compared to the values from fits in the low- region of FLAPW-gBT calculations (see Table 4), which again emphasizes the difficulty of relying on only one data point. The DMI coefficient of a flat supercell-spiral agrees surprisingly well with our reference value (1.60 as compared to 1.43 meV). However, the corresponding energy is considerably lower than the spin-spiral dispersion of FLAPW-gBT calculations for this state [see Fig. 2(d)]. We can trace this discrepancy back to the different treatment of SOC, which is included self-consistently in the supercell calculations but only in first-order perturbation theory in the FLAPW-gBT calculations: if we include SOC in the supercell calculation by a force-theorem approach, i.e. performing a single iteration with SOC on top of a converged scalar-relativistic calculation and compare the sum of single-particle energies, the DM energy agrees within 5% with the gBT calculation (6.5 vs. 6.8 meV). Similar findings have been reported for freestanding Fe/Ir bilayers Meyer et al. (2017). The energies of Yang et al. Yang et al. (2015) (taking their values from Co(1)/Pt(3), where the number in parentheses denotes the number of atomic layers) are 35% higher than our corresponding supercell calculation (compare 2.17 to 1.60 meV), probably due to the different number of substrate layers and different relaxations of atomic positions, and about 50% higher as our reference value from the low- regime. The good agreement between the value of Ref. Yang et al., 2015 and our 1-shot SOC supercell calculation (see Table 4) is due to an accidental cancellation of errors. Note that the authors of Ref. Yang et al., 2015 used a different magnetic volume (corresponding to the volume an atom in fcc-Pt) for the conversion between and (cf. Eq. (9)), which leads to a different ratio of as compared to us.
III.3 Discussion
Overall, our results emphasize the compatibility of the FLAPW-gBT and KKR methods in the relevant micromagnetic limit, as well as the flexibility of the generalized Bloch theorem approach as it can access states beyond the micromagnetic limit more realistically through self-consistent calculations. Our findings regarding the DMI for the present case of a Co monolayer on Pt(111) are in satisfactory agreement to previous studies on this system Simon et al. (2018); Yang et al. (2015); Zimmermann et al. (2018).
Comparing our value for the spin stiffness to the literature, it is considerably higher than most of the experimentally determined values for Co: for bulk-Co a stiffness of about 15 pJ/m (fcc) and 30 pJ/m (hcp-Co) is measured by various methods Hüller (1986), a value of 21 pJ/m was measured for 10 nm thick Co-films Eyrich et al. (2012), and even smaller stiffnesses as low as 14 pJ/m are reported for ultra-thin Co films down to thicknesses of 0.5 nm. The reported spin stiffness by Boulle et al. Boulle et al. (2016) represents an exception to this trend, as they determine 27.5 pJ/m for Co-thicknesses below 1 nm. Still, our results are nearly 50% higher than this value, and they are in line with previous DFT calculations on a Co monolayer on Pt(111) Vida et al. (2016); Simon et al. (2018).
One possibility for the discrepancy between experiment and theory might arise from the fact that experiments are necessarily performed at finite temperatures, whereas the ab-initio calculations neglect any spin-fluctuations. In simplest approximation, the spin-stiffness and DMI spiralization are renormalized by , where is the temperature-dependent magnetization and the saturation magnetization Rózsa et al. (2017). Inserting the values for at K as deduced from Monte Carlo simulations [see Appendix A and Fig. 4(b)], we arrive at a renormalization factor of 0.67 and a room-temperature spin-stiffness of 29.7 pJ/m (using results based on parameters from FLAPW-gBT calculations), which is in much better agreement with the experimentally determined spin-stiffnesses.
The critical temperature at which the magnetization vanishes, , was also extracted from the Monte Carlo simulations. Depending on the parameterization, these values are ranging from 500 K to 620 K for the KKR and the FLAPW-gBT parameterization, respectively. These critical temperatures seem to be quite high for a two-dimensional (2D) system: They are up to 40% higher than one would estimate based on the Curie temperature of bulk Co, K, as calculated by
[TABLE]
with Boltzmann’s constant and uniaxial anisotropy meV. Eq. (15) is based on an renormalization group analysis and assumes unchanged magnetic interaction parameters when going from 3D to 2D (see Ref. Blügel and Bihlmayer, 2007 and references therein for details). But we infer that meV of a monolayer Co/Pt(111) is about 50% higher as compared to hcp-Co, which we determined to meV by the KKR method, in very good agreement to 14.8 meV as determined by LMTO calculations Pajda et al. (2001). This increased explains the discrepancy between Eq. (15) and Monte-Carlo simulations for the monolayer. This difference of the exchange interactions stems from changes in the electronic structure and corrects the estimated value to higher temperatures. Inserting in Eq. (15) the Co-bulk Curie temperature of a hypothetical Co solid with an increased Curie temperature of K corresponding to the nearest neighbor exchange interaction of meV we obtain an estimated of 630 K.
Indeed, this corrected estimate and the Monte-Carlo determined values are well compatible to magneto-optical Kerr effect experiments on a monolayer Co on Pt(111)Shern et al. (1999): a finite magnetization is present up to 623 K. However, the experimental situation is more complex, as the formation of a CoPt surface-alloy is observed above 500 K, and the reported 623 K is an increased value compared to the ideal Co-monolayer caseShern et al. (1999). Consequently, the critical temperature of the ideal Co monolayer on Pt(111) lies in between 500–623 K. All our Monte-Carlo-predicted critical temperatures are compatible with the experiment, whereas estimated from the bulk properties is not. This highlights that a prediction of the Curie temperature of thin films purely from the knowledge of bulk properties is to be taken with caution and possibly fails. Overall, our calculations confirm that the Curie temperature of Co thin-films on Pt(111) lies far above room temperatureChappert et al. (1998); Hashimoto et al. (1990), even for the monolayer.
Having obtained consistent values for the spin-stiffness and DMI spiralization, and considering additionally the uniaxial anisotropy (which converts to ), which we obtained by FLAPW calculations with out-of-plane easy axis, we can determine the ground state. Due to the rather strong magnetocrystalline anisotropy, we obtain a ferromagnetic ground state111The quantity denotes a ferromagnetic ground-state for Zimmermann et al. (2014). With our values for Co/Pt(111), we obtain .. We can in hindsight check the validity of the micromagnetic limit by estimating the domain-wall width
[TABLE]
and with the arguments presented in Appendix B, we are indeed in the micromagnetic limit ( ). Our calculated domain-wall width is in perfect agreement with the experimental value of about 4 nm in monolayer Co/Pt(111) as measured by spin-polarized scanning-tunneling microscopy (SP-STM) at low temperatures Meier et al. (2006). Very recently, a domain-wall width of nm in slightly thicker (1.4 nm) epitaxial Co layers on Pt(111) has been measured by SEMPA at room temperature Corredor et al. (2017).
Another characteristic length scale for spiraling magnetic textures induced by DMI is set by , taking the values from Table 3, which is even an order of magnitude larger than the domain-wall width and approximates the micromagnetic limit even more.
Finally, we have determined the magnetic interaction parameters for the spin-lattice model for the hcp stacking position of Co on Pt(111) and present them in Table 3. The changes are marginal as compared to the fcc stacking position, in agreement to Ref. Vida et al., 2016.
IV Conclusion
We have determined the magnetic interaction parameters for a Co monolayer on Pt(111) by three distinct approaches, (i) performing infinitesimal rotations, (ii) using spin-spirals employing the generalized Bloch theorem for various vectors, (iii) and constraining spin-spirals into a rather small supercell. We obtain consistent results for the spin stiffness and the Dzyaloshinskii-Moriya interaction in the long-wavelength (micromagnetic) limit around the ferromagnetic state using methods (i) and (ii). When going to higher spin-spiral vectors, deviations through differences in the electronic structure play a role for flat spin spirals. In order to still realize the long-wavelength limit in the supercell approach, we propose to use a coned spiraling structure with a small cone angle, which leaves e.g. the magnetic moments unchanged as function of the spin-spiral vector (being inversely proportional to the supercell size). We found that the micromagnetic DMI spiralization might not be accurately inferred by extrapolating from one data point obtained for a large -vector (i.e. small supercell).
Acknowledgments
We acknowledge discussions with Stanislas Rohart. B. Dupé thanks DARPA TEE program (#HR0011727183) and the (DFG) via the project DU 1489/3-1 for financial support. B. Dupé and S. Heinze thank the Deutsche Forschungsgemeinschaft (DFG) for financial support via the project DU 1489/2-1. B. Zimmermann, S. Heinze and S. Blügel thank the European Unions Horizon 2020 research and innovation programme under grant agreement No. 665095 (FET-Open project MAGicSky) and S. Blügel thanks DARPA TEE program (#HR0011831554) from DOI, Deutsche Forschungsgemeinschaft (DFG) through SPP 2137 “Skyrmionics”, and the Collaborative Research Center SFB 1238 (Project C1) for funding. S. Lounis and M. Bouhassoune acknowledge support from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (ERC-consolidator grant 681405 – DYNASORE). We gratefully acknowledge computing time at the HLRN supercomputer as well as JURECA supercomputer (projects jias1a and jias1f) of Jülich Supercomputing Centre and JARA-HPC of RWTH Aachen University.
Appendix A Thermodynamical study of Co/Pt(111)
As discussed in the main text (see Sec. III.1), the KKR as well as FLAPW-gBT methods employing flat or coned spin-spirals all leave the curvature of the spin-spiral dispersion around the -point unchanged. However, large changes occur for larger -vectors corresponding to strongly non-collinear magnetic textures. These states will only be present at high temperature or after an exciting stimulus with high energy. Therefore, to understand the consequences of the different parameterizations of the Heisenberg Hamiltonian on these states, we have performed atomistic Monte Carlo simulations and calculated the critical temperature TC of a Co monolayer in the fcc stacking position on Pt(111).
We took the parameterization of the extended Heisenberg model as determined from KKR and FLAPW-gBT calculations (cf. Table 3), and included an uniaxial anisotropy of meV which was determined by FLAPW calculations. Subsequent Monte Carlo simulations have been performed using a supercell on 192 replicas of the supercell with temperatures geometrically spaced between and Kelvin. During the Metropolis steps, we have averaged over 600 steps the total energy as obtained from Eq. (II.2.1), the total magnetization as well as the specific heat and the magnetization susceptibility .
Fig. 4 shows the results of the Monte Carlo simulations for the three different parameterizations. The black, red and blue curves correspond to the parameterization obtained with via KKR and FLAPW-gBT methods employing coned and flat spin spirals, respectively. Fig. 4(a) shows the averaged total energy as a function of temperature. All energies show a smooth increase with a change of slope at around 500 K which indicates a phase transition. As the total energy of the magnetic lattice increases, the magnetization decreases as displayed in Fig. 4(b): The magnetization decreases linearly up to the critical temperature and then drops abruptly down to zero. The drop occurs at a temperature of 600 K for both FLAPW-gBT parameterizations and at 500 K for the KKR parameterization.
To characterize the phase transition in more detail, we calculate the specific heat and the magnetic susceptibility [cf. Fig. 4(c) and (d), respectively] as explained in Ref. [Böttcher et al., 2018]. When DMI is present in a magnetic system, the phase diagram shows an intermediate region between the paramagnetic phase at high temperature and the long-range ordered phase at low temperature. In this intermediate regime, the number of skyrmions can fluctuate Rózsa et al. (2016); Hou et al. (2017) and may be used to nucleate skyrmions in magnetic multilayers by current pulses Lemesh et al. (2018). This region is present between the critical temperature extracted from the location of the peak of the specific heat, , and the critical temperature extracted from the peak of the magnetization susceptibility, Böttcher et al. (2018).
The specific heat curves show a peak at 600 K for the parameterization from KKR and at 740 K and 800 K for the ones from FLAPW-gBT employing conical and flat spin spirals, respectively. The subsequent phase transition to the long-range ordered phase occurs at 500 K (KKR), 600 K (FLAPT-gBT coned) and 620 (FLAPW-gBT flat). The predicted transition temperatures from the KKR parameterization are lower than the corresponding ones from FLAPW-gBT, with differences of up to 200 K (or about 25%).
Appendix B Arguments for the micromagnetic limit for Co/Pt interfaces
In order to determine the relevant range of spin-spiral vectors for Co/Pt, we perform a power spectrum analysis of the well known domain wall profile. The magnetization of a Néel-type domain wall separating the ferromagnetic domains is given by , with
[TABLE]
where is the domain-wall width.
The Fourier transformation of the component is directly evaluated as
[TABLE]
For the -component, we need to regularize the Fourier transformation and obtain
[TABLE]
Both Fourier transformed components peak around (see Fig. 5), and are significant only for
[TABLE]
This interval has been chosen such that it contains more than 90% of the components . Inserting a typical value for the domain-wall width in Co/Pt based systems () yields
[TABLE]
with the in-plane lattice constant of the Co/Pt(111) film, . These very low values justify the application of the micromagnetic approach for Co/Pt based systems.
Appendix C Treatment of induced Pt moments in non-collinear structures and self-consistent inclusion of SOC
The size of the induced Pt moment can vary strongly in non-collinear structures, as shown by results for flat spin-spirals in Fig. 3(b). To account for such a variation in the KKR calculations, a scheme based on the non-local spin-susceptibility was proposed Polesya et al. (2010), where the size and direction of the induced Pt moment is determined by the orientation of the nearest-neighbor Co-moments,
[TABLE]
The superscript ‘FM’ denoted the moment obtained for the ferromagnetic state.
The induced Pt moment approximated in this way with a values for as calculated by KKR reproduces the results from flat spin-spiral calculations using the FLAPW-gBT method rather well [compare the dashed line and green squares in Fig. 3(b) in the main text]. Remaining deviations are probably due to the variation of the size of the Co moments in FLAPW-gBT calculations (see Fig. 3(a)].
As a next step, we can evaluate the spin-spiral dispersion curves by renormalizing the KKR-obtained parameters and by . However, this treatment does not significantly change the energies, as shown in Fig. 6. In particular, the relevant low- regime around is not affected.
We also tested whether the inclusion of SOC in the self-consistent determination of the charge and spin-density affects the and parameters as determined by infinitesimal rotations. The resulting changes in and (or and are 3% at most, and derived spin-spiral energies are hardly distinguishable as compared to the (1-shot SOC)-procedure employed in the main text (compare dotted to full lines in Fig. 6).
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Fert et al. (2013) A. Fert, V. Cros, and J. a. Sampaio, Nat. Nanotechnol. 8 , 152 (2013).
- 2Dzialoshinskii (1957) I. E. Dzialoshinskii, J. Exp. Theor. Phys. 5 , 1259 (1957).
- 3Moriya (1960) T. Moriya, Phys. Rev. 120 , 91 (1960).
- 4Mühlbauer et al. (2009) S. Mühlbauer, B. Binz, F. Jonietz, C. Pfleiderer, A. Rosch, A. Neubauer, R. Georgii, and P. Böni, Science 323 , 915 (2009).
- 5Münzer et al. (2010) W. Münzer, A. Neubauer, T. Adams, S. Mühlbauer, C. Franz, F. Jonietz, R. Georgii, P. Böni, B. Pedersen, M. Schmidt, A. Rosch, and C. Pfleiderer, Phys. Rev. B 81 , 041203(R) (2010) . · doi ↗
- 6Heinze et al. (2011) S. Heinze, K. von Bergmann, M. Menzel, J. Brede, A. Kubetzka, R. Wiesendanger, G. Bihlmayer, and S. Blügel, Nature Physics 7 , 713 (2011).
- 7Romming et al. (2013) N. Romming, C. Hanneken, M. Menzel, J. E. Bickel, B. Wolter, K. von Bergmann, A. Kubetzka, and R. Wiesendanger, Science 341 , 636 (2013).
- 8Zhang et al. (2015) X. Zhang, M. Ezawa, and Y. Zhou, Scientific Reports 5 , 9400 (2015).
