Ultra-low-energy Electric field-induced Magnetization Switching in Multiferroic Heterostructures

TL;DR

This paper reviews recent advances in electric field-induced magnetization switching in multiferroic heterostructures, highlighting new phenomena, dynamics, and performance metrics relevant for low-energy magnetic memory technologies.

Contribution

It introduces novel insights into magnetization dynamics without potential tilting, and evaluates switching performance metrics in multiferroic heterostructures with various coupling mechanisms.

Findings

Magnetization can switch without potential tilting due to out-of-plane excursions.

Switching dynamics are influenced by interface coupling and exchange interactions.

Performance metrics like delay and energy dissipation are promising for applications.

Abstract

Electric field-induced magnetization switching in multiferroics is intriguing for both fundamental studies and potential technological applications. Here, we review the recent developments on electric field-induced magnetization switching in multiferroic heterostructures. Particularly, we study the dynamics of magnetization switching between the two stable states in a shape-anisotropic single-domain nanomagnet using stochastic Landau-Lifshitz-Gilbert (LLG) equation in the presence of thermal fluctuations. For magnetostrictive nanomagnets in strain-coupled multiferroic composites, such study of magnetization dynamics, contrary to steady-state scenario, revealed intriguing new phenomena on binary switching mechanism. While the traditional method of binary switching requires to tilt the potential profile to the desired state of switching, we show that no such tilting is necessary to switch successfully since the magnetization's excursion out of magnet's plane can generate a built-in asymmetry during switching. We also study the switching dynamics in multiferroic heterostructures having magnetoelectric coupling at the interface and magnetic exchange coupling that can facilitate to maintain the direction of switching with the polarity of the applied electric field. We calculate the performance metrics like switching delay and energy dissipation during switching while simulating LLG dynamics. The performance metrics turn out to be very encouraging for potential technological applications.

Equations

$$E(\theta ,\varphi ,t)=E_{shape}(\theta ,\varphi )+E_{stress}(\theta ,t),$$

$$E_{shape}(\theta ,\varphi )=\frac{{\mu }_0}{2}{M}_s^2\Omega N_d(\theta ,\varphi ),$$

$$\begin{gathered} N_d(\theta ,\varphi )=N_{d-zz}co{s}^2\theta +N_{d-yy}si{n}^2\theta si{n}^2\varphi \\ +N_{d-xx}si{n}^2\theta co{s}^2\varphi \end{gathered}$$

$$E_{barrier}=\frac{{\mu }_0}{2}{M}_s^2\Omega \left(N_{d-yy}-N_{d-zz}\right).$$

$$E_{stress}(\theta ,t)=-(3/2){\lambda }_s\sigma (t)\Omega co{s}^2\theta ,$$

$$E(\theta ,\varphi ,t)=B(\varphi ,t)si{n}^2\theta +C(t),$$

$$B(\varphi ,t)$$

$$B_{shape}(\varphi )$$

$$H_k$$

$$H_d$$

$$B_{stress}(t)$$

$$C(t)$$

$${\mathbf{T}}_{\mathbf{E}}(\theta ,\varphi ,t)$$

$$B_{shape,\varphi }(\varphi )=\frac{{\mu }_0}{2}{M}_s^2\Omega (N_{d-xx}-N_{d-yy})sin(2\varphi ).$$

$$h_i(t)=\sqrt{\frac{2\alpha kT}{|\gamma |M_V\Delta t}}{G}_{(0,1)}(t)(i\in x,y,z),$$

$${\mathbf{H}}_{\mathbf{TH}}(\theta ,\varphi ,t)=P_{\theta }(\theta ,\varphi ,t){\hat{\mathbf{e}}}_{\theta }+P_{\varphi }(\theta ,\varphi ,t){\hat{\mathbf{e}}}_{\varphi },$$

$${\mathbf{T}}_{\mathbf{TH}}(\theta ,\varphi ,t)$$

$$P_{\theta }(\theta ,\varphi ,t)$$

$$P_{\varphi }(\theta ,\varphi ,t)$$

$$\frac{d{\mathbf{n}}_{\mathbf{m}}}{dt}-\alpha \left({\mathbf{n}}_{\mathbf{m}}\times \frac{d{\mathbf{n}}_{\mathbf{m}}}{dt}\right)=-\frac{|\gamma |}{M_V}\left[{\mathbf{T}}_{\mathbf{E}}+{\mathbf{T}}_{\mathbf{TH}}\right],$$

$$\begin{gathered} \left(1+{\alpha }^2\right)\frac{d\theta }{dt}=\frac{|\gamma |}{M_V}[B_{shape,\varphi }(\varphi )sin\theta \\ -2\alpha B(\varphi ,t)sin\theta cos\theta \\ +\left(\alpha P_{\theta }(\theta ,\varphi ,t)+P_{\varphi }(\theta ,\varphi ,t)\right)], \end{gathered}$$

$$\begin{gathered} \left(1+{\alpha }^2\right)\frac{d\varphi }{dt}=\frac{|\gamma |}{M_V}[\alpha B_{shape,\varphi }(\varphi )+2B(\varphi ,t)cos\theta \\ -\{sin\theta {\}}^{-1}\left(P_{\theta }(\theta ,\varphi ,t)-\alpha P_{\varphi }(\theta ,\varphi ,t)\right)](sin\theta \ne 0). \end{gathered}$$

$$\varphi (t)=ta{n}^{-1}\left(\frac{\alpha h_y(t)+h_x(t)}{h_y(t)-\alpha h_x(t)}\right),$$

$$\frac{d\theta }{dt}=-|\gamma |\frac{h_x^2(t)+h_y^2(t)}{\sqrt{(h_y(t)-\alpha h_x(t){)}^2+(\alpha h_y(t)+h_x(t){)}^2}}.$$

$$E_{mag}(t)=-M_{V}Hsin\theta (t)cos\varphi (t),$$

$${\mathbf{T}}_{\mathbf{M}}(t)=-{\mathbf{n}}_{\mathbf{m}}(t)\times \nabla E_{mag}(\theta ,\varphi ).$$

$$\begin{gathered} \left(1+{\alpha }^2\right)\frac{d\theta }{dt}=\frac{|\gamma |}{M_V}[B_{shape,\varphi }(t)sin\theta \\ -2\alpha B(\varphi ,t)sin\theta cos\theta \\ +\alpha M_{V}Hcos\theta cos\varphi -M_{V}Hsin\varphi \\ +\left(\alpha P_{\theta }(\theta ,\varphi ,t)+P_{\varphi }(\theta ,\varphi ,t)\right)], \end{gathered}$$

$$\begin{gathered} \left(1+{\alpha }^2\right)\frac{d\varphi }{dt}=\frac{|\gamma |}{M_V}[\alpha B_{shape,\varphi }(t)+2B(\varphi ,t)cos\theta \\ -[sin\theta {]}^{-1}\left(M_{V}Hcos\theta cos\varphi +\alpha M_{V}Hsin\varphi \right) \\ -[sin\theta {]}^{-1}\left(P_{\theta }(\theta ,\varphi ,t)-\alpha P_{\varphi }(\theta ,\varphi ,t)\right)](sin\theta \ne 0). \end{gathered}$$

$${\varphi }_{min}=co{s}^{-1}\left[\frac{H}{M_s(N_{d-xx}-N_{d-yy})}\right].$$

$$P_d(t)=\frac{\alpha |\gamma |}{(1+{\alpha }^2)M_V}{\left|{\mathbf{T}}_{\mathbf{E}}(t)+{\mathbf{T}}_{\mathbf{M}}(t)\right|}^2.$$

$$E_I(\theta ,t)=-M_V{H}_I(t)cos\theta ,$$

$${\mathbf{T}}_{\mathbf{I}}(\theta ,t)=-{\mathbf{n}}_{\mathbf{m}}\times \nabla E_I(\theta ,t)=-M_V{H}_I(t)sin\theta {\hat{\mathbf{e}}}_{\varphi }.$$

$$\frac{d{\mathbf{n}}_{\mathbf{m}}}{dt}-\alpha \left({\mathbf{n}}_{\mathbf{m}}\times \frac{d{\mathbf{n}}_{\mathbf{m}}}{dt}\right)=-\frac{|\gamma |}{M_V}\left[{\mathbf{T}}_{\mathbf{I}}+{\mathbf{T}}_{\mathbf{TH}}\right].$$

Full text

{history}

1 Introduction

Multiferroics are intriguing materials in which there is intrinsic coupling between one or more ferroic properties, e.g., ferroelectric, ferromagnetic, ferroelastic ordering.[1, 2, 3, 4, 5, 6, 7, 8] In multiferroic magnetoelectrics, application of an electric field can rotate the magnetization via converse magnetoelectric effect, however, such materials in single-phase were thought to be rare.[9] Usually such single-phase multiferroics have issues of weak coupling between polarization and magnetization, and also they usually operate only at low temperatures.[10, 11] Although, new concepts may come along on switching in single-phase materials[11, 12, 13, 14, 15, 16, 17, 18] possibly utilizing Dzyaloshinsky-Moriya (DM) interaction,[19, 20, 21, 22] magnetoelectric coupling in strain-mediated multiferroic composites consisting of a piezoelectric layer coupled to a magnetostrictive nanomagnet (see Fig. 1) is shown to be very effective.[10, 23, 3, 24, 25, 26, 27] Electric field-induced magnetization switching in multiferroics uses a voltage directly thereby eliminating the need to switch magnetization by a cumbersome magnetic field or by a spin-polarized current,[28, 29, 30, 31] although new concepts are being investigated e.g., utilizing giant spin-Hall effect.[32]

The calculation of performance metrics using stochastic Landau-Lifshitz-Gilbert (LLG) equation of magnetization dynamics in strain-mediated multiferroic composites is shown to hold profound promise for computing[7, 33, 34, 35, 36] in beyond Moore’s law era.[37, 38, 39, 40, 41] With a suitable choice of materials and dimensions, when a voltage of few millivolts is applied across such heterostructures, the piezoelectric layer gets strained and the strain is elastically transferred to the magnetostrictive nanomagnet, which can rotate the magnetization. Such switching mechanism dissipates a minuscule amount of energy of $\sim$1 attojoule (aJ) in sub-nanosecond switching delay at room-temperature.[7, 42] This study has opened up a new field called straintronics.[7, 33, 34, 35] Experimental efforts have been undergoing to investigate such device functionality and the induced stress anisotropy in magnetostrictive nanomagnets is demonstrated.[43, 44, 45, 46, 47, 48, 49, 50] The direct experimental demonstration of switching speed (rather than ferromagnetic resonance experiments to get the time-scale of switching) and using low-thickness piezoelectric layers avoiding considerable degradation of the piezoelectric constants [e.g., $<$ 100 nm of lead magnesium niobate-lead titanate (PMN-PT)][51, 52] are still under investigation.

{figurehere}

Schematic diagram of a strain-mediated multiferroic composite (piezoelectric-magnetostrictive heterostructure), and axis assignment. The magnetostrictive nanomagnet is shaped like an elliptical cylinder and it has a single ferromagnetic domain. The mutually anti-parallel magnetization states along the $z$-axis can store a binary bit of information (0 or 1). In standard spherical coordinate system, $\theta$ is the polar angle and $\phi$ is the azimuthal angle. According to the dimensions of the elliptical cylinder in the respective directions, we term the $z$-axis the easy axis, the $y$-axis the in-plane hard axis, and the $x$-axis the out-of-plane hard axis for the nanomagnet.

The binary switching of magnetization between two stable states of a shape-anisotropic single-domain magnetostrictive nanomagnet in a strain-mediated multiferroic composite is very intriguing.[53, 54] The understanding behind binary switching from one stable state to another enables us to design a better switch to address our ever-increasing demand to store, process, and communicate information. From around mid-nineties, the methodology of binary switching for information processing as conceived by famous scientist Landauer and others says that an externally introduced tilt or asymmetry in potential landscape of a bistable element in the desired direction of switching is necessary.[55, 56, 57, 58] This asymmetry in potential landscape can be achieved by utilizing an external magnetic field in a single-domain nanomagnet with two stable states acting as a binary switch.[58] The tilt generates a motion along the direction of switching and the degree of tilt should be sufficient enough to dissuade thermal fluctuations with a tolerable error probability. Such tilt or asymmetry in potential landscape is deemed to be necessary for switching to take place successfully. Fig. 1(a) depicts such traditional methodology of binary switching.

{figurehere}

Traditional methodology versus the proposed methodology of binary switching. (a) In traditional methodology of binary switching, the potential landscape is tilted towards the direction of switching alongwith the lowering of energy barrier separating the two stable states. Note that there are two external agents involved here, one makes the potential landscape monostable and the other one that tilts the potential landscape. At the end, the potential profile is restored back to that of the initial stage to complete the switching process. (b) In the proposed methodology of binary switching, the potential profile remains always ‘symmetric’, i.e., energy barrier is lowered but the potential landscape is not tilted to favor the final state. Switching occurs due to internal dynamics considering complete three-dimensional potential landscape and full three-dimensional motion. Note that in this case only one external agent is involved since it does not require tilting the potential landscape. For a nanomagnet acting as a binary switch, the bold line in its potential landscape corresponds to when magnetization resides on magnet’s plane. Deflection of magnetization from magnet’s plane corresponds to out-of-plane excursion of magnetization, which is the physical mechanism of switching in this case. (Reprinted from Ref. 53.)

Note that the magnetization switching mechanism using spin-transfer-torque,[28, 29, 30] in which a spin-polarized current is passed through a nanomagnet to switch its magnetization, is analogous to the traditional methodology of binary switching. It is well known that the Slonczewski-like spin-transfer-torque that acts in-plane of the nanomagnet cannot be treated as an effective potential,[28] however, the direction of externally applied spin-polarized current induces an equivalent tilt or asymmetry and causes the magnetization switching in the desired direction.

