We propose an exact expression to describe the hysteresis loops of an ensemble of Stoner–Wohlfarth particles undergoing an alternating quasi-static magnetic field. A statistical approach, which treats the quantities characterizing each particle as random variables, is adopted to get the orientation distribution of the local polarizations with respect to the applied field direction and the constitutive equation of the whole particle assembly. The hysteresis loop area gives the energy loss figure, but we have also obtained a straightforward integral expression for this quantity. The analytical relationships for the symmetric loops and the losses are successfully tested against numerical results, and the mathematical method adopted also displayed the ability to reproduce the “elemental loop” associated with any given particle of the system. While having a fundamental character, the proposed approach bears applicative interest, representing a versatile tool as the core of codes that simulate the behavior of devices employing magnetic components.

## I. INTRODUCTION

Magnetic materials exhibit a broad variety of compositions and microstructures and a wide spectrum of responses to a varying magnetizing field $H\u20d7$.^{1–3}

Although such behavior is the result of three-dimensional (3D) magnetization process, a two-dimensional theoretical approach is often regarded as appropriate, in addition to being perfectly adequate for a large class of materials, such as steel sheets, amorphous and nanocrystalline alloys, and thin films. It is also possible, in many cases, to consider a given material as an ensemble of mesoscopic regions (“particles”), either interacting or independent, each characterized by the local polarization $J\u20d7s$ and a magnetic anisotropy {uniaxial [see Fig. 1(a)], biaxial, cubic, etc.}, and even single-domain, if sufficiently small.

The Stoner–Wohlfarth (SW) theory,^{4–8} in its native formulation (which considers single-domain particles as independent, endowed with uniaxial anisotropy, and at *T* = 0 K), fits the simplifications mentioned above and is, thus, widely exploited in modeling. The main fields where the SW approach (often with some modifications) plays an important role are (1) magnetic nanoparticles,^{9–12} notably for magnetic hyperthermia;^{13,14} (2) ferromagnetic thin films^{15–17} and magnetic random access memories;^{18–20} (3) amorphous and nanocrystalline materials;^{21–25} (4) soft magnetic composites;^{26,27} (5) permanent magnets;^{28} and (6) study of the effect of spin-polarized currents.^{29,30}

In the SW model, the response to an applied field $H\u20d7$ of a particle characterized by the local uniaxial anisotropy $K\u2194$ is described by a vector hysteresis transducer providing the constitutive law for the reversible coherent rotations and the irreversible jumps of the local polarization $Js\u20d7$. In other words, one gets the equilibrium orientation *γ* of each $J\u20d7s$ as a function of the anisotropy constant *K*, its orientation *φ*, and $H\u20d7$ [see Fig. 1(a)]. With $J\u20d7$ being the polarization of the whole particle assembly, the $J\u20d7(H\u20d7)$ constitutive equation of the system comes out after integration over the probability density functions (pdfs) of *K* and *φ*. The design of magnetic components needs to implement $J\u20d7(H\u20d7)$ in numeric modeling codes: a step with high computational cost, scaling up with the device size. Indeed, an explicit *γ* = *γ*_{H}(*K*, *φ*) solution for the $J\u20d7s$ orientation is not available so far, being the equilibrium angle obtained either following numeric procedures^{25,31} or using analytic approximations.

In the first case, the problem is commonly faced by adopting the astroid representation for the SW transducer [see Fig. 1(b)] and exploiting the so-called “tangent method,”^{7} as well as searching for the minimum of the Gibbs free energy of each particle utilizing, e.g., the Newton–Raphson method. These are time-consuming approaches, not efficient when implemented in codes used to design magnetic devices.

Researchers following the second way tried, for example, to recover the SW loop branches with *ad hoc* analytic functions,^{32} or proposed an *H* = *H*(*J*) “inverse” constitutive law,^{12} or limited the analytic solution to the low-field regime.^{33} Alternatively, an SW-like hysteron with different kinds of approximations for the anisotropy energy term,^{34,35} or with an equivalent field playing the role of the anisotropy,^{36} has been also adopted.

This work is not focused on the search of the *γ* = *γ*_{H}(*K*, *φ*) single particle solution [see Fig. 1(a)], afterward integrated over the whole system, but it deals at the start from the entire particle assembly by handling each quantity as a random variable (rv) with associated pdf. Statistical methods^{37,38} are then utilized to retrieve *p*_{H}(*γ*): the equilibrium angle pdf, ending with an exact expression for the hysteresis loops (and consequently for the energy losses) of the whole system. The single particle behavior comes out as a particular case. Apart from its importance from the basic viewpoint, this result turns out to be useful for applications (chiefly in the magnetic loss prediction^{39}), shortening the cumbersome and slow procedure needed for the design of devices.

