CC BY
44
УДК 536.21
АСИМПТОТИКА РЕШЕНИЯ МОДЕЛЬНОЙ ЗАДАЧИ ТЕПЛОПРОВОДНОСТИ С ПРОСТЕЙШЕЙ ПОДВИЖНОЙ ГРАНИЦЕЙ И ВНЕШНЕЙ НЕЛИНЕЙНОСТЬЮ
ТИПА СТЕФАНА - БОЛЬЦМАНА
© 2009 г. В.Н. Козлов1, А.А. Аносов2, А.И. Акимов2, Г.М. Гузаиров2
1Санкт-Петербургский государственный политехнический университет, 195251, г. Санкт-Петербург, ул. Политехническая, 29, postbox@stu.neva.ru
' Оренбургский государственный педагогический университет, 460844, г. Оренбург, ул. Советская, 19, root@ospu.esoo.ru
1Saint-Petersburg State Politechnical University, 195251, St-Petersburg, Politechnicheskaya St., 29, postbox@stu.neva.ru
2Orenburg State Pedagogical University, 460844, Orenburg, Sovetskaya St., 19, root@ospu.esoo.ru
Построена асимптотика решения краевой задачи уравнения теплопроводности в точках, лежащих вне пограничного слоя равномерно движущейся границы, теплообмен на которой осуществляется согласно закону Стефана — Больцмана.
Ключевые слова: нелинейная краевая задача, равномерно движущаяся граница, асимптотическое разложение, интегральное представление решения, метод Лапласа, метод последовательных приближений для нелинейного интегрального уравнения.
The outer boundary layer asymptotics of solution of the boundary value problem for heat equation, when on its uniformly moving boundary proceeds heat transfer following the Stefan — Boltzmann law, was built.
Keywords: nonlinear boundary value problem, uniformly moving boundary, asymptotic expansion, integral representation of solution, Laplace method, iteration method for nonlinear integral equation.
Краевые задачи названного типа представляют определенный практический интерес, так как описывают процессы теплообмена в космосе, при отливке массивных тел, в термоэмиссионных преобразователях энергии и т.д. [1]. В настоящей работе исследуется решение T(x,t) следующей краевой задачи:
= Т(х,0) = f (х),
dt 8T(x,t)
ar/
дх
= Y-T\x,t),
(1) (2)
X=k\t
at
где (t,x)eDr = (J,x) | х > kxt, 0 < t <t0 ; Fo = —---
I
критерий Фурье; у, /г, (/г, > 0) - постоянные величины; f(x) - произвольная достаточно гладкая функция. Определим асимптотику функции T (x, t) при условии Fo—>0и вне «пограничного» слоя, т.е. в предположении, что x-kxt = Foc (0 < с < 1).
Задача (1), (2) при Г(х,0) =Т0 = const подробно рассмотрена в [2]. Следуя предложенной в ней методике (которая реализует геометро-оптический асимптотический метод [3]), и, в соответствии с идеей А.Н. Тихонова [4], запишем интегральное представление решения в виде
T (x, t) = J f (x') G( x, x', t-t ' dx'-
- Fo ) G(x, x', t - f) \x,^hf <p(f) df
(3)
где функция Грина II линейной крае-
вой задачи, соответствующей (1), (2), а второе слагаемое в правой части (3) - тепловой потенциал простого слоя с неизвестной плотностью . Функция
Сг(х,х',?-?') может быть записана явно [5]: 0(х, х',) =
1
2^/rFo (t-t')
(x-xf
( x+x'-2k1t' )2 k1 ( x'-k1t: )
e 4Fo (t-t')+e 4Fo(t-t')
Fo
k k1( x' -Щ')
1 e Fo erfc
r , Л x x 2k±t
2Fo
где erfcu = —¡=e "
У/7Г
2^Fo(t-t 1) £ (-1)"Г(« +1/2)
(4)
п-о Г(1/2) и
Обозначив интегралы в правой части (3) Д и /2 соответственно, найдем асимптотику первого из них. С учетом (4) имеем
1
1
,--¡f(x')e TFoTdx'-, -
о 2^Foi
°° (x+x)2 k1x'
jf(x' )e 4Fot+ Fo x 0
1- kj-L seD
Fo „t0( ) Г(1 / 2)
4Fo t ( x + x' )2
n+ 1/2
dx =
применим к его исследованию вышеуказанную схему и получим (здесь и далее ввиду громоздкости возникающих аналитических выражений расчетные формулы приводятся с точностью до О(Ро4)):
Ро кг(х-щ Г
Л2 = -^~е ю (2к^ - х) -
/ • \2к1(~х\_ (/Ш №_х)1Ро+
кг
f111 (2k1t — x) tfIV(2kjt - x)
k2
t2 fV (2kjt-x) 2
k\
Fo2 + 0(Fo3n.
(7)
Подставляя (6) и (7) в (5), окончательно имеем
з
h = IÖ„Fo"+0(FO4)
(8)
где
«о = /(*)>
1 ¿i^-AiQ
fif1=t/"(x) +—е бъ /*(2к^ — х), (9)
К
t2fIV{x) щ=—-— +
2 2
1 (*-■¥) ттт
+ —е--бъ tfni(2k,t-x)~
к\
f"(2klt-x) к,
t3/ww 1 _MfzM
«3 = ' Fo
6 kj
t2 fV (2k^-x)
_ t/1У (2к^~ х) /ш (2к^~ х) к + к2 Прежде чем найти асимптотику интеграла /2, необходимо доопределить неизвестную функцию (р (t). Подставив Т(x,t) из формулы (3) в краевое условие (2), получим для (р нелинейное интегральное уравнение
<p(t) = y
\f (x') G(x, x ', t-t ')
= А1+А2. (5)
Нетрудно проверить, что к интегралам Л1 и Л2 применим метод Лапласа [6]. Так, например, основной вклад в асимптотику первого из них дает точка х = х'. Тогда, разложив функцию /(х1) в ряд по степеням (х-х1) и распространив интегрирование на всю действительную ось (совершив при этом ошибку экспоненциально малого порядка), получим следующий асимптотический ряд:
оо /И г{1п) / \
А1 = X—-—Ро" . (6)
п-0 П\
Выяснив предварительно, что основной вклад в асимптотику интеграла . 12 дает точка х'= 2/г, 1-х .
-Fo jG(x, x',t-t')
t'=o dx
x=k}t 4
x,=kit,(p{f)df X=k\t
(10)
Обозначим первый из интегралов в правой части
(10) через Т*. Заметим, что его асимптотическое разложение может быть найдено из (8) непосредственной подстановкой х = кг1 . В частности, из (9) имеем
Т* = 1>*Рой + 0(Ро4),где
«о = /(VX «Г = tf"(k]t) +
/W)
kx
, tfu\klt) f"(kxt) d-у —--1
k1
k2
o
n-ij
2
o
X
n=0
2
а. _t3fn{kxt) | t2fv{kxt) tfIV(kxt) | /Ш(У) 3 6 2кг к2 к? ' Уравнение (10) с учетом (11) может быть преобра-
зовано к виду <?«) = rt' у
i
T 0
А
2Fo
erfc
k^t-t') 2-JFO (/-/')
k^(i-i') в 4Fo
^ttFOQ -i') 1*
. Fo-1- — 0V!
k\
(12)
Из (12), полагая <px = y, последовательно находим
\2~
<Pl
л ,Fo" 1-4—у+6
К
Fo-
'з= Y
Fo-1-4—у+ 22
к/
к,
/ \2 Fo -
0(Fo ),
+ 0(Fo ) и т.д.
<p(t)= Г\f4(k1t)-4f3(k1t) f \kii)
f//F(Vb ki
-if \kxt )--
ki
Fo+2f 2(kii)
11
r2f(klt)
ki2
-14 rfiv{kJtrikit)+£(M ki \ ki
(13)
+ 2f(klt)
ri2fIV(kit) , if111 (kit) f "(kiOl
--1------ •
2 ki ki2 J
3
/'(V)
Fo2 [■ + 0(Fo3).
Используя (13), перепишем интеграл /2 в следующем виде:
12 =
Fo
2
в 4Fo(i-i')
у/л Fo 0 л/ i-i '
i-
и+1
_ g(-l)"r(„ + l/2)4„Fo„_
„То Г(1 / 2) (х-£/)
(14)
(p (/') dt.
Это уравнение Вольтерра, ядро которого имеет слабую особенность и которое, как было показано в [4], решается классическим методом последовательных приближений. Воспользовавшись тем, что
T* от t не зависит, и введя новые переменные:
cp(t) = <p(t)/T* , у = у - (Т*)3 , получим, аналогично
тому, как это сделано в [2], рекуррентную формулу:
Нетрудно убедиться, что (14) - интеграл Лапласа, т.е. к нему также применим метод Лапласа. Подобно тому, как это делалось для 11, предварительно устанавливаем, что основной вклад в искомую асимпто-— 2кЛ / х
тику дает точка г= —1-, и, разлагая в ряды по
К
степеням с последующим интегрированием по всей действительной оси (допуская при этом ошибку экспоненциально малого порядка), получаем асимптотический ряд для /2:
Ро о>
12=~ге К, 1АРО ,
п-0
где, в частности,
(15)
b0 =7f4(2kit-x), Ьг =
4yf\2kxt-x)
К
i-x) +
Замечание 1. Можно показать, что, вычисляя я-с приближение, мы получаем п-й член разложения ср (I ).
Замечание 2. Как отмечалось в [4], из физического смысла задачи следует, что интеграл в правой части (10) должен принимать значения только из интервала (0;1). Это, в свою очередь, означает выполнение для п е N
неравенства 0 < —(рп < 1, налагающего естественные
К
ограничения на исходные данные задачи (1), (2).
Замечание 3. С сожалением приходится констатировать, что при проведении аналогичных расчетов в [2] была допущена досадная неточность в записи коэффициента при Fo2 : следует читать 22 вместо 10 (что, впрочем, никоим образом не сказалось на основном результате вышеупомянутой статьи).
Итак, принимая во внимание (11), находим решение уравнения (10):
+ 3 (х - кгГ) ф {Щ - х) 3 + к^ /(2к^ - х) /" (2к^ - х) ], ь2 = — (2 V - х) /1У(2 V - х) -
К
-24(х-кхг)/(22кхг-X) Г(2к^-х)1- (16)
- 24 (х - V)/2 (2 V " X) Г (2 V -х)/" (2 V -х) +
+ 6(х - V)2 С(2к^ - х) /'' (2к^ -х)+ С (2к^ - х) 3 ] + + 12 V О - к^)/2 (2V - х) Г (2 V - х) //я (2 V -х) + + 12х (х - V) Г(2к1/ - х) Г' (2к1/ - х) 3Г'' (2к1/ -х)-
-14укх1 /6(2^/ - х)/" (2^/ - х) + 11 / 2/10 (2 V " *) "
- 84 / (х - /5 (2^/ - х) (2^/ - х) 3 -Таким образом, подставляя (15) и (8) в (3), находим асимптотику решения задачи (1), (2):
(17)
п=0
где сп=ап- Ьп (при п = 0,1,2 коэффициенты сп непосредственно выписываются при помощи соотношений (9) и (16)).
Рассмотрим некоторые характерные особенности полученного температурного поля (17). Чтобы не затемнять суть дела длительными аналитическими выкладками, возьмем в качестве /(х) = Т0= const. Кроме того, для определенности ограничимся случаем / > 0 (процесс остывания тела). Разумеется, соответствующие коэффициенты легко получаются из приведенных выше фор-
2
4
мул, однако представляется целесообразным показать, как выглядит решение в данном конкретном случае. Предварительно заметив, что ф(1) - const, из (3) имеем
T(x,t) = Т01G(x,x',t — t')\f=Qdx'-
-Тоср (О} 0(х, , * - О \ М.
о
Воспользовавшись [6, теорема 1.3], можно доказать, что при Ро—>0 первый из интегралов в последнем равенстве асимптотически эквивалентен 1, а вто-
А'1 )
е №
рои - величине
ki
откуда окончательно
Fo h(x-ht) T(x,t)=T0-—e fo 9(f). k,
(18)
Отсюда следует, что в рассматриваемом примере для получения аналитического выражения для температурного поля с заданной степенью точности достаточно продолжить итерационный процесс (12), сходящийся довольно быстро (см. замечание 1). Так, например, с точностью до 0(Ро7) находим
к\(х-к1)
1-Ре Ро х (19)
T(x,t)=T0
(-Ap + llp2 -140/7
3 +969/ -7084/
+ 0(Fo ),
Значение Количество слагаемых
Fo/ к1 1 2 3 4 5 6
0,1 0 - 0,1 0,1 - 0,2 0,2 - 0,3 0,3 - 0,4
0,2 0 - 0,1 0,1 0,1 - 0,2
0,3 0 - 0,1 0,1 0,1
0,4 0 - 0,1 0,1
0,5 0 - 0,1
изменяется в диапазоне от 0,1 до 0,2 , необходимо воспользоваться формулой (19); если, наконец, < 0,2 - шести слагаемых уже не достаточно. Заметим, что аналогичная ситуация возникает и в том
случае, когда — > 0,5. Таким образом, потребность в
К
увеличении количества членов разложения в (18) возникает, очевидно, в том случае, когда смещение границы тела от начального положения весьма незначительно (т.е. границу практически можно считать неподвижной). Изменение температурных кривых с уменьшением кх при фиксированных Ро представлено на рис. 1. Их поведение хорошо согласуется с видом решения задачи (1), (2) для границы х = 0 [4].
В случае подвижной границы имеет место характерная особенность, также отмеченная в [4], - последовательные приближения с точностью до 0(Роп) при четных п аппроксимируют решение задачи снизу, а при нечетных - сверху, что дает возможность оценивать погрешность разностью между двумя последовательными приближениями (рис. 2).
з Ро
где р = уТ0 — и 0 < р < 1 (см. замечание 2). Строго
К
говоря, последнее неравенство следует рассматривать только лишь как своего рода необходимые условия, так как данная оценка хороша, если мы ограничиваемся в формуле (19) двумя слагаемыми. С ростом количества слагаемых неравенства ужесточаются (например, для трех слагаемых 0 < р < 0,25 и т.д.), в связи с чем интересно проследить динамику влияния числа членов разложения на точность результата. Так, если относительная погрешность не превышает 3 %, то соответствующие данные можно оформить в виде таблицы (параметры, входящие в таблицу, берутся с точностью до 0,1).
Таблица
Динамика влияния числа членов разложения на точность результата
Рис. 1. Зависимость температурного поля от ^
Пояснение. Числа в графах означают диапазон из-
з Ро
менения величины уТ0 . Например, если — = 0,2, то
К
трех слагаемых для получения результата (в пределах указанной точности) не достаточно; если уТ$ < 0,1,
достаточно четырех членов разложения; если уТ^
Рис. 2. Последовательные приближения решения @(Х,т) с различной степенью точности 0(Ро") (2 < п < 7)
Наконец, отметим еще один факт, в котором проявляется специфика решения задач с подвижной границей. Из формулы (18) видно, что л: и / входят в решение только в виде разности х - кг1. Значит, температурное
поле существенно зависит от единственного фактора: на каком расстоянии от границы находится точка. Физически это означает, что перепад температур между граничными и внутренними точками сохраняется в течение всего процесса, а геометрически - то, что график Т (х^) может быть получен путем параллельного переноса графика Т (x,tc) , где tc - некоторый фиксированный момент времени (рис. 3), что и было использовано при построении чертежей к настоящей статье.
о
t
1=0.1 т=0.2 т-O.J 1=1).-» т=«.5 1=0.6 т=0Л i=0-( 1=0.4 т=1.0
0 94 0 0.2 0.4 0.6 0.8 X 1
Рис. 3. Эволюция температурного поля с течением времени (точность до О^о7) )
Отметим, что формула (18) хорошо согласуется с результатами расчета температурного поля рассмотренной задачи посредством численных методов [7]. Это позволяет эффективно использовать соотношения (18) и (19), а в общем случае и (17) на практике. Так, на рис. 4 представлено температурное поле тонкого стержня, подвижный конец которого перемещается по закону х = 0,8? (Бо = 0,09).
0.97
Рис. 4. Распределение температуры в тонком полуограниченном стержне с равномерно движущимся концом
При построении графиков (рис. 1-3, =0,6) использованы безразмерные координаты (температура © = 7/70 , время г = ф(]) и пространственный масштаб (X = х/х0 ; х0 = ).
Исследование погранслойной асимптотики решения рассмотренной задачи приводит к изучению эталонного интеграла вида, отличного от описанного выше, к которому метод Лапласа уже не применим. Соответствующие результаты будут изложены в другой работе.
Авторы выражают глубокую признательность Г.А. Несененко за ряд ценных указаний и многочисленные обсуждения.
Литература
1. Беляев Н.М., Рядно А.А. Методы теории теплопроводности. М., 1982. Ч. I. 328 с.
2. Несененко Г.А. Пограничный слой в модельной задаче с подвижной границей при нелинейных краевых условиях типа Стефана - Больцмана // Методы и алгоритмы параметрического анализа линейных и нелинейных моделей переноса. М., 1986. С. 161-169.
3. Несененко Г.А. Геометро-оптический асимптотический метод решения сингулярно возмущенных задач нелинейного тепло- и массопереноса в областях с подвижными границами. М., 1995. 104 с.
4. Тихонов А.Н. Об остывании тел при лучеиспускании, следующем закону Stefan'a - ВоИгтапп'а // Изв. АН СССР. География и геофизика. 1937. № 3. С. 461-479.
5. Карташов Э.М. Аналитические методы в теории теплопроводности твердых тел. М., 1985. 480 с.
6. ФедорюкМ.В. Метод перевала. М., 1977. 368 с.
7. Ин А.Х. Программа расчета температурного поля в подвижных границах // Методы и алгоритмы параметрического анализа линейных и нелинейных моделей переноса. М., 1992. С. 109-114.
Поступила в редакцию_19 июня 2008 г.
