Грант РФФИ 16-07-01037

РФФИ 16-07-01037 «Суперкомпьютерное моделирование динамики жидкости и газа в природных и технических системах на основе лагранжево-эйлерова CSPH-TVD метода высокого порядка точности»

Отчет за 2016 год

  Руководитель проекта
 

Храпов Сергей Сергеевич

  Код и название Конкурса
 

А Конкурс проектов фундаментальных научных исследований

  Аннотация
 

Создан новый лагранжево-эйлеров численный метод (CSPH-TVD) для моделирования течений идеального невязкого газа с произвольными числами Маха в неоднородных и нестационарных внешних полях при наличии в расчетной области сильных разрывов и динамических границ «газ-вакуум». Построено семейство численных схем второго порядка точности, основанных на модифицированных SPH и TVD подходах, и проведено исследование устойчивости, сходимости и эффективности этих схем в 1D, 2D и 3D-моделях. Семейство численных схем CSPH-TVD включает варианты с различными сглаживающими ядрами, TVD-ограничителями и методами решения задачи Римана. Разработанная схема CSPH-TVD является хорошо сбалансированной, что позволяет корректно учитывать как внешние нестационарные силы, так и внутренние, самогравитацию и вязкость газа. 

Проведено детальное сравнение CSPH-TVD-метода с численными схемами SPH, PPM, WENO и MUSCL. Сравнивались устойчивость, точность, сходимость и эффективность численных схем для различных наборов гидродинамических тестов. Результаты численных расчетов с использованием метода CSPH-TVD свидетельствуют о его устойчивости и эффективности при моделировании газодинамических течений с произвольными числами Маха. Численные схемы CSPH-TVD и MUSCL демонстрируют близкое качество решений в областях с большими градиентами давления и в неоднородных гравитационных полях при наличии в расчетной области сильных разрывов и динамических границ типа «газ-вакуум». Преимуществом предложенной схемы CSPH-TVD является хорошая точность моделирования, сравнимая с точностью численных схем годуновского типа при меньших вычислительных затратах. По сравнению с различными модификациями SPH-метода численная схема CSPH-TVD имеет более высокую точность, является менее диссипативной и лучше сбалансированной. Разработанный метод позволяет рассчитывать астрофизические течения в сильно неоднородных гравитационных полях, содержащих нестационарные границы «газ-вакуум». 

CSPH-TVD метод обобщен для интегрирования уравнений Навье-Стокса. Учет вязких сил и теплопроводности осуществляется только на лагранжевом этапе, что позволяет добиться хорошей сбалансированности схемы без существенного снижения производительности вычислений, а также упростить решение задачи Римана на эйлеровом этапе. На основе CSPH-TVD метода созданы последовательные и параллельные OpenMP-CUDA версии программ для двумерных моделей. 

Проведена модификация ранее разработанного CSPH-TVD-метода для интегрирования уравнений мелкой воды: улучшена точность описания гидравлических скачков за счет переменного размера лагранжевых частиц; усовершенствована модель придонного трения. В рамках проекта проведено детальное исследование проблемы граничных условий для уравнений мелкой воды. Модификация CSPH-TVD метода для уравнений мелкой воды позволила получить новые более качественные и точные результаты моделирования динамики затопления территорий в случае зарегулированного стока и чрезвычайных ситуаций (прорыв плотин и ливневые осадки). 

Проведено распараллеливание численных алгоритмов CSPH-TVD метода с использованием современных технологий параллельных вычислений OpenMP и CUDA. На основе модифицированного CSPH-TVD метода для уравнений Сен-Венана созданы CUDA-версия и OpenMP-CUDA-версия программы для численного моделирования гидродинамических течений с учетом реалистичного рельефа земной поверхности. Исследована эффективность распараллеливания численной схемы CSPH-TVD для различных графических процессоров NVIDI TESLA: C2070, K20, K40, K80. 

На основе разработанных подходов при создании численной схемы CSPH-TVD проведена модификация классического SPH метода, которая улучшила свойства численной схемы SPH при описании контактных разрывов и ударных волн. Созданы CUDA-версия и OpenMP-CUDA- версия программ для численного моделирования динамики газовых подсистем астрофизических объектов с учетом самогравитации в рамках трехмерных моделей.

Развернутый научный отчет

  Полученные в 2016 году важнейшие результаты
 

1. Создан новый лагранжево-эйлеров численный метод (CSPH-TVD) для моделирования течений идеального невязкого газа с произвольными числами Маха (M<<1, M<1, M~1, M>1, M>>1) в неоднородных и нестационарных внешних полях при наличии в расчетной области сильных разрывов и динамических границ «газ-вакуум». Построено семейство численных схем второго порядка точности, основанных на модифицированных SPH и TVD подходах, и проведено исследование устойчивости, сходимости и эффективности этих схем в 1D, 2D и 3D-моделях. Семейство численных схем CSPH-TVD включает варианты с различными сглаживающими ядрами (лагранжев этап), TVD-ограничителями и методами решения задачи Римана (эйлеров этап). Разработанная схема CSPH-TVD является хорошо сбалансированной (Well Balanced), что позволяет корректно учитывать как внешние нестационарные силы, так и внутренние, самогравитацию и вязкость газа. Основные результаты анализа численной схемы CSPH-TVD представлены на рис. 1-12 и в табл. 1-4. Часть полученных результатов опубликована в работе [1], а другая их часть подготовлена к публикации.

