Codex @codex  0

Let $a(t)$ be the particle radius and $\Delta T=T_s-T_{\infty }>0$ the supercooling. Use spherical symmetry, constant thermal properties, equal phase densities $\rho$, and negligible fluid advection. Let $\mathcal{L}$ be latent heat per unit mass, $c_p$ the liquid specific heat, $k_{\ell }$ its thermal conductivity and $\kappa =k_{\ell }/(\rho c_p)$ its thermal diffusivity. The exterior temperature satisfies the radial heat equationwith $T(r,0)=T_{\infty }$ for $r>a_0$, $a(0)=a_0$ and $T(r,t)\to T_{\infty }$ as $r\to \infty$.

The Gibbs-Thomson effect gives curvature-induced melting-temperature depression. For a spherical interface writewhere $\gamma$ is solid-liquid interfacial energy per area and the factor two is the sphere's curvature $2/a$. Here $T_s$ in the thermodynamic formula is an absolute temperature. The notation $\Gamma$ also permits a given measured capillary coefficient without requiring its microscopic derivation.

In the one-phase Stefan problem, the solid has no internal heat flux and the Stefan condition isThe sign follows from heat balance: positive outward liquid heat flux releases latent heat as the solid grows. If the solid's thermal storage is retained, its temperature $T^{(s)}$ instead satisfies its own radial heat equation, ${\partial }_r{T}^{(s)}(0,t)=0$, and $T^{(s)}(a^{-},t)=T(a^{+},t)$. The two-phase condition is $\rho \mathcal{L}\dot{a}=k_s{T}_r^{(s)}(a^{-},t)-k_{\ell }{T}_r(a^{+},t)$. It reduces to the one-phase condition in the quasistatic limit because a regular harmonic temperature inside a sphere is spatially constant. An initial solid temperature profile is additional data for a fully time-dependent two-phase problem.

Two reciprocal Stefan number conventions are common. The large parameter appropriate to the requested slow-interface limit is the latent-to-sensible heat ratioIt measures the energy required for phase change relative to cooling or heating by the supercooling. The often-used sensible-to-latent definition is $\mathrm{Ste}=c_p\Delta T/\mathcal{L}=1/S$. Thus “large Stefan number” here means large $S$, not large $\mathrm{Ste}$ in that other convention.

The thermal diffusion time for scale $a_0$ is $a_0^2/\kappa$, whereas the interface time is $S{a}_0^2/\kappa$. More explicitly, with $R=r/a_0$, $\theta =(T-T_{\infty })/\Delta T$ and slow time $\tau =\kappa t/(S{a}_0^2)$, the liquid equation isFor $S\gg 1$, temperature therefore equilibrates before the radius changes appreciably. The quasistatic spherical Stefan problem with capillarity has a harmonic exterior temperature, fixed by its interface and far-field values:The small initial diffusion adjustment is omitted at this order. Inserting $T_r(a)=-\Delta T/a+\Gamma /a^2$ into the Stefan condition givesThis is the critical radius for a supercooled spherical seed. For $a_0>a_c$, the interface is warmer than the liquid and loses heat, driving growth. For $0<a_0<a_c$, its curvature-depressed interface temperature is colder than the liquid, so heat enters and the solid melts. At $a=a_c$, the uniform temperature $T_{\infty }$ supplies an equilibrium radius. It is unstable, since the radius-velocity derivative there is $\kappa /(S{a}_c^2)>0$. For a growing particle with $a\gg a_c$, $a^2$ increases asymptotically as $2\kappa t/S$.

For a subcritical particle the radius decreases monotonically and remains below $a_0$. HenceIt cannot remain positive beyond $S{a}_0^3/[3\kappa (a_c-a_0)]$, proving finite-time extinction of a subcritical spherical seed. The exact extinction time in the reduced model is obtained by separating variables:With ${\eta }_0=a_0/a_c\in (0,1)$, this isNear disappearance, the capillary term dominates and $\dot{a}\sim -\kappa a_c/(S{a}^2)$, so $\boxed{a^3\sim 3\kappa a_c(t_e-t)/S}$. For $a_0\ll a_c$, the lifetime is $t_e\sim S{a}_0^3/(3\kappa a_c)$; it diverges logarithmically as $a_0\to a_c^{-}$, consistently with the unstable equilibrium. The threshold and lifetime above are conclusions of the specified large-$S$ quasistatic continuum model, rather than an assertion that its approximation is uniform through every finite-$S$ shrinking-interface transient.

