On the two-phase fractional Stefan problem

The classical Stefan problem is one of the most studied free boundary problems of evolution type. Recently, there has been interest in treating the corresponding free boundary problem with nonlocal diffusion. We start the paper by reviewing the main properties of the classical problem that are of interest for us. Then we introduce the fractional Stefan problem and develop the basic theory. After that we center our attention on selfsimilar solutions, their properties and consequences. We first discuss the results of the one-phase fractional Stefan problem which have recently been studied by the authors. Finally, we address the theory of the two-phase fractional Stefan problem which contains the main original contributions of this paper. Rigorous numerical studies support our results and claims.



There are no comments yet.


page 1

page 2

page 3

page 4


The one-phase fractional Stefan problem

We study the existence and properties of solutions and free boundaries o...

Performance Analysis of Fractional Learning Algorithms

Fractional learning algorithms are trending in signal processing and ada...

The Algebraic Theory of Fractional Jumps

In this paper we start by briefly surveying the theory of Fractional Jum...

BDF6 SAV schemes for time-fractional Allen-Cahn dissipative systems

Recently, the error analysis of BDFk (1⩽ k⩽5) SAV (scalar auxiliary vari...

On the rate of convergence of a numerical scheme for fractional conservation laws with noise

We consider a semi-discrete finite volume scheme for a degenerate fracti...

Nonlocal diffusion models with consistent local and fractional limits

For some spatially nonlocal diffusion models with a finite range of nonl...

Fractional Covers of Hypergraphs with Bounded Multi-Intersection

Fractional (hyper-)graph theory is concerned with the specific problems ...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1. Introduction

In this paper we will discuss the existence and properties of solutions for the well-known classical Stefan problem and the recently introduced fractional Stefan problem. A main feature of such problems is the existence of a moving free boundary, which has important physical meaning and centers many of the mathematical difficulties of such problems. For a general presentation of Free Boundary problems a classical reference is [16]. For the classical Stefan problem see Section 2 below.

The Stefan problems considered here can be encoded in the following general formulation


where the diffusion operator is chosen as follows:

There is a further choice consisting of considering both types of problems with one and two phases. More precisely, given a constant (latent heat) and (thermal conductivities), we take


for the one-phase problem, and


for the two-phase one. In general, is called the temperature, while the original variable is called the enthalpy. These denominations are made for convenience and have no bearing on the mathematical results.

Formulation (1.1) makes the Stefan problem formally belong to the class of nonlinear degenerate diffusion problems called Generalized Filtration Equations. This class includes the Porous Medium Equation ( with ) and the Fast Diffusion Equation (with ), cf. [24, 25]. Consequently, a part of the abstract theory can be done in common in classes of weak or very weak solutions, both for the standard Laplacian and for the fractional one. However, the strong degeneracy of in (1-Ph) and (2-Ph) (see Figure 1) in the form of a flat interval makes the solutions of (1.1) significantly different than the solutions of the Standard or Fractional Porous Medium Equation.

Figure 1. One-phase and two-phase Stefan nonlinearities

The first work on the fractional Stefan problem that we know of is due to Athanasopoulos and Caffarelli in [2] where it is proved that the temperature is a continuous function in a general setting that includes both the classical and the fractional cases. This is followed up by [15], where detailed properties of the selfsimilar solutions, propagation results for the enthalpy and the temperature, rigorous numerical studies as well as other interesting phenomena are established. Other nonlocal Stefan type models with degeneracies like (1-Ph) and (2-Ph) have also been studied. We mention the recent works [4, 7, 6, 8], where it is always assumed that is a zero order integro-differential operator. There are also some models involving fractional derivatives in time, see e.g. [26].

Organization of the paper. In Section 2 introduce the classical Stefan Problem from a physical point of view, mainly to fix ideas and notations and also to serve as comparison with results on the fractional case. We will also give classical references to the topic and discuss how to deduce the global formulation (1.1).

We address the basic theory of fractional filtration equations in the form of (1.1) in Section 3. We discuss first the existence, uniqueness and properties of bounded very weak solutions. Later we address the basic properties of a class of bounded selfsimilar solutions. Finally, we present the theory of finite difference numerical schemes.

We devote Section 4 to the one-phase fractional Stefan problem, which has been studied in great detail in our paper [15]. We describe there the main results obtained in that article.

Section 5 contains the main original contribution of this article, which regards the two-phase fractional Stefan problem. First, we establish useful comparison properties between the one-phase and the two-phase problems. Later, we move to the study of a selfsimilar solution of particular interest. Thus, in Theorem 5.4 and Theorem 5.6 we construct a solution of the two-phase problem which has a stationary free boundary, a phenomenon that cannot occur in the one-phase problem. Finally, we move to the study of more general bounded selfsimilar solutions. Theorem 5.7 establishes the existence of strictly positive interface points bounding the water region and the ice region . In particular, this shows the existence of a free boundary. Theorem 5.10 establishes the existence of a nonempty mushy region in the case . The last mentioned results are highly nontrivial and require original nonlocal techniques.

Rigorous numerical studies support also the existence of a mushy region in cases when . We comment on this fact in Section 6.

Section 7 is devoted to the study of some propagation properties of general solutions. We also support these results with numerical simulations that present interesting phenomena that were not present in the one-phase problem.

Finally, we close the paper with some comments and open problems.

2. The classical Stefan problem

The Classical Stefan Problem (CSP) is one of the most famous problems in present-day Applied Mathematics, and with no doubt the best-known free boundary problem of evolution type. The mathematical formulation is based on the standard idealization of heat transport in continuous media plus a careful analysis of heat transmission across the change of phase region, the typical example being the melting of ice in water. More generally, by a phase we mean a differentiated state of the substance under consideration, characterized by separate values of the relevant parameters. Actually, any number of phases can be present, but there are no main ideas to compensate for the extra complication, so we will always think about two phases, or even one for simplicity (plus the vacuum state).

The CSP is of interest for mathematicians, because it is a simple free boundary problem, easy to solve today when , but still quite basic problems are open for . It is always interesting for physicists, since there exist several processes of change of phase which can be reduced to the CSP. Finally, it is of interest for Engineers, since many applied problems can be formulated as CSPs, like the problem of continuous casting of steel, crystal growth, and others.

Though understanding change of phase has been and is still a basic concern, the mathematical problem combines PDEs in the phases plus a complicated geometrical movement of the interphase, and as a consequence, the rigorous theory took a long time to develop. A classical origin of the mathematical story are the papers by J. Stefan who around 1890 proposed the mathematical formulation of the later on called Stefan Problem also in dimension when modelling a freezing ground problem in polar regions, [23]. He was motivated by a previous work of Lamé and Clapeyron in 1831 [20] in a problem about solidification. The existence and uniqueness of a solution was published by Kamin as late as 1961 [19] using the concept of weak solution. Progress was then quick and the theory is now very well documented in papers, surveys, conference proceedings, and in a number of books like [22], [21], and the very recent monograph by S. C. Gupta [18].

2.1. The Classical Formulation.

A further assumption which we take for granted in the classical setting is that the transition region between the two phases reduces to an (infinitely thin) surface. It is called the Free Boundary and it is also to be determined as part of the study.

With this in mind let us write the basic equations. First, it is useful to have some notation. We assume that both phases occupy together a fixed spatial domain , and consider the problem in a time interval , for some finite or infinite . On the other hand, the regions occupied by each of the phases evolve with time, so the liquid (water in the standard application) will occupy and the solid (ice) at time . Clearly, for all The initial location of the two phases, and , is also known. Let us introduce some domains in space-time: , and let

The energy balance in the liquid takes the form of the usual linear heat equation


while for the solid we have


Here and are the resp. temperatures in the liquid and solid regions, while and are the resp. thermal conductivities, and the specific heats, and is the density. All of these parameters are usually supposed to be constant (just for the sake of mathematical simplicity). It is however quite natural to assume that they depend on the temperature but then we have to write and . In this paper, we will always consider .

Next we have to describe what happens at the surface separating and , i.e., the free boundary, . A first condition is the equality of temperatures,


This is also an idealization, other conditions have been proposed to describe more accurately the transition dynamics and are currently considered in the mathematical research.

We need a further condition to locate the free boundary separating the phases. In CSP this extra condition on the free boundary is a kinematic condition, describing the movement of the free boundary based on the energy balance taking place on it, in which we have to average the microscopic processes of change of phase. The relevant physical concepts are heat flux and latent heat. The result is as follows: if is the heat flux across , then


being the velocity with which the free boundary moves. The constant is called the latent heat of the phase transition. In the ice/water model it accounts for the work needed to break down the crystalline structure of the ice. Relation (FB2), called the Stefan condition, is not immediate. It is derived from the global physical formulation in the literature. Equivalently, if is the implicit equation for the free boundary in -variables, (FB2) can be written as .

All things considered, we have the complete problem as follows:

Problem about classical solutions. Given a smooth domain and a , we have to:

  1. Find a smooth surface separating two domains in space-time .

  2. Find a function that solves (E1) in and a function that solves (E2) in in a classical sense, Typically we require inside its domain , .

  3. On the free boundary conditions (FB1) and (FB2) hold.

  4. In order to obtain a well-posed problem we add in the standard way initial conditions

  5. Boundary conditions on the exterior boundary of the whole domain for the time interval under consideration. These conditions may be Dirichlet, Neumann, or other type.

The precise details and results can be found in the mentioned literature. Let us remark at this point that it is the Stefan condition (FB2) with that mainly characterizes the Stefan problem, and not the possibly different values of , on both phases.

2.2. The one-phase problem.

Special attention is paid to the simpler case where one of the phases, say the second, is kept at the critical temperature (e.g., in the water-ice example, the ice is at C). Then the classical problem simplifies to:

Finding a subset of bounded by a internal surface and a function such that

( is the fixed lateral boundary ), plus the Stefan condition

where is the normal speed on the advancing free boundary. The theory for the one-phase problem is much more developed, and essentially simpler.

2.3. The global formulation.

In order to get a global formulation we re-derive the model from the general energy balance plus constitutive relations. In an arbitrary volume of material we have

where is the energy contained in at time , is the outcoming energy flux through the boundary , and is the energy created (or spent) inside per unit of time. Therefore, represents an energy density per unit of mass (actually an enthalpy). We need to further describe these quantities by means of constitutive relations. One of them is Fourier’s law, according to which

where is the temperature and is the heat conductivity, in principle a positive constant. Thus, we get the global balance law

It is useful at this stage to include into the function by defining a new enthalpy per unit volume, . Equivalently, we may assume that . Using Gauss’ formula for the first integral in the second member we arrive at the equation

This is the differential form of the global energy balance, usually called enthalpy-temperature formulation. We have now two options:

  1. Either assuming the usual structural hypothesis on the relations between , , and , and performing a partial analysis in each phase, deriving the equations (E1), (E2), plus a free boundary analysis leading to the free boundary conditions (FB1), (FB2), or

  2. trying to continue at the global level, avoiding the splitting into cases. If we take the latter option which allows for a greater generality and conceptual simplicity.

We will take this latter option, which allows us to keep a greater generality and conceptual simplicity. We only need to add a structural relation linking and . This is given by the two statements:

  1. is an increasing function of in the intervals and .

  2. At we have a discontinuity. More precisely,

    jumps from 0 to at .

After some easy manipulations contained in the literature we get the relations stated in the Introduction and the integral formulation in Definition 3.2 with replaced by . If the space domain is bounded, we need boundary conditions on the fixed external boundary of . We see immediately that this is an implicit formulation where the free boundary does not appear in the definition of solution.

3. Common theory for nonlinear fractional problems

The theory of well-posedness and basic properties for fractional Stefan problems can be seen as a part of a more general class of problems that we call Generalized Fractional Filtration Equations (see [24] for the local counterpart). More precisely, one can consider the equation


for , , and


Together with (3.1) one needs to prescribe an initial condition .

Remark 3.1.

Throughout, we always assume , unless otherwise stated. For mathematical simplicity, we also assume in (1-Ph) and (2-Ph).

Our theory is developed in the context of bounded very weak (or distributional) solutions. More precisely:

Definition 3.2 (Very weak solution).

Assume (A). We say that is a very weak solution of (3.1) with initial condition if for all ,

Remark 3.3.

An equivalent alternative for (3.2) is in and

3.1. Well-posedness and basic properties

The following result ensures existence and uniqueness (see [17]).

Theorem 3.4.

Assume (A). Given the initial data , there exists a unique very weak solution of (3.1).

We will also need some extra properties of the solution. For that purpose, we rely on the argument present in Appendix A in [15] where bounded very weak solutions are obtained as a monotone limit of very weak solutions. The general theory for the latter comes from [14, 13]. See also [9, 10] for the theory in the context of weak energy solutions.

Theorem 3.5.

Assume (A). Let be the very weak solutions of (3.1) with respective initial data . Then

  1. (Comparison) If a.e. in , then a.e. in .

  2. (-stability) for a.e. .

  3. -contraction) If , then

  4. (Conservation of mass) If , then

  5. (-regularity) If as , then .

