Electric-field manipulation of magnetization in an insulating dilute ferromagnet through piezoelectromagnetic coupling

Introduction

Repetitive switching of magnetization in ferromagnetic materials is a key process in magnetic recording and information storage. Traditionally, this energy-intensive process is achieved by applying an external magnetic field or spin-polarized currents. An emerging alternative, however, is the use of an electric field to control magnetic anisotropy, which could offer a more energy-efficient solution compared to conventional methods. Different approaches have been proposed in the past to observe the magnetoelectric (ME) effect, in which coupling between magnetism and electricity is present, such as the gate-controlled hole concentration and the Curie temperature in dilute ferromagnetic semiconductors1,2,3,4, voltage induced structural and magnetic phase transition5,6, piezoelectric based modification of magnetism7,8,9 or magnetoelectric effect in multiferroics10,11. The ME effect has been observed in ferromagnetic7, antiferromagnetic12,13, superparamagnetic9, and paramagnetic8 materials.

A first step to accomplish the repetitive switching of magnetization is the demonstration of the control of the key magnetic properties, such as the coercive field, saturation magnetization, or remanent magnetization via electric fields3,14,15. A strong coupling between an electric field and magnetism was shown to exist and theoretically quantified in piezoelectric wurtzite (wz) Ga1−xMnxN in the paramagnetic phase, that is for Mn concentrations x 2.5% and temperatures T ≥ 2 K (ref. 8). Those samples were grown by metalorganic vapour-phase epitaxy (MOVPE)16. Now aided by a carefully elaborated molecular beam epitaxy (MBE) technique, we incorporate substitutionally up to 10% of Mn ions into GaN host17,18. In those semi-insulating layers, Mn3+ ions are coupled by a short-range ferromagnetic (FM) superexchange16,19,20,21, and FM ordering is realized in a percolation-like fashion22,23,24.

Our combined experimental and theoretical studies of the ME effect in the ferromagnetic (Ga,Mn)N address two classes of questions. First, we demonstrate that the sign and magnitude of magnetization generated by an electric field, ME, which we find by superconducting quantum interference device (SQUID) magnetometry, result from a change of uniaxial magnetic anisotropy driven by the inverse piezoelectric effect. Our conclusion is supported by two complementary theoretical approaches: (i) a macrospin model, in which changes of magnetic anisotropy are described in terms of magnetizing or demagnetizing anisotropy fields, depending on the direction of the external magnetic field with respect to the easy axis; (ii) an atomistic Landau-Lifshitz-Gilbert approach for interacting classical spins25,26,27 with input parameters obtained by fitting our magnetization M(H) data. Our results, on the one hand imply that the ME effect can be understood and designed in well characterized magnetic systems. On the other, we find that a short time scale inherent to the LLG approach leads to an overestimation in the theoretical magnitudes of the remanent magnetization and magnetoelectric effect in dilute ferromagnets in the magnetic fields H → 0.

Second, the percolating ferromagnetism means the presence of an infinite FM cluster at T ≤ TC. This cluster encompasses all spins at T = 0 but only about 20% of spins at T approaching TC (ref. 22). Thus, at any non-zero temperature, a large portion of the Mn ions resides in finite clusters of different sizes. Such a magnetic system is expected to exhibit a broad spectrum of relaxation times and, in the hysteretic regime, constitutes a glass far from thermal equilibrium. Thus, a possibility of perturbing ferromagnetic (Ga,Mn)N by an electric field offers a unique opportunity to test theoretical predictions concerning dynamics of driven non-equilibrium glass systems28,29,30.

As we show, our data provide a strong support for persisting theoretical proposals arguing that a non-thermal and non-dissipative drive, bringing a glass beyond the linear response regime, can be described by enlarged effective temperature30. Interestingly, this interpretation appears explaining surprising magnetization behaviour observed in dilute ferromagnetic semiconductor (Bi,Sb,Cr)2Te3 in the regime of the quantum anomalous Hall (QAH) effect31. Tapping of that QHE ferromagnet by gate voltage modulation resulted, according to local SQUID measurements, in magnetization shift towards the direction of the external magnetic field in the hysteretic regime31, as reported here for (Ga,Mn)N. As short-range ferromagnetic superexchange dominates also in topological ferromagnetic insulators32, we expect similar glassy-like relaxation towards thermal equilibrium under a time-dependent drive in both systems.

It is worth noting, however, that our conclusions should not be extended towards experiments, in which magnetization switching in dilute ferromagnetic semiconductors was revealed by resistance measurements33,34, as an important additional ingredient there is Joule heating that may lead to surprising bistable behaviours33. Still another physics may show up in mesoscopic QAH samples, in which resistance switching was found also in magnetic fields far above the coercivity field35,36. In such submicron wires, a highly probable process of backscattering between edge states is controlled by (i) orientation of local magnetization and (ii) spatial arrangement of localized in-gap carriers that are frozen far from thermal equilibrium under the presence of the Hall electric field. Resistance switching dynamics was found to accelerate with the magnetic field in the ferromagnetic QAH case36, as found earlier in mesoscopic antiferromagnetic spin-glasses37. It is to be seen on whether the AC source-drain voltage in QAH nanostructures36 just allowed tracing magnetization dynamics or rather had an effect on the dynamics and effective temperature of the electron glass.

In our work we demonstrate a clear control of magnetic anisotropy of piezoelectric and ferromagnetic (Ga,Mn)N layers by application of an electric field. We observe a reduction of the width of a hysteresis curve from 120 Oe to about 10 Oe as the magnitude of the applied AC voltage VAC increases from 0.2 V to 2 V and a non-reversible switching of magnetization toward the equilibrium direction for magnetic fields close to the coercive field, occurring under a stepped change of the electric field either elongating or shortening the lattice parameter c.

Results and discussion

Samples

Two samples with different thicknesses of (Ga,Mn)N layer d(Ga,Mn)N and similar concentrations of Mn ions x are investigated here, namely, d(Ga,Mn)N = 10 nm, x = 5.6% (Sample A) and d(Ga,Mn)N = 160 nm, x = 6% (Sample B). They were grown along the c axis on a thick GaN buffer followed by an n+-GaN:Si layer (n ~ 3 × 1018 cm−3), which serves as a back-contact contact layer, and characterized according to protocols elaborated previously17,18. The experimental findings19,20 and the theoretical modeling38 agrees that for x < 10% Curie temperature TC of Ga1−xMnxN obeys the power law TCx2.2, which well approximates the percolation exponential formula for small x values. The layers studied here do fall into this dependency exhibiting 7 TC 8 K for x 6%.

In Fig. 1 we compare M(H) data of Sample A acquired at T = 2 K for both in-plane (Hc) and out-of-plane (Hc) configuration for which the demagnetization field was subtracted from the external field value. We see that at that temperature T ≈ TC/4 and the time scale of the experiment, magnetic hysteresis involves only a part of spins, substantiating expectations of the percolation model22,23,24.

Fig. 1: Magnetization, Sample A.
figure 1

a, b Full points: experimental magnetization M of Ga1−xMnxN with x = 5.6% measured at T = 2 K for two perpendicular magnetic field orientations, in-plane (Hc) and out-of-plane (Hc). Shape anisotropy (demagnetization) effects are removed from experimental data. Inset in (b), Close-up on the weak magnetic field regime (H in Oe). Lines in (a) and (b): magnetization calculated for magnetic field sweep direction shown by arrow using atomistic Landau-Lifshitz-Gilbert equations with Langevin dynamics for 9464 classical spins with single-ion magnetocrystalline anisotropy and Heisenberg ferromagnetic interactions between spins (up to 18th neighbours), as described in the theoretical subsection.

