Supercooled Stefan Problem

• The supercooled Stefan problem is a free-boundary model for phase transitions, characterized by finite-time blow-ups, jump phenomena, and fractal interface features.
• Modern approaches leverage classical PDEs, probabilistic reformulations, and obstacle-problem techniques to analyze and overcome singular behaviors and nonuniqueness.
• Regularization via kinetic undercooling and advanced computational methods bridge analytical theory and numerical approximations for reliable solution behavior.

Searching arXiv for the specified paper and closely related work on the supercooled Stefan problem to ground the article in the cited literature. The supercooled Stefan problem is a free-boundary model for the freezing of a liquid that is initially below its equilibrium freezing temperature. In its classical form, heat diffuses in the liquid region while the liquid–solid interface advances according to a Stefan flux law. What distinguishes the supercooled regime from the ordinary melting problem is the onset of singular behavior: the interface velocity can blow up in finite time, the free boundary can jump, new frozen regions can nucleate, and in higher dimensions the freezing set can display fractal features. Modern analysis therefore treats the problem not only as a classical PDE/free-boundary system, but also through probabilistic, obstacle-problem, and optimal-transport formulations that remain meaningful after classical breakdown (Delarue et al., 2019, Engelstein et al., 10 Dec 2025, Chu et al., 10 Dec 2025).

1. Classical free-boundary formulations

In one spatial dimension, the one-phase supercooled Stefan problem is commonly written on the half-line in terms of a temperature variable $u(t,x)$ and an interface $s(t)$, with

$$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$

together with

$$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$

and prescribed initial data $u(0,x)=u_0(x)$, $s(0)=0$ (Delarue et al., 2019). In this normalization, $\alpha>0$ is the latent-heat–inverse-conductivity constant. The same model is often rewritten in enthalpy form, for example as

