Stefan Problem Research Articles

The Stefan Problem With Internal Heat Generation in Spherical Coordinates

Abstract A weakly time-dependent equation for the evolution of the solid–liquid interface in spherical coordinates, driven by internal heat generation, is derived for constant surface temperature boundary conditions. The derivation comes by making an assumption that the interface moves slowly compared to the changes in the temperature so that the technique of separation of variables may be applied for Stefan numbers less than one. Under this approximation, we can separate the nonhomogeneous heat diffusion equation into transient and steady-state terms, and then integrate to get the temperature relations. With the temperature equations in hand, the derivatives are inserted into the interface equation giving a first-order differential equation for the location of the solid–liquid interface as a function of time. The results are compared to a previously derived quasi-static solution and a numerical simulation generated using the method of catching of the front. This method allows for direct tracking of a moving boundary via the calculation of the time it takes to move from node to node in a discretized grid characteristic of classical finite difference methods. Plots of the interface evolution show excellent agreement between the three methods, especially for lower Stefan numbers. The quality of the approximation decreases as the Stefan number increases, but the model is more accurate than the previously studied quasi-static model. For the Stefan numbers St = 1.0 and 10.0, the weakly time-dependent solutions are in better agreement with the numerical results than the quasi-static solutions.

The two-phase Stefan problem with anomalous diffusion

The non-local in space two-phase Stefan problem (a prototype in phase change problems) can be formulated via a singular nonlinear parabolic integro-differential equation which admits a unique weak solution. This formulation makes Stefan problem to be part of the General Filtration Problems; a class which includes the Porous Medium Equation. In this work, we prove that the weak solutions to both Stefan and Porous Media problems are continuous.

Computational assessment of immersed boundary–lattice Boltzmann method for complex moving boundary problems

Computational assessment of immersed boundary–lattice Boltzmann method for complex moving boundary problems

Mesh scheme for a phase transition problem with time-fractional derivative

Abstract The time-fractional phase transition problem, formulated in enthalpy form, is studied. This nonlinear problem with an unknown moving boundary includes, as an example, a mathematical model of one-phase Stefan problem with the latent heat accumulation memory. The posed problem is approximated by the backward Euler mesh scheme. The unique solvability of the mesh scheme is proved and a priori estimates for the solution are obtained. The properties of the mesh problem are studied, in particular, an estimate of movement rate for the mesh phase transition boundary is established. The proved estimate make it possible to localize the phase transition boundary and split the mesh scheme into the sum of a nonlinear problem of small algebraic dimension and a larger linear problem. This information can be used for further construction of efficient algorithms for implementing the mesh scheme. Several algorithms for implementing mesh scheme are briefly discussed.

Non-vanishing sharp-fronted travelling wave solutions of the Fisher-Kolmogorov model.

The Fisher-Kolmogorov-Petrovsky-Piskunov (KPP) model, and generalizations thereof, involves simple reaction-diffusion equations for biological invasion that assume individuals in the population undergo linear diffusion with diffusivity $D$, and logistic proliferation with rate $\lambda $. For the Fisher-KPP model, biologically relevant initial conditions lead to long-time travelling wave solutions that move with speed $c=2\sqrt {\lambda D}$. Despite these attractive features, there are several biological limitations of travelling wave solutions of the Fisher-KPP model. First, these travelling wave solutions do not predict a well-defined invasion front. Second, biologically relevant initial conditions lead to travelling waves that move with speed $c=2\sqrt {\lambda D}> 0$. This means that, for biologically relevant initial data, the Fisher-KPP model cannot be used to study invasion with $c \ne 2\sqrt {\lambda D}$, or retreating travelling waves with $c < 0$. Here, we reformulate the Fisher-KPP model as a moving boundary problem and show that this reformulated model alleviates the key limitations of the Fisher-KPP model. Travelling wave solutions of the moving boundary problem predict a well-defined front that can propagate with any wave speed, $-\infty < c < \infty $. Here, we establish these results using a combination of high-accuracy numerical simulations of the time-dependent partial differential equation, phase plane analysis and perturbation methods. All software required to replicate this work is available on GitHub.

