Diffusion in solids, porous media and unsteady diffusion

Two questions land on Theo's desk in the same week. Maintenance wants to know how long the rubber gaskets on the ethanol loadout arms can stay in service before ethanol swelling makes them leak. And the dryer operators want to know why the grains near the dryer exit sometimes still feel damp inside even when their surfaces are dry. Both are diffusion problems, but in solids rather than in gases or liquids, and both change with time rather than sitting at steady state. This lesson covers diffusion through porous materials, diffusion into dense solids such as polymers, and the unsteady diffusion equation that describes how concentration profiles develop over time.

Two kinds of solid

Diffusion in solids happens in two very different ways. In a dense, non-porous solid, such as a rubber gasket, a polymer liner or a metal, the diffusing molecules dissolve in the material and move through it by jumping between sites or through the free volume between polymer chains. The flux follows Fick's law with a diffusivity of the solute in the solid, and the concentration just inside the surface is set by a solubility or partition coefficient. In a porous solid, such as an adsorbent bead, a catalyst pellet or a grain particle, the molecules move through a network of pores filled with gas or liquid, and the solid itself is just an obstacle course.

For dense solids the useful quantity is often the permeability, the product of the diffusivity and the solubility. A polymer film of thickness \(L\) with the gas partial pressure \(p_1\) on one side and \(p_2\) on the other passes a steady flux \(N_A = \mathcal{P}(p_1 - p_2)/L\), where \(\mathcal{P}\) is the permeability. Ethanol permeates many rubbers and plastics far more easily than water, which is why the fuel ethanol industry pays close attention to elastomer selection for seals and hoses, a topic Theo met in CHE-280.

Diffusion in porous media

In a porous solid, three effects reduce the flux below what the same pore fluid would carry in open space. Only a fraction of the cross-section, the porosity \(\varepsilon\), is open. The pores wind and branch, so the path is longer than the straight-line distance; the tortuosity \(\tau\), typically 2 to 6, accounts for this. And in very small pores, molecules hit the walls more often than they hit each other. The effective diffusivity combines these effects:

In pores much larger than the mean free path of the gas, about 0.1 µm at atmospheric pressure, the pore diffusivity is the ordinary molecular value. In pores much smaller than the mean free path, molecule-wall collisions dominate, and diffusion follows the Knudsen mechanism, in which each molecule moves independently of the others. The Knudsen diffusivity depends on the pore radius \(r\) and the molecular speed but not on pressure or on the other gas:

In the transition range, where both kinds of collision matter, a simple and widely used approximation adds the resistances, \(1/D_{\mathrm{pore}} = 1/D_{AB} + 1/D_K\), known as the Bosanquet formula. Very narrow pores, comparable in size to the molecules themselves, behave differently again: in zeolite crystals such as the 3A molecular sieves that dry the plant's ethanol, molecules hop between adsorption sites in what is called configurational or surface diffusion. CHE-350 treats adsorption in those beds; here the point is that the rate of uptake by a porous particle depends on its internal diffusion as well as on the gas around it.

Unsteady diffusion

When concentrations change with time, the mole balance on a thin slice gains an accumulation term. For one-dimensional diffusion with constant diffusivity and no reaction, the result is Fick's second law, the diffusion equation:

This is exactly the same equation as transient heat conduction, \(\partial T/\partial t = \alpha\,\partial^2T/\partial z^2\), with the diffusivity in place of the thermal diffusivity. Every solution Theo learned for transient conduction in CHE-312 can be borrowed: the semi-infinite solid, the one-term series for slabs, cylinders and spheres, and the charts of Heisler and Gröber. The dimensionless groups translate directly. The mass transfer Fourier number is \(Fo = D_{AB}t/L^2\), and the mass transfer Biot number is \(Bi_m = k_c L/(K D_{AB})\), where \(K\) is the partition coefficient between the solid and the surrounding fluid. A small Biot number means the outside film controls; a large one means internal diffusion controls and the surface quickly reaches equilibrium with the surroundings.

The semi-infinite solid

A thick slab whose surface concentration is suddenly raised to \(c_s\) behaves, for a while, as if it were infinitely deep. With initial concentration \(c_0\), the solution is the complementary error function:

The penetration depth grows as the square root of time. A convenient marker is \(z = 2\sqrt{D_{AB}t}\), where the argument of erfc equals 1 and the concentration has risen to about 16 % of the surface value. The semi-infinite model holds as long as that depth is small compared with the real thickness. The same erfc solution is the basis of the penetration theory of mass transfer coefficients in lesson 5, where a liquid element at a gas interface is treated as a semi-infinite medium for the short time it stays there.

Finite bodies: drying a particle

For bodies of finite size, the solution is an infinite series, but after a short time only the first term matters. For a sphere of radius \(R\) whose surface is held at the equilibrium concentration, the average fraction of the removable solute still inside, the moisture ratio in drying, is