2. Проведено детальное сравнение CSPH-TVD-метода с численными схемами SPH, PPM, WENO и MUSCL. Сравнивались устойчивость, точность, сходимость и эффективность численных схем для различных наборов гидродинамических тестов. Полученные результаты свидетельствуют о том, что относительные ошибки и порядки сходимости всех рассмотренных численных схем довольно близки. Это может быть обусловлено наличием в структуре решения слабого разрыва (разрыва второго рода), который присутствует в решениях практически всех газодинамических задач. Отметим, что наличие слабого разрыва как раз и позволяет понять свойства численных схем в условиях реальных расчетов, а не синтетических тестов с переносом гладких профилей. Вопрос о точности численного решения в окрестности сильного разрыва требует отдельного рассмотрения в 2017 г. Результаты численных расчетов с использованием метода CSPH-TVD свидетельствуют о его устойчивости и эффективности при моделировании газодинамических течений с произвольными числами Маха. Численные схемы CSPH-TVD и MUSCL демонстрируют близкое качество решений в областях с большими градиентами давления и в неоднородных гравитационных полях при наличии в расчетной области сильных разрывов и динамических границ типа «газ-вакуум». Преимуществом предложенной схемы CSPH-TVD является хорошая точность моделирования, сравнимая с точностью численных схем годуновского типа при меньших вычислительных затратах. По сравнению с различными модификациями SPH-метода численная схема CSPH-TVD имеет более высокую точность, является менее диссипативной и лучше сбалансированной. Разработанный метод позволяет рассчитывать астрофизические течения в сильно неоднородных гравитационных полях, содержащих нестационарные границы «газ-вакуум». Основные результаты анализа численной схемы CSPH-TVD представлены на рис. 1-12 и в табл. 1-4. Часть полученных результатов опубликована в работе [1], а другая их часть подготовлена к публикации.

3. CSPH-TVD метод обобщен для интегрирования уравнений Навье-Стокса (сжимаемый и вязкий газ). Учет вязких сил и теплопроводности осуществляется только на лагранжевом этапе (SPH), что позволяет добиться хорошей сбалансированности схемы без существенного снижения производительности вычислений, а также упростить решение задачи Римана на эйлеровом этапе (TVD). На основе CSPH-TVD метода созданы последовательные и параллельные (OpenMP-CUDA) версии программ для двумерных моделей, которые в настоящее время проходят стадию тестирования и сравнения результатов численных расчетов с известными аналитическими решениями и экспериментальными данными. В работе [8], учитывая позитивные свойства схемы CSPH-TVD, рассмотрена возможность применения данного метода для определения кинематических и термодинамических характеристик расплавленного пластика внутри экструдера параллельного 3D-принтера на основе численного гидродинамического моделирования с учетом вязкости и теплопроводности. Применение численного моделирования для исследования характеристик экструдера является менее затратным и более универсальным по сравнению с прямыми экспериментальными методами. 

4. Проведена модификация ранее разработанного CSPH-TVD-метода для интегрирования уравнений мелкой воды (Сен-Венана): улучшена точность описания гидравлических скачков за счет переменного размера лагранжевых частиц; усовершенствована модель придонного трения. В новой модели придонного трения коэффициент шероховатости зависит от скорости течения, уровня турбулизации потока и глубины. Выявления зависимости коэффициент шероховатости от амплитуды и характерных шкал неоднородностей дна, а также типа дна речного русла или затопляемой поверхности требует дополнительного исследования в 2017 г. на основе сравнения результатов моделирования (2D и 3D) и экспериментальных данных. Важным фактором при моделировании динамики затопления территорий с использованием эйлеровых численных схем является корректная постановка граничных условий. В рамках проекта проведено детальное исследование проблемы граничных условий для уравнений мелкой воды. Основные результаты опубликованы в работе [2]. Модификация CSPH-TVD метода для уравнений мелкой воды позволила получить новые более качественные и точные результаты моделирования динамики затопления территорий в случае зарегулированного стока и чрезвычайных ситуаций (прорыв плотин и ливневые осадки). Основные результаты опубликованы в работах [3, 4, 7].

