ON SIMILARITY SOLUTIONS FOR THAWING PROCESSES^{∗}

Ariel L. LOMBARDI  and  Domingo A. TARZIA

Abstract

We review some recent results for a mathematical model for thawing in a saturated semi-infinite porous medium when change of phase induces a density jump and a heat flux condition of the type -q_{0}t^{− \frac{1}{2}} is imposed on the fixed face x=0. Different cases depending on physical parameters are analysed and the explicit solution of the similarity type is obtained when a given condition for the thermal coefficient q_{0} is verified.

Key words. Stefan problem, free boundary problems, phase change process, similarity solution, density jump, thawing processes, freezing, solidification. Resumen. Se realiza una revisión de recientes resultados sobre un modelo matemático de descongelación en un medio poroso semi-infinito saturado cuando un cambio de fase induce un salto de densidad y cuando se impone una condición de flujo de calor del tipo -q_{0}t^{− \frac{1}{2}} en el borde fijo x=0. Se analizan diferentes casos que dependen de diversos parámetros físicos y se obtiene la solución explícita de tipo similaridad cuando una cierta condición sobre el coeficiente térmico q_{0} es satisfecha. Palabras Claves. Problema de Stefan, problemas de frontera libre, procesos de cambio de fase, solución de similaridad, salto de densidad, procesos de descongelación, congelación, solidificación. AMS subject clasification. 35R35, 80A22, 35C05.

1 Introduction

Phase-change problems appear frequently in industrial processes and other problems of technological interest [8]. A large bibliography on the subject was given in [11]. In this paper, we consider the problem of thawing of a partially frozen porous media, saturated with an incompressible liquid, with the aim of constructing similarity solutions.

We have in mind the following physical assumptions (see [2], [4], [5]):

1. A sharp interface between the frozen part and the unfrozen part of the domain exists (sharp, in the macroscopic sense).

2. The frozen part is at rest with respect to the porous skeleton, which will be considered to be indeformable.

3. Due to density jump between the liquid and solid phase, thawing can induce either desaturation or water movement in the unfrozen region. We will consider the latter situation assuming that liquid is continuously supplied to keep the medium saturated.

Although thawing has received less attention than freezing, our investigation is in the same spirit as [3], and [9], with the simplification due to the absence of ice lenses and frozen fringes.

We will study a one-dimensional model of the problem, using the following notation:

∗MAT - Serie A, 3 (2001), 7-12.

ε > 0  : porosity,
ρ > 0  : density; ρ_{w} and ρ_{I}: density of water and ice(g/cm^{3})
c > 0  : specific heat at constant density (\frac{cal}{g^{o}C})
k > 0  : conductivity (\frac{cal}{s cm^{o}C})
u  : temperature of unfrozen zone (^{o}C)
v  : temperature of frozen zone (^{o}C)
u = v = 0 : being the melting point at atmospheric pressure
λ > 0  : latent heat at u = 0 (cal/g)
γ  : coefficient in the Clausius-Clapeyron law (s^{2}cm^{o}C/g)
μ > 0  : viscosity of liquid (g/cm^{3})

and subscripts F, U, I and W refer to the frozen medium, unfrozen medium, pure ice and pure water, respectively, while S refers to the porous skeleton.

The unknowns of the problem are a function x = s(t), representing the free boundary separating Q_{1} = {(x, t) : 0 < x < s(t), t > 0} and Q_{2} = {(x, t) : s(t) < x, t > 0}, and the two functions u(x, t) and v(x, t) defined in Q_{1} and Q_{2}, respectively. Besides standard requirements, s(t), u(x, t) and v(x, t) fulfil the following conditions (we refer to [4] for a detailed explanation of the model):

u_{t} = a_{1}u_{xx} − bρṡ(t) u_{x} ,  in Q_{1}  (1)
v_{t} = a_{2}v_{xx} ,  in Q_{2}  (2)
u(s(t), t) = v(s(t), t) = dρs(t) ṡ(t) ,  t > 0  (3)
k_{F} v_{x}(s(t), t) − k_{U}u_{x}(s(t), t) = αṡ(t) + βρs(t) ṡ^{2}(t) ,  t > 0  (4)
v(x, 0) = v(+∞, t) = −A < 0 ,  x, t > 0  (5)
s(0) = 0  (6)
k_{U}u_{x}(0, t) = −\frac{q_{0}}{\sqrt{t}} ,  t > 0  (7)