In this paper, after a review of the SW model features essential for our purposes (Sec. II A), and some hints about the numerical approach and the adopted statistics (Sec. II B), the core of the work is discussed in Sec. III. Here, the statistical approach is outlined, and the pdfs controlling the equilibrium orientation of the local polarization vs a quasi-static $H\u20d7$ are worked out, separately accounting for the irreversible (Sec. III A) and reversible (Sec. III B) phenomena. In both these subsections, the pdf evolution is described starting from the demagnetized state and following the system magnetization along the branches of the hysteresis loop driven by an alternating field [Fig. 1(c)], separately for the four *Q* quadrants shown in Fig. 1(d). Eventually, the obtained *constitutive equation* (and thus the hysteresis loop and loss) is successfully tested against the numerical procedure, both for the entire assembly (Sec. III B 5) and for the case of a single particle (Sec. III B 6). An approximate analytic expression for the energy loss, not involving the knowledge of the hysteresis loop, is worked out and effectively checked (Sec. IV), and finally, in Sec. V, some possible advancements of the proposed approach are listed.

## II. THE CLASSICAL SW MODEL AND ITS NUMERICAL IMPLEMENTATION

### A. The Stoner–Wohlfarth model

*T*= 0 K, each characterized by a local anisotropy axis with modulus

*K*(and anisotropy field

*H*

_{K}= 2

*K*/

*J*

_{s}, where

*J*

_{s}= $|J\u20d7s|$), forming an angle

*φ*with a reference direction defined by the applied field $H\u20d7$ [Fig. 1(a)]. After introducing the Gibbs free energy density of the particle, together with its first and second derivative (where

*H*= $|H\u20d7|$),

*γ*equilibrium orientations of the local polarization $J\u20d7s$ (corresponding to the minima of

*u*

_{H}) are derived.

^{4–7}This can be done starting from the conditions $uH\u2032=uH\u2033=0$ (identifying the

*u*

_{H}horizontal inflection point generated by the merging of one minimum and one maximum), which allows one to draw, in the (

*H*

_{‖}≔

*H*cos

*φ*,

*H*

_{⊥}≔

*H*sin

*φ*) Cartesian reference frame, a closed curve (astroid) whose contour (with vertices ±

*H*

_{K}) represents the border separating the (

*H*

_{‖},

*H*

_{⊥}) plane region corresponding to two minima of

*u*

_{H}(inside) from the one where one minimum only appears [Fig. 1(b)]. On passing to the polar coordinates (

*φ*,

*H*), after defining the quantity

*H*is instead given, in the (

*φ*,

*K*) plane, where each point represents a particle of the system, one finds the corresponding always positive anisotropy threshold,

*K*is lower (higher) than

*K*

_{c,H}.

The introduction reports the strategies followed in the literature to find the *γ* equilibrium angle for each particle. In the following, the statistics of the system is described by the *ψ*(*K*, *φ*) joint probability density function, with *K* > 0 and *φ* limited to the [0;*π*/2] range for symmetry reasons.

### B. Numerical implementation

*J*(

*H*), worked out in the following, finds the local $J\u20d7s$ equilibrium orientations via an iterative procedure searching for the

*u*

_{H}absolute minimum and follows its time evolution by tracking the memory of the irreversible rotations (IRs), for each particle. On passing to an assembly of SW particles, one gets the behavior of the whole system after integration over the

*ψ*(

*K*,

*φ*) pdf. When we compare the analytic expressions for the hysteresis loop and loss to the numerical simulations,

*φ*and

*K*are assumed to be statistically independent variables. Consequently, we can write

*ψ*(

*K*,

*φ*) =

*f*(

*K*)g(φ), and in particular, we have chosen for the tests, with ⟨

*K*⟩ the average anisotropy value,

## III. THE ANALYTIC-STATISTICAL APPROACH FOR AN ENSEMBLE OF SW PARTICLES

