Diffusion Theory

 Introduction

Diffusion theory is the modeling of photon transport due to photon movement down concentration gradients. Diffusion theory is appropriate in medium dominated by scattering rather than absorption so that each photon undergoes many scattering events before being terminated by an absorption event. The photon has a relatively long residence time which allows the photon to engage in a random walk within the medium.













_{Diffusion theory = Photons falling down a gradient of  photon concentration.}
The instantaneous fluence rate, F(r,t) [J s^-1 cm^-2] or [W cm^-2], is proportional to the concentration of optical energy, C(r,t) [J cm^-3] and the speed of light, c [cm/s], in the medium:

F(r,t) = cC(r,t)

The units of energy concentration [J cm^-3] times the units of velocity [cm/s] yield the units of fluence rate [J s^-1 cm^-2]. This relationship between F and C allows us to discuss the time-resolved spatial distribution of photon concentration, C(r,t), in a light-scattering medium using the same math of diffusion which applies to many things such as solutes in a solution or heat in a material.

Fick's 1st law of diffusion

Diffusion occurs in response to a concentration gradient expressed as the change in concentration due to a change in position, . The local rule for movement or flux J is given by Fick's 1st law of diffusion:











in which the flux J [cm^-2 s^-1] is proportional to the diffusivity  [cm²/s] and the negative gradient of concentration,  [cm^-3 cm^-1] or [cm^-4]. The negative sign indicates that J is positive when movement is down the gradient, i.e., the negative sign cancels the negative gradient along the direction of positive flux.


The flux J is driven by the negative gradient   in the direction of increasing x.

Example arbitrary units

Consider a local concentration of 106 units per cm3 which drops by 10% over a distance of 1 cm. Then the gradient is -105 [units cm-4]. Assume the diffusivity is 103[cm2/s]. Then the flux J equals: J = -(103 [cm2/s])(-105 [units cm-4]) = 108 units cm-2 s-1

For light, the diffusivity is proportional to the diffusion length D [cm] and the speed of light c:

 = cD

where D = 1/(3 µ_s(1-g)). The units of velocity [cm/s] times the units of length [cm] yield the units of diffusivity [cm²/s]. The following example describes the local diffusion of red light in milk.

Example optical energy

Consider a local fluence rate F of 1 W/cm2 in milk (n = 1.33, µs(1-g) = 10 cm-1). The local concentration of optical energy is C = F/c = (1 W/cm2)/(2.2x1010 cm/s) = 4.4x10-11 J cm-3. Assume the concentration drops by 10% over a distance of 3 mm. Then the gradient  is -(0.10)(4.4x10-11 cm-3)/(0.3 cm) = -1.5x10-11 J cm-4. The diffusivity  equals cD where D = 1/(3(10 cm-1)) = 0.033 cm. Therefore  = (2.2x1010 cm/s)(0.033 cm) = 7.5x1010 cm2/s. Then the flux J along the direction x equals - = -(7.5x1010 cm2/s)(-1.5x10-11 J cm-4) = 0.11 J cm-2 s-1 = 0.11 W cm-2.

---

For optical diffusion, Fick's 1st law is expressed as the energy flux J [W cm^-2] proportional to the diffusion constant D [cm] and the negative fluence gradient dF/dx:





which was obtained by substituting cD for  and substituting F/c for C. The factors c and 1/c cancel to yield the above equation.

Fick's 2nd law of diffusion

Consider diffusion at the front and rear surfaces of an incremental planar volume. Fick's 2nd law of diffusion describes the rate of accumulation (or depletion) of concentration within the volume as proportional to the local curvature of the concentration gradient. The local rule for accumulation is given by Fick's 2nd law of diffusion:











in which the accumulation, dC/dt [cm^-3 s^-1], is proportional to the diffusivity  [cm²/s] and the 2nd derivative (or curvature) of the concentration,  [cm^-3 cm^-2] or [cm^-5]. The accumulation is positive when the curvature is positive, i.e., when the concentration gradient is more negative on the front end of the planar volume and less negative on the rear end so that more flux is driven into the volume at the front end than is driven out of the volume at the rear end.


