Abstract

This study analytically analyzes the changes in the temperature profile of a homogenous and isotropic medium having the same thermal parameters as a muscular tissue, due to the heat released by a single magnetic nanoparticle (MNP) to its surroundings when subject to different magnetic field profiles. Exploring the temperature behavior of a heated MNP can be very useful predicting the temperature increment of it immediate surroundings. Therefore, selecting the most effective magnetic field profile (MFP) in order to reach the necessary temperature for cancer therapy is crucial in hyperthermia treatments. In order to find the temperature profile caused by the heated MNP immobilized inside a homogenous medium, the 3D diffusive-heat-flow equation (DHFE) was solved for three different types of boundary conditions (BCs). The change in the BC is caused by the different MF profiles (MFP), which are analyzed in this article. The analytic expressions are suitable for describing the transient temperature response of the medium for each case. The analysis showed that the maximum temperature increment surrounding the MNP can be achieved by radiating periodic magnetic pulses (PMPs) on it, making this MFP more effective than the conventional cosine profile.

1. Introduction

Magnetic Hyperthermia (MH) is one of many approaches currently being tested for cancer therapy [1–3]. The goal of this approach is to specifically heat the regions containing the cancerous cells by means of the magnetic losses caused by the physical properties of the magnetic nanoparticles (MNPs) when being exposed to an external magnetic field (MF).

The MNPs that are often used in MH, are usually made of ferromagnetic or ferrimagnetic materials which strongly react to the externally applied MF [4]. This magnetic reaction is converted by the two dominant relaxation mechanisms, the NΓ©el mechanism and Brownian mechanism, into power dissipation or heat [1, 4].

The eddy currents losses contribution may be neglected due to the low electrical conductivity that characterizes ferro- or ferrimagnetic materials and due to the small particle radius [5–9].

Fannin et al. [10] pointed out that for small enough particles, the anisotropy energy barrier, πΈπ‘Ž, may become so small, that thermal energy fluctuations can overcome it and spontaneous reverse the magnetization of a particle from one easy direction to the other, even in the absence of an applied field. The time it takes for a spontaneously fluctuation to reverse the magnetization after overcoming the energy barrier is characterized by a time constant, referred to as the NΓ©el relaxation time, or πœπ‘. The probability of such a transition is proportional to exp(𝜎), where 𝜎 is the ratio of anisotropy energy to thermal energy or (πΈπ‘Ž/π‘˜π΅π‘‡) [11].

The other distinct mechanism, by which the magnetization of MNPs may relax after an applied field has been removed, is the physical rotational Brownian motion of the particle immobilized inside a medium. When a magnetic field is applied to MNPs, they rotate and progressively align with the magnetic field due to the torque generated by the interaction of the magnetic field with the magnetization [12]. The time associated for an MNP to align with a small external magnetic field is given by the Brownian relaxation time, 𝜏𝐡 [13].

Because these relaxation mechanisms happen simultaneously, they both contribute to the total magnetization and the heat losses and their total influences can be express by an effective relaxation time, πœπ‘’, which is a combination of πœπ‘ and 𝜏𝐡 [14–16].

The two relevant mechanisms to change the magnetization of magnetic particles in an external field are given in Figure 1.

Moreover, our interests in MNPs as heat sources derive from the fact that they are vastly used as MRI agents [17] and their controllable size, ranging from few nanometers to tens of nanometers [18]. This means that the MNP size is smaller than or as same as that of a protein (5–50 nm), a virus (20–450 nm) or a gene (2 nm wide and 10–100 nm long) [11], which enables them to penetrate through the leaky pathological vasculature into the tumor interstitial easily reaching any cell of interest in the body including cancerous cells [19].

In addition, the MNPs can be attached to a specific type of cancerous cells causing a controllable elevation of temperature in them with almost no effects on healthy cells [20, 21]. By selectively heating the cancerous cells to a temperature ranging from 42 to 46Β°C, one can damage the tumors without causing vast harm to the healthy surrounding tissue [17, 19, 22].

Furthermore, in order to ensure that the treatment is biologically noninvasive and thermally tolerated for extended period of time, an experimentally measured tolerable limit of induced heating by an alternative MF was conducted defining a limit to the product of the MF strength (𝐻) and the frequency (𝑓) of the MF (e.g., 𝐻⋅𝑓≀4.85β‹…108 A/mΒ·s [23] or a less rigid criterion 𝐻⋅𝑓≀5β‹…109 A/mΒ·s [24]).

Due to the MNPs submicron length size, the conventional approach to heat conduction problems using macroscopic empirical laws such as Fourier’s law or Joule’s law of heat generation requires justification and even breakdown when the length scale of the system is comparable to the carrier mean free path or when the time scale of the physical process is smaller than the relaxation time of the heat carriers [25]. In this case, transport of heat carriers must be treated using the Boltzmann transport equation as Chen et al. pointed out [26].

Chen [27] suggested that heat is transported by carriers comprising of electrons, phonons, and photons. In dielectric materials, the heat conduction is dominated by phonons, in pure metals predominantly by electrons, and in impure metals or alloys by both phonons and electrons [27, 28]. Therefore, the mean free path of the heat carriers, for an MNP with a Fe core is approximately 0.8 nm [29, 30] and for a biological tissue 0.5 nm [26, 31, 32] allowing the conventional approach to be used for particles having a radius bigger than 10 nm.

Consequently, the temperature gradient caused by the release of the magnetic energy which an MNP absorbs to its immediate surroundings can be found analytically when applying Fourier transforms (FTs) to the DHFE [33] as Shih et al. [34] and Yuan et al. [35] suggested. Using this technique, Liu and Xu [36] analyzed the influences that a sinusoidal heat flux source placed on the skin surface have on the temperature inside it, and Tjahjono [37] analytically analyzed the heating temperature of a slab embedded with gold NPs due to a constant magnetic flux.

By analytically solving the DHFE for different boundary conditions, one can easily describe the dependence of the solution on each parameter composing it, such as the radius of the MNP, the frequency of the MF, and the material properties [38]. This allows us to optimize the solution for better performances, reaching the highest temperature elevation under specific constrains, for example the radius of the MNP. Moreover, when exploring the solution analytically, other parameters and their influence may be observed more clearly, which are usually neglected, or not explored (e.g., the MFP and its effects on the temperature gradient).

Until recently, the MFP was poorly analyzed in context of hyperthermia treatments and how it influences the temperature distribution concerning biological materials and tissues surrounding MNPs. Previous work focused on exploring the influences of different magnetic profiles on biological tissues. These studies were mostly experimental and did not focus on the MNPs contribution to the temperature change when exposed to different types of magnetic field profiles [39–43].

Recently, a numerical simulation model based on the Landau-Lifshitz-Gilbert equation was created for simulating MNPs ensembles when exposed to an incident square wave [44], as opposed to the usual sine wave. This work showed an increase in the normalized heat released by MNPs by at least 50%, as well as a more constant heating efficiency over the spectrum of particle anisotropies due to the infinite number of harmonics contained in an ideal square wave with the possibility of much greater improvement depending on the magnetic anisotropies, volumes, and angles to the incident radiation. However, Morgan and Victora [44] did not elaborate on the temperature dependencies on space and time near the MNP surface, but mostly focused on the dependencies between the angle of the incident wave relative to the anisotropy axis of the MNPs and the magnitude of the normalized output power released from them.

Therefore, the primary aim of this paper is to explore the transient analysis of the changes in temperature (from the steady state temperature 𝑇𝑏=310.15Β°K), as a function of the external MFP applied to a single MNP. By doing so, one can select the most efficient MFP that may improve the efficiency of MH treatments allowing the MF strength and the frequency reductions in order to meet the requirement 𝐻⋅𝑓≀4.85β‹…108 A/mΒ·s [23].

The aims of this paper are as follows:(i) To construct a theoretical model of the magnetic losses for the three different MF profiles studied in this article as follows:(1) Case 1β€”a cosine profile [18, 37],(2) Case 2β€”a PMP [45, 46],(3) Case 3β€”a discontinuous cosine profile.(ii)To explore the maximal temperature elevation and the rate of the temperature change near the MNP’s surface and into the tissue-like surrounding it for each of the above cases. (iii)To investigate who the core radius influences the maximal rate of the temperature change and the maximal temperature value, in order to find the optimal core radius that should be used for each of the above cases [5, 47].(iv)To study the effective confining heat depth (ECHD), symbolized as 𝛿 (see Figure 2), where the temperature elevation has a significant influence for each of the above cases.

2. Methods

In this study, we analytically model the transient temperature field (TTF) produced by a single MNP inside a homogenous and uniform medium having the same thermal parameters as a cancerous muscle cell. The analysis for each of the three MFPs mentioned in Section 1 is presented, after solving the DHFE inside the medium surrounding the MNP with the proper BC corresponding to its specific MFP.

In order to simplify the solution of DHFE that gave us the TTF and the temperature rate change due to the magnetic losses, some assumptions were made:(a)The properties of the surrounding medium are constant and homogeneous having the same thermal properties as the macroscopic-scale muscular tissue [48]. (b)The temperature on the surface of the MNP is uniform.(c)There is a negligible emission and evaporation.(d)There are no β€œthermally significant” blood vessels near the zone of interest therefore the perfusion is negligible.(e)The metabolic-heat generation is neglected.

2.1. The Thermodynamic Analysis

The TTF originating from the surface of a single MNP can be found after solving the 3D DHFE in the homogenous medium surrounding it [39]. The general DHFE can be written in spherical coordinates (due to the problem’s symmetry) as follows: π‘˜π‘šβˆ‡2π‘‡π‘š(π‘Ÿ,𝑑)=πœŒπ‘šπ‘π‘šπœ•π‘‡π‘š(π‘Ÿ,𝑑),πœ•π‘‘(1) where πœŒπ‘š (kg mβˆ’3) is the mass density, π‘π‘š (J kgβˆ’1Β°Cβˆ’1) the specific heat, and π‘˜π‘š (Wmβˆ’1Β°Cβˆ’1) the thermal conductivity of the phantom tissue.

This equation was also used by Keblinski et al. [38] and Govorov et al. [49] for solving nanoscale heat problems.

The general BC for this heat problem is given as follows: βˆ’π‘˜π‘šβ‹…βˆ‡π‘‡π‘š||(π‘Ÿ,𝑑)π‘Ÿ=π‘Ž=π‘žπ‘ ξ…žξ…ž(𝑑),(2) where π‘Ž: is the radius of the MNP in meters and π‘žπ‘ ξ…žξ…ž(𝑑) (Wmβˆ’2) is the heat flux.

The DHFE is valid if the mean free path of the heat carrier phonon or electron is smaller than the characteristic feature size as mentioned in Section 1. For amorphous solids and liquids, due to lack of crystalline, the mean free path is very short and of the order of atomic distances. Consequently, the heat flow in cells can be well described by the diffusive heat equation, even when nanoscopic length scales are involved [31].

Based on the above considerations, we evaluate the temperature field arising from continuous heating of a single particle by solving the heat equation (1) in the region outside the solid sphere surrounding the MNP where there are no heat sources, using a constant heat-flux-boundary condition at the particle surface, caused by the magnetic losses inside the MNP. The constant heat flux from the MNP’s surface becomes heat input to the medium which is then stored within the volumetric penetration depth region as shown in Figure 2.