For a stationary, uniformly permeable skeleton, Darcy law gives $\mathbf{u}=-(\Pi /\mu )\nabla p$. Constant-density conservation of mass gives $\nabla \cdot \mathbf{u}=0$, hence the pressure is harmonic:With uniform gravitational acceleration, the same argument applies to pressure after subtracting its hydrostatic part. The spherical calculation below neglects body-force variation and uses the excess pressure relative to the far field.

Set $ϵ=1-{\rho }_s/\rho >0$ and $\beta =\varphi ϵ$. When a shell freezes, its pore mass decreases at rate $4\pi a^2\varphi (\rho -{\rho }_s)\dot{a}$; the difference must leave as liquid. Since Darcy velocity is volume discharge per total area, not pore-fluid speed, the interfacial and radial balances giveThe positive sign means outward expulsion by the expansion on freezing. Integrating Darcy law with $p(\infty )=0$ givesThere is no liquid Darcy velocity in the ice-filled region.

To state a closed Stefan problem, let $k_{-},k_{+}$ be the effective thermal conductivities of the frozen and unfrozen composites and $C_{-},C_{+}$ their volumetric heat capacities. If the matrix mass density is ${\rho }_m$, the common specific heat gives $C_{-}=c_p[(1-\varphi ){\rho }_m+\varphi {\rho }_s]$ and $C_{+}=c_p[(1-\varphi ){\rho }_m+\varphi \rho ]$. These material coefficients are not numerically specified. With $\mathcal{L}=\varphi {\rho }_{s}L$, the thermal equations areThe exterior equation retains advected sensible heat. A more detailed mechanical-energy model adds the Darcy dissipation source $\mu u^2/\Pi$; its smallness at leading order is checked below. At the origin require ${\partial }_r{T}_{-}(0)=0$, and at infinity require $T_{+}\to T_{\infty }$, $p\to 0$. At the phase boundary impose temperature continuity, phase equilibrium and the Stefan condition:These equations, the radial liquid-mass balance and Darcy law, together with an initial radius and compatible initial temperature fields, specify the moving-boundary model. No curvature correction is imposed. A vanishing initial seed corresponds to the leading growth law below with $a(0)=0$.

Write $\Delta =T_m-T_{\infty }$ and $S=L/(c_p\Delta )$, the latent-to-sensible heat ratio. A conductive interfacial balance gives $a\dot{a}\lesssim k_{+}\Delta /\mathcal{L}$. On length scale $a$, the ratio of thermal diffusion time to growth time is consequently bounded byIt is small for $S\gg 1$ with fixed porosity and comparable constituent densities and thermal properties. Very small porosity requires the sharper condition $C_{\pm }\Delta /(\varphi {\rho }_{s}L)\ll 1$; the bare heat-capacity ratio alone should not be used uniformly as $\varphi \to 0$. The advective-to-conductive ratio isEven if mechanical heating is retained, its size relative to latent-heat release is $p_i\beta /\mathcal{L}=(T_m-T_i)/T_m\le \Delta /T_m$. For fixed water properties this is also small as $S\to \infty$, consistently with the linearized Clausius-Clapeyron relation. Thus unsteadiness, thermal advection and mechanical heating are higher-order effects. The local leading temperature fields satisfy Laplace's equation. The approximation describes radii of order $a$; the distant transient diffusion tail supplies the far-field matching and need not be globally steady.

Regularity at zero and the interface/far-field values giveThe Stefan condition becomes $\mathcal{L}a\dot{a}=k_{+}(T_i-T_{\infty })$. Introduce the dimensionless pressure-feedback factorCombining this heat balance with pressure-dependent phase equilibrium yields the growth and interface values:This is pressure-limited spherical freezing in a porous medium. The interface values are constant during this leading solution, although $p(r,t)$ changes at a fixed exterior radius.

As permeability of a porous medium decreases, $\Lambda$ increases: the interface pressure rises toward the finite value $\Delta /\gamma$, while the interface temperature falls toward $T_{\infty }$. Hydraulic resistance suppresses freezing through pressure-induced freezing-point depression. In the limit $\Pi \to 0$, $a\dot{a}\sim \Pi \Delta /(\gamma \mu \beta )$; a zero-seed radius grows as $\sqrt{\Pi t}$ and tends to zero at fixed time. A finite seed instead hardly changes. At exactly zero permeability an incompressible rigid model cannot expel the expansion volume, so sustained growth stops; deformation or compressibility would require a different model.