Additionally, the results of [2] ensure that for fractional Stefan problems, the temperature is a continuous function. We refer to Appendix A in [15] for an explanation of how the result [2] is applied to our concept of solutions.

Theorem 3.6 (Continuity of temperature).

Assume satisfy either (1-Ph) or (2-Ph). Let be the very weak solution of (3.1) with initial data . Then with a uniform modulus of continuity for . Additionally, if for some open set , then .

3.2. Bounded selfsimilar solutions

The family of equations encoded in (3.1) admits a class of selfsimilar solutions of the form

for any initial data satisfying for all and all . It is standard to check the following result, and we refer the reader to [15] for details.

Theorem 3.7.

Assume (A). Let be the very weak solution of (3.1) with initial data such that for all and all . Then is selfsimilar of the form

where the selfsimilar profile satisfies the stationary equation


When , we can choose a more specific initial data that will lead to a more specific selfsimilar solution from which we will be able to prove several properties for the general solution of (3.1). Indeed, we have the following Theorem which is new in the general context we are treating.

Theorem 3.8.

Under the assumptions of Theorem 3.7, and additionally, that , and that for some ,

Then the corresponding solution is selfsimilar as in Theorem 3.7. Moreover, it has the following properties:

  1. (Monotonicity) If then is nonincreasing while if , then is nondecreasing.

  2. (Boundedness and limits) in , and

  3. (Regularity) If satisfies either (1-Ph) or (2-Ph), then .


