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

Вычислительные технологии

Том 4, № 6, 1999

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ВЗАИМОДЕЙСТВИЙ УЕДИНЕННЫХ ВОЛН С ПРЕПЯТСТВИЯМИ

К. Е. Афанасьев, С. В. Стуколов Кемеровский государственный университет, Россия e-mail: [email protected]

The problems on the interaction of solitary waves in ideal incompressible fluid with obstacles shaped as semi-circular cylinders and an oblique wall are considered. Diagrammes have been constructed on the basis of numerous calculations of the types of emerging wave flows. The problems in a complete nonlinear statement are solved by the methods of boundary elements.

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

В настоящей работе исследуются две задачи. Первая — взаимодействие уединенных волн с погруженным в жидкость препятствием в виде полукругового цилиндрического выступа, вторая — то же, с препятствием в виде наклонной твердой стенки, являющейся боковой границей области течения. Основными определяющими параметрами задач в первом случае являются высота волны и радиус цилиндра, во втором — высота волны и угол наклона твердой стенки. Цель работы — построение некой классификации взаимодействий волн с препятствиями по типу опрокидывания в зависимости от определяющих параметров. Задачи решаются в полной нелинейной постановке. Большое внимание уделяется исследованию кинематических характеристик возникающих волн. Поведение волны в последние перед обрушением моменты времени является существенно нелинейным, что усложняет численное моделирование этого явления. Для решения поставленных задач используются методы граничных элементов: для решения первой задачи — метод на основе интегральной формулы Грина (МГЭ) [6], второй — на основе интегральной формулы Коши (КМГЭ) [7]. Методика решения нестационарных задач гидродинамики идеальной

© К.Е. Афанасьев, С. В. Стуколов, 1999.

жидкости со свободными границами методом граничных элементов на основе интегральной формулы Грина изложена авторами уже в работах [2, 3], а на основе интегральной формулы Коши развивается сравнительно недавно [4]. МГЭ показал высокую эффективность и надежность, однако в работе [1] проведено сравнение МГЭ и КМГЭ и сделан вывод о том, что второй метод является более точным. Поэтому ниже подробно изложена методика решения второй задачи на основе интегральной формулы Коши. Для описания траекторий частиц (точек) свободной границы применяется лагранжев подход. При этом интеграл Коши — Лагранжа и кинематическое условие записываются в виде обыкновенных дифференциальных уравнений первого порядка. Задача является нестационарной, и для ее решения применяется метод Эйлера с автоматическим выбором шага по времени. Начальные условия (положение свободной границы и распределение потенциала на ней), описывающие уединенную волну, получены численно из решения нелинейной стационарной задачи [5]. Работоспособность методов и достоверность получаемых результатов проверяются на тестовых расчетах.

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

В расчетной области течения Б, ограниченной свободной поверхностью С\ и твердыми стенками С2 (рис. 1), решается уравнение Лапласа

=0, ^ = х + гу Е Б

(1)

для функции комплексного потенциала ю(г) = у(х,у) + гф(х,у), где у(х,у) — потенциал скорости и ф(х,у) — функция тока, удовлетворяющие условиям Коши — Римана. На твердых границах выполняется условие непротекания

ф = 0, г Е С2.

(2)

На свободной границе выполняются кинематическое и динамическое условия

d г

ду ду

т.--+ г Е С\,

ох ду

dу 1 |2 ^ _

- _ |Уу|2 + ду = 0, г Е Съ

(3)

(4)

Рис. 1. Схема расчетной области течения с препятствием в виде 1 — полукругового цилиндрического выступа, 2 — наклонной стенки.

где д — ускорение свободного падения. Начальное положение свободной границы и начальное распределение потенциала на ней получено из решения стационарной задачи об уединенной волне [5]. При этом предполагается, что на бесконечности жидкость покоится (для первой задачи у(+то) = у(—то) = 0, г е С1, для второй — у(—то) = 0, г е С1). В численных расчетах область течения ограничивается вертикальными твердыми стенками, на которых задается условие непротекания. Для удобства численной реализации задача приводится к безразмерному виду. В качестве характерных размерных величин выбираются ускорение свободного падения д и глубина бассейна Н. При этом краевая задача (1)-(3) останется без изменений, а уравнение (4) примет вид