In Ref. \refciteroy13_2, it is shown that it is not necessary to tilt the potential landscape by external means to switch successfully, even in the presence of thermal fluctuations. Fig. 1(b) illustrates the basic concept underlying such methodology of binary switching. To understand such switching, the complete three-dimensional potential landscape needs to be considered and thereby it is demonstrated in Ref. \refciteroy13_2 that the intrinsic dynamics can provide an equivalent asymmetry without any requirement of making the potential landscape asymmetric. It needs mention here that for both traditional and proposed methodologies, it is necessary to lower the energy barrier and make the monostable well deep enough to dissuade thermal fluctuations. For a magnetostrictive nanomagnet, stress is the external agent modulating and eventually inverting the potential landscape of the nanomagnet. The voltage-induced stress can be generated on a magnetostrictive nanomagnet by elastically coupling a piezoelectric layer, i.e., using 2-phase multiferroic composites.[10, 3, 24, 59, 60, 61, 33] Such magnetization switching in multiferroic composites can potentially be the basis of ultra-low-energy computing in our future information processing systems.[7, 33, 53, 42]

If we consider only the steady-state scenario, we will come to a conclusion that in a strain-mediated multiferroic composite, the strain transferred by the piezoelectric layer to the magnetostrictive nanomagnet, can only rotate the magnetization $90^{\circ}$ and not the complete $180^{\circ}$, which is the usual perception. However, a complete $180^{\circ}$ switching facilitates us to achieve a higher magnetoresistance while reading the magnetization state using magnetic tunnel junctions (MTJs).[62, 63, 64, 65, 66, 67, 68] Although there are proposals of $90^{\circ}$ switching mechanism,[69, 70, 43] as shown in Ref. \refciteroy13_2 and described above, a complete $180^{\circ}$ switching is possible if we consider the dynamics of magnetization into account rather than assuming steady-state scenario. Basically magnetization’s excursion out of magnet’s plane provides an equivalent asymmetry to cause a complete $180^{\circ}$ switching.[53]

{figurehere}

Schematics of the interface and exchange coupled multiferroic heterostructures. The unique coupling between the polarization in the P-layer and the trilayer $M_{1}$/spacer/$M_{2}$ allows the polarization direction to dictate the magnetic ground state in the trilayer. If the polarization points downward ($P_{\downarrow}$), P-alignment in the trilayer is preferred while an upward polarization ($P_{\uparrow}$) prefers the AP-alignment. Application of a voltage with correct polarity can switch the polarization and hence the magnetization $M_{1}$ gets switched too due to interface and exchange coupling. At the bottom of the figure, the axis assignment for the dynamical motion of magnetization $M_{1}$ in standard spherical coordinate system is shown. (©IOP Publishing. Reproduced by permission of IOP Publishing from Ref. 87. All rights reserved.)

Although the aforesaid switching mechanism in strain-mediated multiferroic composites is intriguing and promising for technological applications, it would be of substantial interest if there exists a strong coupling between polarization and magnetization at the heterostructure interface to harness additional asymmetry during magnetization switching in a specified direction. In Ref. \refciteRefWorks:649, interface and exchange coupled multiferroics are proposed based on density functional theory (DFT) of first-principles calculations. Although this specific case needs to be experimentally demonstrated, the proposed concept is tenable. The first-principles calculations have been proved to be very useful in such respect[9] and with the experimental progress on similar front[72, 73, 74, 75, 76] (also using ferromagnetic oxides[77, 78] rather than ferromagnetic metals), there is a considerable promise on such polarization-magnetization coupling mechanism.[71, 79, 80, 81, 82, 83] In Figure 1, such interface and exchange coupling between polarization and magnetization in a multiferroic heterostructure is depicted.[71] The polarization direction in the P-layer uniquely determines the ground state of the trilayer $M_{1}$/spacer/$M_{2}$, i.e., for downward polarization P↓, parallel alignment (P-alignment) in the trilayer is preferred, while for the upward polarization P↑, antiparallel alignment (AP-alignment) in the trilayer is achieved. This polarization can be switched electrically and such polarization-magnetization coupling mechanism makes the switching of magnetization in the $M_{1}$-layer non-volatile. If a voltage with certain polarity is applied and maintained, the state of the system remains unaltered too. This is advantageous over the strain-mediated switching, which just toggles the magnetization states and therefore requires a read-before-write mechanism. Note that there are other exchanged coupled systems with an insulating spacer layer, however, the interlayer exchange coupling energy is small.[84, 85] There are device structures with non-magnetic spacer layer to preserve large interlayer exchange coupling too, but a high electric field is required and the switching becomes volatile.[86]

In Ref. \refciteroy14_2, magnetization dynamics in the interface and exchange coupled multiferroic heterostructures is studied by solving stochastic Landau-Lifshitz-Gilbert equation in the presence of room-temperature thermal fluctuations. Such phenomenological study of switching has been very useful to understand the performance metrics of magnetic devices.[7, 33, 34, 35, 42, 53] First, an interfacial anisotropy in the interface and exchange coupled multiferroic heterostructures is modeled and subsequently the analysis is performed on the dynamics of magnetization. The results show that switching in sub-nanosecond delay is possible while expending only $\sim$1 aJ of energy at room-temperature. The key point is that the strong interface anisotropy makes the switching error-resilient and fast, allowing to use nanomagnets with very small dimensions (magnetization is stable with $\sim$10 nm lateral dimensions even in the presence of room-temperature thermal fluctuations). Such performance metrics of area, delay, and energy are particularly promising for computing in our future information processing systems.[37, 38, 39, 40, 41]

The rest of the paper is organized as follows. In Section 2, we describe in subsequent subsections two models for electric field-induced magnetization switching: (1) strain-mediated multiferroic composites, and (2) interface and exchange coupled multiferroic heterostructures. For both the models, the stochastic Landau-Lifshitz-Gilbert (LLG) equation of magnetization dynamics in the presence of thermal fluctuations is analytically solved to get a coupled set of equations, which need to be solved numerically onwards. Section 3 presents the simulation results by solving the coupled sets of equations numerically. Then we calculate the performance metrics e.g., switching delay, energy dissipation, which show particularly promising results for technological applications. Finally, Section 4 summarizes this review and provides the outlook on the electric field-induced magnetization switching in multiferroic heterostructures.

2 Model

Here, we will review the models developed for strain-mediated multiferroic composites, and interface and exchange coupled multiferroic heterostructures in subsequent subsections. The emphasis would be on the dynamical motion by solving the Landau-Lifshitz-Gilbert (LLG) equation of magnetization dynamics.[88, 89] We will also consider thermal fluctuations incorporated in the LLG equation making it of stochastic nature.[90, 53, 42, 87]

2.1 Strain-mediated multiferroic composites

Consider a nanomagnet shaped like an elliptical cylinder with its elliptical cross section lying on the $y$-$z$ plane; the major axis is aligned along the $z$-direction and the minor axis along the $y$-direction. (See Fig. 1.) The dimension of the major axis, the minor axis, and the thickness are $a$, $b$, and $l$, respectively. The volume of the nanomagnet is $\Omega=(\pi/4)abl$. In standard spherical coordinate system, as shown in the Fig. 1, $\theta$ is the polar angle and $\phi$ is the azimuthal angle. Note that when $\phi=\pm 90^{\circ}$, the magnetization vector lies on the plane of the nanomagnet ($y$-$z$ plane). Any deviation from $\phi=\pm 90^{\circ}$ is termed as out-of-plane excursion of magnetization.

The total energy of the single-domain, magnetostrictive, polycrystalline (i.e., no net magnetocrystalline anisotropy) nanomagnet, while subjected to uniaxial stress along the easy axis ($z$-axis, the major axis of the ellipse) is the sum of the shape anisotropy energy and the uniaxial stress anisotropy energy:[91]

[TABLE]

where $E_{shape}(\theta,\phi)$ is the shape anisotropy energy and $E_{stress}(\theta,t)$ is the stress anisotropy energy at time $t$. The former is given by[91]

[TABLE]

where $M_{s}$ is the saturation magnetization and $N_{d}(\theta,\phi)$ is the demagnetization factor expressed as[91]

[TABLE]

with $N_{d-zz}$, $N_{d-yy}$, and $N_{d-xx}$ being the components of the demagnetization factor along the $z$-axis, $y$-axis, and $x$-axis, respectively. The parameters $N_{d-zz}$, $N_{d-yy}$, and $N_{d-xx}$ depend on the shape and dimensions of the nanomagnet[91] and are determined from the prescription in Ref. \refciteRefWorks:402.

The in-plane ($\phi=\pm 90^{\circ}$) shape-anisotropic energy-barrier between the two stable states ($\theta=0^{\circ}$ and $180^{\circ}$) can be expressed as

[TABLE]

With $a=100\,nm$, $a=90\,nm$, and $l=6\,nm$, $E_{barrier}\simeq 44\,kT$ at room-temperature ($T=300\,K$), which means that the static error probability due to spontaneous fluctuation of magnetization is $e^{-E_{barrier}/kT}=e^{-44}$.

The uniaxial shape anisotropy favors lining up the magnetization along the major axis ($z$-axis) by minimizing $E_{shape}$, which is why we call the major axis the “easy axis” and the minor axis ($y$-axis) the “in-plane hard axis” of the magnet. The $x$-axis will therefore be the “out-of-plane hard axis” of the magnet and it is “harder” than the in-plane one since the thickness is much smaller than the magnet’s lateral dimensions (i.e., $l<<a,b$).

We assume that an uniaxial stress is generated along the $z$-axis (easy axis) upon application of an electric field along the $x$-axis. It is possible to constrain the expansion along $y$-axis to generate the uniaxial stress along the $z$-direction for piezoelectrics like lead-zirconate-titanate (PZT). We can also use piezoelectrics like lead magnesium niobate-lead titanate (PMN-PT), which generates anisotropic strain in the lateral plane (i.e., the signs of $d_{31}$ and $d_{32}$ coefficients are different) and therefore it can produce more strain for a given voltage or otherwise a lower voltage is required to generate a specified strain. The stress anisotropy energy is given by[91]

[TABLE]

where $(3/2)\lambda_{s}$ is the magnetostriction coefficient of the single-domain nanomagnet and $\sigma(t)$ is the stress at an instant of time $t$. Note that a positive $\lambda_{s}\sigma$ product will favor alignment of the magnetization along the major axis ($z$-axis), while a negative $\lambda_{s}\sigma$ product will favor alignment on the $x$-$y$ plane ($\theta=90^{\circ}$), because that will minimize $E_{stress}$. In our convention, a compressive stress is negative and tensile stress is positive. Therefore, in a material like Terfenol-D that has positive $\lambda_{s}$, a compressive stress will favor alignment on the $x$-$y$ plane ($\theta=90^{\circ}$), and tensile along the major axis. The situation will be opposite with common magnetic materials like iron, nickel, and cobalt, which have negative $\lambda_{s}$.

At any instant of time $t$, the total energy of the nanomagnet can be expressed as

[TABLE]

where

[TABLE]

$H_{k}$ is the Stoner-Wohlfarth switching field,[58] and $H_{d}$ is the out-of-plane demagnetization field.[91] Note that $B_{stress}$ has the same sign as the $\lambda_{s}\sigma$ product. It will be negative if we use stress to rotate the magnetization from the easy axis ($z$-direction) to the plane defined by the in-plane hard axis ($y$-direction) and the out-of-plane hard axis ($x$-direction), i.e., the $x$-$y$ plane ($\theta=90^{\circ}$).

In the macrospin approximation, the magnetization M(t) of the nanomagnet has a constant magnitude but a variable direction, so that we can represent it by the vector of unit norm $\mathbf{n_{m}}(t)=\mathbf{M}(t)/|\mathbf{M}|=\mathbf{\hat{e}_{r}}$ where $\mathbf{\hat{e}_{r}}$ is the unit vector in the radial direction in the standard spherical coordinate system represented by ($r$,$\theta$,$\phi$). The unit vectors for the $\theta$- and $\phi$-rotations are denoted by $\mathbf{\hat{e}_{\theta}}$ and $\mathbf{\hat{e}_{\phi}}$, respectively.

The torque due to shape and stress anisotropy is derived from the gradient of potential profile as

[TABLE]

where

[TABLE]

The effect of random thermal fluctuations is incorporated via a random magnetic field $\mathbf{h}(t)=h_{x}(t)\mathbf{\hat{e}_{x}}+h_{y}(t)\mathbf{\hat{e}_{y}}+h_{z}(t)\mathbf{\hat{e}_{z}}$, where $h_{i}(t)$ ($i=x,y,z$) are the three components of the random thermal field in Cartesian coordinates. We assume the properties of the random field $\mathbf{h}(t)$ as described in Ref. \refciteRefWorks:186. The random thermal field can be written as[90]

[TABLE]

where $\alpha$ is the dimensionless phenomenological Gilbert damping parameter, $\gamma$ is the gyromagnetic ratio for electrons, $1/\Delta t$ is the attempt frequency of thermal fluctuations, $M_{V}=\mu_{0}M_{s}\Omega$, $\Omega$ is the volume, $k$ is the Boltzmann constant, $T$ is temperature, and the quantity $G_{(0,1)}(t)$ is a Gaussian distribution with zero mean and unit variance.

The thermal field and the corresponding torque acting on the magnetization can be written as

[TABLE]

and

[TABLE]

respectively, where

[TABLE]