In this work, we have abandoned the problem of finding the exact analytical solution for the local $J\u20d7s$ equilibrium orientation *γ* of a single SW particle in favor of a statistical approach that considers the whole assembly of SW particles from the beginning. Accordingly, the quantities *K*, *φ*, and *γ* (−*π* to +*π*) are treated as *random variables*, along with their *ψ*(*K*, *φ*) and $pH(b,Q)(\gamma )$ pdfs (the first normalized to 1 for each *φ* angle). In the latter, “*b*” indicates the loop branch considered [see Fig. 1(c)], “*Q*” is the quadrant to which *γ* pertains [see Fig. 1(d)], and *H* plays the role of a parameter. All along the paper, the evolution of the *γ* values, which are bound to increase from −*π* to *π* [see Fig. 1(d)], dictates the sequence adopted to discuss the *Q* quadrants (*III* → *IV* → *I* → *II*).

*H*

_{p}, the corresponding hysteresis loop is obtained as follows: The magnetization process is assumed to start from the demagnetized state (DEM), where

*H*=

*J*= 0 [Fig. 1(c)]—a situation that, for symmetry reasons, can be obtained by formulating the pdf for the $J\u20d7s$ equilibrium orientation still in terms of

*K*and

*φ*instead of

*γ*, the local polarization lying along the easy axis direction in each particle,

*φ*=

*φ*

_{1}and a postulated

*ψ*(

*K*,

*φ*

_{1}) vs

*K*behavior, is reported in Fig. 3(a).] The

*J*vs

*H*evolution is then followed and investigated in sequence [see Fig. 1(c)] along the First Magnetization Curve (FMC), the RECoil curve (REC), and the Negative Descending Branch (NDB) of the hysteresis loop, both in Secs. III A and III B. Note that the ascending branch (↑) of the major loop of vertex (

*H*

_{p},

*J*

_{p}) can be obtained from the descending one (↓), build-up of the REC and NDB curves, being

*J*

_{↑}(

*H*) = −

*J*

_{↓}(−

*H*) for symmetry reasons.

Along each loop branch, the applied field triggers the irreversible rotations (IRs) and drives any local $J\u20d7s$ of the whole particle ensemble to its equilibrium angle (reversible rotations) at the same time. Despite the ensuing complex interplay between reversible and irreversible processes, the technique developed here treats them separately, and the *γ* pdf is consequently worked out in two steps: before finding $pH,IRR(b,Q)(K,\phi )$, i.e., the pdf determined by the IRs only, which bring the local $J\u20d7s$ to their *absolute minimum* (Sec. III A), and subsequently accounting for the coherent reversible rotations (by means of the $J\u20d7s$ relaxation to the *equilibrium orientation*), ending with $pH(b,Q)(\gamma )$ (Sec. III B). This approach is sketched in Fig. 2.

### A. Step 1: Irreversible rotations (IRs)

At this stage, we assume $H\u20d7$ to drive the irreversible processes only, again leaving the local $J\u20d7s$ to point along one of the two sides of the *φ* easy axis directions of the particle to which it pertains. For this reason, the magnetization distribution among the quadrants, rearranged by the occurrence of the IRs, is again described by a pdf stated in terms of *K* and *φ* instead of *γ*: $pH,IRR(b,Q)(K,\phi )$, for the three loop branches of Fig. 1(c), as listed in Sections III A 1 -III A 4, and illustrated in Fig. 3 for a given *φ* = *φ*_{1}.

#### 1. First magnetization curve (FMC)

