Micromagnetic Approach to the Metastability of a Magnetite Nanoparticle and Specific Loss Power As Function of the Easy-Axis Orientation

1. Introduction

It is well-established that nanoparticles with sizes ranging between 1 and 100 nm, exhibit different and appealing characteristics compared with their micro and macro counterparts. In recent years, magnetic nanoparticles (MNPs) have been used in areas of nanotechnology, in particular for biomedical applications [1,2] where precise delivery of anticancer drugs and selective killing of cancer cells is desirable. Magnetic hyperthermia treatment [3] is based on the use of nanoparticles dispersed in a magnetic fluid. When MNPs are under the presence of an external AMF, they are able to heat a specific area of the body causing cells death in tumor tissue [6]. During such process, two relaxation mechanisms for MNPs are known, namely, the Néel [4] and Brown [5] modes. In the former, magnetic moments rotate inside the nanoparticle due to the action of the external AMF. In contrast, in Brown relaxation, magnetic moments are locked to the magnetic easy-axis of the MNP and the AMF causes a mechanical rotation of the particle as a whole. The heat released to the medium is either caused by the friction between the surface of MNPs and the surrounding magnetic-fluid, or by the power dissipated through the hysteresis losses [7]. The heating efficiency of a MNP under an external AMF is expressed in the physical quantity called specific loss power (SLP) [8,9].

The dynamic response of a MNP in an external time-dependent magnetic field of amplitude $H_0$ and period P, is not instantaneous and it exhibits a delay in time. Time delay of the magnetization gives rise to the appearance of two possible regimes [10,11,12]: a dynamically ordered phase (DOP) and a dynamically disordered one (DDP), which depend on the relationship between P, $H_0$ and the metastable lifetime $\tau$[13,14]. In the DOP (see Figure 1a), magnetization is not able to follow the magnetic field, giving rise to a high value of the average magnetization. In DDP (see Figure 1b), magnetization follows the magnetic field even though with a $P-$dependent delay, and it gets to be reversed resulting in a well-defined hysteresis loop. In this case, the metastable lifetime, which is the elapsed time to go from the saturation state and pass through zero magnetization, can be computed.

The universal aspects that allow a dynamic phase transition have been studied in some works [10,11,12,13,14]. B. K. Chakrabarti and M. Acharyya [12] showed that magnetic response in time by an external field depends on the competing time scales, i.e. the time period of oscillation of the external perturbation and the typical relaxation time. When the time period of oscillation is comparable with the effective relaxation time, symmetric hysteresis loops around the origin are obtained. A breaking of the symmetry of the hysteresis appears when the driving frequency of the field increases. This phase transition depends on the field amplitude and temperature. On the other hand, G. Korniss et. al. [11] studied through Monte Carlo simulations the two-dimensional kinetic Ising model with a square-wave oscillating external field of period P. Their results showed that the system undergoes a dynamic phase transition when the half-period $P_{1/2}$ of the field is comparable to the metastable lifetime $\tau$, which in turn depends on temperature T and the field amplitude $H_0$. Several metastable states were also studied by P. A. Rikvold et. al. [14] for an impurity-free kinetic Ising model with nearest neighbor interactions and local dynamics under a magnetic field, in the framework of droplet theory and Monte Carlo simulation. They demonstrated that metastable lifetimes exhibited magnetic field and system-size dependences. Additionally, W. D. Baez and T. Datta [13], established for a two-dimensional ferromagnetic kinetic Ising model with oscillating field and next-nearest neighbor interactions, that metastable lifetimes are determined by the lattice size, amplitude and frequency of the external field, temperature and additional interactions present in the system.

Micromagnetic models [15], based on the solution of the Landau-Lifshitz-Gilbert equation [16,17,18], have also allowed a theoretical description of micro and nanoscale magnetization processes above the exchange length. Micromagnetism integrates classical and quantum mechanical effects, where the spin operators of the Heisenberg model are replaced by classical vectors and also accounts for exchange interaction. The main assumptions of the model are: the distribution of the magnetic moments is considered discrete throughout the volume of the magnetic system by means a set of discretization cells and at the same time it is approximated by a density vector $\mathbf{M}(\mathbf{r},t)$, which is continuous and differentiable with respect to both space $\mathbf{r}$ and time t, and it can be written in terms of a unit vector field $\mathbf{M}(\mathbf{r},t)=M_s\mathbf{m}(\mathbf{r},t)$ where $M_s$ is the saturation value of constant norm.

