Применение методов VOF и Sph для решения задач с развитой свободной поверхностью Текст научной статьи по специальности «Физика»

Известия Института математики и информатики УдГУ

2015. Вып. 2 (46)

УДК 519.63, 532.5

© С. П. Копысов, Л. Е. Тонкое, А. А. Чернова

ПРИМЕНЕНИЕ МЕТОДОВ VOF И SPH ДЛЯ РЕШЕНИЯ ЗАДАЧ С РАЗВИТОЙ СВОБОДНОЙ ПОВЕРХНОСТЬЮ1

В статье рассматриваются методы моделирования течений со свободными поверхностями в рамках эйлерова (метод объема жидкости в ячейке, volume of fluid (VOF)) и лагранжева (метод гидродинамики сглаженных частиц, smoothed particle hydrodynamics (SPH)) описаний. На примере задач о взаимодействии потока несжимаемой жидкости и преграды оцениваются возможности использования методов для описания развитой свободной поверхности и процесса каплеобразования. Численные результаты, полученные с использованием методов VOF и SPH, сопоставляются с экспериментальными данными.

Ключевые слова: математическое моделирование, свободная поверхность, метод гидродинамики сглаженных частиц (SPH), метод объема жидкости в ячейке (VOF).

Введение

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

Наибольшее распространение для моделирования развитой свободной поверхностью получили два подхода. Первый из них основан на явном выделении границы раздела фаз, когда расчетная сетка перестраивается на каждом временном шаге согласно изменению положения и формы поверхности раздела. Такие методы используется, как правило, при моделировании поверхностных волн и малопригодны для описания процессов слияния или дробления подобластей, занятых одной фазой. В моделях многофазных сред используются метод отслеживания фронта (front tracking method) fl], метод объема жидкости в ячейке (VOF) [2], метод частиц в ячейках, метод функции уровня (level of set) [3] и различные их модификации. Данные алгоритмы характеризуются высокой эффективностью расчета движения межфазных границ, но при изменении топологии границ или при обобщении на многомерный случай значительно усложняются [4].

Во втором подходе для моделирования течений жидкости применяется семейство бессеточных методов, в которых частицы используются как численный прием для интегрирования континуальных уравнений сплошной среды. За основу берутся континуальные уравнения сплошной среды, а частицы играют роль дискретных элементов, позволяющих свести уравнения в частных производных к разностной системе обыкновенных дифференциальных уравнений. Особое положение здесь занимает метод гидродинамики сглаженных частиц (smoothed particle hydrodynamics (SPH)) [5].

Выбор данного метода можно обосновать следующими обстоятельствами. Метод допускает решение большого круга рассматриваемых многомерных задач взаимодействия потоков жидкости с твердыми, движущимися или деформируемыми телами (FSI — fluid structure interaction) без каких-либо существенных изменений программной реализации и имеет достаточно простой интерфейс для взаимодействия. Метод очень хорошо распараллеливается и масштабируется на современных гибридных архитектурах. Кроме того, метод предоставляет возможность моделирования развитой свободной поверхности (дробления и коагуляции) жидкости без применения тех или иных специальных алгоритмов.

Работа построена следующим образом. В первом параграфе приводятся формулировка задачи и основные допущения. Во втором параграфе рассматриваются особенности вычисления свободной границы при использовании метода VOF. В третьем параграфе кратко приводится формулировка метода SPH для решения уравнений Навье-Стокса для несжимаемой жидкости. В четвертом параграфе приводится сравнение тестовой задачи с экспериментальными результатами [6].

1 Работа выполнена при финансовой поддержке РФФИ (грант № 14-0Ю0055).

§ 1. Постановка задачи

Рассмотрим задачу о взаимодействии с иод сформируемой преградой потока жидкости, вызванного распадом начального уровня жидкости [6]. Подобные задачи принято называть задачами о разрушении плотины. Решение гидродинами ческой задачи о разрушении плотины, имеющее самостоятельное значение, может быть использовано также при тестировании математических моделей больших волн, образующихся при землетрясениях, оползнях, сильных взрывах, падении метеоритов, быстрой остановке контейнера, частично заполненного жидкостью, захлестывании палубы судна ветровыми волнами и др.