For *H* increasing from DEM [Fig. 1(c)], Eqs. (1b) and (1c) state that the $J\u20d7s$ in quadrant *III* smoothly (i.e., without IRs) crosses the *γ* = −*π*/2 threshold and move to quadrant *IV* when *φ* > *π*/4 and *H* > *H*_{III→IV} = (*K*/*J*_{s})sin 2*φ* [a value always lower or equal to *H*_{c,K}(*φ*): the threshold for IRs]. As stated above, we are not interested at this stage in the equilibrium position of $J\u20d7s$, so we account for these $J\u20d7s$ reversible rotations by assigning the same pdfs to these quadrants: $pH,IRR(FMC,IV)=pH,IRR(FMC,III)$. Note that this “doubling” of the $pH,IRR(FMC,Q)$ does not affect the normalization because each $J\u20d7s$ does not pertain to *Q* = *III* and *Q* = *IV* for the same field, the Step 2 (Sec. III B) taking charge of finding the actual $J\u20d7s$ orientation vs *H*. As *H* increases from zero [where *K*_{c,0} = 0, see Eq. (4)] to the peak value *H*_{p}, the $J\u20d7s$ in quadrant *III* or *IV* switch to quadrant *I* if *K* < *K*_{c,H}(*φ*), and the corresponding $pH,IRR(FMC,Q)(K,\phi )$ is reported in Fig. 3(b). Note that $pH,IRR(FMC,II)$ = 0 because the system geometry [Fig. 1(a)] does not allow irreversible rotations to the second quadrant, a fact in agreement with the detailed discussion reported in the Appendix, which shows that no energy minima of *u*_{H} [Eq. (1a)] can appear in *Q* = *II*.

#### 2. Recoil curve (REC)

When the field is reduced from *H*_{p} to zero, the polarization decreases from its peak value *J*_{p} = *J*(*H*_{p}) to the remanence *J*_{r}. Along this branch, no IRs occur, and thus, $pH,IRR(REC,Q)(K,\phi )$ is frozen to the value reached along FMC at *H* = *H*_{p}: $pHp,IRR(FMC,Q)(K,\phi )$ [see Fig. 3(c)]. Here, again, we have $pH,IRR(REC,IV)$ = $pH,IRR(REC,III)$ because the “doubling” of the pdfs, pointed out in FMC persists, but, as above, it does not affect the normalization.

#### 3. Negative descending branch (NDB)

Eventually, the system is brought from the remanence to the negative tip of the hysteresis loop, lowering *H* from 0 to −*H*_{p}, with the magnetization evolving as described by $pH,IRR(NDB,Q)(K,\phi )$ [Fig. 3(d)]. For what concern the IRs and the smooth $J\u20d7s$ transitions discussed in Sec. III A 1, a “symmetrical” reasoning applies here, with the quadrants *I*, *II*, and *IV* correspondingly playing the role of *Q* = *III*, *IV*, and *I* in FMC, and the *γ* = *π*/2 crossing occurring at the threshold *H*_{I→II} = −*H*_{III→IV} [≤*H*_{c,K}(*φ*)], when *φ* > *π*/4.

#### 4. Thresholds of the $pH,IRR(b,Q)(K,\phi )$ functions

The quantities *K*_{c,H}(*φ*_{1}) and *K*_{c,Hp}(*φ*_{1}) of Fig. 3, signaling a variation of the $pH,IRR(b,Q)$ values for a given *φ* = *φ*_{1}, become the *K*_{c,H}(*φ*) and *K*_{c,Hp}(*φ*) threshold functions when considering the entire 0 ≤ *φ* ≤ *π*/2 range. As *H*, starting from the demagnetized state (DEM), oscillates between the ±*H*_{p} peak values of the hysteresis loop, the *K*_{c,H}(*φ*) threshold follows it, swinging between *K*_{c,0} = 0 and *K*_{c,Hp}(*φ*), with the $pH,IRR(b,Q)(K,\phi )$ function that accordingly changes. An example for *Q* = *I* across NDB is depicted in the (*φ*, *K*) plane of Fig. 4.

### B. Step 2: Coherent reversible rotations and constitutive equation

Given *H* and *γ*, and considering the scheme of Fig. 2, for the three branches of the hysteresis loop indicated in Fig. 1(c), the $pH(b,Q)(\gamma )$ pdf, which accounts for irreversible and reversible processes, is obtained through standard statistical methods.^{37,38}

*γ*′ to satisfy the request $\gamma low\u2032\u2264\gamma \u2032\u2264\gamma $, where the lower limit is −

*π*, −

*π*/2, 0, and

*π*/2 when

*γ*pertains to

*Q*=

*III*,

*IV*,

*I*, and

*II*, respectively [Fig. 1(d)]. Having available the $pH,IRR(b,Q)(K,\phi )$ pdfs, this approach can be stated in terms of

*K*and

*φ*, i.e., identifying in the (

*K*,

*φ*) plane (limited to the 0 ≤

*φ*≤

*π*/2 range for symmetry reasons), a

*D*

_{H}(

*γ*) domain that satisfies the demand above $(\gamma low\u2032\u2264\gamma \u2032\u2264\gamma )$—a result achievable if a

*relationship connecting the rv*γ′

*to the rv’s*K

*and*

*φ*

*:*

*γ ′ = γ*

_{H}

*(K, φ)*exists. In the case of an SW particle, for which it is not possible to get such an explicit

*γ*

_{H}function, the connection between the rv’s can be stated in an implicit way from the request that

*u*

_{H}(

*K*,

*φ*;

*γ*′) [Eq. (1a)] shows a

*minimum,*i.e., from a couple of conditions,

*K*[see Eq. (1b)], provides a relationship valid for

*H*and

*φ*values making it positive,

*K*= $KH*$, becomes [see Eq. (1c)]

*φ*values where it changes sign in the 0 ≤

*φ*≤

*π*/2 interval,

*n*

^{(Q)}= 2, 1, 0, and −1 when

*Q*=

*III*,

*IV*,

*I*, and

*II*, respectively. Furthermore, in the 0 ≤

*φ*≤

*π*/2 range, both the $KH*$ and $uH\u2033*$ functions have a vertical asymptote at

*D*

_{H}(

*γ*) domain (magenta regions in Fig. 5), one has to find before, for each

*γ*′, the $\Delta \phi *(\gamma \u2032)$ interval, i.e., the

*φ*values range making

*u*

_{H}(

*K*,

*φ*;

*γ*′) minimum, which is defined by the overlap region where both $KH*(\phi ;\gamma \u2032)$ and $uH\u2033*(\phi ;\gamma \u2032)$ are positive. This process is described in the Appendix, for the three loop branches [Fig. 1(c)] and the four quadrants [Fig. 1(d)].

*D*

_{H}(

*γ*) domain emerges as the (

*φ*,

*K*) plane region swept by all the $KH*(\phi ;\gamma \u2032)$ functions [each of them limited to its $\Delta \phi *(\gamma \u2032)$ range] for increasing

*γ*′ between $\gamma low\u2032$ to

*γ*(white curves in Fig. 5). This procedure is summarized by the $D\u0302$-operator,

*D*

_{H}(

*γ*) allows one to write the cumulative distribution function, given

*H*and

*γ*, and its derivative—the probability density function,

*K*in Eq. (13a), one must remember that the thresholds

*K*

_{c,H}(

*φ*) and

*K*

_{c,Hp}(

*φ*) [see Eq. (4)] indicate, for each

*φ*, the

*K*values where the $pH,IRR(b,Q)(K,\phi )$ pdf changes, as discussed in Sec. III A 4. As an example, Fig. 4 displays the case

*Q*=

*I*across NDB.

#### 1. First magnetization curve (FMC)

According to Sec. III A 4 and Fig. 3(b), along this branch, the *K*_{c,Hp}(*φ*) curves do not play any role (and thus are not depicted in Fig. 5), whereas the *K*_{c,H}(*φ*) ones not only signal the *K* thresholds where the $pH,IRR(FMC,Q)$ pdfs [retrieved from Fig. 3(b)] change but also represent a portion of the *D*_{H} integration domains (magenta regions in Fig. 5, for *Q* = *III* and *Q* = *IV*) border.

*a)* $\u2212\pi \u2264\gamma \u2264\u2212\pi /2(Q=III)\u0304$