Additionally, there is motion due to Gilbert damping[88, 89] (perpendicular to the precessional motion) through which magnetization relaxes towards the minimum energy position on magnet’s potential landscape. The magnetization dynamics under the action of these two torques $\mathbf{T_{E}}$ and $\mathbf{T_{TH}}$ is described by the stochastic Landau-Lifshitz-Gilbert (LLG) equation as follows.

[TABLE]

where $\alpha$ is the phenomenological Gilbert damping parameter and $\gamma$ is the gyromagnetic ratio of electrons. Solving the above equations, we get the following coupled equations for the dynamics of $\theta$ and $\phi$:

[TABLE]

[TABLE]

We will ignore the random thermal torque due to room-temperature thermal fluctuations while explaining the magnetization dynamics first, following Ref. \refciteroy13_2, however, we will discuss the key consequences of incorporating thermal fluctuations with simulation results later. We assume that the magnetization starts from $\theta\simeq 180^{\circ}$ ($-z$-axis) and the applied stress attempts to switch it to $\theta\simeq 0^{\circ}$ ($+z$-axis).

Out-of-plane excursion of magnetization. We will first intuitively describe how magnetization is deflected from the magnet’s plane ($\phi=\pm 90^{\circ}$, i.e., $y$-$z$ plane), and is stabilized out-of-plane as depicted in Fig. 1(b) due to different torques acting on it. The shape anisotropy energy is in general dependent on azimuthal angle $\phi$ (rather than assuming $\phi=\pm 90^{\circ}$) and it generates additional motions of magnetization in $\hat{\mathbf{e}}_{\theta}$ and $\hat{\mathbf{e}}_{\phi}$ directions [see the motions containing the term $B_{shape,\phi}(\phi)$ in Fig. 1(c) and equations (16) and (17)]. Both of these torques are proportional to $sin(2\phi)$ and vanish when $\phi=\pm 90^{\circ}$. Also, note that the $\phi$-component is proportional to the damping parameter $\alpha$. As shown in the Fig. 1(b), the applied stress generates a torque that attempts to rotate the magnetization anticlockwise and forces the magnetization to deflect from magnet’s plane and stay out of magnet’s plane. As magnetization is deflected from the plane of the magnet ($\phi=\pm 90^{\circ}$), the $\phi$-component of the torque due to shape anisotropy energy as mentioned earlier [$\propto\,\alpha\,sin(2\phi)$] would attempt to bring the magnetization back to magnet’s plane. Because of this counteraction, the magnetization is only stable in the second or the fourth quadrant [i.e., ($90^{\circ}$, $180^{\circ}$) or ($270^{\circ}$, $360^{\circ}$)], among the four possible quadrants of $\phi$. Note that $sin(2\phi)$ is a negative quantity in these two quadrants and hence the motion due to shape anisotropy energy counteracts the precessional motion due to stress. We would term these two quadrants (second or the fourth) as “good” quadrants and the other two (first and third) quadrants as “bad” quadrants, the reasoning behind which would be more prominent onwards, i.e., consideration of the torques due to $\phi$-dependence of shape anisotropy energy is crucial to the magnetization dynamics. The key lesson is that $\phi$ becomes stable only in the “good” quadrants and it facilitates switching of magnetization in the desired direction.

Magnetization’s motion in three-dimensional space. We will now describe the motion of magnetization in the full three-dimensional space intuitively in the presence of various torques originating from the shape and stress anisotropy energy as depicted in the Fig. 1(c). Note that the motion of magnetization needs to be along the $-\hat{\mathbf{e}}_{\theta}$ direction since magnetization is being switched from $\theta\simeq 180^{\circ}$ towards $\theta\simeq 0^{\circ}$. The applied stress generates a precessional motion of magnetization in the +$\hat{\mathbf{e}}_{\phi}$ direction, and the damping of magnetization generates a motion additionally that is perpendicular to both the direction of magnetization ($\hat{\mathbf{e}}_{r}$) and +$\hat{\mathbf{e}}_{\phi}$, i.e., in $-\hat{\mathbf{e}}_{\theta}$ direction. These two motions are shown as $2B(\phi,t)cos\theta\,\hat{\mathbf{e}}_{\phi}$ and $-2\alpha B(\phi,t)sin\theta cos\theta\,\hat{\mathbf{e}}_{\theta}$, respectively in Fig. 1(c), where $\alpha$ is the damping parameter and the quantity $B(\phi,t)$ includes terms both due to the shape anisotropy energy $B_{shape}(\phi)$ and the stress anisotropy energy $B_{stress}(t)$. The quantity $B_{stress}(t)$ is negative and it must overcome the shape anisotropy term $B_{shape}(\phi)$ for switching to get started (mathematically, note that both the quantities $B(\phi,t)$ and $cos\theta$ are negative in the interval $180^{\circ}\geq\theta\geq 90^{\circ}$). Therefore, magnetization starts switching towards its desired direction due to the applied stress. Note that this damped motion in $-\hat{\mathbf{e}}_{\theta}$ direction is considerably weak because of the multiplicative damping parameter $\alpha$, which is usually much less than one (e.g., $\alpha$=0.1 for Terfenol-D).

As described earlier [see Fig. 1(b)], due to applied stress, magnetization rotates out-of-plane and stays in a “good” quadrant for $\phi$ [i.e., ($90^{\circ}$, $180^{\circ}$) or ($270^{\circ}$, $360^{\circ}$)], which generates a motion of magnetization in the $-\hat{\mathbf{e}}_{\theta}$ direction due to $\phi$-dependence of shape anisotropy energy. Thereby a damped motion is generated in the $-\hat{\mathbf{e}}_{\phi}$ direction. These two motions are shown as $-|B_{shape,\phi}(\phi)|sin\theta\,\hat{\mathbf{e}}_{\theta}$ and $-\alpha|B_{shape,\phi}(\phi)|\,\hat{\mathbf{e}}_{\phi}$, respectively in Fig. 1(c), where $B_{shape,\phi}(\phi)\propto sin(2\phi)$. In the “good” quadrants for $\phi$, $B_{shape,\phi}(\phi)$ is negative, therefore, if the magnetization stays out of magnet’s plane in a “good” quadrant, magnetization rotates in its desired direction. Since this motion does not possess any damping factor [note the other damped motion in the $-\hat{\mathbf{e}}_{\theta}$ direction in Fig. 1(c)], it can eventually increase the magnetization switching speed to a couple of orders of magnitude higher. On the other hand, if the magnetization resides out-of-plane but in a “bad” quadrant, the motion of magnetization in its desired direction of switching is hindered. If we apply a higher magnitude of stress, the magnetization is deflected out of magnet’s plane more inside a “good” quadrant (counteracting the random thermal kicks possibly acting in the opposite direction, which will be described later). Note that the damped motion $-\alpha|B_{shape,\phi}(\phi)|\,\hat{\mathbf{e}}_{\phi}$ attempts to bring magnetization back towards the magnet’s plane. As these two motions counteract each other [see Fig. 1(c)], magnetization continues moving in the $-\hat{\mathbf{e}}_{\theta}$ direction and eventually reaches the $x$-$y$ plane ($\theta=90^{\circ}$). Note that without damping, such counteraction does not happen and the magnetization just precesses through “good” and “bad” quadrants consecutively.

Upon reaching the $x$-$y$ plane ($\theta=90^{\circ}$), if magnetization stays in a “good” quadrant for $\phi$ [i.e., $(90^{\circ},180^{\circ})$ or $(270^{\circ},360^{\circ})$], then the torque on the magnetization will be in the correct direction to facilitate magnetization’s traversal towards $\theta\simeq 0^{\circ}$ [see Figs. 2(a) and 2(b)]. Once again this signifies the merit of terminology (“good” or “bad”) used for the four quadrants of $\phi$. At $\theta=90^{\circ}$ (i.e., $cos\theta=0$), the effect of stress on the magnetization rotation has diminished completely [see Fig. 1(c)]. The only two motions that are active at $\theta=90^{\circ}$ are $-|B_{shape,\phi}(\phi)|sin\,\theta\,\hat{\mathbf{e}}_{\theta}$ and $-\alpha|B_{shape,\phi}(\phi)|\,\hat{\mathbf{e}}_{\phi}$ [Fig. 1(c)]. Since $\alpha\ll 1$, magnetization quickly gets out from $\theta=90^{\circ}$ and as the magnetization vector is deflected from $\theta=90^{\circ}$ towards $\theta=0^{\circ}$, the effect of stress again comes into play.

Stress cycle, magnetization directions, and potential profiles. Fig. 3 shows the stress-cycle alongwith the energy profiles and magnetization directions at different instants of time during switching of magnetization. At time $t_{0}$, the magnetization direction is along the easy axis $\theta\simeq 180^{\circ}$ with itse potential landscape unperturbed by stress. Note that the potential profile of the magnet is ‘symmetric’ in both $\theta$- and $\phi$-space with two degenerate minima at $\theta=0^{\circ}$, 180*∘* and a maximum at $\theta=90^{\circ}$ in $\theta$-space, signifying that a binary information can be stored in the nanomagnet corresponding to $\theta=0^{\circ}$, 180*∘*. The anisotropy in the barrier is due to shape anisotropy energy of the nanomagnet only, which is $\sim$44 kT at room-temperature using the nanomagnet’s dimensions and the material parameters used for the magnetostrictive nanomagnet made of Terfenol-D (see Section 3 later). Note that the barrier height separating the two stable states ($\theta$ = $0^{\circ}$ and $180^{\circ}$) is meant when the magnetization stays in-plane (i.e., $\phi=\pm 90^{\circ}$) of the magnet. The barrier becomes higher when the magnetization is deflected from $\phi=90^{\circ}$ as shown in the Fig. 3(c) at time $t_{0}$. The barrier is the highest when the magnetization points along the out-of-plane direction ($\phi$ = $0^{\circ}$ or $180^{\circ}$), which is due to the small thickness of the nanomagnet compared to the lateral dimensions. Note that the magnetization can start from any angle $\phi_{initial}\in(0^{\circ},360^{\circ})$ in the presence of thermal fluctuations [see Fig. 6(b) later].

Fig. 3(c) depicts that as a compressive stress is ramped up on the nanomagnet between time instants $t_{0}$ and $t_{1}$ and a sufficient stress is applied, the potential landscape in $\theta$-space becomes monostable near $\phi=\pm 90^{\circ}$. Since the barrier height is high near $\phi=0^{\circ}$ or $180^{\circ}$, the potential landscape may not become monostable in $\theta$-space therein. However, that is not necessary for switching since application of stress rotates the magnetization in $\phi$-direction and theretofore the magnetization can eventually come near $\phi=\pm 90^{\circ}$, which facilitates switching from $\theta\simeq 180^{\circ}$ towards $\theta=90^{\circ}$ (see Fig. 4). From Fig. 3(c), we can see that the minimum energy position between time instants $t_{1}$ and $t_{2}$ is at ($\theta=90^{\circ}$, $\phi=\pm 90^{\circ}$) and the potential profile at time instant $t_{1}$ is still ‘symmetric’.

Fig. 3 shows that stress is held constant between time instants $t_{1}$ and $t_{2}$ and the magnetization eventually reaches at $x$-$y$ plane ($\theta=90^{\circ}$). For a sufficiently fast ramp rate and high stress, magnetization will reside in “good” quadrants (see Fig. 10 later). This will ensure that magnetization traverses in the correct direction towards $\theta=0^{\circ}$ and switches successfully. A sufficiently fast ramp rate ensures that magnetization would not backtrack towards $\theta=180^{\circ}$ even after crossing $\theta=90^{\circ}$ towards $\theta=0^{\circ}$. If stress is held constant for longer time, magnetization will have higher probability to collapse on magnet’s plane ($\phi=\pm 90^{\circ}$) following which thermal fluctuations will scuttle the magnetization either in a “good” quadrant or in a “bad” quadrant with equal probability, so the success rate would be 50%.

It should be noted from Fig. 3 that the stress is reversed (compression to tensile) between time instants $t_{2}$ and $t_{4}$ rather than just withdrawn, which makes the potential landscape of the nanomagnet more steep in $\theta$-space. However, it does not necessarily mean that the switching will be completed always faster. The reversal of stress can cause magnetization to traverse into “bad” quadrants in $\phi$-space causing magnetization to precess and therefore the switching delay may eventually increase. Particularly for higher stress levels, such increase in switching delay may happen. However, it is noticed that reversing the stress makes the success rate of switching a bit ($<$5%) higher in the presence of thermal fluctuations. The tensile stress is held constant and when $\theta$ becomes $\leq 5^{\circ}$, switching is deemed to have completed. The stress is withdrawn at the end between time instants $t_{5}$ and $t_{6}$ to complete the switching cycle.

Switching failure. As magnetization leaves from $\theta=90^{\circ}$ towards $\theta\simeq 0^{\circ}$ and stress is ramped down, the torque due to stress tries to rotate the azimuthal angle $\phi$ of magnetization clockwise rather than anticlockwise. Mathematically, note that $cos\theta$ is positive for $90^{\circ}\geq\theta\geq 0^{\circ}$ and $B(\phi,t)$ is still negative when stress has not been brought down sufficiently, i.e., still $|B_{stress}(t)|>|B_{shape}(\phi)|$ (see Fig. 1(c) and corresponding discussions in text). For a slow ramp-rate this clockwise rotation may be considerable and magnetization can stray into a “bad” quadrant. Moreover, thermal fluctuations can aggravate the scenario by possibly deflecting the magnetization into the “bad” quadrant further. Switching may impede and magnetization vector may backtrack towards $\theta\simeq 180^{\circ}$ rather than traversing towards $\theta\simeq 0^{\circ}$ causing a switching failure. Therefore, switching failure may happen even after the magnetization has crossed the hard axis ($\theta=90^{\circ}$) towards its destination $\theta\simeq 0^{\circ}$ (see Fig. 5). This is why it does require a sufficiently fast ramp rate during the ramp-down phase of stress. Such switching failure is an intriguing phenomenon, which cannot be conceived if magnetization is always assumed to be on magnet’s plane ($\phi=\pm 90^{\circ}$) and only can be understood if we analyze the magnetization dynamics in a complete three dimensional potential landscape.