Part (a) follows by translation invariance and uniqueness of the equation (i.e. is the solution corresponding to for all ), since by comparison, if then . The bounds in (b) are a consequence of comparison and the fact that any constant is an stationary solution of (3.1). The limits in (b) are obtained by selfsimilarity and the fact that the initial condition is taken in the sense of Remark 3.3 (see Lemma 3.13 in [15] for more details). Finally, (c) follows from Theorem 3.6 and . ∎

Remark 3.9.

By translation invariance, one can obtain selfsimilar solutions not centred at by just considering for any . In this way, one obtains selfsimilar profiles of the form . Moreover, selfsimilar solutions in also provide a family of selfsimilar solutions in by extending the initial data constantly in the remaining directions. See Section 3.1 in [15] for details.

3.3. Numerical schemes

As in [15], we can have a theory of convergent explicit finite-difference schemes (see also [14]). More precisely, we discretize (3.1) by


where is the approximation of the enthalpy defined in the uniform in space and time grid for , i.e

On the other hand, is a monotone finite-difference discretization of (see e.g. [13]). It takes the form:


where are nonnegative weights chosen such that the following consistency assumption hold:


Together with (3.3) one needs to prescribe an initial condition. Since is merely we need to take

or just if has pointwise values everywhere in .

From [15] (see also [14]), we get the following convergence result.

