Demagnetizing fields in chiral magnetic structures
M. A. Moskalenko1'2, I. S. Lobanov1'3, V. M. Uzdin1'3 xITMO University, 197101 St. Petersburg, Russia 2 Science Institute, University of Iceland, 107 Reykjavk, Iceland 3Department of Physics, St. Petersburg State University, St. Petersburg 198504, Russia moskalenko.mary@gmail.com, lobanov.igor@gmail.com, v_uzdin@mail.ru
PACS 75.10.Hk, 75.40.Mg, 75.78.Cd DOI 10.17586/2220-8054-2020-11-4-401-407
A method for calculating the magnetic dipole-dipole interaction in topological magnetic systems has been developed. It can be used to calculate stable states and minimum energy paths that determine the magnetic transition in chiral magnetic structures. Instead of directly summing the dipole interactions between magnetic moments/magnetic elements, we solve a local equation for demagnetizing fields. The states corresponding to the local energy minimum can be found using the Lagrange method for the conditional extrema. The efficiency of the algorithm has been demonstrated by calculating the dependence of the size and shape of magnetic skyrmions and anti-skyrmions on the magnitude of magnetization.
Keywords: topological magnetic structures, skyrmions, demagnetizing fields.
Received: 1 August 2020
Revised: 18 August 2020
1. Introduction
Chiral magnetic structures have attracted much attention in recent years due to the discovery in these systems of localized non-collinear states that can move very rapidly under the action of a spin-polarized electric current [1,2]. In micromagnetic models, these states cannot be destroyed by continuous transformation of magnetization, and it is believed that they are topologically protected from thermal fluctuations [3,4]. In real systems magnetic moments localized on the site of discrete lattice and it is possible to say only about topological stabilization. In lattice model stability of such system can be estimated on the basis of transition state theory (TST) for magnetic degrees of freedom [5,6].
Usually, the properties of skyrmions (Sk), as well as other topological structures, are studied within the framework of the Heisenberg-type Hamiltonian, which includes exchange, anisotropy, external magnetic field, and the Dzyaloshinskii-Moriya interaction (DMI) [5]. All these interactions decay rapidly with distance. This model correctly describes, for example, small Sk in thin PdFe films on the Ir (111) surface, observed experimentally [7]. However, such Sk are stable only at very low temperatures of the order of 10 K. Sk stable at room temperature in ferromagnetic (FM) materials are larger in size, and to describe their behavior it is necessary to take into account the magnetic dipole-dipole interaction which is responsible for the creation of demagnetizing fields. In the first approximation, this interaction can be introduced into the theory by renormalizing the anisotropy parameters [8]. However, more accurate calculations of the energy associated with demagnetizing fields are needed to describe the topological structures of micron sizes [9] and especially three-dimensional magnetic systems [10,11]. Moreover, the stability of Sk with a size of 100 nm and more, according to [12], is always determined by the dipole-dipole interaction. So, if we want to have small stable Sk as future bits of computer memory, we must accurately evaluate the contribution of the demagnetizing fields. The development of methods for manipulating the dipole interaction in magnetic micro- and nanostructures and theoretically estimating its contribution to energy is very important for practical applications, since it allows one to control the size and stability of topological systems. Suppression of the dipole interaction in antiferromagnetic (AF) and ferrimagnetic chiral structures, for example, made it possible to create Sk of about 10 nm in size that are stable at room temperature [13]. Small room temperature Sk in artificial antiferromagnets, also use the reduction of stray fields in multilayer systems with AF interlayer coupling [14].
The lifetime of the topological states of systems without strong dipole interaction at an arbitrary temperature can be evaluated using standard TST approach. It was done for Sk in AF [15] and for small-size Sk in ferromagnetic (FM) materials [6]. Such calculations for FM micromagnetic structures are a challenging problem due to the number of degrees of freedom and need to correctly take into account the dipole-dipole interaction. There are effective methods for calculating the saddle points on the energy surface and the activation energy of Sk annihilation for systems containing millions of magnetic moments [16]. But even the calculating the energy surface for systems with a significant contribution from the dipole-dipole interaction and demagnetizing fields is still a very difficult task.
Below, we present an efficient algorithm for calculating the demagnetizing field and show its efficiency by the example of calculating Sk and antisyrmions (ASk) in a FM system with different saturation magnetization.
2. Micromagnetic and discrete models
The magnetic texture in the micromagnetic approach will be described by the vector field m(r), where r is the coordinate of a point in the magnetic sample, and m is the unit vector in the direction of the local magnetization at this point. The local magnetization value is a fast variable and is assumed to be equal to the fixed saturation magnetization Ms. Below we consider a thin film modeled by the part of the x-y plane Q = [0, d]2 with periodic boundary conditions. The energy E of the magnetic state is expressed in terms of the energy density w as follows:
E = J w(r)dr, W = Wex + WDMI + WK + Wdemag,
n
where we take into account contribution of exchange wex, DMI wDMI, anisotropy wK and demagnetizing field Wdemag. Assuming that the exchange is isotropic:
dm\ 2 / dm\ 2
wex = tA(Vm)2, (Vm)2 = —J + ^ j ,
here A is the exchange stiffness and t is the film thickness. For simplicity, consider uniaxial anisotropy with an easy axis perpendicular to the film (parallel to the z axis):
wK = —tK(m • Z)2,
where parameter of anisotropy K > 0 . The DMI energy is expressed in terms of Lifshitz invariant
t (i) dmk dmk . x , t, T(y) \
L)k = — mk~df ' WDMI = t(DxLXz) + Dy Lyyz)).
If Dx = Dy = D, meta(stable) states of Sks can exist in the system, if Dx = —Dy = D, then ASks can be formed. The exchange stiffness and the anisotropy constant are independent of the magnetic texture and are the same for both Sk and ASk structures. Micromagnetic parameters are selected according to [17] to be
t = 0.6 (nm), A =16 (pJ/m), K = 200 (kJ/m3), D = 2 (mJ/m3). (1)
Saturation magnetization Ms varies from 103 to 4 • 104 (A/m).
The energy of the demagnetizing field is determined by its density:
= _ m H
wdemag 2 Msm- • Hdemag:
where Hdemag is the demagnetizing field, which can be found from the Maxwell equations:
V • B = 0,
V X Hdemag = °
where B = ^0(Hdemag + Msm) is the induction of the demagnetizing field. The second condition can be satisfied introducing the potential $ of the field: Hdemag = —V$. From the first condition we have:
V • Hdemag = —MsV • m.
Therefore the potential $ can be found from the Poisson equation:
V2$ = MsV • m. (2)
It is worth noting that in three dimensions the equation can be solved in terms of the Green's function:
1
/
$(r) = -Msi Vr', ^ m(r') dr' = Ms i m(r') • r r dr'. W sj 4n|r - r'| sj 4n|r - r'|3
Restoring the demagnetizing field from the potential and substituting the result into the energy density, we obtain the standard expression for the energy of the dipole-dipole interaction:
2 f 3(m(r/) • A)(m(r) • A) — m(r/) • m(r) J / r — r
wdemag (r) = —M0Ms J -^ — r/|3-dr , A = ^"TTj . (3)
To compute the dipole-dipole energy, it is necessary to sum over all pairs of magnetic moments in the system, which is computationally expensive. To overcome this difficulty, the Fourier transform is often used, which converts the convolution in the expression for energy to multiplication with a function; however, the approach is applicable only to simple domain such as a square or a cube. Another method for calculating the demagnetizing energy is to obtain
the demagnetizing field numerically, which in many cases is faster than using the Green's function. This approach has been previously used with finite element discretization. Below we develop a finite-difference discretization method that leads to an lattice-type model of the magnetic system.
The state is specified by the directions of the magnetic moment Sn, determined at the sites of the square lattice n. Let's denote by rn the site position n. The energy in the lattice representation contains the same contributions as above:
E = Eex + EDMI + EK + Edemag •
Each energy contribution is an approximation of the integral in the micromagnetic model above, rewritten in terms of the values of the magnetic moments m « Sn localized near points r « rn. For example energy of anisotropy becomes:
Ek = -K ^(Sn ■ 2)2.
n
Anisotropy constant in lattice model K is proportional to K, but it also depends on the size of integration cell. If n are points of a square lattice with lattice constant a, then K = Ka2t.
Derivatives in expressions for wex and wDMI can be estimated by finite differences:
V7 / N I Sn1 — Sn Sn2 — Sn
Vm(r„) = -,-.
aa
where Sn1 and Sn2 are nearest neighbors magnetic moment to Sn along x-axis and y-axis, respectively. Eliminating constant addenda in exchange energy contributions, we get:
"2 (1 — Sn1 ■ Sn) ; • • •
\ dx J a2 the Heisenberg exchange energy is obtained in the following form:
Eex = ] Sn ■ Sk,
<n,k>
where the Heisenberg exchange constant J = 2At, and the sum is taken over all pairs (n, k) of magnetic moments (each pair is taken only once). Similarly, the DMI energy can be expressed as follows:
EdMI = — Dn,k ■ (Sn x Sk),
<n,k>
where all DM vectors Dn k have the same length D > 0 and are directed along y, if k is on the right to n, and along —x, if k is above n; here x, y are basis vectors of x and y axes. The DM vectors are connected to the micromagnetic model by the identity |D| = Da. In simulation we used 100 x 100 square lattice with the following parameters:
a = 4 (nm), J = 1.92 ■ 10-20 (J/bond), D = 4.8 ■ 10-21 (J/bond), K = 1.92 ■ 10-21 (J/spin). (4)
The demagnetizing energy is discretized accordingly to:
E = _ MoM H S
Edemag 2 ' Hn Sn,
n
where m = Msa21 is value of magnetic moment, and Hn « Hdemag (rn) is the demagnetizing field at the point rn. The demagnetizing field is found from a discretization of Poisson equation (2). We use a five-point stencil for the second derivatives with respect to x and y, which gives us a Laplacian of the form A2$ « a-2 where
_ —$x-2,y + 16$x-1,y — 30$x,y + 16$x+1,„ — ^x+2,y
+ 12 ' here (x, y) are coordinates of the moment n on the square lattice. Approximating divergence V ■ m(rn) « a-1GSn in r.h.s of (2) with central finite differences:
-OQ _ Sx+1,y Sx- 1,y . Sx,y+1 Sx,y-1
GSx,y = 2 + 2 ' we obtain the Poisson equation in the discrete form:
= aMs GS.
The last equation is a system of linear algebraic equations, which gives solution upto additive constant. It is convenient to impose additional restriction on value of $ at arbitrary point, e.g. $0,0 = 0. Then, formally, the potential can be
bg ÎH
CD fl H
20 10 0 -10 -20 -30
□
Exchange + Anisotropy
□
40 60 80 100 120 140 160 Radius, nm
fl
œ
Pi
160 140 120 100 80 60 40
>
E(Ms) ■ _
•
• x -
: < »
..................: : • ~
•
V «(Ms ) • —(
1
0 =»
-1 l
i—l
-2
-3 S -
fl
4
0.1
0.2 0.3 Ms, MA/m
0.4
6
FIG . 1. Contributions of different interactions in the total energy (left) and dependencies of the total energy (solid circles) and radius (crosses) (right) for Sk in the Sk lattice on the saturation magnetization Ms Metastable states are computed in 100 x 100 lattice with periodic boundary conditions and lattice constant 4 nm. Film thickness t is set to 0.6 nm, and J = 1.92 • 10-20 J/bond, D = 0.25 J, K = 0.1 J.
/ / ! I
/ III ' / / / //II ■ ' ////II ■■//fill ,,11111
t \ \ V'." ■ - . \ WW • > •. WW \ .
Iu\w -, t w \ \ >. . .. .
■//il I • ■ ; /1 '
m \ \ ■ \\ \ \ \ -, •
I \ WW
' / / m w w w.■
Fig. 2. Single Sk (left) and ASk (right) from the magnetic texture of a square lattice with a distance between Sk (ASk) equal to 400 nm. Simulation is carried out for a 100 x 100 square lattice with a lattice constant a = 4 (nm) for the parameters J = 1.92 • 10-20 J/bond, D = 0.25J, K = 0.1 J, Ms = 0.4MA/m. Only part of moments is shown to make the figure readable. The colors encode the direction of the magnetic moments: green - up, red - down.
found as $ = aMsL 1GS, and the demagnetizing field is given by H = —MsQL 1GS, where Q is gradient expressed in terms of central finite differences:
QJ = ^ Ux+l,y — Ux-l,y Ux,y+1 — Ux,y-1 0
It is worth noting, that the demagnetizing field obtained as a solution of the discretized Poisson equation does not coincide with the common expression for dipole-dipole interaction in lattice model, obtained as discretization of (3):
V V-^ 3(S„ • A)(Sk • A) — Sn • Sk . . rn — rk —— > —-^-75-, where A —
4n nk |rn - rk |3 |r„ - rk
though the energies coincide in the continuous limit. Numeric solution to the discrete version of the Poisson equation can however be found much faster, than direct summation in dipolar interaction or as inversion of matrix L.
Distance to center, nm
Fig. 3. ASk crossection along lattice generator for various values of saturation magnetization Ms. ASk radius grows as Ms increases. The increase is due to a shift of the domain wall on the boundary of ASk rather than scaling of the particle shape without dipole interaction. The core part of the ASk is formed of spins perfectly aligned opposite to FM phase orientation. The domain wall width is also increases for larger magnetization, but slower than radius.
3. Chiral states in demagnetizing field
In addition to the well-studied Sk magnetic structures, in recent years, ASk have attracted attention [9,17]. There are many similarities between these textures, which makes it easy to predict the properties of ASks based on the well-known properties of Sks. The texture of magnetization for Sk and ASk can be matched to each other using the transformation:
mx ^ mx, my ^ — my, mz ^ mz.
Without demagnetizing fields the mapping preserves energy, if DM constants are transformed in the same way as the magnetization. Consequently, a small Sk and an ASk have the same energy surfaces and stability against thermal fluctuations. However, for a large size chiral states, the dipole-dipole interaction plays an important role. The demagnetizing fields for Sk and ASk are significantly different, have different symmetries, leading to distinct shapes of the magnetic structures for a large dipole-dipole interaction.
In this section, we compute (meta)stable Sk and ASk states in the presence of demagnetizing field using the lattice model developed in the previous section. The micromagentic parameters given by (1) are chosen according to [17]. Numerical analysis was carried out for a square lattice of 100 x 100 cells with a lattice constant of 4 (nm). Periodic boundary conditions were applied for both axes, which led to a significant contribution to the energy of dipole-dipole interaction between images of a particle lying on opposite sides of the boundary. Our model is close to modeling the Sk (ASk) lattice with a distance between individual Sk (ASk) 400 (nm), in contrast to the model in [17], where isolated objects were considered. The size of a Sk without dipole interaction is much smaller than the size of a domain, therefore, its energy, radius and shape are practically not affected by an increase in the simulated volume. For a sufficiently large value of Ms, which determines the strength of the dipole interaction, the size of Sk becomes larger, and the contribution of the long-range dipole interaction increases. Therefore, Sk cannot be considered isolated.
Metastable Sk and ASk states were computed using custom optimizer based L-FBGS method. Iterations of L-BFGS method were applied only to magnetic moment directions Sn, while demagnetization field potential $ was computed on each iteration by another solver. In [18] it was shown that L-BFGS method significantly increase convergence speed even for minimum energy path computation. The potential $ was obtained for given moments Sn as solution do discrete Poisson equation by conjugate gradient method, using the potential from previous iteration as an initial approximation. The method demonstrated reasonable convergence rate, providing solution with norm of gradient less than 10-5 J in 1000 iterations for most considered states. Operations on magnetic moments Sn were performed in Cartesian coordinates as stated in the previous section, that differs from commonly used representation of spins in spherical coordinates, stereographic projections [10] or using rotation matrices [19]. Computations in Cartesian coordinates are simple and often faster than usage of other coordinates, see [20]. To take into account constrains on value of magnetic moments, gradients of energy over Sn should be projected to the tangent space of the constrains manifold. Equivalently Lagrange function can be introduced for energy and on Lagrange multiplier for each constrain
Sn = 1, then numerical optimization can be applied to the Lagrange function. To ensure that constrains are satisfied, after each iteration of L-BFGS every direction Sn should be divided by its length.
We computed (meta)stable Sk state for saturation magnetization in range Ms = 0.05 - 0.45 (MA/m). Increasing the magnetization, preserving other exchanges constant, increases the role of demagnetizing field in stabilization of the solitonic state. Increase of Ms decrease total energy of the Sk lattice in non-linear manner, leading to rapid drop in energy about Ms = 0.4 (MA/m) making Sk lattice ground state for larger Ms, as shown in Fig. 1. The radius of Sk grows moderately for Ms below 0.3 (MA/m), after the threshold Sk size grows rapidly being restricted by repulsion between Sk in Sk lattice. Energy (as well as separate contributions) demonstrates perfect linear dependence on the Sk radius for considered values of Ms, see Fig. 1. In contrast to [17] DMI interaction contribute more in absolute value than demagnetizing field, due to partial cancellation of dipole-dipole interaction with other Sk in the Sk lattice.
Sk preserves its circular shape even for large Ms. In contrast ASk changes its shape to a square one for large Ms, see e.g. shape of ASk and Sk in the corresponding lattice for Ms = 0.3 (MA/m) in Fig. 2. The same behaviour was observed in [17]. Shape dependence on saturation magnetization is show in Fig. 3. It can be seen that central part of the ASk is almost precisly in ferromagnetic state, with spins however directed opposite to ferromagnetic phase outside of the Sk. The central "core" part of ASk expands as Ms grows, but domain-wall part of the ASk increases its width moderately. The behaviour is similar to change in anisotropy as expected, since effective anisotropy is commonly used to describe part of the effect of the demagnetizing fields.
4. Conclusion
We have described an approach to computation the demagnetization fields in the framework of the widely used discrete-lattice magnet model. This approach can provide a faster alternative to dipole-dipole interaction computation used e.g. in MuMAX [21] and in OOMMF [18,22]. The approach has been applied to study square lattices ASk and Sk, and the results are in qualitative agreement with the analysis of isolated ASk (Sk) in [17].
Acknowledgements
This work was funded by Russian Science Foundation (Grant 19-42-06302).
References
[1] Fert A., Reyren N., Cros V. Magnetic skyrmions: advances in physics and potential applications. Nature Reviews Materials, 2017, 2, P. 17031.
[2] Finocchio G., Biittner F., Tomasello R., Carpentieri M., Klaui, M. Magnetic skyrmions: from fundamental to applications. Journal of Physics D: Applied Physics, 2016, 49(42), P. 423001.
[3] Nagaosa N., Tokura Y. Topological properties and dynamics of magnetic skyrmions. Nature nanotechnology, 2013, 8(12), P. 899-911.
[4] Bogdanov A. N. and Yablonskii D. A., Thermodynamically stable vortices in magnetically ordered crystals. The mixed state of magnets. Zh. Eksp. Teor. Fiz, 1989, 95(1), P. 178-182.
[5] Bessarab P. F., Uzdin V. M., Jonsson H. Harmonic transition-state theory of thermal spin transitions. Physical Review B, 2012, 85(18), P. 184409.
[6] Uzdin V. M., Potkina M. N., Lobanov I. S., Bessarab P. F., Jonsson H. Energy surface and lifetime of magnetic skyrmions. Journal of Magnetism and Magnetic Materials, 2018, 459, P. 236-240.
[7] Wiesendanger R. Nanoscale magnetic skyrmions in metallic films and multilayers: a new twist for spintronics. Nature Reviews Materials, 2016, 1(7), P. 1-11.
[8] Lobanov I. S., Jonsson H., Uzdin V. M. Mechanism and activation energy of magnetic skyrmion annihilation obtained from minimum energy path calculations. Physical Review B, 2016, 94(17), P. 174418.
[9] Nayak A. K., Kumar V., Ma T., Werner P., Pippel E., Sahoo R., Damay F., Rossler U. K., Felser C., Parkin S. S. P. Magnetic antiskyrmions above room temperature in tetragonal Heusler materials. Nature, 2017, 548(7669), P. 561-566.
[10] Rybakov F. N., Borisov A. B., Bliigel S., Kiselev N. S. New type of stable particlelike states in chiral magnets. Physical review letters, 2015, 115(11), P. 117201.
[11] Vlasov S.M., Uzdin V. M., Leonov A. O. Skyrmion flop transition and congregation of mutually orthogonal skyrmions in cubic helimagnets. Journal of Physics: Condensed Matter, 2020, 32(18), P. 185801.
[12] Buttner F., Lemesh I., Beach G. S. Theory of isolated magnetic skyrmions: From fundamentals to room temperature applications. Scientific reports, 2020, 8(1), P. 1-12.
[13] Caretta L., Mann M., Buttner F., Ueda K., Pfau B., Gunther C. M., Hessing P., Churikova A., Klose C., Schneider M., Engel D. Fast current-driven domain walls and small skyrmions in a compensated ferrimagnet. Nature nanotechnology, 2018, 13(12), P. 1154-1160.
[14] Legrand W., Maccariello D., Ajejas F., Collin S., Vecchiola A., Bouzehouane K., Reyren N., Cros V., Fert A. Room-temperature stabilization of antiferromagnetic skyrmions in synthetic antiferromagnets. Nature materials, 2020, 19(1), P. 34-42.
[15] Potkina M. N., Lobanov I. S., Jonsson H., Uzdin V. M. Skyrmions in antiferromagnets: Thermal stability and the effect of external field and impurities. Journal of Applied Physics, 2020, 127(21), P. 213906.
[16] Lobanov I. S., Potkina M. N., Jonsson H., Uzdin V. M. Truncated minimum energy path method for finding first order saddle points. Nanosys-tems: physics, chemistry, mathematics, 2017, 8(5), P. 586-595.
[17] Camosi L., Rougemaille N., Fruchart O., Vogel J., Rohart S. Micromagnetics of antiskyrmions in ultrathin films. Physical Review B, 2018, 97(13), P. 134404.
[18] Donahue M. J., Porter D. G. OOMMF User's Guide, Version 1.0. Interagency Report NISTIR 6376, National Institute of Standards and Technology, Gaithersburg, MD. 1999.
[19] Ivanov A.V., Dagbartsson D., Tranchida J., Uzdin V. M., Jonsson H. Efficient optimization method for finding minimum energy paths of magnetic transitions Journal of Physics: Condensed Matter, 2020, 32(34), P. 345901.
[20] Lobanov I. S., Uzdin V. M. The lifetime of big size topological chiral magnetic states. Estimation of the pre-exponential factor in the Arrhenius law, 2020. arXiv preprint, arXiv:2008.06754.
[21] Vansteenkiste A., Leliaert J., Dvornik M., Helsen M., Garcia-Sanchez F., and Van Waeyenberge B. The design and verification of MuMax3. AIP Advances, 2014, 97, P. 107133.
[22] Sachdev S. Handbook of magnetism and advanced magnetic materials, 2006.