Fluctuation of magnetization around the easy axis (stable orientation) due to thermal noise. From Equation (8), it can be noted that the torque on the magnetization due to shape and stress anisotropy vanishes when $sin\theta=0$, i.e. when the magnetization vector is exactly aligned along the easy axis ($\theta=0^{\circ}$, 180*∘*).[42, 53] However, thermal fluctuations can deflect the magnetization vector from the easy axis. Considering the case when $\theta=180^{\circ}$, we get the following:

[TABLE]

[TABLE]

From the above equation we can clearly follow that thermal torque can deflect the magnetization from the easy axis since the time rate of change of $\theta(t)$ is non-zero in the presence of thermal agitations. The initial deflection from the easy axis due to the thermal torque does not depend on the component of the random thermal field along the $z$-axis, i.e., $h_{z}(t)$. This is a consequence of having \pm$$z-axis as the easy axes of the nanomagnet. However, once the magnetization direction is even slightly deflected from the easy axis, all the three components of the random thermal field come into play.

Thermal distribution of the initial orientation of the magnetization. We can determine the thermal distributions of $\theta$ and $\phi$ when no stress is applied on the magnetostrictive nanomagnet by solving the Equations (16) and (17) with $B_{stress}$ = 0 (Refs. \refciteroy11_6,roy13_2). This will yield the distribution of the magnetization’s initial orientation when stress is turned on. The $\theta_{initial}$-distribution is Boltzmann peaked at $\theta$ = 0*∘* or 180*∘, while the $\phi_{initial}$-distribution is Gaussian peaked at $\phi=\pm 90^{\circ}$ (see Fig. 6). According to the Boltzmann distribution of $\theta_{initial}$, the most probable value of $\theta$ is either 0∘* or 180*∘*, where stress is ineffective. This will lead to a long tail in the switching delay distribution, which is due to the fact that when magnetization starts out from $\theta=0^{\circ},180^{\circ}$, it needs to wait a while before random thermal fluctuations can set the switching in motion. Thus, switching trajectories starting from an easy axis are very slow and it causes the long tail in the distribution of switching delay.[42, 53]

Application of a bias field to shift the peak of the initial distribution of magnetization from an easy axis. We can eliminate the long tail in the switching delay distribution by applying a small static bias magnetic field that shifts the peak of $\theta_{initial}$ distribution away from the easy axis, so that the most probable starting orientation will no longer be an easy axis.[42] This field is applied along the out-of-plane hard axis (+$x$-direction) of the nanomagnet and thus the potential energy due to the applied magnetic field becomes

[TABLE]

where $H$ is the magnitude of magnetic field. A torque is generated due to this field, which is

[TABLE]

The presence of this torque will modify Equations (16) and (17) to

[TABLE]

[TABLE]

Note that the bias field also makes the potential energy profile of the magnet asymmetric in $\phi$-space and the energy minimum is shifted from $\phi_{min}=\pm 90^{\circ}$ (the plane of the magnet) to

[TABLE]

However, the potential profile will remain symmetric in $\theta$-space, with $\theta=0^{\circ}$ and $\theta=180^{\circ}$ remaining as the minimum energy positions. With the material parameters and the dimensions of the nanomagnet (see Section 3 later), a bias magnetic field of flux density 40 mT would make $\phi_{min}\simeq\pm 87^{\circ}$, i.e. deflect the magnetization $\sim$3*∘* from the magnet’s plane. Application of the bias magnetic field also reduces the in-plane shape anisotropy energy barrier from 44 $kT$ to 36 $kT$ at room temperature. We assume that a permanent magnet will be employed to produce the bias magnetic field and not by a current-carrying coil on-chip. Therefore no additional energy dissipation needs to be considered for this reason.

Energy Dissipation. The energy dissipated during magnetization switching has two components: (1) the energy dissipated in the switching circuit for generating stress in the magnetostrictive nanomagnet with the application of a voltage, and (2) the energy dissipated internally in the nanomagnet due to Gilbert damping. We will term the first component as ‘$CV^{2}$’ dissipation, where $C$ and $V$ denote the capacitance of the piezoelectric layer and the applied voltage, respectively. If the ramp rate is infinite, i.e., the voltage is turned on or off abruptly, the energy dissipated during either turn on or turn off is $(1/2)CV^{2}$. However, if the ramp rate is finite, this energy dissipation is reduced and the exact value will depend on the ramp rate. We calculate it following the procedure described in Ref. \refciteroy11_2.

The second component, which is the internal energy dissipation $E_{d}$ due to Gilbert damping, is given by the expression $\int_{0}^{\tau}P_{d}(t)dt$, where $\tau$ is the switching delay and $P_{d}(t)$ is the power dissipated during switching

[TABLE]

We sum up the two energy dissipations ‘$CV^{2}$’ and $E_{d}$ to get the total dissipation $E_{total}$. The average power dissipation is $E_{total}/\tau$. There is no net dissipation due to random thermal torque, however, thermal fluctuations do affect both $E_{d}$ and ‘$CV^{2}$’ dissipations. It affects $E_{d}$ since it raises the critical stress needed to switch with $\sim$100% probability and ‘$CV^{2}$’ dissipation is also raised since the applied voltage is proportional to stress.

So far the model presented deals with the switching of magnetization (i.e., writing a bit of information), however, the magnetization state needs to be read too, which is usually performed by using a magnetic tunnel junction (MTJ).[62, 63, 64, 65, 66, 67, 68] While reading the bit of information, a material selection issue crops up since magnetostrictive materials in general cannot be utilized as the free layer of the MTJ. This is an important issue since we need to achieve a high magnetoresistance for successful read operation as required by technological applications. In Ref. \refciteroy15_1, it is shown that magnetically coupling the magnetostrictive nanomagnet and the free layer can circumvent this issue. Stochastic Landau-Lifshitz-Gilbert equation of magnetization dynamics in the presence of room-temperature thermal fluctuations is solved and it is shown that such design can eventually lead to a superior energy-delay product.[94]

2.2 Interface and exchange coupled multiferroic heterostructures

We have earlier shown in the Figure 1 the schematic diagram of the interface and exchange coupled multiferroic heterostructure devices and the axis assignment for the orientation of magnetization $M_{1}$. The standard spherical coordinate system with $\theta$ as polar angle and $\phi$ as azimuthal angle is chosen. The magnetization $M_{1}$ orients along $\theta=180^{\circ}$ if polarization points downward ($P_{\downarrow}$), while $M_{1}$ is along $\theta=0^{\circ}$ if polarization points upward ($P_{\uparrow}$). The lateral elliptical cross-section of $M_{1}$ lies on the $y$-$z$ plane ($\phi=\pm 90^{\circ}$) with its major axis pointing to $z$-direction and minor axis in $y$-direction. The dimensions of the nanomagnet $M_{1}$ along the $z$-, $y$-, and $x$-axis are $a$, $b$, and $l$, respectively. So the magnetic easy axis becomes along the $\pm z$-direction and the nanomagnet’s volume $\Omega=(\pi/4)abl$. When magnetization switches between its two stable states ($\theta=0^{\circ},180^{\circ}$), magnetization being a rotational body deflects out of magnet’s plane and any deflection from the magnet’s plane ($\phi=\pm 90^{\circ}$) is termed as out-of-plane excursion.

The interface anisotropy energy in the nanomagnet $M_{1}$ is modeled as[87]

[TABLE]

where $M_{V}=\mu_{0}M_{s}\Omega$, $M_{s}$ is the saturation magnetization, and $H_{I}$ is the interfacial anisotropy field. If $H_{I}=-H_{I,max}$, the magnetization $M_{1}$ points along $\theta=180^{\circ}$ and if we vary $H_{I}$ from $-H_{I,max}$ to $H_{I,max}$, the magnetization orients along $\theta=0^{\circ}$. The total anisotropy of the magnet is the sum of the interface anisotropy alongwith the other anisotropies like magnetocrystalline anisotropy and shape anisotropy,[33, 42]. However, since the interfacial anisotropy is strong compared to the other anisotropies, we consider only the interfacial anisotropy (i.e., $E_{total}\simeq E_{I}$) for brevity.

The magnetization $\mathbf{M}$ of the single-domain nanomagnet $M_{1}$ has a constant magnitude of magnetization but a variable direction. Therefore, the magnetization vector can be represented by the unit vector in the radial direction $\mathbf{\hat{e}_{r}}$ in spherical coordinate system ($r$,$\theta$,$\phi$), i.e., $\mathbf{n_{m}}=\mathbf{M}/|\mathbf{M}|=\mathbf{\hat{e}_{r}}$. The other two unit vectors in the spherical coordinate system are $\mathbf{\hat{e}_{\theta}}$ and $\mathbf{\hat{e}_{\phi}}$ for $\theta$- and $\phi$-rotations, respectively. The torque $\mathbf{T_{I}}$ acting on the magnetization due to interface anisotropy can be derived from the gradient of the energy [see Equation (26)] and is given by

[TABLE]

Note that the torque $\mathbf{T_{I}}$ acts along the out-of-plane direction, so that the magnetization can deflect out of magnet’s plane ($\phi=\pm 90^{\circ}$).

The effect of random thermal fluctuations is incorporated via a random magnetic field and the corresponding torque $\mathbf{T_{TH}}$ is given by the Equation (12) (also see the Equations (13) and (14) for the expressions of $P_{\theta}$ and $P_{\phi}$).

The magnetization dynamics of of the nanomagnet $M_{1}$ under the action of the torques $\mathbf{T_{I}}$ and $\mathbf{T_{TH}}$ is described by the stochastic Landau-Lifshitz-Gilbert (LLG) equation as follows.

[TABLE]

Solving the above equation analytically, we get the following coupled equations of magnetization dynamics for $\theta$ and $\phi$:

[TABLE]

[TABLE]

We solve the above two coupled equations numerically to track the trajectory of magnetization over time, in the presence of room-temperature thermal fluctuations.[87]

From Equations (29) and (30), we see that the torque acting in the $\phi$-direction is much more stronger than the torque exerted in the $\theta$-direction since the damping parameter $\alpha\ll 1$. Although, the nanomagnet has a small thickness (i.e., $l\ll b<a$ and $N_{d-xx}\gg N_{d-yy}>N_{d-zz}$, where $N_{d-zz}$, $N_{d-yy}$, and $N_{d-xx}$ are the components of the demagnetization factor along the $z$-axis, $y$-axis, and $x$-axis, respectively), magnetization cannot remain on the magnet’s plane ($y$-$z$ plane, $\phi=\pm 90^{\circ}$) due to the fact that the interface coupling energy is a few orders of magnitude higher than the shape anisotropy energy. Thus, the magnetization keeps rotating in the $\phi$-direction, but it also traverses towards the anti-parallel direction in $\theta$-space ($\theta\simeq 180^{\circ}$ to $\theta\simeq 0^{\circ}$ or vice-versa) due to damping [see Equation (29)].

Note that exactly at $\theta=180^{\circ}$ or $0^{\circ}$, the torque acting on the magnetization due to interface anisotropy [Equation (27)] is exactly zero, however, as described earlier in the Subsection 2.1, thermal fluctuations can scuttle the magnetization from these points to initiate switching. At the very start of switching, the initial orientation of magnetization is not a fixed value rather a distribution due to thermal agitations (as described in the Subsection 2.1 too). Such distribution is considered during simulations. Thermal fluctuations affect the magnetization switching during the course of switching too.

Energy Dissipation. Due to the application of voltage, we have ‘$CV^{2}$’ energy dissipation, where $C$ and $V$ denote the capacitance of the ferroelectric layer and the applied voltage, respectively. We calculate it following the procedure described in Ref. \refciteroy11_2. The energy dissipated in the nanomagnet due to Gilbert damping can be expressed as $E_{d}=\int_{0}^{\tau}P_{d}(t)dt$, where $\tau$ is the switching delay and $P_{d}(t)$ is the power dissipated at time $t$ given by

[TABLE]

Thermal field with mean zero does not cause any net energy dissipation but it causes variability in the energy dissipation $E_{d}$ by scuttling the trajectory of magnetization.

3 Results and Discussions

Here, we will review the simulation results for strain-mediated multiferroic composites,[42, 53] and interface and exchange coupled multiferroic heterostructures[87] in subsequent subsections. The coupled sets of equations derived in the previous section are numerically solved and the performance metrics e.g., switching delay, energy dissipation during switching are reported.

3.1 Strain-mediated multiferroic composites

Here we present the simulation results for strain mediated multiferroic composites.[53, 42] The magnetostrictive layer is made of polycrystalline Terfenol-D and it has the following material properties – Young’s modulus (Y): 8$\times$1010 Pa, magnetostrictive coefficient ($(3/2)\lambda_{s}$): +90$\times$10*-5*, saturation magnetization ($M_{s}$): 8$\times$105 A/m, and Gilbert’s damping constant ($\alpha$): 0.1 (Refs. \refciteRefWorks:179,RefWorks:176,RefWorks:178,materials). We choose the dimension of the magnetostrictive layer as 100 nm $\times$ 90 nm $\times$ 6 nm, which ensures that the magnet has a single ferromagnetic domain.[42, 92] The tradeoffs between area, switching delay, and energy dissipation have been comprehensively studied for different dimensions of the nanomagnet in Ref. \refciteroy14_5.