Full size image

The data indicate also that magnetocrystalline anisotropy is close to zero for the Mn content and strain in question. In this case, shape anisotropy makes that the easy axis is in-plane in ferromagnetic thin films, as found previously18,20 and here. However, in samples with smaller x values, uniaxial magnetocrystalline anisotropy is clearly non-zero and its sign corresponds to the in-plane easy axis8,39. Those data imply a variation with x of the parameter (xi =c/a-sqrt{8/3}), where c and a are lattice parameters of the wz unit cell. The parameter ξ governs the trigonal deformation and the corresponding uniaxial magnetic anisotropy in magnetically doped wz crystals like Mn-doped GaN8,39,40 or Co-doped ZnO41.

To substantiate those findings, we show in Fig. 2a ξ values obtained by x-ray diffraction for a series of Ga1−xMnxN samples deposited by MOVPE and MBE on a GaN buffer clamping the a value. We also present results of our first-principles calculations with the structure optimization (see “Methods”) by plotting the computed values of ξ and the total energy difference Δε for (Ga,Mn)N with a Mn magnetic moment lying in a plane (along the a-axis) and a Mn magnetic moment aligned with the c-axis as a function of x. As seen in Fig. 2a, the energy difference changes the sign from negative to positive for the Mn concentration of about 2.5%, which correlates but is slightly lower than the computed value corresponding to the sign change of ξ. Those experimental and theoretical data imply that as a function of Mn content, and in the absence of the shape anisotropy, the Mn magnetic moment will have a tendency to reorient from the in-plane to perpendicular direction. Our two samples are close to that reorientation transition with ξ = 0.0029(3) (Sample A) and  − 0.0005(3) (Sample B), in agreement with magnetization data in Fig. 1.

Fig. 2: Wurtzite trigonal deformation.
figure 2

a Left axis: dependence of the parameter (xi =c/a-sqrt{8/3}) on the Mn concentration x measured experimentally at 300 K (stars) essentially identical to evaluated values for 4.2 K (full squares) compared to determined by ab initio computations (open red circles). Right axis: total energy difference between the in plane (along the a-axis) and parallel to the c-axis configurations of the Mn magnetic moment in GaN (open blue squares). Labels indicate ξ values for Samples A, B and MOVPE film from ref. 8. b, c In Ga1−xMnxN, the uniaxial magnetic anisotropy originates primarily from the trigonal deformation of the tetrahedron formed by anion atoms next to Mn. The application of an external electric field E elongates lattice parameter c by the inverse piezoelectric effect and controls uniaxial magnetocrystalline anisotropy.

Full size image

Experimental

As presented in “Methods”, for magnetoelectric investigations Ga1−xMnxN layer resides between the n+-GaN:Si bottom contact and a top HfO2 film covered by a metal gate. Such an arrangement allows for blocking gate currents resulting from a possible conductance via threading dislocations in the (Ga,Mn)N epilayer. In the case of the DC voltage applied to the structure, the electric field on the HfO2 and (Ga,Mn)N films is distributed according to relative resistances of the two layers, so that a rather small electric field is expected in (Ga,Mn)N samples as, presumably, RHfO2R(Ga,Mn)N. In contrast, for resistances and frequencies in question (f = 17 Hz), the AC voltage is distributed according to relative capacitances, meaning that the application of the AC voltage of VAC = 1 V to the investigated structures (see “Methods”) corresponds to the AC electric field EAC of about 0.19 MV/cm and 0.04 MV/cm to the (Ga,Mn)N layer of Sample A and B, respectively. Given the structure parameters and wiring layout, we do not expect temperature changes by capacitor charging currents at experimental temperature of 2 K.

The measurement system, depicted in Methods, utilizes a SQUID magnetometer and phase sensitive detection. An AC voltage VAC is applied between the metal gate and the GaN:Si contact layer. The resulting AC magnetization signal of the sample is picked up by SQUID coils and detected by a lock-in amplifier. Alternatively, we look for magnetization alterations generated by stepped DC voltage changes. As shown in panels b and c of Fig. 2, in our piezoelectric structure with the clamped value of a by the thick buffer layer, an electric field E alters only the lattice parameter c, according to Δc/c = d33E, where d33 = 2.8 pm/V, as obtained for an epitaxial GaN film with the same orientation of the c vector as in our samples, and also with a clamped value of a by a thick GaN buffer layer42.

Experimental results and macrospin model

As mentioned in the introduction, our data reveal the existence of two distinct magnetoelectric phenomena. First of them, discussed in this subsection, corresponds to the regime, where the magnetoelectric response ME is linear in the electric field E. Given a rather unique combination of ferromagnetism and piezoelectricity offered by single crystals of wurtzite (Ga,Mn)N, the data in the linear range makes possible to test quantitatively our understanding of magnetoelectricity in dilute ferromagnets. The second phenomenon, presented in the last subsection of Results, is specified by strong nonlinearities and irreversibilities of ME. As we show, the data in that regime give access to the physics of driven glassy systems far from thermal equilibrium.

We recall that in our experiments the electric field is applied along the c axis of the wurtzite structure, and we measure magnetization changes in the direction of the applied magnetic field H. We have found that in high magnetic fields, ME(H) is linear in E up to the highest electric field supported by our devices. In contrast, in the hysteresis regime, nonlinearties are already apparent in the electric fields lower by a factor of ten. Thus, being interested here in the linear regime, we discuss ME values normalized to VAC = 0.5 V (Sample A). Such a low value of VAC has been employed in low magnetic fields H ≤ 1 kOe but VAC has been risen up to 5 V for H ≥ 10 kOe.

According to the data displayed in Figs. 3 and 4 for Sample A at 2 K, the electric field decreases magnetization for the magnetic field Hc, but increases M for Hc. Interestingly, comparing ME(H) in Fig. 3 to M(H) in Fig. 1, we see that ME is roughly proportional to ∂M/∂H in high magnetic fields but according to Fig. 4, ME becomes proportional to M in weak magnetic fields. The presented data implies that the converse magneto-electric coupling coefficient, obtained from the relation αC = 4πME/E, attains 7.3 Oe cm/MV for wz-Ga0.94Mn0.06N at 2 K and in 0.5 kOe, as E(Ga,Mn)N = 0.095 MV/cm for 0.5 V applied to Sample A. The highest value of ME we have detected is 0.23 emu/cm3 at 5 V, to be compared to the magnetization saturation value Msat = 80 emu/cm3.

Fig. 3: Magnetoelectric effect, Sample A.
figure 3

Areal density of the magnetic moment ME/A normalized to the voltage VAC = 0.5 V in Ga1−xMnxN layer with x = 5.6% at T = 2 K, i.e., in the ferromagnetic phase (blue solid points). The data are plotted as a function of the square root of the magnetic field to show simultaneously weak and strong field regions. Measurements are taken for the in-plane (Hc) and out-of-plane (Hc) configurations. a, c Open red points: theoretical results of magnetoelectric signal ME obtained from the macroscopic spin model (Eq. (1)) with one fitting parameter b = 0.04. b, d Open green stars: the magnetoelectric signal computed by using the atomistic Landau-Lifshitz-Gilbert (LLG) equations with Langevin dynamics.

Full size image
Fig. 4: Magnetoelectric effect in the hysteretic regime.
figure 4