with

a_{1} = α_{1}^{2} = \frac{k_{U}}{ρ_{U}c_{U}} ,  a_{2} = α_{2}^{2} = \frac{k_{F}}{ρ_{F} c_{F}} ,  b = \frac{ερ_{W} c_{W}}{ρ_{U}c_{U}} ,
d = \frac{εγμ}{K},  ρ = \frac{ρ_{W} − ρ_{I}}{ρ_{W}} ,  α = ερ_{I}λ,
β = \frac{ε^{2}ρ_{I} (c_{W} − c_{I}) γμ}{K} = εdρ_{I} (c_{W} − c_{I}) .

Problem I consists of equations (1)-(7), while by problem II we mean the system and (1)-(6) and (8) respectively, where

u(0, t) = B > 0,  t > 0 .  (8)

Problem II was previously studied in [5].

We will look for similarity solutions of Problem I in different cases according to the value of parameters ρ, β y d, following the methods introduced in [5],[10].

First of all, we note that the function u(x, t) = Φ(η), with η = \frac{x}{2α_{1}\sqrt{t}}, is a solution of (1) if and only if Φ satisfies the following equation

\frac{1}{2}Φ′′(η) + (η − \frac{bρ}{α_{1}} \sqrt{t}ṡ(t)) Φ′(η) = 0

and similarly, the function v(x, t) = Ψ(η) is solution of (2) if and only if Ψ satisfies the equation

\frac{1}{2}Ψ′′(η) + ηΨ′(η) = 0.

Therefore, we obtain the following result Theorem 1 The free boundary problem I has the similarity solutions

s(t) = 2ξα_{1}\sqrt{t}

u(x, t) = mξ^{2} + \frac{2q_{0}α_{1}}{K_{U}} g(p, ξ) − \frac{2q_{0}α_{1}}{K_{U}} ∫_{0}^{\frac{x}{2α_{1}\sqrt{t}}}exp(pyr − r^{2}) dr

v(x, t) = \frac{mξ^{2} + A erf (γ_{0}ξ)}{erfc (γ_{0}ξ)} − \frac{mξ^{2}}{erfc (γ_{0}ξ)} erf (\frac{x}{2α_{1}\sqrt{t}})  (9)

if and only if the coefficient ξ satisfies the equation

q_{0} exp ((p − 1)y^{2}) − K_{2}F (m, y) = δy + νy^{3},  y > 0  (10)

where

F (m, y) = (A + my^{2})\frac{exp(−γ_{0}^{2}y^{2})}{erfc(γ_{0}y)},  g(p, ξ) = ∫_{0}^{y} exp(pyr − r^{2}) dr  (11)

erf(x) = \frac{2}{\sqrt{π}} ∫_{0}^{x} exp(−r^{2}) dr , erfc(x) = 1 − erf(x),

and the constants K_{2}, m, δ, ν, γ_{0} are defined as follows

K_{2} = \frac{K_{F}}{α_{2}\sqrt{π}}  m = 2dρα_{1}^{2}  γ_{0} = \frac{α_{1}}{α_{2}} > 0
δ = αα_{1} > 0  ν = 2βρα_{1}^{3}  p = 2bρ .  (12)

In order to analize (10) we need some preliminary results.

Lemma 2 (a) If m > 0, then F grows from A to +∞, when y grows from 0 to +∞. If m < 0, then F has a unique positive maximum, from which it decreases to −∞. In both cases, F (m, y) ∼ \sqrt{π}γ_{0}my^{3} when y → +∞.

(b) For all p > 0, we have

(i) g(p, y) ≥ \frac{1}{py} (exp ((p − 1) y^{2})) − exp (−y^{2}) , y > 0

(ii) g_{y}(p, y) ≥ exp ((p − 1) y^{2}) + \frac{p}{2} (1 − exp (−y^{2})) > 0, y > 0

(iii) g(p, 0) = 0,  g_{y}(p, 0) = 1,  g(p, +∞) = +∞.