For the piezoelectric layer, we consider lead-zirconate-titanate (PZT), which has a dielectric constant of 1000.[42, 33] Since we want any strain generated in the PZT layer is transferred almost completely to the magnetostrictive layer, we assume that the PZT layer is four times thicker than the magnetostrictive layer.[42, 33] We consider that the maximum strain of 500 ppm can be generated in the PZT layer,[100, 101] which would require a voltage of 66.7 mV because $d_{31}$=1.8$\times$10*-10* m/V for PZT.[102] Assuming this strain is transferred completely to the magnetostrictive layer, the corresponding stress in Terfenol-D is the product of the generated strain ($500\times 10^{-6}$) and the Young’s modulus (8$\times$1010 Pa). Therefore, a maximum stress of 40 MPa can be generated in the Terfenol-D nanomagnet. While avoiding considerable degradation of the piezoelectric constants at low-thickness (24 nm as assumed in Ref. \refciteroy11_6) piezoelectric layers is under research,[51, 52] we can use a higher thickness for the piezoelectric layer (e.g., 100 nm) and a concomitant amount of higher voltage, since the energy dissipation due to applying voltage is miniscule. For the piezoelectric layer, if we use materials with high piezoelectric coefficients e.g., PMN-PT, which has $d_{31}$=–3000 pm/V, and $d_{32}$=1000 pm/V (note the anisotropic strain, i.e., the signs of $d_{31}$ and $d_{32}$ are different),[48] we can lower the operating voltage further to a few millivolts.

It should be noted that for the magnetostrictive layer, we need to choose a material that maximizes the product $(3/2)\lambda_{s}\,Y$. Terfenol-D (TbDyFe), which has 30 times higher magnetostriction coefficient in magnitude than the common ferromagnetic materials (e.g., Fe, Co, Ni), has the highest $(3/2)\lambda_{s}\,Y$ (Ref. \refciteroy11_2). If it needs to avoid the rare-earth materials (e.g., Tb and Dy in Terfenol-D), we can also utilize Galfenol (FeGa),[61, 103] which has 6 times less $(3/2)\lambda_{s}$, but twice high $Y$ than that of Terfenol-D.

We consider that magnetization initially situates around $\theta\simeq 180^{\circ}$ and it fluctuates due to random thermal agitations. When a compressive stress is applied to initiate switching, magnetization starts with a certain ($\theta_{initial}$,$\phi_{initial}$) picked from the initial angle distributions. The voltage generated stress is assumed to be ramped up linearly and the stress is kept constant until the magnetization reaches the plane defined by the in-plane and the out-of-plane hard axis (i.e., the $x$-$y$ plane, $\theta=90^{\circ}$). Note that when magnetization reaches at $\theta=90^{\circ}$, the azimuthal angle $\phi$ may not correspond to the magnet’s plane ($y$-$z$ plane, $\phi=\pm 90^{\circ}$) and this has important consequences in switching the magnetization in the correct direction.[53] In any case, the $x$-$y$ plane ($\theta=90^{\circ}$) is always reached sooner or later with thermal fluctuations generating a distribution of time that magnetization takes to reach at $\theta=90^{\circ}$.

After the magnetization reaches the $x$-$y$ plane, the stress is ramped down at the same rate at which it was ramped up, and reversed in magnitude to aid switching. The magnetization dynamics ensures that $\theta$ continues to rotate from $\theta=90^{\circ}$ towards $0^{\circ}$.[53] When $\theta$ becomes $\leq 5^{\circ}$, switching is deemed to have completed. We perform a moderately large number (10,000) of simulations, with their corresponding ($\theta_{initial}$,$\phi_{initial}$) picked from the initial angle distributions, for each value of stress and ramp duration to generate the simulation results. Note that thermal fluctuations affect the magnetization dynamics during the course of switching as well.

Figure 6 shows the distributions of initial angles $\theta_{initial}$ and $\phi_{initial}$ in the presence of thermal fluctuations. No bias magnetic field is applied along the out-of-plane direction (+$x$-axis) here. Note that the peak of the distribution of $\theta_{initial}$ is exactly at $\theta=180^{\circ}$. This means that magnetization is most likely to start from the easy axis $\theta=180^{\circ}$. If magnetization starts near from the easy axis, the torque acting on the magnetization would be vanishingly small. Therefore, only random thermal fluctuations can help magnetization going away from the easy axis and then magnetization can start switching. We have analyzed earlier that a bias field can shift the peak of the distribution $\theta_{initial}$ away from the easy axis and facilitate switching, which we will present next.

Figure 7 shows the distributions of initial angles $\theta_{initial}$ and $\phi_{initial}$ in the presence of thermal fluctuations when a bias magnetic field applied along the out-of-plane direction (+$x$-axis). In Fig. 7(a), note that the bias field has shifted the peak of $\theta_{initial}$ from the easy axis ($\theta=180^{\circ}$) as we desire. The $\phi_{initial}$ distribution in Fig. 7(b) spans mostly within the interval [–90∘,+90∘] since the bias magnetic field is applied in the +$x$-direction. However, the $\phi_{initial}$ distribution is asymmetric in the quadrants (0∘,90∘) and (270∘,360∘), which can be explained as follows. The magnetization is fluctuating around $\theta\simeq 180^{\circ}$ and the precessional motion of the magnetization $[-|\gamma|/\left(1+\alpha^{2}\right)\mathbf{M}\times\mathbf{H}$, where $\mathbf{M}$ is the magnetization and $\mathbf{H}$ is the effective field] due to the +$x$-directed magnetic field is such that the magnetization prefers the $\phi$-quadrant (0∘,90∘) over the $\phi$-quadrant (270∘,360∘). Hence, when the magnetization starts from $\theta\simeq 180^{\circ}$, the initial azimuthal angle $\phi_{initial}$ is more likely to be in the quadrant (0∘,90∘) than the other quadrant (270∘,360∘) in the +$x$-direction.

Figure 8 plots the magnetization dynamics for different values of $\phi_{initial}$ (while keeping fixed values of $\theta_{initial}$, applied stress, and ramp rate) to signify the role of out-of-plane excursion of magnetization. For Figs. 8(a) and 8(b), the magnetization initially lies on the plane of the magnet ($\phi_{initial}=\pm 90^{\circ}$) and the precessional motion of magnetization due to applied stress is in the +$\hat{\mathbf{e}}_{\phi}$ direction, which increases $\phi$ with time. So the magnetization starts out in the “good” quadrants (90∘,180∘) and (270∘,360∘) for Figs. 8(a) and 8(b), respectively. Therefore, both the motions of magnetization [the damped motion due to applied stress and the motion due to out-of-plane excursion shown as $-2\alpha B(\phi,t)sin\theta cos\theta\,\hat{\mathbf{e}}_{\theta}$ and $-|B_{shape,\phi}(\phi)|sin\theta\,\hat{\mathbf{e}}_{\theta}$, respectively in the Fig. 1(c)] are in the $-\hat{\mathbf{e}}_{\theta}$ direction so that $\theta$ decreases with time and the magnetization rotates in the correct direction towards $\theta=90^{\circ}$. The increasing out-of-plane excursion of the magnetization due to $\phi$ with time is eventually opposed by the damped motion due to out-of-plane excursion [depicted as $-\alpha|B_{shape,\phi}(\phi)|\,\hat{\mathbf{e}}_{\phi}$ in Fig. 1(c)], which attempts to bring the magnetization back to the magnet’s plane ($\phi_{initial}=\pm 90^{\circ}$). These two effects balance each other and $\phi$ assumes a stable value in the “good” quadrant which can be clearly observed in the plots (the flat regions in the $\phi$-plots). When $\theta$ reaches $90^{\circ}$, the torque due to stress and shape anisotropy vanishes, however, $\phi$ remains in the respective “good” quadrant for the cases in Figs. 8(a) and 8(b). At this point, stress is reversed with the same ramp rate and the damped motion due to stress and shape anisotropy eventually becomes again in the $-\hat{\mathbf{e}}_{\theta}$ direction. In this way, magnetization continues to rotate in the right direction towards $\theta=0^{\circ}$. Slightly past 0.4 ns, the precessional motion due to stress and shape anisotropy continues to rotate $\phi$ and pushes it into a neighboring “bad” quadrant, but eventually it escapes into the other “good” quadrant. This brief excursion into a “bad” quadrant causes the ripple in $\theta$-plots. However, magnetization ends up switching successfully.

In Figs. 8(c) and 8(d), the magnetization is initially lifted far out of the magnet’s plane ($\phi_{initial}=0^{\circ},180^{\circ}$, respectively), where the huge out-of-plane shape anisotropy energy barrier cannot be overcome by the stress anisotropy and $B(\phi,t)$ remains positive, i.e., $|B_{stress}(t)|<|B_{shape}(\phi)|$. This forces magnetization to precess in the clockwise direction ($-\hat{\mathbf{e}}_{\phi}$) rather than in the anticlockwise direction ($+\hat{\mathbf{e}}_{\phi}$). Therefore $\phi$ decreases with time, which takes magnetization inside the neighboring “good” quadrant and eventually $|B_{stress}(t)|$ becomes greater than $|B_{shape}(\phi)|$. Then $\phi$ assumes a stable value due to the counteraction between the damped motion and the motion due to out-of-plane excursion [depicted as $-\alpha|B_{shape,\phi}(\phi)|\,\hat{\mathbf{e}}_{\phi}$ and $2B(\phi,t)cos\theta\,\hat{\mathbf{e}}_{\phi}$ in Fig. 1(c), respectively]. Thereafter, switching occurs similarly for the cases as in Figs. 8(a) and 8(b). Slightly past 0.3 ns, continuing $\phi$ rotation pushes $\phi$ into a neighboring “bad” quadrant, but eventually it escapes into the other “good” quadrant. Once again, this brief excursion into the “bad” quadrant causes the ripple in $\theta$-plots and eventually successful switching takes place in the end.

So straying into a “bad” quadrant for azimuthal angle $\phi$ during the ramp-down phase does not mean at all that magnetization would fail to switch. Magnetization can rotate in the other “good” quadrant for $\phi$ and complete the switching. Thus, there may be ripples appearing in the magnetization dynamics at the end of switching increasing the switching time, which is because of the transition of azimuthal angle $\phi$ between two “good” quadrants through one “bad” quadrant. This happens particularly near the end of switching, when the precessional motion is strong.

An interesting comparison between the switching delays in the Figs. 8(a) or 8(b), and the switching delays in the Figs. 8(c) or 8(d) reveals that the switching delay decreases by 0.1 ns when magnetization starts from out-of-plane ($\phi_{initial}=0^{\circ},180^{\circ}$) compared to when magnetization starts from in-plane ($\phi_{initial}=90^{\circ},270^{\circ}$). This is a very consequence of the reasoning that out-of-plane excursion of magnetization in the “good” quadrants aids magnetization to move faster in $\theta$-space. When magnetization starts out-of-plane, magnetization spends more time deep inside a “good” quadrant; hence, switching gets faster than that of the case when magnetization starts in-plane of the magnet.

Figure 9 demonstrates a couple of examples when switching fails to occur, i.e., magnetization backtracks to its original position. In Fig. 9(a), when the polar angle $\theta$ reaches $90^{\circ}$, the azimuthal angle $\phi$ has ventured into the “bad” quadrant ($0^{\circ}$,90∘) due to thermal fluctuations. Thus, switching eventually fails. In Fig. 9(b), when the polar angle $\theta$ reaches $90^{\circ}$, the azimuthal angle $\phi$ is greater than $90^{\circ}$ and in the “good” quadrant ($90^{\circ}$,180∘). However, after reaching around $\theta\simeq 50^{\circ}$, the magnetization backtracks to the initial state and therefore switching fails to occur. This happens because of the long ramp duration which forces $\phi$ to decrease over time and eventually $\phi$ enters into the neighboring “bad” quadrant ($0^{\circ}$,90∘). Also during the passage of long duration of ramp, thermal fluctuations have ample opportunity to scuttle the switching. Such switching failure is depicted intuitively in the Fig. 5 earlier.

Figure 10 plots the percentage successful switching rates at room-temperature (300 K) when the magnetostrictive nanomagnet is subjected to stress between 10 MPa and 30 MPa with ramp duration (60 ps, 90 ps, and 120 ps) as a parameter. For each value of stress and ramp duration, a moderately large number (10,000) of simulations were performed to generate these results. Initial angle distributions at 300 K for both $\theta$ and $\phi$ are taken into account during simulations (see Fig. 7). The minimum stress needed to switch the magnetization without considering thermal fluctuations is $\sim$5 MPa, but at room temperature (300 K) this minimum stress is increased to $\sim$14 MPa for 60 ps ramp duration and to $\sim$17 MPa for 90 ps ramp duration. The minimum stress of 5 MPa without considering thermal fluctuations is attributed to the removal of in-plane shape anisotropy energy barrier by the stress anisotropy, while at 300 K an increased magnitude of stress is required since magnetization is scuttled in “bad” quadrants due to thermal fluctuations. When longer ramp duration, a higher stress is required to prevent magnetization traversing into “bad” quadrants. Therefore, it is beneficial to reduce the ramp duration (i.e., having a faster ramp rate) to increase the success rate of switching at a lower stress level. Simulation results show that with 1 ps ramp duration, the critical stress can be reduced by $\sim$2 MPa compared to the case of 60 ps ramp duration.