a Areal density of the in-plane magnetic moment amplitude ME/A of Ga1−xMnxN layer with x = 5.6% induced by AC voltage VAC = 1V at T = 2 K. The measurements are taken in the in-plane (Hc) and out-of-plane (Hc) configuration (full circles and triangles, respectively). These results are compared to directly measured magnetizations  − M(H) and M(H) (grey lines, relative units), measured for Sample B before processing. The virgin magnetization data obtained after cooling in H ≈ 0.15 Oe are shown by empty symbols. b Areal density of the in-plane magnetoelectric signal ME/A at H = 0 as a function of temperature T (blue circles), compared also to experimental value M(T) (grey line, relative units). c Magnetoelectric response to VAC = 1 V measured in the presence of the DC bias VDC = ±28 V applied in the magnetic field beyond the hysteretic regime.

Full size image

We describe ME(H) in this easy-plane piezoelectric ferromagnet in terms of an effective anisotropy magnetic field HE affected by the electric field. If HEH,

$${M}_{E}({H}_{j})={M}_{j}({H}_{j}+{H}_{E})-{M}_{j}({H}_{j})approx {chi }_{jj}({H}_{j}){H}_{E},$$
(1)

where j refers to the two explored configurations of the magnetic field; χjj = ∂Mj/∂Hj and, since uniaxial anisotropy is proportional to ξ(E),

$${H}_{E}({H}_{j})=pm b{d}_{33}E{M}_{j}({H}_{j})c/(axi ),$$
(2)

where  ± corresponds to Hc and Hc, respectively, and b is a dimensionless fitting parameter. The above equations define the form of the magnetoelectric coefficient43αC = ME/E for the case under consideration.

We see in Fig. 3 that our macrospin model with experimental values of χjj(Hj) and b = 0.04 describes quite satisfactorily the field dependence for the two orientations of the magnetic field down to 4 kOe at VAC = 0.5 V.

While in high magnetic fields ME(Hj) χjj(H), in weak magnetic fields, ME(Hj) becomes proportional to Mj(Hj), as shown in Fig. 4. This observation supports experimentally the expectation that αC is non-zero only under time-reversal symmetry breaking43, i.e., αC is odd in Hj in the paramagnetic case, as observed previously8, or odd in Mj at Hj → 0, as found here for the ferromagnetic ordering. This means that for sufficiently small values of Mj, ME should be linear in Mj, the relation fulfilled by our data displayed in Fig. 4a, b as a function of the magnetic field at 2 K and as a function of temperature at H = 0, respectively. Accordingly, as shown in Fig. 4b, the magneto-electric effect vanishes near TC, about 7.5 K for Sample A. Quantitatively, in weak magnetic fields, ME(Hj) = ±βM(Hj), where β = 2 × 10−3 for Sample A and VAC = 1V.

As presented in Fig. 4c, we have also checked that the application of the DC bias of either polarization has no effect on ME(H) generated by VAC. This finding is consistent with the fact that the DC electric field is non-zero only across the HfO2 film.

As far as the origin and magnitude of the parameters b and β are concerned, we note that, in general, crystalline magnetic anisotropy, which we perturb by the electric field via the inverse piezoelectric effect, originates from a single-ion contribution and anisotropy of spin-spin coupling. In magnetic compounds in question, the spin-orbit interaction within anions transmitting exchange between localized spins plays a dominant role44,45, so that this anisotropy should be relatively small for the light nitrogen ions. This fact suggests that single-ion anisotropy, rather sizable for Mn3+ ions with the orbital momentum L = 2, dominates in (Ga,Mn)N. The signs of b and β are consistent with this expectation. However, though we know rather precisely the crystal field parameters of Mn3+ ions in GaN and their dependence on ξ (ref. 8), there are persistent difficulties in relating that information with parameters of the macrospin model46, particularly in dilute magnetic systems47 and in the percolating case. This fact makes a theoretical determination of b and β difficult within the macrospin model. We have, therefore, decided to adapt a stochastic Landau-Lifshitz-Gilbert (LLG) atomistic spin model25,26,27 for the case under consideration.

Landau-Lifshitz-Gilbert atomistic spin model

The way that approach26,27 is adopted to ferromagnetic (Ga,Mn)N is detailed in Methods. Our atomistic LLG model was also employed to describe magnetization of (Ga,Mn)N in the paramagnetic regime, where magnetic properties are determined by single Mn ions, nearest neighbour Mn pairs, and triangles, so that the outcome of the LLG model which treats the Mn spin classically, could be compared to the results of computations describing the spins quantummechanically48.

The magnetic energy of Mn spins is described by a spin Hamiltonian ({{mathcal{H}}}) taking into account the Zeeman energy, Heisenberg spin-spin interactions between Mn spins located up to the 18th nearest neighbour positions, and the magnetocrystalline anisotropy energy, ({{mathcal{H}}}={{{mathcal{H}}}}_{{{rm{Z}}}}+{{{mathcal{H}}}}_{{{rm{Exch}}}}+{{{mathcal{H}}}}_{{{rm{Anizo}}}}), where ({{{mathcal{H}}}}_{{{rm{Anizo}}}}) includes the Jahn-Teller contribution and the trigonal anisotropy term taken as,

$${{{mathcal{H}}}}_{{{rm{Trig}}}}=-frac{1}{2}{K}_{{{rm{Trig}}}}{sum}_{i}frac{1}{2}left[{S}_{iz}^{2}-left({S}_{ix}^{2}+{S}_{iy}^{2}right)right].$$
(3)

Here, KTrig ξ is the anisotropy parameter describing the trigonal deformation that can be controlled by the electric field E, and Si is a unit vector describing the local magnetic moment direction at the site i. The explicit forms of ({{{mathcal{H}}}}_{{{rm{JT}}}}) and ({{{mathcal{H}}}}_{{{rm{Exch}}}}) together with the values of relevant exchange integrals, are given in Methods.

The results on thermal equilibrium values of macroscopic magnetizations Mj(Hj) without and with an electric field E and for two directions j of the external magnetic field are obtained for a simulation box 21 × 21 × 18 nm3 containing 169,000 cation sites, 9464 (5.6%) of which encompass randomly distributed spin vectors. The site number is limited only by the available hardware employed to trace the time evolution towards thermal equilibrium values Mj(Hj) over an acceptable time frame. The simulations start in 50 kOe, in which Mn spins are aligned along the magnetic field, and next the field is decreased in steps specified in “Methods”.

Considering the lack of a clear correspondence between parameters of the quantum and classical spin Hamiltonians, our strategy involves two steps. First, by fitting the LLG model to the experimental data on Mj vs. Hj at E = 0 and for the two orientations of the external magnetic field, as presented in Fig. 1, we have evaluated KJT = 0.85 meV and KTrig = 0.1 meV appearing in system Hamiltonian. A positive sign of KTrig is in agreement with the data on ξ obtained from x-ray diffraction (XRD) measurements, shown in Fig. 2a. However, it should be noted that the simulations point to larger remanence than observed experimentally. Moreover, even rescaling Curie temperature obtained from the numerical modeling for classical spins by a factor (S + 1)/S, we obtain ({T}_{{{rm{C}}}}^{{{rm{sim}}}}cong 4.5) K, which differs significantly from the experimental value TC 7.5 K.

Second, we have performed computations for the determined value of KJT and for KTrig in the absence and presence of the electric field,

