НАУКА и ОБРАЗОВАНИЕ

УДК. 519.63

МГТУ им. Н.Э. Баумана

grishin@power.bmstu.ru

vl.zenkin@gmail.com

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

Опыт вычислительной работы с методом крупных частиц (МКЧ) показал, что для  повышения устойчивости без перехода к неявному алгоритму можно применять приемы предварительного, уточняющего определения параметров газа для использования на последующем эйлеровом этапе.

В базовом алгоритме МКЧ давление и скорость на границах ячеек для использования на эйлеровом этапе при вычислении промежуточного значения энергии определяются расчетом средних арифметических значений между уже определенными давлениями и скоростями в смежных ячейках на  n-ом  шаге по времени [1].

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

Для простоты рассмотрим одномерную постановку задачи без учета трения и теплообмена с окружающей средой. С учетом допущений этапа Эйлера в МКЧ

                                                                                                    (1)

дифференциальное уравнение полной энергии преобразуется к виду [1]:

.                                                                                                                        (2)

Здесь и ниже: p – давление; ρ – плотность; u – скорость; E – полная энергия; e – внутренняя энергия; k – показатель адиабаты.

Выразив величину полной энергии через внутреннюю, получим:

.                                                                                                     (3)

Так как уравнение количества движения для этапа Эйлера имеет вид [1]:

,                                                                                                                                (4)

уравнение энергии можно преобразовать к виду к виду

.                                                                                                                             (5)

Внутренняя энергия идеального калорически совершенного газа выражается в виде

.                                                                                                                      (6)

Подставив это выражение в формулу (5), получим следующее соотношение  для давления:

.                                                                                                                       (7)

Подобные соотношения для двухмерной постановки представлены и в работе [2].

Используя явное разностное представление, с помощью соотношения (7) можно определить предварительные значения давления на следующем, n+1 временном слое  для  i-1/2  и  i+1/2 границ между ячейками  i-1,  i  и  i+1:

.                                                                               (8)

Здесь Δt – шаг по времени; Δx – шаг по координате. Полученные значения давления предлагается использовать для определения промежуточных величин скорости в ячейке на этапе Эйлера:

.

Затем скорости на границах ячейки определяются как средние арифметические от их значений в соседних ячейках:

                                                                                                                (9)

либо вычисляются с использованием разностного представления уравнения (4) по формулам:

                                                                                 (10)

Теперь можно определить промежуточное значение полной энергии с использованием формулы

,                                                                           (11)

которая заменяет соответствующую формулу эйлерова этапа в базовом МКЧ [1]:

.                                                                           (12)

Отметим, что с точки зрения устойчивости вычислений способы определения промежуточных скоростей по формулам (9) и (10) дает практически одинаковые результаты.

Тестирование предлагаемой модификации МКЧ и базового МКЧ выполнялось на задаче распада разрыва в трубе с закрытыми концами. Расчетная схема приведена на рис. 1.

Рис. 1. Расчетная схема тестовой задачи и начальные условия.

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

Очевидно, что с целью ускорения расчетов более эффективным будет метод, позволяющий проводить устойчивые расчеты при возможно более высоких значениях сеточного числа Куранта Cо = aΔt/Δx. Как показал численный эксперимент, базовый МКЧ дает уже при  Со = 0.3 чрезмерные осцилляции давления на фронтах волн, а при Со = 0.9 (см. рис. 3) результат расчета становится неприемлемым, происходит срыв расчета. В то же время описанный выше модифицированный МКЧ дает вполне приемлемые результаты, колебания давления имеют существенно меньшую амплитуду, практически не увеличивающуюся с ростом числа Куранта. Устойчивость вычислительной процедуры не нарушается с увеличением шага повремени вплоть до достижения  сеточного  Cо = 1.1.

Рис. 2. Давление в трубопроводе после 90 шагов по времени. Co=0.3.

Рис. 3. Давление в трубопроводе после 30 шагов по времени. Co=0.9.

Таким образом, предлагаемая модификация МКЧ позволяет значительно повысить устойчивость вычислений и уменьшить амплитуды колебаний параметров газа за фронтами волн без радикального изменения расчетного алгоритма и заметного сокращения производительности. Аналогичный подход может быть применен для решения пространственных задач. Метод может быть рекомендован для решения ресурсоемких задач газовой динамики.

СПИСОК ЛИТЕРАТУРЫ

1.         Белоцерковский О.М., Давыдов Ю.М. Метод крупных частиц в газовой динамике.- М.: Наука, 1982, 392 с.

2.         Гришин Ю.А. Новые схемы метода крупных частиц и их использование для оптимизации газовоздушных трактов двигателей // Математическое моделирование. 2002, т. 14, ╧8, с.51-55.