For 120 ps ramp duration, $\sim$100% success probability is unattainable for any value of stress since thermal agitations have higher latitude to scuttle the magnetization into “bad” quadrants while stress is ramped down. At higher stresses accompanied by a long ramp duration, there occurs higher out-of-plane excursion pushing the magnetization into “bad” quadrants, which further aggravates the error probability. At very long ramp duration, the success (and error) probability becomes 50%, since the magnetization would stay in-plane of the magnet and during the ramp-down phase, random thermal fluctuations may equally scuttle the magnetization either in the “good” quadrants or in the ‘bad” quadrants.

Figure 11(a) shows the distribution of time taken for magnetization polar angle $\theta$ to reach $90^{\circ}$ ($x$-$y$ plane). This wide distribution is caused by: (1) the initial angle distributions in Fig. 7, and (2) thermal fluctuations during the course of transition from some $\theta=\theta_{initial}$ to $90^{\circ}$. We do need to tackle such distribution by keeping the magnetization out-of-plane far enough so that magnetization does not collapse on magnet’s plane. We can use a sensing element to detect when $\theta$ reaches around 90*∘*, so that we can ramp down the stress thereafter. The sensing element can be implemented by measuring the magnetoresistance in a magnetic tunnel junction (MTJ)[62, 63, 64, 65, 66, 67, 68, 53, 34, 35]. We need to get calibrated on the magnetoresistance of the MTJ when magnetization resides at $\theta=90^{\circ}$ ($x$-$y$ plane). Comparing this known signal with the sensed signal of the MTJ, the stress can be ramped down. Such comparator can be implemented with these energy-efficient multiferroic devices, i.e., energy-inefficient charge-based transistors do not need to be utilized.[104, 105, 106] The energy dissipation in this sensing element is not considered here. It should be mentioned that usually it requires several peripheral circuitry in conjunction with the basic device itself, however, energy dissipation considering in the other required circuitry does not change the order of dissipation.[107, 108]

Some tolerance is nonetheless required since the sensing element cannot be perfect. Simulation results show that the internal dynamics works correctly as long as the stress is ramped down when $\theta$ is in the interval [85*∘, 110∘], i.e. it does not have to be exactly 90∘*. This tolerance is due to the motion arising from the out-of-plane excursion of magnetization in a “good” quadrant. If magnetization reaches at $\theta=90^{\circ}$ (even past it towards $\theta=0^{\circ}$) and stress is not withdrawn soon enough, then magnetization can fall on the magnet’s plane ($\phi=\pm 90^{\circ}$, potential energy minima), upon which the success probability would be 50% since thermal fluctuations can put magnetization in either direction of the potential landscape.

Figure 11(b) shows the distribution of azimuthal angle $\phi$ when $\theta$ reaches $90^{\circ}$. Note that $\phi$ always resides in the “good” quadrant [($90^{\circ},180^{\circ}$) or ($270^{\circ},360^{\circ}$)] and has a fairly narrow distribution. As depicted in the Fig. 10, a high stress and fast ramp rate are required to ensure that $\phi$ is in the “good” quadrants, which is conducive to successful switching.

Figure 12 shows the delay and energy distributions for 15 MPa stress and 60 ps ramp duration in the presence of room-temperature (300 K) thermal fluctuations. The high-delay tail in Fig. 12(a) is particularly associated with the switching trajectories that start very close to $\theta=180^{\circ}$. In such trajectories, the torque acting on the magnetization is vanishingly small, which makes the switching sluggish at the beginning. During this time, switching may also becomes susceptible to backtracking due to random thermal fluctuations, which may increas the delay further. Nonetheless, out of 10,000 simulations of switching trajectories peformed, there was not a single one in Fig. 12(a) where the delay exceeded 1 ns meaning that the probability of such happening is less than 0.01%. Since the energy dissipation is the product of the power dissipation and the switching delay, similar behavior is found in Fig. 12(b). Note that the variation in the distribution of energy dissipation happens due to the internal energy dissipation (caused by Gilbert damping) in the presence of random thermal fluctuations since all the trajectories correspond to the same 15 MPa stress and 60 ps ramp duration.

Discussion on magnetization switching. We will now analyze the magnetization switching between the $180^{\circ}$ symmetry equivalent states based on the simulation results as depicted in the Fig. 13. The usual perception is that stress can rotate magnetization of a magnetostrictive nanomagnet only by $90^{\circ}$ from $\pm z$-axis to $\pm y$-axis (see Fig. 13). However, if we determine the expression of torque due to stress from the stress anisotropy energy [see Eq. (5)] as

[TABLE]

we see that the torque due to stress acts along the out of magnet’s plane ($\mathbf{\hat{e}_{\phi}}$ direction) and therefore magnetization lifts out-of-plane (although the out-of-plane demagnetization field due to small thickness of the nanomagnet attempts to keep the magnetization in-plane). This out-of-plane excursion of magnetization generates an intrinsic asymmetry, which can completely switch the magnetization by $180^{\circ}$ (Ref. \refciteroy13_2). Full $180^{\circ}$ switching is desirable since it facilitates having the full tunneling magnetoresistance (TMR) while electrically reading the magnetization state using a magnetic tunnel junction (MTJ).[62, 63, 64, 65, 66, 67, 68]

Such full $180^{\circ}$ switching of magnetization (i.e., memory operation) in piezoelectric-magnetostrictive heterostructures was first addressed by Roy in Ref. \refciteroy11_news. Note that Ref. \refciteRefWorks:154 shows the Bennett clocking operation[57] for logic design, which is not memory operation since using a neighboring nanomagnet to switch another nanomagnet is tantamount to using a magnetic field,[110] the direction of which has to be reversed for switching the memory bit in either direction. Also, Ref. \refciteRefWorks:154 performed an incorrect analysis by assuming that magnetization always resides in-plane, which leads to a very high switching delay and high energy-delay product compared to the traditional transistors.[111] This was addressed and explained by Roy.[7, 112, 113]

Figure 13 illustrates the magnetization switching alongwith the simulation results of different distributions in the presence of room-temperature thermal fluctuations. When magnetization starts switching from $\theta=180^{\circ}$, the initial orientation of magnetization is a distribution due to random thermal fluctuations. When magnetization reaches $\theta=90^{\circ}$, the time taken to reach there is also a distribution due to thermal fluctuations. After magnetization reaches around $\theta=90^{\circ}$, stress needs to be released to switch the magnetization. However, there can be two states magnetization can end up: (1) $\theta=180^{\circ}$ causing a switching failure, and (2) $\theta=0^{\circ}$ leading to successful switching. As shown in the Fig. 13, magnetization’s excursion out of magnet’s plane in a particular direction would lead to successful switching. As presented in the Fig. 10, a sufficiently fast ramp rate of stress accompanied by a sufficiently high stress is conducive to successful switching, i.e., the out-of-plane excursion of magnetization in the specified directions as shown in the Fig. 13 would be maintained. Therefore, the out-of-plane excursion of magnetization generated an equivalent intrinsic asymmetry for magnetization switching in the correct direction. If magnetization starts switching from $\theta=0^{\circ}$ instead of $\theta=180^{\circ}$, the direction of out-of-plane excursions of magnetization as shown in the Fig. 13 would be exactly opposite.

Note that the aforesaid magnetization switching as presented in the Fig. 13 requires to read the magnetization state using spin valve/MTJ to sense when magnetization reaches around $\theta=90^{\circ}$ (since room temperature thermal fluctuations make the traversal time a wide distribution as shown) so that stress can be brought down around that time. This has been discussed in the context of explaining the results in the Fig. 11. The sensing element reads the magnetization state and a comparator can compare the read signal with a pre-calibrated value for the MTJ resistance when magnetization’s orientation is around $\theta=90^{\circ}$ (as explained in Ref. \refciteroy13_2). It is important to note that there is tolerance around $\theta=90^{\circ}$, i.e., stress does not need to be withdrawn exactly at $\theta=90^{\circ}$ since the sensing procedure cannot be perfect.[53] It is shown that the internal magnetization dynamics provides such tolerance.[53] Any additional element for comparison can be built using these energy-efficient multiferroic devices in general.[104, 105, 106] Researchers are trying to replace the traditional switch based on charge-based transistors by a new possible “ultra-low-energy” switch (e.g., using multiferroic composites). Therefore, any circuitry can be built with the energy-efficient switch itself rather than the conventional transistors.[104, 105, 106] Usually, it requires several peripheral circuitry in conjunction with the basic switch in a system.[107, 108] While researchers report on the performance metrics of the basic switch itself, the total energy dissipation considering the other required circuitry does not change the order of energy dissipation, utilizing the respective devices.[107, 108] Therefore energy efficient technologies can be envisaged using such magnetization switching mechanism[7, 33, 93, 53, 42] and computing methodologies[7, 34, 35] based on such switching mechanism.

Also such aforesaid switching is confirmed by others,[114, 115, 116, 117, 118, 119] but they do not perform a detailed analysis using stochastic Landau-Lifshitz-Gilbert equation of magnetization dynamics and therefore do not conceive or mention the requirement of sensing circuitry due to the wide distribution of traversal times in the presence of room temperature thermal fluctuations as explained earlier. While it is desirable that an additional sensing procedure can be avoided, the methodology[53] tackles the variation in the traversal time of magnetization particularly the crucial initial distribution of magnetization effectively for successful switching and it also aids in decreasing the switching delay given a certain error probability. An out-of-plane field (see Ref. \refciteroy11_6) can help breaking the symmetry while switching the magnetization but the variation in traversal times towards the hard axis due to the initial Boltzmann distribution of magnetization is large enough to cause switching failures (note that it requires very low switching error probability e.g., $<10^{-4}$ to be technologically viable). While further research may lower the error probability with a symmetry breaking field or possibly by other means, the critical understandings behind the switching dynamics of magnetization as explained by Roy in Ref. \refciteroy13_2 is an important step forward to the technological applications for building memory and logic using multiferroics.

In Ref. \refciteroy14, information propagation in a chain of nanomagnets using Bennett clocking mechanism[57] in the presence of room-temperature thermal fluctuations is analyzed. It is found that the inherent magnetization dynamics caused by the out-of-plane excursion of magnetization can lead to switching failure even if there is dipole coupling between neighboring nanomagnets that can introduce asymmetry in switching. But such asymmetry is not sufficient and therefore the asymmetry caused by the out-of-plane excursion[53] can be utilized for successful switching of magnetization. Once again, here a sensing element is required to detect when magnetization reaches around $\theta=90^{\circ}$ and multiferroic devices can be utilized to build such functionality, which would not change the order of energy dissipation.

But, one contradictory viewpoint exists that any other circuitry other than the switching device itself needs to be built with traditional transistors,[120] which however undermines the research behind finding a switch replacing the transistor itself[41, 111] and therefore incorrect lacking the understandings over the development of transistor based circuitry.[107, 108] Such viewpoint does not exist anywhere else in literature. A few contradicting facts vis-a-vis the comment made by Bandyopadhyay and Atulasimha (referred as BA onwards) in Refs. \refciteroy13_comment,biswa14_3,munir15 are as follows:[111] (1) BA are coauthors of Roy in Refs. \refciteroy11_news,roy11_2,roy13_2,roy11_5,roy11_6, where energy-efficiency is claimed, requiring the sensing procedure too. In particular, energy efficiency is claimed in the presence of thermal fluctuations in Ref. \refciteroy11_6, which requires the sensing element, explained in details by Roy in Ref. \refciteroy13_2. Also, there is a patent[123] following the research conceived by Roy[33, 93, 54, 42, 53, 7, 112] that claims energy-efficiency, requiring the sensing element therein as well. Therefore, the comment made by BA in Ref. \refciteroy13_comment has no technical basis. (2) Note that Ref. \refcitemunir15, coauthored by BA, uses precisely shaped pulses. Such pulses need to be generated too using some circuitry. According to BA, transistors need to be utilized to build such circuitry and the system would dissipate too much energy, invalidating the claim of energy efficiency in Ref. \refcitemunir15. Note that one additional hardware cannot be shared between many devices distributed on a chip due to interconnect delay and loading effect. Also, note that pulse shaping is an ineffective countermeasure since it is not helping much in reducing the error probability, therefore building and using such circuitry do not make sense either. (3) Ref. \refcitebiswa14_2, in which BA are coauthors, proposed a “toggle” switch (as stated that “a write cycle must be preceded by a read cycle to determine the stored bit”), which would require a similar use of spin-valve or MTJ for reading the known bit, storing it, and then using it for comparison. According to BA, such additional circuitry needs to be constructed with energy-inefficient transistors, invalidating the claim of energy efficiency in Ref. \refcitebiswa14_2.

According to the unsubstantiated and incorrect viewpoint as stated by BA,[120] the switching methodology confirmed by others[114, 115, 116, 117, 118, 119] will require the sensing circuitry too (as already explained by Roy in Ref. \refciteroy13_2) and hence will be energy-inefficient. However, as explained earlier, any circuitry can be built with the energy-efficient switch itself rather than the conventional transistors[104, 105, 106] and furthermore BA contradict themselves as explained above. Therefore, there is no technical reasoning behind such point raised by BA in Ref. \refciteroy13_comment. The other unsubstantiated points raised by BA in Ref. \refciteroy13_comment on Ref. \refciteroy13 are also technically incorrect as explained.[111] Another paper[125] also incorrectly stated that stress has to “to drive the magnetization switch out-of-plane at first”, however, as pointed out by the Eq. (32) that the out-of-plane excursion is inherent to the magnetization dynamics.[53] Such out-of-plane motion also increases the switching speed to the order of GHz.[54, 53] Furthermore, Ref. \refciteRefWorks:888 did not study the consequence of room-temperature thermal fluctuations, which is critical to the magnetization dynamics and error probability of switching.