Incremental planar volume accumulates concentration because the front gradient at x ₁ drives more flux J₁ into the volume than the flux J₂ driven out of the volume by the rear gradient at x₂.

---

The differential equation for optical diffusion is simply Fick's 2nd law with the substitution of the product cD for the diffusivity  and substitution of F/c for concentration C, although the 1/c factors introduced on both sides of the equation cancel:

Connecting optical transport theory and Fick's 1st law of diffusion.

The early work on neutron scattering theory in nuclear reactors first developed the connection between transport theory and Fick's 1st law of diffusion. The resulting statement of diffusion is also applied to optical transport. Consider the following problem:

Light scattered at position r passes through a small aperture with area A contributing to the flux J+.

• given a homogeneous isotropically scattering medium with scattering coefficient µ_s
• given a small aperture with area A at the origin (r = 0)
• given a fluence rate F(0) at the origin
• assume that the fluence rate F(r) at a position r near to the aperture is approximated by
  
• calculate the net flux through the small aperture with area A due to scattering from all the surrounding volume:
  
  
  
  
  • the scattered flux at r is µ_sF(r) = µ_s( F(0) + r(dF(r)/dr_at r=0) ) where r is the distance from r to the origin.
  • the fraction of scattered flux from r that survives a pathlength r without being scattered or absorbed is exp(-µ_tr), where µ_t = µ_s + µ_a. Consider the case of negligible absorption compared to scattering so µ_t = µ_s.
  • the fraction of surviving scattered flux from r that passes through A is Acos/(4r²).
  • integrate the flux from r over all possible  and r for the hemisphere above the aperture:
    
  • a similar integration for the scattering from the hemisphere below the aperture yields a positive value J_- for flux passing through A in the other direction opposite to J₊.
  • the net flux J = J₊ - J_-
• The final expression for the net flux J has the form of Fick's 1st law of diffusion if the diffusion constant D has the following value:
   , where 

---

For photon transport, the term µ_s' = µ_s(1 - g) is substituted the above µ_s which described isotropic scattering used in the above treatment.
In summary, the connection between transport theory and Fick's 1st law of diffusion is based on the linearization of F(r) around the value F(0) and on the approximation of D by the value 1/(3 µ_s').

---

Limits of diffusion theory

There are limits to the applicability of diffusion theory.
It is important to remember that diffusion theory simply attempts to mimic Fick's 1st law of diffusion by a proper choice of D to relate J and -dF/dr. The approximation of F(r) as simply a linear perturbation of the value F(0) neglects higher order terms that depend on d²F/dr². Sometimes the gradient has significant curvature within the spherical regime from which exp(-µ_s' r) allows photons to pass through the aperture. In such regions, the linear approximation to F(r) is inadequate and diffusion theory is inaccurate. Two cases of curvature deserve emphasis:

Near a source

The gradients are very steep near a point source or collection of point sources with both the exponential term exp(-r/(4cDt)) and the term 1/r causing significant curvature in the gradients. Distant from sources, gradients become gradual and diffusion theory is more accurate.

Absorption

Strong absorption prevents photons from engaging in an extended random walk. The approximation µ_t = µ_s' becomes inadequate. In common use is the approximation:



However, we have conducted Monte Carlo simulations of light diffusion from an isotropic point source to yield the spatial distribution for F(r) which suggests a different value for D:
If µ_a << 3µ_s', the effect of absorption can be neglected.

Time-resolved Diffusion Theory

The solution to Fick's 2nd Law of diffusion for the case of a point source of energy, U_o [units], deposited at zero time (t = 0) at the origin (r = 0) is C(r,t) [units cm^-3]:











