Full Text

Несмотря на значительные успехи, которые получены в области разработки численных методов решения трехмерных уравнений теплопроводности, по-прежнему актуальной задачей является разработка численных алгоритмов высокого порядка точности с компактным шаблоном, которые были бы удобны для применения при параллельном программировании. Самым распространенным методом решения уравнений теплопроводности является метод конечного объема, или FVM (Finite Volume Method), который часто применяется в специализированных программах. При этом ищутся средние значения неизвестных параметров в контрольных объемах, а для нахождения значения параметра в конкретной точке контрольного объема применяется полиномиальная аппроксимация с использованием некоторого набора соседних ячеек. Для обеспечения пространственной аппроксимации на границе между элементами необходимо вычислять поток. Метод конечного объема в этом случае показывает хорошую сходимость и требует небольших затрат дополнительной памяти при использовании сеток, близких к равномерным. Однако при решении практических задач часто приходится сталкиваться с сетками, в которых присутствуют сильно вытянутые ячейки. На таких сетках могут возникнуть проблемы с вычислением потоковых переменных, что в итоге может привести к отсутствию сходимости алгоритма при использовании как явной, так и неявной дискретизации по времени. При этом повышение порядка точности метода конечного объема неизменно связано с расширением шаблона аппроксимации. Другим подходом к решению уравнений теплопроводности является использование методов конечного элемента. К достоинствам данных методов стоит отнести тот факт, что высокий порядок аппроксимации достигается за счет использования ограниченного шаблона, к тому же снижаются требования к равномерности сетки. На каждом элементе решение ищется в виде проекции на пространство непрерывных базисных функций. Неизвестные коэффициенты разложения определяются из решения системы дискретных уравнений, которая получается при умножении дифференциальных уравнений на пробные функции. В настоящее время популярным методом является метод Галёркина с разрывными базисными функциями, или DGM (Discontinuous Galerkin Method), который хорошо себя зарекомендовал для решения уравнений Навье-Стокса [1-3]. В разрывном методе Галёркина пространства базисных и пробных функций совпадают. С целью обеспечения устойчивости метода потоки на границах элементов приближенно определяются с использованием стабилизационных добавок [4]. 524 Решение трехмерных уравнений теплопроводности . . . Как правило, в качестве пробных (базисных) функций выбирают полиномы. Порядок точности получаемого решения находится в прямой зависимости от порядка используемых базисных полиномов. В отличие от метода конечного объема, повышение порядка точности не требует расширения шаблона. В данной работе разрывный метод Галёркина применяется для решения трехмерных уравнений теплопроводности на неструктурированных разнесенных сетках. Подобно тому, как это сделано в работах [5, 6] для термодинамических переменных и компонент вектора скорости, в данной работе температура и компоненты теплового потока рассматриваются на двойственных сетках. В отличие от традиционного DGM-подхода [7-9], где решение и потоковые переменные рассматриваются на одной сетке, в настоящей работе температура рассматривается на ячейках основной сетки, а потоковые переменные рассматриваются на ячейках двойственной сетки. Такой подход позволяет избежать проблем с выбором потоковой функции на границе элемента и обеспечивает работоспособность метода на сетках с сильно вытянутыми ячейками [10]. Данная статья является продолжением работ [10, 11] и является их обобщением на трехмерный случай. 1. Реализация метода Галёркина с разрывными базисными функциями при решении трехмерных уравнений теплопроводности. Рассмотрим трехмерное уравнение теплопроводности ρCν ∂u = div(κ · ∂t u) + f, (x, y, z) ∈ G, 0 0) по построению и всегда образует правую тройку векторов, поэтому необходимо проверить смешанное произведение векторов t1 , t2 и t4 = (x4 - x3 , y4 - y3 , z4 - z3 ). Если знак выражения (x4 - x3 ) A + (y4 - y3 ) B + (z4 - z3 ) C положителен, то sign t1 · t2 × n = sign t1 · t2 × t4 ; т.е. данная тройка векторов тоже правая и вектор нормали направлен внутрь тетраэдра. Поэто1 3 2 му необходимо переориентировать грань, например Tk Tk Tk , и вектор внешней нормали будет определяться по формуле n = (-A/l, -B/l, -C/l). Следуя работе [13], возьмем три точки на каноническом треугольнике: ˜ t1 : (ξ1 = 2/3, η1 = 1/6) , ˜ t2 : (ξ2 = 1/6, η2 = 2/3) , ˜ t3 : (ξ3 = 1/6, η3 = 1/6) , pω1 = 1/3, pω2 = 1/3, pω3 = 1/3. 1 2 3 Вычислим интеграл по поверхности {Tk Tk Tk }: 1 f (x, y, z) dS = 1-ξ A2 + B 2 + C 2 1 2 3 {Tk Tk Tk } 0 1 ≈ 2 ˜ f (ξ, η) dξdη ≈ 0 3 A2 + B2 + ˜ f (ξi , ηi ) pωi , C2 i=1 ˜ где f - значение подынтегральной функции в образах квадратурных точек в исходном треугольнике. 3. Вычисление объемного интеграла по тетраэдру Tk и ячейке Dk . Рас1 2 3 4 смотрим тетраэдр Tk с вершинами в точках Tk , Tk , Tk , Tk . Произведем замену переменных так, чтобы тетраэдр Tk отображался на канонический тетраэдр на плоскости {ξ, η, ζ} с координатами: T1 (1, 0, 0), T2 (0, 1, 0), T3 (0, 0, 1), 528 Решение трехмерных уравнений теплопроводности . . . T4 (0, 0, 0). В случае квадратичного приближения необходимо взять 4 точки в каноническом тетраэдре с координатами с соответствующими весами [13]: ˜ t1 : (ξ1 = α, η1 = β, ζ1 = β) , ˜ t2 : (ξ2 = β, η2 = α, ζ2 = β) , ˜ t3 : (ξ3 = β, η3 = β, ζ3 = α) , ˜ t4 : (ξ4 = β, η4 = β, ζ4 = β) , pω1 = 0.25, pω2 = 0.25, pω3 = 0.25, pω4 = 0.25, где α = 0.58541020, β = 0.13819660. Вычислим объемный интеграл по тетраэдру Tk : f (x, y, z) dxdydz = Tk = |J| T 1 ˜ f (ξ, η, ζ) dξdηdζ ≈ |J| 6 4 ˜ f (ξi , ηi , ζi ) pωi , i=1 ˜ где J - якобиан перехода на канонический тетраэдр, f - значение подынтегральной функции в образах квадратурных точек в исходном тетраэдре. 4. Примеры расчетов. В качестве первого примера рассматривалась следующая задача: ∂u ∂2u ∂2u ∂2u = + 2 + 2 , 0 < x < 1, 0 < y < 1, 0 < z < 1, 0 < t 0.02, ∂t ∂x2 ∂y ∂z u(x, y, z, 0) = sin(πx) sin(πy) sin(πz), 0 x 1, 0 y 1, 0 z 1, u(0, y, z, t) = 0, u(1, y, z, t) = 0, 0 y 1, 0 z 1, 0 t T, u(x, 0, z, t) = 0, u(x, 1, z, t) = 0, 0 x 1, 0 z 1, 0 t T, u(x, y, 0, t) = 0, u(x, y, 1, t) = 0, 0 x 1, 0 y 1, 0 t T, которая имеет точное решение 2 uT = e-3π t sin(πx) sin(πy) sin(πz). Порядки точности определены в норме L2 : N uh - uT L2 1/2 2 (uh - uT ) dV = k=1 , Tk где uh - численное решение задачи на сетке с характеристическим размером ячейки h, N - число ячеек в расчетной области. В табл. 1 приведены порядки точности исследуемого метода (DGM) и метода конечных объемов (FVM) для первой задачи на момент времени T = 0.02 с числом Куранта для тепла, равным 2.5 · 10-5 . В следующих задачах порядки точности исследуемого метода r в норме L2 определены по формуле [14] r = log2 uh - uh/2 L2 , uh/2 - uh/4 L2 529 Ж а л н и н Р. В., Л а д о н к и н а М. Е., М а с я г и н В. Ф., Т и ш к и н В. Ф. Таблица 1 N h 48 725 3748 25558 DGM ошибка порядок 0.433 0.226 0.121 0.063 4.212·10-3 3.219·10-4 4.095·10-5 5.369·10-6 - 3.717 2.974 2.936 FVM ошибка порядок 5.085·10-4 3.001·10-4 7.967·10-5 2.261·10-5 - 0.761 1.913 1.817 где uh , uh/2 , uh/4 - численные решения задачи на сетках с характеристическим размером ячеек h, h/2 и h/4 соответственно. Во второй задаче был рассмотрен единичный куб 0 x 1, 0 y 1, 0 z 1, на гранях x = 0 и x = 1 которого задан тепловой поток ω = 1. На остальных четырех гранях происходит свободная теплоотдача с коэффициентом 0.5. Начальное условие имеет вид u(x, y, z, 0) = 0. В табл. 2 приведены порядки сходимости исследуемого метода (DGM) и метода конечных объемов (FVM) для второй задачи на момент времени T = 0.02 с числом Куранта для тепла, равным 2.5 · 10-5 . В третьей задаче определялось температурное поле единичного куба с начальным условием u(x, y, z, 0) = 1 + x2 + y 2 + z 2 . На плоскостях x = 0, y = 0, z = 0 были заданы условия теплоизоляции ω = 0, а на плоскостях x = 1, y = 1, z = 1 фиксировалась температура начального условия. В табл. 3 приведены порядки сходимости исследуемого метода (DGM) и метода конечных объемов (FVM) для третьей задачи на момент времени T = 0.02 с числом Куранта для тепла, равным 2.5 · 10-5 . Таблица 2 Таблица 3 h DGM FVM h DGM FVM 0.433 0.226 3.518 1.619 0.429 1.161 0.433 0.226 1.833 1.639 2.458 1.772 Заключение. Результаты расчетов показывают возможность применения исследуемой методики для решения трехмерных уравнений теплопроводности. Применение описанного метода позволяет получить порядки точности, близкие ко вторым и выше.