One of the critical points to conceive while proposing energy-efficient systems is that it must be area-efficient as well to compete with the traditional transistors e.g., our laptops cannot be 10 times bigger. However, such area-inefficient devices have been proposed in Ref. \refcitebiswa14, which claims a superior design of magnetoelastic memory, compared to an earlier idea.[127, 43, 128, 129] Comparatively, the lateral area chosen by Ref. \refcitebiswa14 is an order higher than that of chosen by Ref. \refciteRefWorks:850. Also, Ref. \refcitebiswa14 chooses to use two pairs of lateral electrodes, used to apply stress at angles, understanding behind which is otherwise known, e.g., Ref. \refcitedsouza11, which has however issues with room-temperature thermal fluctuations. The additional lateral pads cannot be dispensed and particularly they need to be large for application of stress, thereby consuming additional large area for each nanomagnet apart from the area consumption by the nanomagnet itself. The area consumed by the proposal in Ref. \refcitebiswa14 (also Refs. \refcitebiswa14_2,biswa14_3) is so high that the devices become of micro-scale size, therefore, the proposals are untenable for building nanoelectronics.[131, 41] There is an ongoing drive to reduce the area-consumption,[132, 87] but the proposals in Refs. \refcitebiswa14,biswa14_2,biswa14_3 do not bode at all on such crucial front. There is also similar issue regarding area inefficiency in the proposal by Ref. \refciteRefWorks:888 due to consuming large area, which is a crucial issue and therefore untenable for meeting practical standard requirement of area density 1 Tb/in2.

Furthermore, Ref. \refcitebiswa14 performs an incorrect analysis while comparing error probability and switching delay with that of the Refs. \refciteRefWorks:850,RefWorks:559,RefWorks:848,RefWorks:849, and underestimates the MTJ resistance ratio with an incorrect statement “The maximum value of this ratio (assuming $\eta_{1}=\eta_{2}=1$) is 2:1” (for details see Ref. \refcitearXiv:1501.05941v2). Also, there is a very basic issue behind the multi-step switching methodology proposed in Refs. \refcitebiswa14,biswa14_2, which increases the switching delay exponentially while meeting a lowered error probability.[111] It should be noted too that Refs. \refcitebiswa14,biswa14_2 have assumed instantaneous ramp for stress, which is unreasonable and crucial, since ramp that is not fast enough causes switching failures.[53, 42]

Also, switching delay is a major performance metric while proposing energy-efficient systems, since devices with lower switching speed would take more time given a computation task to be performed, therefore making the energy-delay product a tenable performance metric. If we calculate the switching delay of magnetization according to Ref. \refciteRefWorks:154, it will come out $\sim$1000 ns, which is clearly untenable for building nanoelectronics,[131, 41] While comparing with spin-transfer-torque (STT) switching mechanism, if one compares an usual performance metric i.e., switching delay-energy, clearly Ref. \refciteRefWorks:154 performs inferior to STT switching. Also consider the issue that 1 hour of execution would take 100 hours or more using the operation presented in Ref. \refciteRefWorks:154, which however performed an incorrect analysis.[111] If charge based transistors were to operate slow, the energy dissipation would not have been an issue.[107, 131, 41]

In Ref. \refcitemunir15, the authors proposed to reduce error probability of switching in a system of two dipole-coupled magnetostrictive nanomagnets in strain-mediated multiferroic heterostructures using voltage (stress) pulse shaping. The authors conclude that high success probability of switching cannot be achieved at high switching speed ($\sim$1 ns), and therefore their proposed system is only applicable for niche applications (but therein too it is inferior to the traditional transistors). However, such conclusion lacked the critical understandings behind such high error probability for general-purpose logic applications. Fortunately, such analysis and a possible solution along the line of the analysis is present in literature[35] using Bennett clocking[57] for logic design, on which Ref. \refcitemunir15 (and arXiv version Ref. \refcitemunir14) has made some incorrect statements. Earlier, a subset of the authors of Ref. \refcitemunir15 also published a paper[134] on a four magnet system using the Bennett clocking mechanism, where the authors also came up with a similar conclusion of high error probability and the demise of multiferroic nanomagnetic logic, however, without relevant analysis similar to the case as in the Ref. \refcitemunir15. There are several technical issues in Ref. \refcitemunir15 as well.[111] Due to the incorrect analysis performed in the Refs. \refcitefasha13,munir15, in Ref. \refciteRefWorks:894, the authors cast doubt on that front, which is incorrect too.

Rather than using Bennett clocking mechanism[57] for logic design, Ref. \refciteroy13 proposed universal logic gates (NAND/NOR) utilizing a single device with multiple contacts on the device to add up the strains generated in the piezoelectric layer,[136] and a Set input to preset the non-volatile magnetization state and facilitate concatenation. Ref. \refciteroy13 uses the switching methodology presented in Ref. \refciteroy13_2 and it is clearly mentioned while referring the Ref. \refciteroy13_2 that “Computing methodologies utilizing such 180∘ switching mechanism between the two stable states of a shape-anisotropic magnetostrictive nanomagnet have not been proposed so far”. However, BA in Ref. \refciteroy13_comment casts doubt on that, which has no technical basis as explained earlier. Moreover, BA in Ref. \refcitebiswa14_3 follow the same idea presented in Ref. \refciteroy13 by Roy and uses a $\sim 90^{\circ}$ (precisely $86.4^{\circ}$) switching mechanism following an earlier idea by others.[127, 43, 128, 129] Therefore, magnetization does not switch a complete $180^{\circ}$ in Ref. \refcitebiswa14_3 and this leads to a low tunneling magnetoresistance (TMR) while reading the magnetization state and it has serious consequence on read error probability. Ref. \refcitebiswa14_3 states that it requires a sensing circuitry for operation of the proposal in the Ref. \refciteroy13 and indicates as if it is an issue with the computing proposal in Ref. \refciteroy13. First, Ref. \refcitebiswa14_3 does not point out at all that the sensing element is required for complete $180^{\circ}$ switching only and it does not require so if $90^{\circ}$ switching mechanism is used.

Ref. \refcitebiswa14_3 just uses resistors and potential divider (see Fig. 1 in Ref. \refcitebiswa14_3) to accommodate multi-inputs (rather than using intrinsic strain-addition property of piezoelectrics, while following the central idea in Ref. \refciteroy13) and the authors want to use traditional transistors for any other circuit element, which makes the proposal in Ref. \refcitebiswa14_3 energy-inefficient.[111] (Note that a circuit element cannot be shared between many devices distributed on a chip due to interconnect delay and loading effect.) It is stated in Ref. \refcitebiswa14_3 (in the supplementary material) that “The dissipation in the resistance $R$ can be negligible as we can make this resistance arbitrarily high.” This is incorrect (any standard electrical engineering undergraduate textbook can be consulted) since $RC$ delay will be too high. Therefore, the design proposed in Ref. \refcitebiswa14_3 is untenable. There are several technical issues in Ref. \refcitebiswa14_3 as well.[111]

It should be noted that the switching as depicted in the Fig. 13 just toggles the magnetization state upon application and removal of stress. Therefore, it requires to read the magnetization state before writing a bit. If we want to switch the magnetization along a specified direction, it is possible to use spin-transfer-torque with spin-polarized current which was first proposed in Ref. \refciteroy10 by Roy, however, it incurs much higher energy dissipation due to utilizing current-induced spin-torque mechanism. Unless new strategies to lower the energy dissipation in spin-transfer-torque switching can be achieved, a target of 1 aJ energy dissipation per switching a bit of information cannot be reached.

Ref. \refcitebiswa13 follows the idea presented in the Ref. \refciteroy10 with the addition of surface acoustic wave (SAW) and also follows the formulation from Refs. \refciteroy11_6,roy11_3. The basic idea was to switch the magnetization by STT along the desired direction when magnetization comes at the hard-plane ($\theta=90^{\circ}$) upon application of stress.[137] However, thermal fluctuations cause a wide distribution for the time taken by magnetization to reach at $\theta=90^{\circ}$ as explained in Ref. \refciteroy13_2. This is why STT current needs be kept active for almost half of the duration of switching.[138] Figure 3 in Ref. \refcitebiswa13 shows that the energy dissipation with a stress of 6.1 MPa is $5\times 10^{9}$ kT at room-temperature (300 K), which is $\sim$20 pJ. This is 4-5 orders of magnitude higher than that of traditional transistors and therefore untenable for building nanoelectronics.[131, 41] Also, Fig. 3 in Ref. \refcitebiswa13 shows the energy dissipation when no stress is present (only switching in STTRAM), which is $\sim 40\times 10^{9}$ kT = 160 pJ. Note that Ref. \refcitebiswa13 has considered the material Terfenol-D as free layer in STT switching to calculate the energy dissipation. However, the material that is commonly used for the free layer is CoFeB,[139, 132, 67, 140] which has Gilbert damping parameter $\alpha$ an order lower[139, 132] than Terfenol-D ($\alpha$ of Terfenol-D is 0.1 used by Ref. \refcitebiswa13). The critical current of switching is proportional to the damping parameter[139, 132] and the energy dissipation is proportional to the square of the switching current. Therefore, Ref. \refcitebiswa13 calculated incorrectly the switching current an order higher (23 mA) and energy dissipation about two orders higher ($\sim$160 pJ) for STT switching (Ref. \refciteRefWorks:786 correctly determined switching current $\sim$1 mA and energy dissipation $\sim$1 pJ experimentally). Clearly, the comparison with STT switching performed in Ref. \refcitebiswa13 is incorrect and actually the energy dissipation in STT switching is 12.5 times lower (in stead of 8 times higher as incorrectly claimed by Ref. \refcitebiswa13) than the hybrid scheme proposed in Ref. \refcitebiswa13.

Also, the energy dissipation calculated in Ref. \refcitebiswa13 due to SAW is missing a crucial point that will make the energy dissipation exceedingly high. Ref. \refcitebiswa13 says that “SAW is global and affects every memory cell.” This is represented in Fig. 14(a) that with the application of 6.1 MPa stress the barrier height decreases but it is not sufficient enough to make the potential landscape monostable and cause switching. However, when stress is applied the magnetization does rotate. The distribution in Fig. 14(d) upon application of stress is wider than the distribution in Fig. 14(c) when no stress is active. Therefore, the energy dissipation upon application of stress and removal of stress must be calculated, which is missed by Ref. \refcitebiswa13. This energy dissipation turns out to be $\sim$40 kT for one cell. Considering just 1 MB memory, clearly such strategy would dissipate an energy which is quite worse than that of transistors.[131, 41]

Ref. \refciteRefWorks:847 embraces the same concept of using STT (as proposed in Ref. \refciteroy10) and does not consider the sensing procedure (explained in Ref. \refciteroy13_2), which is required for successful switching in strain-mediated multiferroic composites in the presence of thermal fluctuations without STT, as explained earlier. This is why the Ref. \refciteRefWorks:847 incorrectly terms Refs. \refciteroy11_2,roy13_2 as stochastically unstable. Since, using STT leads to energy-inefficiency[137, 138, 141], symmetry breaking by other means rather than using STT is desirable.

In Ref. \refciteroy14_2, it is shown that interface and exchange coupling in multiferroic heterostructures can facilitate switching of magnetization in a particular direction without using spin-torque mechanism, while incurring miniscule energy dissipation and there is no necessity of a sensing element as required for strain-mediated multiferroic composites.[53] The switching methodology presented in the Ref. \refciteroy14_2 can be harnessed for logic design and computing purposes too.[34, 35] In the next subsection, we will present the simulation results for magnetization switching in such interface and exchange coupled multiferroic heterostructures.

Figure 15(a) shows the thermally averaged switching delay versus stress with different ramp duration (60 ps, 90 ps, 120 ps) as a parameter. Only the successful switching events are considered here since the switching delay metric does not make sense for an unsuccessful event. For a certain stress, if the ramp duration is decreased (i.e., the ramp rate is increased), the stress reaches its maximum value quicker and switches the magnetization faster, thereby decreasing the switching delay. For ramp durations of 60 ps and 90 ps, the switching delay decreases with increasing stress since an increasing stress anisotropy rotates the magnetization faster. However, for 120 ps ramp duration, the dependence of switching delay decreases with stress is non-monotonic, due to the exactly the same reasoning that caused the non-monotonicity in Fig. 10. A high stress accompanied by a long ramp duration is harmful since during the ramp-down phase it rotates the magnetization in “bad” quadrants leading to increased switching delay and even switching failures.

Figure 15(b) shows the standard deviation of switching delay distribution with stress for 60 ps ramp duration. At higher values of stress, the increased stress anisotropy energy makes the potential landscape more deep and resist the thermal fluctuations more effectively. Therefore, the spread in switching delay diminishes with a higher stress.

Figure 16(a) shows the thermal mean of the total energy dissipated to switch the magnetization as a function of stress for different ramp durations. For a certain ramp duration, the average power dissipation ($E_{total}/\tau$) increases with stress and for a certain stress it decreases with increasing ramp duration. A higher stress leads to more ‘$CV^{2}$’ dissipation (see Fig. 16(b)) and also more internal dissipation because an increasing stress anisotropy results in a higher torque. A slower switching (i.e. more adiabatic) decreases the power dissipation. But the switching delay curves in the Fig. 15(a) show the opposite trend. At a slower ramp rate (higher ramp duration), the average power dissipation $E_{total}/\tau$ is always smaller than that of a higher ramp rate. However, the switching delay does not decrease as fast as with higher values of stress (in fact switching delay may increase for higher ramp duration, see Fig. 15(a)), which is why the energy dissipation curves in Fig. 16(a) exhibit the cross-overs.