where r [cm] is the distance from the source to the point of observation, t [s] is the time of observation, and  [cm²/s] is the diffusivity. In the denominator, the product t has units of [cm²] which is taken to the 3/2 power to yield units of [cm³]. The exponential in the numerator is dimensionless. The units of U_o can be the units of any extensive variable which is able to diffuse throughout a volume. Hence, C(r,t) has the proper dimensions for concentration, [units cm^-3]. The above expression for C(r,t) is a solution to Fick's 2nd law of diffusion and describes spherically symmetric diffusion from a point impulse source in a homogeneous medium with no boundaries.

---

For the case of optical diffusion, let F = cC and  = cD. Let the point source be an impulse of energy U_o [J]. The expression for fluence rate F(r,t) [W cm^-2] becomes:



The above equation is very useful and the student of tissue optics should know this equation very well.

Time-resolved diffusion Theory:

Examples

The following figures illustrate the time-resolved transport of light from an impulse isotropic point source of energy within a homogeneous unbounded medium with absorption and scattering properties.

The time-resolved fluence rate F(r,t) is plotted versus r at three timepoints: t = 5, 100, and 1000 ps. medium: µa = 1.0 cm-1, µs(1 - g) = 10.0 cm-1, nt = 1.33 impulse energy: Uo = 1 J	Click on figure to enlarge
The same data plotted as a semi-log plot.	Click on figure to enlarge

The C program that generated the data in the above figures is listed here.

Steady-state Diffusion Theory

The time-resolved fluence rate F(r,t) [W/cm²] in response to an impulse of energy U_o [J] and the steady-state fluence rate F_ss(r) [W/cm²] in response to an isotropic point source of continuous power P_o [W] are summarized:

---

---

