Численное моделирование процессов эвтрофирования в нижнем бьефе водохранилища-охладителя

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

Том 2, № 2, 1997

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ЭВТРОФИРОВАНИЯ В НИЖНЕМ БЬЕФЕ ВОДОХРАНИЛИЩА-ОХЛАДИТЕЛЯ**

В. М.Белолипецкий Вычислительный центр СО РАН, Красноярск, Россия

В. Б. Туговиков Красноярский государственный технический университет, Россия

А. А. ЦхАй

Алтайский государственный технический университет Барнаул, Россия

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

1. Введение

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

Наиболее эффективно эта проблема решается с применением технологии вычислительного эксперимента [1, 2]. В настоящей работе рассматриваются следующие стадии матема-

*© В. М. Белолипецкий, В. Б. Туговиков, А. А. Цхай, 1997.

* Работа выполнена при финансовой поддержке Красноярского краевого фонда науки, грант № 6Р0111.

тической технологии: математическое моделирование, численные алгоритмы, программная реализация и ряд численных экспериментов, проведенных для проектируемой Березовской ГРЭС-11.

2. Моделирование гидравлических характеристик потока

В большинстве прикладных гидравлических задач при исследовании неустановившихся течений в открытых водотоках ограничиваются одномерной постановкой. Математические модели основаны на уравнениях Сен-Венана [1]:

дЬ + дх ’

+ А (01) = ш(ді -

дЬ + дх\ш) 9Ш\дх ш2С2я) ,

