Crank J. Free and Moving Boundary Problems
Подождите немного. Документ загружается.
122 Analytical solutions
3.5.3. Embedding technique
A simple example
of
the
meaning of 'embedding' is provided by
the
case of a slab of ice melting under
the
action of
heat
falling on its surface
and with
the
water being instantaneously removed
on
formation.
The
ice
phase, which decreases in size as melting proceeds,
is
considered as
part
of a larger slab which retains
the
dimensions
of
the
original slab of ice
before meltIng started
and
therefore has fixed boundaries.
The
ice phase
is said to
be
"embedded"
in
the
larger slab of constant dimensions in a
mathematical sense.
Oearly,
only those parts
of
the
embedding slab
which coincide with
the
actual ice phase
at
any time have physical
meaning,
and
other
parts
are
fictitious.
In
particular,
the
outer
boundaries
which mark
the
extent of
the
initial slab of ice have fictitious boundary
conditions associated with them. Such conditions
are
not
prescribed
directly by
the
physical problem
but
must
be
constructed in a way such
that "the physical conditions prescribed
on
the
moving, melting interface
are satisfied. This
is
the
example that Boley (1961) used in order
to
introduce his embedding technique.
Consider a slab of thickness
f with surface heating conditions
aT
l
x=O,
O<t~t".,
(3.73)
K-+hT
l
= fl(t),
ax
aT
I
x=f,
O<t~tm'
(3.74)
KI-+hIT
I
=0,
ax
where K,
Kb
h,
hI
are
constants, TI(x, t)
is
the
temperature before
melting commences,
and
fl(t)
is a prescribed function.
The
conditions
(3.73) and (3.74) together with
the
heat flow equation
O<x<f,
(3.75)
where
k
is
heat difiusivity, and some prescribed initial condition, e.g.
O<x<f,
t=O,
(3.76)
suffice
to
define TI(x, t) as long as no melting occurs.
If
the
surface
reaches
the
melting temperature T m at time
tm
then
(3.77)
For
t
~
tm
a new temperature T
2
(x,
t)
is defined which satisfies
the
equation
s(t)<x<f,
(3.78)
Integral-equation formulations
123
and the following conditions
Tz{s(t), t}= T
m
,
aT
z
ds
K-+hTm=pL-
d
+ft{t),
x=s(t),
ax t
aT
z
K1-+h1Tz=0,
x=t,
ax
Tz{x,
t,J
=
T1(x,
t,J,
(3.79)
(3.80)
(3.81)
(3.82)
where
s(t) denotes
the
thickness of the melted (and removed) water i.e.
the
distance
the
melting interface has penetrated the original slab of ice, p
is
density, and L
the
latent heat of melting.
The
solution Tz{x, t), which has physical significance only in the region
s(t)<x<t,
is
now extended over
the
whole region
O<x<t
and this
extended solution
is
denoted by T(x,
t);
thus the solid phase
s(t)<x<t
is
'embedded' in the fictitious space
O<x<t.
Clearly, T(x, t) will satisfy
conditions (3.79-82) and equation (3.78)
will
be
satisfied in
O<x<t.
Condition (3.73) on x = 0, however, must
be
replaced by
aT
K-+
hT=
fl(t) + f(t),
ax
(3.83)
where
f(t)
is
an unknown function which, together with s(t),
is
to
be
determined so that
the
conditions (3.79) and (3.80) are satisfied. Boley
(1961) invokes the appropriate form of Duhamel's theorem (Carslaw and
Jaeger 1959, p. 31), which states:
If
v = F(x, t) represents
the
temperature at the point x and time t in
a solid in which
the
initial temperature
is
zero, while radiation takes
place into a medium at temperature unity, then the solution of
the
problem when radiation takes place into a medium at temperature
~(t)
is
given by
f'
a
v =
~
~(A)
at
F(x, t -
A)
dA.
(3.84)
In
the present problem
satisfies
the
equations:
we write
F(x, t) =
To(x,
t) which therefore
aTo
K
aZTo
-at=---aT'
o<x<t,
t>o,
To(x,
0)
= 0, (3.85)
aTo
x=O, t>O,
(3.86)
K-+hTo=h,
ax
aTo
K1-+h1To=0,
ax
x=t,
t>O.
(3.87)
124
Analytical solutions
The
condition (3.86) expresses radiation into a medium at temperature
unity
on
the
surface x =
O.
A standard solution for To(x, t), satisfying
(3.85-7),
is
easily deduced, for example, from
the
results of Carslaw and
Jaeger (1959, p. 118), and substitution of F by
To
in (3.84) gives
I
tm
aTo
T(x,
t) =
ft(t
1
)
at
(x, t - t
1
)
dtl
I
t
aT.
+
{fl
(t
1
)
+ f(tl)}
_0
(x, t - t
1
)
dtl
(3.88)
t
at
or, on rearranging,
(3.90)
Solution of these two simultaneous integro-diflerential equations, subject
to the initial condition
s(t.J
= 0, yields the two unknown functions
s(t)
and f(t). After this, the temperature
T(x,
t) follows from (3.88)
or
(3.89).
As
an
illustration, Boley (1961) considered
the
problem of a semi-
infinite solid, x
>0,
with a prescribed constant
heat
flux Q
o
at x =
O.
In
this case
h=O,
To = -
~
~(kt)ierfc
(2~(kt))'
(3.92)
and equations (3.90), (3.91) become
L
t-tm
_
dt
~7T
1 • ( S )
-
f(t
-
t)e
.2(t)/(4kt
1) - =
V'T'
--
2Q
(7Tt)'Ilerfc
--
I t! .n...l m k 0
2~(kt)
,
1 (3.93)
-S(t)
f(t-t
)e-·2(t)/(4kt
1
)_1
= 2(k7T)' Q
erf
--
-pL-
.
I
t
-
tm
dt
1{
G
S)
dS}
1
ti
0
~(kt)
dt
(3.94)
Integral-equation formulations 125
The
pre-melting solution,
Tb
in this example
is
_20
0
~.
(~)
T
1
(x, t) - K (kt)
lerlc\2J(kt)
.
(3.95)
In
terms of
the
non-dimensional parameters
1:.(
) = Oos(t)
s Y
1V"'T"
1T2.&.lm
F(Y)=-tg~
,
the
relation between
tm
and T
m'
which from (3.95) is given by
T.
=
20
0
(ktm)!
m K 1T '
becomes
(
pLm)2
tm=k
0
0
'
and
equations (3.93), (3.94) may
be
written
(3.96)
(3.97)
rYF(Y-Y1)e-f2(Y)/Yld~1=2[1-{1T(1+Y)}!ierlc(
~(y)
)],
(3.98)
J
o
Yi
Y'(1
+ y)
~(y)
r
y
F(Y-Y1)e-f2(Y)/Y1
d~1
=Y'1T{erl(
~(y)
)_~
d~(Y)},
(3.99)
J
o
Yi
Y'(1
+ y) m dy
for
y;;:.
0 and the initial condition becomes
~(O)=O.
(3.100)
Boley (1961) obtained series solutions of these equations valid for short
times
and
then proceeded
to
larger times by using numerical methods.
His short-time solutions are
2
1
m 4
(1
2 16m)
~
F(y) =
--yz+-y+-
2-
m
---~
Y"
1T
Y'1T
31T
31T
2
m
{f!6
35)
m 2
1}
2
+
2Y'1T
\3;+8
Y'1T
+m
-2
y +
...
,
(3.101)
2m
! m
2
y2
4m
(1
216m)
s
~(y)
=
31T
Y -
4Y'1T
-
151T
2-
m -
31T~
yz
m
2
(35m
2
1)
3
+
12Y'1T
-
8Y'1T
- m
+2
y +
...
,
, (3.102)
126
Analytical solutions
and the associated temperature
is
given by
~
=
"hr(l
+Y)}ierfc("(l:
y))+
x"y
ierfc(,,:)-~erfc(,,:)
O:;;;;y« 1, (3.103)
where
x = QoX/("7TKTnJ
and
(3.103) has physical meaning only for
x;;;;o~.
In
the same paper, Boley (1961) formulated
the
problem
of
the
solidification of a semi-infinite liquid under arbitrary cooling rates. Both
solid and liquid phases
are
present and so four simultaneous integro-
differential equations are needed to determine two unknown heat-flow
functions, together with
the
movement of the solidification boundary and
the temperature distribution. Series solutions are developed for
the
particular case of constant surface cooling.
It
is
found that
the
first three
terms of the series expressions for this two-phase problem
are
identical
with those of (3.101)
and
(3.102). This means
that
for short times
the
motion of the interface
is
virtually
the
same whether ablation occurs
or
not.
In
a later paper, Boley (1968) developed general short-time solutions
in one space dimension for melting
or
solidifying slabs which include
arbitrary heating conditions, arbitrary rates of liquid removal (including
the
limiting cases of stationary and instantaneously removed liquid),
temperature-dependent properties, and slabs of finite thickness.
Boley and Yagoda (1971) extended their embedding technique
to
obtain three-dimensional starting solutions for a melting slab
or
for
the
companion problem of
the
solidification
of
a liquid. A half-space
is
subjected to a general
heat
input
on
its surface, such that the change
of
phase starts at a point
of
the surface and then
sprel'l.ds
both towards
the
interior of the body
and
along its surface. Mathematically they considered
the half-space z
> 0, initially solid
at
a temperature lower than the
melting temperature T
m,
and heated
on
the
z = 0 surface by a prescribed
heat
input in
the
z-direction Q(x,
y,
t) which has a maximum
at
a
prescribed point of
the
surface, chosen as
the
origin. Thus they write
Q(X,
y,
t) =
a(t)-
b(t)X2-C(t)y2+
.
..
,
b,
c;;;;oO.
(3.104)
The
liquid phase
is
assumed to
be
removed instantaneously
on
formation
but
Boley (1961) and
Wu
and Boley (1966) showed
that
the
starting
solution
is
the same for arbitrary rates of ablation including zero and
infinite. This problem was solved earlier by Boley and Yagoda (1969) for
the
special case of a plane
heat
input (c = 0 in (3.104)) and by Yagoda
and
Boley (1970) for a radially symmetric
heat
input (b
==
c in (3.104)).
Boley and Yagoda (1971) showed
that
for very short times
the
shape of
the phase-change surface, when suitably normalized,
is
independent
of
Integral-equation formulations
127
the distribution and
the
rate
of
the
heat input, and of
the
thermal
properties of the slab.
The
magnitude of the melt depth does, of course,
depend
on
these parameters.
The
initial extents of melting in directions
normal to and along the surface are respectively proportional to
(t-
t,J~
and
(t-
t,J!, where
tm
is
the time
at
which melting starts
at
the origin. In
fact, they established that for short times the melt depth
(X,
Y,
T), where
all quantities are rendered suitably dimensionless,
is
given by
(X,
y,T)=A
2
,J{1-(x;T)Y-(Y;T)Yr.
T«
1, (3.105)
where
Xo(T),
Yo(T)
indicate
the
extent of melting along the surface in the
X, Y coordinate directions,
T
is
the
time after melting commences, and
Xo(T)=AT!,
Yo(T)=BT! subject
to
(XlXO)2=(Y/YO)2~1.
The
constants
A,
B, and
A2
are complicated expressions of
Cl>
Cz,;,
and
Cz
Y
' where
a
(T)
C--
-
1-
at
T '
m
all evaluated
at
X = Y = Z = T =
O.
Boley (1975) summarized investigations of the uniqueness of various
embedding solutions. There
is
a review by Boley (1978) of important
practical problems requiring further study in Wilson
et
al.
(1978).
The
subject of upper and lower bounds for solutions
is
included. Boley (1963)
in discussing an ablation problem proposed a comparison theorem which
allows upper and lower bounds
to
be
constructed by imposing
heat
fluxes
on
the
surface of
the
embedding body which are respectively always
larger
or
smaller than the desired flux.
An
analogous result for a melt that
remains in position was established by Boley (1964).
The
lack of compari-
son theorems and
of
solution bounds for problems with variable heat
parameters
is
noted by Boley (1978) coupled with tentative suggestions
for examining cases of limited variations of heat properties. Bounds and
approximate solutions for linear heat equations with non-linear boundary
conditions
are
established by Glasser and Kern (1978) who insert physi-
cally sensible inequalities into zero and first-order moments of the heat
equation.
Ferriss
and
Hill (1974) used an embedding technique
to
solve
the
oxygen diffusion problem set
out
in equations (1.57-60) of §1.3.10. They
introduced a fictitious condition
ac/ax
= I(t),
x=l,
(3.106)
which replaces conditions (1.59) of
the
original problem and solved the
linear problem defined by (3.106) together with (1.57), (1.58), and (1.60)
in
the
embedding region 0 ~ x ~ 1 by using a cosine integral transform.
Substitution of the two conditions (1.59) into the transform solution
128
Analytical solutions
yields a pair of simultaneous non-linear, integral equations for f(t) and
s(t). Ferriss and
Hill
employed a standard numerical method of solving
the
integral equations starting from a small-time solution
2
-p2n-
2
t
c(x,t)=2:
l:
~cosprrx+t(x2-1),
(3.107)
7T
p=H
P
which is simpler than the series of complementary
error
functions used by
Crank and Gupta
(1972a). They tabulated numerical values of c(x, t) and
of s(t) and made comparisons with a finite-difference method of their own
in the same paper and with the results of Crank
and
Gupta.
Comments by Ferriss and
Hill
(1974) on the embedding technique are
probably true in general. The amount of algebra involved in developing
the embedding formulation
is
by no means trivial and can
be
very
considerable. They found that the computing time required
by
a finite-
difference method was less than for the embedding solution, due mainly
to the tri-diagonal form of the finite-difference equations.
On
the other
hand, the embedding approach permitted larger time-steps for a pre-
scribed accuracy and, furthermore, a variable time-step could
be
conve-
niently introduced such that successive time-steps decreased in geometric
progression.
The
embedding solution also concentrates complete informa-
tion about c and its derivative into the two integral equations, once
f(t) is
known.
It
is worth noting
at
this point that Boley
and
collaborators have
only solved problems in which the solutions in the fixed embedding
domain can
be
written down immediately
and
analytically.
3.5.4. Heat-balance integral method: Goodman
By integrating
the
one-dimensional heat flow equation with respect to
the space variable
x,
and inserting the boundary conditions, Goodman
(1958) produced an integral equation which expresses
the
overall heat
balance of the system. A review article by Goodman (1964)
is
a useful
introduction to integral methods in heat flow problems generally, includ-
ing changes of phase. Earlier papers by Goodman from 1958 onwards
and by other authors are cited.
Successive steps in Goodman's method are:
(i)
Assume a particular form for the dependence of the temperature
on the space variable which
is
consistent with the boundary condi-
tions, e.g. assume a polynomial relationship.
(ii) Integrate
the
heat flow equation with respect to the space variable
over the appropriate interval and substitute
the
assumed tempera-
ture distribution to obtain the heat-balance integral.
(iii) Solve the integral equation to obtain
the
motion of
the
phase-
change boundary and then the time dependence of
the
tempera-
ture distribution.
Integral-equation formulations
129
The
method can conveniently
be
illustrated by solving the one-phase
melting-ice problem in
one
space dimension
as
defined in non-
dimensional terms by equations (1.21-26) with A
= 1. Integration of
(1.21) with respect to
x, from x = 0 to s(t), gives
d
1$(')
(dS
au )
-
udx=
-
-+-(0
t) .
dt dt
ax '
(3.108)
Assume a temperature distribution in the water phase given by
u =
a(t){x-
s(t)} + b(t){X-S(t)}2, (3.109)
which automatically satisfies (1.25).
It
is
convenient
to
modify condition
(1.26)
to
the form
by using the standard formula
du
= au ds + au = 0
dt
ax
dt
at '
x = s(t),
(3.110)
x = s(t)
to
replace ds/dt in (1.26). Remembering that u = 1, x = 0, (3.109) gives
1 =
-as
+
bs
2
and,
on
differentiation, au/ax = a, a
2
u/ax2 =
2b
at x = s(t), so
that
the
substitution into (3.110)
and
some manipulation yield
a =
(1-"/3)/s(t),
b =
(as
+ 1)/s2. (3.111)
Substitution of (3.109) incorporating a(t), b(t) from (3.111), and integrat-
ing
the
left side
of
(3.108) gives
J
ds ds
~(-1
+v3+2)
dt
= -
dt
-(a-2bs).
Further substitution gives
a-2bs=-(1-../3+2)/s
and so finally
ds
6(3-../3)
1 2
s
dt
= 7+../3
=2
a
,
(3.112)
say, from which
it
immediately follows that
s(t) = at1, (3.113)
since
s(O)
= 0 from (1.24).
The
temperature distribution follows readily by
substituting in (3.109) from (3.113) and (3.111). Instead of using the
heat-balance integral (3.108) after finding the temperature profile, Good-
man (1958) suggested use
of
(1.26)
as
an
alternative which leads to the
differential equation
s ds/dt =
(../3-1).
130
Analytical solutions
Goodman and Shea (1960) applied the heat-balance method to the
two-phase problem
of
the melting of a finite slab which
is
initially at a
unifonn temperature below the melting point. Initial sub-cooling was also
considered by Yuen (1980).
The
relevant equations are conveniently
written (Fox 1975)
aU1_ k a
2
U1
at - 1 ax
2
'
aU2_ k a
2
U2
at - 1 ax
2
'
O<x<s(t),
(3.114)
s(t)<x<e,
(3.115)
where the suffixes
1 and 2 denote the liquid
and
solid regions respec-
tively. Certain conditions are also specified
on
x = 0 and
e,
and
on
the
melting boundary we have
K1
aU
1
_
K2
a~=pL
ds }
ax ax
dt
' x = s(t).
U
1
=~=Um'
(3.116)
Goodman and Shea (1960) solved the pre-melting problem by the heat-
balance method though,
of
course, standard solutions are readily availa-
ble (Carslaw and Jaeger 1959) for their boundary conditions
on
x = 0
and
e.
Following the steps outlined in
(i-iii)
above, (3.114) is integrated from
x = 0 to s(t) and (3.115) from
s(t)
to e to obtain
the
integral equations
{
(
au
1
)
(au1)}
1
S
(t)
aU1
d l
S
(t)
ds
k1
- - - =
-dx=-
U1dx-um-,
ax
s(t)
ax 0 at
dt
dt
(3.117)
k2{(au
2
) -
(au
2
) }
=
11
aU2
dx
=~
11
U2
dx + u
ds
(3.118)
ax e ax
.(t)
.(t)
at
dt
s(t)
m
dt
.
Elimination of
aU1/aX
and
a~/ax
at x = s(t) from (3.116-118) now yields
(
K2 u _
K1
Urn
_
pL)
ds +
K1
d6
1
+ K2
d6
2
= K
(au
2
) _
K1
(au
1
)
k2
m
k1
dt
k1
dt
k2
dt
2
ax
e ax
0'
(3.119)
where
1
S
(t)
Ie
6
1
=
U1(X,
t) dx, 6
2
= u
2
(x, t) dx.
s(t)
(3.120)
The
two heat-balance equations (3.117) and (3.118) have been coupled
together through the melting condition (3.116)
to
yield (3.119).
Integral-equation formulations
131
Now suitable fonns for
Ut
and
~
are chosen, e.g.
Ut(x,
t)-
Urn = at(t)(x
-S)+b
t
(t)(X-S)2,}
U2(X,
t)-
U
rn
:=
a2(t)(x
-s)
+ b
2
(t)(X-S)2,
(3.121)
and
the
parameters at(t) etc. detennined by satisfying all the boundary
conditions and
the
melting condition. Proceeding as in the single-phase
problem discussed above
the
parameters at(t), a2(t), etc. are detennined
by satisfying
the
boundary conditions and preserving the heat balance.
Thus, for example, with boundary conditions
K
aUt=_H
tax
'
x=O;
a~=O
ax
'
Urn
= 0,
Goodman
and
Shea (1960) produced
the
differential equations
-pL
ds
+ K
t
d6
t
+ K
2
d6
2
=H
dt
k
t
dt
k2
dt '
(3.122)
(3.123)
(3.124)
(3.125)
Here, (3.123) comes inunediately from (3.119)
on
substituting from
(3.122); (3.124) is derived
by
first evaluating at(t), bt(t) from the condi-
tions (3.122)
and
(3.117) and finally using the first of (3.120). Equation
(3.125) is detennined similarly from the corresponding relations for the
solid phase.
The
initial conditions needed in order
to
solve (3.123-125) are
s(t~:=O,
(3.126)
The
last
of
(3.126) comes by integration with respect
to
x of the
appropriate pre-melting solution at
t:=
tm,
the time when melting com-
mences, which is given by Goodman
and
Shea (1960) as
-He
uz(x,
t~=
2K2
{l-(l-xle)2J.,
(3.127)
and
6z(t~
=
J~
~(x,
t)
dx
from
the
second
of
(3.120) and
the
first of
(3.126).
Goodman
and
Shea obtained solutions in series of equations (3.123-5)
together with (3.126) by a complicated analysis.
But
numerical solutions
are
easily obtainable by standard methods such as the Runge-Kutta
method.