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

УДК 519.6

Н. Н. БЕЛЯЕВ, Л. В. ПОКУТНЕВА (ДИИТ)

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ЗАГРЯЗНЕНИЯ АКВАТОРИИ РЕКИ

На 6a3i моделi потенщально! Te4ii' та моделi переносу домiшки запропонована 2D чисельна модель для прогнозування забруднення рiчок. Наведенi результати розрахунк1в на базi розроблено! моделi.

На базе модели потенциального течения и модели переноса примеси предложена 2D численная модель для прогноза уровня загрязнения в реках. Представлены результаты численного моделирования.

The 2D numerical model to simulate the pollutant dispersion in rivers is offered. The model is based on the equation of potential flow and the transport model. The results of numerical experiment are presented.

Введение

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

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

Математическая модель

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

© Беляев Н. Н., Покутнева Л. В., 2009

д2 P д2 P

дх2 ду2

= 0.

(1)

Решение уравнения (1) проводится при таких граничных условиях [2]:

дР 0

- на твердых стенках — = 0 , где п - единич-

дп

ный вектор внешней нормали;

- на входной границе (границы втекания пото-

дР

ка) — = Vn, где Vn - известное значение ско-

дп

рости втекания;

- на выходной границе

Р = Р * (х = const,у) + const (условие Дирихле).

В рамках модели потенциального течения компоненты вектора скорости потока связаны с величиной потенциала скорости зависимостями

и =-

ср;

дх

v = -

дР_

ду

Определив компоненты вектора скорости водного потока на первом этапе вычислительного эксперимента, можно перейти к решению транспортной задачи. На этом этапе для расчета рассеивания загрязнителя в русле реки применяется уравнение переноса примеси, осред-нённое по глубине реки [1, 2, 4]:

дС диС CvC -+-+-

а

N

сС = div ( gradC) -

дх ду XQ (t)5(х - х,. )5(у - у),

чено на базе методов численного моделирования (СББ модели), например, используя конечно-разностные методы. Для расчета течений в областях сложной геометрической формы в настоящее время используется такие подходы:

• применение процедуры коррекции граничных условий на границах;

• применение неортогональных разностных сеток;

• отображение сложной расчетной области на область более простой геометрической формы;

• применение метода МКЭ;

• применение метода маркирования.

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

где с - концентрация загрязнителя в русле реки; и, V - компоненты вектора скорости течения в русле; ц = (ц,х, цу) - коэффициенты турбулентной диффузии; Q - интенсивность выброса загрязнителя; 5(х - х )5(у - у1) -дельта-функция Дирака; xi, у i - координаты источника выброса; с - коэффициент, учитывающий радиоактивный распад загрязнителя; ( - время.

В работах [2 - 4] рассмотрена постановка краевых условий для уравнения переноса.

Формирование вида расчетной области

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

-л ^

Н ^TfifiS г Н

D.QЭ51'Й£ I —II п

Г >

-Я.^ПГга? 1 R ч:

-41 НЧГ-Й2

-Л..ШИ-Ч2 u

н ¡йи.-ц--ji

b.id^r.ej

—M.iME.ttJ I.VBl-n

Рис. 1. Форма расчетной области

Рис. 2. Изолинии концентрации загрязнителя через 18 с после возникновения чрезвычайной ситуации

,=1

Применение метода маркирования позволяет формировать практически любую форму русла водоема.

Метод решения

Для численного интегрирования уравнения для потенциала скорости используется метод установления решения по времени. В этой связи интегрирование проводится для уравнения вида

- на первом шаге расщепления к = п + — :

Ск_£п 1

ц~ 1 + -(11Ск + 1+Ск) + - Ск =

А 2 х у 4 11

1 ы а

= - (М+„Ск + МхСк + М+Сп + МУСп ) + У-^ 5,;

4 у хх хх уу уу ' 4 ''

, 1

- на втором шаге расщепления к = п +

дР д Р д Р

— = —г + —г

дй дх2 су2

(2)

где й - фиктивное время.

При й ^<х> решение уравнения (2) будет стремиться к «установлению», т.е. к решению уравнения (1).

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

п+— Р 2 _ Рп

1,1_Ц

А й

1

1 п+—

Рп+1 _ Р 2

1,1 1,1

п+—

_ Р ■ 2

1,1

1

п+—

Р 2

2

А й

А х

рп+1 рп+1 Р+1,1 _ Р1, ]

1 1 п +— п+— _Р 2 + Р 2

А У2

А х2

рп+1 рп+1 Р1,1+1 _ Р1, ]

А У2

