Раздел II. Методы построения дискретных математических моделей
УДК 519.6
А.И. Сухинов, Е.Ф. Тимофеева, А.Е. Чистяков
ПОСТРОЕНИЕ И ИССЛЕДОВАНИЕ ДИСКРЕТНОЙ МАТЕМАТИЧЕСКОЙ
МОДЕЛИ РАСЧЕТА ПРИБРЕЖНЫХ ВОЛНОВЫХ ПРОЦЕССОВ
В работе предложена математическая модель движения волн в прибрежной зоне и построена двумерная конечно-объемная модель выхода волны на берег, учитывающая такие физические параметры, как сложная геометрия дна, турбулентный обмен, трение о дно, динамическое изменение уровня возвышения жидкости, и предложена методика построения дискретных математических моделей на сетках с частичной заполненностью ячеек, а также приведен пример расчета заполнености ячеек Для построенной дискретной модели волновой гидродинамики выполнено исследование устойчивости и получены оценки ограничения решения.
Гидродинамика; волновые процессы; метод баланса; принцип максимума.
A.I. Sukhinov, E.F. Timofeeva, A.E. Chistyakov
DEVELOPMENT AND INVESTIGATION OF THE DISCRETE MATHEMATICAL MODEL FOR COASTAL WAVE PROCESS
MODELLING
The paper presents a mathematical model of the waves in the coastal zone and built a twodimensional finite-volume model release waves on the shore, taking into account the physical parameters such as the complex geometry of the bottom, turbulent exchange, the friction on the bottom, the dynamic changes in the level of elevation of the liquid, and the technique of constructing discrete mathematical models on grids with a partially filled cells, as well as an example of the calculation of the cells. For the construction of discrete models used to study the hydrodynamics of wave resistance and estimates of limiting solutions.
Hydrodynamics; wave processes the method of balance; the maximum principle.
Введение. Предметом данной работы является задача выхода волны на берег, имеющая важное значение для исследования волновых процессов у побережья с целью предопределения строительства сооружений и использования конкретного участка береговой линии.
Постановка задачи. Математическая постанов ка задачи может быть сформулирована следующим образом. Слой идеальной несжимаемой жидкости подходит к откосу с изломом на урезе, сопряженному с ровным дном. В этом случае откос приурезовой области нельзя считать плоским, откос и берег расположены под разными углами. Предполагается, что в начальный момент времени жидкость
находится в состоянии покоя. На некотором расстоянии от берега в точке x = 0 .
границе рассматриваемой области жидкости импульс давления. Требуется определить последующее движение воды.
Для описания задачи волновой динамики жидкости исходными уравнениями являются:
♦ уравнен ие Навье-Стокса:
- коэффициенты турбулентного обмена по горизонтальному и вертикальному направлениям соответственно; g - ускорение свободного падения; р - плотность
жидкости; , Ту - составляющие тангенциального напряжения на дне жидкости;
П - поток вектора скорости через боковую поверхность; Ь - расстояние от поверхности жидкости до дна (глубина жидкости с учетом возвышения уровня) на
боковой границе. Система координат выбрана таким образом, что ось Ох со вмещена с поверхностью невозмущенной жидкости и направлена в сторону берега, ось Оу - вертикально вниз.
, ,
-,
(1)
(2)
♦ уравнение неразрывности:
(3)
Уравнение (3) в случае несжимаемой жидкости примет вид
(4)
Уравнения (1)-(4) рассматриваются при следующих граничных условиях: ♦ на дне области:
Р'( х, у, г) = 0, Уп (х, у, г) = 0, рщ'у (х, у, г) = -тх (г ), рцу'х (х, у, г) = -ту (г), ы'х (х, у, г) = 0, чу (х, у, г) = 0;
(5)
♦ на поверхности жидкости:
ы'п (х, у, г) = 0, у'п (х, у, г) = 0, Р
у(х, у,г) =—-, Р'п(х, у,г) = 0; gр
(6)
♦ на боковой границе:
ыП (х, у, г) = 0, у'п (х, у, г) = 0, Р( х, у, г) = —, (7)
тЬ
где V = {и, V} - вектор скорости движения водной среды; Р - давление; м, л
где
с„ (VI )=|
0,0088, |У|<6,6м/с, 0,0026, V^ 6,6м/с
- .
Дискретная модель. Расчетная область представляет собой прямоугольник. Для численной реализации дискретной математической модели поставленной задачи волновой гидродинамики вводится равномерная сетка:
=|г" = пт,х. = ,у. = ]НУ; п = 0,; г = 0,Nх; ] = 0,Ny;
М,т = Т, N^1 = I,, Nyhy = 1у},
где Т - шаг по времени, \, hy - шаги по пространству, Nt - верхняя граница по времени, Nx, Ny - границы по пространству.
Применяется непосредственная аппроксимация. В результате аппроксимации уравнений (1), (2) по временной переменной, вводя промежуточный слой п + О согласно МЛС-методу и расщепляя уравнения по физическим процессам, уравнения приводятся к следующему виду:
ип+О — ип , , , , , ,
• + иих + \иу = (цих), + (ли у )у , (8)
т
;"+а- V’
т
■ + = (^VX )'х + (^У Уу + 5 , (9)
п+1 п+а п'
и - и Р
т Р
vn+l - vn+е = р'
(10)
(11)
Т р
Продифференцировав уравнения (10), (11) и преобразовав с учетом уравнения неразрывности (4), получим:
Р" + Р" =Р
хх уу
т
(и"+°)х +^),
V /
(12)
Расчет задач гидродинамики по данному методу осуществляется в три этапа. На первом этапе на основе уравнений (8), (9) считается поле скоростей. На втором этапе рассчитывается давление по уравнению (12). На третьем этапе из соотноше-(10), (11) .
. -ники, они могут быть заполненными, частично заполненными или пустыми. Центры ячеек и узлы разнесены на hx|2 и hy|2 по координатам X и у соответственно. Обозначим через о{ . заполненность ячейки (г, .) . Поле скоростей и давление рассчитываются в вершинах ячейки, как представлено на рис. 1. Вершинами ячейки (г,.) являются узлы (г,.), (г — 1,.), (г,. — 1), (г — 1,. — 1).
Степень заполненности ячейки определяется давлением столба жидкости внутри данной ячейки. Если среднее давление в узлах, которые относятся к вершинам рассматриваемой ячейки, больше давления столба жидкости внутри ячейки, то ячейка считается заполненной полностью {о1. = 1^ . В общем случае заполненность ячейки можно вычислить по следующей формуле:
0 = РН (р. і)+р-1. н (р.,)+Р., -1Н ( >)+Р<-и-‘ Н Р). (13)
1,1 Р,
где
[1, х > 0,
Н (X I = - функция Хевисайда.
10, х < 0
(г-1, ])
У
Рис. 1. Расположение ячейки относительно прилегающих к ней узлов
В окрестности узла (г,.) лежат ячейки (г,.), (г +1,.), (г,. +1),
(г +1,. +1) , . 2.
У
(і.І-1)
*°и
,°и+1 (І, І) •° 1+1,3+1
а+1,р
(и+1)
Рис. 2. Расположение узлов относительно ячеек
Вводятся коэффициенты к0, к1, к2, к3, к4, описывающие заполненность областей, находящихся в окрестности ячейки. Значение к0 характеризует заполненность области ^0 : XЕ (—1,х{+^ , у Е (у. —1, у.+^ , к1 - Й1: XЕ (,х{+^ ,
у Е(—1,у.+1) , к2 - ^2: хЕ(у —1,х), у Е (у. — Ру.+1) , кз - ^3:
х Е{Xi—1, х+1 ) у Е (у. , у.+1 ) к4 - ^4: х Е (х,-—1> х+1 ) у Е (у.—1, у. ) . Заполненные части областей 0,т будем называть Кт , где т = 0,4. -
кт
(к ) = ^ (ко ) = 0г. + 0г+1,. + 0г+1,. +1 + 0.+1
(*1 )' і = ^7 2 і+1 . (*2 )' і
(*3 . (*4 ) = °'.7 ^7
П.1 2 4 ч/г.1 2
Для рассматриваемой задачи волновой гидродинамики получена аппроксимация уравнений с граничными условиями.
Построение дискретной модели, учитывающей заполненности ячеек. Для
первой подзадачи волновой гидродинамики (8). (9) граничные условия в общем случае примут вид
иХ(х. у. і) = аи .хи + Ри.х. <(х. У. і) = «V.xv + ^,х.
иУ (х. У. і) = аи.уи +Ри.у . ^ (х. У. і) = «.yV + Л.у . (14)
Проинтегрируем по области Оо уравнение (8) и воспользуемся свойством суммы интегралов. в результате чего получим:
и"+& и"
Ц-----------ёхёу + Циих'х'у + Цvuydxdy = Ц(ри'х)'хйхйу + Ц(циу)'уйхйу. (15)
Оо Т Во А, Оо Оо
Вычислим отдельно каждый из полученных интегралов:
п + & п ,,п+а ,,п ]1п+а - ип
Я—т—(*о ).і Л—т—'^у = (*0).іи'1 т г.1 нхну. (16)
Оо ^о
Второй интеграл в выражении (15) запишется так:
Ции'х'у = Циих'х'у + Ции''х'у - (к1 ). і Циих'х'у + (к2 ). і Циих'х'у.
Оо О1 О2 ^1 ^2
Вычисляя интегралы по областям ^1 и ^2. получим:
\\ш',( )(.і и,+1/2.А (м.' - и'.> ) ++'(*2 )'.і ^АП.‘К (.( - Щ~и ) . (17)
Оо 2
. (15):
Л (р,их )'х 'х'у = Л (р,их )'х 'х'у + Л (рих )'х 'х'у.
Оо О1 О2
В последнем выражении для определенности ПОЛОЖИМ. ЧТО > БОг . выделим из области О1 фрагмент О12. смежный с областью О2. причем Б^ = Б^ (рис. 3):
Л(цих)'х'х'у = Л (цих)'х'х'у + Л (Ми'х)'х'х'у -
Оо О1/ О12 О12 иО2
= ((*1)'.і -(*2)'.і) ЛЛ)'Ах'у + (*2)'.і Л(м'х)'Л'у.
В результате получим:
Л(м'х)'Л'у^ (*1 л.і м
*'+1/2.і
и'+1.і и'.і -(* ) .. и'.і и'-1.і
I*2 А'. ;Мч/2.і ,
-((*1 Л'. і -(*2 Л'. і ). і Ли. хи'. і + А. х)
к.
В случае. если БО > Б^. результат будет аналогичным.
(18)
Рис. 3. Схема заполненности ячеек
Подставим в уравнение (15) выражения (16) - (18) и разделим полученное
выражение на площадь ячейки НхНу. в результате чего получим:
(*о)'., + (*1 Л.і U,+mJ-UїUЬі +(*2 Л.і и-1/2.і +
и"+& - и"
и'.і+1-и'.„ и'.і-и'.і-1
і+1/2
■+(*4 Л ^'.і-1
(*1 Л і М+1/2.і
и'+1.і - и'.і
а и . +^
и. х '. і і и. х
'.і г І+1/2.і к2
і+1/2 к2
(*4 Л. і .. і -1/2 '.і , 2'.і -1 - ((*3 Л. і - (*4 Л. і Л.
аи. уии і + Д. у
к
(19)
О
о
Описание граничных условий. С целью упрощения записи уравнений вводятся «маски» граничных условий: т1, т2 . Параметр т1 принимает значение 1 в случае, если узел (г, .) принадлежит множеству граничных узлов, находящихся в придонной области, в противном случае т1 = 0. Параметр т2 = 1 в случае, если узел (г,.) принадлежит боковой границе, в противном случае т2 = 0. С учетом введенных обозначений построена конечно-объемная модель задачи волновой гидродинамики, представленная следующими сеточными урав-:
♦ для расчета поля скоростей:
♦ для составляющей вектора скорости и. :
г, ]
п+о — п п+а/2 — п+О2 п+О2 — п+а/2
—к»). +—к1 ),. „;+1/2, ,-^2.+—к,),. .и»/,, +
п+о/2 п+о/2 п+о/2 п+о/2
, / 7 \ п иг,.+1 — иг,. . /1 \ п I,. — иг,. —1
+ (кз );,. \. +1/2 ^ + (к4 )г,. Vi,.—1/2 2й =
п+о/2 п+о/2 п+о/2 п+о/2
/7 \ иг+1,. — иг,. /1 \ иг,. — иг—1,.
- (к1 )г,. Мг+1/2,. ^ (к2 )г,. М~т,. ^ +
п+о/2 __ п+о/2 п+о/2 __ п+о/2
, /7 \ иг,. +1 иг,. /7 \ иг,. иг,. —1 ,
+ (к3 А',] . +1/2 ^ (к4 А',] Л. —1/2 ^ +
+ ((*3 ),і -(*4 Л. і )) Л Л. і .
Т
рhy
♦ для составляющей вектора скорости V,. ;:
г,./
v"+& _v" v"+a/2 ___v"+a/2 v"+a/2 _v"+a/2
(*°Л'.і ~1~т і +(* Л'.і и'+1/2.і і 2К 1 +(*2 Л'.і и'-1/2.і 1 2Н - +
х х
Vn+&/2 vn+ст/2 vn+ст/2 vn+&/2
+ (*3 )'., + (*4 Л.і VI.,-1/2 ‘ ‘ 2Л ' 1 -1
у у
_.п+&/2 _»п+&/2 _.п+&/2 _»п+&/2
= (*1 Л'.1 М'+1/2.1 ,~Л2 1 (*2 Л'.і М'-1/2.1 1 Л2 ~ +
+(, ( , Л (_Л ,() . V;.-V-
■ ((*1Л J _ (К Л .і ) Ц Лі + (*3Л.»
рЛху 17^ 4 3/^ ‘.;+1'2 Лу2
Vn+&/2 vn+ст/2
-(*4Л.1 Л.1-1/2 '.1 Л2'.1 -1 +(*оЛ'.15;
♦ для расчета поля давления:
Р - Р Р - Р Р - Р
(*),, - (*.),, /г*1+(*з),/ ■ "ч‘ '•'
-(*4 ),•
,2 V 2 /и ,2
Р - Р . . а(ш2).
і,./ г,./ -1 ' 2> г
К
,
Р
т
/7 \ „ п++ (і \ _ п++ /і \ ^шп++ /і \ ^.п+
(*1 ); / иі+1/2,/ - (*2 А, / иг-1/2,/ . (*3 А- / Уї,./+1/2 - (*4 I / Уі,]
• + •
п++ Л
1/2
-““ ^ с1«“ 1,-с. и-.
♦ уравнениями для уточнения поля скоростей по давлению:
(*»),, /
Ы:
Ы
( *.),,,
рп+1 рп+1 рп+1 рп+1
г+.,. - У_ + (*2) 1 •. - г-1,.
( *о), /
уп+1 - уп++
2,Р
2,Р
рп +1 рп +1 рп+1 рп+1
(*3 ),, + (*4 ), ' - '•”
2к,р
2,ур
Для ответа на вопрос, связанный с адекватным описанием сложных физических процессов взаимодействия волны с берегом и построением их численных аналогов, построенная модель исследована на консервативность.
Исследование дискретной модели. Показано выполнение закона сохранения импульса на дискретном уровне для сеточного уравнения и сохранение потока, что согласуется с физическим процессом. Это свидетельствует о том, что дискретные операторы конвективного переноса соответствуют своим непрерывным .
Найдена погрешность аппроксимации модели. В результате получено, что
порядок аппроксимации системы уравнений равен 0(г + Нх2 + Ну2) , т. е. модель
имеет первый порядок погрешности аппроксимации по времени и второй по про.
Устойчивость разностной схемы. При исследовании у стойчивости на основе принципа максимума задачи волновой гидродинамики (8)—(12) используется представление сеточных уравнений в канонической форме [2].
Уравнение вместе с граничными условиями является частным случаем уравнения вида
^у( Р) = Цр)у(р)- X B(p,о)у(о) = F(p), (20)
р)
где L - некоторый сеточный оператор, р = (х, у ■) - центр шаблона,
Ш'(p) = {(xi+l, у} Мх-^ у} ))х, у7+1 ))х, у7-1)} - окрестность центра
2 ЕШ'(р).
(20), -
вующие системе уравнений (8) - (9), могут быть записаны в виде
В( p, дх+)
В( p, дх-) =
В(P, Чу +)
, / +1/2 , Л ,/+1/2 + »Г,
В(P, Чу-)
А( р)
+ В( P, Чх+) + В( P, Чх -) + В( P, Чу+) + В( P, Чу - ) ,
при этом функция правой части примет вид
Г (*0 )■ ■ /
F1( р) = ------— -(В( Р, Чх+ ) + В( Р, Чх-) + В( Р, Чу+ ) + В( Р, Чу - ) )
+В( Р, Чх + )и;+1, ( + В( Р, Чх - )и; -1, ( + В( Р, Чу + )и;, ■+1 + В( Р, Чу - )и;, ■ -1 +
+ (“*3 I;(-(*4 1, ■ )) К )„ ■
Выражение для правой части можно представить в следующем виде:
^ Р) = - £'«;, ( + (“ “ ( - (*4 );,( )) (Ш1 );,(
Аналогично получается представление для уравнения (9) в канонической форме сеточных уравнений (20), если учесть, что уравнение (9) отличаются от (8) . (9)
2 (*0 ), (и;, (
^( P) =
8 .
г ,■ -,(/ рЬх 4 17г,( 4 0/г,(
(12), -
(14), :
'(*1 );, ( +(*2 );, ( .(*3 );, ( +(*4 );, '
А( p) =
к
+
К
, В(P, Чх+ ):
н2
о, ї (*2),,/ , (*3),,/ , (*4),
В( P, Чх-) = —ГГ", В( p, Чу+) = ^т/, В( p, Чу -) =
Н2
к,
где ? (р )=р
(*,),, -(*,),„■ <+£. , + (*з),, <+Гі/2 (*4),, V™,,,;
((*3 )■■ -(*4 )■■ ^ ) / \ а(т2 V
~ Н ((*3 ),, ( -(*4 ),. ■ ) + )■ ((*■ ). ( -(*2 ) ■ )■
Достаточное условие устойчивости схемы для метода «поправок к давлению» определяется на основе принципа максимума [2] при ограничениях на шаг по пространству:
К <
2р н < 2р
и с (Щ) У V с(Щ)
или Ие < 2И . где Ие = и • I / р - числа Рейнольдса. и - скорость распространения водной среды. I - характерный размер области. р - коэффициент турбу-.
Т-Г и П+1 П+1
Для полей и . V получены сл едующие оценки:
п+1 п
и /—\ < Ъ Т
СЩН ) Я=°
р
+(* X-,-(к4 X-, )рр-ы,,
+ 1 *\ (01-.* 1
С( Щн )
С( Щі
|_ъ1 II П
Vn+1\\ (--V < X т
ЩЩН) ^=0
+ (*0 ), к V
С Щн )
С (щ
. -
тематических моделей на сетках с частичной заполненностью ячеек. а также приведен пример расчета заполнености ячеек. Для построенной дискретной модели волновой гидродинамики получены оценки для норм компонентов вектора скорости. гарантирующие устойчивость по начальным условиям и правой части.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Роуч П. Вычислительная гидродинамика. - М.: Мир, 1980.
2. Сам арский АЛ. Численные методы математической физики / А А. Самарский, А.В. Гулин. - 2-е изд. - М.: Научный мир, 2003. - 316 с.
3. Алексеенко Е.В.,Сидоренко Б.В.,Коп^/нова О.В.,Чистяков А£.Сравнительный анализ классических и неклассичнских моделей гидродинамики водоемов с турбулентным обменом // Известия ЮФУ. Технические науки. - 2009. - № 8 (97). - С. 6-18.
4. Шокин ЮМ. Вычислительный эксперимент в проблеме цунами / Ю.И. Шокин, Л.Б. Чуба-ров, Ап. Г. Марчук, КБ. Симонов. - Новосибирск: Наука. Сиб. отд-ние, 1989. - 168 с.
5. Чистяков А.Е. Об аппроксимации граничных условий трехмерной модели движения
// . . - 2010. - 6 (107). - . 66-77.
6. Сух иное А.И.,Чистяков А.Е.,Алексеенко ЕМ. Численная реализация трехмерной модели
// -
тическое моделирование. - 2011. - Т. 23, № 3. - С. 3-21.
7. Гаделыиин В.К.,Любомищенко Д.С.,Сухинов А.И. Математическое моделирование поля ветровых течений и распространения загрязняющих примесей в условиях городского рельефа местности с учетом к-е модели турбулентности // Известия ЮФУ. Технические науки. - 2010. - № 6 (107). - С. 49-65.
8. Тимофеева Е.Ф. Математическая модель движения волн для водоема с нелинейной функцией рельефа дна // Известия ЮФУ. Технические науки. - 2010. - № 6 (107). - С. 95-102.
Статью рекомендовал к опубликованию д.т.н., профессор В.П. Молчанов. Сухинов Александр Иванович
Технологический институт федерального государственного автономного образовательного учреждения высшего профессионального образования «Южный федеральный университет» в г. Таганроге.
E-mail: [email protected].
34792В, г. Таганрог, пер. Некрасовский, ...
Тел.: 88634310599.
; . - . .; .
Чистяков Александр Евгеньевич
E-mail: [email protected].
Тел.: ВВбЗ4З71б0б.
Кафедра высшей математики; ассистент.
Тимофеева Елена Федоровна
Северо-Кавказский государственный технический университет.
E-mail: [email protected].
355029, г. Ставрополь, просп. Кулакова, г.
Тел.: ВВ65272ВВ5В; +79097583970.
Кафедра прикладной математики и компьютерных технологий; старший преподаватель. Sukhinov Alexander Ivanovich
Taganrog Institute of Technology - Federal State-Owned Autonomy Educational Establishment of Higher Vocational Education “Southern Federal University”.
E-mail: [email protected]
44, Nekrasovskiy, Taganrog, 347928, Russia.
Phone: +78бЗ4З10599.
The head of TIT SFedU; Dr. of Phis.-Math. Sc.; Professor.
Chistyakov Alexander Evgenjevich
E-mail: [email protected].
Phone: +78бЗ4З71б0б.
The Department of Higher Mathematics; Assistant.
Timofeeva Elena Fe’dorovna
North Caucasus State Technical University.
E-mail: [email protected].
г, Kulakov Pr., Stavropol, 355029, Russia.
Phone: +7Вб5г7гВВ5В.
The Department of Applied Mathematics and Computer Technology; Senior Lecturer.
УДК 519.6
А.И. Сухинов, А.Е. Чистяков, Е.А. Проценко
ПОСТРОЕНИЕ ДИСКРЕТНОЙ ДВУМЕРНОЙ МАТЕМАТИЧЕСКОЙ МОДЕЛИ ТРАНСПОРТА НАНОСОВ
В работе освещены вопросы построения пространственно-димерной математической модели перемещения наносов в прибрежной зоне мелководных водоемов, которая используется для численного моделирования динамики аккумулятивного берега. Выполнена дискретизация модели на основе метода баланса, при этом предложенные конечноразностные схемы консервативны, устойчивы и имеют второй порядок погрешности аппроксимации по пространственным координатам и первый по времени. Построены сеточные уравнения для задачи транспорта наносов. Предложен алгоритм расчета коэффициентов и правых частей сеточных уравнений.
Математическое моделирование; транспорт наносов; метод баланса; дискретная модель; сеточные уравнения.