Stefan Problem — Phase Change
Calculate interface front positions, phase change rates, and multi-phase temperatures using the Neumann analytical similarity solution.
Analytical Neumann Solution
The classical solution (Franz Neumann, c. 1860) models transient 1D phase change in semi-infinite domains:
- Growth Constant ($\lambda$): Solved iteratively from the heat flux match and latent heat balance at the interface.
- Sensible/Latent Balance: Handles liquid/solid property differences and accounts for superheating/subcooling ($T_i \neq T_m$).
- Front propagation: Proportional to square root of time ($s(t) = 2\lambda\sqrt{\alpha t}$).
Configuration
Results & Visualization
Phase Front & Heat Field Map
Temperature Profile: T(x) across Phases
Live Numerical Summary
Detailed text reports will appear here after clicking "Generate Engineering Report".
Methodology & Theory
1D Stefan Problem (Neumann Solution)
The Stefan problem is a classic moving boundary problem that models transient phase-change processes (solidification or melting) controlled by thermal diffusion. Because heat is absorbed or released during the phase change, the velocity of the front is coupled to the temperature gradients on both sides of the moving interface.
The interface position $s(t)$ is solved using a similarity variable: $$s(t) = 2\lambda\sqrt{\alpha_{grow} t}$$ where $\lambda$ is a dimensionless eigenvalue determined from the interface energy balance (Stefan condition).
Key Reference Formulas:
- Stefan Number (Ste): $$Ste = \frac{C_p |T_m - T_s|}{L}$$ Measures the ratio of sensible heat capacity of the growing phase to the latent heat of phase change.
- Solidification Transcendental Equation: $$\frac{Ste_s}{e^{\lambda^2}\text{erf}(\lambda)} - \frac{k_l \sqrt{\alpha_s}(T_i - T_m)}{k_s \sqrt{\alpha_l}(T_m - T_s) e^{\lambda^2 \alpha_s/\alpha_l} \text{erfc}(\lambda \sqrt{\alpha_s/\alpha_l})} \cdot Ste_s = \lambda \sqrt{\pi}$$
- Melting Transcendental Equation: $$\frac{Ste_l}{e^{\lambda^2}\text{erf}(\lambda)} - \frac{k_s \sqrt{\alpha_l}(T_m - T_i)}{k_l \sqrt{\alpha_s}(T_s - T_m) e^{\lambda^2 \alpha_l/\alpha_s} \text{erfc}(\lambda \sqrt{\alpha_l/\alpha_s})} \cdot Ste_l = \lambda \sqrt{\pi}$$
Academic References:
- Carslaw, H. S., & Jaeger, J. C. (1959). Conduction of Heat in Solids (2nd Edition). Oxford University Press. Chapter 11.
- Alexiades, V., & Solomon, A. D. (1993). Mathematical Modeling of Melting and Freezing Processes. Hemisphere Publishing Corporation.
Worked Engineering Example
A container of liquid water is initially at a uniform temperature of $T_i = 10$°C. The surface is suddenly exposed to a cold plate maintained at $T_s = -15$°C, starting solidification. Calculate the thickness of the ice layer formed after $t = 1$ hour ($3600$ s).
Step-by-step Solution:
1. Identify properties for water/ice:
- Ice (solid): $k_s = 2.22$ W/m·K, $\rho_s = 917$ kg/m³, $C_{p,s} = 2050$ J/kg·K
- Water (liquid): $k_l = 0.57$ W/m·K, $\rho_l = 1000$ kg/m³, $C_{p,l} = 4184$ J/kg·K
- Phase change: $T_m = 0$°C, $L = 333.5$ kJ/kg
2. Compute thermal diffusivities:
$$\alpha_s = \frac{k_s}{\rho_s C_{p,s}} = \frac{2.22}{917 \times 2050} = 1.181 \cdot 10^{-6} \text{ m²/s}$$ $$\alpha_l = \frac{k_l}{\rho_l C_{p,l}} = \frac{0.57}{1000 \times 4184} = 1.362 \cdot 10^{-7} \text{ m²/s}$$ 3. Evaluate Stefan number for solid ice:
$$Ste_s = \frac{C_{p,s}(T_m - T_s)}{L} = \frac{2050 \times (0 - (-15))}{333.5 \times 10^3} = 0.0922$$ 4. Solve the solidification transcendental equation for $\lambda$:
Using iterative root-finding, we find: $$\lambda \approx 0.2014$$ 5. Calculate ice thickness $s(t)$ after 1 hour ($3600$ s):
$$s(t) = 2 \lambda \sqrt{\alpha_s t} = 2 \times 0.2014 \times \sqrt{1.181 \cdot 10^{-6} \times 3600}$$ $$s(3600) = 2 \times 0.2014 \times 0.0652 \approx 0.0262 \text{ m} = 26.2 \text{ mm}$$
Final Result:
The thickness of the ice layer formed in 1 hour is 26.2 mm.