$${K}_{{{rm{Trig}}}}(E)={K}_{{{rm{Trig}}}}(0)[1+{d}_{33}Ec/(axi )],$$
(4)

where ξ = 0.0029(3) for Sample A according to Fig. 2a. Here, in order to reduce uncertainties resulting from statistical fluctuations of computational results, we considered E corresponding to VAC = 5 V, and reduce the numerical values of ME(Hj) by a factor of ten to compare to experimental results depicted in Fig. 3 for VAC = 0.5 V. As seen in Fig. 3, the LLG theory is successful only qualitatively, i.e., predicts properly the sign and shape of ME(Hj). The obtained results bring us to the conclusion that the LLG approach, when applied to dilute ferromagnets, predicts too large magnitudes of remanence in the case of both M(Hj) and ME(Hj). We assign this outcome to time scale of our experiment, which is sufficiently long to allow demagnetization of most Mn spins outside the percolation cluster. Such a demagnetization process does not occur in the necessary much shorter time frame of the LLG simulations, so that most of Mn spins contribute to the theoretical remanence, as shown by solid lines in Fig. 1.

Nonlinearities and irreversibilities

We now discuss surprising switching effects found in the hysteretic regime. The magnetoelectric signal ME(H) in this range and at 2 K is plotted in Fig. 5a as a function of the in-plane magnetic field Hc for three different AC modulation voltages VAC. As seen, the width of the coercive field Hc gets sizably reduced (from 120 Oe to about 15 Oe) by increasing VAC to 5 V. The effect is better visualized in Fig. 5b, c, where we present changes of the coercive field with the applied voltage VAC to the whole structure and with the corresponding electric field EAC across the (Ga,Mn)N layer of the thickness d = 10 and 160 nm for Sample A and B, respectively. We see that the ME appears highly nonlinear in EAC in the hysteretic region. According to Fig. 5c, this nonlinearity is stronger in Sample B in which, judging by the magnitude of Hc, ferromagnetism is softer than in Sample A. Altogether, these results imply that taping of the system by sufficiently high VAC brings the system towards thermal equilibrium.

Fig. 5: Effect of the AC voltage on coercivity.
figure 5

a Areal density of the in-plane magnetic moment amplitude ME/A generated by three different values of the AC voltage VAC as a function of the in-plane magnetic field. b, c Corresponding coercivity field as a function of VAC and the electric field EAC across the (Ga,Mn)N layers of thicknesses 10 and 160 nm for Samples A (green circles) and B (red squares), respectively.

Full size image

While a time-independent DC bias has no effect in our structures (see Fig. 4c), we find that step changes of the DC voltage, if applied in the non-equilibrium range of the hysteresis, moves magnetization toward the direction of the external magnetic field. Our experiment begins by sweeping the external magnetic field to the non-equilibrium range of the hysteresis, and then we change VDC in a step-like manner by ΔVDC at constant H monitoring ME at selected voltage values VDC = 0, 3,…, 18 V or, for the negative sign of ΔVDC, at VDC = 0, −3, …, −18 V. There is at least 1 sec delay between consecutive voltage changes.

Results of such investigations are collected in Fig. 6 which shows that: (i) Magnetization changes generated by voltage alterations occur only in the non-equilibrium regime and those changes are irreversible, i.e. magnetization does not change under removal of voltage. (ii) According to results depicted in Fig. 6a, b, outside the metastable regime, the observed drift in ME magnitudes for consecutive voltage changes is of the order of time-dependent fluctuations in the ME values, shown in Fig. 6c. (iii) Voltage-induced shifts in ME(VDC) values are independent of the DC voltage polarization, and a sufficiently high magnitude of VDC reverses the magnetization direction (Fig. 6a, b). (iv) As shown in Fig. 6d, the consecutive values of ME(VDC) are also virtually independent of the number of steps leading to a given VDC, at least in the range of step heights 1 ≤ ΔVDC ≤ 15 meV.

Fig. 6: Magnetization switching by DC voltage changes.
figure 6

a, b At VDC = 0 V, the magnetic field is swept from −450 to +50 Oe or to −50 Oe (blue circles, full sweep is shown). The VDC is then changed in step-like manner by ΔVDC = 10 mV or −10 meV, and the resulting magnetoelectric signal ME is measured at successive points VDC = 0, 3,…, 18 V (orange points or magenta points in a) or VDC = 0, −3, …, −18 V (orange points or green points in b). Irreversible switching of magnetization is observed independently of the DC voltage polarization if for H < 0, M > 0 (i.e., ME < 0 for Hc configuration), meaning that the system is initially out of thermal equilibrium. c Magnetoelectric response ME as a function of time for ΔVDC = 0. d, Magnetoelectric response ME as a function of VDC for various magnitudes of DC voltage steps ΔVDC.

Full size image

It is remarkable that irreversible shifts of magnetization toward the thermal equilibrium value by voltage alterations, as we report here, were also detected at 0.25 K by local SQUID magnetometry in gated structures of another dilute ferromagnetic semiconductor (Bi,Sb,Cr)2Te3 (ref. 31). As in our case, the effect was independent of gate voltage polarization and the magnetization alteration scaled with the magnitudes of the gate voltage excursion. Also, time-independent gate voltage had little effect on magnetization and hystereses. Furthermore, spatially-resolved magnetization studies pointed to the presence of superparamagnetic regions with noticeable magnetization relaxation over the laboratory time scale.

To explain the origin of that universal behaviour, we first recall that despite magnetic dilution, a well-defined magnetic domain structure was found in (Ga,Mn)As, the fact explained quantitatively in terms of a long-range character of spin-spin interactions in that type of ferromagnets49. In order to address this question in (Ga,Mn)N, we show in Fig. 4a experimental data on virgin magnetization obtained after cooling the sample down to 2 K under near-zero-field conditions (H ≈ 0.15 Oe). The measured form of m(H) corresponds neither to the domain-wall motion nor to the domain-wall pinning in a typical multidomain case for a ferromagnet. In the former case one would expect a rapid increase of m, nearly along the y-axis, in the latter the initial m(H) would be horizontal along the x-axis. Actually, the observed behaviour is consistent with or even supports the view that in the case of short-range interactions, driven here by superexchange, the virgin magnetic properties are dominated more by superparamagnetic regions present in the percolation picture than by a regular domain structure which is formed by dipole interactions in dilute ferromagnets with long-rage exchange couplings49.

More specifically, we expect in both (Bi,Sb,Cr)2Te3 and (Ga,Mn)N, a co-existence, at TC > T > 0, of an infinite ferromagnetic percolating clusters with superparamagnetic regions, characterized by a broad spectrum of blocking temperatures and relaxation times. With no doubts, the presence of long-range dipole interactions, responsible for the domain formation in materials like (Ga,Mn)As, contributes to the complexity of spin dynamics. Therefore, (Bi,Sb,Cr)2Te3 and (Ga,Mn)N in the metastable regime can be regarded as a glass that is far from thermal equilibrium. A question arises: how do such media respond to a time-dependent perturbation? It has been theoretically demonstrated that in driven glassy or granular systems, the degrees of freedom governing slow relaxation thermalize to an effective temperature, whose deviation from bath temperature scales with the drive strength28,29,30. Considering a steady decay of the coercive force with temperature, we see that the data on dilute ferromagnetic semiconductors with a short range of spin-spin interactions provide the experimental corroboration of theoretical predictions concerning tapped glassy systems. However, neither in (Ga,Mn)N, nor in (Bi,Sb,Cr)2Te3 a full magnetization reversal has been attained by time-dependent voltage manipulations, at least in the explored time-domain much longer than the inverse precession frequency. This fact is in accord with our scenario that implies the presence of clusters in which a ferromagnetic spin arrangement is robust to weak perturbations.