Case studies on simulations of flow-induced vibrations of a cooled circular cylinder: Incompressible flow solver for moving mesh problem

The paper presents an in-house code developed for simulating incompressible flow involving moving boundaries. The solver is based on Arbitrarily Lagrangian-Eulerian formulation with Harten Lax and van Leer with Contact for artificial compressibility (ALE-HLLC-AC) method. In the finite volume formulation, the moving boundary problem is taken care of by the arbitrarily Lagrangian-Eulerian formulation, where radial basis function-based interpolation is employed to move the meshes. The Harten Lax and van Leer with -Contact scheme developed to evaluate the convective fluxes in artificial compressibility formulation. High order accuracy is achieved over an unstructured data structure by developing solution-dependent weighted least squares-based gradient calculations. Flow-induced vibrations caused due to the cyclic lift forces over a cooled circular cylinder are simulated and compared with available literature results. For flow-induced vibration of circular cylinder: fluid-structure interaction problem, the elastically mounted cylinder wall temperature is maintained lower than the flowing fluid. The cylinder vibrates transverse to the direction of fluid flow, and its oscillating movement is simulated using a numerical mass-spring-damper system. Both heat transfer scenarios, i.e., without gravity and gravity acting against the flow direction, are simulated. The unsteady, laminar, incompressible flow analysis is carried out at Richardson number (−1, 0) for the cooled circular cylinder with various reduced velocities Ur (3–20). Prandtl number and Reynolds number are Pr = 7.1 and Re = 150, with the fixed mass ratio of 2 and zero damping coefficient. The validation of the results is presented in the literature. It has been noticed that the thermal boundary layer affects the fluid characteristics. The vibration of the cylinder reached the highest amplitude of cylinder dimension, more dominantly in the range of reduced velocities (4–15), showing a galloping kind of flow nature.

A numerical study of a moving boundary problem with variable thermal conductivity and temperature-dependent moving PCM under periodic boundary condition

A numerical study of a moving boundary problem with variable thermal conductivity and temperature-dependent moving PCM under periodic boundary condition

A numerical study of moving boundary problem involving dual phase lag model of heat mass transfer during immersion frying

In this paper, we developed a dual phase lag model of heat mass transfer during immersion frying of foods. This model is a moving boundary problem of a couple system of second order hyperbolic partial differential equations. Modified finite element Legendre wavelet Galerkin method developed for its solution. The discretization in space and then application of Legendre wavelet Galerkin method converts our problem into a coupled system of generalized Sylvester equations. In a particular case, the solution of the present model is compared with exact solution and are approximately the same. L∞ error decreases as Legendre wavelet basis functions increase and strip size decreases. Whole analysis is done in dimensionless form. The effect of variability of moving interface on Fourier number, Luikov number, Posnov number, Kossovitch number, Biot number and relaxation times is discussed in detail.

Numerical approximation of the one-dimensional inverse Cauchy–Stefan problem using heat polynomials methods

Numerical approximation of the one-dimensional inverse Cauchy–Stefan problem using heat polynomials methods

A boundary integral method to investigate pattern formation in a rotating Hele-Shaw cell with time dependent gap

In this article we present a spectrally accurate boundary integral method to investigate the evolution of a Hele Shaw blob in a rotating, radial cell with time dependent gap. The blob of the viscous fluid is affected by the shearing forces induced by a changing gap distance, however the rotation also has significant influence on the shape of the blob. Using the method of complex moments, we characterize the solution of the moving boundary problem that models the evolution of the blob and investigate the morphology using linear analysis. The results of the linear analysis indicate the de-stabilizing influence of the rotation. In the numerical methods, we develop a spectrally accurate boundary integral technique to explore the complete nonlinear dynamics of the problem. Using the boundary integral simulations we demonstrate that, for a given strength of surface tension, the rotational effect causes a more complex pattern on the fluid interface that involves the formation of narrow fingers and more splitting of the interface. The effect is more pronounced if the angular speed increases. We also prescribe a method for controlling the viscous fingering using insights from linear analysis. Our numerical simulations validate the proposed controlling mechanism.