Пусть область О € М3 с внешней границей Г = Г и Г2 (рис. 1) заполнена двухфазной средой и О = иГо и Для определенности будем считать, что в подобласти находится жидкость, ав ^2— газ; Го является границей раздела фаз.

Столб жидкости, имеющий в начальный момент времени форму параллелепипеда, заключен между стенками х = 0и х = 6. В начальный момент времени Ь = 0 перегородка х = Ь мгновенно удаляется и столб жидкости начинает движение под действием силы тяжести. Развивающееся течение ограничено непроницаемыми твердыми стенками Г1, граница Г2 открытая, с параметрами невозмущенной атмосферы. В данной задаче характерной особенностью является процесс разрушения свободной поверхности с интенсивным каплеобразованием при взаимодействия потока жидкости с препятствием. Имеющиеся экспериментальные данные содержат уровни высот столба жидкости на некоторых участках и видеосъемку процесса течения. Проведенное ранее сопоставление экспериментальных и расчетных высот [7] показало корректность результатов, полученных с использованием методов \¥СБРН и УОР. Следует отметить, что при достаточно больших значениях начального перепада уровней часть жидкости отрывается от стенки и поднимается в воздух в виде отдельных струй и капель. Этот процесс при экспериментальном определении высот жидкости не учитывается. Поэтому в данной статье оцениваются возможности корректного определения развитой свободной поверхности рассматриваемыми численными методами с детализацией мгновенных картин течения.

§ 2. Метод объема жидкости в ячейке (УОЕ)

Во многих практически важных случаях поставленная задача может быть сведена к рассмотрению системы из двух неемешивающихея несжимаемых вязких сред, движение каждой из которых описывается уравнениями Навье Стокса и уравнением неразрывности [8]:

где индекс г = {1, 2} — номер среды, р — плотность, и — вектор скорости, р — давление, тг = ^¿(Уи + Уи?) — тензор вязких напряжений, ^ — динамическая вязкость, § — вектор

Го

