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

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

B. А. Кочерыжкин

C.-Петербургский государственный университет, аспирант, vladimir@kocheryzhkin.com

Рассматривается модифицированный метод сглаженных частиц WCSPH (англ. weakly compressible smoothed particle hydrodynamics) для решения нестационарных задач гидродинамики. «Несжимаемость» моделируется при помощи выбора подходящего уравнения состояния, и, тем самым, замыкается система уравнений движения жидкости. Представлены результаты численного моделирования обрушения столба жидкости и сравнение с экспериментальными данными. Расчеты произведены с использованием технологии параллельных вычислений Nvidia CUDA.

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

2. Метод сглаженных частиц. Наиболее подробно метод описан в работах [3, 4]. Ниже кратко представлены основные положения метода.

В основе SPH-интерполяции лежит следующее выражение:

где £(г) —дельта-функция Дирака, г — радиус-вектор, А(г) —макрохарактеристика.

Далее дельта-функция аппроксимируется так называемой сглаживающей функцией, или функцией ядра Ш:

(1)

где (•) —оператор усреднения, h — радиус сглаживания.

© В. А. Кочерыжкин, 2011

Затем интеграл аппроксимируется на ансамбле частиц:

А(*) ~ ^УзА№(г - Ь) или А(г) ~ X т5'—^(г - М, (3)

3 3 Р3

где пробегает по всем соседним частицам, ш- —масса частицы ], г- —положение частицы ], V- — объем частицы ], р- — плотность частицы и А- — соответствующая физическая величина, определенная в частице г- .

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

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

частицы имеют меньше соседних частиц, а значит содержат ошибки интерполяции (меньшую плотность). В работе [1] предлагается использовать уравнение неразрывности для вычисления плотности, однако оно не приводит к устойчивым результатам. Для лучшего приближения используется нормализованное уравнение (3):

ш3 Ш (г - г3 ,к) ...

р{Г) Е^(г-г^)- ()

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

\М.(г) = г — г^, К) и ДА(г) = г — г^, К). (5)

3 р3 3 р3

3. Обрушение столба жидкости. В изотермическом случае движение вязкой ньютоновской жидкости описывается уравнением неразрывности и уравнением Навье—Стокса. В данной постановке система замыкается уравнением состояния в форме Тэта:

р + рУч = 0;

-УР + р% + рДч;

Р = в((^У- 1

где ро —начальная плотность, V — вектор скорости, р — коэффициент динамической вязкости, В = роСо/7 — модуль объемного сжатия, 7 — показатель адиабаты, со — начальная скорость звука.

Рассматривается задача об обрушении столба жидкости. В начальный момент і = 0 столб вязкой жидкости начинает обрушаться под действием силы тяжести. Для

моделирования жесткой границы используются виртуальные частицы первого типа [3]. Используя для расчетов 1024 частиц, получаем хорошее качественное совпадение результатов с расчетами, полученными в работе [2], однако при большей разрешающей способности, с количеством частиц 16384, можно наблюдать качественно новые картины течения (см. рис. 1). На последних двух рисунках отчетливо видно образование полости, что соответствует картинам, недавно полученным в работе [5].

Рис. 1. Задача об обрушении столба жидкости (количество частиц 16385).

Используемые параметры: ^ = 10-3 кг/(м ■ с), р0 = 1000 кг/м 3, △ £ = 10-4с, к = 6г.

Следует отметить, что благодаря технологии многопоточных вычислений Nvidia СИБЛ, время расчета удается сократить многократно. Так, время вычислений в задаче об обрушении столба жидкости с шагом по времени 10-4 и количеством частиц

Рис. 2. Сравнение с экспериментальными данными.

16384 составляет 2.5 минуты на 2 секунды реального времени, в то время, как модель в работе [5] из 5000 частиц и таким же шагом по времени рассчитывалась 71 минуту.

Для сравнения с экспериментальными данными [6] используется конфигурация из 128 х 256 частиц. Правая перегородка отсутствует. Высота столба жидкости H = 2W. На рис. 2 приведены сравнения с данными эксперимента переднего фронта волны и высоты столба жидкости. Листинг доступен по ссылке http://cmag.googlecode.com.

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

Литература

1. Monaghan J. J. Simulating free surface flows with SPH // J. Comp. Phys. 1994. Vol. 110, N 2. P. 399-406.

2. Афанасьев К.Е., Ильясов А.Е., Макарчук Р. С, Попов А. Ю. Численное моделирование течений жидкости со свободными границами методами SPH и MPS // Вычисл. технологии. 2006. Т. 11. Спецвыпуск: Избранные доклады семинара по численным методам и информационным технологиям Кемеровского государственного университета. С. 26-44.

3. Liu G. R., Liu M. B. Smoothed Particle Hydrodynamics: a Meshfree Particle Method. World Scientific, 2003.

4. Li S., Liu W. K. Meshfree Particle Methods. Berlin: Springer Verlag, 2004.

5. Hughes J. P., Graham D. I. Comparison of incompressible and weakly-compressible SPH models for free-surface water flows // J. Hydr. Res. 2010. Vol. 48 (Extra Issue). P. 105-117.

6. Martin J. C., Moyce W. J. And experimental study of the collapse of liquid columns on a rigid horizontal plane // Philos. Trans. Roy. Soc. London. 1952. Ser. A. Vol. 244. P. 312-324.

Статья поступила в редакцию 22 апреля 2011 г.