where the transport factors T(r,t) [cm^-2 s^-1] and T_ss(r) [cm^-2] have been introduced to more carefully distinguish the source, the transport, and the fluence rate. Note on notation: In this class, we use F_ss(r) rather than F(r) to especially emphasize the steady-state fluence rate from the time-resolved fluence rate. However, F(r) should be used outside this class.
The above expression for T_ss(r) can be obtained by integrating T(r,t)exp(-µ_act) over all time to yield the total accumulated amount of photon transport to each position r. The factor exp(-µ_act) accounts for photon absorption (recall that ct = pathlength so this expression is simply Beer's law for photon survival) and causes photon concentration to approach zero as time goes to infinity. The expression for T_ss(r) is derived:



Note that the final expression above has made the substitutions:
which removes the diffusion length D and introduces the optical penetration depth  which is the incremental distance from the source that causes F_ss(r) to decrease to 1/e its initial value. The penetration depth  is a parameter which is very easily understood in experimental measurements and consequently has more intuitive value to some people than D which is important from the perspective of the local step size of the diffusion process. In this class we will often use the following expression for F_ss(r) when we prefer to emphasize the roles of µ_aand :

Steady-state diffusion Theory:

Examples

The following figures illustrate the transport of light from a steady-state isotropic point source of power within a homogeneous unbounded medium with absorption and scattering properties.

The steady-state fluence rate F(r) is plotted versus r. medium: µa = 1.0 cm-1, µs(1 - g) = 10.0 cm-1, nt = 1.33 impulse power: Po = 1 W	Click on figure to enlarge
The same data plotted as a semi-log plot.	Click on figure to enlarge

The C program that generated the data in the above figures is listed here.

Frequency-domain Diffusion Theory:

Isotropic point source in infinite medium

A light source may be modulated sinusoidally to yield a sinusoidally varying fluence rate distribution at a distant observation point within a medium. Such a modulated concentration will propagate in the medium and is often called a photon density wave. Consider a modulated isotropic point source of light within a homogenous turbid medium with no boundaries.
The point source S has a steady-state power S_o [W] which is modulated sinusoidally by a modulation factor M_osin(t), where 0 < M_o < 1:











where  = 2f [radians/s] is the angular frequency of modulation, and the modulation frequency f is in hertz [cycles/s]. The position of observation r is located a distance r from the source. The above equation shows two ways to express S(t), one using a sine function and the other using the equivalent and well-known convention of an exponential with an imaginary exponent.
In response to this modulated source, the modulated fluence rate F(r,t) at r is described:
where T_ss(r) is the steady-state transport, k" is the imaginary wavenumber that describes the attenuation of the photon density wave, and k' is the real wavenumber that describes the phase lag of the observed photon density wave.
The expressions for T_ss(r), k" and k' are:

---

---

where

---

The behavior of k' and k" are shown in the following figure which plots k' and k" as functions of /(µ_ac):



 Click figure to enlarge

---

The above expressions are equivalent to the expressions published by Schmitt et al. 1992 which in turn are equivalent to the expressions published by Fishkin et al. 1991,1993. Link to references..
At the observation point, the persistence of the source modulation is called the modulation, M, and is often described in the literature as (AC_out/DC_out)/(AC_in/DC_in) which equals:


At the observation point, the phase of the signal lags the phase of the source by an angle called the phase,  [radians], which equals:

---

The ratio /(µ_ac) describes the number of radians of modulation cycle that occur in one mean photon lifetime. Only 1/e or 37% of photons survive after a time period of 1/(µ_ac) [s].  is an angle specified by the ratio /(µ_ac), and approaches zero for low modulation frequencies and approaches 90° at the highest modulation frequencies,  >> µ_ac.

• At very low modulation frequencies,  << µ_ac, there is little modulation during the lifetime of a photon. Consequently, photon migration has little impact on the transport of the modulation. The value of k" approaches one so k" approaches 1/, therefore M approaches unity. The value of k' approaches zero, therefore  approaches zero. The observed modulation of F_ss(r,t) mimics the source modulation with no loss of modulation and no phase lag. F(r,t) approaches the behavior S_oT_ss(r)(1 + M_osin(t)).

• At high modulation frequencies, as the modulation frequency approaches the value of µ_ac or greater, there is significant modulation during the lifetime of the photon. Consequently, the photon can diffuse during its lifetime and thereby smear the spatial resolution of the modulation. F(r,t) exhibits significant loss of modulation (M < 1) and an increased phase lag ( > 0): F(r,t) = S_oT_ss(r)(1 + M_oMsin(t - k'r))

Frequency-domain diffusion theory:

Examples

The following figure illustrates the time-resolved signal F(r,t) in blue and the source S(t) in red. Note the greatly decreased F_ss, the slightly decreased modulation M, and the phase lag  at the observation point.

source: f = 400 MHz, So = 1 W, Mo = 1.0 medium: µa = 1.0 cm-1, µs(1 - g) = 10.0 cm-1, nt = 1.33

Click on figure to enlarge

The C program that generated the data in the above figure is listed here.

---

The following figures illustrates the modulation of the photon density wave as it propagates into the medium.

The factor Msin(-) is plotted versus position r. Two values of absorption coefficient µa are used in the calculations, illustrating that lower absorption yields a longer photon lifetime which allows smearing of the spatial resolution of the modulation which attenuates the modulation. Note that very large distances are needed to develop cyclic modulation due to phase lag. Photon density waves do not present multiple cycles in typical experiments. The period of one cycle is rk'/ [s]. source: f = 400 MHz, Mo = 1.0 medium: µa = 1.0 or 0.1 cm-1, µs(1 - g) = 10.0 cm-1, nt = 1.33

Click on figure to enlarge

The fluence rate F(r,t) is plotted versus position r for a specific time (t = 0 [s] or multiples of rk'/). The value Tss(r) falls very quickly causing F(r,t) to fall quickly which makes the modulation difficult to see, especially on a logarithmic scale. The modulation is more apparent for the lower absorption.

Click on figure to enlarge

The C program that generated the data in the above figures is listed here.

Source:http://omlc.ogi.edu/classroom/ece532/class5/fdexamples.html