Об одном методе решения задач теплообмена при течении жидкостей в плоских каналах Текст научной статьи по специальности «Математика»

Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки. 2016. T. 20, № 1. С. 109-120

ISSN: 2310-7081 (online), 1991-8615 (print) doi: http://dx.doi.org/10.14498/vsgtu1449

Математическое моделирование, численные методы и комплексы программ

УДК 517.958:536.24

ОБ ОДНОМ МЕТОДЕ РЕШЕНИЯ ЗАДАЧ ТЕПЛООБМЕНА ПРИ ТЕЧЕНИИ ЖИДКОСТЕЙ В ПЛОСКИХ КАНАЛАХ

А. В. Еремин, И. В. Кудинов, В. В. Жуков

Самарский государственный технический университет, Россия, 443100, Самара, ул. Молодогвардейская, 244.

Аннотация

С использованием дополнительных граничных условий в интегральном методе теплового баланса получено высокой точности приближенное аналитическое решение задачи теплообмена для жидкости, движущейся в плоскопараллельном канале при симметричных граничных условиях первого рода. Ввиду бесконечной скорости распространения теплоты, описываемой параболическим уравнением теплообмена, температура в центре канала изменяется тотчас же после приложения граничного условия первого рода. Путём представления этой температуры в виде дополнительной искомой функции, а также использования дополнительных граничных условий, определяемых так, чтобы искомое решение удовлетворяло исходному дифференциальному уравнению в граничных точках, находится приближенное аналитическое решение краевой задачи. Использование интеграла теплового баланса позволяет свести решение дифференциального уравнения в частных производных к интегрированию обыкновенного дифференциального уравнения относительно дополнительной искомой функции, изменяющейся лишь по продольной переменной. Показано, что выполнение исходного уравнения лишь на границах области с увеличением числа приближений приводит к его выполнению и внутри области. Отсутствие необходимости интегрирования дифференциального уравнения по поперечной пространственной переменной, ограничиваясь лишь выполнением интеграла теплового баланса (осредненного исходного дифференциального уравнения), позволяет применять данный метод к краевым задачам, решения которых не могут быть получены с помощью классических аналитических методов.

© 2016 Самарский государственный технический университет. Образец для цитирования

Еремин А. В., Кудинов И. В., Жуков В. В. Об одном методе решения задач теплообмена при течении жидкостей в плоских каналах // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2016. T. 20, № 1. С. 109-120. doi: 10.14498/vsgtu1449. Сведения об авторах

Антон Владимирович Еремин (к.т.н.; [email protected]; автор, ведущий переписку), доцент, каф. теоретических основ теплотехники и гидромеханики.

Игорь Васильевич Кудинов (к.т.н.; [email protected]), доцент, каф. теоретических основ теплотехники и гидромеханики.

Виталий Владимирович Жуков ([email protected]), аспирант, каф. теоретических основ теплотехники и гидромеханики.

Ключевые слова: теплообмен в жидкости, бесконечная скорость распространения теплоты, интегральный метод теплового баланса, приближенное аналитическое решение, дополнительная искомая функция, дополнительные граничные условия, тригонометрические координатные функции.

Введение. В теории теплопроводности известны методы, основанные на определении глубины термического слоя, называемые интегральными методами теплового баланса [1—6]. Их применение весьма эффективно при решении ряда краевых задач, получение аналитических решений которых классическими методами не представляется возможным. Особенностью данной группы методов является использование допущения о конечной скорости распространения теплоты, несмотря на то, что решению подлежит параболическое уравнение теплопроводности, в основе которого лежит бесконечная скорость перемещения температурного импульса. В случае применения данных методов процесс нагрева (охлаждения) формально разделяется на две стадии, первая из которых характеризуется постепенным перемещением фронта температурного возмущения от поверхности тела к его центру. Дополнительной искомой функцией в первой стадии является величина фронта температурного возмущения (глубина погретого слоя). Во второй стадии изменение температуры происходит по всей толщине тела. Здесь в рассмотрение вводится дополнительная искомая функция, характеризующая изменение температуры в центре тела.

Многочисленные исследования [1-3] показали, что с увеличением числа приближений (п) первой стадии процесса время (^сц) достижения фронтом температурного возмущения центра тела уменьшается и в пределе при п ^ те ^ те. Таким образом, противоречие, связанное с допущением конечного значения скорости распространения теплоты при решении параболического уравнения теплопроводности интегральным методом теплового баланса, снимается ввиду того, что получаемое решение подтверждает бесконечную скорость ее распространения.