To the best of our knowledge, the experiments have revealed novel aspects of driven dynamics, such as the fact that the trajectory toward the thermal equilibrium depends only on the accumulated magnitude of the drive changes but not on the number of those changes. It is worth noting that the standard Langevin dynamics does not capture the effects of time-dependent perturbations, so that the nonlinear and irreversible effects do not show up in our LLG simulations.

Conclusions

In our work, we have exploited the recent progress in the MBE growth of high-quality quantum structures incorporating ferromagnetic (Ga,Mn)N to examine much needed spintronic functionality, i.e., the manipulation of magnetization by an electric field. Our results demonstrate that (Ga,Mn)N with a sufficiently large Mn concentration is one of rather rare piezoelectric ferromagnetic homogeneous compounds, in which the magnetoelectric effect can be detected and, moreover, examined theoretically in a quantitative manner. In another perspective, the magnetization manipulation enlarges a pallet of wide-ranging nitride functionalities50 that have made GaN and related systems the second most important semiconductor family next to Si. In particular, GaN:Mn serves already as a semiinsulating substrate51 in device epitaxy. The spectrum of (Ga,Mn)N functionalities has recently been enriched by the demonstration of efficient generation and detection of spin currents in Pt/(Ga,Mn)N structures up to 50 K (ref. 52) using spin Hall magnetoresistance (SMR) effect53. The same SMR could be then employed to probe the magnetic state of the insulating (Ga,Mn)N in the electrical way, paving the way towards more practical applications.

Our experimental results have revealed the presence of two distinct regimes showing a different pool of phenomena. The first of them concerns a wide range of magnetic and electric fields, in which the magnetoelectric response ME is linear in the AC electric field E and reversible, i.e., the magnetization magnitude returns to its initial value after switching off the electric field. We have developed for that case (i) the macrospin approach and (ii) the atomistic theory with Landau-Lifshitz-Gilbert equations containing the external magnetic field, molecular fields generated by neighbouring Mn spins, and Gaussian magnetic noise—the Langevin dynamics. The macrospin model, containing one fitting parameter, describes ME for both magnetic field configurations, Hc and Hc over a wide field range, including the hysteretic regime. In the case of the LLG theory, isotropic ferromagnetic Mn-Mn exchange integrals come from the theory which describes TC(x), whereas parameters of the single-ion anisotropy, from fitting M(H) for the two configurations. A general agreement between experimental and theoretical results on ME(H) demonstrates that the magnetoelectric response results primarily from the influence of the inverse piezoelectric effect on the wurtzite c lattice parameter, which affects the trigonal uniaxial single-ion magnetocrystalline anisotropy. At the same time, we have identified a limited predictive power of the LLG simulations concerning the magnitudes of remanent magnetization and magnetoelectric effect in dilute ferromagnets at H → 0.

The other magnetoelectric phenomena discussed in our work appear in the metastable range of the low-field hystereses. Our results, and also those found for (Bi,Sb,Cr)2Te331, let us to suggest that dilute ferromagnetic semiconductors with short-range spin-spin couplings show, in the metastable regime, characteristics similar to glasses far from thermal equilibrium. One of pertinent questions is the response of such systems to time-dependent perturbations specified by a time scale long compared to both the inverse Larmor frequency and the transverse spin relaxation time T2, governed by non-scalar spin-spin interactions, but comparable to T1 of superparamagnetic clusters incorporating most of localized spins at T > 0. To explain characteristics of irreversible magnetoelectric responses in (Ga,Mn)N and (Bi,Sb,Cr)2Te3, we refer to extensive numerical simulations of glasses28,29,30, indicating that time-dependent perturbations of such systems lead to fast thermalization to an enhanced temperature magnitude facilitating a unidirectional drift of the magnetization value and direction toward the thermal equilibrium state. Conversely, magnetoelectric data for relatively simple and easy controllable systems, such as (Ga,Mn)N and (Bi,Sb,Cr)2Te3, can motivate further theoretical studies of non-equilibrium physics in disordered systems.

Methods

Material

The (Ga,Mn)N layers investigated here have been grown by plasma assisted molecular beam epitaxy (MBE) method in a Scienta-Omicron Pro-100 MBE at 620 C under near stoichiometric conditions following closely the previously elaborated protocol18. In order to reduce the number of threading dislocations and simultaneously provide the back electrical contact to the final capacitive device we use commercially available conductive low threading dislocation density (LTDD) n-type (0001) GaN:Si templates deposited on a c-plane sapphire, which contain, generally termed, an SiNx mask. During the growth the surface quality of the layers has been investigated in situ by reflection high-energy electron diffraction (RHEED). The post growth [11(bar{2})0] RHEED pattern of Sample B collected after the growth (at a temperature below 200 oC) is shown in Fig. 7. The streaky and sharp RHEED pattern with clearly visible Kikuchi lines indicates that the (Ga,Mn)N surface is smooth and relatively free of oxides or other contaminants.

Fig. 7: Reflection high-energy electron diffraction (RHEED).
figure 7

RHEED pattern along [11(bar{2})0] azimuth of Sample B observed after the growth and cooling the sample down to 200 oC.

Full size image

Two samples have been prepared for this study. Sample A of a low thickness of (Ga,Mn)N layer, d(Ga,Mn)N = 10 nm and a thicker Sample B with d(Ga,Mn)N = 160 nm. The precise values of d(Ga,Mn)N are established upon high resolution x-ray characterization described below. Both are having a similar Mn concentration x, about 6%, established by direct magnetic measurements.

Structural characterization

The crystal structure, epitaxial relationship, in- and out-of-plane lattice parameters, and the strain condition of the samples are established by high resolution x-ray diffraction (HRXRD) using a Philips X’Pert MRD diffractometer operating on wavelength 1.5406 Å of CuKα1 line. Its incident beam is formed by an x-ray mirror and a monochromator [four-bounce (220) Ge asymmetric cut]. Figure 8 collects the most relevant findings. Panel a demonstrates a 2θ/ω scan measured for Sample B over a wide range of detector angle changes, indicating that within a dynamic range of intensity exceeding four orders of magnitude no reflections specific to inclusions of other phases are seen. The reciprocal space map of the (bar{1}bar{1})24 reflection is given in panel b. Two main features can be recognized there: the high intensity peak located at Qz = 4.827 Å−1 – related to the GaN template, and a much weaker one at smaller reciprocal lattice units, at Qz = 4.817 Å−1, from the thinner (Ga,Mn)N layer. The map indicates a nearly perfect pseudomorphic growth, meaning that the (Ga,Mn)N layer adopts the same in-plane lattice parameter as that of the substrate (the GaN LTDD template).

Fig. 8: High-resolution (HR) x-ray characterization of investigated samples.
figure 8

a Diffraction pattern of Sample B in 2θ/ω mode. Three GaN reflections denoted with a star are a result of a multiple reflection and are seen here due to a specific orientation of the sample. The reflections indexed on the grey background are from the sapphire substrate. b The reciprocal space map of the asymmetric (bar{1}bar{1})24 reflection of Sample B. The high-intensity node is related to the GaN template layer. The low-intensity node reflects a practically fully strained (Ga,Mn)N layer, as can be inferred from the location of this node very close to the vertical line, which, together with the inclined dashed line define the triangle of relaxation. c, d HR diffraction pattern for the 0002 Bragg reflection for Sample A and B, respectively (the dark blue lines). The orange lines marks results of simulations upon which the values of the layers thicknesses are established (given in the legends).