Despite all the studies devoted to systems exhibiting dynamic phase transitions, little effort has been paid to the role played by the spatial orientation of the easy-axis of a single nanoparticle relative to the direction and sense of the external applied oscillatory field. It must be stressed that different orientations of the magnetic easy-axis must be interpreted as different ensamble microstates occurring during Brown rotation, so by fixing the magnetic easy-axis in each microstate considered, we can drive our attention to the Néel relaxation at each Brown rotation step. Thus, the purpose of the present work is to calculate the metastable lifetimes for different amplitudes $H_0$ and orientations of the magnetic easy-axis $\mathbf{u}$ as well as the interplay between these two parameters. To do so, micromagnetic simulations were performed using the Ubermag [19] package based on the Object Oriented Micromagnetic Framework (OOMF) [20].

The remainder of the paper is organized as follows. In Section 2 we present an overview of the micromagnetic background of a MNP with a given easy-axis of magnetization, its geometrical and physical aspects and computational details. Numerical results are presented in Section 3 where a proposal of dynamic phase diagram is presented, and in Section 4 we make a summary and a brief conclusion of the main results.

2. Theoretical model and Computational details

2.1. Micromagnetic energy

According to the micromagnetic theory, the energy density $w\left(\mathbf{M}\right)$ is composed by a series of contributions according to the properties of the magnetic material. In our micromagnetic system, consisting of a single-domain magnetic nanoparticle of magnetite, contributions to the energy density are exchange, demagnetization, Zeeman and magnetocrystalline anisotropy [21]:

$$w\left(\mathbf{M}\right)=-\frac{A}{M_s^2}\mathbf{M}·{\nabla }^2\mathbf{M}-\frac{1}{2}{\mu }_0\mathbf{M}·{\mathbf{H}}_d-{\mu }_0\mathbf{M}·\mathbf{H}-\frac{K_1}{M_s^2}{(\mathbf{M}·\mathbf{u})}^2$$

(1)

where A is the stiffness constant accounting for an effective ferromagnetic interaction, ${\mathbf{H}}_d$ is the demagnetizing field caused by magnetostatic interaction, ${\mu }_0=1.2566\times {10}^6$ JA${}^{-2}$m${}^{-1}$ is the permeability of free space, $\mathbf{H}$ is the magnetic field, K is the anisotropy constant and $\mathbf{u}$ denotes the magnetic easy-axis unitary direction. In this work, the magnetic field is a sinusoidal wave applied along the $\mathbf{z}$-direction of the form $\mathbf{H}=H\left(t\right)\widehat{z}$ with $H\left(t\right)=H_{0}sin\left[2\pi f(t-t_0)\right]$, amplitude $H_0$, frequency f and initial displacement $t_0$. The competition between the different micromagnetic energy contributions on minimization, determines the equilibrium distribution of magnetization [22]. At each point or discretization cell of the magnetic sample, an effective field ${\mathbf{H}}_{\mathrm{eff}}$ is calculated by means of:

$${\mathbf{H}}_{\mathrm{eff}}=-\frac{1}{{\mu }_0}\frac{\partial E\left(\mathbf{M}\right)}{\partial \mathbf{M}}$$

(2)

with

$$E\left(\mathbf{M}\right)={\int }_{V}w\left(\mathbf{M}\right)dV$$

(3)

being the energy of the system over the entire volume V of the sample.

2.2. Landau– Lifshitz–Gilbert (LLG) equation

The time evolution of the magnetization of a magnetic nanoparticle, with radius R and saturation magnetization $M_s$, is ruled by the Landau-Lifshitz equation with Gilbert damping, which reads as follows [16,17,18],