Отметим, что с увеличением числа приближений первой стадии процесса диапазон времени, в котором она определена, уменьшается, а диапазон второй стадии — возрастает. Следовательно, при большом числе приближений первой стадии процесса из полученного решения могут быть найдены температуры лишь для малых и сверхмалых значений временной и пространственной переменных. Роль второй стадии в определении температурного состояния при этом возрастает. В настоящей работе рассматривается метод получения приближенного аналитического решения, позволяющий избежать использования первой стадии процесса в интегральном методе теплового баланса.

Постановка задачи. Последовательность применения предлагаемого метода рассмотрим на примере задачи теплообмена для жидкости, движущейся в плоскопараллельном канале (рис. 1) в следующей математической постановке:

*(£, 0) = ¿о; д^|Й = 0; «то,п) = ¿о,

где £ — температура; £ — поперечная координата; п — продольная координата; а — коэффициент температуропроводности; го — половина ширины плоского канала; ¿0 — начальная температура; ¿с —температура стенки; ш(£) = = 1.5шср(1 — £2/го) — распределение скорости по координате £ (0 ^ £ ^ г0); ^ср — средняя скорость.

о

to -Го

5,

¥

V

to

Рис. 1. Схема стабилизированного ламинарного течения жидкости в плоскопараллельном канале [Figure 1. Scheme of stabilization of a laminar flow in a plane-parallel channel]

При решении задачи (1) приняты следующие допущения: 1) течение жидкости ламинарное стабилизированное, процесс теплообмена стационарный;

температура жидкости на входе в канал постоянна по сечению и равна ¿о;

температура внутренней поверхности стенки трубы во всем диапазоне продольной переменной 0 ^ п < постоянна и равна tc, причем tc = ¿о; теплопроводностью жидкости в продольном направлении, а также теплопроводностью стенки канала пренебрегаем; 5) внутренние источники тепла и диссипация энергии не учитываем [7, 8]. Введем следующие безразмерные переменные и параметры:

в =

t - to £ 8 1 п

-; у = —; х =--—,

tc — t0 r0 3 Pe h

(2)

где в — относительная избыточная температура; у — безразмерная поперечная координата; х — безразмерная продольная координата; Н = 2г0 — ширина канала.

С учетом обозначений (2) задача (1) принимает вид (рис. 2)

(1 — y2) ^в^

дх

_ д2в(у,х) ду2 в(у, 0) =0 д в(0, х) ду

в(1,х) = 1.

(х > 0; 0 < у < 1);

0;

(3)

(4)

(5)

(6)

Решение задачи. Введем дополнительную искомую функцию

q(x) = в(0, х),

(7) 111

У

Рис. 2. Схема теплообмена при течении жидкости в плоскопараллельном канале [Figure 2. Scheme of heat transfer of a fluid flow in a plane-parallel channel]

представляющую изменение температуры в центре канала (рис. 2). Температура в точке у = 0 будет изменяться во времени тотчас же после воздействия теплового удара (граничного условия первого рода) на стенку канала (у = 1), что связано с бесконечной скоростью распространения теплоты, заложенной в параболическом уравнении (3). Следовательно, процесс изменения функции д(х) включает весь диапазон изменения переменной 0 ^ х < те и весь диапазон изменения искомой функции 0 ^ 0(0, х) ^ 1. Так как температура в центре пластины является искомой величиной задачи (3)—(6), введение функции д(х) её не изменяет и рассматривается лишь с целью упрощения процесса получения решения данной задачи.

Решение задачи (3)—(6) ищется в виде

п

0(у, х) = 1 + £ Ък(д)сов(^ЛР), г = 2к - 1, (8)

к= 1

где Ък(д) — неизвестные коэффициенты.

Соотношение (8) удовлетворяет граничным условиям (5), (6). Для определения коэффициентов Ък (д) используются соотношение (7) и дополнительные граничные условия, задаваемые в точках у = 0 и у = 1 и определяемые таким образом, чтобы их выполнение решением (8) было эквивалентно выполнению уравнения (3) в этих точках.

Для получения дополнительных граничных условий применительно к точкам у = 0 и у = 1 продифференцируем граничные условия (5), (6) по переменной х:

d20(Q, x)

dydx d 0(1, x)

dx

(9) (10)

Q

Продифференцируем уравнение (3) по переменной у:

-2уд0(у^ + (1 - у2)= . (п)

дх дхду ду3

Сравнивая соотношения (3) и (9), а также (10) и (11), получаем следующие дополнительные граничные условия:

а30(0х) = 0; (12)

ду3 д20(1,х) ду2

0. (13)

В точке у = 1 уравнение (3) приводится к соотношению (13), которое решением (8) выполняется в любом приближении. Следовательно, уравнение (3) в этой точке всегда выполняется, в связи с чем какие-либо другие дополнительные граничные условия здесь больше не рассматриваются.

Для получения следующих дополнительных граничных условий в точке у = 0 продифференцируем условие (12) по х, а уравнение (3) —трижды по переменной у:

д40(0-х) = 0; (14)

дхду3

бд20(у,х) б д30(у,х) 2 д40(у,х) _ д50(у,х)

дхду у дхду2 у дхду3 ду5 .

Сравнивая (14) и (15), применительно к точке у = 0 с учетом (9) находим

д50(0,х) = 0 ду5

Отсюда можно записать общую формулу для дополнительных граничных условий в точке у = 0:

д*0(0,х) , Л

-^^ = 0, г = 3, 5, 7,.... (16)

дуг у 7

Отметим, что благодаря принятой системе координатных функций условия (16) решением (17) выполняются в любом приближении.

Для нахождения неизвестных коэффициентов Ък (д) используется соотношение (7) и некоторые получаемые из него дополнительные условия. Для их определения продифференцируем соотношение (7) по переменной х:

^д(х) д0(0,х) ^х дх '

Сравнивая (17) с уравнением (3), положив у = 0, находим дополнительное условие вида

¿д(х) = д20(0,х) (18)

^х ду2 '

Для получения следующего дополнительного условия продифференцируем соотношение (18) по переменной х:

^2д(х) = д3в(0,х) ^х2 дхду2

Продифференцируем уравнение (3) дважды по переменной у:

,дв(у,х) _ 4 д2в(у, х) + (1 — 2) д3в(у,х) = д4в(у,х) ' дх дхду дхду2 ду4 '

Сравнивая (19) и (20) с учетом (17), находим еще одно дополнительное граничное условие:

^2д(х) 2 ^д(х) = д 4в(0,х) ^х2 ^х ду4

Таким образом, дифференцируя каждое предыдущее дополнительное граничное условие по переменной х, а уравнение (3) — многократно - по переменной у и сравнивая получаемые соотношения, можно получить какое угодно число дополнительных граничных условий. Например, следующие три таких условия будут

^3д(х) ^2д(х) = д6в(0,х); ^х3 ^х2 ду6 '

А(х) 44^3?(х) б0^29(х) _ д8в(0, х) _ ^х4 ^х3 ^х2 ду8 '

<*5?(х)А(х) 844¿39(х) _ д10в(0,х) ^х5 ^х4 ^х3 ду10 '

Для нахождения решения задачи (3)-(6) в первом приближении, подставляя (8) (ограничиваясь одним членом ряда) в соотношение (7), относительно неизвестного коэффициента 61(5) получаем алгебраическое линейное уравнение, из решения которого находим 61(5) = д(х) _ 1. Подставляя найденное значение 61(5) в (8), получаем

в(у,х) = 1 + (д(х) _ 1)сов(пу/2). (22)

Потребуем, чтобы соотношение (21) удовлетворяло не уравнению (3), а некоторому осредненному уравнению — интегралу теплового баланса:

/<1 = Л ^ (23)

Подставляя (22) в (23), после определения интегралов относительно неизвестной функции д(х) получаем обыкновенное дифференциальное уравнение

^о(х) п4 . . п4 .

+ — 5(х)--= 0. (24)

^х 32^ ! 32 у 7

Интегрируя уравнение (24), находим

5(ж) = 1 + С ехр [(—п4/32) ж] , (25)

где С — постоянная интегрирования. Подставляя (25) в (22), получаем

в(у, ж) = 1 + С ехр ((—п4/32)ж) еов(пу/2). (26)

Для определения постоянной интегрирования С составим невязку граничного условия (4) и потребуем ортогональности невязки к первой координатной функции:

/ (1 + Сеов(пу/2)) еов(пу/2)^у = 0. (27)

Jo

Вычисляя интеграл (27), относительно константы интегрирования С получаем алгебраическое линейное уравнение. Его решение С = —4/п.

После определения константы интегрирования решение задачи (3)—(6) в первом приближении находится из (26). Результаты расчетов по формуле (26) в сравнении с численным решением приведены на рис. 3. Из их анализа следует, что в диапазоне 0.2 ^ ж ^ те расхождение результатов не превышает 5%. С уменьшением ж расхождение возрастает.

Для повышения точности решения необходимо увеличивать число членов ряда (8). При этом следует вводить дополнительные граничные условия. Так, во втором приближении (ограничиваясь двумя членами ряда (8)), для определения неизвестных коэффициентов Ь\(ц), 62(5) используются соотношения

(7), (18).

1.0

0.8

3 ®

0.4

0.2

0 0.2 0.4 0.6 0.8 1.0

Рис. 3. (online в цвете) Распределение температуры при x = 0.04 (чёрные кривые), x = = 0.2 (синие кривые), x = 0.4 (зелёные кривые), x = 0.8 (красные кривые); метки 1, 2, 3 — номера приближений; 4 — метод прогонки [Figure 3. (color online) Temperature distribution when x = 0.04 (black curves), x = 0.2 (blue curves), x = 0.4 (green curves), x = 0.8 (red curves); The labels 1, 2, 3 are approximate numbers of solutions for boundary-value problem (3)-(6); the label 4 is a solution by sweep method]

Подставляя (8) в (7), (18), получаем

5(х) _ Ь1(х) _ Ь2(х) _ 1 = 0;

-д(х) п2 . 9п2 . . п (28)

+ т ь1(х) + — (х) =

Из решения системы уравнений (28) находим

, . . 1 -а(х) 9 . . 9 , , . 1 -а(х) 1 , . 1 , ,

С учетом (29) соотношение (8) принимает вид

- I -

2Л2^ _ > + 1)со8(3?)- (30)