^ - 21у^12 + У = 0,г е Сь (5)

2. Алгоритм движения по времени

Краевая задача (1) - (3), (5) является нестационарной, но в отличие от традиционных задач математической физики время явно входит только в граничные условия (3), (5), представляющие собой обыкновенные дифференциальные уравнения первого порядка, для интегрирования которых используется явный метод Эйлера.

Пусть в некоторый момент времени заданы положение свободной границы С\ и распределение потенциала 7°(х,у) на ней. Для нахождения функции тока ф°(х,у) и вектора скорости на С1 необходимо решить уравнение Лапласа (1) с условием 7°(х,у) на С1 и условием (2) на С2. Новое положение свободной границы и распределение потенциала на ней для момента времени + т можно вычислить, используя условия (3) и (5), дискретный аналог которых расписывается по схеме Эйлера следующим образом:

г к+1 = гк + (У^)к т,

7к+1 = 7к + (0, 5 | (Уу>)к |2 - ук)т, (6)

где гк ,ук,7к — значения функций на к-м шаге по времени. Таким образом, получаем смешанную краевую задачу для уравнения Лапласа (1) с граничными условиями (2) и (6), но уже для момента времени + т. Повторное ее решение и использование граничных условий позволяет определить положение свободной границы и распределение потенциала на ней для момента времени + 2т и т. д.

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

Более подробно алгоритм движения и выбора шага по времени изложен в работе [3].

3. Интегральная формула Коши

Для области Д ограниченной кусочно-гладкой замкнутой границей С = 1 и 2, справедлива интегральная формула Коши, которую можно записать с помощью предельных формул Сохоцкого в виде

/ ч 1 [ ^^(г) ,

= 1Ыг] 7-7°,Лг' (7>

С

где e(zo) = 2п для внутренней точки, e(z0) = п для точки на гладкой границе C, e(z0) = ß для угловой точки границы C (ß — угол при вершине). Положительное направление обхода контура C берется таким образом, чтобы область D оставалась слева.

Учитывая, что на свободной поверхности известна действительная часть ip(x,y) функции w(z), а на твердых стенках — ее мнимая часть ф(х,у), мы имеем смешанную краевую задачу для функции w(z). Численное решение этой задачи можно получить, разбив контур C на N линейных элементов Г узлами zj (j = 1, N). Тогда w(z) = lim G(z), где

max|rj | —>0

N N

G(z) — линейная глобальная пробная функция для z £ Г и G(z) = WjЛj (z), где

j=i j=i Wj — значение w(z) в точке zj, Лj• (z) — линейная базисная функция

[ (z — zj)/(zj — zj-l), z £ гj-l, Лj(z) = S (zj+i - z)/(zj+i - zj), z £ rj, 0, z £ Г-iUrj.

После указанного разбиения и линейной аппроксимации функции /ш(г) на границе интеграл Коши можно вычислить аналитически в смысле главного значения при г ^ ^. В результате получим

2niwj = wj+i — wj-i + wj ln

zj+i - zj zj-i — zj

+

N

E

m= 1 m=j,j +1

Im

(8)

где

wm+i — wm +

(zj — zm)wm+i (zj — zm+i)wm

zm+1 zm

zm+1 zm

ln

( zm+i

- zj

\ zm zj

Подставив в это равенство известные действительные или мнимые части функции при ] = 1, N , получим систему N х N линейных алгебраических уравнений

ВХ = Е, (9)

где В — полнозаполненная матрица, X — вектор неизвестных для определения И,е на С2 и 1тна С\, Е — вектор правой части. Для решения системы (9) используется метод Гаусса с выбором ведущего элемента.

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

Пусть двойной узел описывается тождеством гт+\ = гт. Предположим, что гт £ С и в этом узле задан потенциал И,е /ш, а гт+\ £ С2 ив нем задана функция 1т ,ш. В силу непрерывности функции в двойном узле будут выполняться естественные условия

1т 1шт = 1т тт+\, И,е тт+\ = И,е /шт. (10)

В этом случае из системы уравнений (9) следует исключить т- и (т +1)-е строки, заменив их условием (10). Кроме того, в т- и (т + 1)-х элементах столбцов матрицы В системы уравнений (9) отсутствуют вклады интегралов по элементу Гт, имеющему нулевую длину.

m

4. Вычисление давления, кинетической и потенциальной энергий волны

При исследовании взаимодействия солитонов с преградами актуальными являются задачи определения давления Р на преграду, изменения кинетической Ек и потенциальной Ер энергии. Кроме того, сохранение полной энергии Е = Ек + Ер позволяет осуществлять контроль за консервативностью численного метода.

Кинетическая и потенциальная энергии вычисляются на каждом временном шаге по формулам

b

1 f дф 1 х

Ek = - 2 ^~QsdS =- 12 ^ + + + (11)

i=1

Ng

a

1 } 1 Щ