Отметим, что схема является неявной, но расчет значения потенциала скорости Р{ ■ в каждой разностной ячейке на каждом шаге расщепления осуществляется по явной формуле -методу бегущего счета [3]. Это позволяет создать эффективный алгоритм расчета в многосвязных областях и областях сложной геометрической формы, каковыми являются русла рек.

После расчета величины потенциала скорости компоненты вектора скорости рассчитываются по соотношениям:

с=п+—: 4

+ !(ГСк + ГСк) +-Ск = А 2 х У 4 1

1

N

= -(МхСк + М+Сс + МУСк + М+УСс) +У 5,;

4у хх хх уу уу ' 4 1'

, 3

- на третьем шаге расщепления к = п + — ;

с=п+—: 2

С! _ С 1 2

А

11 11 +-(11Ск + ГСк) +-Ск =

2 х у 4 1

1 N

= 1(М_хСс + М+хСк + М_Ск + М+УСс ) + У ^ 5,;

4 4 хх хх УУ УУ > 4 1 '

- на четвертом шаге расщепления к = п +1; 3

с = п +-

4

ск _ с 1

Ай

1 (ГСк

2

у 2(Ь-хСк + Ь+уСк) + ^ Ск =

1 N

= - (МхСк + М+„Сс + МУСс + М+Ск ) + У а 5,.

4 V хх хх уу уу ' 4 1

Пояснение к данным разностным операторам приведено в работе [2]. Из данных выражений можно получить явные формулы для определения неизвестного значения концентрации загрязнителя на каждом шаге расщепления.

ии =

у 1 =

Р _ Р Ах

Р _ Р

Ау

Для численного интегрирования уравнения переноса примеси используется попеременно-треугольная разностная схема [3]. При применении этой разностной схемы процесс решения уравнения переноса расщепляется на четыре шага. Разностные уравнения в операторном виде записываются так [3]:

Практическая реализация модели

На основе рассмотренной численной модели разработан алгоритм расчета переноса примеси в реке. Он состоит в следующем:

• задаётся форма русла реки, положение островов, притоков и т.д.;

• задаётся положение места выброса загрязнителя, интенсивность выброса;

• осуществляется расчет уравнения для потенциала скорости;

• рассчитывается поле скорости водного потока в русле;

1

1

• рассчитывается процесс рассеивания загрязнителя по руслу;

• на печать выводится информация о развитии зоны загрязнения.

На основе построенной численной модели разработан код «RIVER_2» на алгоритмическом языке FORTRAN. Разработанный код был применен для решения следующей модельной задачи. Рассматривается участок акватории реки, где расположены три острова. На самом большом острове имеется отстойник (рис. 4, зона А), где хранятся радиоактивные жидкие отходы, и на этом острове происходит утечка загрязненных сточных вод в акватории реки (рис. 1). Размеры расчетной области: длина -400 м, ширина - 380 м. Скорость потока на входе в расчетную область - 0,2 м/c. Коэффициент диффузии по обоим координатным направлениям равен 0,7 м2/сек. Коэффициент с равен нулю, т.к. рассматривается процесс загрязнения акватории реки за короткий промежуток времени после возникновения аварийной ситуации. Загрязненные сточные воды поступают с острова в реку со скоростью 7 м/с, концентрация загрязнителя в них равна 1 (в безразмерном виде).

Рассмотрим результаты вычислительного эксперимента. На рис. 1 - 3 показана динамика загрязнения акватории реки для различных моментов времени. Из данных рисунков видно, что в акватории формируется зона загрязнения, имеющая сложную геометрическую форму, что вызвано влиянием островов на процесс переноса примеси. Зона загрязнения увеличивается в размерах за счет диффузии и вытягивается в направлении течения реки. Отчетливо видна подзона с высокой концентрацией примеси непосредственно возле места аварийной утечки.

Рис. 3. Изолинии концентрации загрязнителя через 36 с после возникновения чрезвычайной ситуации

В заключение отметим, что для расчета задачи потребовалось около 15 с компьютерного времени.

Рис. 4. Изолинии концентрации загрязнителячерез 63 с после возникновения чрезвычайной ситуации

Выводы

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

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

Численное моделирование процессов загрязнения поверхностных и подземных вод [Текст] / Л. И. Антошкина и др. - Д.: Изд-во ЧП Свидле-ра А. Л., 2004. - 168 с.

Численное моделирование распространения загрязнения в окружающей среде [Текст] / М. З. Згуровский и др. - К.: Наук. думка, 1997. - 368 с.

Марчук, Г. И. Математическое моделирование в проблеме окружающей среды [Текст] / Г. И. Марчук. - М.: Наука, 1982. - 316 с. Основа прогнозирования качества поверхностных вод [Текст]. - М.: Наука, 1982. - 181 с. Савенко, В. Я. Математическая модель механи-

2.

3.

4.

5.

6.

7.

зма поперечной циркуляции в открытых потоках при неизотропных коэффициентах турбулентной вязкости [Текст] / В. Я. Савенко, Е. С. Славинская // Вестник ХГАДТУ. -Вып. 7. - Х., 1998. - С. 50-53. Самарский, А. А. Теория разностных схем [Текст] / А. А. Самарский. - М.: Наука, 1983. -616 с.

Шеренков, И. Л. Прикладные плановые задачи гидравлики спокойных потоков [Текст] / И. Л. Шеренков. - М.: Энергия, 1978. - 240 с.

Поступила в редколлегию 26.03.2009. Принята к печати 03.04.2009.