rus
Issue #7-8/2018
O.Posnansky
Dynamic diffusion in binary composites measured by nuclear magnetic resonance
Dynamic diffusion in binary composites measured by nuclear magnetic resonance
Self-diffusion of water molecules measured with nuclear magnetic resonance (NMR) explores composite soft matter non-invasively. Detecting structural properties from NMR data is a challenging inverse problem and has to account for variability in shapes, sizes and concentrations of the compartments of soft matter. In this paper we develop a numerical method, providing a fundamental understanding of contrast diffusion mechanism in NMR signal closely linked to both the compartment hierarchy and the temporal spectrum. We model a composite with self-similar random structure and demonstrate that we retrieve the effective dynamic diffusion spectrum as function of geometries of different shapes, sizes and concentrations of disordered compartments of material.
Теги: continues fractions nuclear magnetic resonance renormalization group theory ядерный магнитный резонанс signal attenuation temporal diffusion velocity autocorrelation автокорреляционная функция скоростей временная диффузия затухание сигнала непрерывная дробь ренормализационная группа
INTRODUCTION
Since introduction of diffusion attenuated nuclear magnetic resonance (NMR) it has become a popular marker for detection of subtle differences between composites with important application in material classification and prediction of their physical properties [1]. NMR data from conventional pulsed field gradient spin echo (PGSE) experiments are commonly providing mean diffusivity reflecting packing and size of components of heterogeneous structures. NMR data interpretation must consider a number of effects operating on microscopic hierarchical length scales [2]. These effects range from direct influences on the water mobility such as sizes and shapes of inhomogeneities, to ensemble effects from orientation dispersion or partial volumes of differing types of comprising components.
Sensitizing the data acquisition itself to correlations between multiple parameters is one way to achieve a richer information content from the sample and could be a way forward to improve specificity to individual signal components. Parameter dependencies of such materials usually characterized with various nonlinearities, which especially distinct near the percolation threshold [3]. A number of classical models of the percolation theory are based on random impedance networks (Fig.1a) which have been widely applied to study the transport properties and relaxations in disordered composites. Such models demonstrate the presence of structural transition [4] and possess hierarchical properties of magnetic field distributions (Fig.1b).
In the present work, we propose to apply the renormalization group theory in combination with method of continuous fractions for a unified description of the temporal diffusion in random composites, which are widely used as a lattice model of physical properties in disordered materials [4, 5].
DIFFUSION WEIGHTED SIGNAL ATTENUATION AND MAGNETIC FIELD GRADIENT WAVEFORM
In the magnetic field gradient waveform, G(t), the random and time dependent location of particles modulates the spin phase. As a consequence spin-echo (SE) attenuation, E(t), includes the correlation between velocities of the particle [1,6] and results in an exponential decay:
E(t) = e– α(t), (1)
if a Gaussian approximation of the spin-phase modulation is assumed. The attenuation factor, α(t), in Eq.(1) is yielded as
form_all_eng.ai, (2)
where form_all_eng.ai, and form_all_eng.ai is a velocity autocorrelation function after Fourier transform. Representing magnetic field gradient waveform in a frequency domain, i.e. form_all_eng.ai, it is possible to obtain the attenuation in a form:
form_all_eng.ai, (3)
with gradient spectra, form_all_eng.ai.
The usual spin-echo pulse sequence has time varying magnetic field gradients, G(t), of amplitude g and width δ separated by time Δ (Fig.1d). Spins are excited with 90° radio-frequency (RF) magnetic field pulse, inverted by a 180° RF magnetic field pulse and distinguished by opposite signs of the phase |k| after refocusing. The gradient spectra in this case is yielded [6]
form_all_eng.ai. (4)
Autocorrelation function of velocities (dynamic diffusivity) in general case can be modelled with
form_all_eng.ai, (5a)
which corresponds to the static case in the limit ω → ∞:
form_all_eng.ai. (5b)
For example, if we combine Eq.(4) with Eq.(3) and accept a static diffusion from Eq.(5b) then a classical Stejskal-Tanner result [6],
form_all_eng.ai, (6)
is valid. In the limit Δ >> δ and correlation times, τc, much smaller than time of strong attenuation, and Eq.(6) approaches
form_all_eng.ai . (7)
A HIERARCHICAL MODEL OF RANDOM COMPOSITE
Assume the balance of all fluxes at any location within the composite. Then its discrete version for conservation law is represented as:
form_all_eng.ai, (8)
if ψ1 is a concentration of particles on l th node of lattice. The distribution of local properties, D(lm), in Eq.(8) is governed by binary function:
form_all_eng.ai, (9)
where indicator δ (lm) = {0,1} represents a set E of all possible lattice configurations. Using Eq.(8) we map composite on a lattice of S inner and SГ contact boundary Г nodes. Every bond of this lattice is bounded with the closest neighbour nodes (lm). The bond can either belong to phase (1) or (2) with probabilities p and q satisfied p + q = 1 and characterized by diffusivities D(1) and D(2) correspondingly. The Dirichlet boundary condition in Eq.(8) corresponds Ψm|Г+ = 0 and Ψm|Г– = 1 on the left, Г+, and right, Г–, sides of the lattice respectively. Using Eq.(8) and boundary conditions it is possible to estimate an equivalent diffusivity, form_all_eng.ai, averaged for all realizations of Eq.(9) on a lattice.
We model global (or effective) properties of a composite with randomly distributed local properties according to a real space renormalization group method [9]. The binary function of Eq.(9) in this case is an iterative transformation dependent on index n. Thus the self-similar recursion for Eq.(9) is defined as
form_10_eng.ai, (10)
where the function of equivalent diffusivity is separated into components (k) = {(1), (2)} depending on the connectedness of bonds in the set E.
The transformation of p in Eq.(10) corresponds to connectedness function, R(p), and can be expressed by the polynomial equations
form_11_eng.ai (11)
characterized with unstable points p* = {0,52D, 0,213D}
(Fig.1c). The number of iterations for building a hierarchical lattice is defined by the constraint |pn – pn + 1 | ≤ ε, where ε > 0 is an infinitely small number. Due to this condition, the two components of a composite are indistinguishable after multi-scale hierarchical averaging leading to a convergence of binary properties to solitary with effective properties corresponding to component (1) or (2) depending on selection of initial probability p0 at the first iteration step n = 0.
The averaged equivalent diffusivity can be separated depending on the connectivity of realizations in a set E on every step of iteration. The transfer of the diffusivity from n to the n + 1 iteration level can be carried out according the rules of continued fraction:
form_all_eng_1.ai (12)
(notation of Pringsheim [10]) if p0 < p* .
RESULTS OF NUMERICAL MODELLING USING RENORMALIZATION DIFFUSION COEFFICIENTS
We carried out numerical tests of our hierarchical model for 2D and 3D cases according to Eqs.(8–11). Diffusivities of components (1) and (2) were distributed on H-cell of renormalization composing all possible configurations. In our tests the ratios of initial values of static diffusivities and correlation times were D(1) / D(2) = 103 and τc(1) / τc(2) = 10–1 (Table 1). Extracting frequency dependent dynamic part of velocity autocorrelation function according Eq.(5a), we can verify our iterative procedure at points p → 1 and p → 1 corresponding pristine low- or high-diffusivity networks respectively for 2D and 3D cases (Fig.2a,b). In the double logarithmic scale, log10(D(ω)) ↔ log10(ω), diffusivity increases linearly at low frequency regime and approach a steady state at high frequency regime. The transition to the steady state depends on the correlation time, τc, and different for components (1), (2) at p → 0 and p → 1. With increasing of concentration p of component (1) (p > 0 and p < p*) a variability of configurations provides more contribution, and, as a consequence, effective diffusivity increases nearly linearly in semi-logarithmic scale, log10(D(ω)) ↔ p.
The significant nonlinear jump of effective diffusivity is observed in the vicinity of the unstable point p* and high frequency limit, ω → ∞. On the other hand, behaviour of effective diffusivity changes from increasing to decreasing near the unstable point p* and low frequency limit, ω → 0. As soon as p → 1 effective diffusivity increases linearly in both low and high frequency limits. The behaviour of diffusivity is similar for 2D and 3D cases except a point p*, which depends on the Euclidean space. Because effective diffusivity was transferred according rules of continuous fraction (Eq.(12)) at p0 < p* we can observe a non-integer power law in the interval log10(ω) = (–1,5, –1), i.e. D (ω) ~ ω1,2 (see insert of Fig.2c). At low and high frequency limit the power approaches 2, i.e.D (ω) ~ ω2, and 0, i.e. D (ω) ~ ω0, corresponding Eq.(5) (see insert of Fig.2c). The first, d(log10(D(ω))) / d(log10(ω)), and second, d2(log10(D(ω))) /
d(log10(ω))2, derivatives of dynamic diffusivity in the logarithmic scale allows to estimate the transition between different power laws quantitatively at p → 0 (Fig.2c). It is independent on the Euclidean space according our modelling. Non-integer power law of dynamic diffusivity is absent if p → 1 (Fig.2d).
We can perform inverse Fourier transform of the frequency dependent diffusivity and investigate the time dependence of attenuation assuming Eq.(7). Figs.3a,b demonstrate the time dependent diffusivity for 2D and 3D cases correspondently in semi-logarithmic scale, log (D(t)) ↔ t. In case p → 0 diffusion is nearly time independent. It can be explained by dominating of component (2) characterized with low static diffusion (Table 1). Slight nonlinearity is quantified by the first derivative of the diffusion (Fig.3c). Increasing of the concentration of the component (1) leads to the linear increase of diffusion in the whole range of time interval. The same behavior is observed in case p → 1. Near the unstable point p* diffusion increases nonlinearly if t → 0. In the opposite situation, i.e. t → ∞, diffusion approaches its maximum in p* and then drops. If p → 1 diffusion is linear dependent which can be easily detected due to the high value of static diffusivity of component (1) and low correlation time (Table 1). Qualitatively time dependent diffusion is similar in 2D and 3D except of behaviour at p*.
CONCLUSIONS
In this article, we have developed the effective medium framework for the attenuation of diffusion weighted signal measured with NMR. We introduced dispersive characteristic of diffusion, which incorporates all measurable properties of the heterogeneous random composite. We also demonstrated the use of the real space renormalization method in combination with continued fraction for the heterogeneous diffusion, finding the frequency and time-dependent diffusion coefficient, and showed that the diffusion demonstrates a power-law form in a general case. This form, however, differs qualitatively for different concentrations of composite components. In general, we showed that the relevant parameters of a complex medium are the variance of local diffusivity, as well as the distinct length scales in the spatial diffusivity distribution. ■
Since introduction of diffusion attenuated nuclear magnetic resonance (NMR) it has become a popular marker for detection of subtle differences between composites with important application in material classification and prediction of their physical properties [1]. NMR data from conventional pulsed field gradient spin echo (PGSE) experiments are commonly providing mean diffusivity reflecting packing and size of components of heterogeneous structures. NMR data interpretation must consider a number of effects operating on microscopic hierarchical length scales [2]. These effects range from direct influences on the water mobility such as sizes and shapes of inhomogeneities, to ensemble effects from orientation dispersion or partial volumes of differing types of comprising components.
Sensitizing the data acquisition itself to correlations between multiple parameters is one way to achieve a richer information content from the sample and could be a way forward to improve specificity to individual signal components. Parameter dependencies of such materials usually characterized with various nonlinearities, which especially distinct near the percolation threshold [3]. A number of classical models of the percolation theory are based on random impedance networks (Fig.1a) which have been widely applied to study the transport properties and relaxations in disordered composites. Such models demonstrate the presence of structural transition [4] and possess hierarchical properties of magnetic field distributions (Fig.1b).
In the present work, we propose to apply the renormalization group theory in combination with method of continuous fractions for a unified description of the temporal diffusion in random composites, which are widely used as a lattice model of physical properties in disordered materials [4, 5].
DIFFUSION WEIGHTED SIGNAL ATTENUATION AND MAGNETIC FIELD GRADIENT WAVEFORM
In the magnetic field gradient waveform, G(t), the random and time dependent location of particles modulates the spin phase. As a consequence spin-echo (SE) attenuation, E(t), includes the correlation between velocities of the particle [1,6] and results in an exponential decay:
E(t) = e– α(t), (1)
if a Gaussian approximation of the spin-phase modulation is assumed. The attenuation factor, α(t), in Eq.(1) is yielded as
form_all_eng.ai, (2)
where form_all_eng.ai, and form_all_eng.ai is a velocity autocorrelation function after Fourier transform. Representing magnetic field gradient waveform in a frequency domain, i.e. form_all_eng.ai, it is possible to obtain the attenuation in a form:
form_all_eng.ai, (3)
with gradient spectra, form_all_eng.ai.
The usual spin-echo pulse sequence has time varying magnetic field gradients, G(t), of amplitude g and width δ separated by time Δ (Fig.1d). Spins are excited with 90° radio-frequency (RF) magnetic field pulse, inverted by a 180° RF magnetic field pulse and distinguished by opposite signs of the phase |k| after refocusing. The gradient spectra in this case is yielded [6]
form_all_eng.ai. (4)
Autocorrelation function of velocities (dynamic diffusivity) in general case can be modelled with
form_all_eng.ai, (5a)
which corresponds to the static case in the limit ω → ∞:
form_all_eng.ai. (5b)
For example, if we combine Eq.(4) with Eq.(3) and accept a static diffusion from Eq.(5b) then a classical Stejskal-Tanner result [6],
form_all_eng.ai, (6)
is valid. In the limit Δ >> δ and correlation times, τc, much smaller than time of strong attenuation, and Eq.(6) approaches
form_all_eng.ai . (7)
A HIERARCHICAL MODEL OF RANDOM COMPOSITE
Assume the balance of all fluxes at any location within the composite. Then its discrete version for conservation law is represented as:
form_all_eng.ai, (8)
if ψ1 is a concentration of particles on l th node of lattice. The distribution of local properties, D(lm), in Eq.(8) is governed by binary function:
form_all_eng.ai, (9)
where indicator δ (lm) = {0,1} represents a set E of all possible lattice configurations. Using Eq.(8) we map composite on a lattice of S inner and SГ contact boundary Г nodes. Every bond of this lattice is bounded with the closest neighbour nodes (lm). The bond can either belong to phase (1) or (2) with probabilities p and q satisfied p + q = 1 and characterized by diffusivities D(1) and D(2) correspondingly. The Dirichlet boundary condition in Eq.(8) corresponds Ψm|Г+ = 0 and Ψm|Г– = 1 on the left, Г+, and right, Г–, sides of the lattice respectively. Using Eq.(8) and boundary conditions it is possible to estimate an equivalent diffusivity, form_all_eng.ai, averaged for all realizations of Eq.(9) on a lattice.
We model global (or effective) properties of a composite with randomly distributed local properties according to a real space renormalization group method [9]. The binary function of Eq.(9) in this case is an iterative transformation dependent on index n. Thus the self-similar recursion for Eq.(9) is defined as
form_10_eng.ai, (10)
where the function of equivalent diffusivity is separated into components (k) = {(1), (2)} depending on the connectedness of bonds in the set E.
The transformation of p in Eq.(10) corresponds to connectedness function, R(p), and can be expressed by the polynomial equations
form_11_eng.ai (11)
characterized with unstable points p* = {0,52D, 0,213D}
(Fig.1c). The number of iterations for building a hierarchical lattice is defined by the constraint |pn – pn + 1 | ≤ ε, where ε > 0 is an infinitely small number. Due to this condition, the two components of a composite are indistinguishable after multi-scale hierarchical averaging leading to a convergence of binary properties to solitary with effective properties corresponding to component (1) or (2) depending on selection of initial probability p0 at the first iteration step n = 0.
The averaged equivalent diffusivity can be separated depending on the connectivity of realizations in a set E on every step of iteration. The transfer of the diffusivity from n to the n + 1 iteration level can be carried out according the rules of continued fraction:
form_all_eng_1.ai (12)
(notation of Pringsheim [10]) if p0 < p* .
RESULTS OF NUMERICAL MODELLING USING RENORMALIZATION DIFFUSION COEFFICIENTS
We carried out numerical tests of our hierarchical model for 2D and 3D cases according to Eqs.(8–11). Diffusivities of components (1) and (2) were distributed on H-cell of renormalization composing all possible configurations. In our tests the ratios of initial values of static diffusivities and correlation times were D(1) / D(2) = 103 and τc(1) / τc(2) = 10–1 (Table 1). Extracting frequency dependent dynamic part of velocity autocorrelation function according Eq.(5a), we can verify our iterative procedure at points p → 1 and p → 1 corresponding pristine low- or high-diffusivity networks respectively for 2D and 3D cases (Fig.2a,b). In the double logarithmic scale, log10(D(ω)) ↔ log10(ω), diffusivity increases linearly at low frequency regime and approach a steady state at high frequency regime. The transition to the steady state depends on the correlation time, τc, and different for components (1), (2) at p → 0 and p → 1. With increasing of concentration p of component (1) (p > 0 and p < p*) a variability of configurations provides more contribution, and, as a consequence, effective diffusivity increases nearly linearly in semi-logarithmic scale, log10(D(ω)) ↔ p.
The significant nonlinear jump of effective diffusivity is observed in the vicinity of the unstable point p* and high frequency limit, ω → ∞. On the other hand, behaviour of effective diffusivity changes from increasing to decreasing near the unstable point p* and low frequency limit, ω → 0. As soon as p → 1 effective diffusivity increases linearly in both low and high frequency limits. The behaviour of diffusivity is similar for 2D and 3D cases except a point p*, which depends on the Euclidean space. Because effective diffusivity was transferred according rules of continuous fraction (Eq.(12)) at p0 < p* we can observe a non-integer power law in the interval log10(ω) = (–1,5, –1), i.e. D (ω) ~ ω1,2 (see insert of Fig.2c). At low and high frequency limit the power approaches 2, i.e.D (ω) ~ ω2, and 0, i.e. D (ω) ~ ω0, corresponding Eq.(5) (see insert of Fig.2c). The first, d(log10(D(ω))) / d(log10(ω)), and second, d2(log10(D(ω))) /
d(log10(ω))2, derivatives of dynamic diffusivity in the logarithmic scale allows to estimate the transition between different power laws quantitatively at p → 0 (Fig.2c). It is independent on the Euclidean space according our modelling. Non-integer power law of dynamic diffusivity is absent if p → 1 (Fig.2d).
We can perform inverse Fourier transform of the frequency dependent diffusivity and investigate the time dependence of attenuation assuming Eq.(7). Figs.3a,b demonstrate the time dependent diffusivity for 2D and 3D cases correspondently in semi-logarithmic scale, log (D(t)) ↔ t. In case p → 0 diffusion is nearly time independent. It can be explained by dominating of component (2) characterized with low static diffusion (Table 1). Slight nonlinearity is quantified by the first derivative of the diffusion (Fig.3c). Increasing of the concentration of the component (1) leads to the linear increase of diffusion in the whole range of time interval. The same behavior is observed in case p → 1. Near the unstable point p* diffusion increases nonlinearly if t → 0. In the opposite situation, i.e. t → ∞, diffusion approaches its maximum in p* and then drops. If p → 1 diffusion is linear dependent which can be easily detected due to the high value of static diffusivity of component (1) and low correlation time (Table 1). Qualitatively time dependent diffusion is similar in 2D and 3D except of behaviour at p*.
CONCLUSIONS
In this article, we have developed the effective medium framework for the attenuation of diffusion weighted signal measured with NMR. We introduced dispersive characteristic of diffusion, which incorporates all measurable properties of the heterogeneous random composite. We also demonstrated the use of the real space renormalization method in combination with continued fraction for the heterogeneous diffusion, finding the frequency and time-dependent diffusion coefficient, and showed that the diffusion demonstrates a power-law form in a general case. This form, however, differs qualitatively for different concentrations of composite components. In general, we showed that the relevant parameters of a complex medium are the variance of local diffusivity, as well as the distinct length scales in the spatial diffusivity distribution. ■
Readers feedback