$$\frac{d\mathbf{M}}{d\mathbf{t}}=-\frac{\gamma }{1+{\alpha }^2}\left(\mathbf{M}\times {\mathbf{H}}_{\mathrm{eff}}\right)-\frac{\gamma }{1+{\alpha }^2}\frac{\alpha }{M_s}\mathbf{M}\times \left(\mathbf{M}\times {\mathbf{H}}_{\mathrm{eff}}\right)$$

(4)

where $\alpha$ denotes the dimensionless damping parameter and $\gamma =2.2128\times {10}^5$ mA${}^{-1}$s${}^{-1}$ is the gyromagnetic ratio. LLG equation is strictly valid only at absolute zero temperature and the normalized magnetization vector is $\mathbf{m}=\mathbf{M}/M_s$. The system transfers energy and angular momentum from the movement of the magnetic moments to other degrees of freedom through the $\alpha$ parameter[23]. The first term of LLG equation is associated with the gyromagnetic ratio that promotes a uniform precession of the magnetization around the effective field. The second term is a phenomenological damping parameter that accounts for the energy loss of the system before reaching the equilibrium state.

2.3. Metastable lifetime and Specific Loss Power (SLP)

Ferromagnetic systems attached to a periodically oscillating magnetic field can exhibit order-disorder DPT. A significant aspect in DPT is the metastable lifetime. To study this DPT behavior, first the system must be fully magnetized with all the magnetic moments aligned along the direction of the maximum external magnetic field. Once the field is allowed to change in time, we analyze the conditions under which the magnetization becomes reversed and it passes through the zero value of magnetization. The elapsed time between these two magnetic states corresponds to the metastable lifetime and it is schematically shown in Figure 2. This is a measure of the time it takes for the system to escape from the metastable saturation state in the free-energy landscape.

A complete decay of the metastable phase occurs only in a dynamically disordered phase, as shown in Figure 1b. In this phase, the metastable lifetime is comparable with the period of the external magnetic field, the magnetization follows the magnetic field and its average over one period is close to zero. Consequently, the hysteresis loop is symmetric. The area enclosed by the hysteresis loop is the irreversible work dissipated in form of heat, which is quantified by the SLP energy. This heat is given by the following relationship [24,25]:

$$SLP=\frac{{\mu }_{0}f}{\rho }\oint M\left(t\right)dH\left(t\right)$$

(5)

where $\rho$ is is the mass density of the particle. The integration is done over one period of the oscillating magnetic field. The thermal energy quantified by SLP is an important magnitude for magnetic hyperthermia studies[3]. SLP depends on several parameters, such as $H_0$ and f of the external field, composition and system size, shape and magnetic state of the MNP [8].

2.4. Computational model

The simulation and the dynamic response of the magnetization of our micromagnetic system is carried out by means of the Ubermag micromagnetic software package [19,20]. The discrete representation of the magnetization distribution is obtained through the finite-difference method as is schematically represented in Figure 3a, where we have considered a magnetic nanoparticle with pseudo-spherical geometry and a crystal structure corresponding to magnetite ${\mathrm{Fe}}_3{\mathrm{O}}_4$. The physical properties used in our simulation for magnetite are shown in Table 1.

A discretization cubic cell of 2 nm edge was employed. This linear dimension is much smaller than the exchange length ($l_{ex}$) as required, which in the case of magnetite is given in Table 1. A typical value of $\alpha =0.07$ for a bulk ferrimagnetic material like magnetite was employed [28]. The frequency of the AMF was $f=1.0$ MHz, which is within the range of the radio-frequency currently used in magnetic hyperthermia [29]. Additionally, the initial displacement considered was $t_0=-{\left(4f\right)}^{-1}$.

In this work, initially the system is prepared in a fully positively magnetized state along the $\widehat{z}$-direction and different orientations of the magnetic easy-axis $\mathbf{u}$ of anisotropy are considered. In Cartesian coordinates, we have $\mathbf{u}(\theta ,\phi )=sin\left(\theta \right)cos\left(\phi \right)\mathbf{x}+sin\left(\theta \right)sin\left(\phi \right)\mathbf{y}+cos\left(\theta \right)\mathbf{z}$, where the azimuthal angle was fixed at $\phi ={90}^{\circ }$, and the different values of the polar angle $\theta$ were ranged between $2.5^{\circ }$ and $87.5^{\circ }$ .

