Identification of Curie temperature distributions in magnetic particulate systems

This paper develops a methodology for extracting the Curie temperature distribution from magnetisation versus temperature measurements which are realizable by standard laboratory magnetometry. The method is integral in nature, robust against various sources of measurement noise, and can be adopted to a wide range of granular magnetic materials and magnetic particle systems. The validity and practicality of the method is demonstrated using large-scale Monte–Carlo simulations of an Ising-like model as a proof of concept, and general conclusions are drawn about its applicability to different classes of systems and experimental conditions.


Introduction
Development of magnetic nanotechnologies, such as nanostructured high temperature permanent magnets [1], heat assisted magnetic recording technologies (HAMR) [2], or biomedical applications [3,4], relies upon the availability of methodologies for the accurate large-scale characterisation of magnetic nanoparticles and granular materials at elevated temperatures. The presence of non-uniformities leads to a broadening of the phase transition region and the consequent difficulties in determining a single Curie temperature T c . Instead, identification of a distribution of T c is often required.
For example, in HAMR, the extent of the broadening of the T c distribution determines the noise performance and the quality of the recording process, and its accurate quantification is essential for the optimisation and quality control [2]. In heat assisted cancer therapy, developing self-regulated magnetic hyperthermia requires optimisation of T c distributions in assemblies of low-T c nanoparticles [4]. Identifying T c distributions in such systems, dominated by disorder and complex spatial inhomogeneities, is in general a challenge which hinders the optimisation of their performance at high temperatures.
Direct techniques, such as those based upon specialised laser systems, have been developed for the experimental determination of the switching temperature distribution in HAMR media films [5][6][7]. Such techniques become difficult to implement in the case of irregular and highly disordered magnetic particle distributions. Indirect techniques have instead heat-assisted magnetic recording (Some figures may appear in colour only in the online journal)