где х, і — пространственная вдоль потока и временная координаты, г — координата свободной поверхности, 9 — ускорение свободного падения, 0 — расход, ш — площадь поперечного сечения, Я — гидравлический радиус: Я = ш/х, здесь ш = /с ь(х,£№, X = Ь(х, 0) + 2 /о у/1 + 0.5(дЬ/д£)2й£, где х — смоченный периметр, Ь(х, £) — ширина реки на расстоянии £ от дна, Н — глубина реки, С — коэффициент Шези — вычисляется по

формуле Маннинга: С = — Я1/6, где п — шероховатость русла.

п

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

Если изменением 0 и г во времени можно пренебречь, то из (1) получим уравнения для установившегося течения в неравномерных руслах:

д0

дх Ъ

(1-г')дх+2иг=

где д — боковая приточность внутри участков, и = 0 — средняя скорость течения,

ш

и2В

г' =--- — число Фруда, и — средняя по сечению горизонтальная скорость потока,

дш 9дшш

—— = -7— , В — ширина свободной поверхности потока, К2 = ш2С2Я. При решении

дхн дх |Ь=сопеі

уравнений Сен-Венана для докритических течений (Г' < 1) в области хн < х < х^, і > 0

задают по одному граничному условию на левой (х = хн) и правой (х = х^) границах:

0|х=хн = /

Н|х=хь = Н

Для нижних бьефов гидроузлов / (і) задается таблично и соответствует проектному или реальному режиму попусков. Функция Н(0) в наблюдаемых створах задается по опорному графику расходов.

3. Моделирование температурного режима нижнего бьефа

Температурный режим нижнего бьефа определяется из упрощенного уравнения переноса для температуры [1]:

dT + dT Sn + 5Дн + Sv B + q

Ж + uax =----------------CP--------Ш + Ш

где x — пространственная координата, направленная вдоль потока, t — время, T — температура воды, Tq — температура путевого притока, с — удельная теплоемкость, р — плотность воды, Sn — суммарный тепловой поток через свободную водную поверхность, и — средняя по сечению скорость течения, Sv = Svw/B, Sv — объемные источники тепла, к которым относится поток, обусловленный переходом части механической энергии в тепловую, S№ — приток тепла от русла.

Для уравнения (5) необходимо задать начальное и граничное условия:

T|t=t0 = T0

T |x=xo = ^

Одним из основных факторов, влияющих на изменение термического режима, является суммарный тепловой поток. Составными частями в Sn входят: тепловой поток, обусловленный поглощением солнечной радиации Sr; поглощение длинноволнового излучения атмосферы Sa; встречное излучение водной поверхности Sw; теплообмен с атмосферой за счет конвекции Sp; теплообмен за счет испарения и конденсации Sis. С учетом направления передачи тепла имеем [3]:

Sn = Sr + Sa — Sw ± Sp ± Sis.

Для расчета поглощенной солнечной радиации используется следующая зависимость [1, 3]:

Sr = 0.94Sr

где N — общая облачность, доли единицы, Sr — солнечная радиация при ясном небе, рассчитываемая по формуле [4]:

sr = 0.66+0.34Y п 0:+0;Л „2Л,s .he^.

(7 — 0.9 + 0.4sin hc\ к sin2 hc

0.1 + 0.4sinhc ) p2(sinhc + 0.107)"

Здесь к — функция от фа, фа = (1.55 + 0.046c)e1'075 — влагосодержание атмосферы. y — склонение солнца, вычисляется по формуле [5]:

где t — номер дня, для которого проводится расчет, отсчитывается от некоторого дня с номером d (d отсчитывается от начала года).

Sa = 4.46  10-13(Ta + 273.15)6(1 + 0.17N2), Sw = 4.47  10-8(T + 273.15)4,

Sp = 0.459/(w)(T — Та),

Sis = / (w)(es — ба),

/ 5278 \

где es = 25.4 exp I 17.62 — —------------- I — давление водяного пара на высоте 2 м, ба =

у Т + 273.15 J

5278

25.4 exp ( 17.62 — —----- I — давление водяного пара на уровне водной поверхности,

V Td + 273.15/

т (Та + 273.15)5278 ,

Td = ------------------- ----- — 237.15 — температура конденсации (точка росы), oC, ф — от-

5278 — 1пф

носительная влажность воздуха, %, w - скорость ветра на высоте 2 м, м/с, /

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

4.3w, если A0v < 0.0148w3 (“не подогретая”

/ (w) = { поверхность),

3.54w + 3.O9(A0v)1, если A0v > 0.0148w3,

Т + 273.15 Та + 273.15 где A0' = 1 — 0.378 6s/P — 1 — 0.378 ба/P (Р — давление' рт ст^>'

4. Моделирование процессов эвтрофирования нижнего бьефа

Для прогнозирования хода процессов эвтрофирования в исследуемом объекте применяется метод имитационного моделирования [6]. При этом динамика содержания планктонного звена водной экосистемы оценивается путем описания достаточно полного цикла биохимической трансформации соединений фосфора (P), как правило, лимитирующего развитие гидробионтов. Моделируется изменение под влиянием комплекса природных и антропогенных факторов концентраций компонентов экосистемы в единичном объеме воды, перемещающемся с русловым потоком.

Модель рассматривает взаимодействие в водной толще растворенных минерального (DIP) и органического (DOP) соединений фосфора, его фракций в составе детрита (PD), биомасс фитопланктона (FP) и бактерий (BP). В целом имитационная модель учитывает те фракции фосфора и процессы, которые имеют первостепенное значение в динамике экосистемы, в развитии фитопланктона и при эвтрофировании водных объектов, а именно — продукцию фитопланктона и потребление водорослями DIP, бактериальную продукцию и минерализацию DOP, метаболические выделения фитопланктона и бактерий, их отмирание с образованием детрита и его последующим разложением, а также обменные процессы биогенным веществом в слое вода — дно.

Уравнения биохимической трансформации соединений фосфора в водной толще имеют вид [11, 12]:

+ dwUC

dt dx

где С « — содержание компонентов водной экосистемы, гР/м3, г = 1 для биомассы бактерий, г = 2 для биомассы фитопланктона, г = 3 для растворенного органического фосфора, г = 4 для растворенного минерального фосфора, г = 5 для взвешенного детритного фосфора; Д(г) — скорость трансформации каждой фракции фосфора, гР/(м3-сут); 7(г) — поток минерального фосфора в водную толщу из донных отложений, гР/(м2•сут); С(г) — скорость поступления компонентов с боковой нагрузкой внутри участков, связанной со смывом и сбросами, гР/(м-сут).

Начальные и граничные условия для уравнения (11) задаются в виде (6), (7).

Скорости изменения концентраций каждой фракции фосфора Д(г) заданы следующими выражениями:

д(1) = (и(1) - £(1) - М(1))С(1), д(2) = (и(2) - Ь(2) - М(2))С(2), д(3) = к(3)С(5) + Ь(2)С(2) - и(1)С(1),

Я(4) = Ь(1)С(1) - и(2)С(2),

д(5) = М(1)С(1) + М(2)С(2) - К(3)С(5) - К8еЙС(5) + Кь^2,

Л2

где и(1) и и(2) — удельные скорости потребления биогенных веществ бактериями и фитопланктоном соответственно, сут -1; ь(2) — удельные скорости метаболических выде-

лений бактерий и фитопланктона соответственно, сут-1; М(1) и

М

— удельные скорости

отмирания бактерий и фитопланктона соответственно, сут-1; К(3) — скорость разложения детрита, сут-1. Скорость убыли детрита из-за седиментации в водохранилище-регуляторе К8еа = К0Л0/Л, где Л — средняя глубина водохранилища в текущий момент времени. Для водотока нижнего бьефа используется постоянное значение К8е^ = 0.1, сут-1 [7]. Величина К характеризует процесс взмучивания и оценивается на этапе калибровки.

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

Рост бактерий и фитопланктона:

и

С

1 + 1 +

С

Метаболические выделения бактерий и фитопланктона:

Ь(1) = пи(1), Ь(2) = г 2 и(2),

где

и

а4 | /1 а2 , /1 а1\

Г1 = ^-4-----------------+ (1--------), Г2 = -+ (1-).

_____+ и(1) а4 ______+ и(2) а2

а4 а2

Смертность бактерий и фитопланктона:

С

М

Скорость разложения детрита:

К

1.2  10-4

1 + 3.0  10-4(е°-351Т)'

Использованы следующие зависимости трофических функций от условий внешней среды — температуры воды и освещенности [10]:

= КТ(е-гх - е-г),

Ке = Ка + Кь  СЫа, СЫс

Г е-КЛ, Г

С

а

I = ^ах 1 + со8 ^ , I

2 ^ 1 ^ \ I ? тах ^ 5

О = 0 3 , 3.68  10-3

Отв = 03 + 1 + 5.25  10-3е°-4°3Т , ТР = ^2 + 1 + 2.8  10-2е°-21Т .

Здесь I — функция освещенности, кал/(см2-сут), Iopt — оптимальное значение освещенности, кал/(см2-сут), 1тах — максимальное значение освещенности, кал/(см2-сут), 1ср — среднесуточное значение освещенности, *п — текущее время суток, — время стояния солнца в зените (12 час).

Для оценки выноса минерального фосфора из дна в водную толщу реализован специальный модельный блок, описывающий трансформацию соединений фосфора в донных отложениях [8]:

дС(6) = 3 - к С(6) д* ктС ,

= ктаС(6) - ^(Кр(С(7) - С(4)) + УфС(4)) - ^,

С(8) = 7рС(7).

Здесь С(6), С(7) и С(8) — концентрации соответственно обменного, минерального и сорбированного фосфора в илах в пересчете на единицу объема порового раствора; а — стехиометрический кэффициент; f — фотопериод; а и Ь1 — пористость и толщина донных отложений; кт — коэффициент минерализации органики в донных отложениях; Кр

— коэффициент переноса фосфатов в илах; Vф — скорость фильтрации; 7р — константа обратимой линейной сорбции фосфора.

В работе использована температурная зависимость параметров процессов в донных отложениях типа Аррениуса:

Вг = 1.09т-2°В°г,

где под г-ми компонентами вектора В подразумеваются соответственно Кр и 7р . Упомянутые массовые потоки фосфора на границе имеют вид:

3(5) = -квсаС (5)П + Кь ^

п

3(4) = кр(С(7) - С(4)) + Уф С(4).

I

г

X

5. Численные алгоритмы

Численные алгоритмы решения рассмотренных задач основаны на методе конечных разностей. Для аппроксимации дифференциальных уравнений разностными вводится пространственно-временная сетка (tn, Xj): tn+i = tn + т (n = 1, ... , N), xi+1 = xi + Ai (i = 1, ... , где т = const — шаг по времени, Ai — величина шага по пространственной переменной, постоянная внутри каждого гидравлического участка. Значения параметров в узлах сетки f (xi,tn) обозначим через f™.

Расчет длины шага по x осуществляется таким образом, чтобы значение Ai укладывалось по длине гидроучастка целое число раз, но было бы близко к заданному для всей расчетной области значению A. Таким образом, чем меньше отношение A к длине наименьшего участка, тем равномернее сетка по x.

Алгоритм расчета гидравлических характеристик потока. Значения гидравлических характеристик в стационарном приближении определяются из решения системы (2) согласно алгоритму:

Qi+1 Qi ___ Qnp

Ai

Ai

zi+1 _ zi

^0i+1

+

1

Bi

Ai

Qi+1 |Qi+1 1 Wi+1Ci+1 Ri+1 £Ші+1

= q

дш \ dxhJ,

Fr_

H+1 -^n/ i+l

2ui+1 Qi+1 Qi

Ai

1 _ Fr

Сначала по (12) рассчитывается значение Q в узлах расчетной сетки по известному значению ^н в створе плотины гидроузла. На следующем шаге по известному значению Q(жi) определяются Л-(жі) и 2^ = гдн + Л,£. Затем, решая уравнение (13) против потока, определяем начальное распределение г, по которому восстанавливается Л.

Алгоритм расчета температурного режима реки. Для определения температурного режима применялся алгоритм, основанный на решении задачи (5)-(7) по неявной схеме бегущего счета первого порядка [1, 9]:

т

тn+1 тn+1 + un+1 i _ i-1

+ ui-1 A

A,- і

Sn + S + Sv B cp ш

+

ш

in (тп _

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

Результаты расчетов температурных режимов рек отличаются от натурных данных не более чем на 10% [1, 9].

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

C(k)n+1 _ C(k)n

т

+

u

n+lC(k)"+1 _ un+1C(fc)n+1

n+^(k) _ ui-1 Ci-1

A

i- 1

дГ+1 +

J <k>-

ш

n+1

1

n

n

В случае наличия притока в каком-либо створе реки средняя по сечению температура вычисляется по формуле:

Т _ Т-^~ + Тпр^пр

Q- + Q

пр

где Тпр и Qпр — температура и расход притока, знак (+), (—) означает величину параметра соответственно ниже и выше створа впадения притока.

Концентрация каждого моделируемого вещества в этом случае также пересчитывается как средневзвешенная, аналогично (16):

См _ С-У + Спр^пР

Q- + Q

пр

где Спр) — концентрация к-го вещества в притоке.

6. Структура комплекса программ

Описанные выше алгоритмы реализованы в виде комплекса программ на 1ВЫ ГО. Тексты программ написаны на алгоримическом языке РОКТИА^77 с использованием технологии модульного программирования. Это позволяет компановать программы в соответствии с особенностями решаемых задач, используя готовые блоки, и дополнять их блоками для моделирования ранее не учитываемых физических процессов [9].

Программы, входящие в структуру комплекса, можно разбить на две основные группы: сервисные и функциональные [9]. К первым относятся программы предварительной обработки исходных данных. Они преобразуют исходную информацию из текстовых файлов и записывают в файлы базы данных, которая функционирует в рамках комплекса программ. Информация в текстовые файлы вводится пользователями комплекса при помощи текстовых редакторов из таблиц. К группе функциональных модулей относятся а) программы начала и продолжения счета для проведения вычислительных экспериментов, требующих больших времен счета на ЭВМ; б) программы для проведения расчетов на системе русел; в) программы краткосрочных расчетов; г) программы для проведения расчетов отдельных физических параметров.

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

Проектнлш вариант 1

Проехтнвариант £

Обводной канал

Рис. 1. Схемы проектных вариантов Верхне-Урюпского гидроузла.

7. Вычислительный эксперимент для

Верхне-Урюпского водохранилища-охладителя

С помощью описанного выше комплекса программ был проведен ряд вычислительных экспериментов для проектируемого водохранилища-охладителя Березовской ГРЭС-11. Необходимость проведения расчетов вызвана тем, что в результате эксплуатации водохранилища-охладителя в нем будет интенсивно развиваться процесс эвтрофирования [10, 11], и теплая, сильно загрязненная биостоком (фито- и бактериопланктоном) и биогенами вода станет сбрасываться из охладителя в нижний бьеф (р. Урюп).

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

Оба проектных варианта схематически изображены на рис. 1.

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

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

а) определение гидравлических параметров руслового потока;

б) расчет изменения вдоль по потоку температуры воды в реке;