(iv) lim_{y→+∞} \frac{g(p,y)}{y^{2}} = 0 if p ≤ 0 and lim_{y→+∞} \frac{g(p,y)}{y^{2}} = +∞ if p > 0.

(v) lim_{y→+∞} \frac{y g(p,y)}{exp((p−1)y^{2})} = \{ \begin{matrix}\frac{1}{p−2} if p > 2 \\ +∞ if p ≤ 2 \end{matrix}

(c) If m > 0 and ν ≥ −m\sqrt{π}γ_{0}K_{2}, then the function

H(y) = F (m, y) + \frac{ν}{K_{2}} y^{3}

is strictly increasing.

Proof. The assertions (a) and (b) (i)-(iii) were proved in [5], and (b) (iv) and (v) easily follow from the identity (see [1])

\frac{2}{\sqrt{π}} g(p, y) = exp (\frac{p^{2}y^{2}}{4} ) ( erf (\frac{py}{2} ) + erf (\frac{(2 − p)y}{2} )).

Owing to function

Φ_{2}(y) = \frac{y}{\sqrt{π}} \begin{matrix}exp & (−y^{2}) \\ erfc & (y) \end{matrix} − y^{2}

is increasing (see [6],[7]) the assertion (c) holds.

In Section 2 we prove the existence and the uniqueness of the similarity solution for different values of the physical parameters ρ, β and d. In Section 3, we dicuss the equivalence of the problems I and II, and we extend some existence results for problem II obtained in [5].

2 Existence and Uniqueness of Similarity Solutions

In order to solve equation (10) we introduce the following function

Q_{0}(y) = \frac{δy + νy^{3} + K_{2}F (m, y)}{exp ((p − 1) y^{2})}  (13)

defined for y > 0 which verifies Q_{0}(0) = K_{2}A > 0.

Theorem 3 Let m be a positive real number. We define the following sets in the plane ν, p:

R_{1} = \{(ν, p) ∈ ℝ^{2} : −\sqrt{π}K_{2}γ_{0}m ≤ ν, p ≤ 1\} , R_{2} = ℝ^{2} − R_{1}.

We have:

(a) If (ν, p) ∈ R_{1} then the problem I has a unique similarity solution if and only if q_{0} > \frac{K_{F}}{α_{2}\sqrt{π}} A .

(b) If (ν, p) ∈ R_{2} then the problem I has a similarity solution if and only if 0 < q_{0} ≤ max_{y≥0}Q_{0}(y) .

Proof. To prove the existence (and uniqueness) of similarity solution to problem I, it is necessary and sufficient to verify that the equation (10) has a (unique) solution. The equation (10) has a solution ξ if and only if q_{0} = Q_{0}(ξ). The proof is splitted in four cases [7]:

(i) m > 0, ν ≥ 0 and p ≤ 1; (ii) m > 0, ν ≥ 0 and p > 1;
(iii) ν < −\sqrt{π}K_{2}γ_{0}m;  (iv)m > 0, −\sqrt{π}K_{2}γ_{0}m ≤ ν < 0 and p ≤ 1.

Remark 4 For m = 0, i.e. dρ = 0, there exist a unique solution of equation (10) if and only if the inequality q_{0} > K_{2}A is verified. This result has already been found in [10].

Remark 5 We note that in the case p > 1, if max_{y≥0} Q_{0}(y) > q_{0} > K_{2}A there exist at least two solutions. On the other hand, if q_{0} is sufficiently small, then there exists a unique solution. The situation is a bit different in the problem II, studied in [5], where it was proved the existence and uniqueness of similarity solutions in the case m > 0, ν ≥ 0, p ≤ 2.

Similarly, we can obtain the following results.

Theorem 6 Let m < 0. We define the sets

R_{3} = \{(ν, p) ∈ ℝ^{2} : ν > −\sqrt{π}K_{2}γ_{0}m, p ≤ 1\} , R_{4} = ℝ^{2} − R_{3}.

Then

(a) If (ν, p) ∈ R_{3}, there exists a solution when q_{0} > K_{2}A.

(b) If (ν, p) ∈ R_{4}, there exists a solution when 0 < q_{0} ≤ max_{y>0} Q_{0}(y) .

Proof. By using Proposition 2 we have, δy +νy^{3} +K_{2}F (m, y) ∼ ν +K_{2}\sqrt{π}γ_{0}my^{3} when y → +∞ and it follows that if (ν, p) ∈ R_{3} then lim_{y→+∞} Q_{0}(y) = +∞, from which we have [K_{2}A, +∞) ⊂ Range(Q_{0}). This proves part (a).

If (ν, p) ∈ R_{4}, it is easy to see that the function Q_{0} has a positive finite maximum, and then the second part is also proved.

Remark 7 In [7] it was studied the physical acceptability of the similarity solution as a function of the thermal coefficients.

3 Relationship between Problems I and II

Let (s, u, v) be given by (9), for some constant ξ > 0. Then u(0, t) is a constant given by

u(0, t) = mξ^{2} + \frac{2q_{0}}{K_{U}} g(p, ξ) > 0.  (14)

Then, we can consider the problem II, by imposing this new temperature as u(0, t) at the fixed face x = 0.

Theorem 8 Let m > 0, ν > 0, p ≤ 1 and q_{0} > K_{2}A. If (s, u, v) is the unique similarity solution of Problem I, then (s, u, v) is the unique similarity solution of Problem II, provided the constant B in the condition (8) is given by

B = mξ^{2} + \frac{2q_{0}α_{1}}{K_{U}} g(p, ξ)  (15)

where ξ is the unique solution of equation (10).

Proof. We know that (s, u, v) is given by (9) where ξ is the unique solution of (10), which can be written as

Q_{0}(y) = q_{0},  y > 0  (16)

By the results obtained in [5], there exists a unique solution to Problem II, with B defined by (15), given by

s̄(t) = 2ξ̄α_{1}\sqrt{t}

ū(x, t) = B − \frac{mξ̄− B}{g(p, ξ̄)} ∫_{0}^{\frac{x}{2α_{1}\sqrt{t}}}exp(pξ̄r − r^{2}) dr  (17)

v̄(x, t) = \frac{mξ̄^{2} erfc (\frac{x}{2α_{2}\sqrt{t}}) + A ( erf (γ_{0}ξ̄) − erf (\frac{x}{2α_{2}\sqrt{t}}))}{erfc (γ_{0}ξ̄)}

where ξ̄ is the unique solution of the equation

\frac{\sqrt{π}}{2} K_{1} (B − my^{2}) \frac{exp ((p − 1) y^{2})}{g(p, y)} − K_{2}F (m, y) = δy + νy^{3} ,  y > 0  (18)

and K_{1} = \frac{K_{U}}{α_{1}π} . It is easy to see that the solutions given by (9) and (17) are coincident if and only if ξ = ξ̄. Then, it is sufficient to see that ξ is a solution of (18) [7].

Suppose that (s, u, v) is a solution to problem I, with the boundary condition (7). By the results of Section 1, we know that (s, u, v) are given by (9), where ξ must satisfy the equation (10). For this solution, the temperature in the fixed boundary is constant and equal to B = u(0, t) = T_{0}(q_{0}, ξ), where T_{0} is the real function defined by

T_{0}(q, y) = my^{2} + \frac{2}{\sqrt{π}} \frac{q}{K_{1}} g(p, y),  q > 0, y > 0 .  (19)

Assuming that q > 0, we will describe some properties of function T_{0}. First of all, we note that T_{0}(q, 0) = 0. Besides, it follows from the proposition 2 that if m > 0 and p > 0, then T_{0}(q, y) is an increasing function in both of its arguments, with T_{0}(q, +∞) = +∞. If m < 0 and p > 0, then T_{0}(q, +∞) = +∞, and if m < 0 and p ≤ 0, then T_{0}(q, +∞) = −∞. Finally, if m > 0 and p < 0 then T_{0}(q, +∞) = +∞.

Suppose that m > 0. For each ξ̄ > 0 let

q̄_{0} = Q_{0}(ξ̄)  (20)

where Q_{0} is the function defined by (13). Let B̄ = T_{0}(q̄_{0}, ξ̄) = mξ̄^{2} + \frac{2}{\sqrt{π}} \frac{q̄_{0}}{K_{1}} g(p, ξ̄), then a solution to problem I with q_{0} = q̄_{0}, which is given by (9) with ξ = ξ̄ because of (20), corresponds to a solution to Problem II with B = B̄. Then, given B > 0, we can show the existence of solution to Problem II, by proving that B belongs to the image set of the function J(·) = T_{0} (Q_{0} (·) , ·). For different values of the physical parameters, we study function J and we obtain the following results [7].

Theorem 9 Let m > 0. If ν ≥ −\sqrt{π}K_{2}mγ_{0}, then there exists a similarity solution to Problem II. If, in addition, 0≤ p ≤ 1, then the similarity solution is unique. For the case ν < −\sqrt{π}K_{2}mγ_{0}, a sufficient condition in order to have the existence of a solution to Problem II, is that B verifies the inequality

B < m\frac{δ}{|ν|} + \frac{2}{\sqrt{π}} \frac{M}{K_{1}} g (p, \sqrt{\frac{δ}{|ν|}}) .

Remark 10 The last Theorem, extends a result of [5], where it was proved that if m > 0, ν < 0, then there exists a solution to Problem II when

B < \frac{mδ}{|ν|} .

Proposition 11 Suppose m < 0 and ν > \sqrt{π}K_{2} |m| γ_{0}, then:

(a) If p ≤ 2 the problem II has a similarity solution if

B > max (0, J(y_{0}))  (21)

(b) If p > 2 the problem II has a similarity solution if m + \frac{2(ν+\sqrt{pi}K_{2}mγ_{0})}{\sqrt{pi}K_{1}(p−2)} > 0 and (21) are verified.

Acknowledgements The first author has been supported by Universidad de Buenos Aires under grant TX048 and by ANPCyT under grant PICT 03-00000-00137 (Argentina). The second author has been supported by PIP # 4798/96 ”Free Boundary Problems for the Heat-Diffusion Equation” from CONICET-UA, Rosario (Argentina).

References

[1] M. Abramowitz, I. Stegun, ”Handbook of mathematical functions”, Dover, New York (1972).

[2] Ch. Charach, I. Rubinstein, ”Pressure-temperature effects in planar Stefan problems with density change”, J. Appl. Phys., 71 (1992), 1128-1137.

[3] Z. Ding, F. Talamucci, ”Primary frost heave: a quasi steady approach ”, Canad. Appl. Math. Quart., 4 (1996), 31-63.

[4] A. Fasano, Z. Guan, M. Primicerio, I. Rubinstein, ”Thawing in saturated porous media”, Meccanica, 28 (1993), 103-109.

[5] A. Fasano, M. Primicerio, D. A. Tarzia, ”Similarity solutions in a class of thawing processes”, Math. Models Meth. Appl. Sci., 9 (1999), 1-10.

[6] G. H. Hardy, J. E. Littlewood, G. Polya, ”Inequalities ”, Cambridge Univ. Press, London (1934).

[7] A. L. Lombardi, D. A. Tarzia, ”Similarity solutions for thawing processes with a heat flux condition at the fixed boundary”, Meccanica, to appear.

[8] V. J. Lunardini, ”Heat transfer with freezing and thawing”, Elsevier, Amsterdam (1991).

[9] Y. Nakano, ”Quasi-steady problems in freezing soils. I. Analysis on the steady growth of an ice layer”, Cold. Reg. Sci. Tech. 17 (1990), 207-226.

[10] D. A. Tarzia, ”An inequality for the coefficient σ of the free boundary s(t) = 2σ\sqrt{t} of the Neumann solution for the two-phase Stefan problem”, Quart. Appl. Math., 39 (1981-1982), 491-497.

[11] D. A. Tarzia, ”A bibliography on moving-free boundary problems for the heat-diffusion equation: The Stefan problem and related problems”, MAT, Serie A, Rosario, 2 (2000) (with 5869 titles on the subject). See www.austral.edu.ar/MAT-SerieA/2(2000)/

Ariel L. Lombardi, Departamento de Matemática, Facultad de Ciencias Exactas y Naturales, Universidad de Buenos Aires, 1428 Buenos Aires, Argentina. E-mail: aldoc7@dm.uba.ar

Domingo A. Tarzia, Departamento de Matemática - CONICET, FCE, Universidad Austral, Paraguay 1950, S2000FZF Rosario, Argentina. E-mail: Domingo.Tarzia@fce.austral.edu.ar ; tarzia@uaufce.edu.ar