*D*

_{H}(

*γ*) domain is obtained employing the $D\u0302$-operator [Eq. (12)] with $\gamma low\u2032$ = −

*π*. Being $\phi \u221e(III)(\gamma )$ =

*γ*+

*π*[Eq. (11) with

*n*

^{(Q)}= 2], Eq. (13a) becomes

*γ*′ =

*γ*and

*n*

^{(Q)}= 2, and the integrand is defined as follows:

*b)* $\u2212\pi /2\u2264\gamma \u22640(Q=IV)$

*D*

_{H}(

*γ*) domain is obtained by utilizing the $D\u0302$-operator [Eq. (12)] with $\gamma low\u2032$ = −

*π*/2. The top-left border of

*D*

_{H}(

*γ*) is given [Eq. (8)] by the curve $KH*(\phi ;\u2212\pi /2)=HJs/sin2\phi $, with asymptote

*φ*=

*π*/2—a value not accounted for by Eq. (11). The cumulative and density functions [the second one remembering again Eq. (A2)] become

*γ*′ =

*γ*and

*n*

^{(Q)}= 1, and the integrand defined similarly as in

*Q*=

*III*[Eqs. (16) and (17)]

*c)* $0\u2264\gamma \u2264\pi /2(Q=I)$

*D*

_{H}(

*γ*) domain is obtained using the $D\u0302$-operator [Eq. (12)] with

*γ*′ increasing from $\gamma low\u2032$ = 0. The cumulative function [being

*φ*

_{∞}(

*γ*) =

*γ*from Eq. (11) with

*n*

^{(Q)}= 0], and the pdf associated [observing that $KH*(\gamma ;\gamma )\u2192\u221e$, from Eq. (8)] become