5. Проведено распараллеливание численных алгоритмов CSPH-TVD метода с использованием современных технологий параллельных вычислений OpenMP и CUDA. На основе модифицированного CSPH-TVD метода для уравнений Сен-Венана созданы CUDA-версия (для одного GPU - получено свидетельство о регистрации [9]) и OpenMP-CUDA-версия (для 2 и 4 GPU - готовится к регистрации) программы для численного моделирования гидродинамических течений с учетом реалистичного рельефа земной поверхности. Исследована эффективность распараллеливания численной схемы CSPH-TVD для различных графических процессоров NVIDI TESLA: C2070, K20, K40, K80. На рис. 13-14 приведены доля времени выполнения различных этапов численной схемы CSPH-TVD на одном временном слое и общее время расчета динамики затопления ВАП с учетом эффективности загрузки различных GPU. Показано, что параллельная реализация численной модели мелкой воды приводит к ускорению расчетов по сравнению с последовательной версией в 65 - 1200 раз в зависимости от типа и количества GPU. Для параллельной CUDA- и OpenMP-CUDA реализации модели мелкой воды на основе СSPH-TVD метода разработан zoom-in (расчет внутри расчета) подход для критически важных регионов затопления. На рис. 15 представлены результаты моделирования динамики затопления ВАП с использованием zoom-in подхода. Основные результаты опубликованы в работах [5, 6].

6. На основе разработанных подходов при создании численной схемы CSPH-TVD проведена модификация классического SPH метода, которая улучшила свойства численной схемы SPH при описании контактных разрывов и ударных волн. Созданы CUDA-версия (для 1 GPU - готовится к регистрации) и OpenMP-CUDA- версия (для 2 и 4 GPU - готовится к регистрации) программ для численного моделирования динамики газовых подсистем астрофизических объектов с учетом самогравитации (прямой расчет гравитационного взаимодействия частиц - каждая с каждой) в рамках трехмерных моделей. Переход на параллельную реализацию для графических процессоров привел к 40-кратному ускорению расчетов по сравнению с параллельной версией для CPU (сравнивались расчеты на 2 x NVIDIA TESLA K80 и 2 x Intel Xeon E5-2687W v3 3.5GHz).

  Сопоставление полученных результатов с мировым уровнем
 

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