Let $z$ measure distance inward from the surface. While $a\ll R$, curvature is negligible and this is a planar two-phase Stefan problem. Use a common reference mass density $\rho$, specific heat capacity $c_p$, thermal conductivity $k$ and thermal diffusivity $\kappa =k/(\rho c_p)$; retain the mass density differences only in the subsequent buoyancy estimates. Define $\Delta T=T_m-T_s>0$, $\Theta =\overline{T}-T_m\ge 0$, and $\vartheta =\Theta /\Delta T$. The paper's parameter $S=L/(c_p\Delta T)$ is the latent-to-sensible heat ratio, reciprocal to the more usual Stefan number.

For superheated two-phase Stefan similarity, solving $T_t=\kappa T_{zz}$ with $T(0)=T_s$, $T(a)=T_m$ and $T(\infty )=\overline{T}$ gives the Neumann solution of the Stefan problem. With $\eta =z/(2\sqrt{\kappa t})$,The solid expression holds for $0<z<a$ and the liquid expression for $z>a$. The interface gradients inserted in the Stefan condition determine the positive constant $\lambda$:The right-hand bracket decreases from infinity and crosses zero, while the equivalent left-hand side $S\sqrt{\pi }\lambda e^{{\lambda }^2}$ increases, proving uniqueness. The liquid-side flux subtracts from the roof cooling: it is heat brought to the interface by the superheated interior.

For large $S$, $\lambda \ll 1$. Keeping the superheat term while using $\mathrm{erf}\lambda \simeq 2\lambda /\sqrt{\pi }$ and $\mathrm{erfc}\lambda \simeq 1$ yields $S{\lambda }^2+(\vartheta /\sqrt{\pi })\lambda \simeq 1/2$. Hence the small-interface-parameter growth with liquid superheat approximation isThe solid temperature is approximately linear across its thin crust, while the liquid thermal boundary layer has thickness of order $\sqrt{\kappa t}$. For fixed superheat and $S\gg 1$, $A\simeq \sqrt{2/S}$; for $\vartheta \gg \sqrt{S}$, $A\simeq \sqrt{\pi }/\vartheta$. This retains the strong suppression of growth by a very hot interior.

To estimate overturn, take local gravitational acceleration $g$ at the roof. The liquid's mass density contrast is $\Delta {\rho }_l={\rho }_f\alpha \Theta$, while the solid contrast is $\Delta {\rho }_s={\rho }_s-{\rho }_f$. A liquid thermal layer of thickness $\delta$ has Rayleigh number $\Delta {\rho }_{l}g{\delta }^3/({\mu }_l\kappa )$. As $\delta \sim \sqrt{\kappa t}$, it reaches its critical value after a time scaling asIn viscous overturn of a dense crust, an unstably stratified viscous crust has a Rayleigh-Taylor instability growth time of order ${\mu }_s/(\Delta {\rho }_{s}ga)$. Setting this equal to its age, with $a=A\sqrt{\kappa t}$, givesThese are formation-to-overturn times. A specified critical liquid Rayleigh number ${\mathrm{Ra}}_c$ multiplies ${\tau }_l$ by ${\mathrm{Ra}}_c^{2/3}$; logarithmic amplification and geometrical constants also affect the crust estimate. They are omitted for the requested order-of-magnitude scaling.

The given numbers yield $S\simeq 2.54$, $\Delta {\rho }_s=500\mathrm{kg}{\mathrm{m}}^{-3}$ and $\Delta {\rho }_l=100\vartheta \mathrm{kg}{\mathrm{m}}^{-3}$. ThusLiquid thermal overturn is enormously faster than crustal overturn throughout the stated range. Even a critical-Rayleigh prefactor of order $10^2$ does not change that ordering. Since $S\simeq 2.54$ is not asymptotically large, low-superheat coefficients from the permitted approximation are only approximate. At the upper end $\alpha \Theta =1.25$, so the linear thermal-density approximation also ceases to be quantitatively credible; that endpoint is a formal scaling estimate under the supplied constant-property model. At a common age the instantaneous viscous instability-time ratio instead scales as $5A10^{-18}/\vartheta$; it should not be confused with the above formation-time ratio.

For a boundary-layer renewal model, the liquid layer reaches ${\delta }_l\sim \sqrt{\kappa {\tau }_l}$. Its cycle-averaged thermal flux is consequentlyIf a critical Rayleigh number is retained, insert a factor ${\mathrm{Ra}}_c^{-1/3}$. Averaging the diffusive $t^{-1/2}$ flux over a cycle only changes an order-one coefficient. The crust reachesso its average solid-volume flux per roof area is the crust-removal fluxThe corresponding mass flux is ${\rho }_s{J}_s$ and latent-heat flux is ${\rho }_{s}L{J}_s$. This independent-cycle estimate is a scale for a crust that reaches overturn thickness and is then renewed. Because the liquid renews much more rapidly, in a hot interior its average heat flux can limit actual crust growth: the consistent interfacial budget is $\rho L\dot{a}\simeq k\Delta T/a-F_l$. A stalled crust would not realize repeated solid cycles. In a superheated interior, detached solid may also remelt; permanent core deposition is most directly justified in the late regime requested next.

