Industrial solids drying simulator

H. K. D. H. Bhadeshia and H. Yan, Queen Mary University of London

For two semi-infinite bars of compositions $ c_1 $ and $ c_2 $ joined at $ x = 0 $, the analytical solution (assuming a constant diffusion coefficient $ D $) is:

c(x, t) = \frac{c_1 + c_2}{2} - \frac{c_1 - c_2}{2}\operatorname{erf}\left(\frac{x}{2\sqrt{Dt}}\right)

At $ x = 0 $, $ \operatorname{erf}(0) = 0 $, meaning $ c(0, t) = \frac{c_1 + c_2}{2} = \text{constant} $ for all $ t > 0 $.

Diffusion couple (left) and a semi-infinite slab where the concentration is maintained constant at the flat surface (right).

The same solution applies to a semi-infinite slab where the concentration is maintained constant at a flat surface. Consequently, the error function describes the early stages of drying for a flat slab, where the surface is maintained at the equilibrium moisture content and the drying front has not yet reached the centreline.

Simulation of non-linear moisture diffusion, falling-rate kinetics, and differential shrinkage stress in porous ceramics, polymers, and food materials. The exponential, power law and constant buttons below are different representations of how the diffusion coefficient varies with moisture content, because the material shrinks as it loses moisture. More details in the "mathematical modelling" link at the bottom of this page.

Notice that stress develops as a part of the drying process. The saturated interior resists the contraction of the rapidly drying exterior skin, generating high surface tensile stresses, and core compression. This can lead to skin peeling and general distortion of the slab.

The theory applies to: (a) cast ceramics (tape-case electronic substrates, concrete); (b) extruded cylindrical pellets, for example, alumina pellets which in their pores support CoMo catalysts used in petroleum refineries. The pellets at first have 60% moisture but then are dried slowly to avoid cracking due to the development of differential shrinkage causing corresponding stresses. After charging with the aqueous catalyst, they have to be dried again. For cylindrical pellets, the method presented here would need to be modified because diffusion occurs in all radial directions. (c) solution-cast polymer films or coatings; (d) slurry cast battery-electrodes.

The image below shows extruded porous γ-Al₂O₃ catalyst support pellets. Each pellet is about 3 mm in length. Compare the solid cylinders versus the trilobe extrudates, that maximise surface area.

Diffusivity model

Simulation time ($t$) 0.0 h (0 s)

Animation playback speed 1.0×

Slab half-thickness ($L$) 10.0 mm

Initial moisture concentration ($C_0$) 0.35 kg/kg

Equilibrium moisture concentration ($C_\mathrm{e}$) 0.04 kg/kg

Base diffusivity ($D_0$) 2.0 × 10⁻¹⁰ m²/s

Concentration coefficient ($\beta$) 4.0

Mass transfer Biot number ($\mathrm{Bi}_\mathrm{m}$) 50.0 (near-equilibrium)

Shrinkage stiffness factor $\left[\frac{E \beta_\mathrm{s}}{1 - \nu}\right]$ 45.0 MPa

Critical tensile strength ($\sigma_\mathrm{crit}$) 6.0 MPa

Mechanical integrity status:
Initialised at zero time. All stresses zero.

Spatial profile: moisture $C\{z, t\}$ and tensile stress $\sigma\{z, t\}$ across half-slab

Moisture concentration $C\{z, t\}$

Equivalent constant $D$

Stress $\sigma\{z, t\}$

Kinetics: drying rate $(-\mathrm{d}\bar{C}/\mathrm{d}t\text{ vs }\bar{C})$ and surface stress history

Drying rate $(-\mathrm{d}\bar{C}/\mathrm{d}t)$

Surface stress $\sigma\{L, t\}$

$\sigma_\mathrm{crit}$ limit

Mathematical modelling, falling-rate physics, and differential shrinkage mechanics

During industrial drying of dense capillary-porous solids (structural ceramics, polymer sheets, food aggregates, and catalyst supports), the drying rate exhibits two distinct regimes: a constant-rate period governed by surface water evaporation, followed by a falling-rate period once the critical surface moisture concentration is reached. In the falling-rate regime, liquid transport from the solid interior to the evaporative boundary is strictly rate-limiting.

Transient non-linear 1D conservation law:

$$\frac{\partial C}{\partial t} = \frac{\partial}{\partial z} \left[ D_\mathrm{eff}\{C, T\} \frac{\partial C}{\partial z} \right]$$

Boundary conditions:

Centreline symmetry at $z = 0$: $$\left.\frac{\partial C}{\partial z}\right|_{z=0} = 0$$
Convective surface flux at $z = L$: $$-D_\mathrm{eff} \left.\frac{\partial C}{\partial z}\right|_{z=L} = k_\mathrm{m} \left[ C\{L, t\} - C_\mathrm{e} \right]$$

Moisture-dependent diffusivity models: As liquid evaporates and pores empty, liquid continuity breaks down into vapour diffusion and capillary meniscus flow, causing effective moisture diffusivity to drop by orders of magnitude:

Exponential model:

$$D_\mathrm{eff}\{C\} = D_0 \exp\left[ \beta (C - C_0) \right]$$

Power-law percolation model:

$$D_\mathrm{eff}\{C\} = D_0 \left( \frac{C}{C_0} \right)^\alpha$$

Differential shrinkage and transient stress generation: Because the exterior surface loses moisture rapidly while the core remains saturated, uneven volumetric shrinkage induces mechanical deformation. When neighbouring elements are constrained by the rigid core, a severe stress gradient develops:

Local free shrinkage strain:

$$\varepsilon\{z, t\} = \beta_\mathrm{s} \left[ C_0 - C\{z, t\} \right]$$

Average slab strain:

$$\bar{\varepsilon}\{t\} = \frac{1}{L} \int_0^L \varepsilon\{z, t\} \,\mathrm{d}z$$

Elastic stress profile under plane-strain constraint:

$$\sigma\{z, t\} = \left( \frac{E}{1 - \nu} \right) \left[ \bar{\varepsilon}\{t\} - \varepsilon\{z, t\} \right] = \left( \frac{E \beta_\mathrm{s}}{1 - \nu} \right) \left[ \bar{C}\{t\} - C\{z, t\} \right]$$

At the exposed boundary ($z = L$), $\varepsilon\{L, t\} > \bar{\varepsilon}\{t\}$, causing severe surface tensile stress ($\sigma > 0$). Conversely, the moist centre ($z = 0$) is forced into compression ($\sigma < 0$). If the surface tensile stress exceeds the critical tensile strength $\sigma_\mathrm{crit}$, catastrophic surface checking, micro-cracking, or structural warping occurs.