ер = ^ у2 ^ = 6 + + У*+1)(х - )• (12)

а

Здесь а и Ь — абсциссы точек пересечения границ Сх и С2 (см. рис. 1), Ь — длина ¿-го элемента, N — количество точек на свободной поверхности, / = (дф/дз)г — производная по касательной.

Для вычисления давления необходимо решить дополнительную краевую задачу:

) = 0, г = х + гу е Б, (13)

1 2

Ъ = -2 1^|2 - У, г е Сь (14)

фг = 0, г е С2, (15)

где = —д(-) = Ъ + гфг. Решая систему уравнений (9), но уже для функции при

граничных условиях (14)-(15) находим щ на С2. После этого давление вычисляется по формуле

Р(х, у) = - (Ъ + 1 |^|2 + у], г е С2. (16)

5. Тестовые расчеты

Для подтверждения работоспособности выбранного метода проводится ряд тестов.

Первый тест заключается в том, чтобы найти решение уравнения Лапласа в области D = {0 < x < 2п; —1 < у < 0, 5sin(x)}, в которой на дне и вертикальных стенках ставится условие непротекания ф(х,у) = 0, а на верхней границе — условие <^(x,y) = — cos(x) ch (у + 1), правая часть в котором является гармонической функцией. Численные значения функции тока ф(х,у), найденные комплексным методом граничных элементов, сравниваются с точным решением: гфЬехЬ(х,у) = sin(x) sh (у + 1). В таблице при-