Подставляя (30) в (23), после определения интегралов относительно неизвестной функции д(х) получаем обыкновенное дифференциальное уравнение:

896+ (1952п2 + 108п4) ^^ + 81п6д(х) _ 81п6 = 0. (31) -х2 -х

Общее решение уравнения (31) имеет вид

д(х) = 1 + С ехр(_30.38233х) + С ехр(_2.86059х),

где С1, С2 — постоянные интегрирования.

Для определения постоянных интегрирования С1, С2 составим невязку начального условия (4) и потребуем ортогональности невязки к первой и второй координатным функциям:

[ в(у, 0)сов(пу/2)-у = 0; [ в(у, 0) сов(3пу/2)-у = 0. (32) 00

После подстановки (30) в (32) и определения интегралов соотношения (32) относительно С1 и С2 приводятся к системе линейных уравнений

( _0.20629С1 + 0.48812С2 + 0.63662 = 0; |0.70707С1 + 0.009959С2 _ 0.21221 = 0,

из решения которой получаем С1 = 0.31661, С2 = —1.17042.

После определения неизвестной функции д(х) и констант интегрирования С1, С2 решение задачи (3)-(6) во втором приближении находится из (30). Результаты расчетов безразмерной температуры по формуле (30) в сравнении с точным решением представлены на рис. 3. Их анализ позволяет заключить, что в диапазоне изменения продольной координаты 0.1 ^ х < те расхождение результатов не превышает 1 %.

На рис. 3 приведены также результаты расчетов температуры в третьем приближении (три члена ряда (8)). Из их анализа следует, что в диапазоне 0.1 ^ х < те полученное решение практически совпадает с численным.

На рис. 4, 5 приведены графики изменения невязок уравнения (3) и граничного условия (4) в первом, втором и третьем приближениях. Их анализ позволяет заключить, что с увеличением числа приближений невязки уменьшаются, что свидетельствует о сходимости рассматриваемого метода решения краевой задачи (3)-(6).

1.0

О f-i

В j

2 g

3 g.

S ®

щ J3 Л

H "S

-1.0

\ у- з 4-Х. ч •• \ ч ч. ч

- - *- " -

ч.......v......... Ч Ч ч s V Ч ■-..Л ""V..... Ч ...... ч ч^ y'/f . • ' ' '>' v /1 ..*-<■" t / t / /

0.2

0.4

0.6

0.8

1.0

У

Рис. 4. Невязка уравнения (3) при x = 0.05: 1, 2, 3 — первое, второе и третье приближения [Figure 4. The residual error of the equation (3) when x = 0.05; The labels 1,2,3 are approximate numbers of solutions for boundary-value problem (3)—(6)] 1.0

с

8:1

g-g 0.5

-a 8

■Э >>

8 iS

oj ^

£ I 0

H

о

-0.5

."У/

7/

/ /

/

..............

0.2

0.4

0.6

0.8

1.0

У

Рис. 5. Невязка граничного условия (4): 1, 2, 3 — первое, второе и третье приближения [Figure 5. The residual error of the boundary condition (4); The labels 1, 2, 3 are approximate numbers of solutions for boundary-value problem (3)—(6)]

Заключение. Рассмотренный метод получения решения, основанный на использовании дополнительной искомой функции в интегральном методе теплового баланса, позволяет свести решение дифференциального уравнения в частных производных к интегрированию обыкновенного дифференциального уравнения. При этом требуется выполнение не исходного дифференциального уравнения, а этого же уравнения, осредненного по поперечной пространственной переменной. Такой метод позволяет находить решения краевых задач со сложными дифференциальными операторами в уравнениях (нелинейные, с переменными коэффициентами и др.), получение решений которых с помощью классических аналитических методов либо затруднительно, либо вообще

не представляется возможным.

Показано, что выполнение дифференциального уравнения на границах области с увеличением числа приближений приводит к его выполнению и внутри нее с точностью, зависящей от числа приближений.

Благодарности. Работа выполнена при поддержке Российского фонда фундаментальных исследований (проект № 16-38-00059-мол_а).

ОКСГО8

Антон Владимирович Еремин: http://orcid.org/0000-0002-2614-6329 Игорь Васильевич Кудинов: http://orcid.org/0000-0002-9422-0367 Виталий Владимирович Жуков: http://orcid.org/0000-0003-4960-1416

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. Кудинов В. А., Стефанюк Е. В. Аналитический метод решения задач теплопроводности на основе введения фронта температурного возмущения и дополнительных граничных условий // Инженерно-физический журнал, 2009. Т. 82, №3. С. 540-558.

2. Стефанюк Е. В., Кудинов В. А. Получение приближенных аналитических решений при рассогласовании начальных и граничных условий в задачах теории теплопроводности// Изв. вузов. Матем., 2010. №4. С. 63-71.

3. Кудинов В. А., Кудинов И. В., Скворцова М. П. Обобщенные функции и дополнительные граничные условия в задачах теплопроводности для многослойных тел // Ж. вычисл. матем. и матем. физ., 2015. Т. 55, №4. С. 669-680. doi: 10.7868/ 80044466915040080.

4. Тимошпольский В. И., Постольник Ю. С., Андрианов Д. Н. Теоретические основы теплофизики и термомеханики в металлургии. Минск: Бел. навука, 2006. 560 с.

5. Федоров Ф. М. Граничный метод решения прикладных задач математической физики. Новосибирск: Наука, 2000. 220 с.

6. Глазунов Ю. Т. Вариационные методы. Москва, Ижевск: Регулярная и хаотическая динамика; Институт компьютерных исследований, 2006. 470 с.

7. Петухов Б. С. Теплообмен и сопротивление при ламинарном течении жидкости в трубах. М.: Энергия, 1967. 412 с.

8. Цой П. В. Системные методы расчета краевых задач тепломассопереноса. М.: МЭИ, 2005. 568 с.

Поступила в редакцию 21/1Х/2015; в окончательном варианте — 18/Х1/2015; принята в печать — 24/Х1/2015.

Vestn. Samar. Gos. Techn. Un-ta. Ser. Fiz.-mat. nauki

[J. Samara State Tech. Univ., Ser. Phys. & Math. Sci.], 2016, vol. 20, no. 1, pp. 109-120

ISSN: 2310-7081 (online), 1991-8615 (print) doi: http://dx.doi.org/10.14498/vsgtu1449

MSC: 35K05, 80A20

A METHOD FOR SOLVING PROBLEMS OF HEAT TRANSFER DURING THE FLOW OF FLUIDS IN A PLANE CHANNEL

A. V. Eremin, I. V. Kudinov, V. V. Zhukov

Samara State Technical University,

244, Molodogvardeyskaya st., Samara, 443100, Russian Federation.

Abstract

Using the integral method of heat-transfer with the additional boundary conditions we obtain the high precision approximate analytical solution of heat-transfer for a fluid, moving in plate-parallel channel with symmetric boundary conditions of the first kind. Because of the infinite speed of heat propagation described by a parabolic equation of heat-conduction, the temperature in the centre of channel would change immediately after the boundary conditions (of the first kind) application. We receive the approximate analytical solution of boundary value problem using the representation of this temperature in the form of additional required function and introducing the additional boundary conditions to satisfy the original differential equation in boundary points by the desired function. Using of the integral of heat balance we reduce the solving of differential equation in partial derivatives to integration of ordinary differential equation with respect to additional required function, that changes depending on longitudinal variable. We note that fulfillment of the original equation at the boundaries of the area with increasing number of approximations leads to the fulfillment of that equation inside the area. No need to integrate the differential equation on the transverse spatial variable, so we are limited only by the implementation of the integral of heat-transfer (averaged original differential equation), that allows to apply this method to boundary value problems, unsolvable using classic analytical methods.

Keywords: heat conduction in fluid, infinite speed of heat propagation, integral method of thermal balance, approximate analytical solution, additional required function, additional boundary conditions, trigonometric coordinate functions.

© 2016 Samara State Technical University. Please cite this article in press as:

Eremin A. V., KudinovI.V., Zhukov V. V. A method for solving problems of heat transfer during the flow of fluids in a plane channel, Vestn. Samar. Gos. Tekhn. Univ., Ser. Fiz.-Mat. Nauki [J. Samara State Tech. Univ., Ser. Phys. & Math. Sci.], 2016, vol. 20, no. 1, pp. 109-120. doi: 10.14498/vsgtu1449. (In Russian) Authors Details:

Anton V. Eremin (Cand. Tech. Sci.; [email protected]; Corresponding Author), Associate Professor, Dept. of Theoretical Fundamentals of Heat-Engineering and Hydromechanics. Igor V. Kudinov (Cand. Tech. Sci.; [email protected]), Associate Professor, Dept. of Theoretical Fundamentals of Heat-Engineering and Hydromechanics.

Vitaliy V. Zhukov ([email protected]), Postgraduate Student, Dept. of Theoretical Fundamentals of Heat-Engineering and Hydromechanics.

EpeMUH А. В., KyguHOB M. B., }KyKOB B. B.

ORCIDs

Anton V. Eremin: http://orcid.org/0000-0002-2614-6329 Igor V. Kudinov: http://orcid.org/0000-0002-9422-0367 Vitaliy V. Zhukov: http://orcid.org/0000-0003-4960-1416

Acknowledgments. This work was supported by the Russian Foundation for Basic Research (project no. 16-38-00059-moLa).

REFERENCES

1. Kudinov V. A., Stefanyuk E. V. Analytical solution method for heat conduction problems based on the introduction of the temperature perturbation front and additional boundary conditions, Journal of Engineering Physics and Thermophysics, vol.82, no. 3, pp. 537-555. doi: 10.1007/s10891-009-0223-8.

2. Stefanyuk E. V., Kudinov V. A. Approximate analytic solution of heat conduction problems with a mismatch between initial and boundary conditions, Russian Math. (Iz. VUZ), 2010, vol.54, no. 4, pp. 55-61. doi: 10.3103/S1066369X10040079.

3. Kudinov V. A., Kudinov I. V., Skvortsova M. P. Generalized functions and additional boundary conditions in heat conduction problems for multilayered bodies, Comput. Math. Math. Phys., 2015, vol.55, no. 4, pp. 666-676. doi: 10.1134/S0965542515040089.

4. Timoshpol'skii V. I., Postol'nik Yu. S., Andrianov D. N. Teoreticheskie osnovy teplofiziki i termomekhaniki v metallurgii [Theoretical Foundations of Thermophysics and Thermomechanics in Metallurgy]. Minsk, Bel. navuka, 2006, 560 pp. (In Russian)

5. Fedorov F. M. Granichnyi metod resheniia prikladnykh zadach matematicheskoi fiziki [A Boundary Method for Solving Applied Problems of Mathematical Physics]. Novosibirsk, Nauka, 2000, 220 pp. (In Russian)

6. Glazunov Yu. T. Variatsionnye metody [Variational Methods]. Moscow, Izhevsk, Reguliarnaia i khaoticheskaia dinamika; Institut komp'iuternykh issledovanii, 2006, 470 pp. (In Russian)

7. Petukhov B. S. Teploobmen i soprotivlenie pri laminarnom techenii zhidkosti v trubakh [Heat transfer and resistance at laminar flow of fluids in pipes]. Moscow, Energiia, 1967, 412 pp. (In Russian)

8. Tsoi P. V. Sistemnye metody rascheta kraevykh zadach teplomassoperenosa [System methods for calculating the heat and mass transfer boundary value problems]. Moscow, Moscow Power Engineering Institute Publ., 2005, 568 pp. (In Russian)

Received 21/IX/2015;

received in revised form 18/XI/2015;

accepted 24/XI/2015.