Figure 16(b) plots the ‘$CV^{2}$’ energy dissipation due to the application of voltage-induced stress for different ramp durations. Stress is proportional to the applied voltage $V$, and therefore the ‘$CV^{2}$’ energy dissipation increases with stress for a certain ramp duration. This ‘$CV^{2}$’ energy dissipation however is a small fraction of the total energy dissipation ($<$ 15%) particularly because a miniscule voltage-generated stress is required to switch the magnetization of a magnetostrictive nanomagnet in a piezoelectric-magnetostrictive multiferroic. The ‘$CV^{2}$’ dissipation decreases with an increasing ramp duration (i.e., slower ramp rate) for a certain stress since the switching becomes more adiabatic. This component of the energy dissipation would come as several orders of magnitude higher if we switch the magnetization with an external magnetic field[142] or spin-transfer-torque.[30]

Note that we would require an adiabatic circuit to take advantage of the energy saving due to using an adiabatic pulse. In any case, using a piezoelectric layer of higher piezoelectric coefficients e.g., utilizing lead magnesium niobate-lead titanate (PMN-PT) instead of lead-zirconate-titanate (PZT), the energy dissipation can be reduced further. Also, PMN-PT can generate anisotropic strain, which allows us to work with lower voltage for a required strain, thereby reducing the energy dissipation even further. PMN-PT layer has a dielectric constant of 1000, $d_{31}$=–3000 pm/V, and $d_{32}$=1000 pm/V (Ref. \refciteRefWorks:790). With the piezoelectric layer’s thickness $t_{piezo}$=24 nm (Ref. \refciteroy11_6), $V=1.9$ mVs (2.9 mVs) of voltages would generate 20 MPa (30 MPa) compressive stress [$\sigma=Y\,d_{eff}\,(V/t_{piezo})$, where $d_{eff}=(d_{31}-d_{32})/(1+\nu)$] in the magnetostrictive Terfenol-D layer, which has $Y=80$ GPa,[42] and Poisson’s ratio $\nu=0.3$ (Ref. \refciteRefWorks:821). Modeling the piezoelectric layer as a parallel plate capacitor ($\sim$100 nm lateral dimensions), the capacitance C=2.6 fF and thus $CV^{2}$ energy dissipation turns out to be $<$ 0.1 aJ. This is without considering any energy saving with the help of an adiabatic circuit. Such miniscule energy dissipation is the basis of ultra-low-energy computing using these multiferroic devices.[33, 7, 34, 35]

Figure 17 plots the switching delay versus energy dissipation for different ramp durations. This plot can be extracted from the Fig. 15(a) (stress versus switching delay) and Fig. 16(a) (stress versus energy dissipation). This plot signifies that as the switching delay increases, energy dissipation deceases. This points out the usual delay-energy trade-off. For 120 ps ramp duration, there is a opposite trend at higher stress values, the reasoning behind which has been already described while explaining the results in the Fig. 10.

3.2 Interface and exchange coupled multiferroic heterostructures

Here we present the simulation results for interface and exchange coupled multiferroic heterostructures.[87] Fig. 1 shows that the nanomagnets ($M_{1}$ and $M_{2}$ layers) and the ferroelectric $P$-layer are made of Fe and $PbTiO_{3}$, respectively, while the spacer layer is made of $Au$ and the thicknesses of the trilayer $M_{1}$/spacer/$M_{2}$ are 1/4/1 monolayers.[71, 79, 144, 145] The Fe layer has a unit cell length of 0.287 nm and it possess the following material parameters: saturation magnetization ($M_{s}$) = 1e6 A/m, and damping parameter ($\alpha$) = 0.01.[146, 147, 148] The elliptical lateral cross-section ($y$-$z$ plane, Fig. 1) of the vertical stack has a dimension of $15\,nm\times 7\,nm$. The $P$-layer has a unit cell length of 0.388 nm and it has 5 layers in the vertical direction ($x$-direction, Fig. 1)[71]. The energy difference between the P-alignment and AP-alignment is 10 meV/atom,[71] and the absolute value of energy is calculated to be about 10 eV or 385 kT at room-temperature. This huge interface coupling makes potential landscape of $M_{1}$ monostable at $\theta=180^{\circ}$ or at $\theta=0^{\circ}$ depending on the $P$-alignment or the $AP$-alignment in the trilayer, respectively. Due to the monostable energy landscape and huge energy barrier, spontaneous switching of magnetization between $\theta=180^{\circ}$ and $\theta=0^{\circ}$ cannot take place. The interface coupling energy is a few orders of magnitude higher than the shape anisotropic energy and hence consideration of shape anisotropy does not make any significant difference, however, it is included during the simulations.

If we consider the $P$-layer as a parallel plate capacitor and use a relative dielectric constant of 1000,[149] the capacitance $C$ of the layer becomes $\sim$0.4 fF. If the $P$-layer is accessed with a 10 $\mu m$ long silver wire with resistivity $\sim$2.6 $\mu\Omega$-$cm$,[150] the resistance $R$ can be calculated as $\sim$3 $k\Omega$. Therefore, the $RC$ time constant turns out to be the order of 1 ps. We assume that the ferroelectric $PbTiO_{3}$ has a coercive voltage of 20 MV/m[151] and hence a voltage $V\simeq 40$ mv is required to switch its polarization. Note that polarization switching is possible in less than 100 ps[152] and a voltage ramp with period $T=100$ ps or more is considered to enforce the quasistatic (adiabatic) assumption ($T\gg RC$). Without any adiabatic assumption, the metric $CV^{2}$ can be calculated as $0.5$ aJ and hence the energy dissipation due to the application of the voltage is miniscule. With 100 ps ramp period, the “$CV^{2}$” dissipation is determined as a negligible value of 0.01 aJ.[93, 42] Note that we do not calculate any standby leakage current through the thin ferroelectric since the device operation is non-volatile, i.e., it is possible to turn off the voltage without loosing the information in long term. However, during the active mode of operation, the leakage needs to be considered, however, the tunneling current is small ($<$ 1 nA, Ref. \refciteRefWorks:792) leading to negligible energy dissipation. Another issue that needs to be considered is ferroelectric fatigue, which may make the coercive field higher over time,[153, 154] i.e., it would require a higher voltage to switch the polarization. However, the energy dissipation due to applied voltage is miniscule, and therefore it does not appear to be a bottleneck provided the polarization switches and the interface coupling between the $P$-layer and the trilayer persists. In any case, further progress in the experimental front can handle such issues better.

Figure 18 presents a sample magnetization dynamics in the presence of room-temperature (300 K) thermal fluctuations. The ramp period is considered to be 100 ps and the it turns out from the simulation that switching has completed in less than 175 ps. Note that during the course of switching, random thermal kicks have forced magnetization to backtrack temporarily, however, the strong interface anisotropy has eventually enforced magnetization to switch from $\theta\simeq 180^{\circ}$ to $\theta\simeq 0^{\circ}$.

Figure 19 plots the average switching delay versus average energy dissipation for a range of ramp periods (0.1 ns – 1 ns). A moderately large number of simulations (10,000) in the presence of room-temperature (300 K) thermal fluctuations are performed and they are averaged to generate each point in the curve. For each trajectory of the 10,000 simulations, when the magnetization reaches $\theta\leq 5^{\circ}$, the switching is deemed to have completed. Note that as we increase the ramp period of applied voltage across the heterostructure, the switching delay also increases and less energy is dissipated in the switching process, elucidating the well-established delay-energy trade-off for a device in general. The results clearly demonstrate that that switching in sub-nanosecond delay is possible while dissipating a miniscule amount of energy of $\sim$1 aJ. Note that the “$CV^{2}$” energy dissipation is a couple of orders of magnitude lower than the energy dissipation due to Gilbert damping and it decreases with the increase of ramp period since the switching becomes more adiabatic.[93, 42] While the Figure 19 provides the mean of the switching delay distribution, the standard deviation in switching delay is also an important performance metric. the standard deviation for ramp period of 0.1 ns is about 22 ps and it increases about twice when the ramp period is increased to 1 ns. At higher ramp period, thermal fluctuations have more time to scuttle the magnetization and cause variability in switching time, increasing the standard deviation.

To understand the effect of temperature on the performance metrics, simulations have been performed at an elevated temperature (400 K) and the metrics switching delay and energy dissipation turn out to be similar (within 5%) compared to that of room-temperature (300 K) case.[87] Interestingly, the mean switching delay at a higher temperature T = 400 K decreases compared to the case at T = 300 K, which can be traced out from the reasoning that the initial deflection of magnetization due to thermal fluctuations increases at a higher temperature. Hence, magnetization is likely to start far away from the easy axis at a higher temperature for different trajectories, leading to the decrease in the mean switching delay. It turns out that this decrease in mean switching delay at T = 400 K is very small (less than 2%) compared to that of T = 300 K. However, the trend of standard deviation in switching delay with the increase in temperature shows an opposite trend to that of the mean. It can be understood by considering that the standard deviation of random thermal field at higher temperature increases with temperature [see Eq. (10)]. The mean energy dissipation decreases with increasing temperature and this decrease at T = 400 K is quite small (less than 3%) compared to the case of T = 300 K. This once again signifies the switching delay-energy trade-off for this device.

So far we have considered the writing of a bit of information by switching the magnetization from one state to another, however, the magnetization state needs to be read too. In this interface and exchange coupled structure, the giant magnetoresistance (GMR)[155, 156] of the trilayer is calculated to be of the order of 30%,[71, 157] which provides a way to read the magnetization states (P-alignment or AP-alignment). Although this GMR is not that high compared to tunneling magnetoresistance (TMR),[62, 63, 64, 65, 66, 67, 68] suitable design strategies can be possibly be devised to work with this moderate value of GMR and also it may be possible to increase the magnetoresistance by suitable material choice and design. It is also argued in Ref. \refciteRefWorks:649 that even with the variance in the smaller thicknesses of the layers, it is still possible to interface-couple the polarization and magnetization in the proposed structures. It should be noted that the modeling of interface anisotropy is not limited to the way that is performed here, however, any strong interface-coupled system would facilitate switching of magnetization from one state to the another.

The important consequence of having such a strong interface anisotropy is that we can achieve devices of very small lateral area, and it is therefore possible to cram an enormous amount of devices on a single chip. If we use an area density of $10^{-12}$ cm*-2*, the dissipated power would be 10 mW/cm2 considering 1 aJ energy dissipation in a single nanomagnet with 1 ns switching delay and 10% switching activity (i.e., 10% of the magnets switch at a given time). Such extremely dense and ultra-low-energy non-volatile computing systems can be powered by energy harvesting systems without the need of an external battery.[158, 159, 160, 161]

4 Summary and outlook

We have reviewed the dynamical systems study for electrical field-induced magnetization switching in strain-mediated multiferroic composites, and interface and exchange coupled multiferroic heterostructures. The magnetization switching dynamics in strain-mediated multiferroic composites using Landau-Lifshitz-Gilbert equation revealed intriguing phenomena in binary switching mechanism. It is shown that binary switching in a ‘symmetric’ potential landscape can be successful in the presence of room-temperature thermal fluctuations. To achieve such symmetry-breaking, we require the following two criteria: (1) a sufficiently high stress that keeps the magnetization more out of magnet’s plane inside the so-called “good” quadrants; and (2) a sufficiently fast ramp rate that reduces the possibility of backtracking causing switching failure. A high stress and a fast ramp rate also increase the switching speed and counters the detrimental effects of thermal fluctuations. This can potentially open up a new methodology of binary switching since tilting the potential landscape would not be necessary and such findings would stimulate experimental research to establish the proposed methodology of binary switching. As stated, it requires to sense when magnetization reaches around hard-plane so that stress can be brought down thereafter to achieve successful switching by breaking the symmetry of having equal probability of successful switching and switching failure in a‘symmetric’ potential landscape. Other ways to break this symmetry may be harnessed and is subject to further research. Note that the aforesaid switching in strain-mediated multiferroic composites just toggles the magnetization state without being able to maintain the direction of switching. It is shown that magnetization switching in interface and exchange coupled multiferroic heterostructures can maintain the direction of switching.

The calculated performance metrics from LLG simulations e.g., switching delay and energy dissipation show a profound promise for technological applications. The results show that switching can take place in sub-nanosecond delay while expending a miniscule amount of energy of $\sim$1 attojoule. This energy dissipation is at least 2-3 orders of magnitude lower than that of the other emerging devices. So multiferroic magnetoelectrics are intriguing in respect to both basic physics of binary switching and applied physics. Also, a strong interface and exchange coupling in multiferroic heterostructures enforces error-resiliency during the switching process and facilitates to scale down the lateral dimensions to unprecedented dimensions of $\sim$10 nm even in the presence of room-temperature thermal fluctuations. This is very crucial since it will help competing with the traditional charge-based electronics and consequently for building future nanoelectronics. Due to these superior performance characteristics of multiferroic magnetoelectrics as described here, currently it is of immense interest to analyze different possible theoretical designs followed by experimental demonstrations. Successful experimental implementations must tackle the issue of process variation at low dimensions, which traditional transistors face too. Processors built on such devices can lead to unprecedented applications that can work by harvesting energy from the environment without the need of an external battery e.g., medically implanted device to warn an impending epileptic seizure by monitoring the brain signals, wearable computers powered by body movements etc. Moreover, the basic building blocks are non-volatile permitting instant turn-on computer, facilitating computational designs, and avoiding energy dissipation for bit storage elements in a system. This can potentially perpetuate Moore’s law beyond the year 2020 and turn out to be an unprecedented opportunity for ultra-low-energy computing in our future information processing systems. Experimental efforts are emerging and successful experimental implementation can potentially pave the way for our future nanoelectronics.