( max 1фПх — фПхЛ1\

ведена относительная погрешность оф =-——-- точного и численного зна-

\ max |фп | J

чений функции ф от числа узлов N по границе с указанием числа узлов на свободной

границе Ng. Вектор скорости (Vx,Vy) в каждом узле на свободной поверхности вычислялся методом, изложенным в работе [3], и сравнивался с точным значением V!cext = sin(x) ch (y + 1), Vyiexi = — cos(x) sh (y + 1). В третьей и четвертой колонках приведены

гт, max |Vxtexi — Vxtexil max | Vytext — Vytext|

значения oVx =-:—т—--, oVy =-1-1-. В пятой колонке приводятся

x max |Vxtext| ' y max |Vytext|

числа обусловленности K(B) матрицы B системы линейных алгебраических уравнений (9). Сравнительно небольшое число обусловленности объясняется тем, что в основе численного метода лежит уравнение Фредгольма второго рода.

Таблица

N N 6ф OVx OVy K (B)

72/30 5.5E-03 2.0E-03 1.6E-02 3.4

145/60 1.3E-03 6.5E-04 9.3E-03 3.6

290/120 3.1E-04 6.2E-04 8.7E-03 4.3

580/240 7.6E-05 6.0E-04 5.7E-03 5.7

Второй тест проводится на решении нестационарной задачи о движении уединенной волны амплитуды = 0.5 по бассейну постоянной глубины Н = 1. В этом тесте важным является то, что уединенные волны в процессе движения не изменяют амплитуду и скорость, сохраняют форму и полную энергию. Для расчета была рассмотрена область О = {-15 < х < 15; —1 < у < у0}, где у0 описывает уединенную стационарную волну. На границе области взято 350 элементов, из них 200 — на свободной поверхности. Вершина волны при Ь = 0 находилась в точке х = —5, у = 0.5. Расчет проводился до момента безразмерного времени Ь = 8.28, когда вершина волны перешла в точку с абциссой х = 5. К этому моменту времени волна прошла путь, равный 1.3 длины волны, определяемой длиной отрезка по оси х, на котором выполняется условие 1т г(¿) > 0.01А(Ь), г € С [9].

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

Данный тест дополняется расчетом взаимодействия уединенной волны амплитуды А = 0.4 с вертикальной стенкой, имеющей абциссу х = 5. Остальные параметры задачи были взяты из предшествующего теста. На рис. 2, б приведены профили свободной поверхности для нескольких моментов безразмерного времени (Ь = 0, утах = 0.4 — первоначальная форма солитона; Ь = 8.91, утах = 0.95 — форма свободной поверхности в момент максимального заплеска; Ь = 19.12, утах = 0.392 — форма восстановленного солитона, вершина которого в этот момент времени находилась в точке с абциссой х = —7.5. Отражение

у, 0.0 < t < 8.28

0.6 -

0.4 -

0.2 -

-0.2 --о.4: // = = 1 А= 0.5 Е = 0.61323 АЕ = 1 .1 -Ю-2

-8.0 -6.0 -4.0 -2.0 0.0 2.0 4.0 6.0 8.0 X

у, t = 8.91

0.5 - í = = 19.121 4 = 0.0 ¿ = 14.121 г = 4.538 /

0.0

Н= 1 ,4 = 0.4 Е = 0.4426 /\Е = 1 4-Ю-2

0.5

-8 .0 -6.0 -4.о -г1о о!о г'о 4.'о \

Р 1 — А = 0.2 2 — ,4 = 0.3

0.8 3 — А = 0.4

0.6 4 — ,4 = 0.5

5 — А = 0.6

0.4

0.2

6.0 8.0 10.0 12.0 1

Рис. 2. Результаты тестовых расчетов.

волны от стенки приводит к изменению амплитуды волны и формированию хвоста из вторичных волн малой амплитуды, что находится в полном соответствии с результатами численных исследований [9].

Кроме кинематических особенностей возникающих течений интерес представляет задача определения динамического воздействия волн на вертикальную стенку. На рис. 2, в для волн различных амплитуд показаны хронограммы волнового давления в точке стенки, совпадающей с начальным урезом жидкости. Точками * и + на график нанесены результаты расчетов для волн амплитуды А = 0.4 и А = 0.5 соответственно из работы [11]. При амплитудах А > 0.3 расчетные хронограммы имеют два локальных максимума. Это явление можно объяснить следствием действия сил инерции. Кроме того, первый максимум имеет большее значение, чем второй, при этом момент максимального заплеска волны на стенку не совпадает с моментами локальных максимумов давления. Эти особенности при накате солитонов на вертикальную стенку подтверждаются экспериментами [8].

6. Численные результаты

6.1. Движение над дном с полукруговым выступом

Задача о накате волны на полукруговой цилиндрический выступ подробно описана в работе [12], где показано, что взаимодействие волны и погруженного цилиндра порождает различные волновые картины. Авторы провели классификацию волновых движений и определили пять зон, соответствующих различным волновым режимам: В-Ц — волновая цепь, О-В — опрокидывание вперед, О-Н — опрокидывание назад, О-Г — обмен

гребнями, Н-Т — неустойчивость Танаки. По данной задаче нами также проводились расчеты независимо от вышеуказанной работы. После знакомства с приведенной в работе диаграммой было проведено детальное сопоставление результатов. Оказалось, что наши расчеты полностью вкладываются в предлагаемую волновую классификацию и хорошо ложатся на диаграмму (рис. 3), взятую из [12] (точками отмечены выполненные численные расчеты). На рис. 4 приведены результаты расчетов, показывающие характерные для соответствующих зон волновые процессы. Геометрия области течения схематически показана на рис. 1 (случай 1). Глубина бассейна Н = 1, центр цилиндра всегда находился в точке х = 0, у = -1, вершина волны располагалась в начальный момент времени в точке х = -8, у = А (А — амплитуда волны). На границе области бралось 370 элементов, из них 197 — на свободной поверхности.

Волновая цепь — это режим, при котором солитон, проходящий над препятствием, практически не изменяет своей формы, а сзади него формируется последовательность затухающих волн малой амплитуды. Расчет по данному волновому режиму представлен на рис. 4, а для волны амплитуды = 0.7 и радиуса Я = 0.4. Под режимом "опрокидывание вперед" понимается такой режим, при котором волна, проходя над препятствием, трансформируется так, что на ее переднем фронте зарождается "вторая" волна. При этом

Рис. 3. Диаграмма волновых режимов при взаимодействии уединенной волны с полукруговым цилиндрическим препятствием. — амплитуда солитона, К — радиус цилиндра.

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

На рис. 4, б представлен расчет по накату уединенной волны амплитуды = 0.462 на полукруговой цилиндр радиуса R = 0.8. Для данного случая режим опрокидывания напоминает скользящий бурун. Этот конкретный случай просчитан в работе [12] и на рис. 5, а проводится сравнение формы свободной поверхности в моменты, близкие к обрушению. Другой расчет, также приведенный в работе [12] для амплитуды = 0.462 и радиуса R = 0.7, но уже относящийся к режиму "опрокидывание назад", представлен на рис. 4, в. Данный режим характеризуется тем, что солитон в процессе прохождения над препятствием принимает форму волны, имеющей амплитуду, меньшую первоначальной, но на заднем фронте прошедшей волны образуется впадина, в которой формируется всплеск, опрокидывающийся "против движения" основной волны. Момент, близкий к опрокидыва-

Рис. 4. Волновые картины течений при взаимодействии солитона с полукруговым цилиндрическим выступом. а: В —Ц — волновая цепь, б: О-В — опрокидывание вперед, в: О —Н — опрокидывание назад; г: Н —Т — неустойчивость Танаки.

нию, а также сравнение с расчетом уже цитируемых авторов представлены на рис. 5, б.

Режим "обмен гребнями" соответствует режиму, при котором происходит плавный "переток" гребня уединенной волны через препятствие. При подходе волны к препятствию ее амплитуда падает, а за телом волна вновь зарождается в виде другого солитона. В режимах О-В и О-Н также наблюдается обмен гребнями, но здесь более важными являются процессы, связанные с обрушением волн при трансформации солитонов. Характерный для этой зоны режим распространения солитона показан на рис. 4, г, где приведена форма свободных границ в различные моменты времени для солитона амплитуды = 0.25 при Я = 0.65.

Режим "неустойчивость Танаки" характеризуется тем, что волны большой амплитуды > 0.78 при взаимодействии с препятствиями радиуса, меньшего 0.5, могут вести себя двояко. С одной стороны, происходит опрокидывание гребня вперед по движению, а с другой — может наблюдаться трансформация волны в волну меньшей амплитуды, но имеющую ту же полную энергию, что и первоначальная волна. Относительная ошибка изменения полной энергии для всех расчетов находилась в пределах 2 %.

Рис. 5. Сравнение с расчетом [12] — кривая 1, расчет авторов — кривая 2.

6.2. Взаимодействие уединенных волн с наклонной твердой стенкой

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

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

Геометрия области течения схематически показана на рис. 1 (случай 2). Глубина бассейна Н = 1, дно начинает возвышаться в точке х = 0, у = —1 под углом а, вершина волны в начальный момент находилась в точке х = —5, у = А (А — вершина волны). На границе области бралось 350 элементов, из них 200 — на свободной поверхности. Контроль точности решения осуществлялся по закону сохранения полной энергии. Относительная ошибка изменения полной энергии для всех расчетов не превышала 2%. Меняющимися параметрами задачи были угол наклона правой стенки, который принимал значения от 5 до 90° с шагом 5и амплитуда 0.2 < А < 0.6.

В результате расчетов были выявлены четыре зоны течений по типу опрокидывания в зависимости от угла наклона стенки (рис. 6).

Первая зона (О-В) соответствует малым углам наклона стенки (пологое дно) и характеризуется тем, что волна опрокидывается вперед (по ходу движения в начальный момент времени) во время наката на берег. На рис. 7, а представлены профили свободной поверхности в последние перед обрушением моменты времени для угла наклона стенки а = 5° и А = 0.5.

Вторая зона (О-В-О) характеризуется тем, что волна опрокидывается по-прежнему вперед, но во время отката от берега. На рис. 7, б демонстрируются результаты расчета для угла наклона стенки а = 20 ° и А = 0.5.

Третья зона (О - Н) характеризуется тем, что волна опрокидывается во время отката от берега, но в противоположном направлении относительно движения в начальный момент времени. Пример такого поведения волны при накате на наклонную стенку (а = 50 °) приведен на рис. 7, в.

В диапазоне углов наклона стенки, близких к вертикальным, волна откатывается без опрокидывания, порождая за собой след из затухающих волн малой амплитуды (В-Ц) (рис. 7, г).

Заштрихованные (см. рис. 7) переходные зоны соответствуют параметрам задачи, при которых возникающие волновые картины имеют признаки соседствующих зон. Строились они следующим образом: отмечались параметры задачи, для которых возникающую волновую картину можно с уверенностью отнести к одному типу опрокидывания, а затем

А

ОБО 0 н

0.5 о * ® ®-5- -

- В Ц

0.4 - —

0.3

10 20 30 40 50 60 70 80 а, град

Рис. 6. Диаграмма волновых режимов при взаимодействии уединенной волны с наклонной стенкой. — амплитуда солитона, а — угол наклона боковой стенки.

Рис. 7. Волновые картины течений при взаимодействии солитона с боковой наклонной стенкой. а: О-В — опрокидывание вперед, б: О-В-О — опрокидывание вперед во время отката, в: О-Н — опрокидывание назад, г: В-Ц — волновая цепь (1-4 — £ = 0; 5.37; 8.48; 13.07).

Рис. 8. Хронограммы волнового давления при накате волны амплитуды А = 0.5 на наклонную стенку в точке х = 0, где дно начинает возвышаться (а), и в точке стенки, совпадающей с начальным урезом жидкости (б). 1-4 — а = 5, 20, 50 и 60°.

параметры, соответствующие другому типу опрокидывания. Область параметров, лежащих между ними, штриховалась.

На рис. 8 показаны хронограммы волнового давления при накате волны амплитуды А = 0.5 на наклонную стенку в точках х = 0, где дно начинает возвышаться, и в точке стенки, совпадающей с начальным урезом жидкости. Кривая 1 на рис. 8, а соответствует расчету наката волны на пологое дно (а = 5°) и имеет одно максимальное значение. Этот факт объясняется тем, что волна успевает полностью миновать точку излома дна до своего опрокидывания. Кривая 2 соответствует углу а = 20° и также имеет один максимум, но с некоторого момента времени давление начинает увеличиваться. Это объясняется тем, что к моменту опрокидывания волны ее некоторая часть уже отразилась от стенки и проходит над точкой излома дна. Наличие двух максимумов кривой 3 объясняется тем, что на момент опрокидывания волны ее большая часть уже отразилась от стенки и прошла точку излома. Кривая 4 соответствует режиму, когда опрокидывание волны не наблюдается. Данная хронограмма похожа на кривую, полученную при накате на вертикальную стенку. Аналогичные выводы можно сделать и для кривых на рис. 8, б. Здесь кривая 1, соответствующая пологому дну, отсутствует в связи с тем, что при накате волны амплитуды А = 0.5 на наклонную стенку под углом а = 5° волна опрокидывается до подхода к начальной точки уреза жидкости.

На рис. 9 представлены графики максимального заплеска волн в зависимости от угла наклона стенки. Заметно качественное совпадение поведения графиков. Количественное же различие можно объяснить с одной стороны тем, что используемая для решения задачи модель не учитывает вязкости жидкости, с другой стороны тем, что при проведении

у

0.5

Н -75---,--

Ч-1-1-8 | | | |

20 зо 40 50 60 70 80 а, град

А=0.6

0.5

0.4

0.3

0.2

Рис. 9. Графики максимального заплеска волн на наклонную стенку. □ и о — сравнение с результатами расчетов [10], х и * — сравнение с экспериментальными данными [8], ♦ и И--то же,

экспериментов трудно получить уединенную волну заданной амплитуды. В результате движения волнопродуктора порождается волна, которая при последующем движении по ровному дну принимает форму солитона, но при этом может уменьшиться ее амплитуда.

7. Заключение

Расчеты задачи о взаимодействии уединенных волн с полукруговым цилиндрическим выступом полностью согласуются с данными работы [12]. Решение же этой задачи КМГЭ представляется нецелесообразным, так как подробные результаты отражены в работе [12], где в основе решения лежал метод, близкий к КМГЭ.

Расчеты задачи о накате солитона на наклонный берег выявили четыре зоны волновых течений, позволяющие в зависимости от амплитуды волны и наклона дна определить тип волновой картины течения в прибрежной зоне и предсказать гидродинамическое воздействие волны на преграду. Данное исследование проводилось в диапазоне амплитуд волн от 0.2 до 0.6 и в этом смысле не может претендовать на полную классификацию возникающих течений.

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

[1] Афанасьев К. Е. Решение нелинейных задач гидродинамики идеальной жидкости со свободными границами методами конечных и граничных элементов. Дис. ... д-ра физ.-мат. наук. Кемерово, 1997.

[2] Афанасьев К.Е., Афанасьева М.М., Терентьев А. Г. Исследование эволюции свободных границ методами конечных и граничных элементов при нестационарном движении тел в идеальной несжимаемой жидкости. Изв. АН СССР. МЖГ, №5, 1986,

[3] Афанасьев К. Е., Самойлова Т. И. Техника использования метода граничных элементов в задачах со свободными границами. В "Вычислительные технологии". ИВТ СО РАН, Новосибирск, 4, №11, 1995, 19-37.

[4] Афанасьев К.Е., Стуколов С. В. Накат уединенной волны на наклонный берег. Вестник Омского ун-та, №3, 1998, 9-12.

[13].

8-13.

[5] Афанасьев К. Е., Стуколов С. В. О наличии трех решений при обтекании препятствий сверхкритическим установившимся потоком тяжелой жидкости. ПМТФ, №1, 1999, 27-35.

[6] Бреббия К., Теллес Ж., Вроубел Л. Методы граничных элементов. Мир, M., 1987.

[7] ГРОМАДКА Т., Лей Ч. Комплексный метод граничных элементов. Мир, M., 1990.

[8] МаноЙлин С. В. Некоторые экспериментально-теоретические методы определения воздействия волн цунами на гидротехнические сооружения и акватории морских портов. Препринт №5. ВЦ СО АН СССР, Красноярск, 1989.

[9] Протопопов Б. E. Численный анализ трансформации уединенной волны при отражении от вертикальной преграды. Изв. АН СССР. МЖГ, №5, 1990, 115-123.

[10] ФРАНК А. М. Дискретная нелинейно-дисперсионная модель мелкой воды. ПМТФ, №5, 1993, 15-24.

[11] Шокин Ю.И., РУЗИЕВ Р. А., ХАКИМЗЯНов Г. С. Численное моделирование плоских потенциальных течений жидкости с поверхностными волнами. Препринт №12. ВЦ СО АН СССР, Красноярск, 1990.

[12] Cooker M. J., Peregrine D.H., Vidal C., Dold J.W. The interaction between a solitary wave and a submerged semicircular cylinder. J. Fluid Mech., 215, 1990, 1-22.

[13] Synolakis C. E. The runup of solitary waves. J. Fluid Mech., 185, 1987, 523-545.

Поступила в редакцию 1 декабря 1998 г.