Letter
Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. proven practical, based on identifying the T c distribution from temper ature dependent measurement of the magnetic moment m(T) or the AC susceptibility using inverse problem solving techniques [8][9][10]. Such techniques have been applied predominantly to thin film samples thus far and did not include the finite size effects, nor did they relate the T c distribution to the intrinsic properties of elementary particles in a systematic manner.
In this work, an indirect approach is developed to extract the T c distribution from m(T) data in assemblies of finite size magnetic grains or particles. The method systematically incorporates knowledge of the finite grain size distribution, critical exponents, and bulk Curie temperature T b c , and as a consequence of universality in the phase transition region, the method can be adopted to a broad class of different material systems. This universality also allows validation of the methodology using simplified models without sacrificing its generality, such as the two-dimensional Ising model used in this work for which the T b c and critical exponents are known from analytical calculations [11,12].
The key notion adopted in this work is that of the finite system size critical temperature T c = T c (D), to be distinguished from the bulk Curie temperature T b c relevant in the thermodynamic limit [6,[11][12][13][14][15][16][17]. Strictly speaking, the phase trans ition temperature is defined only in the thermodynamic limit as a critical point which marks non-analytic and divergent behaviour of the thermodynamic state functions [11]. In finite size systems the divergent point becomes a rounded peak and the notion of the finite size Curie temperature T c becomes indefinite. It is standard to adopt the definition that T c is related or indeed equal to the peak of a rounded state function, such as the derivative of magnetisation M(T) of a grain, dM/dT, or the temperature dependent magnetic susceptibility [14,17]. In practice, both definitions give similar results and any differences diminish as D increases towards the bulk, as illustrated in figures 1(a)-(c). The trend obeys the well known scaling law [11][12][13][14][15][16][17]: where ν is a material dependent critical exponent associated with the the correlation length in the phase transition region.
The study outlined in this paper will consider assemblies of independent particles of variable size, each viewed as an elementary thermodynamic system characterised by a unique value of size-dependent T c which obeys relation (1), and develop an approach for extracting the finite size T c distribution from a typical m(T) measurement of such a particle assembly. The application domain includes the analysis of dilute systems of magnetic nanoparticles, or granular magnetic materials for heat assisted magnetic recording where the intergranular exchange interactions are optimised to produce weak correlations.
It is worthwhile mentioning that although the peak of dM/dT can be used to define T c of a single grain, it is not true that the derivative of m(T) of a non-interacting granular assembly will similarly define the distribution of T c of grains.
The dm/dT , being a superposition of dM/dT of all grains within an assembly, will also contain convoluted contributions from the finite widths of functional dependences of dM/dT of grains, expected to statistically vary from grain-to-grain. This will lead to broadening of dm/dT and add a fictitious contribution to the intrinsic distribution of T c of grains. The essence of identification methodologies such as developed in this article is to systematically deconvolve these contributions and extract the genuine T c distribution.

Overview and application of the framework
The key result derived in this article is the integral expression for temperature dependent magnetic moment of an assembly of independent magnetic grains: where m 0 is a constant, d is the dimension of the smallest characteristic region of a magnetic particle, β is the magnetisation universal critical exponent, and μ is the universal scaling function to be discussed in detail. Expressions analogous to equation (2) can also be derived for other thermodynamic variables, such as magnetic susceptibility or specific heat. (2) is represented in terms of the reduced Curie temperature t c : If the particle size distribution is chosen to be the experimentally relevant lognormal distribution, f t (t c ) will be shown to also take lognormal form: where t c and σ 2 tc are respectively the logarithmic mean and variance. The corresponding arithmetic mean and variance follow from the standard properties of lognormal distribution as t c = exp(t c +σ 2 tc /2) and σ 2 tc = t c 2 (exp(σ 2 tc ) − 1). According to (3), the mean non-reduced Curie temperature is The parameters t c and σ tc in (4) can be identified by least-square fitting equation (2) to experimental m(T) data. This in principle also allows extracting the values of the bulk Curie temperature T b c and the critical exponents β, ν as fit parameters, although any knowledge of these parameters from independent experiments or simulations aids in reducing the fit parameter correlation. For example, T b c can be estimated from finite lattice size simulations using the Binder cumulant expansion method [19,20], and the critical exponents can be estimated using the modified Arrott (Kouvel-Fisher) technique [21], or the finite size scaling analysis used below.
The relation (2) can be adapted to different classes of material systems in a straightforward manner by choosing appropriate critical exponents and scaling functions (table 1). This generality is a result of universality arising from scale invariance in the phase transition region to be discussed in more detail next.

Derivation of the fitting function
Consider a system of magnetic grains of variable size ( figure  1(a)). The volume of a grain is V = CD d , where for example for cylinder D is the base diameter, d = 2, and C = Lπ/4 with L being the height, assuming L > D. For a sphere D is the diameter, d = 3 and C = π/6. For convenience, the dimensionless size D = D/D 0 is introduced, where D 0 is some reference length such as the atomic lattice spacing. Each grain in the ensemble is assumed to be an elementary thermodynamic system with magnetic moment m g having length m g = | m g | dependent on T and D, which can be expressed using the definitions above as: where m 0 = CD d 0 M s , M s is the saturation magnetisation at T = 0, and the temperature and size dependent magnetisation M(T, D) is normalised to be dimensionless. The fact that M depends upon D is a manifestation of the finite size effect.
where near the phase transition the M(T, D) dependence is known to take universal form [11,16,17]: The scaling function μ(x) can be established by matching the magnetisation data M(T, D) for grains of different size D using scaled c and critical exponents β and ν are known (figures 1(b) and (d)). Alternatively, if unknown, β, ν, and T b c can be found by obtaining the collapse of the M(T, D) data for different D. 3 This collapse then gives the functional form of μ.
Inserting equation (7) into (6) and arranging using (1) and (3), leads to the expression for ensemble magnetic moment as given by equation (2). The explicit relationship between f t (t c ) and f D (D) can be obtained by substituting equation (1) into (6): which is the standard transformation between two probability distributions of random variables related by a functional relation [18]. For example, for the experimentally relevant lognormal grain size distribution: where D and σ 2 D are the logarithmic mean and variance, combining (1), (3) and (8) gives the f t (t c ) also in the form of lognormal distribution as given in (4) with t c = −D/ν and σ tc =σ D /ν.
It is worthwhile noting that the validity of the scaling relation (7) is a consequence of the universality and emergent scale invariance near phase transitions. The universal behaviour of different classes of materials can represented by the same sets of critical exponents and scaling functions [11,23]. This universality also means that complex materials can be studied by using simplified model systems with the same critical exponents, i.e. models in the same universality class, making simulations based on simplified model systems especially powerful tool for investigating phase transitions.

Proof of concept
A zero field two-dimensional Ising model is used as a proof of concept model system for validating the methodology. The zero field Ising Hamiltonian reads H = −J ij s i s j , where the double-sum extends over the nearest neighbour spin pairs, and J > 0. Metropolis Monte-Carlo simulation technique will be used to study the T dependent behaviour.
The advantage of choice of the Ising model is in the availability of analytical results. Moreover, given the low value of the critical exponent ratio β/ν, equation (7) implies only a weak dependence of M on D, which turns the Ising model into a convenient test-case for studying the resolution of the methodology. However, it should be emphasised that this model choice does not relate to any inherent limitations of the identification methodology and is made mostly for computational convenience. The identification approach can be adopted to different materials outside of the Ising class simply by selecting the appropriate critical exponents and scaling functions (table 1).
The values of the critical exponents β, ν, scaling function μ, and T b c required in equation (2) are available in the literature for a broad range of materials. Alternatively, they can be obtained from independent simulations of the critical behaviour. To demonstrate the procedure, the Monte-Carlo simulations of the zero-field Ising model were used here starting from the paramagnetic phase at high temperature T and proceeding to ferromagnetic phase at low temperatures in small decrements. At each temperature, the system was equilibrated during 10 000 Monte-Carlo sweeps and magnetisation M(T) computed by averaging the spins over the lattice. Figure 1(b) shows the computed data for three Ising lattices of different size D. The form of the scaling function μ can be found by collapsing these data through equation (7) using established techniques [24]. The resulting μ is shown in figure 1(d) and its empirical expression given in table 1. The scaling procedure also generated optimised numerical values β = 0.150 600 ± 8 × 10 −6 , ν = 1.0256 ± 10 −4 , and T b c = 2.264 7583 ± 5 × 10 −7 , which compare well with the theoretical values in table 1. The procedure can also be applied to identify the critical exponents and scaling functions of previously unidentified materials. Alternatively, equation (2) can be fitted to m(T) data of a granular assembly directly without a priori knowledge of the critical parameters. However, this may lead to increased fit parameter correlation.

Discussion
As benchmark granular model system for validating equation (2), four ensembles of independent planar Ising spin lattices of circular shape imitating grains were considered (figure 1(a)). The diameters of grains D were drawn from the lognormal grain size distribution (9) with D = 100 spins and σ D / D = 10%, 20%, 30%, 40%. Each ensemble consisted of 95 000 grains of varying size D. The temperature dependence of magnetic moment, m(T), for any given ensemble with specific σ D was evaluated by using the Monte-Carlo method and by superimposing the contributions from each grain in the ensemble. Figure 2(a) shows examples of m(T) for granular ensembles with σ D / D = 10% and 40%.
Equations (2) and (4) were least square fitted to the generated m(T) data. During the fitting procedure T b c , the critical exponents β and ν, and the scaling function μ of the Ising model were used (table 1). The parameters t c and σ tc of f t (t c ) were considered to be the only fit parameters varied during the fitting procedure, which were then converted to T c and σ Tc . Figure 2(b) shows the fitted σ Tc / T c plotted as a function of reference σ Tc / T c obtained by binning the values of D and finding the average value of T c (D) of the susceptibility peaks associated with each bin. The resulting trend deviates from the expected result, represented by the diagonal dashed line, which suggests strong parameter correlation between the t c and σ tc during the fitting.
Combining t c = −D/ν and σ tc =σ D /ν gives σ tc /t c = −σ D /D, which after using the standard properties of lognormal distributions (see the discussion of equations (4) and (9)) leads to the arithmetic form: Applying formula (10) simultaneously with equations (2) and (4) during the fitting procedure now allows unambiguous identification of the T c distribution as shown in figure 2(b). Given that information about the granular or particle size distributions such as D and σ D is typically available from independent experiments [15,[25][26][27], the introduction of the constraint in (10) does not pose any conceptual challenge. The identification approach developed here has been based strictly on the zero field analysis, which is when genuine second order phase transitions can be observed in bulk systems.
In recent experiments the T c distribution has been associated with measurements of the so-called switching temperature distribution in HAMR media [5,6]. The measurement principle, based on spatially resolved tracking of temperatures at which grains reverse, required significant reversal magnetic field to bring the frequency of sub-critical thermal fluctuations into the time window of the measurement apparatus. This brings about the question of the effect of external magnetic fields in extracting the genuine T c distributions.
The effect of the applied field H can be studied by considering simple mean-field theory, where the magnetisation of a single grain is represented as M(H, T) = tanh((VM s H + J mf M)/k B T), with J mf being the mean-field interaction (Weiss molecular field) in the units of energy. Setting H = 0, the critical temperature of a grain can be found to be T c = J mf /k B , and |M| 1 for T < T c , whereas M = 0 for T > T c . Taking the peak of the derivative dM(H, T)/dT as an estimate of the fictitious H-dependent Curie temperature T h c of a grain, differentiating and arranging gives T h c = T c + (J −1 mf VM s M peak H)T c , where M peak is the magnetisation associated with the peak at T h c . Solving this equation together with the mean-field formula for M(H, T h c ), assuming small magnetisation |M| << 1, and expanding the solution to the first order gives T h c ≈ T c + q|H|T c , where the constant q = J −1 mf VM s . Then, averaging over a statistical ensemble of grains allows to express the mean and standard deviation of the T c distribution as T h c = T c + q|H| T c and σ h c = σ c + q|H|σ c . Although these simple calculations based on the mean-field theory inherently do not incorporate thermal fluctuations and finite size effects included in equation (7), they qualitatively suggest that the presence of external field leads to apparent shift and broadening of the T c distribution. This is consistent with experimental reports overestimating the zero-field theoretical calculations [7][8][9]13].

Conclusions
The developed approach based on least-square fitting a typical experimentally measurable m(T) dependence by expression (2), possibly with constraint (10), has been demonstrated to allow for the extraction of the T c distribution in noninteracting magnetic particle assemblies. The constraining relation (10) might prove unnecessary in general materials with more pronounced dependence of m(T) on σ D / D .
Most of the similar previously developed T c distribution identification techniques assume bulk relations for describing thermodynamic state functions near the T c [7][8][9][10], and thereby cannot systematically incorporate the finite size effects of grains. Validating these techniques against a theoretically consistent physical picture, such as that which underlies the approach presented in this article, goes beyond the scope of the present work and will require consideration in the future.
The present approach can be systematically adapted to different materials by specifying relevant critical exponents and scaling functions. These can be found in broad literature or, alternatively, obtained through independent experiments or simulations of simplified models in the same universality class by using the finite size scaling procedure illustrated here based on the Ising model. The methodology can be consistently extended to account for the anisotropy distributions and effects of surface disorder by incorporating nonuniversal corrections in equation (7). It can also be extended to include external magnetic field through generalised forms of the field-dependent scaling functions. In addition, the given derivation of key formulas presents a recipe for including the effects of inter-granular interactions, such as demagnetizing fields. This thus opens prospects for developing identification methodologies for broad class of granular and particulate systems.