Full size image

HRXRD patterns for the 0002 Bragg reflection are marked as the thicker blue line in Fig. 8c, d for Samples A and B, respectively. Clearly resolved x-ray interference fringes imply a high structural perfection of the layers and good quality of the interfaces. The peaks corresponding to the (Ga, Mn)N epitaxial layers are shifted to smaller angles with respect to that of the GaN buffer, a result of a larger perpendicular lattice parameter of (Ga,Mn)N. The orange lines mark results of simulations performed in the frame of the dynamical theory of x-ray diffraction using commercially available PANalytical EPITAXY software. In these computations we put the elastic stiffness constants in (Ga,Mn)N as those in bulk GaN and c = 5.1850(5) Å, and a = 3.1885(5) Å, as established in ref. 54. The best fit to the experimental data is obtained for d(Ga,Mn)N = 9.5(1) nm for Sample A and d(Ga,Mn)N = 161(1) nm for Sample B, from which simplified values of 10 and 160 nm values are adopted to characterize both layers. We can parenthetically add that allowing for 5% of lattice relaxation (to account for a small horizontal shift of the (Ga,Mn)N node seen in Fig. 8b) does not call for a greater modification of d(Ga,Mn)N than the modeling uncertainty specified above.

The XRD measurements are also used to establish the values of the trigonal distortion parameter (xi =c/a-sqrt{8/3}). Generally, the value of c is obtain from the symmetrical 0002 reflection, a parameter is calculated based on data from both 0002 one and asymmetrical −1–124 reflections. Bragg angle of (bar{1}bar{1})24 is usually determined from a relevant reciprocal space map, as the one presented in Fig. 8b. The bare results for ξ obtained from XRD in Samples A and B are: ξA = 0.0029(3) and ξB = −0.0005(3), respectively. However, these values quantify the distortion at room temperature, whereas PEME measurements are performed at 2 K. So one can expect an additional deformation in response to the specific thermal contractions of all the constituent layers embedded within the structure. In order to assess the magnitude of the temperature induced modification of the trigonal distortion of the (Ga,Mn)N unit cell upon cooling we transform ξ(300 K) into ξ(4 K) by taking into account55,56 relevant elastic constants57 and thermal expansion coefficients of a and c in sapphire58 and GaN, ref. 59. We find that for this particular material combination [the layer stack is presented in a later section (Fig. 9)]. The magnitude of the expected change of ξ upon cooling to liquid helium temperatures is marginally small. The calculated change does not exceed +0.0001, which is only a third of the typical uncertainty of establishing ξ from HRXRD. We report all these values in Table 1.

Fig. 9: Sample A right after the lithographic processing.
figure 9

Gold pads W at the corners are the contacts to the back gate layer of n+-GaN:Si, the topmost layer of the substrates on which the 10 nm-thick (Ga,Mn)N layers have been deposited. Small contact pads C are deposited on 60 nm HfO2-dielectric oxide deposited by atomic layer deposition onto the surface of (Ga,Mn)N to prevent electrically shorts via the threading dislocations. Each pad forms a parallel plate capacitor with the back gate conductive layer of GaN:Si. Upon further electrical testing only the non-leaking small capacitor (up to 12) are connected in parallel to form the final metal-oxide-semiconductor capacitor structure with capacitance of up to 1.2 nF.

Full size image
Table 1 Parameters of sample A and B
Full size table

Importantly for this project, the atomic force microscopy AFM investigations yield a low surface roughness. As reported in ref. 60 for the same Sample B and in ref. 18 for a range of sistering samples, the root mean square of the roughness of the surface does not exceed 2 nm on a 5 × 5 μm2 field of view. The main reason of the surface roughening is the formation of spiral hillocks known to form at the locations where the dislocations piercing the (Ga, Mn)N layer reach the free surface, otherwise clear monoatomic steps are observed. The hillocks’ tops are these locations on the surface of GaN, which provide short-circuiting electrical passes for the vertical current across the layer61, so they have to be blocked for achieving effective gating of the material. It was previously estimated60 that the surface density of hillocks in our samples, grown on LTTD substrates, does not exceed 108 cm−2, which represents an improvement exceeding two orders of magnitude compared to typical growth on plain GaN-templated sapphire substrates.

The growth has been performed on quarters of 2” wafers. In order to avoid a high lateral gradient of the Mn concentration, and so of the other magnetic characteristics caused by the inhomogeneous substrate temperature during the MBE growth18, for the current study the specimens are cut from the center of the quarter-wafers substrates, where the homogeneity of Mn concentration and related magnetic characteristics is the best18,62.

Magnetic characterization

All the magnetic measurements are performed using a Quantum Design MPMS XL superconducting quantum interference device (SQUID) magnetometer following a procedure elaborated earlier8,63,64 and adjusted for measurements of ultra-thin (Ga,Mn)N films65 deposited on sapphire substrates66. As shown previously for MBE-grown epilaers17,18,20, and checked for the present samples, our high-resolution magnetic data do not reveal the presence of any high-TC ferromagnetic precipitation.

SQUID magnetoelectric measurements

For the gate-dependent measurements Ga1−xMnxN layers need to be embedded in metal-insulator-semiconductor (MIS) structures. As the back gate electrode we use the conductive GaN:Si top buffer layer. To get an access to it we open by reactive ion etching (RIE) at least two 1 × 1 mm2 windows in the corners of our 4 × 5 mm2 specimens (marked “W” in Fig. 9). After placing photoresistive masking over these windows the whole surface is covered with 60 nm gate oxide (HfO2) grown by atomic layer deposition. Next, after re-opening of the windows by a lift-off process, we lithographically pattern between 20 to 60 of separated gold pads over the surface covered by the gate dielectric. At the same process we deposit Au in the window areas to make the electrical contacts the GaN:Si layer—the back gate electrode(s). The gold pads are of an area about 0.2 × 0.2 mm2 each and they form top capacitor plates to the back gate GaN:Si conductive layer. Figure 9 shows one of the investigated samples right after completing the patterning process. The capacitors are then electrically tested individually at room temperature. Only capacitors characterized by a high breakdown voltage VG ≥ 20 V are connected in parallel using silver-filled electrically conductive epoxy. A resulting total capacitance of the investigated structures varies between 0.5 to 1.2 nF. A simplified general layout of the investigated structures is given in the central part of Fig. 10.

Fig. 10: Schematic diagram of the experimental setup to perform magnetoelectric measurements in SQUID magnetometer.
figure 10

The sample is located at the center of the magnetometer’s pick up coils and the electrical gate induced changes of magnetic moment of (Ga,Mn)N are detected by phase sensitive fashion using a Lock-in amplifier.

Full size image

The completed capacitor-like structure is mounted on a modified sample holder provided by Quantum Design for performing electrical measurements in their MPMS magnetometers. In particular the modifications allows for both in plane and perpendicular orientation of the sample for electrical measurements. We perform the magnetoelectric measurements following the protocol described in ref. 8 with a number of small practical improvements.