в) расчет изменения вдоль по потоку концентраций компонент водной экосистемы.

Если известно подробное распределение граничных условий во времени, то все три

шага расчетной схемы выполняются на каждом шаге по времени. В описываемых вычислительных экспериментах определение гидравлических параметров руслового потока и температуры воды р. Урюп на моделируемом участке нижнего бьефа водохранилища-охладителя велось по среднемесячным данным расходов и метеорологических параметров (кроме солнечной радиации, которая рассчитывается по формулам (9), (10) в зависимости от времени суток), что позволило вынести из цикла по времени расчет гидравлических характеристик. Тестовые расчеты показали, что это не существенно меняет значения температуры и параметров эвтрофирования. Выяснилось, что выведение из временного цикла расчета температуры воды и замена ее значений на среднемесячные приводит к существенным отличиям в результатах расчетов.

Одной из важных проблем при использовании математической модели эвтрофирова-ния явилось определение внутренних параметров (идентификация). Ввиду того что исследуемый объект — Верхне-Урюпский гидроузел — не существует в природе, о процедуре идентификации в полном смысле этого слова говорить нельзя. Поэтому для детализации констант и коэффициентов модели использованы характеристики водоема-аналога. В качестве последнего было взято оз. Балатон, которое испытывает интенсивное антропогенное эвтрофирование, и по нему существует достаточно полный набор необходимых данных [13]. На рис. 2 приведены результаты натурных наблюдений и расчетов по предложенной в работе фосфорной модели водной экосистемы. Как видно, модель достаточно точно воспроизводит режим трансформации компонентов эвтрофной водной экосистемы.