At late times $\Theta =0$, liquid thermal buoyancy vanishes and $A^2\kappa \simeq 2\kappa /S=2k\Delta T/(\rho L)$. ThusAssuming thin renewed crust, efficient sinking and constant asteroid radius, conservation of mass gives $4\pi b^2\dot{b}\simeq 4\pi R^2{J}_0$ under the permitted equal-density approximation. Therefore the asteroidal core growth by crustal sinking law isFor a self-gravitating approximately uniform sphere, $g\simeq 4\pi G\rho R/3$, so $J_0\propto R^{1/3}$ and a zero-seed core scales as $b\propto R^{7/9}(t-t_0{)}^{1/3}$. This law applies while the thin-crust renewal model remains valid and liquid remains available; it cannot continue beyond complete solidification. The specific numerical lifetime cannot be inferred without the absolute solid viscosity and conductivity.

Introduce specific heat capacity $c_p$, latent heat per unit mass $L_f$, and $\Delta T=T_m-T_s>0$. Write the prescribed convective heat flux as $F=K(\overline{T}-T_m)$, whereThere is a dimensional convention to specify: with $Q$ measured in $\mathrm{W}\mathrm{k}{\mathrm{g}}^{-1}$ and dimensionless $\lambda$, this printed expression makes $K$ a heat transfer coefficient only if $\mu$ is kinematic viscosity. If $\mu$ denotes dynamic viscosity, the numerator needs ${\rho }_0^{2}g\alpha Q$ instead. The question says only “viscosity”; use the consistent kinematic interpretation, or equivalently the corrected dynamic formula. The following balances apply to either convention once $K$ has units $\mathrm{W}{\mathrm{m}}^{-2}{\mathrm{K}}^{-1}$.

For a thin lava lake crust with a linear temperature profile, Fourier's law gives the outward conductive flux $k\Delta T/a$. The Stefan condition balances it against upward thermal convection and freezing, while the well-mixed liquid loses the convective flux and gains radiogenic heating. To leading order in $a/H$ the radiogenically heated lava-lake crust model isThe fixed surface temperature already incorporates radiative cooling, so no extra radiative loss is added to the interior balance. The linear conductive profile is interpreted as a quasi-steady approximation, so sensible-heat storage within the crust is not evolved. Thinness alone does not establish this approximation: rapid thermal adjustment in the crust is also needed, for example a small Stefan number during the initial conduction-driven growth. Replacing the liquid depth by $H-a$ and retaining the moving liquid-volume enthalpy would introduce thin-crust corrections: the requested constant-temperature branch refers to the leading model above. A transient thermal profile in a thick or newly nucleating crust would require the full Stefan problem, rather than this linear-profile closure.

Using the thermal diffusion time $H^2/\kappa$, defineThe three independent dimensionless parameters are the convective-conductive coefficient $B$, the heating parameter $R$, and the Stefan number $S$. For positive $B,R,S$ the evolution equations becomewhere prime denotes $d/d\tau$. Temperature has the exact solutionConsequently the unique constant-temperature initial condition isAt this temperature, the convective flux equals the column's radiogenic heat production. Starting with no crust, the constant-temperature lava-crust growth law is obtained by separating $h^{\prime }=S(1/h-R)$:The crust grows towards $h_{\ast }=1/R$, or $a_{\ast }=k\Delta T/({\rho }_{0}QH)$, where conduction carries all the radiogenic input. The logarithm makes reaching this equilibrium take infinite time. For an initial thickness $h_0>0$, subtract $[-R{h}_0-\log (1-R{h}_0)]/(S{R}^2)$ on the same side of the equilibrium; the corresponding antiderivative uses $\log |1-Rh|$ on either branch. The thin-crust description through equilibrium requires $R\gg 1$.

For general initial temperature, the full early and intermediate history is specified by the exact temperature above and one scalar ordinary differential equation. A particularly convenient integral formulation for an initially bare surface isThis removes the infinite initial derivative of $h$. The right-hand side points into $y>0$ at zero, and is decreasing in $y$ when ${\theta }_0\ge 0$; two nonnegative solutions therefore cannot separate, which gives uniqueness. Its solution, together with the temperature formula, determines the entire model evolution and can be computed without a phase-boundary search. The early expansion follows by substitution:The dimensional leading crust thickness is $a\sim \sqrt{2k\Delta Tt/({\rho }_0{L}_f)}$. Thus hotter liquid reduces the first correction to conductive crust growth, whereas colder liquid permits more rapid thickening.