The schematic diagram of the experimental set up is given in Fig. 10. Contrary to the standard method of SQUID magnetometry in which the sample is transported up and down across the SQUID pick-up coils, in our approach the sample is kept stationary at the central location of the SQUID’s 2nd order gradiometer. Having the sample stationary in SQUID a time-dependent response is provided by the application of an AC modulation voltage to the gate (capacitor) over the (Ga,Mn)N layer. The AC voltage VAC of an acoustic frequency f (3 ≤ f ≤ 17 Hz) via the inverse piezoelectric effect modulates the magnetic properties of the (Ga,Mn)N layer, and these changes are picked up by SQUID coils. After amplification by the SQUID electronics, these signals are fed to the lock-in amplifier, where the amplitude of the changes is established in the phase sensitive fashion. We apply both the DC and AC external voltages: VDC to change the magnetic anisotropy and, usually much smaller, VAC to detect the resulting ME. The method allows for detection changes down to 10−10 emu, however on an expense of progressively longer acquisition times. At f = 17 Hz of AC modulation it takes about 30 mins to complete each experimental point, including 5 min of initial delay to stabilize the field flux in the bore of the magnetometer’s superconducting magnet. Using about 0.1 s delay time during the data acquisition we collect about 12,000 single readouts. By plotting a histogram of the sample data, we obtain the magnetoelectric response, ME, as the mean value of a fitted normal distribution. It should be emphasized that this technique does not sense time-independent magnetic properties of the substrate or other elements of the whole set-up. The in-phase AC magnetic response is generated only in (Ga,Mn)N located between the plates of the selected individual capacitors.

Density functional theory calculations

The first-principles spin-polarized calculations are performed using a plane-wave basis set and the Perdew-Burke-Ernzerhof exchange-correlation functional, as implemented in the Quantum-ESPRESSO package67. The non-collinear magnetism is taken into account by performing calculations that include spin-orbit coupling and fully relativistic norm-conserving pseudopotentials. The plane-wave cutoff is set to be 50 Ry. Five different unit cells are used, namely 2 × 2 × 1, 2 × 2 × 2, 2 × 2 × 3, 3 × 3 × 2, and (3times 2sqrt{3}times 2) with 16, 32, 48, 72, and 96 atoms, respectively. In each case one Ga ion is replaced with Mn and the positions of all ions are fully relaxed. For each unit cell, the lattice parameter a is fixed to the theoretical value for GaN (i.e., to the substrate), and the c lattice parameter is fully optimized. The structural optimization is done using ultrasoft pseudopotentials and a 30 Ry kinetic energy cut-off (240 Ry for the electronic density). For the Brillouin zone sampling 8 × 8 × 8 and 3 × 3 × 3 k-point meshes are used for the smallest and largest unit cells, respectively. The calculated equilibrium values of lattice parameters for GaN are a = 0.3217 nm and c = 0.5241 nm, and agrees well with earlier theoretical reports68.

Atomistic spin model simulations

The magnetic energy of localized spin vectors Si is described by the Hamiltonian,

$${{mathcal{H}}}={{{mathcal{H}}}}_{{{rm{Z}}}}+{{{mathcal{H}}}}_{{{rm{Exch}}}}+{{{mathcal{H}}}}_{{{rm{Aniso}}}},$$
(5)

with the three terms relating to the Zeeman energy ({{{mathcal{H}}}}_{{{rm{Z}}}}=-{mu }_{S}{sum }_{i}{{{bf{S}}}}_{i}{{bf{H}}}), the exchange interaction energy in the Heisenberg form ({{{mathcal{H}}}}_{{{rm{Exch}}}}=-{sum }_{ine j}{J}_{i,j}{{{bf{S}}}}_{i}{{{bf{S}}}}_{j}), and the magnetocrystalline anisotropy energy (MAE)48. Here, Si is a unit vector denoting the local magnetic moment direction at the site i, μS = gμBS is an actual value of magnetic atomic moment with a spin S = 2, the Landé g-factor g = 2.0, and H is the external magnetic field. Following ref. 38, we parameterize the exchange energies Jijvs. Mn-Mn distances Rij by the exponential function ({J}_{ij}={J}_{0}{e}^{-{R}_{ij}/b}), with b = 1.11 Å and J0 52.9 meV.

The anisotropy term ({{{mathcal{H}}}}_{{{rm{Aniso}}}}) contains both the trigonal (uniaxial along [001]c) ({{{mathcal{H}}}}_{{{rm{Krig}}}}) and triaxial ({{{mathcal{H}}}}_{{{rm{JT}}}}) components (refs. 8,40,48,69,70). The uniaxial term includes the single ion magnetic anisotropy of the Mn3+ ions, driven by the trigonal distortion specific to the deviation from the perfect wurtzite atoms arrangement in the two-atom basis of GaN lattice. The triaxial term describes the influence of the Jahn-Teller deformation, which in wurtzite crystals grown along c-axis is directed along three inclined out-of-plane directions. These three Jahn-Teller centers are distorted along ({e}_{{{rm{JT}}}}^{A}=[sqrt{2/3},0,sqrt{1/3}]), ({e}_{{{rm{JT}}}}^{B}=[-sqrt{1/6},-sqrt{1/2},sqrt{1/3}]) and ({e}_{{{rm{JT}}}}^{C}=[-sqrt{1/6},sqrt{1/2},sqrt{1/3}]). It is important to remark here that it is precisely the presence of the three energy barriers that impedes in plane reversal of M and allows formation of hysteresis loops. The anisotropy terms take the form,

$${{{mathcal{H}}}}_{{{rm{Trig}}}}=-frac{1}{2}{K}_{{{rm{Trig}}}}{sum}_{i}frac{1}{2}left[{S}_{iz}^{2}-left({S}_{ix}^{2}+{S}_{iy}^{2}right)right]$$
(6)
$${{{mathcal{H}}}}_{{{rm{JT}}}}=-frac{1}{2}{K}_{{{rm{JT}}}}{sum}_{i}{sum}_{j=A,B,C}{left({{{bf{S}}}}_{i}cdot {{{bf{e}}}}_{{{rm{JT}}}}^{j}right)}^{4},$$
(7)

where KTrig and KJT are anisotropy parameters for uniaxial and Jahn-Teller deformations, respectively. The trigonal anisotropy energy can be expressed in more than one way due to the normalization condition ({S}_{x}^{2}+{S}_{y}^{2}+{S}_{z}^{2}=1). The lowest-order term is (-frac{1}{2}{K}_{{{rm{Trig}}}}{S}_{z}^{2}=-frac{1}{2}{K}_{{{rm{Trig}}}}(1-{S}_{x}^{2}-{S}_{y}^{2})). For KTrig > 0 the zc axis is magnetic easy axis and the xy plane is the hard one. Since constant term can be omitted, we have opted for the form presented in Eq. (6).

For the modelling of the magnetoelectric effect, the trigonal term is of a particular importance. It is exactly this term of ({{mathcal{H}}}) which is affected by the application of the electric field E and we compute the magnitude of gate driven effects through the relation KTrig(E) ξ(E).

Simulations of the dynamical properties of the system are based on the stochastic Landau-Lifshitz-Gilbert (LLG) equation applied at the atomistic level26,27,48, which is given by

$$frac{partial {{{bf{S}}}}_{i}}{partial t}=-frac{gamma }{1+{alpha }_{G}^{2}}left[{{{bf{S}}}}_{i}times {{{bf{H}}}}_{{{rm{eff}}}}^{i,tot}+{alpha }_{G}{{{bf{S}}}}_{i}times left({{{bf{S}}}}_{i}times {{{bf{H}}}}_{{{rm{eff}}}}^{i,tot}right)right],$$
(8)