Establishment and analysis of numerical model of deuterium fuel ice migration in indirect-drive cryogenic target

A cryogenic target is the core component in inertial confinement fusion. For meeting the experimental requirements, the uniformity of the fuel ice layer in the capsule needs to reach 99% at least. The capsule has a nonuniform temperature distribution natively resulting from the columnar hohlraum geometry, leading to a nonuniform fuel ice layer profile. A proper temperature distribution along the hohlraum can be used to improve the uniformity of the fuel ice layer. A numerical model that can predict the kinetics of the fuel solid-vapor interface during the migration process is proposed by coupling the heat transfer and flow equations with the Stefan condition. This numerical model has good calculational accuracy verified with the pure-ice sublimation experiment. By studying different hohlraum temperatures and fuel ice layer initial eccentricities, the basic characteristic of ice immigration has been figured out. For the manufacture process, a temperature control scheme has been proposed by which the uniformity of the fuel ice layer can reach 99.57%. For the retention process, a perfectly uniform fuel ice layer can be maintained for 2000s before it cannot meet the ignition standard. These results are of significant guidance on cryogenic target experiments.

A new look in heat balance integral method to a two-dimensional Stefan problem with convection

In the current work, we developed a new approximation function for temperature profile with the help of Legendre wavelet in heat-balance integral method (HBIM) to solve a two-dimensional moving boundary problem with moving phase change material (PCM). It is assumed that PCM moves with induced velocity u along x and y direction. In heat transfer mechanism conduction and convection driven by fluid flow in liquid region is considered. To validate the current approximate method, we compared our numerical results with a previous work and found in strong acceptance. In particular, to show the accuracy of the present approximate method, we compared our numerical results against exact solution by converting present problem into a one-dimensional standard melting problem and found in good acceptance. The effect of Péclet number on temperature profile and moving melting front are analyzed in detail. Furthermore, it is shown that with a moving phase change material (PCM) the liquid/solid interface get accelerated and hence, the melting process becomes fast. This study may be applicable in thermal management and energy storage system.

Maximal regularity for evolution equations and application to the Stefan problem

Maximal regularity is a useful tool for solving abstract parabolic evolution equations. A variant of the one-phase quasistationary Stefan problem can be reduced to a single evolution equation. We approach this problem by maximal regularity and then apply an existence theorem for this type of evolution equation. We use as an assumption that one particular result on the solvability of a degenerate oblique derivative problem extends in an appropriate way.

A semi-analytical benchmark for the Stefan problem in arbitrary dimension – assessing accuracy of enthalpy-based methods

