CC BY
59
УДК 536.21; 27.35; 25
ВЗАИМОДЕЙСТВИЕ ИМПУЛЬСНОГО ИСТОЧНИКА ТЕПЛОВОЙ ЭНЕРГИИ С АНИЗОТРОПНЫМ ПРОСТРАНСТВОМ,
ТЕПЛОФИЗИЧЕСКИЕ ХАРАКТЕРИСТИКИ КОТОРОГО ЗАВИСЯТ ОТ ТЕМПЕРАТУРЫ
В.Ф. Формалев, Е.Л. Кузнецова, Л.Н. Рабинский
В статье впервые получено аналитическое решение задачи теплопереноса в анизотропном пространстве, компоненты тензора теплопроводности которого зависят от температуры. Теплоперенос инициирован импульсным точечным источником тепловой энергии. Оказалось, что теплоперенос носит волновой характер с конечной скоростью распространения фронта тепловой волны, хотя уравнение теплопроводности имеет параболический тип. Исследованы случаи различных степенных зависимостей компонентов тензора теплопроводности от температуры, показано, что фронты тепловых волн имеют вид эллипсоидов в пространстве и эллипсов на плоскости. Результаты могут быть использованы при исследовании теплообмена в композиционных материалах при их лазерном облучении.
Ключевые слова: импульс тепловой энергии, анизотропное пространство, нелинейные компоненты тензора теплопроводности, тепловая волна, геометрия фронта тепловой волны, скорость распространения тепловой волны.
При взаимодействии мощных излучений с поверхностями, ограничивающими анизотропные тела, теплофизические характеристики которых зависят от температуры, тепловые потоки и температурные поля могут распространяться с конечной скоростью, даже если теплоперенос описывается на основе гипотезы Фурье [1]. Этот анализ актуален при взаимодействии импульсных излучений с тепловой защитой гиперзвуковых летательных аппаратов, изготовленной из слоистых композиционных материалов, поскольку уровень температур в зоне взаимодействия может быть существенно занижен, а на периферии - завышен. Аналогичную проблему можно наблюдать при возникновении импульсного источника тепловой энергии внутри тела за счет различных физико-химических превращений.
Для изучения теплопереноса в анизотропных средах, теплофизиче-ские характеристики которых являются функциями в виде однородных многочленов температуры, в данной работе поставлена и аналитически решена задача о теплопереносе в нелинейном анизотропном пространстве при воздействии на него импульсного источника тепловой энергии.
Предположение о волновом характере теплопереноса с конечной скоростью распространения фронта теплового потока и температуры из-за нелинейности теплофизических характеристик полностью подтвердились, причем аналитический характер решения позволяет не только качественно, но и количественно исследовать волновой теплоперенос.
Аналогичные явления исследовались Самарским А. А. [1] с соавторами, Лыковым А.В. [2], Формалевым В.Ф. [3], Соболевым С.Л. [4] и др., однако аналитическое решение нелинейной задачи теплопроводности в анизотропном пространстве получено впервые.
Постановка задачи
Рассматривается задача о нестационарном распределении температур T (х, y, z, t) в трехмерном анизотропном пространстве V от мгновенного
точечного источника с энергией E0, приложенного в начале координат х = 0, y = 0, z = 0 в начальный момент времени t = 0, т.е. следующая задача Коши:
дт д (п ,^дтЛ дт
( дТ Л
cp— = — \lx (Т)— I+— l (Т)— + —14 (т) — I + 2— Я„ (т) — И дг дх { ххК дх I ду { ^ }ду I дz { zzK !дг) дх I ^ }ду
д ( Л , ^ч дт Л , д
( дт Л
+
+2 i 1 f]+2 |1 , y,z}e (-¥;¥), ' >0, (1)
E (х, y, z,0) = E0d( x - 0)d( y - 0)d( z - 0), (2)
где d(x- 0), d(y-0), d(z-0) - дельта-функции Дирака, E0 - импульсная
энергия, причем интеграл по пространственным переменным от распределения температур, инициированного этой энергией, есть величина постоянная
cpjjj T (х, y, z, t) dxdydz = E0 = const, (3 )
V
откуда начальное условие можно определить следующим образом:
jjj T (х, y, z, t) dxdydz = E = const. (4)
Таким образом, вместо начального условия (2) рассматривается условие (4), поэтому решение задачи (1), (4) будет зависеть от E0 / cp.
Компоненты тензора теплопроводности в уравнении (1) определяются соотношениями [5]
1 (T) = 1 (T )ax +1 (T )а +1 (T )aZx, 1 (T) = 1 (T ) = 1 (T +1 (T +1 (T )a(xa(y, 1 (T ) = 1 (T )=1 (T +1 (T )avxavz (T )aH, (5)
1 (T) = 1 (T H +1 (T )< +1 (T )aZy,
1 (Т)=1у (Т)=1 (ТК"* +1 (ТКК + 1 (Т)КК,
1 (т )=1 (т к+1 (т к+1 (т к,
где главные компоненты тензора теплопроводности зависят от температуры следующим образом:
1= кхТа, 1= кцТа, Л( = (6)
а к, I = {Х,ц,С}, j ={х, у, г} - направляющие косинусы углов между глав-
ными осями OX, Oh, OZ тензора теплопроводности и осями Ox, Oy, Oz прямоугольной декартовой системы координат; kX = const, kh = const, к^ = const.
Необходимо найти нестационарное распределение температур Т (х, y, z, t) под действием точечного источника (4), приложенного в точке с
координатами x = 0, y = 0, z = 0 в начальный момент t = 0, причем в остальных точках пространства в соответствии с условием (2) температура равна нулю.
Метод решения
Используем линейное преобразование поворота вокруг начала координат декартовой прямоугольной системы координат до совпадения с главными осями тензора теплопроводности, определяемое соотношениями
Х = а^хх + а^уу + , Л = апхх + апуу + апг2,
С = а?хх + а^у + (7)
Поскольку матрица линейного преобразования (7) невырожденна и ортогональна, то обратная матрица совпадает с транспонированной, вследствие чего из соотношений (7) получаем обратное преобразование
х = аххХ + а11хп + а?хС, У = аху% + а11у1 + а?у£, г = ахЛ + апЛ + а;Л.
Подстановка преобразования (7) в задачу (1), (4) приводит к задаче для уравнения, не содержащего смешанных производных
дТ д cp— = —
И dt дХ
f к? °дТ 1'
X дХ) дЛ
к Ts — khT dh
Л
+
)
дС
г
дТ 1 кТ ° —
С дС
[X,h,C}^(-¥-;¥), t > 0;
j]J Т (X, h, С, t) dXdhdZ = = const.
cP
Перейдем к новой системе координат
с/ \1/2 , ч 1/2 \ 1/2
X] =X(L / kx) , Х2 =h(L / kh) , Хз =Z(L / k() .
(8) (9)
(10)
где L - любое (например, L = 1), получим из (8), (9)
дТ д — = a— Эt дх1
Ts дТ . дх1.
+ a-
д_
Эх,,
Тс
дТ_ Эх
+ a-
2 )
Эх,
Г
дТ_ дх
3
[хг,хз]е (-¥;¥), t > 0; \l XkhkС ссс
L3/2 JJJ Т ( X1, х2 , X3, t) ^х2^3
El cp
(11) (12)
Здесь a = L / cp.
Будем искать решение задачи (11), (12) в автомодельном виде
T (^ ^ ^t ) = г<х'в( Ч» Ч^ Чз), (13)
где
х.
х,
_ х1 _ х2 _ х3
Ч1 = , Я = , Я3 = .
(14)
В выражениях (13), (14) показатели степеней а и р определим подстановкой (13), (14) в задачу (11), (12), получим
ав( д^ Я3 )-р
(
дв + Я 2
дЧх
дв
дв
Л'
+ я3
дЧ2 дЧ3 у
= &+1)-2р а
д
дч1 дч у дч.
+-
д
в°
дв
Л
д
(
+
дЧ2 У дЧ3
дв
\
г
,+^ЩК ^в( я, Я2, Я3) = Г
ва — дЯ3 у
Еп
(15)
(16)
V2 V -------- ' "" "" СР
Из уравнений (15), (16) формируются следующие два соотношения для определения а и р:
а-1 = а(&+1)-2р, а + 3р = 0,
откуда
1
(17)
а =
"3 :, Р= 1
3& + 2 3^ + 2
С учетом (17) задача (15), (16) преобразуется к следующему стационарному виду:
а-
+-
д_ 1
вс
дв
л
+ а-
д
дЧ У дЧ2
вс
дв
Л
+ а-
д
V
3^+2
дв
дв
+Ч2Т"
дЧ дЧ.
дЧ2 У дЧ3 дв
вс
дв
Л
Л
+Ч3 Л
2 дЧ3 у
+
3
3^+ 2
дЧ3 у в = 0.
+
Ш в (Ч^ Я2 , Я3 ) ^3
Е^
(18)
(19)
Пусть функция в( ч1, ч2, ч3) центрально симметрична, то есть зависит от одной координаты г сферической системы координат
Г ■■
Ч1 = г соб у- соб г = ч2 + 42 + Ч32,
Ч2 = г соб у- Бт ^ = аг^ —,
Ч1
(20)
Ч3 = г бш у, ^ = ягсбш
Ч3
л/Чо2
2 , 2 , 2 Ч1 + Ч2 + Ч3
Тогда задача (18), (19) трансформируется в следующую задачу Коши для нелинейного обыкновенного дифференциального уравнения:
а-1
г
4 d {r qв=0, (21)
r dr ^ dr) 37 + 2 dr 37 + 2
jr2 q(r)dr=p. (22)
При переходе к сферической системе координат (20) интеграл в (19)
по переменной d равен 2p, по переменной g - двум, а якобиан равен
2
r cosg.
Поскольку функция q(q1,q2,q3) является центрально симметричной и трансформируется в функцию в(r), то для уравнения (21) должно выполняться условие симметрии
(0) = 0. (23)
Таким образом, задача (21), (22) для нелинейного обыкновенного дифференциального уравнения (21) представляется в виде
(r вв)' + —-—- (в- r3)' = 0, (24)
v ; a (37 + 2) ' У J
j r 2-q( r) dr = Р (25)
в7в (0) = 0. (26) Первый интеграл уравнения (24) будет
r в7в+—1-- r в = С,
a (3s + 2) 1
причем в силу условия симметрии (26) при r = 0 постоянная интегрирования C1 равна нулю. Следовательно,
r в7в+—1-- r в = 0; (27)
a (3s + 2) v 7
pJMKjre^dr = p (28)
L r cP
Уравнение (27) - нелинейное обыкновенное дифференциальное уравнение с разделяющимися переменными. Его общее решение имеет вид
{ 7r2 7r2 Y7
в( r )= 7 ч--7—- . (29)
v ; ^ 2a (37 + 2) 2a (37 + 2))
В нем постоянная интегрирования определяется выражением
r 2
' г
0 о „2
а r02 может быть определена из условия постоянства энергии
2a (37 + 2)
(28). Из него с учетом (29) находим
4р
1/7
L
■3/2
f \ 77
7 1 jr2 (r02 -r2 )1/7 dr = EP, (30)
2a (37 + 2)] I v 0 ' cp
где re [0;r0].
Интеграл в левой части выражения (30) вычисляется в квадратурах,
Г
для чего преобразуем его к виду (— изменяется от 0/ г0 = 0 до г0 / г0 = 1)
4р
4кхКК ^
ч1/а
Е
3/2
а
2а (3а+2)
'0
3а2 1 / г ~\2 Г а _
'0
0 V г0 J
1 \1/а
1 -
V ^' J J
а
г г \
V г0 J
Ел ер
(31)
В соответствии с [6] имеем
1 т{ г\р , г(Р + 1)-Г(т +1)/2 г | wm (1 - )ра* = ^ , {Р +1, т +1}> 0,
(32)
2Г( р + (т + 3)/2) где Г(^) - гамма-функция, р = 1/ а, т = 2.
Упрощая выражение (32) и подставляя его в (31), получим постоян
ную г02
'0 (Е ):
(Е0 /ер) Еш (2а(3а+2)уа (3а+2)(а+2)Г((а+2)/2а)
2а 3а+2
. (33)
4р ^якХкцк£ V а J а Г(1/а)
Выражения (29) и (33) определяют решение задачи Коши (24)-(26) и через (13) - решение (11), (12). Возвращаясь к декартовым координатам, получаем решение исходной задачи (1), (2)
Г У/а
Т (х, у, г, X):
1
X
3/(3а+2)
а
2а (3а+2)
х
х
(аххх + ахуУ + ах2г )2 Е + (ацхх + (ХцуУ + ц )2 Е + (а?хх + а?уУ + аСгг )2 Е
2 _кХ_кЦ_кС
Г -
2Х 3а+2
1/а
(34)
где г02 определяется выражением (33).
Из решения (34) видно, что если выражение в квадратных скобках равно нулю, то поверхность второго порядка
'0 (Е):
(аХхХ + аХуУ + аХгг )2 (ацхх +ацуУ + ацгг )2 (а(хх +а^уУ + )2
(
— 2кх
3а+2 _X
V Е
2
- + -
(
2
- + -
(
2к£
3а+2 _£
V Е
2
(35)
определяет подвижный фронт, разделяющий область с ненулевой температурой и остальную часть пространства с нулевой температурой. Из выражения (35) видно, что этой областью является эллипсоид с полуосями
г X
0
— 2кх
3а+2 I_X
г X
0
30* В?
Г X
0
— 2к£
3а+2 I_£
(36)
V Е ' " V Е ' и \ Е вдоль главных осей ОХ, Оц, 0£ тензора теплопроводности. Таким образом, (35) определяет подвижный фронт тепловой волны, разделяющей область ненулевого решения и область, не затронутую тепловыми возмущениями.
2
2
2
Из (33) и (34) видно, что решение не существует при а = -2/3 и
а = 0. Кроме этого, из условия р +1 =—+1 > 0 в (32) следует, что
а
сте(-¥;-1)и(-1; -2/3)и( 0; , а также аргументы гамма-функций в (33) не
должны быть отрицательными целыми и нуле, т.е.
1 а+2
— Ф-п и -Ф-п, п = 0,1,2,....
а 2а
В задаче с нулевой начальной температурой отрицательные значения а приводят к бесконечной теплопроводности. Поэтому из указанных промежутков а остается интервал ае(0;.
Анализ результатов
Анализ решения (34) в трехмерном случае проводить громоздко и неудобно. Поэтому по изложенному методу получено решение двухмерной задачи на плоскости. Оно имеет вид
Т (X, у, г ) = -
Д/(а+1)
а
1
2а (2а + 2)
Ь
(х 008 (р + у ътф) — + (—х 8т( + у со8(('
Ь
1/(а+1)
,(37)
где ( - угол между осью ОХ и декартовой осью Ох, а г02 определяется вы-
ражением
() =
Г(Е0/ер) а+1 У+1 (4а(а+1)'
а
(а+1)
а
(38)
Из анализа условия —+1 > 0 в (32) и а ф -1, аФ 0 следуют интервалы
а
а, при которых формально существует решение (37)
ае(-¥;-1)и( 0; ~). (39)
Из равенства нулю в (37) выражения в квадратных скобках следует геометрическое место точек на плоскости
Г2 (Е)
(х соб(+у б1п()2 Ь + (- х Бт(+у соб()2 —
_кХ_кп
г
1/(а+1)
(40)
являющееся кривой второго порядка, которая определяет подвижный фронт волны, разграничивающей область ненулевого решения, и остальную область с начальным распределением температуры (нулевым). Таким образом, фронты бегущей тепловой волны по холодному пространству на плоскости представляют собой эллипсы с полуосями
1 = гг 2а+2 • 1 = гг 2а+2
1х - V ^ ь ' П 01 \ Ь '
а
Г
1
2
Л
Г
0
2
Г,
Из анализа выражения (37) следует, что на фронте бегущей волны температура непрерывна, производные первого порядка по пространственным переменным разрывны, плотности тепловых потоков за счет наличия сомножителя Та непрерывны, а вторые производные температуры по пространственным переменным разрывны.
На рис. 1-3 показаны распределения температур вдоль каждой из координатных осей при трех фиксированных моментах времени, полученные по формулам (37), (38) со следующими входными данными: а = 1; а = 2; а = 3; ( = 0; кХ = 5Вт/(м-Ка+1); кп = 1Вт/(м-Ка+1);
ер = 2000 Дж/(м К), Е0=1000 Дж. Результаты приведены в виде сечений по осям у = 0 и х = 0. При этом величина Е0/ ер в соответствии с формулой (38) входит в выражение для г02, характеризующего координаты границы фронта тепловой волны.
Рисунки подтверждают волновой характер распределения температур с фронтом тепловой волны на изотерме Т = 0, описываемой эллипсами с полуосями (41). При г ® 0 + 0 решение (37) стремится к дельта-функции.
Видно, что областями ненулевого решения в различные моменты времени являются области, ограниченные эллипсами с разными скоростями распространения тепловых волн в различных направлениях.
Для любого угла ( из выражения (37) находим координаты х фронта волны при у = 0 и координаты у при х = 0, когда Т (х, у, г ) = 0
х = г0г
2а+2
Ь 2 Ь • 2
-С08 --81П (
V кХ К )
-1/2 , ^ ч-1/2
у = Г0г
2а+2
Ь 2 Ь 2 -81П --С08 (
V кХ К )
откуда определяются скорости движения точек фронта тепловой волны
dх г.
-1/2
V =■
dг 2а+2
у = Ф = Г0
у dг 2а+2
Ь 2 Ь 2а+1
С08 --81П (
V кХ К )
г 2а+2, (42)
(г т ^-1/2 2а+1 Ь • 2 Ь 2 -81П --С08 (
V кХ К
г 2а+2. (43)
Результаты показывают, что чем ближе время к начальному моменту, тем выше температура, а с увеличением времени она резко падает. Так при а = 1 в момент времени г = 10-6 е температура в начале координат составляет примерно 3000К, а при г = 10-5е и а = 1 □ 1000К. При наносекунд-ных длительностях температура в начале координат может значительно превышать температуру фазовых переходов любых материалов. При возрастании а температура резко падает, так как с увеличением теплопроводности тепловая энергия интенсивнее отводится на периферию от места приложения источника энергии.
а б
Рис. 1. Температура и температурный фронт тепловой волны для а = 1 в моменты времени 1 - х = 10-6 е, 2 - х = 2 10-6 е, 3 - х = 5 10-6 е: а - сечение
у = 0, б - сечение х = 0.
а б
Рис. 2. Температура и температурный фронт тепловой волны для а= 1 в моменты времени 1 - х = 10-5 е, 2 - х = 2 10-5 е, 3 - х = 5 10-5 е: (а) - сечение
У = 0, (б) - сечение х = 0.
а
б
Рис. 3. Температура и температурный фронт тепловой волны для а = 2 в моменты времени 1 - х = 2, 2 - х = 4е, 3 - х = 8с: (а) - сечение у = 0,
(б) - сечение х = 0.
Заключение
1. Изложен метод получения аналитического решения задач теплопроводности с нелинейными теплофизическими характеристиками в анизотропном пространстве от импульсного источника тепловой энергии.
2. Оказалось, что процесс теплопереноса носит волновой характер, определяемый теплопроводностью в виде однородного многочлена от температуры с четко выраженным фронтом тепловой волны, нестационарно продвигающемся по холодному пространству. При этом, в соответствии с анизотропным характером теплопереноса, фронты тепловых волн являются эллипсоидами в трехмерном случае и эллипсами - в двухмерном.
3. Установлено, что на фронтах тепловых волн температура непрерывна, первые и вторые производные температуры по пространственным переменным разрывны, а плотности тепловых потоков непрерывны вследствие наличия сомножителя Та, что является неожиданным, так как производные первого порядка разрывны, а тепловые потоки непрерывны.
Работа выполнена при финансовой поддержке грантов Президента Российской Федерации НШ-1387.2014.8, МД-3297.2013.8 и РФФИ 12-0133095, 14-01-00479.
Список литературы
1. Самарский А. А., Галактионов В. А., Курдюмов С.П., Михайлов А.П. Режимы с обострением в задачах для квазилинейных параболических уравнений. М.: Наука, 1987. 480с.
2. Лыков А.В. Теория теплопроводности. М.: Высшая школа, 1967.
600с.
3. Формалев В.Ф. О тепловых ударных волнах в нелинейных твердых средах // ТВТ. 2012. Т. 50. № 6. С. 799-803.
4. Соболев С. Л. Процессы переноса и бегущие волны в неравновесных системах // Успехи физических наук. 1991. Т. 161. № 3. С. 5-29.
5. Формалев В.Ф., Кузнецова Е.Л. Тепломассоперенос в анизотропных телах при аэрогазодинамическом нагреве. М.: Изд-во МАИ, 2010. 308с.
6. Двайт Г.Б. Таблицы интегралов и другие математические формулы. М.: Наука, 1964. 228с.
7. Кузнецова Е.Л., Формалев В.Ф. О тепловых волнах в нелинейном анизотропном пространстве // Инженерная физика. Серия «Теплофизика и тепломеханика». 2010. № 5. С. 43-47.
8. Кузнецова Е.Л. Математическое моделирование тепломассопере-носа в композиционных материалах при высокотемпературном нагреве в элементах ракетно-космической техники. М.: Изд-во МАИ-ПРИНТ. 2010. 158с.
9. Лурье С.А., Рабинский Л.Н., Соляев Ю.О., Нгуен Д.К., Исследование локальных эффектов в распределении температурных напряжений на контактных границах слоистых сред. Электронный журнал «Труды МАИ», выпуск №71, М., МАИ, 2013
Формалев Владимир Федорович, д-р ф.-м. нпук, проф., Россия, Москва, Московский авиационный институт (национальный исследовательский университет),
Кузнецова Екатерина Львовна, д-р ф.-м.наук, проф., ведущий научный сотрудник, Россия, Москва, Московский авиационный институт (национальный исследовательский университет),
Рабинский Лев Наумович, д-р ф.-м.наук, проф., декан, Россия, Москва, Московский авиационный институт (национальный исследовательский университет)
INTERACTION OF THE PULSE SOURCE OF THERMAL ENERGY WITH ANISOTROPIC SPACE WHICHHEATPHYSICAL CHARACTERISTICS DEPEND ON TEMPERATURE
V.F. Formalev, E.L. Kuznetsova, L.N. Rabinskiy
In article the analytical solution of a problem of heat transfer in anisotropic space is for the first time received, components of which tensor of heat conductivity depend on temperature. Heat transfer is initiated by a pulse dot source of thermal energy. It appeared that heat transfer has wave character with a final speed of distribution of the front of a thermal wave though the equation of heat conductivity has parabolic type. Cases of various sedate dependences of components of a tensor of heat conductivity on temperature are investigated, is shown that fronts of thermal waves have an appearance of ellipsoids in space and ellipses on the plane. Results can be used at research of heat exchange in composite materials at their laser radiation.
Key words: impulse of thermal energy, anisotropic space, nonlinear components of a tensor of heat conductivity, thermal wave, geometry of the front of a thermal wave, speed of distribution of a thermal wave.
Formalev Vladimir Fedorovich, doctor of physical and mathematical sciences, professor, formalev38 a mail.ru, Russia, Moscow, Moscow aviation institute (national research university)
Kuznetsova Ekaterina Lvovna, doctor of physical and mathematical sciences, professor, the leading researcher, lareynaa mail.ru, Russia, Moscow, Moscow aviation institute (national research university)
Rabinsky Lev Naumovich, doctor of physical and mathematical sciences, professor, the dean, rabinskiy@,mail. ru, Russia, Moscow, Moscow aviation institute (national research university)