where γ is the gyromagnetic ratio and αG = 1 is the precession damping term. The total effective magnetic field acting on spin i is composed of two terms ({{{bf{H}}}}_{{{rm{eff}}}}^{i,tot}={{{bf{H}}}}_{{{rm{eff}}}}^{i}+{{{bf{H}}}}_{th}^{i}) (refs. 25,26,27). The net effective field ({{{bf{H}}}}_{{{rm{eff}}}}^{i}=-(1/{mu }_{S})partial {{mathcal{H}}}/partial {{{bf{S}}}}_{i}) is obtained from the derivative of the total Hamiltonian with respect to the atomic moment μSSi. The thermal field in each spatial dimension is represented by a normal distribution Γ(t) with a standard deviation of 1 and mean of zero, and is given by

$${{{bf{H}}}}_{th}^{i}={{mathbf{Gamma }}}(t)sqrt{frac{2{alpha }_{G}{k}_{B}T}{gamma {mu }_{S}Delta t}},$$
(9)

where Δt = 5  10−6 ns is the integration time step and T is temperature.

In atomistic spin model simulations, the parallelization of the code is implemented. We divide the simulated region into sections, where each processor simulates specified part of the complete system. The computations are accelerated using graphics processing units (GPUs). We use Euler’s method as the integration scheme. In order to speed up the calculations, strict synchronization of the processors (cores) is applied every fifth step of the LLG integration. The scaling and parallelization of the code is found to be very good. We obtain practically the same results independent of number of nodes or cores used in simulations as well as independent of used computer architecture (based on GPU or CPU). As the velocity of magnetization precession depends on the strength of the magnetic field, we assume that the initialization and averaging steps of the LLG simulations depend on the field value, as given in Table 2.

Table 2 Parameters of atomistic spin model and LLG simulations
Full size table

Related Articles

Iron homeostasis and ferroptosis in muscle diseases and disorders: mechanisms and therapeutic prospects

The muscular system plays a critical role in the human body by governing skeletal movement, cardiovascular function, and the activities of digestive organs. Additionally, muscle tissues serve an endocrine function by secreting myogenic cytokines, thereby regulating metabolism throughout the entire body. Maintaining muscle function requires iron homeostasis. Recent studies suggest that disruptions in iron metabolism and ferroptosis, a form of iron-dependent cell death, are essential contributors to the progression of a wide range of muscle diseases and disorders, including sarcopenia, cardiomyopathy, and amyotrophic lateral sclerosis. Thus, a comprehensive overview of the mechanisms regulating iron metabolism and ferroptosis in these conditions is crucial for identifying potential therapeutic targets and developing new strategies for disease treatment and/or prevention. This review aims to summarize recent advances in understanding the molecular mechanisms underlying ferroptosis in the context of muscle injury, as well as associated muscle diseases and disorders. Moreover, we discuss potential targets within the ferroptosis pathway and possible strategies for managing muscle disorders. Finally, we shed new light on current limitations and future prospects for therapeutic interventions targeting ferroptosis.

Type 2 immunity in allergic diseases

Significant advancements have been made in understanding the cellular and molecular mechanisms of type 2 immunity in allergic diseases such as asthma, allergic rhinitis, chronic rhinosinusitis, eosinophilic esophagitis (EoE), food and drug allergies, and atopic dermatitis (AD). Type 2 immunity has evolved to protect against parasitic diseases and toxins, plays a role in the expulsion of parasites and larvae from inner tissues to the lumen and outside the body, maintains microbe-rich skin and mucosal epithelial barriers and counterbalances the type 1 immune response and its destructive effects. During the development of a type 2 immune response, an innate immune response initiates starting from epithelial cells and innate lymphoid cells (ILCs), including dendritic cells and macrophages, and translates to adaptive T and B-cell immunity, particularly IgE antibody production. Eosinophils, mast cells and basophils have effects on effector functions. Cytokines from ILC2s and CD4+ helper type 2 (Th2) cells, CD8 + T cells, and NK-T cells, along with myeloid cells, including IL-4, IL-5, IL-9, and IL-13, initiate and sustain allergic inflammation via T cell cells, eosinophils, and ILC2s; promote IgE class switching; and open the epithelial barrier. Epithelial cell activation, alarmin release and barrier dysfunction are key in the development of not only allergic diseases but also many other systemic diseases. Recent biologics targeting the pathways and effector functions of IL4/IL13, IL-5, and IgE have shown promising results for almost all ages, although some patients with severe allergic diseases do not respond to these therapies, highlighting the unmet need for a more detailed and personalized approach.

Management practices and manufacturing firm responses to a randomized energy audit

Increasing the efficiency of industrial energy use is widely considered important for mitigating climate change. We randomize assignment of an energy audit intervention aimed at improving energy efficiency and reducing energy expenditures of small- and medium-sized metal processing firms in Shandong Province, China, and examine impacts on energy outcomes and interactions with firms’ management practices. We find that the intervention reduced firms’ unit cost of electricity by 8% on average. Firms with more developed structured management practices showed higher rates of recommendation adoption. However, the post-intervention electricity unit cost reduction is larger in firms with less developed practices, primarily driven by a single recommendation that corrected managers’ inaccurate reporting of transformer usage at baseline, lowering their electricity costs. By closing management-associated gaps in awareness of energy expenditures, energy audit programmes may reduce a firm’s unit cost of energy but have an ambiguous impact on energy use and climate change.

Optical sorting: past, present and future

Optical sorting combines optical tweezers with diverse techniques, including optical spectrum, artificial intelligence (AI) and immunoassay, to endow unprecedented capabilities in particle sorting. In comparison to other methods such as microfluidics, acoustics and electrophoresis, optical sorting offers appreciable advantages in nanoscale precision, high resolution, non-invasiveness, and is becoming increasingly indispensable in fields of biophysics, chemistry, and materials science. This review aims to offer a comprehensive overview of the history, development, and perspectives of various optical sorting techniques, categorised as passive and active sorting methods. To begin, we elucidate the fundamental physics and attributes of both conventional and exotic optical forces. We then explore sorting capabilities of active optical sorting, which fuses optical tweezers with a diversity of techniques, including Raman spectroscopy and machine learning. Afterwards, we reveal the essential roles played by deterministic light fields, configured with lens systems or metasurfaces, in the passive sorting of particles based on their varying sizes and shapes, sorting resolutions and speeds. We conclude with our vision of the most promising and futuristic directions, including AI-facilitated ultrafast and bio-morphology-selective sorting. It can be envisioned that optical sorting will inevitably become a revolutionary tool in scientific research and practical biomedical applications.

Power price stability and the insurance value of renewable technologies

To understand if renewables stabilize or destabilize electricity prices, we simulate European power markets as projected by the National Energy and Climate Plans for 2030 but replicating the historical variability in electricity demand, the prices of fossil fuels and weather. We propose a β-sensitivity metric, defined as the projected increase in the average annual price of electricity when the price of natural gas increases by 1 euro. We show that annual power prices spikes would be more moderate because the β-sensitivity would fall from 1.4 euros to 1 euro. Deployment of solar photovoltaic and wind technologies exceeding 30% of the 2030 target would lower it further, below 0.5 euros. Our framework shows that this stabilization of prices would produce social welfare gains, that is, we find an insurance value of renewables. Because market mechanisms do not internalize this value, we argue that it should be explicitly considered in energy policy decisions.

Responses

Your email address will not be published. Required fields are marked *