(.'кого равновесия:

1.С

Рис. 1. Расчетная область

^^ + V • (рии)г = -У'Рг + V • П + У • иг = 0,

(2.1) (2.2)

(т 1 - т2)е = (р1 - р2 + аК)е, и1 = и2,

(2.3)

(2.4)

Где е — единичный нормальный вектор, внешний по отношению к Q1; K — кривизна поверхности Го; и — коэффициент поверхностного натяжения. Вводя в рассмотрение скалярную функцию

íi, X е Q1, a(t, x) = < V 7 [0, x е Q U Го,

имеющую в данном случае смысл объемной концентрации жидкости, и определив функции плотности р = р2 + (р1 — р2)а, давления p = p2 + (p1 — p2)a и вязкости ц = ц2 + (ц1 — ц2)а, перепишем систему (2.1), (2.2) с учетом (2.4), (2.3) следующим образом:

^ + V • (рии) = -Vp + V • т - aKVa + pg, (2.5)

V- u = 0. (2.6)

Система (2.5), (2.6) дополняется уравнением переноса для а

да „

— = и ■ Va. (2.7)

Численное решение системы (2.5)^(2.7) строится применением метода контрольного объема для дискретизации исходных уравнений. Для согласования полей скоростей и давлений в данной работе была использована процедура расщепления операторов PISO [9].

В рамках рассматриваемой задачи на непроницаемых стенках ставятся условия прилипания: u|rx = 0, а значит и все промежуточные значения u*, u**, u*** на Г1 также должны быть нулевыми. Тогда

V(p** — р*)|гх = V(p* — pn)|r, =0.

Уравнение (2.7) решается на каждом шаге по времени после отыскания решения системы (2.5), (2.6). Для этого (2.7) удобно переписать в консервативном виде

да

— + V • (ua) + V • (a(l - a)ur) = 0, (2.8)

dt

где ur = u2 — u1 — скорость относительного движения фаз на границе раздела Го- В исходной модели Го является бесконечно тонкой и ur =0, однако из-за численной диффузии при вычислениях неизбежно появляется переходная область. Следуя методу контрольного объема, запишем аппроксимацию (2.8):

^м'-^м + ^м1 + £ F/K+1} + £ (1_ аТмЩШ = 0, (2.9)

Í&Nm Í&NM

где

|qurS| Va ,„ ...

F(a)=ca шах q, q = |^|. (2.10)

Скалярный коэффициент ca служит для управления искусственным сжатием решения в области разрыва для компенсации эффекта численной диффузии.

Метод VOF является наиболее надежным для решения широкого круга практических задач отслеживания свободной поверхности. В рассматриваемом подходе маркировка ячеек выполняется следующим образом: если a < 0.5, то ячейки содержат газ, если a > 0.5 — то жидкость. Все, что между ними, определяет поверхностные ячейки свободной границы. Проблема маркировки осложняется и тем, что при при моделировании некоторых течений малых объемов (капля, тонкая струя и т.д.), характерный размер которых меньше ячейки, метод VOF не позволяет разрешить поверхность и приводит к появлению ячеек, которые невозможно отнести к жидкости, газу или поверхности.

§ 3. Метод гидродинамики сглаженных частиц (БРН)

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

Пусть в пространстве П фупкция /(г) известна в N точках:

/(г) =/ /(г')Ш(г - г',Н)йг', г, г' € П, (3.1)

где Ш(г — г',Н) — сглаживающее ядро, удовлетворяющее следующим условиям:

Jw(г,Н)йг = 1; Ш(г — г', Н) ^ 0; Ш(г — г',Н) = Ш(г' — г,Н); Ш(г — г',Н) € Ср, (3.2)

где Н — радиус сглаживания. При Н ^ 0 выражение (3.1) сходится к

/(г) = У /(г'Жг — г')^г'. (3.3)

Проводя приближенно замену интеграла на конечную сумму [5], получаем

N

т* " " - (3.4)

ЯГ) = Е ,=1 р( .?)

Для уравнений гидродинамики в объеме жидкости выделим частицу г, имеющую ш^, рi — массу и плотность. Определим некоторую область сглаживания = {г : |г — г | < кН}, тогда аппроксимации функции и ее градиента примут вид

Л=£ V/, = £ (3.5)

П = П(г4) и = Ш(гi — гj,Н).

Для решения рассматриваемой задачи использовалась функция ядра в виде кубического сплайна:

Ш{т,Н) = к1\{2-Г)\ 1 ^ ^ < 2, (3.6)

4 V Ь' ' ^ ^ /I ^

^ > 2 /г ^ ^

где к = 1/(пН3) для трехмерного пространства. Из уравнения (3.5) можно вычислить плотность заменой произвольной функции / на р:

Pi = £ ш, Wij (г — г,-, Н). (3.7) j

г

^ = - Е -• ->*)• (3-8)

БРН-аппроксимация в обозначениях, введенных выше для уравнений движения сплошной сжимаемой жидкости, будет иметь вид

^ = - Е ^ (| +1 + ^ + ^ (3-9)

где Tij — диссипативныц член в форме Монагана [5]:

Uij • rij < 0 J

Uij • rij > О,

(3.10)

где ßij = (huijrij)/(r2j + n2), Cij = (ci + Ci)/2, Pij = (pi + Pj)/2 nj = 0.01h2, Z — коэффициент искусственной вязкости. Анализ записи вязкого члена, предложенного Монаганом [5], позволяет получить связь коэффициента искусственной вязкости с кинематической вязкостью:

V = (Z • h • Cij)/p.

Таким образом, кинематической вязкости (v = 10-6 м2/с) воды соответствуют значения коэффициента Z = 0.02 ^ 0.03.

Новое положение частицы вычисляется по скоростям окружающих ее частиц:

(з-п>

где £ — константа.

Для расчета течений несжимаемой жидкости использовалась слабосжимаемая формулировка метода WCSPH (weakly compressible SPH), предполагающая замыкание системы уравнений (3.7)^(3.11) уравнением состояния Тейта:

ft = (B[(S),-1l' (312)

10, pi < po,

где B = pocO/Y — модуль объемного сжатия; po — плотность жидкости на свободной поверхности (p = 0); Co — константа, имеющая физический смысл скорости звука; y — константа, для рассматриваемых задач y = 7.

В методе гидродинамики сглаженных частиц условие на свободной границе не выставляется, а давление в частицах, в том числе и в тех, которые принадлежат свободной границе, рассчитывается явно из уравнения состояния (3.12). Тем самым условие постоянства давления на свободной границе удовлетворяется приближенно. Для стабилизации описания свободных поверхностей и нефизичных осцилляций плотности и давления жидкости вблизи твердой стенки применялся фильтр Шепарда:

р* = Е pjWiiT- = Е (3-13)

j pj j

с корректировкой функции ядра WFij = Wij/ j(Wijmj/pj)•

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

§ 4. Численное моделирование

В методе гидродинамики сглаженных частиц жидкость моделировалась 2500000 частицами, расположение которых в начальный момент времени задавалось гексаэдрической решеткой. При решении поставленной задачи методом VOF, в открытом программном комплексе OpenFOAM 2.3.0, использовалась ортогональная сетка, состоящая из 1972608 шестигранников. В вычислительных экспериментах для метода WCSPH использовался GPU-ускоритель GeFORCE GTX 680, 4 Гб, а расчеты методом VOF в OpenFOAM выполнялись на 8 ядрах CPU (два 4-ядерных процессора Intel Xeon Е5420, 2.5 GHz). Общее время решения задачи составило 4.25 и 32.45 часа для методов WCSPH и VOF соответственно. Жидкость (вода) определялась

■ГШ

Рис. 2. Обрушение столба жидкости: а) \VCSPH; б) УОР; в) эксперимент [6]

следующими параметрами: начальная плотность ро = 1000 кг/м3, коэффициент искусственной вязкости для \VCSPH ( = 0.021 (V = 0.0007 м2/с), коэффициент кинематической вязкости для УОР V = 1 • 10-6 м2/с. Величина длины сглаживания Н = 0.003 м.

Образованный вследствие распада начального уровня жидкости поток достигает преграды, и в результате взаимодействия жидкости с ней происходит разделение потока (рис. 2, в, £ = 0.64 с). Часть жидкости обтекает преграду и направляется в боковые зазоры, а центральный

фронт при взаимодействии с преградой образует водяной столб (рис. 2), в котором наблюдается распад непрерывного потока жидкости на отдельные капли и струи. Отметим, что водяной столб, полученный методом \¥СБРН (рис. 2, а), согласуется с экспериментальным (рис. 2, в) как качественно, так и по высоте. При этом водяной столб, полученный в рамках УОР-метода (рис. 2, б), не только ниже экспериментального на 22%, но и не содержит отдельных капель и струй.

Другой характерной особенностью течения в момент времени Ь = 0.64 с является формирование тонкого слоя жидкости, включающего и отдельные капли (рис. 2, в). Данный слой имеет веерную ориентацию, направленную по диагонали от водяного столба к левым углам емкости. В рамках проведенных вычислительных экспериментов веерный слой получен только с использованием метода \¥СБРН (рис. 2, а). Разрешить тонкий слой жидкости методом УОР на используемых расчетных сетках не удалось. Так как толщина веерного слоя не превышает 2 мм, для его разрешения требуется не только более подробная расчетная сетка, но и дополнительная реконструкция свободной поверхности. Отметим, что в других исследованиях, например [12], получить методом \¥СБРН веерную структуру потока не удавалось.

Другим характерным моментом времени является Ь = 0.96 с, когда наблюдается не только рост водяного столба, обусловленный продолжающимся взаимодействием потока с преградой, но и заплеск жидкости на левую стенку емкости с последующим образованием второго водяного столба (рис. 2, Ь = 0.96с). Высота водяных столбов превышает высоту емкости как в эксперименте (рис. 2, в), так и в расчетах (рис. 2, а, б). При этом водяные столбы, полученные методом УОР (рис. 2, б), ниже экспериментальных (рис. 2, в). В момент времени Ь = 0.96с также отмечено взаимодействие отраженного от левой стенки потока с набегающим. На рисунке 2 видно формирование перекрестного центрального водяного столба и образование обратного течения вблизи левой стенки емкости.

Рассмотрим характерные особенности потока жидкости в момент времени Ь = 1.24 с. Отраженный поток жидкости достигает преграды и центрального водяного столба, расположенного у ее подножья. Взаимодействие потока с преградой приводит к образованию обратного течения перед преградой (рис. 2, в, Ь = 1.24с). Как видно из рисунка 2, а, б оба используемых численных метода позволяют отследить как формирование вторичного возвратного течения перед преградой, так и взаимодействие отраженного потока с водяным столбом. Однако метод УОР при данном разрешении стенки не позволяет детально описать развитую свободную поверхность. А полученная с использованием метода \¥СБРН мгновенная картина течения (соответствующая времени Ь = 1.24 с) не согласуется с экспериментальной в месте локализации водного столба (рис. 2, а), что является следствием численного замедления потока. Дополнительно проведенные исследования с различным числом частиц показали, что увеличение числа частиц приводит к уменьшению отставания фронта жидкости от экспериментального, что подтверждается и результатами [12].

Дальнейшее замедление потока в \¥СБРН приводит к тому, что в момент времени Ь = 1.76 с гребень волны, качественно соответствующий экспериментальному, фиксируется в середине участка емкости до перегородки х = Ь (рис. 2, а, Ь = 1.76 с), а не вблизи перегородки, как в эксперименте (рис. 2, в). Метод УОР не разрешает гребень волны (рис. 2, б), но позволяет отразить профиль колебаний свободной поверхности, соответствующий экспериментальному.

Момент времени Ь = 1.96 с также является характерным, так как бурный отраженный поток достигает местоположения перегородки х = Ь (рис. 2, в, Ь = 1.96 с). Отмеченное ранее замедление потока в \¥СБРН приводит к тому, что отраженная волна не достигает участка х=Ь

и эволюция развитой свободной поверхности, хорошо согласуются с данными эксперимента (рис. 2, в). Отраженная волна жидкости, полученная УОР-методом (рис. 2, б, Ь = 1.96с), хорошо согласуется с экспериментальной, за исключением разрешения эволюции развитой свободной поверхности.

Таким образом, сравнение методов УОР и \¥СБРН показало, что метод УОР для разрешения развитой свободной поверхности требует существенно более детальную дискретизацию расчетной области по пространству, чем метод БРН. Для более корректного выделения меж-

фазных границ требуются дополнительные затраты на реконструкцию свободной поверхности. Несмотря на то что полученные VOF-методом картины течений качественно соответствуют экспериментальным, говорить о точном представлении свободной поверхности в виде отдельных капель, брызг не имеет смысла. Поэтому применение метода VOF для моделирования процесса каплеобразования, подразумевающего подробное разрешение нестационарной во времени развитой свободной поверхности, ограничивается сверхвысокими вычислительными затратами и является недостаточно эффективным.

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

Список литературы

1. Osher S., Sethian J.A. Fronts propagating with curvature-dependent speed: Algorithms based on Hamilton-Jacobi formulations // Journal of Computational Physics. 1988. Vol. 79. № 1. P. 12-49. DOI: 10.1016/0021-9991 (88)90002-2

2. Hirt C.W., Nichols B.D. Volume of fluid (VOF) method for the dynamics of free boundaries // Journal of Computational Physics. 1981. Vol. 39. №1. P. 201-225. DOI: 10.1016/0021-9991(81)90145-5

3. Sethian J.A. Level set methods and fast marching methods. Cambridge: Cambridge University Press, 1999. 404 p.

4. Тонкое Л.Е. Численное моделирование динамики капли вязкой жидкости методом функции уровня // Вестник Удмуртского университета. Математика. Механика. Компьютерные науки. 2010. Вып. 3. С. 134-140.

5. Monaghan J.J., Thompson M.C., Hourigan К. Simulation of free surface flows with SPH // Journal of Computational Physics. 1994. Vol. 110. №2. P. 399-406. DOI: 10.1006/jcph.l994.1034

6. Kleefsman K.M.T., Fekken G., Veldman A.E.P., Iwanowski В., Buchner B. A Volume-of-Fluid based simulation method for wave impact problems // Journal of Computational Physics. 2005. Vol. 206. № 1. P. 363-393. DOI: 10.1016/j.jcp.2004.12.007

7. Копысов С.П., Тонков Л.Е., Чернова А.А. Численное моделирование взаимодействия с преградой потока несжимаемой жидкости методом SPH // Сеточные методы для краевых задач и приложения: материалы X междунар. конф. Казан, ун-т. Казань, 2014. С. 412-417.

8. Лойцянский Л.Г. Механика жидкости и газа. М.: Наука, 1978. 736 с.

9. Issa R.I. Solution of implicitly discretised fluid flow equations by operator-splitting // Journal of Computational Physics. 1986. Vol. 62. №1. P. 40-65. DOI: 10.1016/0021-9991(86)90099-9

10. Violeau D. Fluid mechanics and the SPH Method. Oxford: Oxford University Press, 2012. 616 p.

11. Macia F., Leo M., Gonzalez L.M., Jose L., Cercos-Pita J.L., Souto-Iglesisas A. A boundary integral SPH formulation: consistency and applications to ISPH and WCSPH // Progress of Theoretical Physics. 2012. Vol. 128 (3). P. 439-462.

12. Crespo A.C., Dominguez J.M., Barreiro A., Gomez-Gesteira M., Rogers B.D. GPUs, a new tool of acceleration in CFD: efficiency and reliability on smoothed particle hydrodynamics methods // PLoS ONE. 2011. Vol. 6. №6. P. 1-13. DOI: 10.1371/journal.pone.0020685

URL: http://journals.plos.org/plosone/article?id=10.1371/journal.pone.0020685

Поступила в редакцию 30.09.2015

Копысов Сергей Петрович, д. ф.-м.н., заведующий лабораторией, Институт механики УрО РАН, 426067, Россия, г. Ижевск, ул. Т. Барамзиной, 34. E-mail: s.kopysov@gmail.com

Тонков Леонид Евгеньевич, к.ф.-м.н., старший научный сотрудник, Институт механики УрО РАН, 426067, Россия, г. Ижевск, ул. Т. Барамзиной, 34. E-mail: tnk@udman.ru

Чернова Алена Алексеевна, к.т.н., научный сотрудник, Институт механики УрО РАН, 426067, Россия, г. Ижевск, ул. Т. Барамзиной, 34. E-mail: alicaaa@gmail.com

S.P. Kopysov, L. E. Tonkov, A. A. Chernova

Application of VOF and SPH to solve problems with a developed free surface

Keywords: mathematical modeling, the free surface, smoothed particle hydrodynamics (SPH), volume of fluid (VOF). MSC: 76D27, 76M25

This article discusses methods for simulation of flows with free surfaces within the Euler (Volume of fluid — VOF ) and Lagrangian (Smoothed Particle Hydrodynamics — SPH) description. On the example of problems of interaction between the flow of an incompressible fluid and obstacles, we evaluate the possibility of using the methods to describe the developed free surface and the process of dropping. Numerical results obtained by the methods VOF and SPH are compared with experimental data.

REFERENCES

1. Osher S., Sethian J.A. Fronts propagating with curvature-dependent speed: Algorithms based on Hamilton-Jacobi formulations, Journal of Computational Physics, 1988, vol. 79, no. 1, pp. 12-49. DOI: 10.1016/0021-9991(88)90002-2

2. Hirt C.W., Nichols B.D. Volume of fluid (VOF) method for the dynamics of free boundaries, Journal of Computational Physics, 1981, vol. 39, no. 1, pp. 201-225. DOI: 10.1016/0021-9991(81)90145-5

3. Sethian J.A. Level set methods and fast marching methods, Cambridge: Cambridge University Press, 1999, 404 p.

4. Tonkov L.E. Computation of viscous drop dynamics with level set method, Vestn. Udmurt. Univ. Mat. Mekh. Komp'yut. Nauki, 2010, no. 3, pp. 134-140 (in Russian).

5. Monaghan J.J., Thompson M.C., Hourigan K. Simulation of free surface flows with SPH, Journal of Computational Physics, 1994, vol. 110, no. 2, pp. 399-406. DOI: 10.1006/jcph.1994.1034

6. Kleefsman K.M.T., Fekken G., Veldman A.E.P., Iwanowski B., Buchner B. A Volume-of-Fluid based simulation method for wave impact problems, Journal of Computational Physics, 2005, vol. 206, no. 1, pp. 363-393. DOI: 10.1016/j.jcp.2004.12.007

7. Kopysov S.P., Tonkov L.E., Chernova A.A. Numerical simulation of the interaction with the barrier of an incompressible fluid by SPH, Mesh methods for boundary-value problems and applications: Proc. of the 10th Int. Conf., Kazan University, Kazan, 2014, pp. 412-417 (in Russian).

8. Loitsyanskii L.G. Mekhanika zhidkosti i gaza (Mechanics of liquid and gas), Moscow: Nauka, 1978, 736 p.

9. Issa R.I. Solution of implicitly discretised fluid flow equations by operator-splitting, Journal of Computational Physics, 1986, vol. 62, no. 1, pp. 40-65. DOI: 10.1016/0021-9991(86)90099-9

10. Violeau D. Fluid mechanics and the SPH method, Oxford: Oxford University Press, 2012, 616 p.

11. Macia F., Leo M., Gonzalez L.M., Jose L., Cercos-Pita J.L., Souto-Iglesisas A. A boundary integral SPH formulation: consistency and applications to ISPH and WCSPH, Progress of Theoretical Physics, 2012, vol. 128, no. 3, pp. 439-462.

12. Crespo A.C., Dominguez J.M., Barreiro A., Gomez-Gesteira M., Rogers B.D. GPUs, a new tool of acceleration in CFD: efficiency and reliability on smoothed particle hydrodynamics methods, PLoS ONE, 2011, vol. 6, no. 6, pp. 1-13. DOI: 10.1371/journal.pone.0020685

http://j ournals.plos.org/plosone/article?id=10.1371/journal.pone.0020685

Received 30.09.2015

Kopysov Sergei Petrovich, Doctor of Physics and Mathematics, Head of Laboratory, Institute of Mechanics, Ural Branch of the Russian Academy of Sciences, ul. T. Baramzinoi, 34, Izhevsk, 426067, Russia. E-mail: s.kopysov@gmail.com

Tonkov Leonid Evgen'evich, Candidate of Physics and Mathematics, Senior Researcher, Institute of Mechanics, Ural Branch of the Russian Academy of Sciences, ul. T. Baramzinoi, 34, Izhevsk, 426067, Russia. E-mail: tnk@udman.ru

Chernova Alena Alekseevna, Candidate of Engineering, Researcher, Institute of Mechanics, Ural Branch of the Russian Academy of Sciences, ul. T. Baramzinoi, 34, Izhevsk, 426067, Russia. E-mail: alicaaa@gmail.com