На этом основании при прогнозном расчете процесса эвтрофирования в объектах Верхне-Урюпского гидроузла использовались найденные в результате оптимизационных расчетов значения коэффициентов модели, характеризующие зависимость интенсивности биохимических процессов: потребления, выделения, смертности, метаболизма планктона и др. — в эвтрофной водной экосистеме от задаваемых внешних условий. Расчеты проводились для условий очень маловодного года [15], что позволяет оценить наиболее неблагоприятный случай динамики параметров эвтрофирования в проектируемом водном объекте.

По данным натурных наблюдений и прогнозным оценкам о химическом составе поверхностных вод рассматриваемого участка речного бассейна [10, 11, 15, 16] были выбраны сценарии содержания компонентов фосфорной модели в притоке, поступающем в водохранилище-регулятор, сбросах из водохранилища-охладителя и р. Береш, являющейся притоком р. Урюп.

Сценарий среднемесячных концентраций компонентов в стоке из водохранилища-охладителя выбирался из соображений, заведомо ухудшавших ситуацию. Считалось, что концентрация бактерий увеличивается более чем в 100 раз по сравнению с притоком — р. Урюп до впадения в водохранилище-регулятор [10]. Концентрация фитопланктона в сбрасываемой воде в мае — сентябре считалась равной 40 мгР/л [11]. Концентрация фосфатов в эвтрофных водных экосистемах в период интенсивного развития водорослей обычно бывает невысока ввиду их выедания фитопланктоном [13, 14] ив расчетах полагалась близкой к 0.01 мгР/л.