*Q*=

*III*[Eqs. (16) and (17)].

*d)* $\pi /2\u2264\gamma \u2264\pi (Q=II)$

*D*

_{H}, and agrees with the fact that, in this region, we have $pH,IRR(FMC,II)$ = 0 (Sec. III A and Fig. 2). Therefore,

#### 2. Recoil curve (REC)

*K*

_{c,H}threshold decreasing with

*H*, from $Kc,Hp$ to zero. However, as shown in Fig. 3(c), due to the absence of IRs along this branch, the “irreversible” pdfs are frozen to the last value reached across FMC, i.e., when

*H*=

*H*

_{p},

*Q*quadrants. Accordingly, the integrands [Eq. (16)] turn into

#### 3. Negative descending branch (NDB)

*D*

_{H}shapes in

*Q*=

*III*,

*IV*,

*I*, and

*II*recovering the ones found along FMC (or REC) in

*Q*=

*I*,

*II*,

*III*, and

*IV*, respectively, but with different

*φ*integration limits (see Fig. 5). Now, the top-left border of

*D*

_{H}in

*Q*=

*II*is given [Eq. (8)] by the curve $KH*(\phi ;\pi /2)=\u2212HJs/sin2\phi $, with asymptote

*φ*=

*π*/2—a value not accounted for by Eq. (11). Thus, a procedure like the one described above supplies the

*γ*pdfs for the four quadrants,

*K*

_{c,H}(

*φ*) and

*K*

_{c,Hp}(

*φ*), as shown in Fig. 3(d).

#### 4. Compact expression for the $pH(b,Q)(\gamma )$ pdfs

*γ*with $H\u20d7$ in a very compact integral form, depending on the “

*b*” branch (FMC, REC, or NDB) and the

*Q*(=

*III*,

*IV*,

*I*, and

*II*) quadrant occupied by the equilibrium angle

*γ*,

*b*” and “

*Q*,” the integration limits of

*Q*=

*II*in FMC and REC and of

*Q*=

*IV*in NDB are “artificially” put equal to zero, thus obtaining $pH(FMC,II)$ = $pH(REC,II)$ = $pH(NDB,IV)$ = 0.

#### 5. Constitutive law and hysteresis loop

*constitutive equation*, calculated over the whole

*γ*domain,

*W*is given by its area. Moreover, as described below (Sec. IV), it is also possible to work out a semi-empirical but very accurate integral expression for

*W*not entailing the knowledge of the hysteresis loop.

Figure 7 shows, for the whole −*π* ≤ *γ* ≤ *π* range, an example of $pH(b,Q)(\gamma )$ vs *H* evolution, calculated with Eq. (26), and the corresponding hysteresis loop [Eq. (28)] with *J*_{s} = 1.61 T. We have assumed *ψ*(*K*, *φ*) =*f*(*K*)*g*(φ), with the marginal pdfs supplied by Eqs. (5a) and (5b) and $K$ = 3000 J/m^{3}.

The effectiveness of the statistical method is apparent in Fig. 8, where the quasi-static hysteresis loops, computed at different peak polarization *J*_{p} using Eq. (28) and the corresponding loss figures (inset), perfectly recover the ones obtained via the numerical procedure presented in Sec. II B [again it is assumed *ψ*(*K*, *φ*) = *f*(*K*)*g*(φ), with the marginal pdfs furnished by Eqs. (5a) and (5b)].

#### 6. Elemental hysteresis loops

*γ*equilibrium angle for the local polarization $J\u20d7s$ [see Fig. 1(a)]—a goal that would directly supply the analytic expression for the elemental loops associated with a single particle. This difficulty has been circumvented by the mathematical technique outlined above, which supplies, in a sense, a “

*statistical solution*” of Eq. (8), working out

*p*

_{H}(

*γ*), the probability to find $J\u20d7s$ laid along the

*γ*direction. In this framework, the elemental hysteresis loops corresponding to a particle characterized by

*φ*=

*φ*

_{1}and

*K*=

*K*

_{1}come out defining the

*K*and

*φ*joint statistics as follows:

*δ*

_{D}is the Dirac delta function.

*γ*, for any

*K*

_{1},

*γ*= 0 [in particular, the $J\u20d7s$ belonging to the third quadrant when

*H*= 0 remains confined in it if

*φ*≤

*π*/4 or, when

*φ*>

*π*/4, as long as

*H*≤

*H*

_{III→IV}(see Sec. III A 1)], and no loop appears. When the alternating

