CC BY
34
УДК 669.041:536.2.001.573
DOI: 10.14529/power150402
ЧИСЛЕННАЯ АППРОКСИМАЦИЯ КОНВЕКТИВНОГО ГРАНИЧНОГО УСЛОВИЯ ДЛЯ СЕТОК С ПОДВИЖНЫМИ УЗЛАМИ
С.В. Панферов, В.И. Панферов
Южно-Уральский государственный университет, г. Челябинск
Обычно для решения уравнения теплопроводности в областях с переменными во времени границами применяют метод ловли границы в узел пространственной сетки, что обуславливает необходимость использования при расчетах переменного шага по времени, кроме того, переменным будет и число пространственных узлов. Однако во многих случаях более предпочтительным может быть метод сеток с подвижными узлами, в этом случае нет необходимости в изменении числа пространственных узлов и шага по времени. В данной работе для сеток с подвижными узлами рассмотрена задача аппроксимации конвективного граничного условия. Непосредственная замена производных в граничном условии конечными разностями приводит к большой погрешности вычисления температуры поверхности и, вследствие этого, и всего температурного поля тела. При использовании сетки с постоянным шагом по пространству с целью повышения точности расчетов для конечно-разностной замены граничного условия можно использовать формулу Бека. В литературе для сеток с подвижными узлами формулы, аналогичной формуле Бека, нет, поэтому возникает задача по определению такой формулы. Для решения поставленной задачи аппроксимации применен метод теплового баланса для элементарной ячейки у поверхности тела. Выполнена апробация полученной конечно-разностной формулы, в том числе и с помощью вычислительного эксперимента. Полученные результаты могут быть использованы при построении вычислительных схем метода сеток с подвижными узлами.
Ключевые слова: конечно-разностная схема, конвективное граничное условие, метод сеток с подвижными узлами, расчетная область с подвижными границами, температурное поле, аппроксимация.
Введение
При плавлении, затвердевании и окислении металла, при аварийном замерзании воды в теплопроводах и в системах водоснабжения, при расчете процесса охлаждения и замерзания подвижного фронта теплоносителя во время заполнения пустого трубопровода при пуске в зимних условиях, при промерзании влажного грунта, в том числе и при наличии снежного покрова с переменной толщиной и в других случаях возникает задача расчета температурных полей в областях с переменными во времени границами [1-7]. При конечно-разностном решении такой задачи, как правило, применяют метод ловли границы в узел пространственной сетки [7, 8], что обусловливает необходимость использования при расчетах переменного шага по времени, кроме того, переменным будет и число пространственных узлов. Здесь более предпочтительным может быть метод сеток с подвижными узлами [9-11]. Это позволяет, в частности, избежать изменения количества узлов расчетной сетки и в связи с этим и размерности используемых информационных массивов, а также и шага по времени, что может быть достаточно привлекательным, например, при разработке программного обеспечения. Однако при этом следует иметь в виду, что в любом случае - как при численных расчетах с постоянными размерами шагов по пространственным координатам, так и при переменных размерах таких шагов необходима конечно-разностная аппроксимация граничных условий, описывающих особенности теплообмена иссле-
дуемого твердого тела с окружающей средой. При этом известно, что решение этой задачи непосредственной заменой производных конечными разностями может привести к большим погрешностям вычислений. Поэтому необходима разработка специальных подходов и решений этого вопроса.
Актуальность исследуемого вопроса
Согласно данным работ [1, 5-7] актуальность задачи расчета температурных полей в областях с переменными во времени границами для настоящего времени весьма значительна, в литературе отмечается недостаточная развитость и обоснованность некоторых подходов и приемов, в частности, метода сеток с подвижными узлами. Поэтому изучение и выявление всех сторон и особенностей метода сеток с подвижными узлами имеет достаточно большое значение.
Постановка задачи исследования
Чаще всего теплообмен на границе описывается граничным условием III рода, которое имеет вид:
-X
dt_
dN
=a(tc -i|г),
(1)
где t = t (Ы, т) - температура тела в точке Ы в момент времени т ; N - нормаль к границе Г (поверхности тела); X, а - соответственно коэффициенты теплопроводности и теплоотдачи; Г -соответственно температура окружающей среды и
г
температура поверхности тела (температура тела на границе).
Известно, что уравнение (1) обычно аппроксимируется следующей конечно-разностной схе-
мой:
X tn tn-1 _ h
=«(4 - tkn ),
(2)
где tn - температура поверхности тела (в узле п ) в момент времени к -Дх, а - температура тела в соседнем узле п -1, удаленном от поверхностного узла на величину шага по пространству h , в тот же момент времени к -Дх; ^ - температура среды в момент времени к -Дх; Дх - размер расчетного шага по времени. Здесь предполагается, что размер расчетной области по направлению нормали N делится на п частей (шагов h).
Известно также [12], что аппроксимация граничного условия (1) конечно-разностной схемой (2) дает заметную ошибку при определении температуры поверхности тела, если
а h л -< 1.
X
(3)
В этом случае для повышения точности определения температуры поверхности и в связи с этим и всего температурного поля тела при расчете методом сеток с постоянным шагом по пространственной координате h можно использовать формулу аппроксимации, предложенную Беком [13]. Согласно [12] формула Бека прошла успешную апробацию в практике вычислений.
Если же для расчета температурного поля используется метод сеток с подвижными узлами [9-11], то, естественно, возникает вопрос: какой вид будет иметь формула, аналогичная формуле Бека, но для сеток с подвижными узлами. В данной работе дается ответ на этот вопрос.
Теоретическая часть исследования
При выводе формулы аппроксимации в работе [13] использовался достаточно известный прием: для получения разностного решения, хорошо описывающего реальное температурное поле, целесообразно выполнение закона сохранения энергии и для самой разностной схемы [14, 15]. Данный метод часто называют методом конечного контрольного объема [15] или методом теплового баланса для элементарных объемов [8, 14]. Следует заметить, что в отличие от [13] при выводе формулы аппроксимации граничного условия для сеток с подвижными узлами будем использовать усреднение не температур на интервале времени Дх , а плотностей тепловых потоков.
Обозначим плотность теплового потока теплоотдачей в начале интервала времени Дх через q1k,
а в конце - через +1, среднее его значение - че-
рез q1 =( q1 + +1 , причем нетрудно видеть, что
qi =с
tC + tc+1.
Л +1k+1a
ln + tn 2
(4)
Далее, плотность теплового потока теплопроводностью в поверхностном слое тела обозна-
а к
чим соответственно в начале Дх через q2, а в конце - через q2+1, среднее его значение - через =(+ +1 )/2 , причем нетрудно видеть, что
к Л
q2 =x
с -1
k+1 Jt+1 Л
n—1
с -1
n—1
2hk
2h
+1
(5)
Здесь к и к + - расстояние между узлами расчетной сетки в моменты времени - Дх и (к +1) - Дх соответственно.
Следуя [8, 13-15], примем, что размер контрольного объема (элементарной ячейки) для поверхности равен полуслою и оценим теплосодержание (энтальпию) полуслоя на поверхности в начале и в конце временного шага Дх :
ср^кк /2 и соответственно ср^+1 кк+х ¡2 . Здесь с, р - соответственно удельная теплоемкость и плотность вещества, кроме того, также как и в [13] считалось, что средняя температура полуслоя равна температуре поверхности тела.
Разность между количеством теплоты, подведенным к полуслою на поверхности теплоотдачей, и количеством теплоты, отведенным от него за время Дх теплопроводностью во внутрь тела, согласно закону сохранения энергии представляет запасенное количество теплоты, идущее на изменение теплосодержания (энтальпии) полуслоя. Математически это запишется так:
(q1- q2 )-at = ср
+1 +1 tn h
2
tkhk Л 2
(6)
Преобразовывая это уравнение соответствующим образом, получим искомую формулу аппроксимации:
tk+1=-
Xhk a -At
+a t.
(4+tC+1)+x
(tk
.k+1 Л
n-1 1 n-1 hk hk+1
xh
k+1
(7)
a -At
-+а + -
h
k+1
Здесь а - коэффициент температуропроводности.
Заметим также, что в случае, если кк = кк+1 = к , т. е. если граница раздела сред неподвижна, то из (7) вытекает формула аппроксимации Бека [13], которая в данном случае будет иметь вид:
tk
( h2
k+1 __V
ah ---1
a -At X
Л
« (tC +t
k+1
)+ th
+1
k +1
С =
h2 ah ,
-+— + 1
a -At X
(8)
t
a-
k
h
X
Панферов С.В., Панферов В.И.
Численная аппроксимация конвективного граничного условия для сеток с подвижными узлами
Апробация формулы аппроксимации
В работах [12, 13] показано, что предложенная Беком формула аппроксимации конвективного граничного условия совместно с известными методами [7, 14, 16] конечно-разностной замены дифференциального уравнения теплопроводности обеспечивает достаточно точное описание нагрева (охлаждения) массивных тел. Это в известной мере подтверждает ее адекватность реальным физическим процессам.
На адекватность полученной для сеток с подвижными узлами формулы аппроксимации (7), учитывая вышеотмеченное, в определенной мере
указывает тот факт, что при ^ = ^+1 из нее частным случаем получается формула аппроксимации Бека. Кроме того, с целью апробации формулы (7) применялся следующий прием: сравнивались результаты расчета нагрева стальной пластины без учета окисления и с учетом его, но при предположении, что металл сразу после окисления осыпается и не влияет на теплообмен, т. е. окисление приводит только к уменьшению толщины пластины. Такое сравнение обусловлено тем, что в литературе отсутствуют точные данные о распределении температуры в окалине и металле с учетом переноса теплоты через поверхностный слой окисленного металла. Кроме того, натурный эксперимент по определению, например, температуры подвижной границы весьма затруднителен, в частности, из-за того, что в этом случае необходим подвижный датчик температуры.
При выполнении расчетов в первом случае граничное условие (1) аппроксимировалось формулой Бэка, а во втором - формулой (7). Уравнение теплопроводности всех случаях аппроксимировалось неявной разностной схемой, которая решалась методом прогонки. В случае учета окисления металла использовалась схема с подвижными узлами приведенная в [10].
Совершенно ясно, что уменьшение геометрических размеров (толщины) пластины из-за окисления при расчетах должно приводить к лучшему
прогреву металла. Однако температурные распределения в обоих случаях должны незначительно отличаться друг от друга вследствие относительно незначительного уменьшения толщины стальной пластины в результате окисления металла. Это и подтверждается при сравнении результатов расчета.
Представляет также интерес сравнение результатов расчета температурных полей при различных способах аппроксимации конвективного граничного условия для сеток с постоянным шагом по пространству: формулой (2), полученной простой заменой производной конечной разностью, и формулой Бека.
В таблице приведены результаты расчетов симметричного нагрева стальной пластины толщиной 0,1м при а = 0,02м2/ч , Х = 30,24Вт/(м-°С)
, а = 348,9 Вт/(м2 -°С). При этом полагалось, что в начальный момент времени температура во всех точках по толщине пластины одинакова и равна 700 °С (так называемый горячий посад), температура греющей среды tС = 1300 °С , а окисление металла описывается следующим соотношением: дhoк| дт = -39, 4//*ок (т) • ехр{-7580/ [%м (т)+273]} х
х10-4,м/ч, полученным аппроксимацией экспериментальных данных. Здесь ИОК - толщина окалины, а %М - температура поверхности неокис-ленного металла (в реальных условиях под слоем окалины).
В таблице tПМ■ и ^ - значения температур
поверхности и центра нагреваемой пластины при следующих условиях: ■ = 1 - нагрев без учета окисления при аппроксимации граничного условия формулой (2); ■ = 2 - то же самое, но при аппроксимации граничного условия формулой Бека; ■ = 3 - нагрев с учетом окисления, приводящего только к уменьшению толщины пластины (считается, что окисленный металл сразу осыпается), и аппроксимации граничного условия формулой (7).
Температура поверхности и центра пластины
Время, мин Температура,°С Половина толщины неокисленного металла, м
tпм1 tпм2 tпм3 ^2
0 700 700 700 700 700 700 0,05
3 903,6 903,53 904,76 777,73 780,46 783,55 0,04897
6 976,77 967,26 969,27 870,46 867,48 872,20 0,04892
9 1035,45 1024,55 1027,67 948,22 941,22 946,99 0,04884
12 1083,4 1071,54 1075,43 1011,98 1002,45 1008,97 0,04872
15 1122,67 1110,54 1115,02 1064,19 1053,24 1060,29 0,04858
18 1154,81 1142,88 1147,79 1106,94 1095,36 1102,78 0,04842
21 1181,13 1169,70 1174,92 1141,93 1130,29 1137,94 0,04825
24 1202,68 1191,97 1197,36 1170,58 1159,26 1167,02 0,04807
27 1220,32 1210,38 1215,92 1194,04 1183,28 1191,04 0,04788
30 1234,76 1225,68 1231,24 1213,25 1203,20 1210,87 0,04770
Значения половины толщины пластины, приведенные в таблице, относятся только к случаю i = 3 .
Как видно из таблицы, расхождение значений tПМ2 и ?ПМ3, а также и составляет относительно незначительную величину, что позволяет сделать вывод о том, что формула (7) обеспечивает удовлетворительное описание процесса и может быть рекомендована для использования при расчете температурных полей в областях с подвижными границами.
Практическая значимость результатов
Полученная формула аппроксимации конвективного граничного условия для сеток с подвижными узлами, как нам представляется, может быть неким дополнением к теоретическим основам используемого в практике вычислений температурных полей в областях с переменными границами метода сеток с подвижными узлами.
Выводы
Рассмотрена задача аппроксимации конвективного граничного условия для сеток с подвижными узлами. Используя закон сохранения энергии для элементарной ячейки у поверхности тела, получили формулу численной аппроксимации граничного условия, аналогичную известной в литературе формуле Бека. Формула аппроксимации может быть использована для повышения точности расчета температурного поля тела с подвижными границами.
Литература
1. Цаплин, А.И. Моделирование теплофизиче-ских процессов и объектов в металлургии: учеб. пособие / А.И. Цаплин, И.Л. Никулин. - Пермь: Изд-во ПГТУ, 2011. - 299 с.
2. Панферов В.И. К вопросу об оптимальном управлении процессами нагрева (охлаждения) и затвердевания металла / В.И. Панферов // Известия вузов. Черная металлургия. - 1982. - № 4. -С. 129-132.
3. Панферов, В.И. Об оптимальном управлении нагревом окисляющихся массивных тел при теплообмене со средой через поверхностный слой окалины / В.И. Панферов // Известия вузов. Черная металлургия. - 1984. - № 2. - С. 87-90.
4. Панферов, В.И. Идентификация тепловых режимов трубопроводных систем /В.И. Панферов // Вестник ЮУрГУ. Серия «Строительство и архитектура». - 2005. - Вып. 3, № 13 (53). - С. 85-90.
5. Сосновский, А.В. Математическое моделирование влияния толщины снежного покрова на деградацию мерзлоты при потеплении климата / А.В. Сосновский //Криосфера Земли. - 2006. - Т. X, № 3. - С. 83-88.
6. Горелик, Я.Б. Особенности расчета теплового состояния мерзлых грунтов в основании факельной установки / Я.Б. Горелик, С.Н. Романюк,
A.А. Селезнев // Криосфера Земли. - 2014. -Т. XVIII, № 1. - С. 57-64.
7. Кузнецов, Г.В. Разностные методы решения задач теплопроводности: учеб. пособие / Г.В. Кузнецов, М.А. Шеремет - Томск: Изд-во ТПУ, 2007. - 172 с.
8. Арутюнов, В.А. Математическое моделирование тепловой работы промышленных печей /
B.А. Арутюнов, В.В. Бухмиров, С.А. Крупенников. -М.: Металлургия, 1990. - 239 с.
9. Соловьев, А.Е. Решение задачи о движении границы раздела двух сред условия / А.Е. Соловьев, Н.М. Ященко // Инженерно-физический журнал. -
1981. - Т. X, № 2. - С. 370-371.
10. Панферов, В.И. Моделирование нагрева окисляющихся массивных тел методом сеток с «подвижными» узлами / В.И. Панферов, Б.Н. Пар-сункин // Известия вузов. Черная металлургия. -
1982. - № 4. - С. 105-109.
11. Панферов, В.И. Решение задачи Стефана для отключенного теплопровода / В.И. Панферов, Ю. О. Миханькова // Теплофизика и информатика в металлургии: достижения и проблемы: материалы междунар. конф. - Екатеринбург: УГТУ-УПИ, 2000. - С. 284-288.
12. Жеребятьев, И.Ф. Математическое моделирование уравнений типа теплопроводности с разрывными коэффициентами / И. Ф. Жеребятьев, А.Т. Лукьянов. -М.: Энергия, 1968. - 56 с.
13. Бек, Дж. Численная аппроксимация конвективного граничного условия / Дж. Бек // Труды американского общества инженеров-механиков. Теплопередача (русский перевод). - 1962. - № 1. -
C. 109-110.
14. Дульнев, Г.Н. Применение ЭВМ для решения задач теплобмена / Г.Н. Дульнев, В.Г. Парфенов, А.В. Сигалов. -М.: Высш. шк., 1990. - 207 с.
15. Бек, Дж. Некорректные обратные задачи теплопроводности: пер. с англ. / Дж. Бек, Б. Бла-куэлл, Ч. Сент-Клэр, мл. - М.: Мир, 1989. - 312 с.
16. Рябенький, В.С. Введение в вычислительную математику: учеб. пособие /В.С. Рябенький. -М.: Физматлит, 2000. - 296 с.
Панферов Сергей Владимирович, канд. техн. наук, доцент кафедры «Теплогазоснабжение и вентиляция», Южно-Уральский государственный университет, г. Челябинск; [email protected].
Панферов Владимир Иванович, д-р техн. наук, профессор, зав. кафедрой «Теплогазоснабжение и вентиляция», Южно-Уральский государственный университет, г. Челябинск; [email protected].
Поступила в редакцию 23 апреля 2015 г.
Панферов С.В., Панферов В.И.
Численная аппроксимация конвективного граничного условия для сеток с подвижными узлами
DOI: 10.14529/power150402
NUMERICAL APPROXIMATION OF CONVECTIVE BOUNDARY CONDITIONS FOR GRIDS WITH MOBILE NODES
S.V. Panferov, [email protected], V.I. Panferov, [email protected]
South Ural State University, Chelyabinsk, Russian Federation
Usually, to solve the equation of heat conduction in the areas with variable time boundaries are used in the method of catching a boundary node spatial grid, which necessitates the use in the calculations of the step of alternating in time, moreover, be variable, and the number of spatial nodes. However, in many cases, may be more preferable method meshes with mobile nodes, in this case there is no need to change the number of spatial nodes and a time step. In this paper, for meshes with mobile nodes consider the problem of approximating the convective boundary condition. Direct replacement of derivatives in the boundary condition by finite differences leads to large error calculating surface temperature and, therefore, the whole temperature field of the body. When using a grid with a constant pitch in space in order to increase the accuracy of calculations for finite-difference replace the boundary condition formula can be used Beck. In the literature for meshes with mobile nodes formulas similar to Beck, is not, so there is the problem of determining such a formula. To solve the problem of approximation of the method of heat balance of the unit cell in the body surface. Performed testing of the resulting finite-difference formulas, including using a computational experiment. The results obtained can be used in the construction scheme for computing with mobile nodes.
Keywords: finite-difference scheme, convective boundary condition, the grid method with mobile nodes, computational domain with moving boundaries, temperature field approximation.
References
1. Tsaplin A.I., Nikulin I.L. Modelirovanie teplofizicheskikh protsessov i ob"ektov v metallurgii: ucheb. poso-bie [Modeling of Thermal Processes and Objects in Metallurgy: a Tutorial]. Perm', Perm' State Technical University Publ., 2011. 299 p.
2. Panferov V.I. [On the Optimal Control of the Processes of Heating (Cooling) and the Solidification of the Metal]. Izvestiya vuzov. Chernaya metallurgiya [News of High Schools. Ferrous Metallurgy], 1982, no. 4, pp. 129-132. (in Russ.)
3. Panferov V.I. [The Optimal Control of Heating Oxidation of Massive Bodies in Heat Exchange with the Environment Through the Surface Layer of Scale]. Izvestiya vuzov. Chernaya metallurgiya [News of High Schools. Ferrous Metallurgy], 1982, no. 2, pp. 87-90. (in Russ.)
4. Panferov V.I. [Identification of Thermal Conditions of Piping Systems]. Bulletin of the South Ural State University. Ser. Construction Engineering and Architecture, 2005, iss. 3, no. 13 (53), pp. 85-90. (in Russ.)
5. Sosnovskiy A.V. [Mathematical Modeling of the Influence of the Thickness of the Snow Cover on Permafrost Degradation Under Climate Warming]. Kriosfera Zemli [Earths Cryosphere], 2006, vol. X, no. 3, pp. 83-88. (in Russ.)
6. Gorelik Ya.B., Romanyuk S.N., Seleznev A.A. [Mathematical Modeling of the Influence of the Thickness of the Snow Cover on Permafrost Degradation Under Climate Warming]. Kriosfera Zemli [Earths Cryosphere], 2014, vol. XVIII, no 1. pp. 57-64. (in Russ.)
7. Kuznetsov G.V., Sheremet M.A. Raznostnye metody resheniya zadach teploprovodnosti: ucheb. posobie [Difference Methods for Solving Heat Conduction: a Tutorial]. Tomsk, TPU Publ., 2007. 172 p.
8. Arutyunov V.A., Bukhmirov V.V., Krupennikov S.A. Matematicheskoe modelirovanie teplovoy raboty promyshlennykh pechey [Mathematical Modeling of the Thermal Performance of Industrial Furnaces]. Moscow, Metallurgy, 1990. 239 p.
9. Solov'ev A.E., Yashchenko N.M. [Solution of the Problem of the Motion of the Interface Between Two Media Conditions]. Inzhenerno-fizicheskiy zhurnal [Journal of Engineering Physics], 1981, vol. X, no. 2, pp. 370-371. (in Russ.)
10. Panferov V.I., Parsunkin B.N. [Modeling of Heating Oxidation of Massive Bodies with the Method of Nets "mobile" Sites]. Izvestiya vuzov. Chernaya metallurgiya [News of High Schools. Ferrous Metallurgy], 1982, no. 4, pp. 105-109. (in Russ.)
11. Panferov V.I., Mikhan'kova Yu.O.[Solution of the Stefan Problem for a Disconnected Heating Pipeline]. Teplofizika i informatika v metallurgii: dostizheniya i problemy: materialy mezhdunar. konf. Ekaterinburg, UGTU-UPI [Thermal Physics and Computer Science in Metallurgy: Achievements and Challenges: Proceedings
of the International Conference. Ekaterinburg, Ural State Technical University]. Ekaterinburg, Ural State Technical University, 2000, pp. 284-288. (in Russ.)
12. Zherebyat'ev I.F., Luk'yanov A.T. Matematicheskoe modelirovanie uravneniy tipa teploprovodnosti s razryvnymi koeffitsientami [Mathematical Modeling of the Thermal Conductivity Type Equations with Discontinuous Coefficients]. Moscow, Energy, 1968. 56 p.
13. Beck J. [Numerical Approximation of the Convective Boundary Condition]. Trudy amerikanskogo obshche-stva inzhenerov-mekhanikov. Teploperedacha (russkiy perevod) [Proceedings of the American Society of Mechanical Engineers. Heat Transfer (Russian Translation)], 1962, no. 1. pp. 109-110. (in Russ.)
14. Dul'nev G.N., Parfenov V.G., Sigalov A.V. Primenenie EVMdlya resheniya zadach teplobmena [The Use of Computers for Solving Problems of Heat Transfer]. Moscow, Higher School, 1990. 207 p.
15. Beck J., Blackwell B., St. Clair C., Jr. Nekorrektnye obratnye zadachi teploprovodnosti: per. s angl. [Incorrect Inverse Heat Conduction Problems: Translation From English]. Moscow, Higher School, 1990. 207 p.
16. Ryaben'kiy, V.S. Vvedenie v vychislitel'nuyu matematiku: ucheb. posobie [Introduction to Computational Mathematics: a Tutorial]. Moscow, Fizmatlit, 2000. 296 p.
Received 23 April 2015
ОБРАЗЕЦ ЦИТИРОВАНИЯ
FOR CITATION
Панферов, С.В. Численная аппроксимация конвективного граничного условия для сеток с подвижными узлами / С.В. Панферов, В.И. Панферов // Вестник ЮУрГУ. Серия «Энергетика». - 2015. -Т. 15, № 4. - С. 13-18. DOI: 10.14529/power150402
Panferov S.V., Panferov V.I. Numerical Approximation of Convective Boundary Conditions for Grids with Mobile Nodes. Bulletin of the South Ural State University. Ser. Power Engineering, 2015, vol. 15, no. 4, pp. 13-18. (in Russ.) DOI: 10.14529/power150402