Расчеты в нижнем бьефе проведены для следующих вариантов (см. рис. 1):

Вариант 1: нижний бьеф начинается от водохранилища-охладителя, из которого в течение всего года производится сброс воды.

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

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

В течение года учитывается разбавляющее влияние притока — р. Береш. Заканчивается нижний бьеф впадением в р. Чулым.

Гидравлические параметры рассчитывались по уравнениям Сен-Венана в квазистаци-онарном приближении для условий р. Урюп и пересчитывались на начало каждого месяца в связи с изменением среднемесячных значений расходов. Изменение температуры и концентрации гидрохимических компонентов определялось по нестационарным моделям с шагом 1 час. Влияние притока р. Береш рассчитывается по формулам (16), (17).

В вычислительных экспериментах по проектным вариантам 2а и 2б наличие водохранилища-регулятора и обводного канала привело к необходимости провести серию расчетов для трех смежных русел. Расчеты выполнялись в следующей последовательности:

определяются гидротермический режим водохранилища-регулятора по стационарному варианту и режим эвтрофирования [14];

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

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

Некоторые из результатов расчетов приведены на рис. 3. Осцилляции на графиках связаны с суточными колебаниями значений исходных параметров.

Основные выводы, полученные из результатов вычислительных экспериментов.

1. При реализации первого проектного варианта биомасса водорослей в замыкающем створе нижнего бьефа в мае — августе может составлять 30-32 мгР/л в сыром весе, снижаясь до 15 мгР/л в сентябре. Ниже будет содержание фитопланктона при реализации вариантов 2а и 2б — соответственно 14-16 мгР/л в мае — августе и до 5 мгР/л в сентябре. Варианты 2а и 2б по степени эвтрофирования нижнего бьефа практически не отличают-ся.Таким образом, полученные в предположении о заданном (40 мгР/л) содержании фитопланктона в сбросах из охладителя значения концентраций водорослей в нижнем бьефе в летний период в соответствии с эколого-санитарной классификацией качества поверхностных вод суши [17] характеризуют водоток нижнего бьефа как сильно загрязненный для случая реализации любого из трех проектных вариантов. Концентрация растворенного органического фосфора будет колебаться в пределах 50-100 мгР/л, минерального — в

Рис. 2. Результаты верификации фосфорной модели по данным эвтрофирования оз. Балатон. • — опытные,------расчетные данные.