PurposeThis paper aims to propose a semi-analytical benchmarking framework for enthalpy-based methods used in problems involving phase change with latent heat. The benchmark is based on a class of semi-analytical solutions of spatially symmetric Stefan problems in an arbitrary spatial dimension. Via a public repository this study provides a finite element numerical code based on the FEniCS computational platform, which can be used to test and compare any method of choice with the (semi-)analytical solutions. As a particular demonstration, this paper uses the benchmark to test several standard temperature-based implementations of the enthalpy method and assesses their accuracy and stability with respect to the discretization parameters.Design/methodology/approachThe class of spatially symmetric semi-analytical self-similar solutions to the Stefan problem is found for an arbitrary spatial dimension, connecting some of the known results in a unified manner, while providing the solutions’ existence and uniqueness. For two chosen standard semi-implicit temperature-based enthalpy methods, the numerical error assessment of the implementations is carried out in the finite element formulation of the problem. This paper compares the numerical approximations to the semi-analytical solutions and analyzes the influence of discretization parameters, as well as their interdependence. This study also compares accuracy of these methods with other traditional approach based on time-explicit treatment of the effective heat capacity with and without iterative correction.FindingsThis study shows that the quantitative comparison between the semi-analytical and numerical solutions of the symmetric Stefan problems can serve as a robust tool for identifying the optimal values of discretization parameters, both in terms of accuracy and stability. Moreover, this study concludes that, from the performance point of view, both of the semi-implicit implementations studied are equivalent, for optimal choice of discretization parameters, they outperform the effective heat capacity method with iterative correction in terms of accuracy, but, by contrast, they lose stability for subcritical thickness of the mushy region.Practical implicationsThe proposed benchmark provides a versatile, accessible test bed for computational methods approximating multidimensional phase change problems. The supplemented numerical code can be directly used to test any method of choice against the semi-analytical solutions.Originality/valueWhile the solutions of the symmetric Stefan problems for individual spatial dimensions can be found scattered across the literature, the unifying perspective on their derivation presented here has, to the best of the authors’ knowledge, been missing. The unified formulation in a general dimension can be used for the systematic construction of well-posed, reliable and genuinely multidimensional benchmark experiments.

The similarity method and explicit solutions for the fractional space one-phase Stefan problems

In this paper we obtain self-similarity solutions for a one-phase one-dimensional fractional space Stefan problem in terms of the three parametric Mittag-Leffler function \(E_{\alpha ,m,l}(z)\). We consider Dirichlet and Neumann conditions at the fixed face, involving Caputo fractional space derivatives of order \(0<\alpha <1\). We recover the solution for the classical one-phase Stefan problem when the order of the Caputo derivatives approaches one.

Global solutions to the supercooled Stefan problem with blow-ups: regularity and uniqueness

We consider the supercooled Stefan problem, which captures the freezing of a supercooled liquid, in one space dimension. A probabilistic reformulation of the problem allows to define global solutions, even in the presence of blow-ups of the freezing rate. We provide a complete description of such solutions, by relating the temperature distribution in the liquid to the regularity of the ice growth process. The latter is shown to transition between (i) continuous differentiability, (ii) Holder continuity, and (iii) discontinuity. In particular, in the second regime we rediscover the square root behavior of the growth process pointed out by Stefan in his seminal paper [Ste89] from 1889 for the ordinary Stefan problem. In our second main theorem, we establish the uniqueness of the global solutions, a first result of this kind in the context of growth processes with singular self-excitation when blow-ups are present.

On a moving boundary problem associated with the swelling drug release platforms

This paper studies a parabolic-type system of time-dependent partial differential equations coupled with Robin boundary conditions to model the solvent penetration and subsequently a controlled drug release from a swellable polymeric spherical platform. The model's key feature is that the solvent penetration and the polymer swelling result in moving boundaries. From an analytical perspective, some properties of the solutions, such as the positivity and the uniqueness of the solutions, are established. A second-order iterative procedure based on the backward finite difference method is carried out to solve the problem numerically. Finally, the numerical results and the asymptotic solutions are compared, intending to illustrate the validity and applicability of the numerical scheme.

Pore-Scale Study on Convective Drying of Porous Media.

In this work, a numerical model for isothermal liquid-vapor phase change (evaporation) of the two-component air-water system is proposed based on the pseudopotential lattice Boltzmann method. Through the Chapman-Enskog multiscale analysis, we show that the model can correctly recover the macroscopic governing equations of the multicomponent multiphase system with a built-in binary diffusion mechanism. The model is verified based on the two-component Stefan problem where the measured binary diffusivity is consistent with theoretical analysis. The model is then applied to convective drying of a dual-porosity porous medium at the pore scale. The simulation captures a classical transition in the drying process of porous media, from the constant rate period (CRP, first phase) showing significant capillary pumping from large to small pores, to the falling rate period (FRP, second phase) with the liquid front receding in small pores. It is found that, in the CRP, the evaporation rate increases with the inflow Reynolds number (Re), while in the FRP, the evaporation curves almost collapse at different Res. The underlying mechanism is elucidated by introducing an effective Péclet number (Pe). It is shown that convection is dominant in the CRP and diffusion in the FRP, as evidenced by Pe > 1 and Pe < 1, respectively. We also find a log-law dependence of the average evaporation rate on the inflow Re in the CRP regime. The present work provides new insights into the drying physics of porous media and its direct modeling at the pore scale.