For \(Fo\) above about 0.1, the first term alone is accurate to within about 2 %, which gives a direct estimate of the time to reach a target moisture. The result explains the dryer operators' observation. A particle loses moisture from its surface quickly, but the inside lags. When the dryer residence time corresponds to a Fourier number near 0.1, the average moisture can be on target while the center is still noticeably wetter. Larger particles are much worse: the time to reach a given moisture ratio scales as \(R^2\), so doubling the particle size quadruples the drying time in the diffusion-controlled period. Lesson 11 returns to this as the falling-rate period of drying.

# Diffusion in porous solids and unsteady diffusion
import math
# Knudsen versus molecular diffusion of ethanol in a porous adsorbent or catalyst pellet
T, M = 393.15, 46.07            # K, g/mol
for r_nm in (2.0, 10.0, 100.0, 1000.0):
    r = r_nm * 1e-9
    DK = 97.0 * r * math.sqrt(T / M)            # m2/s (Knudsen, r in m)
    DAB = 2.0e-5                                 # m2/s molecular, illustrative at 120 C
    D = 1 / (1 / DK + 1 / DAB)                   # Bosanquet combination
    print(f"pore radius {r_nm:7.1f} nm: DK = {DK:.2e}, combined D = {D:.2e} m2/s")
eps, tau = 0.40, 3.0
print(f"effective diffusivity for 10 nm pores = {0.40 / 3.0 * 1 / (1 / (97 * 10e-9 * math.sqrt(T / M)) + 1 / 2.0e-5):.2e} m2/s (eps = {eps}, tau = {tau})")

# Semi-infinite solid: ethanol diffusing into a polymer gasket after a step change at the surface
D_poly = 5e-13                  # m2/s, illustrative value for ethanol in a rubber gasket
def erfc(x): return 1 - math.erf(x)
for t_days in (1, 7, 30):
    t = t_days * 86400
    depth = 2 * math.sqrt(D_poly * t)            # depth where c/cs = erfc(1) = 0.16
    print(f"after {t_days:2d} days the ethanol front (c/cs = 0.16) is {depth * 1000:.2f} mm deep")
x = 1e-3
for t_days in (7, 30, 90):
    t = t_days * 86400
    print(f"c/cs at 1 mm after {t_days:2d} days = {erfc(x / (2 * math.sqrt(D_poly * t))):.3f}")

# Drying of a wet grain particle: one-term solution for a sphere, surface held at equilibrium
D_w = 2.0e-10                   # m2/s, moisture diffusivity (illustrative)
Rp = 1.0e-3                     # m, particle radius
def frac_left(Fo, terms=50):
    # average moisture ratio for a sphere: 6/pi^2 sum exp(-n^2 pi^2 Fo)/n^2
    return 6 / math.pi ** 2 * sum(math.exp(-(n * math.pi) ** 2 * Fo) / n ** 2 for n in range(1, terms + 1))
for t_min in (2, 5, 10, 20):
    Fo = D_w * t_min * 60 / Rp ** 2
    print(f"t = {t_min:3d} min, Fo = {Fo:.3f}, moisture ratio = {frac_left(Fo):.3f}")
# time to reach a moisture ratio of 0.2 with the one-term form
Fo20 = -math.log(0.2 * math.pi ** 2 / 6) / math.pi ** 2
print(f"Fo for moisture ratio 0.2 = {Fo20:.3f}, time = {Fo20 * Rp ** 2 / D_w / 60:.0f} min")

pore radius     2.0 nm: DK = 5.67e-07, combined D = 5.51e-07 m2/s
pore radius    10.0 nm: DK = 2.83e-06, combined D = 2.48e-06 m2/s
pore radius   100.0 nm: DK = 2.83e-05, combined D = 1.17e-05 m2/s
pore radius  1000.0 nm: DK = 2.83e-04, combined D = 1.87e-05 m2/s
effective diffusivity for 10 nm pores = 3.31e-07 m2/s (eps = 0.4, tau = 3.0)
after  1 days the ethanol front (c/cs = 0.16) is 0.42 mm deep
after  7 days the ethanol front (c/cs = 0.16) is 1.10 mm deep
after 30 days the ethanol front (c/cs = 0.16) is 2.28 mm deep
c/cs at 1 mm after  7 days = 0.198
c/cs at 1 mm after 30 days = 0.535
c/cs at 1 mm after 90 days = 0.720
t =   2 min, Fo = 0.024, moisture ratio = 0.548
t =   5 min, Fo = 0.060, moisture ratio = 0.351
t =  10 min, Fo = 0.120, moisture ratio = 0.187
t =  20 min, Fo = 0.240, moisture ratio = 0.057
Fo for moisture ratio 0.2 = 0.113, time = 9 min

When the surface is not at equilibrium

The solutions above assume the surface reaches its equilibrium value instantly, which is the case when the external mass transfer coefficient is large, or the Biot number is large. When the Biot number is small, the outside film controls, the concentration inside the body stays nearly uniform, and the lumped model from heat transfer applies: \(V\,d\bar{c}/dt = -k_c A(c_s - c_\infty)\), giving an exponential approach with time constant \(V/(k_cA)\). Between the two limits, the one-term solutions use eigenvalues that depend on the Biot number, exactly as in CHE-312, and the Heisler charts can be read with mass transfer groups in place of heat transfer ones.