Вычислительные технологии
Том 5, № 3, 2000
НЕЛИНЕЙНЫЕ ПОВЕРХНОСТНЫЕ И ВНУТРЕННИЕ ВОЛНЫ, ВЫЗВАННЫЕ КОЛЕБАНИЯМИ КРУГОВОГО ЦИЛИНДРА В МНОГОСЛОЙНОЙ ЖИДКОСТИ *
С. И. Горлов Омский филиал Института математики им. С. Л. Соболева СО РАН, Россия e-mail: gorlov@iitam.omsk.net.ru
Within the framework of the nonlinear theory the initially-boundary problem about periodic oscillations of the circular cylinder in three-layer fluid is considered. The liquid in each stratum is ideal, incompressible, heavy and homogeneous. The cylinder is completely located in a low layer. The system of integro-differential equations of the problem is obtained expressing kinematic and dynamic conditions on the interface of liquid media as well as the condition of the zero velocity normal component at the contour points. These equations contain an unknown function, inherent in the intensities of vortex sheets located along the liquid boundaries and sources, simulating the contour. The solution of the given system is based on two iterative processes, one of which is connected with the integration over time under the fifth-order Runge-Kutta-Fehlberg scheme, and the other involves the solution of a system of the linear algebraic equations obtained by the discretization of integral equations at each time step with the help of high-order panel method. In order to reduce the computation costs due to the numerous collocation points on the boundaries, the method of areas decomposition is used. The horizontal and vertical oscillations of circular cylinder in three-layer fluid (salty, fresh water and air medium) and in two-layer fluid (water-air medium) are considered. The obtained form of the fluids interfaces and total hydrodynamic performances of circular cylinder are indicated. For parameters considered in work, the essential influence of an additional stratum on fluid characteristics is not revealed.
Введение
В последнее время большое внимание уделяется исследованию влияния поверхностных и внутренних волн на гидродинамические характеристики контура. Задача о колебаниях контура в многослойной весомой жидкости принадлежит к этому классу и имеет несомненный практический интерес, связанный с проектированием транспортных средств и подводных сооружений. Значительные успехи в решении задач подобного рода получены в рамках линейной теории, позволяющей значительно упростить кинематические и динамические условия на границах раздела сред [3, 4]. Настоящая работа посвящена разработке
* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант №96-01-00093.
© С. И. Горлов, 2000.
численного метода решения нелинейных начально-краевых задач о движении контура в многослойной весомой жидкости, состоящей из слоя соленой, пресной воды и воздушной среды. Следует отметить, что реальные океанские течения имеют более сложную структуру с неоднородным распределением плотности внутри каждого слоя, но многочисленные теоретические и экспериментальные исследования показывают, что в рамках такой модели можно качественно оценить влияние генерируемых волн на распределенные и суммарные гидродинамические нагрузки погруженного тела.
1. Система интегродифференциальных уравнений задачи
Рассмотрим жидкость, состоящую из слоев D1, D2 и D3 (нумерация начинается снизу, слои D1 и D3 полубесконечны). Жидкость в каждом слое является идеальной, несжимаемой, тяжелой и однородной. Обозначим через t время, рд. — плотность жидкости в слое Dk (k = 1, 2, 3), Li(t) — граница раздела областей Di(t) и Dl+1(t) (/ = 1, 2), L0(t) — контур, описывающий круговой цилиндр радиуса R, g — ускорение силы тяжести, H — расстояние между L1(t) и L2(t) в начальный момент времени. Ось x введенной системы координат совпадает с невозмущенным уровнем границы раздела L2(t). Цилиндр полностью расположен в нижнем слое, а его центр имеет координаты:
— Vo
xc(t) = — cos (ut)(1 — S), yc(t) = — h--cos (ut)S,
и и
где S =1 соответствует вертикальным колебаниям, а S = 0 — горизонтальным, и — частота колебаний, h — ордината осредненного положения центра цилиндра.
Введем в рассмотрение интенсивности вихревых слоев Yi(si,t) и Y2(s2,t), расположенных вдоль границ раздела L1(t) и L2(t) соответственно, а также слой источников q(s0,t), моделирующих контур L0(t) (sj — дуговая координата точки z(sj) G Lj(t) (j = 0,1, 2)). Предположим, что в бесконечно удаленных точках границ раздела L1(t) и L2(t) затухают возмущения функций Y1(s1,t) и Y2(s2,t), т.е. Yi(±œ,t) = 0 (/ = 1, 2). Тогда движение жидкости в областях Dk (t) (k = 1, 2, 3) будет описываться функцией
—( = Г Y1(s1,t) ds1 + 1 f Y2(s2,t) ds2 + ^ f g(so,t) dso m
(Z, ) = z — Z(s1) + z — Z(s2) + 2n J z — Z(so) ' ( )
Li(t) L2(t) Lo(t)
На основе метода, описанного в [1], получена система интегродифференциальных уравнений, соответствующих кинематическому и динамическому условиям на границах раздела L1(t) и L2(t) и условию непротекания в точках L0(t), которая имеет вид:
= —i(z(si ),t), z(si ) G Li (t) (l = 1, 2), (2)
ÔGi(si,t) i(z(si),t)|2 , T Y2(si,t)\ i o\ ^
~Pi*-Ô--+ gImz(si)--3- , z(si) G Li(t) (l = 1, 2), (3)
dt ^V 2 8
si
Gi(si,t)= / fl^ + 2pi,V^i,t)) dai, = Pl — Pl+1 (l = 1, 2), J V 2 y Pi + pi+1
Vjs(Sj,t) = Re (Vj(z(Sj),t)ei0j(Sj'4)) , z(s,) G Lj(t) (j = 0,1, 2), = Im ((Vo(z(so),t) - VLo(t)^0'^) , z(s0) G L0(t),
Vl0(t) = Xc(t) - iyyc(t). (4)
Здесь Vj(z(sj),t) определяется (1) для z(sj) G Lj(t) (j = 0,1, 2), при этом несобственные интегралы, входящие в выражение (1) для комплексной скорости, следует понимать в смысле главного значения по Коши. Обозначения: 6j (sj ,t) — угол между касательной в точке z(sj) G Lj(t) и осью x (j = 0,1, 2), Vjs(sj, t) — касательная составляющая скорости в точках контура Lj (t) (j = 0,1, 2), Vl0 (t) — комплексная скорость колебаний кругового цилиндра в выбранной системе координат.
В бесконечно удаленных точках областей Di(t), D2(t) и D3(t) затухают возмущения скоростей и границ раздела сред:
lim V(z,t) = 0, lim Imz(s) = (/ - 2)H при z(s) G L(t) (l = 1, 2). (5)
В начальный момент времени отсутствуют возмущения скоростей и границ раздела сред:
Imz(si) = (/ - 2)H при z(si) G L(0), y(s, 0) = q(so, 0) = 0 (l = 1, 2). (6)
Распределенные и суммарные гидродинамические характеристики определяются при помощи интеграла Коши — Лагранжа, записанного в подвижной системе координат:
p(so,t) - f (t) = -Pl ^dt / Vos(ao,t) dao - Re (Vlo(t)Vo(z(so), t)) + ¡VoizM^2j , (7)
Rx - iR = i J (p(so,t) - f (t)) е-*°(в°'4) dso, (8)
Lo(t)
где f (t) — некоторая функция, зависящая только от времени.
2. Решение системы интегродифференциальных уравнений
Полученная система интегродифференциальных соотношений (1)-(6) является нелинейной. Эта нелинейность обусловлена двумя факторами: интенсивности особенностей 71^1,¿), 72(з2,и д(з0,£) входят в граничные условия нелинейным образом и неизвестны границы раздела сред (¿) и Ь2(^), которые определяются в процессе решения. Этот факт свидетельствует о наличии определенных трудностей при решении системы.
Для интегрирования системы уравнений (2), (3) по времени будем использовать метод Рунге — Кутта — Фельберга пятого порядка точности. При этом в каждый момент времени ¿п (п = 1, 2,...) вычисляются значения функций СП(^П) и точки границ раздела zn(sn) С ¿П (1 = 1, 2, верхний индекс служит для обозначения функции, вычисляемой на
п-м временном шаге). Задача нахождения значений 72 и 5п(50) на каждом шаге
по времени сводится к решению системы линейных интегральных уравнений:
+ 2^(0 = -^, ^п(^п) е ¿п (1 = 1,2), (9)
qn(so)
Im ((V0(zn(sо)) - V2o)е^) , zn(so) E Ц. (10)
2 ---- и\~ ■ ) 1 - и/ — —'о ■
Систему интегральных уравнений (9), (10) будем решать методом панелей высокого порядка [6]. Для этого введем разбиение контуров ¿п и ¿п на интервалы [зп_ 1, з^] (I =1, 2; г = 1,... , I) и ] (3 = 1,... , J) соответственно. На этих интервалах вы-
берем точки коллокации ^^ € ¿п К** е и ^) € ¿п (з0а е N ^-ь^.
Уравнения (9) и (10) будем рассматривать в точках гп(Зп*) (1 = 1, 2; г = 1,... , I) и гп(з0 ^) (3 = 1,... , J) соответственно. Контуры ¿п на г-м интервале [Зч^, и контур ¿п на 3-м интервале ] аппроксимируются параболой, а функции 7гп(зп) и дп(з0) на этих
же интервалах — линейной функцией. После дискретизации интегральных уравнений (9), (10) с учетом (1) получается система линейных алгебраических уравнений относительно значений функций 7гп(зп) и ^п(з0) на концах интервалов [Зч^ ,5^] и ] (I = 1, 2;
г = 1,..., I; 3 = 1,...^). Решив эту систему, из (1) найдем значения V (г) в точках контура гп(з0 ^) € ¿п (3 = 1,..., J), а из (7), (8) — распределенные и суммарные гидродинамические характеристики.
3. Решение конкретных задач
При помощи описанного численного метода рассмотрены задачи о вертикальных и горизонтальных колебаниях кругового цилиндра в трехслойной (соленая, пресная вода и воздушная среда, pi* = 0.0148, р2* = 1) и двухслойной (pi* = 1) жидкости. Безразмерными параметрами задач являются следующие величины: Fr = VO/VgR, а = ^R/VO, h/R, H/R. Вычислялись коэффициенты суммарных гидродинамических нагрузок (Cx,Cy) = 2(Rx,Ry)/(piV)2R) и форма границ раздела сред.
При численном решении задачи расчетная область рассматривалась в интервале |x/R| < 25. Число узлов на границах раздела и контуре выбиралось равным 500 и 80 соответственно. Шаг интегрирования менялся динамически от Ат = 0.05 до Ат = 0.01 (т = tVO/R — безразмерное время). В интервале 20 < |x/R| < 25 был введен демпфирующий слой для подавления отраженных волн от границ вычислительной области по методике, описанной в [5]. Развитие неустойчивости Кельвина — Гельмгольца предотвращалось при помощи фильтрационной процедуры, разработанной в [7]. Значения производных dGi (si, t)/dsi и углов 9i (si, t) (/ = 1, 2) вычислялись на основе кубических сплайнов [2]. Большое количество узлов на границах раздела приводит к огромным вычислительным затратам при решении систем линейных алгебраических уравнений, полученных дискретизацией интегральных уравнений (9), (10). Существенно уменьшить машинное время позволил метод декомпозиции областей, который приобретает в последнее время все большую популярность в вычислительной волновой гидродинамике [8].
Вычисления контролировались на основе интегрального закона сохранения энергии [5], распространенного на случай многослойной жидкости:
Рис. 1. Волны, вызванные вертикальными колебаниями кругового цилиндра в трехслойной (а, pi* = 0.0148, р2* = 1, H/R = 0.3) и двухслойной (б, pi* = 1) жидкости при Fr = 0.5, а = 2,
h/R = 2.
Рис Pi*
. 2. Волны, вызванные горизонтальными колебаниями кругового цилиндра в трехслойной (а, = 0.0148, Р2* = 1, H/R = 0.3) и двухслойной (б, pi* = 1) жидкости при Fr = 0.5, а = 2,
h/R = 2.
Рис. 3. Гидродинамическая нагрузка кругового цилиндра, совершающего вертикальные (а) и горизонтальные (б) колебания в многослойной жидкости при h/R = 2, H/R = 0.3, pi* = 0.0148,
p2* = 1, Fr = 0.2; 0.5; 1 (кривые 1-3).
E = Pl +/1+1 1 I (S ,t)ji (ßi,t) dsi+ i=i
2
+ 2pi* J (si, t)^oiSi (si,t) dsi + g J y2(si,t) (si, t) dsi
—те /
(x(si,t), y(si,t)) G Li(t), Ф01 (si,t) и ^0i(si,t) (l = 1, 2 — полусумма функций тока и потен-
2
циалов при подходе к границе раздела сверху и снизу). При всех указанных предположениях относительно числа узлов и величины временного шага изменение энергии E в процессе вычислений не превышало 0.5 %.
На рис. 1, 2 (см. выше) представлены поверхностные и внутренние волны, вызванные вертикальными и горизонтальными колебаниями кругового цилиндра. Для рассмотренных параметров поверхностные волны в случаях колебаний в двухслойной и трехслойной жидкости мало отличаются друг от друга. Этот факт объясняется тем, что основной вклад в возмущения на границах раздела вносит контур, а взаимодействие самих границ между собой довольно слабое. Следует отметить, что на границе раздела соленой и пресной воды в случае горизонтальных колебаний новые гребни волн образуются раньше, чем на свободной поверхности. Соответствующие гидродинамические нагрузки имеют осциллирующий характер и приведены на рис. 3.
Заключение
Разработанный ранее численный метод решения нелинейных начально-краевых задач о нестационарном движении контура вблизи границы раздела двух жидких сред распространен на случай движения контура в многослойной жидкости. Рассмотрены горизонтальные и вертикальные колебания кругового цилиндра в трехслойной (соленая, пресная вода и воздушная среда) и в двухслойной (водно-воздушная среда) жидкости. Для параметров, при которых проводился численный эксперимент, существенного влияния дополнительного слоя на гидродинамические нагрузки контура и профили поверхностных волн не обнаружено.
Список литературы
[1] Горлов С. И. Нелинейная задача о волнах, возникающих на границе раздела сред, при одновременных разгонных и колебательных движениях кругового цилиндра. Вычисл. технологии, 3, №6, 1998, 21-29.
[2] Русаков С. В. Методы сплайн-функций в вычислительной гидродинамике. Изд-во Пермского гос. ун-та, 1987.
[3] Степанянц Ю. А., Стурова И. В., ТЕодорович Э. В. Линейная теория генерации поверхностных и внутренних волн. Итоги науки и техники. Сеp. МЖГ, ВИНИТИ, М., 21, 1987, 93-179.
[4] Стурова И. В. Линейная теория генерации поверхностных и внутренних волн локальными возмущениями. Дис. .. . докт. физ.-мат. наук. Новосибирск, 1994.
[5] BAKER G.R., MeirüN D.I., Orszag S.A. Generalized vortex methods for free-surface flow problems. J. Fluid Mech., 123, 1982, 477-501.
[6] Hess J. L. Higher-order numerical solution of the integral equation for the two-dimensional Neumann problem. Comput. Meth. Appl. Mech. and Engng., 2, No. 1, 1973, 1-15.
[7] Longuet-Higgins M.S., Cokelet E. D. The deformation of steep surface waves on water. I. A numerical method of computations. Proc. R. Soc., London, A350, 1976, 1-26.
[8] Wu G. X., Eatock Taylor R. Time stepping solutions of the two-dimensional nonlinear wave radiation problem. Ocean Engng., 22, No. 8, 1995, 785-798.
Поступила в редакцию 1 ноября 1998 г.