*H*overcomes the inversion field

*H*

_{p}= $Hc,K1(\phi 1)$, the $J\u20d7s$ particles in

*Q*=

*III*or

*Q*=

*IV*jump to the first quadrant in the 0 ≤

*γ*≤

*φ*

_{1}range [with

*p*

_{H}(

*γ*) →

*δ*

_{D}(

*γ*) when

*H*→ +

*∞*]. A similar and, to some extent, specular behavior is found following the recoil curve and the negative descending branch.

## IV. ENERGY LOSS: AN APPROXIMATE INTEGRAL EXPRESSION

*H*

_{c,K}(

*φ*) [Eq. (3)], peak polarization

*J*

_{p}(

*φ*), remanence

*J*

_{r}(

*φ*), and loss

*W*. Now, we take as peak value the polarization reached immediately after the irreversible switch occurring at

*H*

_{c,K}under increasing field:

*J*

_{p}(

*φ*) ≔ $J[Hc,K+]$ (a larger

*H*moves

*J*along a reversible branch, non-modifying the loop area). An example obtained for

*φ*= 35° (

*J*

_{s}= 1 T,

*K*= 1 J/m

^{3}) is given by the green solid line in Fig. 10(a), where instead the orange dashed rectangle corresponds to the case where

*irreversible rotations*(180° switches between

*φ*and

*φ*–

*π*of the local $Js\u20d7$) only are present. The area (giving the energy loss in J/m

^{3}) of the latter, with peak irreversible polarization

*J*

_{IRR,p}≡

*J*

_{r}=

*J*

_{s}cos

*φ*, is supplied by the following relationship:

*A*(

*φ*) given by Eq. (2).

*φ*, the difference between

*W*and

*W*

_{IRR}is concentrated in the gray pseudo-triangles with vertices ±

*J*

_{r}, ±

*J*

_{c}, and ±

*J*

_{p}[Fig. 10(a)]: all linear quantities dependent on

*J*

_{s}and independent of

*K*[

*J*

_{r}=

*J*

_{s}cos

*φ*, and for

*J*

_{c}and

*J*

_{p}, see Eqs. (8.34), (8.38), and (8.41) of Ref. 5]. Being

*H*

_{c,K}[Eq. (3)] inversely proportional to

*J*

_{s}, as this last varies the loops aspect ratio modifies, but

*W*and

*W*

_{IRR}[Eq. (31)] do not. Moreover,

*H*

_{c,K}is the same for both loops, and one gets

*W*

_{IRR}(

*aK*,

*φ*) =

*a*

*W*

_{IRR}(

*K*,

*φ*) and

*W*(

*aK*,

*φ*) =

*a*

*W*(

*K*,

*φ*). For the arguments above, the ratio between the two losses [these last plotted in Fig. 10(b)] turns out to be

*φ*dependent only, as reported in Fig. 10(c), where it is interpolated by the following

*general empirical law*:

*W*/

*W*

_{IRR}vs $log10[1/(1\u22122\pi \phi )]$ [Fig. 10(c) inset].

*H*=

*H*

_{p}, the energy loss of the whole particle ensemble is obtained after integration over the anisotropy axes values and orientations [from Eqs. (31) and (32), and with

*K*

_{c,Hp}given by Eq. (4)],

*H*

_{p}, calculated on the First Magnetization Curve with Eq. (28),

## V. CONCLUSIONS

The developed mathematical tool analytically predicts the response of an assembly of non-interacting SW particles driven by a quasi-static alternating magnetic field. The exact integral relationships for the hysteresis loops and losses allow one to bypass the computationally slow numeric procedure necessary to implement the SW model in software used to predict the behavior of devices containing magnetic components.

The versatility of this approach lies in its skills listed below.

It works with any anisotropy distribution, permitting one to easily predict the role played by the sample texture (e.g., the effect of a macroscopic easy axis induced either by annealing or by applied stress).

The same analytic-statistical procedure, although more complicated, can be exploited with any field history (e.g., rotating fields,

^{31,34}loops with bias, and asymmetrical minor loops).Marginal modifications of the approach allow it to account for the presence of domain walls inside the particles, following the idea suggested in Ref. 22.

The whole mathematical structure is prone to a generalization for systems of particles endowed with biaxial anisotropy.