Рис. 3. Концентрации компонентов экоситемы в замыкающем створе нижнего бьефа водохранилища-охладителя.

пределах 20-120 мгР/л, детритного (взвешенного) — в пределах 2-14 мгР/л. Реализация вариантов 2 а и 2 б приведет к примерно одинаковым последствиям по степени влияния на процесс эвтрофирования нижнего бьефа и более предпочтительна, чем вариант 1, так как при этом биомасса фитопланктона в замыкающем створе снизится вдвое по сравнению с водохранилищем-охладителем, но тем не менее на два и более порядка может превышать современный уровень эвтрофирования р. Урюп.

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

3. Содержание абиотических форм фосфора (растворенной и взвешенной органики и фосфатов) в нижнем бьефе при реализации всех трех проектных вариантов в течение года будет практически близким и станет изменяться в пределах, обычных для поверхностных вод Сибири.

Авторы выражают благодарность В. Ю. Агейкову за помощь в выполнении рачетов.

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

[1] Белолипецкий В. М., Генова С. Н., Туговиков В. Б., Шокин Ю. И. Численное моделирование задач гидроледотермики водотоков. ИВТ СО РАН, Новосибирск, ВЦ СО РАН, Красноярск, 1993.

[2] Самарский А. А. Математическое моделирование и вычислительный эксперимент. Вест. АН СССР, 5, 1979, 38-49.

[3] HARLEMAN D., Brocard D., NAJARIAN T. A predictive model for transient temperature distributions in unsteady flows. Cambridge, 1973 (MIT, Report, No175).

[4] Казаков А.Л., Лыкосов В.Н. О параметризации взаимодействия атмосферы с подстилающей поверхностью при численном моделировании атмосферных процессов. Труды Западно-Сибирского Регионального НИИ: “Численные методы прогноза погоды”, вып. 55, М., 1982, 3-20.

[5] Абраменков Н. М. Моделирование процесса замерзания шугоносных рек. Труды САНИИ им. В. А. Бугаева. Гидрометеоиздат, вып. 101(182), М., 1984, 3-100.

[6] ЦхАй А. А., АгЕйков В. Ю. Математическое моделирование экосистемы проектируемого водохранилища. В “Применение компьютера в гидротехнике и охрана водных ресурсов”, София, 1990, 428-439.

[7] АйзАтуллин Т. А., Лебедев Ю. М. Моделирование трансформации органических загрязнений в экосистемах и самоочищения водотоков и водоемов. В “Итоги науки и техники. Сер. Общая экология. Биоценология. Гидробиология”, т. 4, ВИНИТИ, М., 1977, 43.

[8] ЦхАй А. А., Тушев А. Н., Щербинина Л. Ю. Модельная оценка выделения биогенов из ложа водохранилища при изменении кислородного режима. В “Математическое моделирование в проблемах рационального природопользования”, Ростов-на-Дону, 1991, 39.

[9] Туговиков В. Б. Численное моделирование гидроледотермических режимов нижних бьефов ГЭС. Дис. ... канд. физ.-мат. наук. Новосибирск, 1991.

[10] Разработка прогноза химического состава воды водоема-охладителя Березовской ГРЭС-2: Отчет о НИР ГХИ Госкомгидромета (научн. рук. Л. Н. Назарова). Ростов-на-Дону, 1990.

[11] Оценка современного состояния и разработка прогноза качества воды водных объектов р. Чулым по гидрохимическим и гидробиологическим показателям: Отчет о НИР ГХИ Госкомгидромета (научн. рук. Л. Н. Назарова). Ростов-на-Дону, 1990.

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

[13] Леонов А. В. Математическое моделирование трансформации соединений фосфора в пресноводных экосистемах на примере оз. Балатон. М., 1986.

[14] Эволюция круговорота фосфора и эвтрофирование природных вод (Под. ред. К. Я. Кондратьева и И. С. Коплан-Дикса). Наука, Ленинград, 1988.

[15] Разработка модели эвтрофирования нижнего бьефа Верхне-Урюпского водохранилища: Гидрологическая записка. Красноярск, 1991.

[16] Разработка прогноза качества вод реки Урюп с учетом вклада в их загрязнение сточных вод предприятий первой очереди строительства КАТЭКа: Отчет о НИР СКТБ “Наука"СО АН СССР (научн. рук. А. М. Мартынова). Красноярск, 1990.

[17] Романенко В. Д. и др. Экологическая оценка воздействия гидротехнического строительства на водные объекты. Киев, 1990.

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