About the authors

Ruslan V Zhalnin

Ogarev Mordovia State University

Email: zhrv@hpc.mrsu.ru
68, Bol'shevistskaya st., Saransk, 430005, Russian Federation
(Cand. Phys. & Math. Sci.; zhrv@hpc.mrsu.ru; Corresponding Author), Head of Dept., Dept. of Applied Mathematics, Differential Equations & Theoretical Mechanics

Marina E Ladonkina

Keldysh Institute of Applied Mathematics of Russian Academy of Sciences

Email: ladonkina@imamod.ru
4, Miusskaya pl., Moscow, 125047, Russian Federation
(Cand. Phys. & Math. Sci.; ladonkina@imamod.ru), Senior Researcher

Victor F Masyagin

Ogarev Mordovia State University

Email: vmasyagin@gmail.com
68, Bol'shevistskaya st., Saransk, 430005, Russian Federation
Asistant Lecturer, Dept. of Applied Mathematics, Differential Equations & Theoretical Mechanics

Vladimir F Tishkin

Keldysh Institute of Applied Mathematics of Russian Academy of Sciences

Email: v.f.tishkin@mail.ru
4, Miusskaya pl., Moscow, 125047, Russian Federation
(Dr. Phys. & Math. Sci.; v.f.tishkin@mail.ru), Deputy Director

References

Statistics

Article Metrics

Refbacks

• There are currently no refbacks.