Theory and PDE

PenguinStefan.jl advances one-phase and two-phase Stefan problems on Cartesian cut-cell meshes. The solver couples heat diffusion with interface motion.

For each phase, temperature obeys a diffusion equation

\[\partial_t T = \nabla \cdot (D \nabla T) + s,\]

with boundary/interface conditions provided through PenguinBCs.jl and the cut-cell operators from CartesianGeometry.jl and CartesianOperators.jl.

The interface speed is computed from Stefan balance (normal heat-flux jump) and used to update the geometric representation (AbstractInterfaceRep, e.g. LevelSetRep).

For Gibbs-Thomson thermodynamics, the trace law supports isotropic and 2D anisotropic variants:

\[T_\Gamma = T_m - \sigma_{\mathrm{eff}}\,\kappa_\Gamma - \mu_{\mathrm{eff}}\,V_\Gamma\]

with optional harmonic anisotropy factors from PenguinBCs.HarmonicAnisotropy: 1 + \epsilon\cos(m(\theta-\theta_0)), and optional capillary stiffness factor 1 + \epsilon(1-m^2)\cos(m(\theta-\theta_0)).

Geometric representations

The Stefan solver is coupled to an AbstractInterfaceRep backend. The current public backends are:

  • LevelSetRep (general-purpose geometry handling),
  • GlobalHFRep for graph-like interfaces,

For GlobalHFRep, the promised 3D scope is graph interfaces representable as x = x_f(y, z, t) after selecting the graph axis (equivalently y = y_f(...) or z = z_f(...) for other axis choices).

Implemented model families

  • Monophasic Stefan: StefanMonoProblem, StefanMonoState, step!.
  • Diphasic Stefan: StefanDiphProblem, StefanDiphState, step!.
  • Similarity/manufactured references for validation in 1D, 2D, and 3D graph cases.

LevelSet inner iteration

This table shows a sample LevelSet inner-iteration history (absolute interface error and final coupling residual) from a :ls_newton run.

Itersnum - srefFinal stefan residual
00.06801644907936871-
10.00117295687753796170.0014828223090695009
20.000200611704756936689.652345108024374e-5
40.00032499574184899354.584366535692608e-7