Asymptotic Analysis on the Sharp Interface Limit of the Time-Fractional Cahn--Hilliard Equation

In this paper, we aim to study the motions of interfaces and coarsening rates governed by the time-fractional Cahn--Hilliard equation (TFCHE). It is observed by many numerical experiments that the microstructure evolution described by the TFCHE displays quite different dynamical processes compared with the classical Cahn--Hilliard equation, in particular, regarding motions of interfaces and coarsening rates. By using the method of matched asymptotic expansions, we first derive the sharp interface limit models. Then we can theoretically analyze the motions of interfaces with respect to different timescales. For instance, for the TFCHE with the constant diffusion mobility, the sharp interface limit model is a fractional Stefan problem at the timescale $t=O(1)$. However, on the timescale $t=O(\varepsilon^{-\frac1\alpha})$, the sharp interface limit model is a fractional Mullins--Sekerka model. Similar asymptotic regime results are also obtained for the case with one-sided degenerated mobility. Moreover, the scaling invariant property of the sharp interface models suggests that the TFCHE with constant mobility preserves an $\alpha/3$ coarsening rate, and a crossover of the coarsening rates from $\frac{\alpha}{3}$ to $\frac\alpha4$ is obtained for the case with one-sided degenerated mobility, in good agreement with the numerical experiments.

Cartesian grid method with a consistent direct discretization approach for numerical simulation of non-Newtonian fluid flow

This paper proposes an extended version of the discretization scheme for a consistent direct discretization approach in a Cartesian grid method for flow simulations of incompressible non-Newtonian fluids. Based on the concept of the original method, which was developed for Newtonian fluid flow, the Navier–Stokes and the pressure Poisson equations are discretized directly in the vicinity of a solid boundary while maintaining the consistency between the velocity and pressure fields. The main difficulty in the discretization arises from the evaluation of the velocity gradients on the fluid–solid interface, which is necessary for handling the viscous term for spatially distributed shear-dependent viscosity. In the extended method, a new procedure for the calculation of the velocity gradient using the relation that holds on a rigidly moving solid surface is proposed. The validity of the method is ascertained in fundamental test problems of non-Newtonian fluids: laminar Poiseuille/Couette flows in a concentric annulus and flows around a circular cylinder. It is found that a second-order spatial accuracy with the desired numerical stability is achieved in both the shear-thinning and shear-thickening viscous properties, and the numerical error is considerably decreased compared with some conventional methods. As a realistic application of the extended method, a polymer mixing simulation inside a co-rotating twin-screw extruder is demonstrated by applying the method to the moving boundary problem. In a screw arrangement with a sequential combination of a kneading disk (KD), a cylindrical ring (CR), and a backward mixing screw (BMS), large- and small-scale structures of bifurcating and converging flows are respectively formed in the KD and BMS zones, whereas elongational flow structures and velocity fluctuations are still maintained in the CR zone due to flow interaction with the upstream and downstream screw elements, which presumably contributes to moderate mixing at low shear stress. • A Cartesian grid method for non-Newtonian fluid flow is developed. • The present method extends a consistent direct discretization approach. • A procedure for evaluating velocity gradients on an immersed boundary is proposed. • Second order spatial accuracies for velocity and its derivatives are ascertained. • A polymer mixing simulation in a co-rotating twin-screw extruder is demonstrated.