3. Results and Discussion

The magnetization in the micromagnetic model is calculated by solving the LLG equation at zero temperature. First, we calculate the $z-$component of magnetization for amplitudes of $H_0=5$ mT and $H_0=200$ mT, as is shown in Figure 4a–c and Figure 4c–e respectively. Simulations were performed for ${\mathbf{u}}_{\theta }={\mathbf{u}}_{7.5^{\circ }}$, ${\mathbf{u}}_{\theta }={\mathbf{u}}_{{45}^{\circ }}$ and ${\mathbf{u}}_{\theta }={\mathbf{u}}_{82.5^{\circ }}$. First, note that for an amplitude of the external AMF of $H_0=5$ mT, and for each orientation of the easy-axis, the system is found in a DOP. In contrast, a DDP is evidenced for $H_0=200$ mT. Hence, the results shown in Figure 4 demonstrate that both DOP and DDP depend on the $H_0$ of AMF and the easy-axis orientation. A sharp reversal of magnetization takes place for $\theta$ values small, i.e. $\theta =7.5^{\circ }$.

Now that it is known the dependence of DOP and DDP with the amplitude of the external AMF, we can calculate the critical amplitudes $H_{0C}$ for which the DPT of the type DOP-DDP occurs for a wide range of values of the polar angle $\theta$ as shown in Figure 5. To achieve this dynamic magnetic phase diagram, it was necessary to carry out systematically several simulations like those shown in Figure 4, in order to determine the critical amplitudes $H_{0C}$ of the AMF for which the magnetization was not able to be switched. As can be observed in this figure, the lower region corresponds to the DOP whereas the upper region indicates the DDP. Thus, the DPT is strongly dependent on both the amplitude and the easy-axis orientation of the particle, and both quantities follow a non-linear relationship. Moreover, as the easy-axis is oriented nearly-parallel to the magnetic field direction, DPT takes place at higher values of amplitude. In contrast, for an orientation nearly-perpendicular of ${\mathbf{u}}_{\theta }$ DPT occurs at lower amplitude values. Thus, a nearly perpendicular orientation of ${\mathbf{u}}_{\theta }$ makes the magnetic system less rigid, giving rise to a DPT at small magnetic fields. This fact constitutes a way of tuning the magnetically hard or soft character of the system.

As larger amplitudes of the external AMF are considered, the magnetization quickly escapes from its metastable state and $\tau$ slowly decreases as ${\mathbf{u}}_{\theta }$ is oriented nearly-perpendicular to the external AMF direction. In contrast, the small the amplitude, more time is required to escape from the metastable state, and $\tau$ quickly decreases with the increase of $\theta$. Thus, we conclude that the metastable lifetime $\tau$ not only depends on the amplitude of the external AMF at a given frequency, but also on the easy-axis orientation. Several numerical results for different orientations of ${\mathbf{u}}_{\theta }$ were obtained as is shown in Figure 6.

In addition, we also obtained the dynamic hysteresis loops for different easy-axis orientations of at a given field amplitude. The respective results are shown in Figure 7 for $H_0=200$ mT. Hence, as already unveiled, it is clear that the angle $\theta$ of the easy-axis strongly affects magnetization reversal. As long as the easy-axis is oriented nearly perpendicular to the external AMF main direction, coercivity is significantly reduced. This presumably implies that as the external field is applied in a direction away from that of the easy-axis, it is easier to "pull" the magnetization out of the easy direction of anisotropy through different paths in phase space. This fact is also the responsible for the crossing of branches observed in the descending and ascending branches of the hysteresis loops at the switching fields, which is more noticeable for large values of $\theta$. This peculiar and surprisingly behavior of the crossing of hysteresis branches has been already demonstrated to occur in similar systems [30,31]. The crossing of the ascending and descending branches of the magnetization is observable when the applied external AMF is considerably far of the easy-axis, and this fact is associated with the energy minimum around the metastable state. Mathews et al. in [30] demonstrated that such scenario can be occur for a Stoner–Wohlfarth particle with a unique anisotropy. Hence, we have been able to show that the crossing of branches is significantly affected by the orientation of the particle.