After solving (1) and (2) (see the detailed formulations in Appendix A (A.1)–(A.12)), the temperature elevation inside the medium surrounding the MNP can be expressed using (A.12) as follows: Ξ”π‘‡π‘š(π‘Ÿ,𝑑)=πœƒ(𝑅,𝑑)π‘˜π‘šπ‘Ÿβˆšπ›Όπ‘š=πœ™(𝑅,𝑑)βˆ—π‘žπ‘ ξ…žξ…ž(𝑑)π‘˜π‘šπ‘Ÿβˆšπ›Όπ‘š,(3) where βˆ— symbolizes the convolution between two functions, βˆšπ‘…=1/π›Όπ‘šπ‘Ÿ and π›Όπ‘š=π‘˜π‘š/πœŒπ‘šπ‘π‘š.

In order to analytically calculate (3), the general expression of π‘žπ‘ ξ…žξ…ž(𝑑) must be found for each case, which depends on 𝐻(𝑑) (Amβˆ’1), the magnetic field, and on 𝑀(𝑑) (Amβˆ’1), the magnetization.

When a linear and isotropy material is assumed, the relation between 𝑀(𝑑) and 𝐻(𝑑) in the frequency domain (using FTs) may be described by the magnetic susceptibility [5] ξ‚‹ξ€œπ‘€(πœ”)=π‘‘π‘‘ξ…žπœ’ξ€·π‘‘ξ…žξ€Έπ‘’βˆ’π‘–πœ”π‘‘β€²ξ€œπ‘‘π‘‘π»(𝑑)π‘’βˆ’π‘–πœ”π‘‘ξ‚=ξ‚πœ’(πœ”)𝐻(πœ”).(4) The magnetic susceptibility ξ‚πœ’(πœ”) in the frequency domain can be expressed as [4] πœ’ξ‚πœ’(πœ”)=01+π‘–πœ”πœπ‘’=πœ™πœ‡0𝑀𝑠2𝑣𝑝3π‘˜π΅π‘‡11+π‘–πœ”πœπ‘’,(5) where πœ’0 is the static susceptibility, πœπ‘’=𝜏Neel||𝜏Brown is the effective relaxation time given by Fannin [14], πœ™ is the volume fraction solid [4], 𝑣𝑝 the particle volume, πœ‡0 the vacuum permeability, π‘˜π΅ is the Boltzmann constant, and 𝑀𝑠 is the magnetic saturation.

Moreover, in order to calculate the total heat generated by a single MNP caused by the conversion of the absorbed magnetic energy to heat inside a linear ferromagnetic medium, we must introduce Poynting’s theorem for small electric fields and neglecting ohmic losses [5, 6] as follows: βˆ‡β‹…π‘†ξ…žξ…žoutξ€œξ€œ(𝑑)=βˆ’π‘‘πœ”π‘‘πœ”ξ…žπ»(πœ”)⋅𝐻(πœ”)πœ”πœ‡0Im(πœ‡(πœ”))𝑒𝑖(πœ”βˆ’πœ”β€²)π‘‘βˆ’πœ•π‘ˆ(𝑑)πœ•π‘‘=βˆ’π‘ƒLoss(𝑑)βˆ’πœ•π‘ˆ(𝑑),πœ•π‘‘(6) where π‘†ξ…žξ…žout(𝑑) represents the energy flowing out through the boundary surfaces of the volume per unit time, 𝐻(πœ”) is the conjugate of 𝐻(πœ”),πœ‡0=4πœ‹10βˆ’7 (Vs/Am) is the vacuum permeability, πœ‡(πœ”)=πœ‡π‘Ÿ(πœ”)βˆ’π‘–πœ‡im(πœ”)=πœ‡0(1+ξ‚πœ’(πœ”)) is the complex magnetic permeability [5, 6], Im() is the imaginary part of πœ‡(πœ”), πœ•π‘ˆ(𝑑)/πœ•π‘‘ is the time rate change of the effective electromagnetic energy density given by (7), and 𝑃Loss(𝑑) represents the conversion of the magnetic energy into heat not counting conductive losses. It is worth mentioning that only the imaginary part of the complex permeability is causing heat losses, and πœ•π‘ˆ(𝑑)/πœ•π‘‘ can be found using [5, 6] as follows: π»β‹…πœ•π΅=ξ€œξ€œπœ•π‘‘π‘‘πœ”π‘‘πœ”ξ…žπ»(πœ”)⋅𝐻(πœ”)πœ”πœ‡0Imπœ‡(πœ”)𝑒𝑖(πœ”βˆ’πœ”β€²)𝑑+πœ‡0πœ•ξ€œξ€œ2πœ•π‘‘π‘‘πœ”π‘‘πœ”ξ…žπ‘‘π»(πœ”)⋅𝐻(πœ”)ξ€Ίπœ”π‘‘πœ”ξ€»πœ‡(πœ”)×𝑒𝑖(πœ”βˆ’πœ”β€²)𝑑=𝑃Loss(𝑑)+πœ•π‘ˆ(𝑑),πœ•π‘‘(7) where πœ‡(πœ”) is the conjugate of πœ‡(πœ”).

Next, the explicit analytic expressions for the temperature gradient profile from the equilibrium tempertaure, Ξ”π‘‡π‘š(π‘Ÿ,𝑑) and πœ•π‘‡π‘š(π‘Ÿ,𝑑)/πœ•π‘‘ are deduced for three different BC derived from the three MFPs mentioned earlier in Section 1.

2.2. The Analytic Expressions of the TTF for Three Different MFPs

Case 1 (a cosine profile). The magnetic field has a cosine profile so π»ξ€·πœ”(𝑑)=𝐴cos0𝑑.(8) Taking the inverse FT of (A.19) deduced in Appendix A using (A.13)–(A.19), one can find that πœƒ(𝑅,𝑑) can be written in this case as πœƒ(𝑅,𝑑)=π‘Žπœ‡0𝐴2πœ”06πœ”0πœ’0πœξ€·πœ”1+0πœξ€Έ2β‹…βŽ›βŽœβŽœβŽπ‘Žπ‘…0𝑒𝑖2πœ”0𝑑2𝑅0βˆšπ‘–2πœ”0ξ‚ξ‚€βˆ’βˆš+1exp𝑖2πœ”0𝑅+π‘Žπ‘…0π‘’βˆ’π‘–2πœ”0𝑑2𝑅0βˆšβˆ’2πœ”0ξ‚ξ‚€βˆ’βˆšπ‘–+1expβˆ’π‘–2πœ”0π‘…ξ‚βŽžβŽŸβŽŸβŽ ,+1(9) where πœƒ(𝑅,𝑑) is a function in the complex domain, therefore the temperature profile has a magnitude and phase as often occurs in many problems of physics or engineering such as theory of heat conduction, particularly when nonsteady heat conduction is concerned [50, 51]. Moreover, (9) is related to the TTF by (3).
Sometimes, the derivative of the temperature profile, or the rate of the change in the temperature surrounding the MNP is taken in consideration in order to verify that the treatment is safe for inducing controlled MH [47, 52]. For Case 1, this equals πœ•πœƒ(𝑅,𝑑)πœ•π‘‘=π‘–πœ”0π‘Žπœ‡0𝐴2πœ”06πœ”0πœ’0πœξ€·πœ”1+0πœξ€Έ2β‹…βŽ›βŽœβŽœβŽπ‘Žπ‘…0𝑒𝑖2πœ”0𝑑𝑅0βˆšπ‘–2πœ”0ξ‚ξ‚€βˆ’βˆš+1exp𝑖2πœ”0π‘…ξ‚βˆ’π‘Žπ‘…0π‘’βˆ’π‘–2πœ”0𝑑𝑅0βˆšβˆ’2πœ”0ξ‚ξ‚€βˆ’βˆšπ‘–+1expβˆ’π‘–2πœ”0π‘…ξ‚βŽžβŽŸβŽŸβŽ .(10)

Case 2 (a PMP). The magnetic field has a rectangular pulse shape profile with a period of 𝑇𝑠=2πœ‹/πœ”0 and a pulse width of Ξ”, so 𝐻(𝑑)=2𝐴0≀𝑑≀Δ,𝐻(𝑑)=0Δ≀𝑑≀𝑇𝑠.(11) The amplitude of the pulse wave was chosen to be twice the amplitude of the cosine single in order to maintain the same peak-to-peak value, for this case and the previous one. For this case, the temperature elevation as a function of time can be expressed using the inverse FT of (A.25) found in Appendix A that was deduced using (A.20)–(A.25) to receive the following: =π‘Žπœƒ(𝑅,𝑑)34πœ‡0𝐴2πœ”0πœ‹2sinπ‘šπœ‹Ξ”/π‘‡π‘ ξ€Έπ‘šξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπ‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2β‹…ξƒ©π‘Žπ‘…0𝑒𝑖(𝑛+π‘š)πœ”0𝑑𝑅0βˆšπ‘–(𝑛+π‘š)πœ”0ξ‚€βˆ’βˆš+1exp𝑖(𝑛+π‘š)πœ”0𝑅+π‘Žπ‘…0π‘’βˆ’π‘–(𝑛+π‘š)πœ”0𝑑𝑅0βˆšβˆ’π‘–(𝑛+π‘š)πœ”0ξ‚€βˆ’βˆš+1expβˆ’π‘–(𝑛+π‘š)πœ”0𝑅+π‘Žπ‘…0𝑒𝑖(π‘šβˆ’π‘›)πœ”0𝑑𝑅0βˆšπ‘–(π‘šβˆ’π‘›)πœ”0ξ‚€βˆ’βˆš+1exp𝑖(π‘šβˆ’π‘›)πœ”0𝑅+π‘Žπ‘…0π‘’βˆ’π‘–(π‘šβˆ’π‘›)πœ”0𝑑𝑅0βˆšβˆ’π‘–(π‘šβˆ’π‘›)πœ”0ξ‚€βˆ’βˆš+1expβˆ’π‘–(π‘šβˆ’π‘›)πœ”0𝑅ξƒͺ+π‘Ž34πœ”0πœ‡0𝐴2β‹…Ξ”πœ‹π‘‡π‘ ξ“ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπ‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2Γ—ξƒ©π‘Žπ‘…0π‘’π‘–π‘›πœ”0𝑑𝑅0βˆšπ‘–π‘›πœ”0ξ‚€βˆ’βˆš+1expπ‘–π‘›πœ”0𝑅+π‘Žπ‘…0π‘’βˆ’π‘–π‘›πœ”0𝑑𝑅0βˆšβˆ’π‘–π‘›πœ”0ξ‚€βˆ’βˆš+1expβˆ’π‘–π‘›πœ”0𝑅ξƒͺ.(12) Again, (12) is related to the TTF by (3).

