ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 1, с. 113-120
УДК 519.634
ПРИМЕНЕНИЕ КООРДИНАТНЫХ ПРЕОБРАЗОВАНИЙ ДЛЯ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ НАКАТА ВОЛН ЦУНАМИ
МЕТОДОМ КРУПНЫХ ЧАСТИЦ1
© 2015 г. А. В. Кофанов, В. Д. Лисейкин, А. Д. Рычков
(630090 Новосибирск, просп. Лаврентьева, 6, Ин-т вычисл. технологий СО РАН;
630090 Новосибирск, ул. Пирогова, 2, Новосибирский гос. ун-т) e-mail: avkof87@gmail.com, lvd@ict.nsc.ru, rych@ict.nsc.ru Поступила в редакцию 30.09.2013 г.
Переработанный вариант 31.03.2014 г.
Рассматривается численный алгоритм решения задачи о накате уединенной волны цунами на берег со сложным рельефом береговой кромки. Алгоритм включает применение технологии построения координатных отображений, преобразующих равномерную прямоугольную сетку в эталонной вычислительной области в сетку в физической области с ячейками, сгущающимися в зоне уреза воды. Применение таких координатных отображений позволяет существенно уменьшить число точек сетки и, соответственно, время счета задачи. В качестве математической модели используются уравнения мелкой воды, а для решения задачи применяется метод крупных частиц. Приведены иллюстрации результатов расчета криволинейной сетки и конфигурация зоны затопления на реальном примере. Библ. 17. Фиг. 4.
Ключевые слова: цунами, уравнения мелкой воды, метод частиц в ячейках, координатные преобразования, адаптивные разностные сетки, обращенные уравнения Бельтрами и диффузии.
БО1: 10.7868/80044466915010147
1. ВВЕДЕНИЕ
Японский термин "цунами" означает набегающую на побережье большую волну, вызванную подводным землетрясением, извержением подводного вулкана или подводным оползнем. Это опасное по своим последствиям природное явление часто приводит к гигантским разрушениям и гибели большого числа людей. При разработке средств прогнозирования последствий такого явления ключевую роль играет численное моделирование процессов возникновения и распространения волн цунами, а также их наката на побережье. Такая задача рассматривается в данной работе в рамках модели мелкой воды, которая достаточно адекватно описывает подобные процессы. Основной целью моделирования является определение зоны затопления суши при накате на побережье уединенной волны. Такого класса задачи относятся к задачам волновой гидродинамики со свободной поверхностью, для решения которых разработаны достаточно эффективные численные методы. Прежде всего, это конечно-разностные методы, в том числе основанные на применении метода дробных шагов, методы конечных элементов (см. [1], [2]), неконсервативные методы характеристик (МОС) (см. [3]) и методы конечных объемов (см. [4], [5]). Однако следует заметить, что в случае сложных нерегулярных рельефов дна, характерных для реальных прибрежных районов, большинство из этих подходов оказываются малоэффективными и требуют применения особых численных алгоритмов, как это, например, предложено в [4] для случая отсутствия сухого дна. В последнее время большое внимание уделяется развитию бессеточных методов типа SPH (Smoothed particle hydrodynamics) (см. [6]), являющегося дальнейшим развитием метода частиц в ячейках (PIC). В отличие от метода PIC, в SPH отсутствует эйлеров этап. Это позволило, с одной стороны, сделать метод более робастным, применимым для расчета течений в сложных областях, но с другой стороны, потребовалось разработать ряд сложных процедур сглаживания для описания взаимодействия между частицами, что не обеспечивает полную
1) Работа выполнена при финансовой поддержке РФФИ (код проекта 13-01-00231) и Президентской программы поддержки ведущих научных школ РФ № НШ-5006.2014.9.
8
113
консервативность метода. В [7] предложен оригинальный комбинированный подход, основанный на сочетании конечно-разностного метода и метода SPH, который позволяет решать уравнения мелкой воды в водоемах со сложным нерегулярным рельефом дна. В данной работе для моделирования наката волны цунами на сухой берег используется достаточно хорошо известный метод крупных частиц (см. [8]), в котором, как и в методе SPH, используются частицы в виде "элементарного столбика воды", называемые в дальнейшем "крупными частицами". Метод крупных частиц является развитием метода "частиц в ячейках" Харлоу, широко применяется для решения эволюционных систем дифференциальных уравнений, связанных с описаниями течений сплошной среды, и основывается на расщеплении исходной системы уравнений по физическим процессам. Процесс решения состоит из двух этапов — эйлерова и лагранжева. На эйлеровом этапе расчета изменяются лишь величины, относящиеся к ячейке в целом, а жидкость предполагается моментально заторможенной. На лагранжевом этапе определяются потоки массы и импульса через границы ячеек в предположении, что они переносятся только за счет нормальной к границе составляющей скорости. Это позволяет существенно упростить процедуру учета взаимодействия между частицами, обеспечивая в то же время достаточно точное воспроизведение наката (отката) воды на сухой берег и строгое выполнение законов сохранения. Заметим, что число частиц, необходимых для решения двумерных задач в методе SPH с необходимой точностью, составляет примерно 100000 (см. [9]). Это число сравнимо с числом узлов разностной сетки в методе крупных частиц, однако программная реализация последнего существенно проще, а число операций, затрачиваемое на расчет движения одной крупной частицы, примерно на два порядка меньше, чем в SPH. Кроме того, разностная схема, основанная на данном методе, является хорошо сбалансированной, обеспечивающей выполнение условия гидростатического равновесия жидкости в покоящейся воде.
Эффективность и экономичность численного моделирования также существенно зависят от введения новых координатных преобразований, с помощью которых конструируются разностные сетки, адаптирующиеся к быстрым изменениям физических характеристик (скорость движения среды), а также к рельефу дна и кромке береговой линии. За счет технологии адаптации (сгущение только в нужных зонах) удается существенно уменьшить число ячеек разностной сетки при приемлемой точности расчетов. Сам процесс построения координатного преобразования или разностной сетки должен быть автоматизирован для проведения расчетов, обеспечивающих оперативные предсказания возможных последствий воздействия цунами на прибрежную зону. В данной работе для решения задачи наката волны цунами применяется современная технология построения координатных преобразований и разностных сеток, базирующаяся на решении обращенных уравнений Бельтрами и диффузии (см. [10]).
2. ЧИСЛЕННЫЙ МЕТОД
В рамках модели мелкой воды рассматривается задача о накате уединенной волны на берег со сложным рельефом береговой кромки. С учетом силы трения, осредненной по глубине, система уравнений, описывающая движение такой волны в декартовой системе координат, координатная плоскость которой г = 0 совпадает с поверхностью "спокойной" (неподвижной) воды, записывается в виде
дН + дНи + дНу = 0 dt дх ду
дНи + dHu¿ + Жиу + g Н = нdh _ uhcr , (1)
дt дх ду 2 дх дх
дНу + cHuv + дНу2 + g дН2 = gH _ vHC дt дх ду 2 ду ду где H(x, y, t), h(x, y) — глубина воды и профиль дна соответственно, отсчитываемые от поверхности воды, u(x, y, t), v(x, y, t) — декартовы координаты вектора скорости u(x, y, t), CR — коэффициент трения, g — ускорение силы тяжести. Для определения коэффициента CR использовалась известная формула
2
Cr = Нт\U,
где n — эмпирический коэффициент (коэффициент Шези). Областью решения системы (1) являлся прямоугольник D{0 < x < Lx, 0 < y < Ly} с профилем дна, взятым из соответствующей батиметрии. На верхней и нижней границах в области волны задавались неотражающие граничные
(б)
5000
10000
15000
20000 X
10000
12000
14000
шШш&ш
жмгт ШтшШ :ттш
шшт
Г'"!:»»:!"'!
иш
16000 X, м
Фиг. 1. Расчетная криволинейная сетка со сгущением ячеек к линии уреза волны (а) и фрагмент сетки (б).
условия, в области суши все параметры жидкости в (1) полагались равными нулю. На правой границе задавалось изменяющееся во времени значение величины И(Ьх, у, ?), взятое из анализа подхода к побережью океанской волны.
Для численного решения системы (1) применялся метод крупных частиц (см. [8]) на адаптивной сетке. Использование "классического" метода частиц в ячейках оказалось невозможным в силу статистического характера получаемого с его помощью решения. Такое решение не удовлетворяет условию гидростатического равновесия жидкости в покоящейся воде и значительно искажает процесс распространения волны. Метод крупных частиц лишен этого недостатка и позволяет получить сбалансированную разностную схему. При реализации данного метода определенную проблему представляет выбор разностной сетки. Как показал наш опыт, использование прямоугольных равномерных сеток возможно, но для получения приемлемой точности требуются сетки с очень мелкими ячейками во всей области и, значит, с большим числом узлов, что в силу ограничения на шаг по времени, присущего явным разностным схемам, требует достаточно большого времени счета. Поскольку в таких задачах наибольший интерес представляет расчет зон затопления суши, можно уменьшить число узлов разностной сетки, если использовать сетки, сгущающиеся вблизи линии уреза, например, можно строить блочные прямоугольные сетки. Однако из-за значительной искривленности линии уреза (подвижной линии, разделяющей воду и сушу) получается весьма сложная нерегулярная структура блоков с несогласованными сетками на их границах. Поэтому в данной работе предложено использовать адаптивные криволинейные разностные сетки, которые строятся с помощью зависящего от времени координатного преобразования
(хЦ,%),Ж%)) :Н2 ^ Д % = (2)
из вычислительной области, являющейся единичным квадратом Е2 = {0 < £,1,£,2 < 1}, в физическую область решения Б{0 < х < Ьх, 0 < у < Ьу}. Отображение с помощью этого преобразования
узлов равномерной прямоугольной сетки в 22 задает узлы неравномерной сетки в физической области Б. Координатное преобразование (2) определяется таким образом, чтобы ячейки полученной криволинейной сетки были
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.