Under whatever field history and accounting for interactions, in Ref. 40, a mathematical technique capable to track the irreversible switches in a particle assembly described by the SW model has been worked out (so preserving the non-local memory of the system), together with a graphic representation that turns out to be an extension of the one operating in the Preisach model. One can then envisage a coupling between such a tool and the one developed here, with the aim to assess a more general analytical approach to the particle systems description.

## ACKNOWLEDGMENTS

This research work was carried out in the framework of the Grant No. 19ENG06 HEFMAG project, which was funded by the EMPIR program, and co-financed by the Participating States and the European Union’s Horizon 2020 research and innovation program.

## AUTHOR DECLARATIONS

### Conflict of Interest

The author has no conflicts to disclose.

### Author Contributions

**C. Appino**: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Investigation (lead); Methodology (lead); Software (lead); Supervision (lead); Validation (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead).

## DATA AVAILABILITY

The data that support the findings of this study are available within the article.

### APPENDIX: MINIMA OF THE GIBBS FREE ENERGY DENSITY *u*_{H} VS THE ANISOTROPY AXIS ORIENTATION *φ*

For a general *γ*′ equilibrium orientation of $J\u20d7s$, with $\gamma low\u2032\u2264\gamma \u2032\u2264\gamma $ (where *γ* depends on the *Q* quadrant considered, and $\gamma low\u2032$ is its lower limit), the $\Delta \phi *(\gamma \u2032)$ range, corresponding to a minimum of the particle Gibbs energy *u*_{H}(*K*, *φ*; *γ*′) [Eq. (1a)], is determined by the overlap between two *φ* regions: the first one, where $KH*(\phi ;\gamma \u2032)$ [Eq. (8)] is positive, and the second one satisfying Eq. (9), $uH\u2033*(\phi ;\gamma \u2032)>$ 0, as illustrated in Fig. 11 (remember that *φ* is limited between 0 and *π*/2 for symmetry reasons).

The $\Delta \phi *(\gamma \u2032)$ interval is found before for the field values 0 ≤ *H* ≤ *H*_{p}, i.e., when the First Magnetization Curve (FMC) and the RECoil curve (REC) are covered (upper strip of Fig. 11). In this case, for *Q* = *III* and *Q* = *IV*, the two regions partially overlap, and the $\Delta \phi *(\gamma \u2032)$ limits are [$\phi 0(III)$; $\phi \u221e(III)$ = *γ*′ + *π*] and [$\phi 0(IV)$; *π*/2], respectively [see Eqs. (10) and (11)]. In *Q* = *I*, the entire region where $KH*>$ 0 falls in the range where $uH\u2033*>$ 0; $\phi 0(I)$, which drops outside, does not play any role, and the $\Delta \phi *(\gamma \u2032)$ limits are [$\phi \u221e(I)$; *π*/2]. Eventually, in *Q* = *II*, the two regions do not overlap, signaling the fact that, in this quadrant, no stable equilibrium position for $J\u20d7s$ exists. This fact agrees with the result $pH,IRR(FMC,II)$ = 0, which is reported at the end of Sec. III A, and is derived simply by observing the system geometry [Fig. 1(a)].

When *H* < 0 (lower strip of Fig. 11), i.e., along the Negative Descending Branch (NDB), a very similar and, to some extent, specular scenario appears—the situations in *Q* = *III*, *IV*, *I*, and *II* now recovering the ones found along FMC and REC in *Q* = *I*, *II*, *III*, and *IV*, respectively.

*Q*=

*III*and

*Q*=

*IV*for FMC and REC, and

*Q*=

*I*and

*Q*=

*II*for NDB). In this circumstance, inserting the $\phi 0(Q)$ value [Eq. (10)] into Eq. (8), one gets (considering that

*K*> 0)

*γ*′ parameter) identifying, in the (

*φ*,

*K*) plane, a curve that coincides with

*K*

_{c,H}, being this last obtained under the conditions $uH\u2032$ = $uH\u2033$ = 0, as well. Consequently, $\phi 0(Q)$ represents the abscissa of the tangent point between $KH*$ and

*K*

_{c,H}, and one gets

*p*

_{H}(

*γ*), in Sec. III B 1.

## REFERENCES

*Modern Magnetic Materials: Principles and Applications*

*Wiley Encyclopedia of Electrical and Electronics Engineering*

*The Science of Hysteresis*

*γ*-Fe

_{2}O

_{3}nanoparticles for magnetic hyperthermia

*Probability, Random Variables, and Stochastic Processes*