$(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$

in one dimension, or with a different diffusion normalization in higher-dimensional formulations (Chau et al., 1 Jan 2026, Choi et al., 2024).

These formulations encode the same phase-transition mechanism with different state variables. In the classical PDE picture, the interface condition $u=0$ fixes the freezing temperature and the Stefan law identifies interface motion with boundary flux. In the enthalpy picture, the latent heat is absorbed into the indicator $\chi_{\{\eta>0\}}$, so the free boundary is carried implicitly by the positivity set. This duality is central in later developments, because obstacle-type methods and weak formulations are far more robust once blow-ups or jumps occur.

In higher dimensions, the problem is posed in a container $s(t)$0 with liquid region $s(t)$1 and free boundary $s(t)$2. In the liquid one has $s(t)$3, while the interface velocity satisfies

$s(t)$4

Under mild hypotheses, the free boundary may be written as a space-time graph $s(t)$5, where $s(t)$6 is the freezing-time function (Engelstein et al., 10 Dec 2025). The higher-dimensional theory makes explicit that the supercooled problem is not merely a one-dimensional instability: it allows infinite-speed propagation, nucleation, and jump phenomena even when a global graph structure persists.

2. Probabilistic reformulation and solution concepts

A decisive development was the reformulation of the one-dimensional problem as a McKean–Vlasov hitting-time equation. Given $s(t)$7 with density $s(t)$8 and a Brownian motion $s(t)$9, one seeks a nondecreasing right-continuous process $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$0 such that

$$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$1

and

$$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$2

On intervals where $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$3 is differentiable, the density $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$4 of $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$5 solves

$$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$6

and then $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$7, $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$8 recover the Stefan problem (Delarue et al., 2019).

This probabilistic representation permits several distinct solution notions.

Notion	Defining feature	Source
Classical solution	Smooth free boundary; PDE and Stefan law hold pointwise	(Delarue et al., 2019)
Physical solution	At each jump, the boundary takes the minimal jump compatible with preservation of mass	(Munoz, 23 Jun 2025)
Minimal solution	Smallest fixed point of the McKean–Vlasov map $$\partial_t u(t,x)=\tfrac12\,\partial_{xx}u(t,x), \qquad x>s(t),\ t>0,$$9	(Cuchiero et al., 2020)
Maximal solution	Selected by maximizing a suitable notion of average freezing time via free-target transport	(Chu et al., 10 Dec 2025, Chau et al., 1 Jan 2026)

For minimal solutions, one defines $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$0 for the process driven by a candidate loss profile $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$1, and iterates from $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$2. The increasing limit yields a global minimal solution. Under the integrability condition $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$3, the minimal solution is physical, so minimality and the minimal-jump rule coincide (Cuchiero et al., 2020). This equivalence closes a conceptual gap between order-theoretic and physically selected continuations.

A different variational selection appears in the free-target optimal Brownian transport approach. There the target measure is not prescribed in advance but optimized under density constraints, producing maximal weak solutions in a stochastic order. In one dimension, the maximal weak solution is unique; in higher dimensions, existence is available but universality can fail (Choi et al., 2024, Chau et al., 1 Jan 2026, Chu et al., 10 Dec 2025). The coexistence of physical, minimal, and maximal notions is not a terminological redundancy: it reflects genuinely different selection principles for an intrinsically singular free-boundary problem.

3. Blow-up, jumps, nucleation, and regularity

The signature analytical difficulty is finite-time blow-up of the freezing rate. In the McKean–Vlasov formulation, blow-up occurs when the boundary density reaches a critical threshold. More precisely, a discontinuity time is characterized by the condition that the jump $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$4 is the minimal $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$5 such that

$$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$6

while in PDE variables the threshold is detected by

$$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$7

for instantaneous blow-up (Delarue et al., 2019).

This leads to a tripartite regularity classification in one dimension. In the sub-critical regime, $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$8, the interface is $$u(t,s(t))=0, \qquad s'(t)=\tfrac{\alpha}{2}\,\partial_x u\bigl(t,s(t)^+\bigr),$$9 to the right of $u(0,x)=u_0(x)$0. In the critical regime, $u(0,x)=u_0(x)$1 but $u(0,x)=u_0(x)$2, the interface is $u(0,x)=u_0(x)$3-Hölder and obeys the square-root law

$u(0,x)=u_0(x)$4

In the super-critical regime, $u(0,x)=u_0(x)$5, the interface has a positive jump (Delarue et al., 2019). The square-root behavior recovers the law already noted by Stefan in 1889 for the ordinary Stefan problem.

Recent one-dimensional regularity theory sharpens this picture for physical solutions. Assuming only that the initial temperature is integrable, the free boundary, although allowed to jump as a function of time, is $u(0,x)=u_0(x)$6 as a function of space, and $u(0,x)=u_0(x)$7 outside of a closed, countable set. The set of positive jump times is locally finite, and therefore cannot accumulate. The associated freezing time

$u(0,x)=u_0(x)$8

admits a regular/singular classification: the singular set is countable, while the regular set is an open $u(0,x)=u_0(x)$9 manifold (Munoz, 23 Jun 2025).

In arbitrary dimensions, the free boundary also has a robust structure. It decomposes into

$s(0)=0$0

where the regular part is a $s(0)=0$1 hypersurface moving with finite speed, the singular part consists of points with infinite speed or nucleation and satisfies $s(0)=0$2, and the jump set occurs only at a finite or countable collection of isolated times and lies on a space-time smooth graph. Moreover, the freezing-time map satisfies $s(0)=0$3, and singular points coincide with the critical set $s(0)=0$4 (Engelstein et al., 10 Dec 2025).

These structural theorems coexist with strong nonuniqueness phenomena for general weak solutions. In arbitrary dimensions one can force fractal freezing on prescribed closed null sets, and one may achieve $s(0)=0$5 at a prescribed time. By contrast, maximal solutions have a transition zone that is open modulo a low-dimensional set and admit a much finer regularity analysis through obstacle-problem theory (Chu et al., 10 Dec 2025). The central lesson is therefore two-sided: the supercooled Stefan problem is highly unstable at the level of unrestricted weak solutions, yet selected solution classes recover substantial geometric regularity.

4. Well-posedness and uniqueness

Global existence in one dimension was established through the probabilistic formulation under hypotheses that the initial density is bounded and changes monotonicity only finitely often on compacts. Under these assumptions there exists a physical global solution through any countable sequence of blow-ups, and such a solution is unique under the minimal-jump convention (Delarue et al., 2019). This was the first uniqueness result covering growth processes with singular self-excitation in the presence of blow-ups.

Subsequent work clarified the relation between local and global uniqueness. For physical solutions, short-time uniqueness implies global uniqueness: if any two physical solutions agree on some interval $s(0)=0$6, then they agree for all $s(0)=0$7. This yields global uniqueness for broad classes of data, including “flat-near-0” profiles and the oscillatory regime studied by Shkolnikov (Munoz, 23 Jun 2025). A complementary sensitivity theory shows that the physical solution map is continuous under weak perturbations of the initial law, provided the underlying datum admits a unique physical solution. For minimal solutions, right-shift perturbations are continuous, whereas left continuity fails unless uniqueness of physical solutions holds; equivalently, two-sided continuity of the minimal-solution map is equivalent to uniqueness of physical solutions (Baker, 2023).

In higher dimensions, Choi–Kim–Kim proved global-time existence of weak solutions for a general class of initial data. Their solution is maximal in a stochastic order among comparable weak solutions with the same initial data, and the construction is based on a free-target optimization problem for Brownian stopping times with a superharmonic cost (Choi et al., 2024). Kim–Kim obtained global-time existence and weak-strong uniqueness for the freezing problem for a well-prepared class of initial domains generated from the initial data, again through a Brownian transport formulation (Kim et al., 2021).

Uniqueness, however, remains strongly dependent on the chosen selection principle. In one dimension, maximal weak solutions are unique, but the same work shows that the supercooled Stefan problem lacks monotonicity and $s(0)=0$8-Lipschitz stability, although it does have stability under weak convergence of measures (Chau et al., 1 Jan 2026). In arbitrary dimensions, maximal solutions are in general non-universal: different strictly superharmonic costs can select distinct maximal solutions (Chu et al., 10 Dec 2025). A persistent source of confusion is therefore the phrase “the solution.” For the supercooled Stefan problem, that phrase is precise only after the admissible class—classical, physical, minimal, or maximal—has been fixed.

5. Regularization, zero kinetic undercooling, and computation

A natural regularization introduces kinetic undercooling. In one dimension, the regularized unknowns $s(0)=0$9 satisfy

$\alpha>0$0

Here the interface temperature is depressed below equilibrium by an amount proportional to its speed. A maximum principle gives $\alpha>0$1, so $\alpha>0$2 remains bounded and no blow-up occurs (Baker et al., 2020).

The zero-undercooling limit rigorously connects this regularized problem to the singular classical problem. Assuming $\alpha>0$3, $\alpha>0$4, $\alpha>0$5, and $\alpha>0$6, the kinetic-undercooling problem has a unique Lipschitz free boundary $\alpha>0$7, the family $\alpha>0$8 increases pointwise as $\alpha>0$9, and

$(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$0

uniformly on compact time intervals, where $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$1 is the unique global probabilistic solution of the classical supercooled Stefan problem without kinetics (Baker et al., 2020). The key tools are a Feynman–Kac formula expressing the free boundary through the local time of a reflected process and a comparison principle in the parameter $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$2. This establishes kinetic undercooling as a physically relevant regularization and provides a bridge between PDE and stochastic-hitting formulations.

On the computational side, the natural explicit Euler time-stepping scheme for the probabilistic formulation converges globally in time to the physical free boundary in the Skorokhod $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$3 topology, even when jumps occur. The analysis also yields an explicit local error bound, while numerical experiments show why $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$4 rather than uniform or $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$5 convergence is appropriate in the blow-up regime (Kaushansky et al., 2020). Implicit approximations improve the jump resolution further: an implicit time-stepping scheme and a Donsker-type fully discrete approximation converge even in the presence of blow-ups, and under stronger assumptions one obtains a convergence rate arbitrarily close to $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$6 (Cuchiero et al., 2022). A different computational direction is the deep level-set method, which parameterizes the interface by a neural-network level-set function trained against the probabilistic Stefan condition; it handles supercooling, topology changes, and surface tension effects (Shkolnikov et al., 2023).

6. Generalizations, interacting particles, and related models

Several recent works place the supercooled Stefan problem inside broader stochastic and continuum frameworks. One line links it to interacting particle systems and aggregation models. In dimension one, solutions to the McKean–Vlasov equation arise as mean-field limits of particle systems interacting through hitting times, yielding propagation of chaos and a rigorous link to systemic-risk and integrate-and-fire models (Cuchiero et al., 2020). In the plane, the scaling limit of external multi-particle DLA satisfies the probabilistic growth rule only as a subsolution inequality in general, while still exhibiting a rigorous connection to classical and weak Stefan solutions (Nadtochiy et al., 2021). For non-integrable initial data, infinite particle systems starting from Poisson point processes produce barrier limits that represent the supercooled Stefan free boundary and allow polynomial growth rates, finite-time blow-up criteria, and critical scaling laws (Blore et al., 22 Jul 2025).

Another line studies stochastic perturbations and control. With transport noise on the half-line, one obtains two weak formulations: one for continuous evolution and one allowing jumps. In the subcritical regime the solution evolves continuously, whereas if part of the initial profile exceeds the critical level $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$7, there is finite-time blow-up with positive probability. The global continuation is again selected by a minimal principle, now phrased as minimal temperature increase (Ledger et al., 2024). In systemic-risk modeling, a drift-controlled version of the supercooled Stefan problem arises as the mean-field limit of bailout policies for defaultable institutions, and the optimal control is computed numerically through a policy-gradient method on a regularized PDE system (Cuchiero et al., 2021).

The model has also been extended beyond pure heat flow. Thermodynamically consistent two-phase Stefan problems with variable surface tension and optional kinetic undercooling generate local semiflows, admit Lyapunov structure through entropy production, and exhibit precise stability criteria for spherical equilibria (Pruess et al., 2011). In a thermomechanical setting with fracture and fluid flow, the Stefan relation is coupled to viscoelasticity, damage, and an $(\eta-\chi_{\{\eta>0\}})_t-\tfrac12\,\eta_{xx}=0$8-theory for the heat equation, with kinetic superheating/supercooling introduced as a relaxation of the sharp Stefan graph (Roubíček, 2020). Variable thermophysical properties and velocity-dependent phase-change temperatures lead to one-phase reductions and finite-difference schemes that remain relevant for nanoparticle melting and dendrite formation (Myers et al., 2019).

Taken together, these developments show that the supercooled Stefan problem is not a single well-posed model in the classical sense, but a family of closely related free-boundary theories organized around a common instability. Classical PDE, McKean–Vlasov hitting times, obstacle problems, free-target transport, and particle-system limits each capture a different aspect of that instability. The current theory is strongest in one dimension for physical solutions and in arbitrary dimensions for selected maximal or graph-type solutions, while nonuniqueness, universality, and the structure of admissible continuations remain central themes in the ongoing analysis of supercooled freezing (Munoz, 23 Jun 2025, Chu et al., 10 Dec 2025, Choi et al., 2024).