CC BY
20
УДК 536.2
Применение метода наименьших квадратов к задаче о переносе излучения в шаровой полости
© В.С. Зарубин, О.В. Пугачев, И.Ю. Савельева МГТУ им. Н.Э. Баумана, Москва, 105005, Россия
Многие используемые в технике теплозащитные материалы имеют пористую структуру. При интенсивном тепловом воздействии возникает необходимость учитывать перенос тепловой энергии путем излучения в порах таких материалов. Построена математическая модель, описывающая теплообмен излучением в шаровой полости, форму которой можно рассматривать как среднюю статистическую по отношению к формам замкнутых пор в твердых телах. Для количественного анализа этой модели использован метод наименьших квадратов. Введен эквивалентный коэффициент теплопроводности условной сплошной среды, заполняющей пору, что позволяет рассматривать материал с пористой структурой как сплошное неоднородное твердое тело.
Ключевые слова: шаровая полость, математическая модель переноса излучения, метод наименьших квадратов.
Введение. Перенос тепловой энергии в пористых теплозащитных материалах наряду с теплопроводностью происходит и путем теплового излучения [1, 2]. К таким материалам можно отнести, например композиты с сотовым заполнителем [3, 4], конструктивно-орто-тропные системы [5] и сферопластики [6]. Пористая структура теплозащитного материала может формироваться в процессе теплового воздействия вследствие термодеструкции [7]. Для получения количественной оценки интенсивности теплопереноса излучением в таких материалах целесообразно использовать методы математического моделирования [8, 9].
Существуют различные подходы к построению математических моделей процесса переноса тепловой энергии излучением через замкнутые поры в твердых телах. При этом пору условно заменяют включением, материал которого имеет некоторый эквивалентный коэффициент теплопроводности Xк , подлежащий определению с использованием той или иной математической модели. В книгах [10, 11] рассмотрены работы, посвященные определению величины Xк для замкнутых пор различных форм, среди которых в качестве средней статистической часто принимают форму шаровой полости. В большинстве этих работ использованы математические модели, в которых не учитывается градиент температуры в окрестности поры. Однако влияние этого фактора возрастает по мере интенсификации теплового воздействия на
пористые теплозащитные материалы [12, 13]. В статье представлено решение задачи о переносе излучения в шаровой полости с учетом неоднородного распределения температуры в окрестности этой полости с применением метода наименьших квадратов.
Модель теплопереноса излучением в шаровой полости. Свойства непрозрачной поверхности шаровой полости радиусом r0 примем соответствующими свойствам диффузно-серой поверхности с коэффициентом излучения 8 = const, а среду в этой полости будем считать диатермичной, т. е. не поглощающей и не рассеивающей излучение [14]. Для определения теплового потока, проходящего через эту полость при неоднородном распределении температуры T(P) на ее поверхности S (P е S ), введем плотности потоков падающего q* (P) и собственного 8 q0(P) излучений, где, согласно закону Стефана —
4 —8
Больцмана [15], q0(P) = a0T (P), а а0 = 5,67-10
Вт/(м2- К4). Пада-
ющее на непрозрачную поверхность излучение частично поглощается (Ад*) и частично отражается (Яд*), для диффузно-серой поверхности коэффициент поглощения А = 8 и коэффициент поглощения Я = 1 - А = = 1 - 8. Собственное излучение после сложения с отраженным излучением составляет эффективное излучение с плотностью потока
q* (P) = 8q0 (P) + Rq* (P), P е S.
(1)
Рассмотрим на сферической поверхности $ полости две произвольные точки N и Р (рис. 1). Элементарная площадка ^) в окрестности точки N е $ посылает на единичную площадку в окрестности точки Р е $ поток эффективного излучения
dq* (P) = q*(N)dФNp,
(2)
где й— элементарный угловой коэффициент [10], равный, в соответствии с законом Ламберта [11] для распределения диффузного излучения по направлениям, определяемым углом РN между отрезком ^ и нормалью к площадке ^),
а ф ш = С05 ^ 2С05 Рр а$ (N ) =
'кlpN
Рис. 1. Схема к определению элементарного углового коэффициента в шаровой полости
dS (N) 4лг02
(3)
где РР = РN = arccos(lPN /2г0) — угол между нормалью в точке Р и отрезком PN длиной 1т. Формулы (1)-(3) справедливы для еди-
ничной площадки в окрестности любой точки, принадлежащей поверхности S. Поэтому с учетом данных формул после интегрирования по поверхности полости приходим к интегральному уравнению Фредгольма второго рода в виде [16]
q*(P) — |У(N)dS(N) = sqo(P), P, N е S. (4)
W S
Отсюда с учетом равенства (1) следует
q. (P) = -Ц- fq*(N)dS(N) = const, (5)
W S
т. е. плотность потока падающего излучения в шаровой полости одинакова для всех точек ее поверхности и равна средней плотности потока эффективного излучения.
Умножив уравнение (4) на
dS(P)/ (4toq ), после интегрирования по сферической поверхности S при условии R = 1 — s = const вместо равенства (5) получим
1
q. =——
4ro) S
Количество энергии, теряемое путем излучения единицей площади поверхности твердого тела в единицу времени, называют плотностью потока результирующего излучения [17], равной q° =sq0 —Aq., или, для шаровой полости с учетом формулы (6) для каждой точки P е S,
q (P) = sqo(P) —jqo(P)dS(P) = «si^P —|T4(P)dS(P). (7) 4лГо S 4лг02 S
jqo(M )dS (M ). (6)
Отсюда в соответствии с законом сохранения тепловой энергии при ее переносе излучением следует
jq° (P ) dS ( P ) = 0.
При этом суммарный тепловой поток, передаваемый через полость, можно вычислить по формуле
ео = г д (Р) (р)\^(р) (8)
» 9
Температурное поле в окрестности шаровой полости. Примем, что шаровая полость радиусом r0 находится в неограниченной области, заполненной однородным материалом с коэффициентом теплопроводности X. В центре полости поместим начало сферической системы координат r, 9, %, а на большом расстоянии r от центра полости зададим вектор градиента температурного поля в этом материале, имеющий модуль G и направленный вдоль оси, от которой происходит отсчет угловой координаты 9. В этом случае при r ^ да установившееся осесимметричное (не зависящее от угловой координаты % ) распределение температуры в этом материале будет описывать функция (r, 9) = T0 + Gr cos 9, где T0 — температура в плоскости при 9 = л /2. Эта функция удовлетворяет уравнению Лапласа, которое в сферических координатах с учетом осевой симметрии имеет вид [17]
1 d
i
r
2 dr
дТ Л 1 д
i
-2 -1 +
dr J r sin 9д9
sin 9—| = 0. (9)
d9 J
Наличие в однородном материале шаровой полости вызовет возмущение температурного поля, описываемое также удовлетворяющими уравнению (1) дополнительными слагаемыми [18]
Б,
P, (cos 9), Бj = const, j = 1,2,..., (10)
rj+1
где Pj (cos 9) — ортогональные полиномы Лежандра степени j: P (cos 9) = cos 9, P2 (cos 9) = (3cos2 9-1)/2 и для j >2
jPj (cos 9) = (2 j -1)P,-1 (cos 9)cos 9 - (j - 1)Р,-2(cos 9). (11)
Среднее по поверхности сферы значение любого из ортогональных полиномов Лежандра при j >0 равно нулю, поскольку
л
Jp, (cos 9)sin 9d 9 = 0. (12)
Это свойство обеспечивает выполнение в полости условия теплового баланса, необходимого для существования установившегося распределения температуры на ее поверхности. Следует отметить, что уравнению (9) удовлетворяет также частное решение, пропорциональное 1/ г, но оно противоречит физическому смыслу рассматриваемой задачи, поскольку описывает либо только подвод теплоты путем теплопроводности к полости, либо только ее отвод из полости.
Температурное поле в окрестности шаровой полости должно удовлетворять на ее поверхности нелинейному граничному условию
F (0) = АдТ (''-в)
дг
- q (в) = 0, 0е[0,я]. (13)
г=г0
Нелинейность этого условия вызвана тем, что плотность ({ (9) результирующего излучения, согласно формуле (7), зависит от четвертой степени значений температуры Т*(9) = Т (г0,9) на поверхности полости. В силу нелинейности условия (13) не удается в аналитическом виде найти коэффициенты В^ в дополнительных слагаемых,
определяемых равенствами (10). Один из путей вычисления этих коэффициентов состоит в применении метода наименьших квадратов, потребовав минимума интегральной квадратичной невязки условия (13), которую можно представить в виде
л
А* = | ( ^(9))281И 9d9, (14)
где F (0) = F (0)/(M0/ го).
Запишем с учетом дополнительных слагаемых распределение на поверхности полости безразмерной температуры
T*(0) = = 1 + G cos 0 + YjBjPj (cos 0), (15)
То j=i
где G = Gr0 / T0 и Bj = Bj /(r0J+1T0). После перехода в равенстве (7) от
T(P) к T*(0), подстановки q (P) = q (0) в формулу (14) и выполнения
*
интегрирования невязка А будет функцией лишь безразмерных неизвестных коэффициентов Bj . В соответствии с необходимым условием минимума этой невязки коэффициенты Bj должны удовлетворять бесконечной системе нелинейных алгебраических уравнений * ^
дА / dBj = 0, J = 1,2,.. . Ясно, что полное решение такой системы
недостижимо. Последовательно увеличивая выбор наибольшего значения j и рассматривая таким образом усеченные системы уравнений, можно вычислить конечное число безразмерных коэффициентов — *
Bj , контролируя при этом убывание значения А . В дальнейшем ограничимся рассмотрением первого приближения.
Первое приближение. Распределение температуры в окружающей шаровую полость среде в первом приближении с учетом соотношения (10) представим равенством
Т (r, 9) = Т0 + (Gr + Б1/ r 2)cos 9. (16)
Прежде всего оценим границы возможных значений коэффициента Б1 , рассматривая предельные случаи как отсутствия теплопереноса в шаровой полости, так и весьма интенсивного теплопереноса. Если ввести эквивалентный коэффициент XR теплопроводности условной среды в полости, то отсутствию теплопереноса будет соответствовать значение XR =0, а предельной интенсификации теплопереноса XR ^ да.
В силу ограниченности температуры условной среды в полости при r = 0 уравнению (9) будет удовлетворять в этой среде распределение температуры [19]
Т° (r, 9) = Т0 + B°r cos 9, Б° = const. (17)
В соотношения (16) и (17) входят два неизвестных коэффициента Б1 и Б°, которые можно найти из условий непрерывности при r = r0 распределения температуры и плотности теплового потока
Т(r, 9) = Т° (r0, 9); ХдТ / dr |r=r0 = Xr дТ° / dr |r=r0.
Отсюда с использованием формул (16) и (17) найдем
G + Б1/ r03 = Б°; G - 2Б1 / r03 = (XR / X)Б°,
или Б° = 3GX / (2X + XR ) и Б1 / r03 = G(X-XR )/(2X + XR ). Таким образом, при XR = 0 получим Б1 / r0 = G / 2, а при XR ^ да —
3 — 2
Б1 / r0 ^-G. В итоге для безразмерной величины Б1 = Б1/(г0Т0) промежуток возможных значений соответствует полуинтервалу (-G,G / 2]. С помощью эквивалентного коэффициента XR теплопроводности условной среды можно представить проходящий через шаровую полость суммарный тепловой поток
е = 2rcr02Xr J |rr0 sin 9d9^r02XrB° = ^r02:3XR+XL. (18)
В рассматриваемом первом приближении равенство (15) для безразмерной формы распределения температуры на поверхности полости примет вид
T* = 1 + (G + B^cos 0. (19)
Эта формула в сочетании с равенством (7) после вычисления интеграла в правой части этого равенства позволяет записать
q0 (0)=^Гг1=ё(1+(G+Bi)cos 0)4-8(1+2(G+Bi)2+!(G+Д)4), (20)
где 8=8G0T)3r0/ X. Тогда вместо формулы (13) с учетом равенства (16) получим
F^G) = G-2B1cos0 +1 + 2(G+B1)2 + ^(G+B1)4 - (1 + (G+B1)cos0)4. (21) 8 8 5
Использовав формулу (21) при вычислении интеграла в соотно-
* —
шении (14), из условия ЗА / 6B1 = 0 получим алгебраическое уравнение седьмой степени относительно коэффициента B1 в виде
¿oB = 0, (22)
n=0
где a0 =4G7/25 + 12G5/5 + 8G3 + G(1 +1/8-2/82)/3; a1=28G6/25 + + 12G4 + 24G2 + (1+2/8)2/3; a2 =84G5/25 + 24G3 + 24G; a3 = 28G4/5 + + 24G2 + 8; a4 =28G3 /5 + 12G; a5=84G2/25 +12/5; a6 =28G/25; a7 = 4/ 25. Отметим, что от параметра г зависят лишь коэффициенты a0 и a1, причем a0 =0 при G = 0, т. е. отсутствие при r ^го градиента температуры приводит к имеющему физический смысл решению уравнения (22) B1 =0.
Анализ первого приближения. При достаточно малых значениях радиуса r0 шаровых пор и допустимых шкалой абсолютных температур значений G получим, что G = Gr0 / T0 ^ 1, а из установленного промежутка (—G, G /2] возможных значений B1 следует | B11. Эти неравенства при количественном анализе первого приближения дают основание сначала в уравнении (22) и в его коэффициентах удержать лишь первые степени B1 и G, а затем провести сравнение с результатами расчетов с использованием уравнения (22).
В линейном приближении из уравнения (22) найдем
-2
-2-8+8' (2 + 8 )2
B1 = G ^ _ч2 . (23)
в = 0,02
Параметр в можно рассматривать (подобно эквивалентному коэффициенту Хк теплопроводности условной среды в шаровой полости) как величину, характеризующую интенсивность теплопереноса излучением. Если в качестве ориентира принять в = 0,86, Т0 = 1600 К,
—3 ^
Г0 =10 м и X = 0,1 Вт/(м • К), то получим в « 2. В предельных случаях при в =0 из равенства (23) следует В1 ^ О /2, а при В^да имеем В1 ^ —О, т. е. приходим к ранее установленным границам возможных значений безразмерного коэффициента В1.
На рис. 2 в полулогарифмических координатах приведены результаты расчетов по уравнению (22) в виде зависимостей коэффициента В1 от параметра в при различных значениях параметра О (сплошные, штриховые и штрихпунктирные кривые). Штрихпунктирная линия построена по формуле (23) при значении О = 0,2 и мало отличается от соответствующей штриховой кривой. При О <0,2 по мере
убывания параметра О отличие результатов расчетов по формулам (22) и (23) быстро уменьшается (кривая, построенная по формуле (23) при О = 1, в масштабе рис. 2 полностью совпадает с соответствующей штрихпунктирной кривой).
Представленные на этом рисунке зависимости позволяют с помощью полученного ранее равенства В1 / г0 = О(X — Хк) / (2Х + Хк) вычислить для условной среды в полости относительный эквивалентный коэффициент теплопроводности
X = Хк = О — 2 В1 к х О+В1 .
При изменении угла 0 в промежутке от 0 до л плотность результирующего излучения монотонно убывает, принимая в этом промежутке нулевое значение. Из соотношения (20) следует, что для рассматриваемого первого приближения значению ({ =0 соответствует угол
Рис. 2. Зависимость коэффициента В1 от параметра в при различных значениях параметра О
0* = arccos
(1 + 2(G + B1)2 + (G + B1)4 / 5)1/4 -1 G + B1
e= 10,0
0,01 0,02 0,05 0,1 0,2
1,0 G
Рис. 3. Зависимость угла 0* от параметра О при различных значениях параметра г
Поскольку B1 — 0 при G — 0, предельным переходом из этой формулы можно установить, что 0* — п /2 при G — 0. На рис. 3 по этой формуле с использованием результатов количественного анализа уравнения (22) построены в полулогарифмических координатах зависимости угла 0* от параметра
G для различных значений параметра s. Формально при 8 —> 0 и G = 1 предельное значение B1 =0,5.
В этом случае 0* «1,1611, что соответствует углу ~ 66,53°.
Таким образом, в полуинтервале (0*; п /2] изменения угла 0 на сферической поверхности рассматриваемой шаровой полости в данном приближении возникает несоответствие между подводимым к этой поверхности путем теплопроводности и отводимым путем излучения тепловыми потоками, т. е. происходит локальное нарушение закона сохранения энергии. На рис. 4 в логарифмических координатах для различных значений параметра s приведены зависимости от параметра
G отношения ^ = cos2 0* теплового потока, подводимого к полости от окружающей среды путем теплопроводности через круговой пояс сферической поверхности при 0е(0*; п /2), ко всему подводимому к полости тепловому потоку. При s — 0 и G = 1 значению B1 =0,5 соответствует предельное значение ^^0,15866; его следует рассматривать в качестве оценки наибольшей возможной интегральной погрешности, которая может возникнуть при использовании рассматриваемого пер-
Е = 0,01
0,1 0,2
л 0,1
0,05
0,02 0,01 0,005
0,002 0,001 0,0005
0,0002 0,0001
\
S> / -7 -
XI Г /
/// / /
л S/. /
7-
-Л / / /
/ / /
А 4 % / / /
т У —г t-
-1» /
W / / /
V, / /
0,5 2,0
5,0
10,0
0,01 0,02 0,05 0,1 0,2
0,5 1,0 G
Рис. 4. Зависимость параметра ^ от G при различных значениях s
вого приближения. При фиксированном значении параметра О значение ^ достаточно быстро убывает по мере возрастания параметра
в. Если к тому же учесть, что реальные значения параметра О существенно меньше единицы, то интегральную погрешность первого приближения можно считать сопоставимой с погрешностью используемых исходных данных.
Заключение. Методом наименьших квадратов построено приближенное решение нелинейной задачи теплопереноса путем излучения в шаровой полости. Проведена оценка наибольшей возможной интегральной погрешности этого решения. Введено понятие эквивалентного коэффициента теплопроводности условной сплошной среды, заполняющей пору, что позволяет рассматривать материал с пористой структурой как сплошное неоднородное твердое тело. Получена формула для вычисления этого коэффициента в первом приближении.
Работа выполнена по грантам НШ-1432.2014.8 и МК-6573.2015.8 программ Президента РФ государственной поддержки ведущих научных школ и молодых кандидатов наук, а также в рамках проекта 1712 в сфере научной деятельности в части государственного задания № 2014/104 Минобр-науки РФ и государственного задания по проекту № 1.2640.2014.
ЛИТЕРАТУРА
[1] Елисеев В.Н., Товстоног В.А. Теплообмен и тепловые испытания материалов и конструкций аэрокосмической техники при радиационном нагреве. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2014, 396 с.
[2] Димитриенко Ю.И. Механика композиционных материалов при высоких температурах. Москва, Машиностроение, 1997, 368 с.
[3] Димитриенко Ю.И., Федонюк Н.Н., Губарева Е.А., Сборщиков С.В., Прозоровский А.А. Многомасштабное конечно-элементное моделирование трехслойных сотовых композитных конструкций. Наука и образование: электронное научно-техническое издание, 2014, № 7. doi: 10.7463/0714.0717805
[4] Димитриенко Ю.И., Федонюк Н.Н., Губарева Е.А., Сборщиков С.В., Прозоровский А.А., Ерасов В.С., Яковлев Н.О. Моделирование и разработка трехслойных композиционных материалов с сотовым заполнителем. Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки, 2014, № 5, с. 66-81.
[5] Димитриенко Ю.И., Губарева Е.А., Сборщиков С.В. Асимптотическая теория конструктивно-ортотропных пластин с двухпериодической структурой. Математическое моделирование и численные методы, 2014, № 1, с. 36-57.
[6] Димитриенко Ю.И., Сборщиков С.В., Соколов А.П., Садовничий Д.Н., Гафаров Б.Р. Численнное и экспериментальное моделирование прочностных характеристик сферопластиков. Композиты и наноструктуры, 2013, № 3, с. 35-51.
[7] Димитриенко Ю.И., Минин В.В., Сыздыков Е.К. Численное моделирование процессов тепломассопереноса и кинетики напряжений в термоде-структирующих композитных оболочках. Вычислительные технологии, 2012, т. 17, № 2, с. 43-59.
[8] Зарубин В.С. Математическое моделирование в технике. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2010, 496 с.
[9] Зарубин В.С., Кувыркин Г.Н. Особенности математического моделирования технических устройств. Математическое моделирование и численные методы, 2014, № 1 (1), с. 5-17.
[10] Дульнев Г.Н., Заричняк Ю.П. Теплопроводность смесей и композиционных материалов. Ленинград, Энергия, 1974, 264 с.
[11] Чудновский А.Ф. Теплофизические характеристики дисперсных материалов. Москва, Физматгиз, 1962, 456 с.
[12] Зарубин В.С., Кувыркин Г.Н. Математическое моделирование термомеханических процессов при интенсивном тепловом воздействии. Теплофизика высоких температур, 2003, т. 41, № 2, с. 300-309.
[13] Зарубин В.С., Кувыркин Г.Н., Савельева И.Ю. Радиационно-кондук-тивный теплоперенос в шаровой полости. Теплофизика высоких температур, 2015, т. 53, № 2, с. 243-249.
[14] Зигель Р., Хауэлл Дж. Теплообмен излучением. Пер. с англ. Москва, Мир, 1975, 936 с.
[15] Исаев С.И., Кожинов И.А., Кофанов В.И. и др.; Леонтьев Ф.И., ред. Теория тепломассообмена. Москва, Изд-во МГТУ им. Н.Э. Баумана, 1997, 684 с.
[16] Зарубин В. С. Оценка эквивалентного коэффициента теплопроводности при переносе излучения в шаровой полости. Инженерный журнал: наука и инновации, 2013, вып. 8. URL: http://engjournal.ru/catalog/mathmodel/ hidden/887.html
[17] Зарубин В.С., Кувыркин Г.Н. Математические модели механики и электродинамики сплошной среды. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2008, 512 с.
[18] Корн Г., Корн Т. Справочник по математике для научных работников и инженеров. Пер. с англ. Москва, Наука, 1968, 720 с.
[19] Зарубин В.С., Кувыркин Г.Н., Савельева И.Ю. Эффективный коэффициент теплопроводности композита с шаровыми включениями. Тепловые процессы в технике, 2012, № 10, с. 470-474.
Статья поступила в редакцию 18.05.2015
Ссылку на эту статью просим оформлять следующим образом:
Зарубин В.С., Пугачев О.В., Савельева И.Ю. Применение метода наименьших квадратов к задаче о переносе излучения в шаровой полости.
Математическое моделирование и численные методы, 2015, № 4 (8), с. 53-65.
Зарубин Владимир Степанович — д-р техн. наук, профессор кафедры «Прикладная математика» МГТУ им. Н.Э. Баумана. e-mail: fn2@bmstu.ru
Пугачев Олег Всеволодович — д-р физ.-мат. наук, профессор кафедры «Прикладная математика» МГТУ им. Н.Э. Баумана. e-mail: fn2@bmstu.ru
Савельева Инга Юрьевна — канд. физ.-мат. наук, доцент кафедры «Прикладная математика» МГТУ им. Н.Э. Баумана. e-mail: fn2@bmstu.ru
Application of the least squares method to the problem of radiation transfer in a spherical cavity
© V.S. Zarubin, O.V. Pugachev, I.Yu. Savelyeva
Bauman Moscow State Technical University, Moscow, 105005, Russia
A lot of heat-shielding materials used in engineering have porous structure. When there is an intensive thermal exposure, there occurs a necessity to consider thermal energy transfer by means of radiation in pores of such materials. We contructed a mathematical model describing heat exchange by radiation in a spherical cavity. Its form can be considered as an average statistical form in relation to forms of closed pores in solid bodies. For the quantitative analysis of this model we used the method of the least squares and introduced an equivalent coefficient of thermal conductivity in the conditional continuous environment filling a pore. This allows to regard the material with porous structure as a continuous non-uniform solid body.
Keywords: spherical cavity, mathematical model of radiative transfer, least squares method.
REFERENCES
[1] Eliseev V.N., Tovstonog V.A. Teploobmen i teplovye ispytaniya materialov i konstruktsiy aierokosmicheskoy tekhniki pri radiatsionnom nagreve [Heat exchange and thermal tests of materials and designs of space equipment at radiation heating]. Moscow, BMSTU Publ., 2014, 396 p.
[2] Dimitrienko Yu.I. Mekhanika kompozitsionnykh materialov pri vysokikh temperaturakh [Mechanics of composite materials at high temperatures]. Moscow, Mashinostroenie Publ., 1997, 368 p.
[3] Dimitrienko Yu.I., Fedonyuk N.N., Gubareva E.A., Sborshchikov S.V., Prozorovskiy A.A. Nauka i obrazovanie: elektronnoe nauchno-tekhnicheskoe izdanie — Science and Education: Electronic Scientific and Technical Edition, 2014, no. 7. doi: 10.7463/0714.0717805
[4] Dimitrienko Yu.I., Fedonyuk N.N., Gubareva E.A., Sborshchikov S.V., Prozorovskiy A.A., Erasov V.S., Yakovlev N.O. VestnikMGTU im. N.E. Baumana. Ser. Estestvennye nauki — Herald of Bauman Moscow State Technical University. Ser. Natural Sciences, 2014, no. 5, pp. 66-81.
[5] Dimitrienko Yu.I., Gubareva E.A., Sborshchikov S.V. Matematicheskoe modelirovanie i chislennye metody — Mathematical Modeling and Computational Methods, 2014, no. 1, pp. 36-56.
[6] Dimitrienko Yu.I., Sborshchikov S.V., Sokolov A.P., Sadovnichiy D.N., Gafarov B.R. Kompozity i nanostruktury — Composites and Nanostructures, 2013, no. 3, pp. 35-51.
[7] Dimitrienko Yu.I., Minin V.V., Syzdykov E.K. Vychislitelnye tehnologii — Computing Technologies, 2012, vol. 17, no. 2, pp. 43-59.
[8] Zarubin V.S. Matematicheskoe modelirovanie v tekhnike [Mathematical modeling in engineering]. Moscow, BMSTU Publ., 2010, 496 p.
[9] Zarubin V.S., Kuvyrkin G.N. Matematicheskoe modelirovanie i chislennye metody — Mathematical Modeling and Computational Methods, 2014, no. 1 (1), pp. 5-17.
[10] Dulnev G.N., Zarichnyak Yu.P. Teploprovodnost smesey i kompozitsionnykh materialov [Heat conductivity of mixes and composite materials]. Leningrad, Energiya Publ., 1974, 264 p.
[11] Chudnovskiy A.F. Teplofizicheskie kharakteristiki dispersnykh materialov [Heat and physical characteristics of disperse materials]. Moscow, Fizmatgiz Publ., 1962, 456 p.
[12]Zarubin V.S., Kuvyrkin G.N. Teplofizika vysokikh temperatur — High Temperature, 2003, vol. 41, no. 2, pp. 300-309.
[13]Zarubin V.S., Kuvyrkin G.N., Saveleva I.Yu. Teplofizika vysokikh temperatur — High Temperature, 2015, vol. 53, no. 2, pp. 243-249.
[14]Zigel R., Hauell Dzh. Teploobmen izlucheniem [Heat exchange by radiation]. Moscow, Mir Publ., 1975, 936 p. (in Russian).
[15]Isaev S.I., Kozhinov I.A., Kofanov V.I., et al. Teoriya teplomassoobmena [Theory of a heat and mass exchange]. Leontev A.I., ed. Moscow, BMSTU Publ., 1997, 684 p.
[16] Zarubin V.S. Inzhenernyiy zhurnal: nauka i innovatsii — Engineering Journal: Science and Innovation, 2013, no. 8. Available at: http://engjournal.ru/ catalog/mathmodel/hidden/887.html
[17] Zarubin V.S., Kuvyrkin G.N. Matematicheskie modeli mekhaniki i elektrodinamiki sploshnoy sredy [Mathematical models of mechanics and electrodynamics of the continuous environment]. Moscow, BMSTU Publ., 2008, 512 p.
[18]Korn G., Korn T. Spravochnik po matematike dlya nauchnykh rabotnikov i inzhenerov [The reference book on mathematics for scientists and engineers]. Moscow, Nauka Publ., 1968, 720 p. (in Russian).
[19]Zarubin V.S., Kuvyrkin G.N., Savelyeva I.Yu. Teplovyeprotsessy v tekhnike — Thermal Processes in Engineering, 2012, no. 10, pp. 470-474.
Zarubin V.S., Dr. Sci. (Eng.), Professor of the Applied Mathematics Department, Bauman Moscow State Technical University. e-mail: fn2@bmstu.ru
Pugachev O.V., Dr. Sci. (Phys. & Math.), Professor of the Applied Mathematics Department at Bauman Moscow State Technical University. e-mail: fn2@bmstu.ru
Savelyeva I.Yu., Cand. Sci. (Phys. & Math.), Assoc. Professor of the Applied Mathematics Department at Bauman Moscow State Technical University. e-mail: fn2@bmstu.ru