As for Case 1, we can calculate the rate of the change in the temperature surrounding the MNP and receive the following: πœ•πœƒ(𝑅,𝑑)πœ•π‘‘=π‘–πœ”20π‘Ž34πœ‡0𝐴2πœ‹2sinπ‘šπœ‹Ξ”/π‘‡π‘ ξ€Έπ‘šξ‚΅Γ—sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπ‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2⋅(𝑛+π‘š)π‘Žπ‘…0𝑒𝑖(𝑛+π‘š)πœ”0𝑑𝑅0βˆšπ‘–(𝑛+π‘š)πœ”0ξ‚€βˆ’βˆš+1Γ—exp𝑖(𝑛+π‘š)πœ”0π‘…ξ‚βˆ’(𝑛+π‘š)π‘Žπ‘…0π‘’βˆ’π‘–(𝑛+π‘š)πœ”0𝑑𝑅0βˆšβˆ’π‘–(𝑛+π‘š)πœ”0ξ‚€βˆ’βˆš+1Γ—expβˆ’π‘–(𝑛+π‘š)πœ”0𝑅+(π‘šβˆ’π‘›)π‘Žπ‘…0𝑒𝑖(π‘šβˆ’π‘›)πœ”0𝑑𝑅0βˆšπ‘–(π‘šβˆ’π‘›)πœ”0ξ‚€βˆ’βˆš+1Γ—exp𝑖(π‘šβˆ’π‘›)πœ”0π‘…ξ‚βˆ’(π‘šβˆ’π‘›)π‘Žπ‘…0π‘’βˆ’π‘–(π‘šβˆ’π‘›)πœ”0𝑑𝑅0βˆšβˆ’π‘–(π‘šβˆ’π‘›)πœ”0ξ‚€βˆ’βˆš+1Γ—expβˆ’π‘–(π‘šβˆ’π‘›)πœ”0𝑅ξƒͺ+π‘–πœ”20π‘Ž34πœ‡0𝐴2Ξ”πœ‹π‘‡π‘ ξ“ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άβ‹…π‘›2πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2Γ—ξƒ©π‘Žπ‘…0π‘’π‘–π‘›πœ”0𝑑𝑅0βˆšπ‘–π‘›πœ”0ξ‚€βˆ’βˆš+1expπ‘–π‘›πœ”0π‘…ξ‚βˆ’π‘Žπ‘…0π‘’βˆ’π‘–π‘›πœ”0𝑑𝑅0βˆšβˆ’π‘–π‘›πœ”0ξ‚€βˆ’βˆš+1expβˆ’π‘–π‘›πœ”0𝑅ξƒͺ.(13)

Case 3 (a discontinuous cosine profile). The magnetic field has a periodic discontinuous cosine profile with a time constant of 𝑇𝑠=2πœ‹/πœ”1 and a pulse width of Ξ”, so π»ξ€·πœ”(𝑑)=𝐴cos0𝑑0≀𝑑≀Δ,𝐻(𝑑)=0Ξ”β‰€π‘‘β‰€π‘‡π‘ πœ”0β‰ πœ”1.(14) For this third case, the temperature elevation as a function of time can be expressed using the inverse FT of (A.35) found in Appendix A that was deduced using (A.30)–(A.35) to receive: =π‘Žπœƒ(𝑅,𝑑)3πœ‡0𝐴𝑇𝑠2π‘Žπ‘…0Γ—ξ€·ξ€·πœ”ξ“ξ“ξƒ¬ξƒ©sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒͺΓ—ξƒ©ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒͺβ‹…ξ€·π‘›πœ”1ξ€Έ2πœ’0πœξ€·1+π‘›πœ”1ξ€Έ2𝜏2×𝑒𝑖(𝑛+π‘š)πœ”1𝑑𝑅0βˆšπ‘–(𝑛+π‘š)πœ”1ξ‚€βˆ’βˆš+1Γ—exp𝑖(𝑛+π‘š)πœ”1𝑅+π‘’βˆ’π‘–(𝑛+π‘š)πœ”1𝑑𝑅0βˆšβˆ’π‘–(𝑛+π‘š)πœ”1ξ‚€βˆ’βˆš+1Γ—expβˆ’π‘–(𝑛+π‘š)πœ”1𝑅+𝑒𝑖(π‘šβˆ’π‘›)πœ”1𝑑𝑅0βˆšπ‘–(π‘šβˆ’π‘›)πœ”1ξ‚€βˆ’βˆš+1Γ—exp𝑖(π‘šβˆ’π‘›)πœ”1𝑅+π‘’βˆ’π‘–(π‘šβˆ’π‘›)πœ”1𝑑𝑅0βˆšβˆ’π‘–(π‘šβˆ’π‘›)πœ”1ξ‚€βˆ’βˆš+1Γ—expβˆ’π‘–(π‘šβˆ’π‘›)πœ”1𝑅.ξƒͺξƒ­(15) Again, (15) is related to the TTF by (3).