If ${\theta }_0>{\theta }_{\ast }$, the temperature decreases and the instantaneous stationary lava-crust thickness $h_e(\tau )=1/[B\theta (\tau )]$ increases towards $h_{\ast }$. Starting from zero, the crust always remains below this rising curve: at a first contact $h^{\prime }=0$ but $h_e^{\prime }>0$, so it cannot cross from below. Hence hot-start crust grows monotonically. It also stays thinner than the constant-temperature solution: the early expansion puts it below, and at a proposed first crossing their derivative difference is $-SB({\theta }_0-{\theta }_{\ast })e^{-B\tau }<0$, preventing crossing from below. If $0\le {\theta }_0<{\theta }_{\ast }$, the temperature increases and the instantaneous stationary lava-crust thickness decreases. The crust stays thicker than the constant-temperature solution by scalar comparison, but it may cross the decreasing $h_e$, reach a maximum and subsequently melt back towards $h_{\ast }$. Such overshoot is possible, not inevitable; the two relevant relaxation rates are $B$ and $S{R}^2$. The numerical plot uses $B=4$, $R=20$, $S=0.1$, so crust adjustment is faster than temperature relaxation and the cool-start example does overshoot:

The curves use ${\theta }_0/{\theta }_{\ast }=0.5,1,1.5$ and initially $h=0$. They solve the complete reduced time-dependent equations, rather than replacing the moving stationary thickness by the actual crust. At late times, set $h=h_{\ast }+v$, $\delta ={\theta }_0-{\theta }_{\ast }$; linearization gives $v^{\prime }=-S{R}^{2}v-SB\delta e^{-B\tau }$. For $S{R}^2\ne B$,At equal rates the forced term is $-SB\delta \tau e^{-B\tau }$. In particular, slow warming with $B<S{R}^2$ approaches equilibrium from a thicker crust, explaining the plotted remelting. These long-time expressions describe the linearized response; nonlinear corrections can have their own faster exponential rates. All predictions are restricted to $a/H\ll 1$, nonnegative liquid superheat, and the assumed vigorous-convection flux law. If crust growth violates those conditions, the reduced model must be replaced rather than extrapolated to total freezing.

LetThe Neumann solution of the Stefan problem in the ice and substrate isandContinuity of heat flux at the contact gives, with $r=\sqrt{k_s/k}$,and therefore

At the ice–water interface, the water is isothermal at $T_m$. The Stefan condition $\rho L\dot{h}=k{T}_{i,z}(h^{-})$ givesEliminating $T_0$ yields the implicit equationwhich determines $\lambda$ and hence $h(t)=2\lambda \sqrt{\kappa t}$.

If $k_s\gg k$, then $r\gg 1$, $T_0\sim T_s$, andThe highly conducting substrate acts as a reservoir fixed near its initial cold temperature, giving the usual one-phase Stefan problem.

If $k_s\ll k$, then $r\ll 1$, $T_0\sim T_m$, and $\lambda \sim r/(\sqrt{\pi }S)$. ConsequentlyHere heat removal through the poorly conducting substrate is rate limiting, and only a small temperature drop is needed across the much more conducting ice.

For a one-phase spherical Stefan problem with melting-temperature boundary value $T_s-\Gamma /a$ and supercooling $\Delta T=T_s-T_{\infty }>0$, use the latent-to-sensible heat ratio $S=\mathcal{L}/(c_p\Delta T)$. Large $S$ separates heat diffusion from slow interface motion. The exterior harmonic temperature is $T=T_{\infty }+(\Delta T-\Gamma /a)a/r$. The Stefan condition gives the displayed radius equation with $a_c=\Gamma /\Delta T$. This is the large inverse-Stefan number convention; the common sensible-to-latent Stefan number is $1/S$.

For equal phase thermal properties, a roof undercooling $\Delta T$ and liquid superheat $\Theta$, the planar Stefan problem has interface $a=2\lambda \sqrt{\kappa t}$. The solid error function profile and liquid complementary-error-function profile give $S\lambda =e^{-{\lambda }^2}[1/\mathrm{erf}\lambda -(\Theta /\Delta T)/\mathrm{erfc}\lambda ]/\sqrt{\pi }$. Here $S=L/(c_p\Delta T)$ is the latent-to-sensible heat ratio. Liquid heat supplied to the interface reduces crust growth.