CC BY
93
ТЕХНИЧЕСКИЕ НАУКИ
УДК 519.63 М. В. Васильева, Д. А. Стальное
ЧИСЛЕННОЕ УСРЕДНЕНИЕ ДЛЯ ЗАДАЧИ ТЕПЛОПРОВОДНОСТИ В НЕОДНОРОДНЫХ И ПЕРФОРИРОВАННЫХ СРЕДАХ
Рассматривается нестационарное уравнение теплопроводности в перфорированных и неоднородных средах с высокопроводящими и низкопроводящими включениями. Такие области встречаются при рассмотрении современных строительных материалов, которые обладают существенной неоднородностью, что позволяет улучшить тепловые и механические характеристики. Среди таких материалов можно выделить материалы с пористой (перфорированной) структурой или с неоднородными включениями. Такая неоднородная структура может существенно сказываться на возникающих тепловых потоках. Для численного решения уравнения необходимо построить расчетную сетку, элементы которой могут сеточно разрешить имеющиеся неоднородности. Получаемая таким образом расчетная сетка является мелкомасштабной и может приводить к большому количеству неизвестных. Для численного решения рассматриваемой задачи мы строим аппроксимацию уравнения на крупной сетке с использованием процедуры численного усреднения, которая основана на вычислении эффективных коэффициентов для неоднородной области. Численное сравнение результатов решения модельной задачи проводится для двумерной области для случая перфорированных и неоднородных сред с высокопроводящими и низкопроводящими включениями. В качестве эталонного решения было взято решение методом конечных элементов на мелкой сетке. Вычислительная реализация построена с использованием библиотеки FEniCS. Для построения геометрической области используется программа Gmsh. Визуализация полученных результатов происходит с использованием программы Paraview.
Ключевые слова: уравнение теплопроводности, параболическое уравнение, уравнение с частными производными, уравнение математической физики, перфорированная область, неоднородная область, численные методы, метод численного усреднения, метод конечных элементов, Fenics.
ВАСИЛЬЕВА Мария Васильевна - доцент-исследователь научно-исследовательской кафедры вычислительных технологий ИМИ СВФУ им. М. К. Аммосова.
E-mail: [email protected]
VASILYEVA Maria Vasilyevna - Associate Professor-Researcher of the Research Department of Computational Technologies, Institute of Mathematics and Informatics, M. K. Ammosov North-Eastern Federal University.
СТАЛЬНОВ Денис Алексеевич - магистрант научно-исследовательской кафедры вычислительных технологий ИМИ СВФУ им. М. К. Аммосова.
E-mail: [email protected]
STALNOV Denis Alekseevich - Master's Student of the Research Department of Computational Technologies, Institute of Mathematics and Informatics, M. K. Ammosov North-Eastern Federal University.
M. V. Vasilyeva, D. A. Stalnov
Numerical Averaging for Heat Conduction Problems in Heterogeneous and Perforated Domains
In this paper we consider the time-dependent heat equation in perforated and heterogeneous domains by high- and low- conducting inclusions. Such areas are encountered when considering the modern building materials that have a significant heterogeneity that can significantly improve the thermal and mechanical characteristics. Among such materials we can identify materials with the porous (perforated) structure or with the heterogeneous inclusions. Such heterogeneous structure can significantly affect the resulting heat fluxes. For the numerical solution of the equation is necessary to construct a computational grid, elements of which can allow the existing grid heterogeneity. The computational grid obtained in this way is a small-scale, and can lead to a lot of variables. For the numerical solution of the problem we construct an approximation of the equation on a coarse grid using a numerical averaging procedure, which is based on the calculation of the effective rates for the heterogeneous field. Numerical comparison of the results of the solution of model task carried out for the two-dimensional field for the case of perforated and heterogeneous environments by high- and low- conducting inclusions. The reference solution was taken decision by finite element method on a small grid. Computational implementation built using FEniCS library. To construct a geometric area used Gmsh program. Visualization of the results is done using Paraview program.
Keywords: heat equation, parabolic equation, partial differential equation, the equation of mathematical physics, perforated media, heterogeneous media, numerical methods, numerical averaging method, finite element method, Fenics.
Введение
Современные строительные материалы обладают существенной неоднородностью свойств, которые позволяют улучшить их различные характеристики, например, прочность или теплопроводность. Среди таких материалов можно выделить материалы с пористой (перфорированной) структурой или с неоднородными включениями. При рассмотрении тепловых процессов в таких средах необходимо учитывать неоднородную структуру, которая может существенно сказываться на возникающих тепловых потоках. При учете мелкомасштабных неоднородностей при аппроксимации задачи получаемые расчетные сетки могут приводить к задачам с большим количеством неизвестных [1, 2]. При решении таких задач принято использовать методы численного усреднения [3-6] или различные многомасштабные методы [8-12], которые позволят существенно снизить количество неизвестных посредством построения крупномасштабной модели.
Методы усреднения используются для построения аппроксимации задачи на грубой сетке и позволяют вычислить эффективные характеристики материала. Методы численного усреднения дают макроскопические законы и параметры на основе локальных вычислений, однако такие подходы обычно основываются на априорно сформулированных допущениях. В отличие от этих подходов в многомасштабном методе мы имеем двусторонний обмен информацией между микро- и макромасштабами [11].
В работе рассмотрен процесс теплопереноса в областях с неоднородными коэффициентами и перфорациями. В первой части приводится постановка задачи. Во второй части рассмотрена аппроксимация задачи на мелкой сетке, которая разрешает все неоднородности на сеточном уровне. В третьей части представлен метод численного усреднения, который позволяет вычислять эффективные характеристики и проводить расчет на грубой сетке. Результаты численного исследования представлены в последней части работы.
Постановка задачи
Тепловое состояние твердого тела описывается следующей системой уравнений (закон сохранения энергии и закон Фурье для потока тепла):
с — + q = /, х еО и 0 < I < I
д1
q = — grad Т, х еО,
тах'
(1)
где k - коэффициент теплопроводности, С - коэффициент объемной теплоемкости и ^ - мощность внутренних источников/стоков тепла.
Подстановка выражения для потока тепла в закон сохранения энергии позволяет получить параболическое уравнение теплопроводности [1, 2]
дТ
с-— (к ^ Т ) = |, х ей,0 < I < 1тах.
дг
Дополним уравнение теплопроводности следующим начальным условием
Т (х,0 ) = 0, х еП, t = 0
и граничными условиями
т = g(х), х е Г п,
-к — = 0, х е Г дп
(2)
(3)
(4)
В качестве расчетной области ^ возьмем следующую неоднородную область с неоднородными включениями (рис. 1). Рассмотрим далее два типа неоднородностей в В первом случае мы будем рассматривать случай неоднородной среды О = О! ^ 02 , где коэффициенты уравнения зададим следующим образом:
к =
кх, х ей] к?, х ей
, с =
С], х ей] С2, х е Й^
(5)
Рис. 1. Расчетная область. Слева: неоднородная среда, содержащая включения. Справа: перфорированная область
где индексы 1, 2 соответствуют нашей основной среде и включениям. Во втором случае будем рассматривать перфорированную область, т. е. П = / П2. На границе перфораций зададим однородное граничное условие Неймана
-к 1Т = 0'х е гр' (6)
дп
где дП = Г0 и Гм и Гр.
Аппроксимация задачи на мелкой сетке
Для численного решения задачи на мелком масштабе построим расчетную сетку с треугольными элементами, которая сеточно разрешает все неоднородности. Для распределения тепла в расчетной области О проведем аппроксимацию задачи (2)-(4) с использованием метода конечных элементов. Для задания вариационной постановки задачи умножим уравнение (2) на тестовую функцию q и проинтегрируем по области О
¡а — qdx + |У(А:УГ}qdx = |fqdx, Vq еК, (7)
а д
где
V = { е Я1 (О): V = q, х е Г0= { е Я1 (П): V = 0, х е Г
Используя формулу интегрирования по частям и неявную разностную схему для аппроксимации по времени, получим следующую вариационную постановку задачи: найти Т е V такую, что
_
|с-qdx + , Vq)dx = J/qах, Vq е V. (8)
о т о о
Отметим, что в качестве базисных функций будем использовать стандартные линейные базисные функции. Далее для простоты будем предполагать, что / = 0 .
Метод численного усреднения для решения задачи на крупной сетке
Рассмотрим решение уравнения теплопроводности в неоднородных и перфорированных областях с использованием метода численного усреднения для решения задачи на грубой сетке (рис. 2). Метод заключается в нахождении эффективных коэффициентов на грубой сетке посредством решения независимых локальных задач для каждого участка области (ячейки грубой сетки) [4, 5, 7].
Пусть расчетная область разбивается на треугольные элементы, где Мн и Ми -крупномасштабная и мелкомасштабная расчетные сетки. Для вычисления эффективных коэффициентов, используемых при решении задачи на грубой сетке, необходимо численно решить локальные задачи в подобластях К., где j - номер треугольного элемента грубой сетки, Мн = К^ (j = 1, N, NH - количество элементов грубой сетки).
Запишем аппроксимацию уравнения методом конечных элементов на грубой сетке Мн
Т-Т
л т
п 1 п
*
qdx + К*ЧТ, Vq )dx = 0, Vq е К
(9)
Здесь коэффициенты Си k являются эффективными коэффициентами, которые определены в ячейках К. , у = 1, N.
Для вычисления эффективных характеристик будем решать следующие локальные задачи в подобласти К. для случая двумерной области
Рис. 2. Крупномасштабная расчетная сетка для метода численного усреднения (слева) и мелкомасштабная расчетная сетка (справа)
-div(к grad Тт) = 0, х = (х1, х2) е Kj,
Тт = Хт, х е дKj / Гр, т = 1,2,
дТт
-к-= 0, х е Г п
дп р
для случая перфорированной области и
(к §гай Тт) = 0, х = (х1, х2) е Kj, Тт = хт, х едК., т = 1,2
(10)
(11)
для случая неоднородных коэффициентов k, определенные в (5).
/ *
Эффективный коэффициент теплопроводности k будем вычислять следующим обра-
k * =
кп кх 2 k21 к22
(12)
где компоненты тензора имеют вид
к* — | | I к йх, к-, 2 — | | I к йх ¡К]1'^ дх, I дх
К . г К1 дх2
к*21 — ¿11 Кк5 — I111
к йх.
(13)
КУК) дх2
3
Для вычисления усредненного коэффициента С будем использовать среднее значение по локальной области
*
С =1
ЛI
К *
\КАк
с (х.
(14)
* *
Вычисленные коэффициенты Сил будем использовать при решении задачи (9) на крупномасштабной сетке.
Таким образом, вычислительный алгоритм решения задачи можно представить следующим образом:
1) строим крупномасштабную расчетную сетку Мн ;
# 2) в каждой локальной подобласти К'. вычисляем эффективные характеристики к и С посредством решения задачи (10) или (11), а также выражений (12)-(14);
3) решаем задачу (9) на грубой сетке Мн с использованием вычисленных коэффициен-
Отметим, что в данной конечно-элементной аппроксимации мы используем линейные базисные функции для аппроксимации на мелкой и грубой расчетных сетках.
Численные результаты
Рассмотрим численное решение задачи теплопроводности (2)-(4) в неоднородной среде и в перфорированной области с однородными граничными условиями Неймана (8) на границе перфораций. Геометрические области с мелкомасштабной расчетной сеткой представлены на рис. 3 для случая неоднородных коэффициентов и перфорированной области. Расчетная область является единичным квадратом. Мелкомасштабная расчетная сетка содержит 30529 узлов и 62767 ячеек для неоднородной среды и 29968 узлов и 61309 ячеек для случая перфорированной области (крупномасштабная сетка имеет 16 узлов и 9 ячеек и представлена на рис. 3).
Для проведения численного исследования погрешности рассмотренного метода были использованы следующие параметры:
- коэффициент теплопроводности: к = 0.1 и к2 =1000 ;
- коэффициент теплоемкости: С =1 и С2 = 0.1;
- правая часть / = 0 ;
- моделирование проводилось при £тах = 15 и Т = 0.1;
- начальное условие 7 = 0 ;
- на левой границе задавалось граничное условие Дирихле 7^=1 и на остальных - одно -родные граничные условия Неймана.
Рис. 3. Мелкомасштабные расчетные сетки для неоднородной области (слева) и для перфорированной области (справа). Слева: 30529 узлов и 62767 ячеек. Справа: 29968 узлов и 61309 ячеек
1 = 2 \ = 1 1 = 15
Рис. 4. Распределение температуры в неоднородной среде в различные моменты времени ? = 2,7 и 15. Сверху: решение на мелкой сетке с неоднородными коэффициентами (количество неизвестных N = 30529). Снизу: решение задачи с использованием метода численного усреднения на грубой сетке (количество неизвестных N = 16)
1 = 2 1 = 7 1 = 15
1 = 2 \ = 7 1 = 15
Рис. 5. Распределение температуры в перфорированной области в различные моменты времени 1 = 2,7 и 15. Сверху: решение на мелкой сетке (количество неизвестных N = 29968). Снизу: решение задачи с использованием метода численного усреднения на грубой сетке (количество неизвестных N = 16)
1 = 2 1 = 7 1= 15
1 = 2 1 = 7 1= 15
Рис. 6. Распределение температуры в области с низкопроводящими включениями в различные моменты времени t = 2,7 и 15. Сверху: решение на мелкой сетке с неоднородными коэффициентами (количество неизвестных N = 30529). Снизу: решение задачи с использованием метода численного усреднения на грубой сетке (количество неизвестных N = 16)
Рис. 7. Распределение температуры на конечный момент времени вдоль линий у = 0.5. (слева) и х = 1.0 (справа). Расчеты на мелкой расчетной сетке
Вычислительная реализация построена с использованием библиотеки FEniCS [13, 14]. Для построения геометрической области используется программа Gmsh [15]. Визуализация полученных результатов происходит с использованием программы Paraview [16].
Результаты расчетов на мелкой и грубой сетках представлены на рис. 4-5. Распределение температуры представлено в различные моменты времени для случая неоднородной среды (рис. 4) и для случая перфорированной области (рис. 5). На рисунках сверху изображено рас -пределение температуры при расчете на мелкой сетке, где количество неизвестных системы очень велико, а на рисунках снизу представлено решение на грубой сетке с использованием процедуры численного усреднения. Отметим, что количество неизвестных становится существенно ниже. Результаты расчетов в области с низкопроводящими включениями при ^ = 0.1, k2 = 0.0001 и ех = 1, 02 = 10.0 на мелкой и грубой сетках для неоднородной среды в различные моменты времени представлены на рис. 6.
На рис. 7 представлено распределение температуры на конечный момент времени вдоль
Рис. 8. Относительная Ь. погрешность для неоднородной области (красный), для перфорированной области (зеленый) и для неоднородной области с низкопроводящими включениями (синий) в различные моменты времени
линий у = 0.5 и х = 1.0. Сравнение было проведено на мелкой расчетной сетке и иллюстрирует влияние включений и перфораций на поле температур. Наличие высокопроводя-щих включений повышает температуру на правой границе, а низкопроводящие включения понижают ее.
Распределение Ь2 погрешности на различные моменты времени представлено на рис. 8. Относительная Ь2 погрешность вычислялась следующим образом:
е
\ис - иг
2
dx
4 'Ч |2 1 \иЛ dx
где ис и и^ - решения на грубой и мелкой сетках (эталонное решение). В качестве эталонных решений были взяты решения на мелкой сетке.
Из приведенных результатов можно сделать вывод, что решение задачи теплопроводности методом численного усреднения позволяет существенно уменьшить размерность задачи, при этом погрешность решения очень мала по сравнению с решением на мелкой сетке. Погрешность решения на конечный момент времени для перфорированной области около 2 %, для неоднородной среды около 2 % в случае высокопроводящих включений и около 7 % для случая низкопроводящих включений.
Отметим, что предлагаемый метод дает малую погрешность для случая, когда удается построить грубую расчетную сетку, в которой неоднородности и/или перфорации не касаются границ ячеек грубой сетки. Для случая же, когда неоднородности могут касаться границ ячеек, необходимо строить специальные модификации метода, например, вычисление локальных задач в расширенных областях для предотвращения возникновения высоких градиентов решения на гранях расчетной сетки. Такие методы будут рассмотрены в дальнейших работах для случая двумерных и трехмерных областей.
Заключение
В данной статье было проведено исследование метода численного усреднения для нестационарного уравнения теплопроводности. Рассмотренные неоднородности иллюстрируют строительные материалы с высокопроводящими или низкопроводящими включениями. Проведено сравнение метода численного усреднения для решения задачи на грубой сетке, содержащей 16 узлов, с методом конечных элементов на мелкой сетке, содержащей около 30 000 узлов. Полученные результаты иллюстрируют малую погрешность при моделировании с использованием метода усреднения как для областей с неоднородными свойствами, так и для перфорированной области. Однако было замечено, что решение с малой погрешностью удается получить только в случае отсутствия неоднородностей или перфораций на границе ячеек крупной сетки. В таком случае рекомендуется использовать модификации метода или другой метод, который будет рассмотрен авторами в следующих работах.
Работа выполнена при финансовой поддержке РФФИ (проект №17-01-00732), Правительства РФ (договор №14.Y26.31.0013).
Л и т е р а т у р а
1. Вабищевич П. Н., Самарский А. А. Вычислительная теплопередача. - М.: Едиториал УРСС, 2003.
2. Вабищевич П. Н., Васильева М. В., Горнов В. Ф., Павлова Н. В. Математическое моделирование искусственного замораживания грунтов. Вычислительные технологии. - 2014. - 19 (4), 19-31.
3. Bakhvalov N. S., Panasenko G. Homogenisation: averaging processes in periodic media: mathematical problems in the mechanics of composite materials. - Springer Science \& Business Media, 2012. - Т. 36.
4. Talonov A., Vasilyeva M. On numerical homogenization of shale gas transport // Journal of Computational and Applied Mathematics. - 2016.
5. Hales J. D. et al. Asymptotic expansion homogenization for multiscale nuclear fuel analysis // Computational Materials Science. - 2015. - Т. 99. - С. 290-297.
6. L. J. Durlofsky, Numerical calculation of equivalent grid block permeability tensors for heterogeneous porous media, Water Resour. Res., 27. - 1991, 699-708.
7. X. H. Wu, Y. Efendiev, and T. Y. Hou, Analysis of upscaling absolute permeability, Discrete and Continuous Dynamical Systems, Series B, 2. - 2002, 185-204.
8. Eric T. Chung, Yalchin Efendiev, Guanglian Li, Maria Vasilyeva. Generalized multiscale finite element methods for problems in perforated heterogeneous domains // Applicable Analysis. Volume 95, Issue 10. -2016. - Pp. 2254-2279
9. Chung, E.T., Leung, W.T., Vasilyeva, M. Mixed GMsFEM for second order elliptic problem in perforated domains // Journal of Computational and Applied Mathematics. - 2016. - 304.
- pp. 84-99.
10. G. Allaire and R. Brizzi, A multiscale finite element method for numerical homogenization, SIAM J. Multiscale Modeling and Simulation, 4(3). - 2005, 790-812.
11. Y. Efendiev and T. Hou, Multiscale finite element methods, Theory and applications, Springer, New York, 2009.
12. Федоренко Р. П., Страховская Л. Г Об одном варианте метода конечных элементов // Журнал вычислительной математики и математической физики. - 1979. - Т. 19, №. 4.
- С. 950-960.
13. Logg A., Mardal K. A., Wells G. (ed.). Automated solution of differential equations by the finite element method: The FEniCS book. - Springer Science and Business Media, 2012. - Т. 84.
14. http://fenicsproject.org Вычислительный пакет FEniCS.
15. http://geuz.org/gmsh/. Программный пакет GMSH.
16. http://www.paraview.org Программный пакет Paraview.
R e f e r e n c e s
1. Vabishhevich P. N., Samarskij A. A. Vychislitel'naja teploperedacha. - M.: Editorial URSS, 2003.
2. Vabishhevich P. N., Vasil'eva M. V., Gomov V. F., Pavlova N. V. Matematicheskoe modelirovanie iskusstvennogo zamorazhivanija gruntov. Vychislitel'nye tehnologii. - 2014. 19 (4), 19-31.
3. Bakhvalov N. S., Panasenko G. Homogenisation: averaging processes in periodic media: mathematical problems in the mechanics of composite materials. - Springer Science \& Business Media, 2012. - T. 36.
4. Talonov A., Vasilyeva M. On numerical homogenization of shale gas transport // Journal of Computational and Applied Mathematics. - 2016.
5. Hales J. D. et al. Asymptotic expansion homogenization for multiscale nuclear fuel analysis // Computational Materials Science. - 2015. - T. 99. - S. 290-297.
6. L. J. Durlofsky, Numerical calculation of equivalent grid block permeability tensors for heterogeneous porous media, Water Resour. Res., 27. - 1991, 699-708.
7. X. H. Wu, Y. Efendiev, and T. Y. Hou, Analysis of upscaling absolute permeability, Discrete and Continuous Dynamical Systems, Series B, 2. - 2002, 185-204.
8. Eric T. Chung, Yalchin Efendiev, Guanglian Li, Maria Vasilyeva. Generalized multiscale finite element methods for problems in perforated heterogeneous domains // Applicable Analysis. Volume 95, Issue 10. -2016. - Pp. 2254-2279.
9. Shung, E.T., Leung, W.T., Vasilyeva, M. Mixed GMsFEM for second order elliptic problem in perforated domains // Journal of Computational and Applied Mathematics. - 2016. - 304. - pp. 84-99.
10. G. Allaire and R. Brizzi, A multiscale finite element method for numerical homogenization, SIAM J. Multiscale Modeling and Simulation, 4(3). - 2005, 790-812.
11. Y. Efendiev and T. Hou, Multiscale finite element methods, Theory and applications, Springer, New York, 2009.
12. Fedorenko R. P., Strahovskaja L. G. Ob odnom variante metoda konechnyh jelementov // Zhurnal vychislitel'noj matematiki i matematicheskoj fiziki. - 1979. - T. 19, №. 4. - S. 950-960.
13. Logg A., Mardal K. A., Wells G. (ed.). Automated solution of differential equations by the finite element method: The FEniCS book. - Springer Science and Business Media, 2012. - T. 84.
14. http://fenicsproject.org Software package FEniCS.
15. http://geuz.org/gmsh/. Software package GMSH.
16. http://www.paraview.org Software package Paraview.
МИП ООО «Адгезия-металлоконструкции»
Оказывает полный комплекс услуг по проектированию, производству и комплектации быстровозводимых и энергоэффективных зданий из легких стальных тонкостенных конструкций (ЛСТК), а также комплекс работ по изготовлению термопрофилей.
Телефон: +7 (4112) 25-81-68, +7 (4112) 47-32-86.
E-mail: [email protected].
Сайт: http://www.adgesia.com/.
Адрес: Республика Саха (Якутия), г. Якутск, ул. Чернышевского, 115/1