As for Cases 1 and 2, we can calculate the rate of the change in the temperature surrounding the MNP to receive the following: πœ•πœƒ(𝑅,𝑑)πœ•π‘‘=π‘–πœ”1π‘Ž3πœ‡0𝐴𝑇𝑠2Γ—ξ€·ξ€·πœ”ξ“ξ“ξƒ¬ξƒ©sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒͺΓ—ξƒ©ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·(ξ€·πœ”sinΞ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒͺβ‹…ξ€·π‘›πœ”1ξ€Έ2πœ’0πœξ€·1+π‘›πœ”1ξ€Έ2𝜏2×(𝑛+π‘š)π‘Žπ‘…0𝑒𝑖(𝑛+π‘š)πœ”1𝑑𝑅0βˆšπ‘–(𝑛+π‘š)πœ”1ξ‚€βˆ’βˆš+1Γ—exp𝑖(𝑛+π‘š)πœ”1π‘…ξ‚βˆ’(𝑛+π‘š)π‘Žπ‘…0π‘’βˆ’π‘–(𝑛+π‘š)πœ”1𝑑𝑅0βˆšβˆ’π‘–(𝑛+π‘š)πœ”1ξ‚€βˆ’βˆš+1Γ—expβˆ’π‘–(𝑛+π‘š)πœ”1𝑅+(π‘šβˆ’π‘›)π‘Žπ‘…0𝑒𝑖(π‘šβˆ’π‘›)πœ”1𝑑𝑅0βˆšπ‘–(π‘šβˆ’π‘›)πœ”1ξ‚€βˆ’βˆš+1Γ—exp𝑖(π‘šβˆ’π‘›)πœ”1π‘…ξ‚βˆ’(π‘šβˆ’π‘›)π‘Žπ‘…0π‘’βˆ’π‘–(π‘šβˆ’π‘›)πœ”1𝑑𝑅0βˆšβˆ’π‘–(π‘šβˆ’π‘›)πœ”1ξ‚€βˆ’βˆš+1Γ—expβˆ’π‘–(π‘šβˆ’π‘›)πœ”1𝑅.ξƒͺξƒ­(16) For Cases 2 and 3, there is a limitation regarding the MFS and the frequency in order for the MH treatment to be safe (see (A.30) and (A.39)). Moreover, for frequencies lower than 10 MHz, there is essentially no attenuation of the MFS within muscle-equivalent materials, limiting the maximal harmonic frequency to 10 MHz [16].

In conclusion, (9)–(16) can be used to predict the TTP and the special temperature profile for a single-MNP subject to three different magnetic field profiles and using the same equations we can also explore the influence that the core radius has on the temperature profile estimating the ECHD for each case.

2.3. The Simulations Parameters

The mathematical expressions of the TTP were simulated using MATLAB and COMSOL (COMSOL results can be seen in Appendix B) for a single MNP immobilized inside a uniform and isotropic phantom medium having the same biological thermal properties as a muscular tissue [48] and are summarized in Table 1. These assumptions were made in order to simplify the theoretical calculations.

The thermal parameters are considered to be constant with temperature and space as will be latter proven. Moreover, the magnetic parameters of the MNP were measured at 𝑇𝑏=310.15Β°K based on the findings of Fannin [14] and are summed up in Table 2.

The external magnetic field strength (MFS) for all three cases, was chosen as 8.8 kAmβˆ’1 and the MF frequency as 𝑓0=400KHz. These values are based on previous works made by Kettering et al. [52], Hergt et al. [53], and Hilger et al. [54].

For all three profiles mentioned in this section, the simulations were plotted for 0β‰€π‘Ÿβˆ’π‘Žβ‰€10 nm and 0≀𝑑≀5 μs, where π‘Ÿ is the distance from the center of the MNP and π‘Ž its radius. The upper value for distance simulation was chosen accordantly to the thickness of the cell membrane thatis about 5–10 nm [55–57] and damaging it can cause the destruction of the cell [58]. The upper time value was chosen so several cycles of the magnetic field could be simulated and plotted.

For all the simulations the volume fraction solid was defined as πœ™=0.032. This value is been justified in Section 4.

In Section 3, as already mentioned in Section 1, the maximal temperature elevation and the temperature change rate near the MNP’s surface and into the tissue surrounding it, are simulated. Moreover, the influence the core radius has on the maximal temperature change rate and on the maximal temperature elevation was also explored, in order to find the optimal value that must be used for each case as suggested by Rosensweig [4] and Kappiyoor et al. [47]. Furthermore, the ECHD that is defined by the point the temperature reaches π‘’βˆ’1 of the maximal value, has also been explored for each case, defining the confining heat region and can be compared with the thickness of a cell membrane which varies between 3–10 nm [53–56].

3. Results

The mathematical expressions of the TTP were simulated in this section, using MATLAB, for a single MNP immobilized inside a uniform and isotropic phantom medium having the same biological thermal properties as a muscular tissue (Table 1). Moreover, the MNP’s magnetic parameters are summarized in Table 2.

For Case 1, the mathematical expression of the temporal and spatial temperature increment, (9), is presented in Figure 3 for 𝑇𝑠=2.5 μs.

It can be seen from Figure 3 that the temperature changes periodically with a time period of 1.25 μs that is equivalent to a frequency oscillation of 800 kHz, which is twice the frequency of the external applied MF, as predicted by (9). This can be explained by the multiplication of the magnetic field and the magnetic induction, both being a function of 𝑓0 or πœ”0.

Moreover, the temperature increment reached its maximum value after 0 μs, reaching Δ𝑇max=2.1 nΒ°K on the surface of the MNP. As expected, the hottest spots are on the surface of the MNP, and as the point of view gets further from the surface the temperature declines as (9) predicted. This value causes only a low-temperature gradient in the thermal properties of the surrounding medium therefore, the thermal parameters of the phantom cell can be considered constants, as assumed.

According to Figure 3, the temperature profile has a β€œDC” level that can be found from calculating the first term of (9) making the temperature increment to be always higher than the initial temperature as expected, because the magnetic losses inside the MNP irradiate heat to the medium surrounding it, at all times [59, 60].

Furthermore, at a distance of 12 nm apart from the MNP’s surface the temperature maximal value equals 0.8 nΒ°K that is equivalent to π‘’βˆ’1 of the absolute maximal value, defining the ECHD or 𝛿=12 nm.

In order to have a unique quantity to be compared with each case and does not depend on time, we averaged the TTP over one cycle. In Case 1, the averaged value over one time period equals 0.55 nΒ°K.

Next, we explored the maximal temperature rise rate as a function of the core radius using (9) and (10) receiving the data in Figure 4.

As can be seen from Figure 4, the maximal temperature rate change equals 0.011Β°Ksβˆ’1 and the maximal temperature rise equals 4.7 nΒ°K, both received for a core radius of 9.3 nm. The temperature rate rise and the maximal temperature are considerably small, due to the relaxation time that depends on the volume of the particle making this MFP to be safe to use for MH treatments [47]. For radii larger or smaller than 9.3 nm, the magnetic heat dissipation start to decrease as the magnetic relaxation time gets bigger or smaller, respectively, reducing the denominator or numerator in (9) and (10).

Equations (9) and (10) enable us to understand that the changes in the temperature depend on many parameters such as, the magnetic field strength, the magnetic field frequency, the magnetic properties of the material, and the core radius. Consequently, in order to optimize the heat losses we must select the most effective radius for a specific type of MNP.

For Case 2, the mathematical expression of the temporal and spatial temperature increment are plotted in Figure 5 for the summation of 25 indexes (not to exceed 10 MHz [16]) and Ξ”=0.2𝑇𝑠, where 𝑇𝑠=2.5 μs.

For convenience, Figure 5 describes the temperature profile for the first two cycles as given by (12). This equation shows that the characteristic behavior of the temperature repeats itself every 𝑇𝑠=2.5 μs, that is, the cycle of the magnetic field, therefore one can limit the study to only a final number of cycles.

As can be seen from Figure 5, the temperature builds up very fast due to the steep elevation of the magnetic field caused by the Heaviside-shaped MP and reaches a maximal value of Δ𝑇max=8.8 nΒ°K on the MNP’s surface after 0.45 μs. Then the temperature begins to drop after 0.1 μs from the time the MF was turned on reaching a minimal value of 6 nΒ°K near the surface of the MNP. From that point on, the temperature profile temporal behavior is defined by the summation of the total numbers of harmonics composing the MF until the MF is turned off again as can be found from (12). Furthermore, the temperature reaches its maximum value close to the surface of the MNP and decreases with distance reaching a maximal value of 3 nΒ°K, 12 nm apart from the MNP’s surface.

For this case, the maximal value is 4 times higher than the one received in Case 1, making it a preferable MFP to be used in MH as Morgan and Victora suggested [44].

Again, the thermal parameters can be considered constant and not dependent on temperature near the MNP’s surface because the temperature rise is less than 1Β°K.

From Figure 5, the ECHD can be found as 𝛿=12nm, which is the same as the value received in Case 1, meaning that the temperature decreases as fast as the cosine case and is confined to a specific area near the MNP’s surface.

Moreover, in order to have a unique quantity to be compared in each case that does not depend on time, we averaged the total temperature rise over one cycle. In this case the averaged temperature elevation was 1.3 nΒ°K after been normalized to the time period. This value is about 2.4 times higher than the value received in Case 1, making this MFP a better candidate for MH treatments.

Next, we explored the maximal temperature rise rate as a function of the core radius. For Case 2, we can use (11) and (12) receiving the data in Figure 6.

As seen from Figure 6, the absolute maximal temperature elevation equals 0.032 μ°K received for a core radius of 8.3 nm, and the maximal temperature derivative 1.01 °Ksβˆ’1 is received for a core radius of 8.2 nm. Because this MFP produces temperature changes that are too rapid to be safe for inducing MH [47], the radius that we chosen for a safer treatment is in consistence with Case 1 and equals 10 nm. Consequently, the NP size plays an important role in determining the amount of heating that an MFH treatment can provide as Kappiyoor et al. [47] already mentioned.

Again, the maximal temperature rate rise and the maximal temperature are considerably small due to the relaxation time that depends on the volume of the particle. For radii larger or smaller than 8.4 nm, the magnetic heat dissipation starts to decrease as the magnetic relaxation time gets bigger or smaller, respectively, due to its affect on the relaxation time, reducing the denominator or numerator in (11) and (12).

By comparing Case 2 to Case 1 we see that for the same MNP radius (10 nm) having the same magnetic material properties (given by Table 2), the maximal temperature rise received is about 4 times higher for Case 2 than in Case 1 and the maximal temperature derivative for this case is 40 times higher than Case 1, making the periodic pulse-shaped MFP a better magnetic field source for MH treatments.

For Case 3 the analytic expression for the TTP can be plotted for Ξ”=0.2𝑇𝑠 and 𝑇𝑠=15 μs and are shown in Figure 7 for the summation of 25 indexes (not to exceed 10 MHz [16]). The cosine MF time period that multiplies the Heaviside function equals 2.5 μs, and is equivalent to a frequency of 400 KHz.

As Figure 7 shows, the changes in the temperature profile result from the MNP reaction to two different MFPs, the cosine profile and the periodic rectangular pulse profile. The last is responsible for switching on and off the MF.

The influence that the periodic rectangular-pulse-shaped MF has on the temperature gradient can be seen by the steep temperature elevation at the beginning and at the end of every cycle caused by the derivative of Heaviside function composing the magnetic flux density, 𝐡(𝑑), and the influence that the cosine MFP has on the temperature gradient can be seen as the cosine β€œripple” that is added. This β€œripple” has 3 peaks that are separated 1.25 μs apart, which is twice the cosine MF frequency received in Case 1. On the MNP surface, the maximal temperature gradient reaches Δ𝑇max=2.3 nΒ°K after 2 ns from the time the MF was turned on, and repeats itself every 15 μs, which is equivalent to the time period of the signal. This value is higher than the value received for the cosine-shaped MF, but lower than the one received in Case 2. However, after the highest peak the maximal value of the cosine β€œripple” reaches the same one as in the cosine case or 2.1 nΒ°K as expected.

For this case, the ECHD equals 𝛿=12nm that is the same as for the other two cases where the temperature change reaches a value of 0.8 nΒ°K. After 0.2 𝑇𝑠, the temperature elevation becomes insignificant as the MF is turned off.

Again, in order to have a unique quantity to be compared in each case that does not depend on time, we can average the total temperature rise over one cycle. In this case the average temperature elevation equals 0.235 nΒ°K after we normalize it to the time period. This value is about 2.5 times lower than the value received in Case 1, making this MFP a less preferable candidate for MH treatments.

Next we explored the maximal temperature rise rate as a function of the core radius. For Case 3, we used (15) and (16) receiving the data in Figure 8.

As seen from Figure 8, the absolute maximal temperature rate elevation equals 5.4 nΒ°K received for a core radius of 9.2 nm. The maximal temperature derivative 0.024Β°Ksβˆ’1 is received for a core radius of 9.1 nm. Again, the maximal temperature rate rise and the maximal temperature are considerably small due to the relaxation time that depends on the volume of the particle. For radii larger or smaller than 9.2 nm, the magnetic heat dissipation starts to decrease as the magnetic relaxation time gets bigger or smaller, respectively, reducing the denominator or the numerator in (15) and (16).

By comparing Case 3 to Case 1 for a core radius of 10 nm and the same magnetic material properties (given by Table 2), the maximal temperature rise received for Case 3 is about two times higher than Case 1, and the maximal temperature derivative for this case is 2.2 times higher than Case 1. However, due to the long period for which the MF is turned off, and consequently the lower heat released from the MNP over one cycle, this MFP is less preferable than the cosine MFP, for MH treatments.

In order to make it easier to understand the differences between the three cases analyzed in this paper, Table 3 is added that summarizes the most significant parameters.

Moreover, a summarizing figure, Figure 9, describing the temperature rise as a function of time is also added for a particle with a core radius of 10 nm when the observation point is on the MNP surface, in order to easily evaluate the differences in the three cases studied in this article.

4. Discussion

An analytical analysis of the TTP was preformed for three MFPs. The mathematical models were received by solving the DHFE for different BC matching each MFP, using the FTs.

Major work have been done in the past to solve the DHFE equation for a cosine-MF source as can be found in [36, 37, 49, 59]. Keblinski et al. [38] found that a laser source having a constant power of 1.4β‹…10βˆ’8 W heating a single MNP with a radius of 65 nm can cause a temperature change of 0.06 K at the particle surface. Moreover, for a cosine-MF heat source the local temperature was found to be even lower, causing a maximum change in temperature of 0.1 mΒ°K for a particle having a radius of 50 nm at a frequency of 2 MHZ [48]. Both results are negligible from the point of view of biological applications as expected.

However, Keblinski et al. [38] and others [4, 20] solved the DHFE equation only for a constant heat flux having the average power of a cosine-MF, without exploring the temperature temporal behavior. In addition, until now there has not been an explicit mathematical formulation that solves the DHFE equation for other periodic MFPs that can be used as radiation sources for MH treatments. Morgan and Victora [44] showed that the use of an incident square wave, as opposed to the usual sine wave, increases the normalized power heat by at least 50%, however this conclusion was based on calculating only the Poynting vector and not based the solving the DHFE in order to find the explicit temperature change.

In consequence to the above, we should explore the influences that different magnetic irradiation profiles have on the induced temperature gradients inside tumor cells, for the same physical and thermal MNP’s parameters in order to verify what Morgan and Victora [44] suggested. Furthermore, optimizing the heat power is of great importance from biological point of view. A typical cell having a diameter ranging from 2–10 μm [61] can uptake a maximal quantity of anionic MNPs, that varies between 2.8β‹…105 and 7.2β‹…106 per cell, consequently limiting the total amount of magnetic material per cell. Moreover, high concentration of MNPs with different types of coatings can cause a toxic reaction to the central nervous system [62] or may cause cellular perturbations [63] therefore, it is important to reduce the MNP’s concentration. Nevertheless, reaching these quantities in vivo proves to be a very difficult task all types of cancerous cells [64, 65]. Hence, one must optimize other parameters such as, the profile of the MF in order to use lower magnetic concentration in order to reach the same temperature gradient values.

Consequently, this paper focuses on the influences that three different MFPs have on the temperature surrounding a single MNP, as mentioned in Section 1, when being exposed to it, analytically proving to be the most effective one in causing the highest temperature rise, using the same magnetic and thermal parameters.

For all three cases, the MATLAB (in Section 3) and COMSOL (in Appendix B) simulations results showed, that the maximum temperature rise for a given core radius of 10 nm ranges between 2.1 nΒ°K and 8.8 nΒ°K depending on the MFP.

Similar results were received by Keblinski et al. [38] and Rabin [31] for a constant heat flux and an MNP’s having approximately the same physical and magnetic properties. The very low absolute change in temperature caused by a single MNP can be explained by its low magnetic susceptibility πœ’0, and by the effective relaxation time that changes drastically with the MNP’s volume [15]. Therefore, a single MNP can release only a small amount of heat causing a very small change in the temperature surrounding it.

However, for larger magnetic concentration occupying a single cell such as 1 ng of Fe3O4 per human cell, that is, equivalent to 5β‹…107 MNPs per cell (for a particle radius 10 nm and cell radius 5 μm), Linh et al. [66] and Balivada et al. [67] showed that a local temperature elevation of several degrees can be reached making MH treatments effective. In addition, these quantities were also proposed by Vera and Bayazitoglu [59], Chan et al. [61], Huang [64], and Melancon et al. [65] who proved their efficiency in inducing MH. Consequently, one can produce a significant global temperature increment inside a cell even if the local temperature increment of each particle is negligible as long as we heat many particles in the same volume of interest.

For the multiparticle case, Rabin [31], and Keblinski et al. [38] calculated the temperature rise inside a spherical region with radius 𝑅 (m) consisting of many randomly dispersed heat sources using Δ𝑇global(𝑑)=(𝑅2Δ𝑇nano(𝑑)/π‘Ÿπ‘2)(4πœ‹/3)π‘Ÿπ‘3πœŒπ‘, where πœŒπ‘ (mβˆ’3) is the number of MNPs per unit volume, π‘˜ (Wmβˆ’1Β°Cβˆ’1) is the thermal conductivity of the medium, π‘Ÿπ‘ (m) is the radius of MNPs, and Δ𝑇nano(𝑑) is the temperature gradient caused by a single MNP.

For πœŒπ‘=5β‹…1021 (mβˆ’3) and an average tumor cell radius of π‘Ÿcell=7 μm [52], the number of MNPs inside a single cell can be calculated to equal 8β‹…106 that fits the concentrations found by Linh et al. [66] and Balivada et al. [67]. From πœŒπ‘, we calculated the distances between two neighboring particles that is approximately 58 nm. This means that the volume fraction of the MNPs inside the cell is about 0.02.

By choosing a solid volume fraction of πœ™=0.032, the calculated distance between two neighboring particles is about 50 nm fitting a concentration of πœŒπ‘=8β‹…1021 mβˆ’3, that is, in the toxicity safety range, for a tumor cell having an average radius of 7 μm [52, 66].

Due to the large distances between the particles, we assumed that the interparticle interactions are negligible, so the relaxation time and magnetic susceptibility can be calculated using the same expression as (5) [68].

The total temperature increment for the three cases analyzed in this paper can be found by substituting the received values for the single-MNP case, (9), (12), and (15) into Δ𝑇global(𝑑) when average tumor radius of 𝑅=4 mm was assumed in consistence with magnitudes of cancer tumors [31, 38, 67].

For the cosine MFP an average value of Δ𝑇global_cos(𝑑)=2.9∘K over one cycle is received near the MNP’s surface. This means that the MNP’s concentration is not sufficient to give increment to a dramatic temperature gradient under the parameters summarized in Tables 1 and 2.

In order to receive a 6Β° increment, that is needed for MH, in the temperature near the MNP’s surface, a larger amount than the proposed of particles is required.

For the PMP, an average value of Δ𝑇global_pulsed(𝑑)=7.2∘K over one cycle is received for the same parameters summarized in Tables 1 and 2 that is sufficient to induce MH.

For the discontinuous pulse-shaped MF, a maximum peak of Δ𝑇global_pulsed_cos(𝑑)=1.23∘K over one cycle is received, meaning that the MNP’s concentration in this case as in Case 1 is not sufficient to ensure that MH can occur.

The comparison between the three temperature gradients received for each case, shows that the preferable MFP for MH is the PMP one, compared to Case 1 and Case 3. For Case 2, the temperature gradient at the surface of the MNP is sufficient to cause damage to biologic cells [58, 69, 70]. Therefore, using a periodic pulse MFP can reduce the necessary amount of MNPs by a factor or even more, allowing a wider range of markers to be used for hyperthermia treatments, and simplifying the biological processes to conjugate them to a cell.

In addition, we also explored the influence that the MNP’s radius has on the maximal temperature gradient an on its rate rise. As seen from Figures 4, 6, and 8, the NP size has a great influence on determining the amount of heat released from the MFP’s surface effecting both the temperature gradient as well as the temperature rise rate, as previous works showed [50, 51, 70].

For the first Case 1 studied, the optimal core radius was found as 9.3 nm where the maximal temperature reaches 4.7 nΒ°K and the temperature change rate equals 0.011Β°Ksβˆ’1 (Figure 4). This optimal radius was also received by Kappiyoor et al. [47] for almost the same MF properties and magnetic material properties. However, because the equation solved by Kappiyoor et al. [47] is different than (1), the maximal value is slightly lower that the values received by Rosensweig [4] and Kappiyoor et al. [47]. Moreover, the maximal value is also affected by the parameters chosen to describe the magnetic properties of the MNP, as demonstrated by Kappiyoor et al. [47]. Our magnetic parameters are slightly different than the ones used by Rosensweig [4] and Kappiyoor et al. [47], which may account for the differences in the maximal values in this study as Kappiyoor et al. [47] showed.

For Case 2 studied, the optimal core radius was found as 8.3 nm where the maximal temperature gradient reaches 32 nΒ°K, and the temperature change rate equals 1.1Β°Ksβˆ’1 for a summation of 25 indexes (Figure 6). As can be seen by comparison, there is a benefit in using a PMP rather than a cosine MF, due to the higher temperature gradient received in the MNP’s surrounding and the sufficient average temperature gradient received per cycle that is about 2.5 times higher in Case 2 than in Case 1.

Although for a total summation of 25 indexes, the temperature change rate is approximately 1Β°Ksβˆ’1, (suggested to be less safe [47]), one can reduce the received value by limiting the number of the summation indexes composing the MF to a lower number such as 𝑁=10 instead of 𝑁=25, making the treatment safer but also maintaining higher temperature values that in Case 1 (Figure 14). Furthermore, when looking at the results of multiplying each coefficient’s amplitude with its matched harmonic, the limitation for the treatment to be biologically noninvasive remains valid as long as 𝐴eff𝑓0≀5β‹…109 A/mΒ·s as mathematically justified in (A.29) and (A.30) limiting the total summation index to a value lower than 𝑁=25.

For the Case 3 studied, the optimal core radius was found as 9 nm where the maximal temperature reaches 0.64 μ°K and the temperature change rate equals 27.3Β°Ksβˆ’1 (Figure 8). Although the values are higher than the ones received for the first case, the average temperature elevation was found to be lower after normalizing it to the time period, making this type of MFP less preferable.

As we can see for all the three cases analyzed in this paper, the optimal radius depends very much on the magnetic material properties [47] and on the profile of the magnetic field as we have proven in Section 2. Therefore, for each case studied and for each magnetic material the equations, developed for Cases 1–3 must be solved separately in order to optimize the MH treatment.

Another interesting finding driven from the mathematical equations, is the confinement of the temperature to an area having an average radius of 12 nm from the MNP’s surface for all three cases. This means that most of the heat dissipation occurs in the vicinity of the heat sources confining the temperature increment in the proximity of tumor cells alone unaffecting the healthy cells.

The importance of this paper lies in the fact that, until now there was no explicit mathematical formulation that solves the DHFE equation for other types of periodic MFPs used as excitation sources for MH treatments. As we found out, changing the profile of the MF radiation can induce higher temperature gradients in tumor cells, for the same physical and thermal parameters enabling reduction of the MNPs concentration per cell. This is of great importance because, a typical cell has a maximal quantity of MNPs that it can uptake, and because high concentration of MNPs with different types of coatings can cause a toxic reaction to the central nervous system [62]. Therefore, lowering the magnetic concentration per cell, but still receiving the same temperature gradients may be of great use.

With the outcome of this paper, we are moving forward to in vitro studies in order to verify the theoretical results received in this paper experimentally.

5. Conclusions

This study investigates the effects of different heat-flux profiles on a single MNP immobilized inside a phantom having the same thermal properties as a muscle tissue. The exact solution of DHFE was solved for different boundaries conditions using FTs. According to the analytic solutions, the PMP profile was found to be the more effective in rising the temperature of the medium surrounding the MNP than the cosine profile, making it a better candidate for hyperthermia treatments rather than the conventional cosine MP.

Moreover, in order to reach a significant temperature gradient for all cases studied (a) a cosine profile, (b) a PMP profile, and (c) a discontinuous cosine profile, there is a need for a larger number of MNPs to be immobilized inside the cell medium as Rabin [31] and Keblinski et al. [38] previously suggested. Using their techniques, a significant temperature rise was achieved for the periodic pulse-shaped MF in comparison to the other two cases studies.

In order to understand the influences that a denser cluster has on the temperature gradient, other studies should be done investigating the interparticle interactions affecting the temperature increment and its derivative.

Appendices

A. Methods

In this appendix, we are deducing the equations for the temperature profiles introduced in Section 2, step by step. For simplicity, new variables are used to solve (1), where: ξƒŽπ‘…=1π›Όπ‘šπ‘Ÿ,π›Όπ‘š=π‘˜π‘šπœŒπ‘šπ‘π‘š,πœƒ(𝑅,𝑑)=π‘˜π‘šξ€·π‘‡π‘š(𝑅,𝑑)βˆ’π‘‡π‘ξ€Έπ‘…,𝑅0=ξƒŽ1π›Όπ‘šπ‘Ž.(A.1) Therefore, by substituting the new variables from (A.1) into the left part of (1) we receive the following: π‘˜π‘šβˆ‡2π‘‡π‘šπ‘˜(π‘Ÿ,𝑑)=π‘šπ‘Ÿπœ•2π‘Ÿπ‘‡π‘š(π‘Ÿ,𝑑)πœ•2π‘Ÿ||||π‘‡π‘šπ‘šπ‘…βˆš(π‘Ÿ,𝑑)β†’πœƒ(𝑅,𝑑)/π‘˜π‘Ÿβ†’π‘…π›Όπ‘š=1π›Όπ‘šπœ•2πœƒ(𝑅,𝑑)π‘…πœ•2𝑅.(A.2) And by substituting the new variables from (A.1) into the right part of (1) we receive that: πœŒπ‘šπ‘π‘šπœ•π‘‡π‘š(π‘Ÿ,𝑑)=πœŒπœ•π‘‘π‘šπ‘π‘šπ‘˜π‘šπœ•πœƒ(𝑅,𝑑)|||π‘…πœ•π‘‘π›Όπ‘šβ‰‘π‘˜π‘š/(πœŒπ‘šπ‘π‘š)=1π›Όπ‘šπœ•πœƒ(𝑅,𝑑).π‘…πœ•π‘‘(A.3) So (1) can be rewritten as follows: πœ•2πœƒ(𝑅,𝑑)πœ•2𝑅=πœ•πœƒ(𝑅,𝑑).πœ•π‘‘(A.4) The same procedure can be done to the BC, substituting the new variables from (A.1) into the left part of (2) to receive the following: βˆ’π‘˜π‘šβˆ‡ξ‚΅π‘‡π‘š(π‘Ÿ,𝑑)π‘Ÿξ‚Ά||||π‘Ÿ=π‘Ž=βˆ’π‘˜π‘šπœ•π‘‡π‘š(π‘Ÿ,𝑑)||||πœ•π‘Ÿπ‘‡π‘šπ‘šπ‘…βˆš(π‘Ÿ,𝑑)β†’πœƒ(𝑅,𝑑)/π‘˜π‘Ÿβ†’π‘…π›Όπ‘š1=βˆ’βˆšπ›Όπ‘šπœ•ξ‚΅πœ•π‘…πœƒ(𝑅,𝑑)𝑅1=βˆ’βˆšπ›Όπ‘šβˆ‡π‘…ξ‚΅πœƒ(𝑅,𝑑)𝑅|||||𝑅=𝑅0,(A.5) And (2) can be rewritten as follows: ξ‚΅βˆ’βˆ‡πœƒ(𝑅,𝑑)𝑅||||𝑅=𝑅0=π‘žπ‘ ξ…žξ…žβˆš(𝑑)π›Όπ‘š.(A.6) By taking the FT of (A.4) (defined as in (4)) one receives the transformation in the frequency domain, so Μƒπœ•0=π‘–πœ”πœƒ(𝑅,πœ”)βˆ’2Μƒπœƒ(𝑅,πœ”)πœ•2𝑅.(A.7) The general solution of (A.7) can be found as follows: Μƒπœƒ(𝑅,πœ”)=𝑐2(πœ”)π‘’βˆ’βˆšπ‘–πœ”π‘….(A.8) Substituting (A.8) into the LT of (A.6), the BC can be written as follows: 𝑐2π‘Žβˆšβˆšπ‘–πœ”+π›Όπ‘šπ‘Ž2=Μƒπ‘žπ‘ ξ…žξ…ž(πœ”)π‘’βˆšπ‘–πœ”π‘…0.(A.9) So Μƒπ‘Žπœƒ(𝑅,πœ”)=2βˆšπ›Όπ‘šβˆšξ‚€ξ‚€π‘Ž/π›Όπ‘šξ‚βˆšξ‚π‘–πœ”+1Μƒπ‘žπ‘ ξ…žξ…ž(πœ”)π‘’βˆ’βˆšπ‘–πœ”(π‘…βˆ’π‘…0)=π‘Žπ‘…0𝑅0βˆšπ‘–πœ”+1Μƒπ‘žπ‘ ξ…žξ…ž(πœ”)π‘’βˆ’βˆšπ‘–πœ”(π‘…βˆ’π‘…0).(A.10) Using technical computing software (Maple or/and Wolfram Mathematica) the inverse FT of Μƒπœƒ(𝑅,πœ”)/Μƒπ‘žπ‘ ξ…žξ…žξ‚(πœ”)=πœ™(𝑅,πœ”) for 𝑑>0, can be found by substituting π‘–πœ”β†’π‘  in (A.10) and taking the inverse Laplace transform of the received equation so that βŽ›βŽœβŽœβŽπ‘’πœ™(𝑅,𝑑)=π‘Žβˆ’(π‘…βˆ’π‘…0)2/4π‘‘βˆšβˆ’ξ‚€ξ€·πœ‹π‘‘erfcπ‘…βˆ’π‘…0ξ€Έβˆš/2βˆšπ‘‘+𝑑/𝑅0𝑒(π‘…βˆ’π‘…0)/𝑅0+𝑑/𝑅02𝑅0⎞⎟⎟⎠.(A.11) This function converges to 0 for π‘‘β†’βˆž or/and for π‘…βˆ’π‘…0β†’βˆž.

So the changes in the temperature can be found using (A.1) and (A.10) as follows. Ξ”π‘‡π‘š(π‘Ÿ,𝑑)=πœƒ(𝑅,𝑑)π‘˜π‘šπ‘Ÿβˆšπ›Όπ‘š=πœ™(𝑅,𝑑)βˆ—π‘žπ‘ ξ…žξ…ž(𝑑)π‘˜π‘šπ‘Ÿβˆšπ›Όπ‘š.(A.12) Equation (A.12) slightly differs than the one received by Keblinski et al. [38] due to the BC that define the heat flux coming from the surface of the MNP defining the heat created by the magnetic losses inside it, whereas Keblinski et al. [38] suggested that the heat sources are inside the medium of interest and that the heat-power density is constant in time. In order to analytically calculate (3) or (A.12), the general expression of π‘žπ‘ ξ…žξ…ž(𝑑) (Wmβˆ’2) must be found for each case.

Case 2 (a cosine MFP). For Case 1, the magnetization, 𝑀(𝑑), can be found in the time domain after substituting (5) and the MF in (4) and taking the inverse FT of it, that results in 𝑀(𝑑)=πœ’(𝑑)βˆ—=𝐻(𝑑)π΄πœ’0πœξƒ©ξ€·πœ”cos0π‘‘ξ€Έπœ+πœ”0ξ€·πœ”sin0𝑑ξƒͺ1(1/𝜏)2+πœ”02.(A.13) Substituting (A.13) into the magnetic induction [5] results in 𝐡(𝑑)=πœ‡0𝐻(𝑑)+πœ‡0𝑀(𝑑).(A.14) Further substituting the received magnetic induction described in (A.14) into (7), one can calculate the conversion of the magnetic energy into heat losses, resulting in 𝐇(𝐭)β‹…πœ•π(𝐭)=ξ€πœ•π‘‘π‘‘πœ”π‘‘πœ”ξ…žπ»ξ€·π‘Ÿ,πœ”ξ…žξ€Έξ€·π‘–πœ”πœ‡0πœ‡ξ€Έ(πœ”)⋅𝐻(π‘Ÿ,πœ”)𝑒𝑖(πœ”βˆ’πœ”β€²)𝑑=πœ‡0𝐴2ξ€·πœ”cos0π‘‘ξ€Έβ‹…ξ€œξ‚€ξ‚€πœ’π‘‘πœ”π‘–πœ”1+0×𝛿1+π‘–πœ”πœξ‚ξ‚πœ”βˆ’πœ”0ξ€Έξ€·+π›Ώπœ”+πœ”0ξ€Έ2π‘’π‘–πœ”π‘‘=πœ‡0𝐴2ξ€·πœ”cos0π‘‘ξ€Έπœ”0Γ—ξ‚΅πœ”0πœ’0𝜏1+πœ”02𝜏2ξ€·πœ”cos0π‘‘ξ€Έξ€·πœ”βˆ’sin0π‘‘ξ€Έξ‚Έπœ’1+01+πœ”02𝜏2ξ‚Ήξ‚Ά=𝑃Loss(𝑑)+πœ•π‘ˆ.πœ•π‘‘(A.15) Or 𝑃Loss(𝑑)=πœ‡0𝐴2ξ€·πœ”cos0π‘‘ξ€Έπœ”0πœ”0πœ’0𝜏1+πœ”02𝜏2ξ€·πœ”cos0𝑑.(A.16) Because 𝑃Loss(𝑑) is only a function of time between 0<π‘Ÿ<π‘Ž (isotropic and homogeneous material), then the outward heat flux at π‘Ÿ=π‘Ž can be calculated as follows: π‘žπ‘ ξ…žξ…ž(π‘Ÿ=π‘Ž,𝑑)4πœ‹π‘Ž2=4πœ‹π‘Ž33𝑃Loss(𝑑).(A.17) Or π‘žπ‘ ξ…žξ…ž(π‘Ÿ=π‘Ž,𝑑)=π‘Žπœ‡0𝐴2πœ”06ξƒ©πœ”0πœ’0πœξ€·πœ”1+0πœξ€Έ2ξ€·ξ€·cos2πœ”0𝑑ξƒͺ.+1(A.18) Taking the FT of (A.18) and substituting it in (A.10) one can calculate the FT of πœƒ(𝑅,𝑑) to receive the following: Μƒπœƒ(𝑅,πœ”)=π‘Žπœ‡0𝐴2πœ”06Γ—ξƒ©πœ”0πœ’0πœξ€·πœ”1+0πœξ€Έ2ξƒ©π›Ώξ€·πœ”βˆ’2πœ”0ξ€Έξ€·+π›Ώπœ”+2πœ”0ξ€Έ2β‹…+𝛿(πœ”)ξƒͺξƒͺπ‘Žπ‘…0𝑅0βˆšξ‚€βˆ’βˆšπ‘–πœ”+1exp.π‘–πœ”π‘…(A.19)

Case 3 (a PMP profile). The PMP, (11), can be decomposed using the theory of Fourier’s series into its harmonics to receive [71, 72] the following: Δ𝐻(𝑑)=2𝐴⋅𝑇𝑠+βˆžξ“π‘›=14π΄ξ‚΅πœ‹π‘›sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άξ€·cosπ‘›πœ”0𝑑.(A.20) By substituting (A.20) into (A.14) and then using them in (7) we can calculate the total heat dissipation for this case as follows: 𝐇(𝐭)β‹…πœ•π(𝐭)=ξ€πœ•π‘‘π‘‘πœ”π‘‘πœ”ξ…žπ»ξ€·π‘Ÿ,πœ”ξ…žξ€Έξ€·π‘–πœ”πœ‡0πœ‡ξ€Έ(πœ”)⋅𝐻(π‘Ÿ,πœ”)𝑒𝑖(πœ”βˆ’πœ”β€²)𝑑=πœ‡0ξ‚΅Ξ”2𝐴⋅𝑇𝑠+4π΄ξ‚΅πœ‹π‘šsinπ‘šπœ‹Ξ”π‘‡π‘ ξ‚Άξ€·cosπ‘šπœ”0𝑑⋅4π΄ξ‚΅πœ‹π‘›sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άξ€œξ‚€ξ‚€πœ’π‘‘πœ”π‘–πœ”1+0×𝛿1+π‘–πœ”πœξ‚ξ‚πœ”βˆ’π‘›πœ”0ξ€Έξ€·+π›Ώπœ”+π‘›πœ”0ξ€Έ2π‘’π‘–πœ”π‘‘=πœ‡0ξ‚΅Ξ”2𝐴⋅𝑇𝑠+4π΄ξ‚΅πœ‹π‘šsinπ‘šπœ‹Ξ”π‘‡π‘ ξ‚Άξ€·cosπ‘šπœ”0𝑑⋅4π΄πœ‹ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπœ”0Γ—ξƒ©π‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2ξ€·cosπ‘›πœ”0π‘‘ξ€Έξ€·βˆ’sinπ‘›πœ”0π‘‘ξ€Έξƒ¬πœ’1+0ξ€·1+π‘›πœ”0ξ€Έ2𝜏2ξƒ­ξƒͺ=𝑃Loss(𝑑)+πœ•π‘ˆ.πœ•π‘‘(A.21) Therefore, we can find that the heat losses equal to 𝑃Loss(𝑑)=πœ‡0ξ‚΅Ξ”2A⋅𝑇𝑠+4π΄ξ‚΅πœ‹π‘šsinπ‘šπœ‹Ξ”π‘‡π‘ ξ‚Άξ€·cosπ‘šπœ”0𝑑⋅4π΄πœ‹ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπœ”0ξƒ©π‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2ξ€·cosπ‘›πœ”0𝑑ξƒͺ.(A.22) And for Ξ”<𝑑<𝑇𝑠: 𝐇(𝐭)β‹…πœ•π(𝐭)πœ•π‘‘=0.(A.23) Using (A.17) and (A.22) we can calculate the heat flux at the surface of the MNP. π‘žπ‘ ξ…žξ…ž(π‘Ÿ=π‘Ž,𝑑)=𝑃Lossπ‘Ž(𝑑)3=π‘Ž3πœ‡0ξ‚΅Ξ”2𝐴⋅𝑇𝑠+4π΄ξ‚΅πœ‹π‘šsinπ‘šπœ‹Ξ”π‘‡π‘ ξ‚Άξ€·cosπ‘šπœ”0𝑑⋅4π΄πœ‹ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπœ”0ξƒ©π‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2ξ€·cosπ‘›πœ”0𝑑ξƒͺ.(A.24) By taking the FT of the resulted heat flux and substituting it in (A.10) one can receive Μƒ=π‘Žπœƒ(𝑅,πœ”)38πœ”0πœ‡0𝐴2β‹…Ξ”πœ‹π‘‡π‘ ξ“ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άπ‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2Γ—ξƒ©π›Ώξ€·πœ”βˆ’π‘›πœ”0ξ€Έξ€·+π›Ώπœ”+π‘›πœ”0ξ€Έ2ξƒͺπ‘Žπ‘…0𝑅0βˆšξ‚€βˆ’βˆšπ‘–πœ”+1exp+π‘Žπ‘–πœ”π‘…38πœ‡0𝐴2πœ”0πœ‹2sinπ‘šπœ‹Ξ”/π‘‡π‘ ξ€Έπ‘šξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚ΆΓ—π‘›πœ”0πœ’0πœξ€·1+π‘›πœ”0ξ€Έ2𝜏2β‹…ξƒ©π›Ώξ€·πœ”βˆ’(𝑛+π‘š)πœ”0ξ€Έξ€·+π›Ώπœ”+(𝑛+π‘š)πœ”0ξ€Έ2+π›Ώξ€·πœ”βˆ’(π‘šβˆ’π‘›)πœ”0ξ€Έξ€·+π›Ώπœ”+(π‘šβˆ’π‘›)πœ”0ξ€Έ2ξƒͺΓ—π‘Žπ‘…0𝑅0βˆšξ‚€βˆ’βˆšπ‘–πœ”+1exp.π‘–πœ”π‘…(A.25)

When looking at (A.25), the multiplication of each coefficient’s amplitude with its matched harmonic must meet the biologically noninvasive limitation π΄π‘šβ‹…π‘“π‘šβ‰€5β‹…109A/mΒ·s. The mathematical justification is deduced next.

Looking at the eddy currents that evolve in the body [73] 𝐸(πœ”)=βˆ’π‘–πœ”π‘Ÿ2𝐡𝑧,𝐽(πœ”)=βˆ’π‘–πœ”π‘ŸπœŽ2π΅π‘§π‘Ÿβ†’πΈ(𝑑)=βˆ’2πœ•π΅π‘§(𝑑)πœ•π‘‘=βˆ’πœ‡0π‘Ÿ2πœ•π»(𝑑),πœ•π‘‘π½(𝑑)=βˆ’π‘ŸπœŽ2πœ•π΅π‘§(𝑑)πœ‡πœ•π‘‘=βˆ’0π‘ŸπœŽ2πœ•π»(𝑑).πœ•π‘‘(A.26) They can be written using (A.20) as follows: 𝐸(𝑑)=βˆ’πœ‡0π‘Ÿ2πœ•π»(𝑑)πœ•π‘‘=πœ‡0π‘Ÿ2πœ”04π΄πœ‹ξ‚΅sinπ‘šπœ‹Ξ”π‘‡π‘ ξ‚Άξ€·sinπ‘šπœ”0𝑑,πœ‡π½(𝑑)=βˆ’0π‘ŸπœŽ2πœ•π»(𝑑)=πœ‡πœ•π‘‘0π‘ŸπœŽ2πœ”04π΄πœ‹ξ‚΅sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άξ€·sinπ‘›πœ”0𝑑.(A.27) So the eddy losses inside the body can be found using Poynting theory [5] as follows: 𝑃average=1π‘‡π‘ ξ€œ1𝑃(𝑑)𝑑𝑑=π‘‡π‘ ξ€œ=𝜎𝐸(𝑑)𝐽(𝑑)𝑑𝑑𝑇𝑠2π‘Ÿπœ”0πœ‡0π΄πœ‹ξ‚Ά2ξ€œξ“ξ‚΅sinπ‘šπœ‹Ξ”π‘‡π‘ ξ‚Άξ€·sinπ‘šπœ”0𝑑×sinπ‘›πœ‹Ξ”π‘‡π‘ ξ‚Άξ€·sinπ‘›πœ”0𝑑=πœŽπ‘‘π‘‘π‘‡π‘ ξ€·πœ‹π‘Ÿπ‘“0πœ‡0ξ€Έ2ξ‚€4π΄πœ‹ξ‚2sin2ξ‚΅π‘šπœ‹Ξ”π‘‡π‘ ξ‚ΆΓ—ξ€œsin2ξ€·π‘šπœ”0𝑑𝑑𝑑=πœŽπœ‹π‘Ÿπ‘“0πœ‡0ξ€Έ2ξ‚€4π΄πœ‹ξ‚212sin2ξ‚΅π‘šπœ‹Ξ”π‘‡π‘ ξ‚Ά=𝐴eff2πœŽξ€·πœ‹π‘Ÿπ‘“0πœ‡0ξ€Έ2.(A.28) The last expression is the same as the one received by Atkinson et al. [16]. Therefore the limitation on the MFS and the frequency can be summarized as follows [11, 16, 24]: 𝐴eff𝑓0=𝑓04π΄πœ‹βˆš2ξƒͺξƒŽξ“sin2ξ‚΅π‘šπœ‹Ξ”π‘‡π‘ ξ‚Άβ‰€5β‹…109A/mβ‹…s.(A.29) For 𝑁=25 and a duty cycle of 𝑑=Ξ”/𝑇𝑠=0.2, the treatment is safe as long as 𝐴eff𝑓0=𝑓04π΄πœ‹βˆš2ξƒͺξƒŽξ“sin2ξ‚΅π‘šπœ‹Ξ”π‘‡π‘ ξ‚Ά=𝐴𝑓03.1≀5β‹…109A/mβ‹…s.(A.30) Consequently, as long as (A.30) is valid the treatment is safe. Choosing other maximal summation index values such as 𝑁=20 will result in a new constraint over the frequency and the MFS that must fulfill 𝐴eff𝑓0=𝑓0√(4𝐴/πœ‹βˆš2)βˆ‘sin2(π‘šπœ‹Ξ”/𝑇𝑠)β‰ˆπ΄π‘“02.8≀5β‹…109 A/mΒ·s and so on.

Moreover, for frequencies lower than 10 MHz there is essentially no attenuation of the MFS within cylinders of muscle-equivalent material, therefore the maximal harmonic frequency should not exceed 10 MHz [16].

Case 4 (a discontinuous cosine MFP). As for the previous case, we decompose the MF using the theory of Fourier’s series into its harmonics for 0≀𝑑≀Δ to receive [71, 72] =𝐻(𝑑)βˆžξ“π‘›=12𝐴𝑇𝑠(ξ€·πœ”sinΞ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·(ξ€·πœ”sinΞ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­ξ€·Γ—cosπ‘›πœ”1𝑑.(A.31)

For Δ≀𝑑≀𝑇𝑠, the magnetic power losses are zero, because the MF dose not exists.

By substituting (A.31) into (A.14) and then using them in (7) we can calculate the total heat dissipation for this case as follows:𝐇(𝐭)β‹…πœ•π(𝐭)=ξ€πœ•π‘‘π‘‘πœ”π‘‘πœ”ξ…žπ»ξ€·π‘Ÿ,πœ”ξ…žξ€Έξ€·π‘–πœ”πœ‡0πœ‡ξ€Έ(πœ”)⋅𝐻(π‘Ÿ,πœ”)𝑒𝑖(πœ”βˆ’πœ”β€²)𝑑=πœ‡0ξ‚΅2𝐴𝑇𝑠2⋅(ξ€·πœ”sinΞ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒ­ξ€·Γ—cosπ‘šπœ”0𝑑ξƒͺβ‹…ξ“ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­Γ—ξ€œξ‚€ξ‚€πœ’π‘‘πœ”π‘–πœ”1+0×𝛿1+π‘–πœ”πœξ‚ξ‚πœ”βˆ’π‘›πœ”1ξ€Έξ€·+π›Ώπœ”+π‘›πœ”1ξ€Έ2π‘’π‘–πœ”π‘‘=πœ‡02𝐴𝑇𝑠2Γ—ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒ­ξ€·Γ—cosπ‘šπœ”1𝑑ξƒͺβ‹…ξ“ξƒ©ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’nπœ”1ξƒͺΓ—ξƒ©ξ€·π‘›πœ”1ξ€Έ2πœ’0πœξ€·1+π‘›πœ”1ξ€Έ2𝜏2ξ€·cosπ‘›πœ”1π‘‘ξ€Έξ€·βˆ’sinπ‘›πœ”1π‘‘ξ€Έξƒ¬πœ’1+0ξ€·1+π‘›πœ”1ξ€Έ2𝜏2ξƒ­ξƒͺ=𝑃Loss(𝑑)+πœ•π‘ˆ.πœ•π‘‘(A.32) Therefore, we can find that the heat losses equal to 𝑃Loss(𝑑)=πœ‡02𝐴𝑇𝑠2Γ—ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒ­ξ€·Γ—cosπ‘šπœ”1𝑑ξƒͺβ‹…ξ“ξƒ©ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒͺΓ—ξƒ©ξ€·π‘›πœ”1ξ€Έ2πœ’0πœξ€·1+π‘›πœ”1ξ€Έ2𝜏2ξ€·cosπ‘›πœ”1𝑑ξƒͺ.(A.33) Using (A.17) and (A.33) we can calculate the heat flux at the surface of the MNP. π‘žπ‘ ξ…žξ…ž(π‘Ÿ=π‘Ž,𝑑)=𝑃Lossπ‘Ž(𝑑)3=π‘Ž3πœ‡02𝐴𝑇𝑠2ξ“ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒ­ξ€·cosπ‘šπœ”1𝑑ξƒͺβ‹…ξ“ξƒ©ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒͺΓ—ξƒ©ξ€·π‘›πœ”1ξ€Έ2πœ’0πœξ€·1+π‘›πœ”1ξ€Έ2𝜏2ξ€·cosπ‘›πœ”1𝑑ξƒͺ.(A.34) By taking the FT of the resulted heat flux and substituting it in (A.10) one can receive Μƒ=πœƒ(𝑅,πœ”)2π‘Ž3πœ‡0𝐴𝑇𝑠2π‘Žπ‘…0𝑅0βˆšξ‚€βˆ’βˆšπ‘–πœ”+1expξ‚Γ—ξ€·ξ€·πœ”π‘–πœ”π‘…ξ“ξ“ξƒ¬ξƒ©sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒͺΓ—ξƒ©ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒͺβ‹…ξ€·π‘›πœ”1ξ€Έ2πœ’0πœξ€·1+π‘›πœ”1ξ€Έ2𝜏2Γ—ξƒ©π›Ώξ€·πœ”βˆ’(𝑛+π‘š)πœ”1ξ€Έξ€·+π›Ώπœ”+(𝑛+π‘š)πœ”1ξ€Έ2+π›Ώξ€·πœ”βˆ’(π‘šβˆ’π‘›)πœ”1ξ€Έξ€·+π›Ώπœ”+(π‘šβˆ’π‘›)πœ”1ξ€Έ2.ξƒͺξƒ­(A.35) When looking at (A.35) the multiplication of each coefficient’s amplitude with its matched harmonic must meet the biologically noninvasive limitation π΄π‘šβ‹…π‘“π‘šβ‰€5β‹…109A/mΒ·s. The mathematical justification is deduced next.

For 𝑓0>𝑓1 and 𝑇𝑠=1/𝑓1, we find that (A.26) becomes 𝐸(𝑑)=βˆ’πœ‡0π‘Ÿ2πœ•π»(𝑑)πœ•π‘‘=πœ‡0π‘Ÿ2πœ”1𝑛2π΄π‘‡π‘ ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·(ξ€·πœ”sinΞ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­ξ€·sinπ‘›πœ”1𝑑,πœ‡π½(𝑑)=βˆ’0π‘ŸπœŽ2πœ•π»(𝑑)=πœ‡πœ•π‘‘0π‘ŸπœŽ2πœ”1Γ—ξ“π‘š2π΄π‘‡π‘ ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘šπœ”1ξ€Έξ€Έπœ”0+π‘šπœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘šπœ”1ξ€Έξ€Έπœ”0βˆ’π‘šπœ”1ξƒ­ξ€·sinπ‘šπœ”1𝑑.(A.36) Consequently, (A.28) becomes 𝑃avr=1π‘‡π‘ ξ€œ1𝑃(𝑑)𝑑𝑑=π‘‡π‘ ξ€œ=𝜎𝐸(𝑑)𝐽(𝑑)𝑑𝑑𝑇𝑠2π‘Ÿπœ”1πœ‡0𝐴𝑇𝑠2β‹…ξ€œξ“π‘›ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­ξ€·sinπ‘›πœ”1π‘‘ξ€ΈΓ—ξ“π‘šξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­ξ€·sinπ‘›πœ”1𝑑=πœŽπ‘‘π‘‘π‘‡π‘ ξ‚΅2π‘Ÿπœ”1πœ‡0𝐴𝑇𝑠2×𝑛2ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­2Γ—ξ€œsin2ξ€·π‘›πœ”1𝑑𝑑𝑑=πœŽπœ‹π‘Ÿπ‘“1πœ‡0ξ€Έ2ξ‚΅4𝐴𝑇𝑠2Γ—12𝑛2ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­2=𝐴eff2πœŽξ€·πœ‹π‘Ÿπ‘“1πœ‡0ξ€Έ2.(A.37) The last expression is the same as the one received by Atkinson et al. [16]. Therefore, the limitation on the MFS and the frequency can be summarized as follows [11, 16, 24]: 𝐴eff𝑓1=𝑓14π΄π‘‡π‘ βˆš2ξƒͺΓ—ξ„Άξ„΅ξ„΅βŽ·ξ“π‘›2ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­2≀5β‹…109A/mβ‹…s.(A.38) For 𝑁=25, a duty cycle of 𝑑=Ξ”/𝑇𝑠=0.2 and πœ”0=6πœ”1, the treatment is safe as long as 𝐴eff𝑓1=𝑓14π΄π‘‡π‘ βˆš2ξƒͺΓ—ξ„Άξ„΅ξ„΅βŽ·ξ“π‘›2ξƒ¬ξ€·ξ€·πœ”sin(Ξ”/2)0+π‘›πœ”1ξ€Έξ€Έπœ”0+π‘›πœ”1+ξ€·ξ€·πœ”sin(Ξ”/2)0βˆ’π‘›πœ”1ξ€Έξ€Έπœ”0βˆ’π‘›πœ”1ξƒ­2=𝐴𝑓14.6<𝐴𝑓0≀5β‹…109A/mβ‹…s.(A.39)

Consequently, as long as (A.39), the treatment will be safe. Moreover, for frequencies lower than 10 MHz, there is essentially no attenuation of the MFS within cylinders of muscle-equivalent material, therefore the maximal harmonic frequency should not exceed 10 MHz [16].

B. COMSOL: Results

In order to validate the analytic solutions and the MATLAB simulations, a numerical simulation was performed using COMSOL, for the same thermal and magnetic properties given in Tables 1 and 2. The simulation results can be seen for each case studied in Methods and Results parts in this Appendix.

For Case 1, the mathematical expression of the temperature increment, (9), was plotted as a function of time and space, where the results are given in Figure 10 for 𝑇𝑠=2.5 μs.

The maximal temperature elevation in Figure 10 reached a value of 2.25 nΒ°K on the surface of the MNP which is 0.15 nΒ°K higher than the one received for the analytic simulation, Figure 3.

At 2 nm apart from the surface of the MNP surface, the temperature elevation reached a value of 2.05 nΒ°K that is 0.2 nΒ°K higher than the one receive in Figure 3. Again, there is a small difference between both simulations results. As the observation point gets further from the surface of the MNP, the temperature differences get bigger; reaching a value of 0.4 nΒ°K at an observation point located 10 nm apart from the surface.

This may be caused by the triangles constructing the COMSOL’s numeric mesh, which are used to solve numerically the heat problem, that get larger and bigger, as the observation point gets further from the MNP surface contributing to the error.

Comparing between Figures 3 and 10 we conclude that the numerical simulation fits the analytic solution.

For Case 2, as in Case 1, in order to validate the analytic solution a numerical simulation was also performed using COMSOL for the same thermal and magnetic properties (Tables 1 and 2). The simulation result can be seen in Figure 11.

The maximal temperature elevation in Figure 11 reaches a value of 8.5 nΒ°K on the surface of the MNP, which is 0.4 nΒ°K higher than the one receive in Figure 5.

At 2 nm apart from the surface of the MNP surface, the temperature elevation reached a value of 7.5 nΒ°K, that is 0.2 nΒ°K higher than the one receive in Figure 5. Again, it seems that there exists a small difference between the simulations results. As the observation point gets further from the surface of the MNP, the differences gets bigger reaching a value of 0.8 nΒ°K at an observation point located 10 nm a part from the surface.

This may be caused by the bigger triangles in the mesh that are formed in the COMSOL software as the observation point gets further from the MNP surface contributing to the error.

As can be seen from Figure 11, there are 5 peaks during the time that the MF is tuned on that fit the number of peaks in Figure 5, these peaks evolve due to the final number of harmonics that form the PMP MF, as given by (11). However, there is a slightly difference in the temperature profiles between Figures 11 and 5; in Figure 11 the first peak is lower than the others in comparison to Figure 5 where the first peak is about the same high as the last peak.

Again, there are some small changes between both software simulations as expected; however the results for both simulations conclude that there is a benefit in using the PMPs instead of the cosine MFP due to the higher temperature rise values received for the same magnetic parameters.

For Case 3, we used again the numerical simulation COMSOL in order to validate the analytic solution for the same thermal and magnetic properties. The simulation result can be seen in Figure 12.

The maximal temperature elevation in Figure 12 reached a value of 2.3 nΒ°K on the surface of the MNP which is the same as the one receive in Figure 7.

At 2 nm apart from the surface of the MNP surface the temperature elevation reached a value of 2 nΒ°K, that is, 0.2 nΒ°K higher than the one received in Figure 7. Again, there is a small difference between the simulations results. As the observation point get further from the surface of the MNP, the differences gets bigger reaching a value of 0.3 nΒ°K at an observation point located 10 nm a part from the surface.

As explain before, this may be caused by the bigger triangles in the mesh that are formed in the COMSOL software as the observation point gets further from the MNP surface contributing to the error. Although, there are some small changes between both simulations as expected, the maximal temperature rise is almost the same as the cosine MFP.

C. The Effects the Maximal Number of Indexes Has on Cases 2 and 3: Results

In Appendix C, we examined the influences that the maximal numbers of indexes composing the MF signal have on the temperature rise and on the temperature rate rise for Case 2 and Case 3.

The maximal index numbers for summation were chosen as 𝑁=100, 𝑁=25, 𝑁=15, 𝑁=10, and 𝑁=1. Above 𝑁=25 the MF is practically absorbed in the tissue [14], but this fact was not taken in consideration in the simulations results.

Case 3. The temperature rise for Case 2, as a function of the maximal summation indexes can be seen in Figure 13.

From Figure 13, we concluded that the maximal temperature rise depends on the number of harmonics composing the MF signal. For 𝑁=100, the maximal temperature rise reaches a value of 50 nΒ°K for a core radius of 7.7 nm, that is, 1.6 times higher than the maximal value received for 𝑁=25. The summation of 100 indexes can be seen as the ideal PMPs shaped MFP.

For 𝑁=15, we receive a temperature rate of 28 nΒ°K for a core radius of 8.8 nm and for 𝑁=10 we received a value of 26 nΒ°K for a core radius of 9 nm. Furthermore, we can see that the number of indexes composing the MF changes the optimal radius as it gets smaller as the index number gets bigger.

Now, we examined the influences that the number of maximal summation indexes composing the MF signal has on the temperature rate rise. The chosen numbers were 𝑁=100, 𝑁=25, 𝑁=15, 𝑁=10, and 𝑁=1.

From Figure 14, we concluded that the maximal temperature rate rise depends on the number of harmonics composing the MF signal. For 𝑁=100, the maximal temperature rate reaches a value of 3.9Β°Ksβˆ’1 is received for a core radius 7.6 nm, and is 3.9 times higher than the value received for 𝑁=25. For 𝑁=15 we receive a temperature rate of 0.5Β°Ksβˆ’1 for a core radius of 8.2 nm, that is, half the value received for 𝑁=25, and for 𝑁=10, we received a value of 0.3Β°Ksβˆ’1 for a core radius of 8.5 nm.

Case 4. Now, we examined the influences that the maximal number of summation indexes composing the MF signal has on the temperature rise. The chosen numbers were 𝑁=100, 𝑁=25, 𝑁=15, 𝑁=10, and 𝑁=1.
From Figure 15, we concluded that the maximal temperature rise depends on the number of harmonics composing the MF signal. For 𝑁=100, the maximal temperature rise reaches a value of 5 nΒ°K for a core radius of 9.3 nm, that is, 1.1 times higher than the maximal value received for 𝑁=25. For 𝑁=15, we receive a temperature rate of 5.5 nΒ°K for a core radius of 9.2 nm, and for 𝑁=10, we received a value of 6.1Β°K for a core radius of 9.3 nm. As already mentioned, there is a limitation to the highest frequency that can be used for MH and should not exceed 10 MHz [16], in our case this limits the summation to 25 indexes that compose the MF signal. Moreover, we can see that the number of indexes composing the MF changes the optimal radius; it gets smaller as the index number gets higher. Furthermore, we can see that the number of indexes composing the MF changes the optimal radius by getting smaller as the index number gets bigger.

Now, we examined the influences that the number of indexes composing the MF signal has on the temperature rate rise. The chosen numbers were 𝑁=100, 𝑁=25, 𝑁=15, 𝑁=10, and 𝑁=1.

From Figure 16, we concluded that the maximal temperature rise rate depends on the number of harmonics composing the MF signal. For 𝑁=100, the maximal temperature rate reaches a value of 0.09Β°Ksβˆ’1, is received for a core radius 8.4 nm, and is 4.5 times higher than the value received for 𝑁=25. For 𝑁=15, we receive a temperature rate of 0.016Β°Ksβˆ’1 for a core radius of 9.3 nm, that is, half the value received for 𝑁=25, and for 𝑁=10 we received a value of 0.015Β°Ksβˆ’1 for a core radius of 9.3 nm.

As already mentioned there is a limitation to the highest frequency that can be used for MH and should not exceed 10 MHz [16], in our case this limits the summation to 25 indexes that compose the MF signal. Moreover, we can see that the number of indexes composing the MF changes the optimal radius and it gets smaller as the index number gets higher.