On the other hand, the micromagnetic energy contributions for each orientation of the easy-axis are shown in Figure 8. Figure 8a–c show the exchange energy for ${\mathbf{u}}_{7.5^{\circ }}$, ${\mathbf{u}}_{{45}^{\circ }}$ and ${\mathbf{u}}_{82.5^{\circ }}$, respectively. Figure 8d–f stand for the demagnetization energy, Figure 9a–c for the anisotropy energy and Figure 9d–f for Zeeman energy, and for the same angles. We observe that the main micromagnetic energy contribution is given by the Zeeman energy, whereas the exchange energy contribution is the smallest one.

For the case of ${\mathbf{u}}_{82.5^{\circ }}$, a zoom of the anisotropy and Zeeman energies is shown for clarity in Figure 9. Some jumps at $H\left(t\right)/H_0=\pm 0.403$, $H\left(t\right)/H_0=\pm 0.409$ and $H\left(t\right)/H_0=\pm 0.414$ are observed for both the ascending and descending field branches, which is consistent with a sort of transition zone. In particular, a significant increase in the value of the anisotropy energy is observed for $H\left(t\right)/H_0=\pm 0.403$, which implies that the magnetic moments are not aligned with the easy-axis of the MNP. At the same value of $H\left(t\right)/H_0$, a minimum in the Zeeman energy is observed, which implies an alignment of the magnetic moments with the field direction. In contrast, for $H\left(t\right)/H_0=\pm 0.409$, the opposite scenario takes place. This means, that during the transition zone, magnetic moments become magnetically frustrated, as a consequence of the competition between these two energies.

Finally, all of the above results lead to ask about the amount of heat that can be dissipated in terms of the specific loss power, which is proportional to the hysteresis loop area. Such quantity was computed using Equation 5, and the respective angular dependence is shown in Figure 10. The SLP is maximized when the easy-axis of the MNP is oriented nearly parallel to the magnetic field direction, that is, for small angle values there is a greater efficiency in the amount of heat released to the environment. This is a novel result, because in the literature it is known that the specific loss power depends on parameters such as the size of the particles, the frequency and amplitude of the external AMF [8,9], but to the best of our knowledge, no studies on this specific regard have been addressed up to now. This result may shed lights on new ways to maximize the SLP for hyperthermia purposes.

4. Conclusions

In this study, we have presented numerical calculations of the micromagnetic properties and dynamic critical behavior of the magnetization for a single-domain MNP of magnetite under the influence of an oscillating magnetic field. Exchange, demagnetization, anisotropy and Zeeman energies were considered and the time evolution of the magnetization was obtained by solving the LLG equation at zero temperature. The numerical simulation was performed by means of the micromagnetic simulator Ubermag package based on the Object Oriented Micromagnetic Framework (OOMMF). We studied the effect of the orientation of the easy-axis of the MNP upon the critical values of the amplitude of the external AMF dealing with a dynamic phase transition. Thus, results were summarized in a proposal of a dynamic magnetic phase diagram. In addition, we calculated the time it takes for the magnetization to achieve the first-passage through the zero magnetization value from a saturated state, called the metastable lifetime, for different angles between the easy-axis and the direction of the external applied field. For a fixed amplitude $H_0=200$ mT, we obtained the hysteresis loops and calculated the specific loss power to know the amount of heat that can be potentially released to the medium. Our results show that the dynamic phase transition, the metastable lifetime and the specific loss power, strongly depend on the amplitude of the external field, and the orientation of the particle. Furthermore, our results show that the crossing of hysteresis branches are affected also by such orientation. Finally, SLP results may shed lights on new ways to maximize this quantity for hyperthermia purposes by tuning the orientation of the particles in the system.