О МАТЕМАТИЧЕСКОЙ МОДЕЛИ ОДНОФАЗНОЙ ЗАДАЧИ СТЕФАНА
К. А. Калиева
Казахский национальный педагогический университет им. Абая, 050010, Алма-Ата, Республика Казахстан
УДК 517.946
Рассмотрена математическая модель двумерной однофазной задачи Стефана. Метод нахождения зависимости теплового поля и определение скорости движения границы фазового перехода основаны на использовании функции Грина для уравнения теплопроводности в однослойной среде с нерегулярной границей области. Функция Грина теплофизического процесса построена на основе интегральных преобразований Лапласа, Ханкеля и конечного sin-преобразования Фурье, что позволяет получить явное аналитическое представление зависимости теплового поля и определить область фазового перехода однофазной двумерной задачи Стефана.
Ключевые слова: однофазная задача Стефана, функция Грина в плоском двухгранном угле, теплофизические процессы в однослойной среде, скорость движения фазового перехода.
This paper is devoted to constriction of mathematical model for the one-phase Stefan problem with the irregular diving boundary of the region. We consider a two-dimensional heat equation with the irregular boundary of the division between phases. Firstly, we construct the Green function for heat equation with fixed boundary. Then, using the constructed Green function we obtain analytic solution for the one-phase Stefan problem. Uniqueness and regularity of the constructed solution with the irregular diving boundary have been proved in the weighted Sobolev spaces. This result gives some general properties and features for the analytic solution of the two-phase Stefan problem.
Key words: one-phase Stefan problem, Green’s function in a flat dithered angle, heat physical processes, velocity of diving boundary motion.
Введение. Интерес к задаче Стефана, являющейся одной из наиболее трудных краевых задач нестационарной теплопроводности, обусловлен сложностью ее математической постановки. К задачам Стефана приводят математические модели изменения энергии при кристаллизации или плавлении в случае наличия границы раздела фаз, при этом основная техническая трудность заключается в том, что характерные скорости конвективного движения в расплаве (измеряются в миллиметрах в секунду) и перемещения фронта кристаллизации (измеряются в миллиметрах в час) существенно различаются. Это затрудняет исследование нестационарной задачи. С точки зрения вычислений эта трудность выражается в том, что классические исследования сводятся к решению интегральных уравнений Вольтер-ры — Фредгольма, что не позволяет применять принцип сжатых отображений. В настоящей работе предлагается решение задачи Стефана с использованием функции Грина модельной задачи сопряжения для уравнения теплопроводности с разрывными коэффициентами [1], а также предпринята попытка применить метод функции Грина при построении решения задач со свободными границами. Метод функции Грина позволил преодолеть технические
трудности, возникающие при решении задач для уравнения теплопроводности с разрывными коэффициентами, выходящего на границу области в угловую точку. С помощью методики решения подобных задач можно использовать построенную функцию Грина для решения двумерной двухфазной задачи Стефана, в которой температура фазового перехода и скорость ее движения получены в явном виде.
Успехи, достигнутые в вопросах аналитической теории исследования задач Стефана, весьма незначительны, что обусловлено трудностями вычислительного характера, возникающими при нахождении точного аналитического решения. По-видимому, необходимо конструирование нового математического аппарата для задач типа задач Стефана.
Математические аспекты классической задачи Стефана представлены в работе [2], в которой доказана теорема существования классического решения задачи Стефана для параболического уравнения на малом промежутке времени. Решение задачи получается как предел при е ^ 0 решений вспомогательных “регуляризованных” задач. Для решений вспомогательных задач удается получить оценки, позволяющие говорить о компактности семейств решений в пространстве C(2,1).
Краевые задачи для уравнения теплопроводности в плоском угле D = D {r> 0, 0 <p< в} при 0 < в < 2п исследованы В. А. Солонниковым, Е. В. Фроловой, Е. В. Радкевичем, В. А. Кондратьевым и О. А. Олейник. Е. В. Радкевич рассмотрел начально-краевую задачу (без разрыва коэффициента) для уравнения теплопроводности в двухгранном угле, удовлетворяющего условию
0 < 1 + k
у <
V2e
1
1
jj(k+2,k/2+1) / s~\ \ /л
в весовом пространстве H_k_i (Qt), где в — плоский угол раствора.
В [1, 3] построена функция Грина в двухгранном угле для уравнения теплопроводности без разрыва коэффициента и получены коэрцитивные оценки в гельдеровских нормах. В. А. Солонниковым и Е. В. Фроловой с помощью преобразования Меллина построена функция Грина для уравнения Лапласа, полученные результаты использованы для доказательства разрешимости начально-краевых задач для уравнения теплопроводности в плоском угле. Задача для уравнения теплопроводности в однослойной среде сведена к разностному уравнению на комплексной плоскости. Получены априорные оценки решений в весовых пространствах Соболева Hj^_lk_k1/‘1 2+1'> (Qt) при выполнении условий 0 < 1 + k — у < П, где в — плоский угол раствора.
1. Основные обозначения и определения функционального весового пространства Соболева. Обозначим через Li1(Qt) банахово пространство, состоящее из всех измеримых на Qt функций с конечной нормой
u )
T
dt
r2^ \u(r, <y, t)\2 dx
1
2
Если l — целое число, то пространство ство с конечной нормой
H (г+1) n^_l
(П) определяется как гильбертово простран-
где
ИЦ+фп)
с\\ / к(р)г2(м 21+Ы) \Diu(r,P)\2dn j\=i+i П
1
2
1
j = (juh), \j \ = ji + j2,
5j
dxi1 dxi? ’
k(x)
Ki, x G П i, K2, x G ^2-
Пусть a — некоторая каноническая область, в которую отображается граница дП области П. Положим ат = a х [0,T]. Определим пространство Hij-+l--^2l//‘+1/^> (ат) следов по норме Ф(г, t):
|Ф|
н
(1+1/2,1/2+1/4) Я-1
т
(ат)
5] в(<e,s,j)r2(P-2(l+1)+2*+j.) |D,>D\1 Ф(гьт)|2daT ) +
2s+ji=l
ат
r
+ ( Idt Iв(V’sj)dax
\2s+j1=l 0 ax
(T T
^ [ dt f dT f в(p,s,j)rl(
2s+\1=i о о a
2(^-2(l+l/2)+2s+j1)
\r — r1\‘
DSDj1
Ф(г, t)
ф
-D!Dj■ Ф(г1Л
V/r1
Л 2(M-2(l+1)+2s+j1) \ЩЩ Ф(г1>^ — ЩЧ\ ^(rb t)
T T
+ | ^ f dt f dT f в(p,s,jy2^-
2s+n=i-iJo 0 a
\t — T\3/2
-2(l+1)+2s+j1) \DtDj\ Ф(г1 , t) — DSDr1 Ф(г1,Т )
da ) +
-da ) +
\t — T |5/2
2 \ 2 da)
2
2
1
e(^,S,j)
р(Ус(Ус\\, \j\< l + 1,
C\k(v), \j \ = l + 2-
Обозначим через нечную норму
(Qt) пространство измеримых функций u(x,t), имеющих ко-
u
н
(l+2,l/2+1) /л — l — 1
(QT )
E / e(T,s,j)r2(^-2(+1)+2s+a)\DtDi u{r,p,t}\2 dQT
2s+ \j\=l+2QT
1
2
При таком определении функциональных весовых пространств сохраняются основные свойства невесовых пространств, связанные со следами функции и их соответствующими вложениями.
2. Постановка задачи. Рассмотрим область DT = D х (0,T), где D = {r> 0, 0 < p < £(r, t)}, с границей dD = 7 U Г, где 7 = {r > 0,p = 0}; Г = {r > 0,p = £(r, t)} — граница области движения фазы. При этом 0 < p0 < п (p0 — тангенс угла наклона касательной к кривой £(r,t)). В области DT = D х (0,T) 7Т = 7 х [0,T] введем обозначение u(r,p,t) для точки M(r,p,t) G DT.
Требуется определить распределение температуры u(r,p,t) в среде с фазовыми переходами, удовлетворяющее уравнению
du
~8t
a2(
д 2u dr2
+
1 du r dr
+
1 d 2u r2 dp2
) + f (r,P,t)
в области Dt , начальному условию
u lt=o = uo(r, P)
в области D, граничным условиям на границе области
u IyT = P(r,t)
и условиям Стефана на фронте фазового перехода
(1)
(2)
(3)
u ^^(r.^+o u (r,t);
(4)
1 dui _ * dj (r,t)
Kr dP =Pq dt
' dp <p=$(r,t)+0 dt
(5)
Здесь j(r, t) — граница области движения, которую необходимо определить; q*, р — соответственно скрытая теплота и плотность кристаллизации или плавления, отнесенная к единице массы твердого тела.
Предположим, что начальная температура удовлетворяет условиям
u0(r,v) = O(r1 a), duo('’ Р) = O(r ф при r — 0,
uo(r, p) = o(rX a)
dr
duo(r, p) dr
o(r a) при r —- <X).
При i = 1, 2 граничные условия имеют вид
PjrJ! = 0{ri-,% 1 МД) = 0(r-n при r — 0,
r r dr
P(r,t) ( 1-UN 1 dP(r,t) ( -as
= o(r a),-------------= o(r a) при r —— <X),
r r dr
кроме того, выполняются условия согласования
(-1 )i+1^ra-2(m+1)+n+2s+1DjlDlmDst f (r,v,t)\t=o + a2mra-2(m+1)+jl+2s+1 A("Puo(r, p) = 0. (6)
3. Построение функции Грина для задачи сопряжения уравнения теплопроводности. Рассмотрим область D = {r > 0, 0 < p < po }, (0 < po < п), где po — тангенс
угла наклона касательной к кривой j(r,t); Dt = D x (0,T). В указанной области построим функцию Грина аналогично тому, как это сделано в работе [4]. Используя интегральные преобразования Лапласа, Ханкеля и конечного sin-преобразования Фурье, нетрудно получить функцию Грина в области Dt = {r > 0, 0 < p < po } x (0,T):
2 , 2 r +ri
g(r,p,ri,pl,t) = Е
Фи(<р)
ж е 4а2 £
-/
ГгЛ Ф(^)Фп(^ 1)
n=i a2t 2aft) | Фп |2 ’
sin \п>£ „2
--------, ||Фп|
К^0
sin AppO 11 2 sin2 Ап^0
Как показано в работе [2], для функции
Ап
2 \ ’ Х'п
n
2n +1 2p0
(7)
-п.
g(r,p,ri,pl,t) =
22 ж — r +/1 e 4а2 £
-1л
п=1
ГГi \ Ф(^)Фп(^ 1)
a2t Лп V 2a2^ II Фп 112
справедливы представления
g(r,p,ri,p ^ =
2 2 2 2 — r—2Г1 + A1 cos(^-^l) — r—2r1 + A1 cos(^+^l)
e 4а2£ 2а2£ ^ ^ ’ e 4а2£ 2а2£ VY^ ^ ’
4na2t
+
+
4na2t
r2+r2 rr, / 7+£« „
e—^cos m , { 1 7 Фт(р,^)Фm(p, p 1 )e0mP
2na2t
'Jn 1 2ni
p th(ppo)
dp I +
у—гж
r 2 + r 2 П
e -Ztr Г rr, cos 0 rr1 sin 9 | 1
+ 1------^T- I e 2а‘2£ --1—— Ja
Y+КЖ
Ф т(гр,р)Ф m(p,p 1 )e'
0p
2na2t
2a2t l 2ni
p cth(ppo)
-dp I d9-
r2 + r2 ж
у—гж
4+гж „ „
PZPi t ^ chz rrishzja IJ_ Г ф m(p,V)4 m(p, p rfp | dz. (8)
2na2t I 2a2t a 2ni I ~ л 1 1 [ J
p cth(ppo)
у—гж
ж (—i)r
Здесь Ja(x) = E
2m+a
x
функция Бесселя первого рода; Г — гамма-
г=0 т!Г(т + а + 1) \ 2 функция Эйлера.
Функция ФУ, p) под знаком контурного интеграла определена по формуле
ФУ,р)
sh(pp) sh(ppo)'
4. Построение температурного поля и границы фазового перехода для задачи Стефана. Полагая p0 = £(r, t), распределение температурного поля для задачи Стефана (1)-(5) будем искать в виде суммы тепловых потенциалов с неизвестной плотностью u(r, p, t) на неизвестной границе:
ж pri,t)
u(r,p,t) = J dri j Uo(ri,pl)g(r,p,ri,pl,t)dpl +
00
t ж pri,T)
+ dT dri f (ri,pl,T)g(r,p,ri,pl,t - T)dpl+
t СЮ
+ 2a2J drj u* (ri, т) 0 0
дд(г,<р,г1у<р l,t) 3y\
t (Ю
dri + 2a2dT P(ri, т)
- V1 =£(г1,т) J J
* ^ ' 00
t Ю 5(ri,r)
+ J dT J dri j u(ri,v1, т)g(r,p,ri,pl,t — т)dp 1 + 0 0 0
( t (
dg(r,v,ri,vl,t)
др 1
dri +
- ip\ =0
+ 2a2 dW Kxi [C(rijт)]g(r,p,rijpl,t)\p=i{riT)dri■ (9)
00
Здесь K\1 [•] = dpi + idr — pi; u(r,t) — плотность объемного потенциала, которую необходимо определить; £(r,t) — граница области движения фазы температурного поля. Непосредственно проверяется, что интегральные представления решения (9) задачи (1)-(5) в виде тепловых потенциалов удовлетворяют неоднородным начально-краевым условиям (2), (3) в соответствующих областях и первому условию сопряжения. С использованием второго условия сопряжения для определения границы раздела фаз получаем уравнение
d£(rjt) в2 (дЧ(r,t) + 19j(r,t) А|
dt I dr2 r dr r2 j
0,
(10)
К
где в =—.
pq*
Из условия соответствия интегрального представления (7) оператору теплопроводности
L [•] = d-a2
д2 1 д 1 д2
+ -ж- +
д^ ^ \St2 r дr r2 др2у
следует интегральное уравнение Вольтерры — Фредгольма второго рода
t (Ю S(ri,r)
u(r,p,t) + J d^y dri J riu(ri,pI ,т)Н (r,p,rijp\t — т)dp1 + F(r,t) = 0, 0 0 0
где
H (r,P,rijPl,t) = Y^ h{n)(rjrijt)$n(p ljrijt)$n(p,rjt)j
n=i
h(n) (r, ri,t) = d(n) (r, ri, t)g(n) + c(n) (r, ri, t)
дд(
n)
d(n)(r,rijt) = L0 [C(r,t)] + a2^
дr
2
(n) 2 дФп д£ д 2( д2 1 д
сС \r-jTijt) = 2a2^^ — , L0D = — — a2 1 — + - —
д£2 \St
2 (sp
дt a V дт2 ' r дr
те
F(r,p,t) = L {j dri j U0(ri,pI)g(r,p,ri,pI,t)dpI +
0 0
t те Z(rx,T)
+ / dT dri f (ri,p p,T)g(r,p,ri,pU — T)dp 1+
0 0 0 t те
,2
+ 2a2 dr P (ri}r)
00 t OO
dg(r,p,ri,pI,t)
dp 1
dri+
- (£>\ =0
+ 2a2 j dT j KXl [£(ri,r)] g(r,p,ri,pl,t)\vl=i{riT)dri+
~dg(r,p,ri,pI,t)
00
t OO
+2a2 J dr J u*(ri,T) 00
Рассмотрим ядро системы интегральных уравнений
dp \
dr-\\ .
- v\=pri,T)
Н(r,p,ri,pI,t) = Y^ h(%)(r,ript)^n(p
n=i
Используя формулы [4] те i
xY-ie-axI2v(2c^xdx) = Г(2~ , c-ie2a—-YM-Yv (—^ , Ref7 + v + 2 ) > 0,
r(2v + 1)
0
где (z) — функция Уиттекера, определяемая формулой
1
M^ v (z)
zv+2
22vB {v + p + i ,v — p + 2^
(1 + x)u-11-2 (1 — x)v+^-1 e ^ dx,
i
B [v + p + i ,v — p + 2) — в -функция, для ядер системы интегральных уравнений получаем
\Н(r,p,ri,pU)\ < Mo
__ (^-^1 Г
e 8a0t
4a0t 2
Введем обозначение
и + Н * и = F.
i
Резольвента интегрального уравнения (8) имеет оценку
\R\ < MQeboVl, (11)
где Mo,bo — определенные константы, зависящие от координат области движения фазы перехода.
Определяя резольвенту на основе оценки (11), найдем решение матричного интегрального уравнения Вольтерры — Фредгольма второго рода (8)
и = R * F — F
5. Формулировка основного результата. Доказана следующая
Теорема. Пусть Ai = ^, l > 0 — целое число, вещественные числа ц > 0, Ао = min {1, А^ удовлетворяют неравенству 1/2 < 1 + l — ц < А0, тогда при А1 > 1 и любых функциях ио(г,р) е Я^1(П), Р(r,t) е н{--2’1/2+3/4)(yt), fi (r,p,t) е ),
удовлетворяющих условиям согласования (6), можно утверждать, что
1) однофазная задача Стефана (1)-(5) имеет единственное решение u(r,p,t) е
Н,
(/+2,//2+1)
^—1—1
(DT), определяемое по формуле (9);
2) граница области фазового движения £(r,t) определяется по формуле
2 , 2 r +ri
£ (r,t)= ri£0(ri)
e 4,32t
в 2t
Ai
rr1
2рЧ
dr\
2 к
где Со(r) = r sin po; 0 < Po <n; в2 = —- ;
pq*
3) для решения u(r,p,t) справедлива оценка нормы
U || H(l + 2,l/2+i) ( ) < C
(Dt )-
Uo
1н(—-1р)
+ II f
H
(1,1/2)
-l-i
(DT )
+ || p || H(l+3/2,l/2+3/4)
l-i
(YT )
+
+ 1 KAi [£(r,t)]
\h (1 + 1/2,1/2 + 1/4)(Гт )
li
Здесь KAi [•] = +1 jL — Pl — оператор, из которого следует, что граница области фазового
движения £(r,t) определяется по формуле
£ (r,t)
J r1£o(r1)e
о
22 r +ri
4@21
/34 lAl
rr1
2p4
dr1,
2 к
£o(r) = r sin po; о < po <n; в2 = —-.
pq*
Заключение. Рассматриваемая задача о распространении тепла в средах с изменяющимся фазовым состоянием используется в различных задачах теплофизики (кристаллизация расплава; образование непрерывного слитка; растворение пузырька газа в жидкости; изучение динамики одномерного неизотермического испарения идеальной смеси; плавка цилиндрического стержня; прямые и обратные задачи теории фильтрации (такие, как задача нагнетания гидромеханического раствора в грунт; задача продвижения водонефтяного контакта в условиях упругого режима и т. п.); конвекция, вызываемая кристаллизацией переохлажденного расплава; обтекание твердого тела высокоскоростным потоком вязкой
несжимаемой жидкости; вращение жесткого цилиндра в вязкопластической среде; неустановившееся вязкопластическое течение между плоскопараллельными стенками и др.) Экспериментальный выбор охлаждения или нагрева теплового режима является дорогостоящим и не всегда реализуемым процессом, поэтому наиболее подходящим способом исследования задачи Стефана является нахождение аналитического решения, описывающего изменение температуры фазового перехода и границу области ее движения. Для каждой задачи тепломассообмена в соответствующих областях с подвижными границами требуется построение своей функции Грина с учетом начальных и граничных условий.
Список литературы
1. Solonnikov V. A., Frolova E. V. Lp-theory for the Stefan problem // J. Math. Sci. 2000. V. 99, iss. 1. P. 989-1006.
2. МЕйРМАНОВ А. М. О классическом решении многомерной задачи Стефана для квазилинейных параболических уравнений // Мат. сб. 1980. Т. 112, № 2. С. 170-192.
3. Солонников В. А. О разрешимости классических начально-краевых задач для уравнения теплопроводности в двухгранном угле // Краевые задачи математической физики и смежные вопросы теории функции. № 16. Зап. научн. сем. ЛОМИ. Л.: Наука. Ленингр. отд-ние, 1984. Т. 138. С. 3-165.
4. Калиева К. А. Задачи сопряжения для уравнения теплопроводности в плоском угле: Дис. ... канд. физ.-мат. наук. Алма-Ата, 2004.
5. Градштейн И. С. Таблицы интегралов, сумм, рядов и произведений. 4-е изд. / И. С. Град-штейн, И. М. Рыжик. М.: Физматгиз, 1963.
Калиева Куляш Абиловна — канд. физ.-мат. наук, доц. Казахского национального педагогического университета им. Абая;
e-mail: kklara 09@mail.ru
Дата поступления — 19.02.13