Theorem 3.10.

Assume (A). Let be the very weak solution of (3.1) with as initial data such that for all , be such that , be such that (3.4) and (3.5) hold, and be the solution of (3.3). Then, for all compact sets , we have that

The above convergence is the discrete version of convergence in .

Remark 3.11.

We would like to mention that all the results of Section 3 apply also in the local case, i.e., replacing by in (3.1). More precisely,

  1. The existence part as in Theorem 3.4 is a classical matter (see [24]). We also refer to Appendix A in [15] for a modern reference in a more general local-nonlocal context.

  2. Properties as in Theorem 3.5 follow from the results in Appendix A in [15]. See also [12, 11].

  3. Regularity of as in Theorem 3.6 is the classical result of Caffarelli and Evans in [5].

  4. Convergence of numerical schemes as in Theorem 3.10 follows from the results of [14] replacing in (3.4) by the standard monotone finite-difference discretization of the Laplacian:

4. The one-phase fractional Stefan problem

Here we list a series of important results regarding the one-phase fractional Stefan problem


where is given by (1-Ph) with . Since such is locally Lipschitz, all the results listed in Section 3 apply for the very weak solution of (4.1). Moreover, we list (without proofs) a series of interesting results recently obtained in [15].

We start by stating the fine properties of the selfsimilar profile.

Theorem 4.1.

Assume that is given by (1-Ph) with , and let the assumptions of Theorem 3.8 hold with and for . The profile has the following additional properties:

  1. (Free boundary) There exists a unique finite such that . This means that the free boundary of the space-time solution at the level is given by the curve

    Moreover, depends only on and the ratio (but not on ).

  2. (Improved monotonicity) is strictly decreasing in .

  3. (Improved regularity) . Moreover, , for some , and (SSS) is satisfied in the classical sense in .

  4. (Behaviour near the free boundary) For close to and ,

  5. (Fine behaviour at ) For all , we have and for ,

  6. (Mass transfer) If , then

    If both integrals above are infinite.