Результаты численных расчетов с использованием метода CSPH-TVD свидетельствуют о его устойчивости и эффективности при моделировании газодинамических течений с произвольными числами Маха (M << 1, M < 1, M ~ 1, M > 1, M >> 1). Численные схемы CSPH-TVD и MUSCL демонстрируют близкое качество решений в областях с большими градиентами давления и в неоднородных гравитационных полях при наличии в расчетной области сильных
разрывов и динамических границ типа «газ-вакуум». Преимуществом предложенной схемы CSPH-TVD является хорошая точность моделирования, сравнимая с точностью численных схем годуновского типа при меньших вычислительных затратах. По сравнению с различными модификациями SPH-метода численная схема CSPH-TVD имеет более высокую точность, является менее диссипативной и лучше сбалансированной [Хоперсков А.В. и др. // Астроном. ж., 2012; Khoperskov A.V., Khrapov S.S., Nedugova E.A. // Astronomy Letters, 2003]. Разработанный метод позволяет рассчитывать астрофизические течения в сильно неоднородных гравитационных полях, содержащих нестационарные границы «газ-вакуум» [Еремин М.А. и др. // Изв. Волгогр. гос. техн. ун-та, 2010]. Отметим, что возможны различные варианты реализации алгоритма CSPH-TVD. В частности, можно использовать различные сглаживающие ядра, TVD-ограничители (сохраняющие TVD-свойство разностной схемы), методы решения задачи Римана, а так же возможен учет нестационарных источников/стоков, вязких сил и самогравитации газа (нестационарных внешних сил).
Проведенное исследование показывает, что метод CSPH-TVD имеет сопоставимые с популярными эйлеровыми численными схемами (MUSCL, PPM, WENO) точность и сходимость. Отличительной его чертой является возможность расчетов на границе с вакуумом без необходимости дополнительной регуляризации, что успешно применяется при решении различных газодинамических задач [Хоперсков А.В. Храпов С.С. и др.// Проблемы управления, 2012; Кузьмин Н.М. и др. // Вестник ВолГУ, Серия1, 2015; Khrapov S.S. et al. // Advances in Mechanical Engineering, 2013]. Отметим также, что данный метод может быть применен при решении задач для различных предметных областей, подобных описанным в [Белослудцев А.А. и др.// Вестник ВолГУ, Серия 1, 2009; Кузьмин Н. М. и др. // Астрономический журнал, 2007; Кузьмин Н.М. // Вестник ВолГУ, Серия 1, 2015].

При анализе свойств численной схемы CSPH-TVD проводилось сравнение с мировыми аналогами, представленными в открытой печати.
- ADER [Toro E.F., Titarev V.A. (2005) ADER schemes for scalar non-linear hyperbolic conservation laws with source terms in three-space dimensions, J. Comput. Phys., 202, No. 1, 196-215; Castro C.E., Toro E.F., Kaeser M. (2012) ADER scheme on unstructured meshes for shallow water: simulation of tsunami waves, Geophysics Journal International, 189, 1505-1520]. Отметим также: 
- Метод прыжкового переноса [Головизнин В.М., Карабасов С.А. (2000) Метод прыжкового переноса для численного решения гиперболических уравнений. Точный алгоритм для моделирования конвекции на эйлеровых сетках, Препринт ИБРАЭ № IBRAE-2000-04, 40 с.] примечателен тем, что позволяет получить точное решение задачи о переносе профиля в рамках линейных уравнений адвекции. Но распространение его на систему нелинейных уравнений гидродинамики или мелкой воды далеко нетривиально. Кроме того, сама идеология метода предполагает однородность интегрируемых уравнений, что далеко не всегда соответствует реальным условиям в практических задачах. Отметим также, что метод достаточно чувствителен к значениям числа Куранта-Фридрихса-Леви, а полученные с его помощью результаты существенно зависят от начальных условий [Т.А. Александрикова, М.П. Галанин, Т.Г. Еленина, Нелинейная монотонизация схемы К.И. Бабенко для численного решения уравнения переноса, Матем. моделирование, 16:6 (2004), 44-47].
- Схема «КАБАРЕ» [Головизнин В.М., Самарский А.А. (1998) Разностная аппроксимация конвективного переноса с пространственным расщеплением временной производной, Журн. Мат. Моделирования, 10, № 1.86–100; Карабасов С.А. (2010) О возможностях методов второго порядка аппроксимации на примере модельных задач газо- и гидродинамики, Матем. моделирование, 22, № 7, 93-120; Глотов В.Ю., Головизнин В.М. (2013) Схема КАБАРЕ для двумерной несжимаемой жидкости в переменных «скорость-давление», Журн. Выч. Матем. и Матем. физики. № 6. С.898; Остапенко В.В. (2012) О сильной монотонности схемы “КАБАРЕ”, Журн. Выч. Матем. и Матем. Физики, № 3. С. 447]. Основная идея схемы «кабаре» заключается в расщеплении временной производной. Данная схема, по сути, является суммой схем «явный уголок» и «неявный уголок». Т.е., она представляет собой модификацию сеточного шаблона. Такой подход используется достаточно давно и довольно успешно. В качестве примера можно привести такие методы, как схема Лакса-Вендроффа (и ее монотонизированная TVD-версия), схема Leap-frog (на наш взгляд, в силу схожести основной идеи и сеточного шаблона, являющаяся предтечей схемы «кабаре»). Однако, несмотря на все достоинства подхода, основанного на расщеплении временной производной, он не позволяет построить схему, обладающую свойством сбалансированности (well balanced), без существенных дополнительных модификаций, которые неизбежно достаточно сильно усложнят исходный алгоритм. А указанное свойство сбалансированности численной схемы крайне важно при численном интегрировании уравнений мелкой воды на нерегулярном рельефе.
- Схемы MUSTA [Castro M.J., Pares C., Pardo A., Toro E.F. (2008) Well-balanced high-order MUSTA schemes for nonconservative hyperbolic systems, Numerical Mathematics and Advanced Applications, Proceedings of ENUMATH 2007, the 7th European Conference on Numerical Mathematics and Advanced Applications, Graz, Austria, September 2007, Springer Berlin Heidelberg, Germany (2008)];
- FORCE [Toro E.F., Hidalgo A., Dumbser M. (2009) FORCE schemes on unstructured meshes I: Conservative hyperbolic systems, J. Comput. Phys., 228, 3368–3389].
- ADER-WAF [Titarev V.A., Toro E.F. (2007) Analysis of ADER and ADER-WAF schemes, Journal of Numerical Analysis, 27 616–630]. 

В основном, исследователи в данной области сосредоточены на решении вопросов, касающихся создания новых приближенных решений задачи Римана [Tian B., Toro E.F., Castro C.E. (2011) A path-conservative method for a ?ve-equation model of two-phase ?ow with an HLLC-type Riemann solver, Computers and Fluids, 46, 122-132; Castro C.E., Toro E.F. (2008) Solvers for the high-order Riemann problem for hyperbolic balance laws, J. Comput. Phys., 227, 2481-2513]; создание схем для хоть и очень важных, но частных задач [Colella P., Dorr M.R., Hittinger J.A.F., Martin D.F. (2011) High-Order Finite-Volume Methods in Mapped Coordinates, J. Comput. Phys., 230, No. 8, 2952-2976; Kallemov B., Miller G.H., Trebotich D. (2010) Numerical Simulation of Polymer Flow in Microfluidic Devices, Proceedings of the Fourth SIAM Conference on Mathematics for Industry (MI09), 93–98]; повышении формального порядка точности [Shi J., Zhang Y.-T., Shu C.-W. (2003) Resolution of high order WENO schemes for complicated ?ow structures, J. Comput. Phys., 186, 690–696]. Среди прочих подходов, разрабатываемых для решения систем гиперболических уравнений, отметим работу [Stone J.M., Gardiner T.A. (2008) Athena: A New Code for Astrophysical MHD, Astrophys. J. Suppl., 178, 137], посвященную созданию численного кода для МГД-задач в астрофизике, применение метода Галеркина [Tryoen J., Ern A., Le Maitre O., Ndjinga M. (2010) Intrusive Galerkin methods with upwinding for uncertain nonlinear hyperbolic systems, J. Comput. Phys., 229, No. 18, 6485-6511; Burman E., Quarteroni A., Stamm B. (2010) Interior penalty continuous and discontinuous finite element approximations of hyperbolic equations, J. Sci. Comput., 43, No. 3, 293-312], использование неструктурированных [Ivanov M.S., Bonfiglioli A., Paciorri R., Sabetta F. (2010) Computation of weak steady shock reflections by means of an unstructured shock-fitting solver, Shock Waves, 20, No. 4, 271-284; Lipnikov K., Svyatskiy D., Vassilevski Y. (2010) A monotone finite volume method for advection-diffusion equations on unstructured polygonal meshes, J. Comput. Phys., 229, No. 11, 4017-4032] и адаптивных сеток [Guzik S., McCorquodale P., Colella P. (2012) A Freestream-Preserving High-Order Finite-Volume Method for Mapped Grids with Adaptive-Mesh Refinement, 50th AIAA Aerospace Sciences Meeting Nashville, TN, January 9-12, 2012; Zhang Q. (2011) High-order, multidimensional, and conservative coarse-fine interpolation for adaptive mesh refinement, Comput. Methods Appl. Mech. Engrg., 200(45-46), 3159-3168]. Метод SPH [Monaghan J.J. Particle methods for hydrodynamics // Computer Physics reports. 1985. 3, N2. 71-124]) в гидродинамике получил особое развитие для моделирования астрофизического газа, поскольку наличие многосвязных областей раздела «газ–вакуум» (очень больших перепадов плотности) является типичным для межзвездной и межгалактической среды [Monaghan J.J. Smoothed Particle Hydrodynamics // Ann. Rev. Astron and Astrophysics. 1992. 30. 543-574; Hubber D.A., Batty C.P., McLeod A., Whitworth A.P. SEREN — a new SPH code for star and planet formation simulations. Algorithms and tests // Astronomy & Astrophysics, 2011. 529, A27; Алиев А.В. Применение метода сглаженных частиц для решения задач физической газовой динамики // Вычислительные методы и программирование. 2008. 9. 40-47]. Метод SPH используется и для решения уравнений мелкой воды. Одним из вариантов является GPH-метод (Godunov-type Particle Hydrodynamics), отличающийся использованием решения задачи Римана. Важной явилась работа [Monaghan J.J. Simulating free surface flows with SPH // Journal of Computational Physics. 1994. 110, N2. 399-406], в которой SPH-подход был применен для свободной поверхности несжимаемой жидкости. Авторы работы [Ata R., Soulaimani A. A stabilized SPH method for inviscid shallow water flows // Intern. Journal for Numerical Methods in Fluids. 2005. 47. 139-159] модифицировали метод SPH для модели мелкой воды, также учетом Riemann solver.
Современная вычислительная техника позволяет решать самые различные задачи гидрологии с учетом реального рельефа местности и всех значимых физических факторов для больших территорий [Писарев А.В. и др. // Вестник Удмуртского университета. Математика. Механика. Компьютерные науки. 2013]. Выделим проблему управления гидрологическим режимом пойменных ландшафтов во время весеннего паводка на крупных реках [Воронин А.А. и др. // Проблемы управления. 2012]. Решение этой задачи требует максимально эффективных численных моделей на основе параллельных технологий [Glinskiy B. et al. // Supercomputing frontiers and innovations. 2014]. Для эколого-экономического управления территорией, например, Волго-Ахтубинской поймы, необходимо решать задачу оптимизации гидрографа для конкретных условий в каждом году, и делать экспертные заключения для различных режимов работы десятков гидросооружений в пойме и по строительству новых [Воронин А.А. и др. // Управление большими системами. 2015]. Каждая такая задача требует проведения сотен численных экспериментов по прямому гидродинамическому моделированию весеннего гидрологического режима для территории площадью от 2 тыс. до 20 тыс. кв. км. Практика использования больших суперкомпьютеров (в частности, на вычислительных ресурсах Суперкомпьютерного комплекса МГУ [Sadovnichy V. et al. // Chapman & Hall/CRC Computational Science. Boca Raton, USA, CRC Press. 2013]) для проведения массовых гидродинамических расчетов выявила ряд проблем, связанных с необходимостью проводить за короткий период большое число расчетов и передачей через сеть больших объемов данных для последующей обработки и анализа. Оперативность проведения расчетов и обработки данных моделирования являются важнейшим фактором применения такого рода моделей в практической деятельности. Дополнительной проблемой является визуализация проведенных расчетов, что представляется общей трудностью для очень высокопроизводительных машин [Moreland K. et al. // Supercomputing Frontiers and Innovations. 2015]. Переход на вычислительные ресурсы класса "персональные суперкомпьютеры" на основе GPU с использованием разработанных в рамках проекта пакетов программ для параллельных расчетов на вычислительных узлах NVIDIA TESLA C2070, K20, K40, K80 позволяет эффективно решать указанные проблемы.

  Методы и подходы, использованные в ходе выполнения Проекта
 

В CSPH-TVD-методе расчетная область покрывается неподвижной (эйлеровой) сеткой в центры ячеек которой в начальный момент времени помещаются «жидкие» частицы. Данная схема кардинально отличается от идеологии классического метода частиц в ячейках [Харлоу Ф.Х. М.: Мир, 1967, с.316], а также различных его модификаций [Петренко В.Е., и др. // Числ. методы мех. сплош. среды. 1976. Т.7. № 4. С.130-148; Ким В.В., и др. // Математ. моделир. 2006. Т. 18. № 8. С.5-11].
Основная идея предлагаемого метода заключается в том, чтобы на первом лагранжевом (частицы) этапе учитывать силы (градиент давления, внешние и вязкие) и источники/стоки массы и энергии, а на втором эйлеровом (сетка) этапе рассчитывать потоки массы, импульса и энергии, обусловленные перемещением жидкости через границы ячеек (адвекция).

На первом лагранжевом этапе CSPH-TVD-метода используется SPH подход [Monaghan J.J. // Comp. Phys. reports. 1985. V.3, N2. P.71-124; Monaghan J.J. // Ann. Rev. Astron and Astroph. 1992. V.30. P.543-574], модифицированный нами для обеспечения хорошей сбалансированности и уменьшения численной диссипации алгоритма.
Классический SPH-метод основан на дискретизации сплошной среды конечным набором лагранжевых частиц, движущихся со скоростью потока и допускающих произвольную связность между собой. В данном подходе частицы представляются не в виде точечных объектов, а как «жидкие» частицы, непрерывно распределенные в пространстве и имеющие конечный объем. Пространственное распределение величин, характеризующих жидкие частицы (плотность, энергия), представляется с использованием сглаживающей функции W. В качестве сглаживающего ядра W применяются сплайн-функции различных порядков или Гауссово ядро.

Повышение порядка аппроксимации (по времени и пространству) CSPH-TVD-метода на лагранжевом этапе обеспечивается посредством:
1) применения различных схем интегрирования по времени (Рунге-Кутта, прогноза-коррекции) сохраняющих TVD-свойство метода;
2) увеличения порядка полиномиального представления сглаживающего ядра W.

На втором эйлеровом этапе CSPH-TVD-алгоритма используются сеточные методы удовлетворяющие принципу TVD [Harten A. // J. Comput. Phys. 1983. V.49. N3. P.357-393]. Выполнение условия TVD осуществляется посредством применения функций-ограничителей, не позволяющих антидиффузионным слагаемым разностной схемы становиться слишком большими [Sweby P.K. // SIAM J. Numer. An. 1984. V.21. N5. P.995-1011]. Таким образом, можно строить схемы второго порядка точности и выше, у которых отсутствуют или существенно уменьшены паразитные осцилляции вблизи разрывов.

Поскольку решение задачи Римана (Riemann solver) является важным моментом при построении такого рода схем, то на данном этапе используются приближенные методы (повышающие производительность вычислений) решения задачи Римана (Куранта-Изаксона-Риса, Роу, Ошера, Хартена-Лакса-ван Лира и др.) [Roe P.L. // J. Comput. Phys. 1981. V.43. N.2. P.357–372; Harten A., Lax P.D., van Leer B. (1983) // SIAM Review. 1983. V.255. N.1. P.32-61; Toro E.F., Chakroborthy A. // Aeronaut. J. 1994. V.98. N.979. P.325-332; Toro E.F. // Phylos. Trans. Royal Soc. London. 1989. Ser A341. N.1662. P.499-530]. Отметим также, что для сохранения свойства монотонности или удовлетворения требованию TVD для наклонов распределения параметров в ячейках используются различные ограничители.

Повышение порядка аппроксимации (по времени и пространству) CSPH-TVD-метода на эйлеровом этапе обеспечивается посредством:
1) использования для расчетов потоков усредненных по времени величин, рассчитанных на первом лагранжевом этапе;
2) увеличения порядка полиномиальной реконструкции распределений вычисляемых величин внутри расчетных ячеек с использованием стандартных процедур или с помощью сглаживающего ядра W.
3) использования либо точного решение задачи Римана (GRP), либо приближенного решения, являющегося нелинейным и полным (в решении учитываются ударные волны, волны разряжения, контактные и тангенциальные разрывы).

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

Определение точности разрабатываемых в рамках проекта схем CSPH-TVD-метода осуществляется посредством сравнения численного решения с точным на сжимающихся сетках для различных наборов гидродинамических тестов: распад произвольного разрыва давления, плотности и скорости; длительное распространение профилей разных (гладких и разрывных) плотности; стационарные и медленно распространяющиеся (против набегающего потока) ударные волны; динамика развития гидродинамических неустойчивостей (Кельвина-Гельмгольца, Релея-Тейлора, сверхотражения). На основе этого подхода оцениваются диссипативные и дисперсионные свойства численных схем.

Для оценки скорости сходимости разрабатываемых в рамках проекта схем CSPH-TVD-алгоритма использован метод Рунге [Остапенко В.В. // Вычисл. технол. Т.1997. Т.2, №5. С.57-65; Чирков Д.В., Черный С.Г. // Вычисл. технол. 2000. Т.5. №5. С.86-107; Сафронов А.В. // Вычисл. методы и программ. 2010. Т.11. С.137-143.], основанный на сравнении численных решений получаемых на последовательно сжимающихся сетках для различных наборов гидродинамических тестов.
Важной характеристикой численных схем является эффективность - точность по отношению к вычислительным затратам. Более эффективной будет схема, для которой при заданной точности время вычислений меньше [Titarev V.A., Toro E.F. // J. of Numer. Analysis. 2007. V.27. P.616–630]. Оценка эффективности схем типа CSPH-TVD проводится по отношению к порядку аппроксимации и различным методам решения задачи Римана. Кроме этого, проведено сравнение эффективности CSPH-TVD-метода с другими лагранжевыми и эйлеровыми численными схемами.

Практическое использование численных схем связано с проведением гидродинамических расчетов с высоким разрешением в двумерных и трехмерных моделях, требующих применения технологий параллельных вычислений. Для оценки эффективности распараллеливания численных алгоритмов CSPH-TVD-метода использованы технологии параллельных вычислений OpenMP (системы с общей памятью CPU), и CUDA (графические ускорители NVIDIA GPU) [Антонов А.С. – М.: Изд-во МГУ. 2002; Букатов А.А. и др. – Ростов-на-Дону: Издательство ООО «ЦВВР». 2003; Антонов А.С. – М.: Издательство МГУ. 2009; Gorobets A. at all // Computer & Fluids. 2011. V. 49. №1. P. 101-109; Деги Д.В., Старченко А.В. // Вестник Томского гос. университета. 2012. №2 (18). С.88-98; Горобец А.В. и др. // Вестник Южно-Уральского государственного университета. Серия: Математическое моделирование и программирование. 2011. №25(242). С.76-86; Боресков А.В., Харламов А.А. – М.: ДМК Пресс, 2011. – 232 с.; Сандерс Дж., Кэндрот Э. – М.: ДМК Пресс. 2011. – 232 с.; Казённов А.М. // Комп. исслед. и моделир. 2010. Т.2. №3. С.295-308].

Наиболее перспективным направлением, обеспечивающим высокую производительность при моделировании динамики жидкости и газа, является адаптация вычислительного параллельного кода для гибридных суперкомпьютеров на базе графических процессоров (GPU) NVIDIA, используя технологию CUDA. Архитектура графических процессоров NVIDIA, начиная с Fermi, и вычислители на их основе обеспечили значительно более высокую производительность на операциях с вещественными числами двойной точности (по сравнению с предыдущим поколением), ECC-защиту памяти и привнесли ряд серьезных усовершенствований, позволяющих использовать GPU для широкого круга научных и прикладных задач и промышленных вычислений. Развитием архитектуры Fermi является новая вычислительная архитектура NVIDIA Kepler, которая втрое более экономична.

Для моделирования сложных двумерных и трехмерных гидродинамических расчетов использована двухуровневая гибридная схема распараллеливания численных алгоритмов OpenMP-CUDA [Khrapov S.S., at all. Advances in Mechanical Engineering. 2013; Богданов П.Б., Горобец А.В., Суков С.А. Журнал вычислительной математики и математической физики. 2013]. Двухуровневая гибридная схема OpenMP-CUDA предназначена для работы на гибридных вычислительных системах с общей памятью, в составе которых имеется несколько GPU. Для проведения расчетов на вычислительном кластере, состоящем из нескольких гибридных систем с общей памятью, целесообразно использовать трехуровневую гибридную схему распараллеливания MPI-OpenMP-CUDA. При использовании нескольких GPU необходимо применять технологию NVIDIA GPUDirectTM (Direct Access, Direct Transfers) для увеличения производительности гидродинамических расчетов.

Для детального описания мелкомасштабных структур применяются адаптивные сетки и технология Zoom-in Simulation (расчет внутри расчета). При моделировании в рамках Zoom-in подхода на первом этапе расчет ведется глобально для всей области. На следующем этапе рассчитываются заданные локальные области, с учетом граничных и начальных условий согласованных с предыдущим этапом. В зависимости от конкретной задачи и учета физических механизмов может использоваться несколько таких вложенных расчетов. Данный подход, например, используется в космологических расчетах [например, Semelin, Combes 2005, A&A, 441, 55; Governato et al 2009, MNRAS, 398, 312; Scannapieco et al 2009, MNRAS, 396, 696], имеются примеры использования zoom-in для локальных моделей динамики молекулярной фазы [Haugboelle et al., 2011, LPI Contr. N1639, 9116].

  Вклад каждого члена коллектива в выполнение Проекта в 2016 году
 

Храпов С.С. (руководитель проекта): Общее руководство проектом и организация работ основных исполнителей проекта. Разработка нового лагранжево-эйлерова численного метода CSPH-TVD для моделирования динамики жидкости и газа в природных и технических системах. Построение семейства численных схем различного порядка точности. Разработка численных алгоритмов на основе гибридных технологий OpenMP-CUDA. Подготовка статей для публикации результатов НИР в рецензируемых журналах индексируемых. Подготовка отчетов по проекту.

Хоперсков А.В.: Постановка модельных задач для моделирования динамики жидкости и газа в природных и технических системах. Оценка адекватности и физическая интерпретация полученных результатов численного моделирования. Подготовка статей для публикации результатов НИР в рецензируемых журналах индексируемых.
Кузьмин Н.М.: Компьютерная реализация семейства численных схем CSPH-TVD. Тестирование и сравнение CSPH-TVD-метода с численными схемами PPM, WENO и MUSCL. 

Писарев А.В.: Разработка алгоритмов и комплекса программ для обработки и визуализации результатов гидродинамического моделирования. Проведение гидродинамических расчетов и обработка результатов моделирования.

Агафонникова Е.О.: Разработка программных комплексов для обработки и визуализации результатов гидродинамического моделирования. Проведение гидродинамических расчетов и обработка результатов моделирования.

Дьяконова Т.А.: Разработка CUDA-, OpenMP-CUDA-версий параллельных программ. Исследование эффективности распараллеливания численной схемы CSPH-TVD по сравнению с другими численными алгоритмами. Проведение гидродинамических расчетов и обработка результатов моделирования.

Бочкарева Е.В.: Проведение газодинамических расчетов и обработка результатов моделирования.

Белоусов А.В.: Компьютерная реализация семейства численных схем CSPH-TVD. Тестирование и детальное сравнение CSPH-TVD-метода с другими численными схемами. Тестирование CUDA-, OpenMP-CUDA-версий параллельных программ. Проведение газодинамических расчетов и обработка результатов моделирования.

  Количество научных работ по Проекту за 2016 год
 

Всего: 9

Из них в изданиях, включенных в перечень ВАК: 6

Из них в изданиях, включенных в библиографическую базу данных РИНЦ: 7

Из них в изданиях, включенных в международные системы цитирования: 1

Количество научных работ, подготовленных в ходе выполнения Проекта и принятых к печати в 2016 году: 1

  Участие в 2016 году в научных мероприятиях по тематике Проекта
 

Международная конференция «Суперкомпьютерные дни в России» Москва, 26-27 сентября 2016 г. (устный доклад)

  Библиографический список всех публикаций по Проекту
 

1. Храпов С.С., Кузьмин Н.М., Бутенко М.А. Сравнение точности и сходимости для метода CSPH-TVD и некоторых эйлеровых схем для решения уравнений газодинамики // Вестник Волгоградского государственного университета. Серия 1: Математика. Физика. 2016. № 6 (37). С. 155-163.

2. Дьяконова Т.А., Храпов С.С., Хоперсков А.В. Проблема граничных условий для уравнений мелкой воды // Вестник Удмуртского университета. Математика. Механика. Компьютерные науки. 2016. Т. 26. № 3. С. 401-417.

3. Дьяконова Т.А., Хоперсков А.В., Храпов С.С. Компьютерное моделирование динамики затопления территорий в случае чрезвычайных ситуаций с использованием технологий параллельных вычислений // Кибернетика и программирование. 2016. - №3. - С. 17-34.

4. Агафонникова Е.О., Хоперсков А.В., Храпов С.С. Проблема прогноза и управления гидрологическим режимом на горной территории в период ливневого паводка на основе гидродинамических численных экспериментов // Кибернетика и программирование. 2016. - № 3. - С. 35-53.

5. Дьяконова Т.А., Хоперсков А.В., Храпов С.С. Численная модель мелкой воды: использование графических процессоров NVIDIA CUDA / Суперкомпьютерные дни в России: Труды международной конференции (26-27 сентября 2016 г., г. Москва). 2016. С. 741-752.

6. Dyakonova T.A., Khoperskov A.V., Khrapov S.S. Numerical model of shallow water: the use of NVIDIA CUDA graphics processors // Communications in Computer and Information Science. V. 687. 14 p.

7. Воронин А.А., Васильченко А.А., Писарев А.В., Храпов С.С., Радченко Ю.Е. Проектирование механизмов управления гидрологическим режимом Волго-Ахтубинской поймы на основе геоинформационного и гидродинамического моделирования // Вестник Волгоградского государственного университета. Серия 1: Математика. Физика. 2016. № 1 (32). С. 24-37.

8. Белоусов А.В., Храпов С.С., Тен А.В., Садчиков Н.В., Болдырева Ю.А. Параллельный FDM принтер // Вестник Волгоградского государственного университета. Серия 1: Математика. Физика. 2016. № 4 (35). С. 116-131.

9. Дьяконова Т.А., Хоперсков А.В., Храпов С.С. Параллельная CUDA-версия программы для численного моделирования гидродинамических течений на основе СSPH-TVD метода / Свидетельство о государственной регистрации программы для ЭВМ N 2016610820 от 19.01.2016.

Аннотации публикаций

Форма 503

CSPH-TVD схема 2016

Springer CCIS