Again, we remind the reader that selfsimilar solutions in also provide a family of selfsimilar solutions in by extending the initial data constantly in the remaining directions. Once the above properties are established in that case as well, one can prove that the temperature has the property of finite speed of propagation under very mild assumptions on the initial data.

Theorem 4.2 (Finite speed of propagation for the temperature).

Let be the very weak solution of (4.1) with as initial data and . If for some , , and , then:

  1. (Growth of the support) for some and all .

  2. (Maximal support) for all with

Moreover, the temperature not only propagates with finite speed, but it also preserves the positivity sets, an important qualitative aspect of the solution.

Theorem 4.3 (Conservation of positivity for the temperature ).

Let be the very weak solution of (4.1) with as initial data and . If in an open set for a given time , then

The same result holds for if is either or strictly positive in .

Finally, we have that the enthalpy

has infinite speed of propagation, with precise estimates on the tail. For simplicity, we state it only for positive solutions.

Theorem 4.4 (Infinite speed of propagation and tail behaviour for the enthalpy ).

Let be the very weak solution of (4.1) with as initial data.

  1. If in for and , then for all .

  2. If additionally for and , then

The question of asymptotic behaviour is still under study, but we refer to the preliminary results of the one-phase work [15].

5. The two-phase fractional Stefan problem

In this section we treat the two-phase fractional Stefan problem, i.e.,


where is given by the graph (2-Ph). Again, we make the choice .

5.1. Relations between one-phase and two-phase Stefan problems

Here we will see that any solution of the two-phase Stefan problem is essentially bounded from above and from below by solutions of the one-phase Stefan problem.

Proposition 5.1.

Let be the very weak solution of (5.1) with as initial data; the very weak solution of (4.1) with ; the very weak solution of (4.1) with ; and define .

Then in .

We need two lemmas to prove this result.

Lemma 5.2.

Let . Then is a very weak solution of (5.1) if and only if is a very weak solution of (4.1).


By comparison, implies that . Thus,

which concludes the proof. ∎

Lemma 5.3.

Let . Then is a very weak solution of (5.1) with initial data if and only if is a very weak solution of (4.1) with initial data .


By comparison, implies that . Moreover,

Now take . Then, in , and

That is,

Finally, the initial data relation follows from Remark 3.3. ∎

Proof of Proposition 5.1.

Lemmas 5.2 and 5.3 ensure that and are solutions of (5.1) with initial data and respectively. By the relation and comparison for problem (5.1), we have that . ∎

5.2. A selfsimilar solution with antisymmetric temperature data

We continue to address properties of the same type as in Theorem 4.1 for the two-phase Stefan problem. In the case where the initial temperature is an antisymmetric function, the solution has a unique interphase point between the water and the ice region, and it lies at for all times (stationary interphase). As a consequence, the enthalpy is continuous for all and discontinuous at for all times . This is the precise result:

Theorem 5.4.

Assume , , and

Let be selfsimilar solution (given by Theorem 3.8) of (5.1) with initial data . Let also and be the corresponding profiles. Then, additionally to the properties given in Theorem 3.8, we have that:

  1. (Antisymmetry) for all (hence, for ).

  2. (Interphase and discontinuity) is discontinuous at , where it has a jump of size . More precisely: if , if , and if and only if .

To prove this theorem, we need a simple lemma.

Lemma 5.5.

is a very weak solution weak solution of (5.1) with initial data if and only if is a very weak solution of

and initial data .


Clearly in and in . The initial data relation follows from Remark 3.3. ∎

Proof of Theorem 5.4.

1) Antisymmetry. We consider the translated problem as in Lemma 5.5 and prove that and are antisymmetric. We recall that the initial datum is

and . We avoid the superscript tilde on and in the rest of the proof for convenience.

To prove that is antisymmetric define . Note that

Then and in which ensures that

Note also that and thus is a very weak solution with initial data . By uniqueness, this implies that , which proves the antisymmetry result. The antisymmetry of follows. Note that the translation did not affect the .

2) Interphase points. We go back to the original notation without translation. Since is antisymmetric and also continuous we have . Moreover, since is nonincreasing, if and if . Define

We already know that , and moreover, (since and is continuous and nonincreasing). By antisymmetry of we also have that

3) Conclusion. Assume that , i.e., for all . Take any . Then, by antisymmetry of , we get

Note that for all , and