CC BY
68
Вычислительные технологии
Том 2, № 2, 1997
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ЭВТРОФИРОВАНИЯ В НИЖНЕМ БЬЕФЕ ВОДОХРАНИЛИЩА-ОХЛАДИТЕЛЯ**
В. М.Белолипецкий Вычислительный центр СО РАН, Красноярск, Россия
В. Б. Туговиков Красноярский государственный технический университет, Россия
А. А. ЦхАй
Алтайский государственный технический университет Барнаул, Россия
Рассматривается задача моделирования процессов эвтрофирования и гидротерми-ки в нижнем бьефе водохранилища-охладителя. Приводятся метематические модели, численные алгоритмы и результаты численных экспериментов для проектируемой Березовской ГРЭС-П.
1. Введение
Температурный режим водохранилищ-охладителей, создаваемых при тепловых электростанциях, существенно отличается от естественного. Высокая температура воды вызывает изменение биохимических параметров водоема. В воде повышается концентрация фито- и биопланктона, изменяется содержание кислорода, затем наступает стадия массовой гибели живых организмов, в результате чего образуются продукты их жизнедеятельности — биогены. Совокупность данных процессов, называемая эвтрофированием, в условиях гидротермического режима водохранилищ-охладителей приводит к существенному загрязнению водоемов. В незамкнутых циклах охлаждения и очистки загрязненная в результате нарушения естественого функционирования водной экосистемы вода сбрасывается в реки или озера, на берегах которых располагаются населенные пункты. Для более эффективного решения задач технического и хозяйственно-питьевого водоснабжения необходимо иметь возможность тщательной проработки различных вариантов экологической защиты уже на стадии проектирования.
Наиболее эффективно эта проблема решается с применением технологии вычислительного эксперимента [1, 2]. В настоящей работе рассматриваются следующие стадии матема-
*© В. М. Белолипецкий, В. Б. Туговиков, А. А. Цхай, 1997.
* Работа выполнена при финансовой поддержке Красноярского краевого фонда науки, грант № 6Р0111.
тической технологии: математическое моделирование, численные алгоритмы, программная реализация и ряд численных экспериментов, проведенных для проектируемой Березовской ГРЭС-11.
2. Моделирование гидравлических характеристик потока
В большинстве прикладных гидравлических задач при исследовании неустановившихся течений в открытых водотоках ограничиваются одномерной постановкой. Математические модели основаны на уравнениях Сен-Венана [1]:
дЬ + дх ’
+ А (01) = ш(ді -
дЬ + дх\ш) 9Ш\дх ш2С2я) , (1)
где х, і — пространственная вдоль потока и временная координаты, г — координата свободной поверхности, 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|х=хн = / (і) (3)
Н|х=хь = Н(0). (4)
Для нижних бьефов гидроузлов / (і) задается таблично и соответствует проектному или реальному режиму попусков. Функция Н(0) в наблюдаемых створах задается по опорному графику расходов.
3. Моделирование температурного режима нижнего бьефа
Температурный режим нижнего бьефа определяется из упрощенного уравнения переноса для температуры [1]:
dT + dT Sn + 5Дн + Sv B + q (T T^ (Г)
Ж + uax =----------------CP--------Ш + Ш (Tq- T >■ (Г)
где x — пространственная координата, направленная вдоль потока, t — время, T — температура воды, Tq — температура путевого притока, с — удельная теплоемкость, р — плотность воды, Sn — суммарный тепловой поток через свободную водную поверхность, и — средняя по сечению скорость течения, Sv = Svw/B, Sv — объемные источники тепла, к которым относится поток, обусловленный переходом части механической энергии в тепловую, S№ — приток тепла от русла.
Для уравнения (5) необходимо задать начальное и граничное условия:
T|t=t0 = T0(x) при Хо < x < xL, (6)
T |x=xo = ^(t)- (7)
Одним из основных факторов, влияющих на изменение термического режима, является суммарный тепловой поток. Составными частями в Sn входят: тепловой поток, обусловленный поглощением солнечной радиации Sr; поглощение длинноволнового излучения атмосферы Sa; встречное излучение водной поверхности Sw; теплообмен с атмосферой за счет конвекции Sp; теплообмен за счет испарения и конденсации Sis. С учетом направления передачи тепла имеем [3]:
Sn = Sr + Sa — Sw ± Sp ± Sis. (8)
Для расчета поглощенной солнечной радиации используется следующая зависимость [1, 3]:
Sr = 0.94Sr(1 — 0.65N2), (9)
где N — общая облачность, доли единицы, Sr — солнечная радиация при ясном небе, рассчитываемая по формуле [4]:
sr = 0.66+0.34Y п 0:+0;Л „2Л,s .he^. (І0)
(7 — 0.9 + 0.4sin hc\ к sin2 hc
0.1 + 0.4sinhc ) p2(sinhc + 0.107)"
Здесь к — функция от фа, фа = (1.55 + 0.046c)e1'075 — влагосодержание атмосферы. y — склонение солнца, вычисляется по формуле [5]:
( 2п Y = 23.5sin ( —(t — d)
где 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пф(Та + 273.15)
носительная влажность воздуха, %, w - скорость ветра на высоте 2 м, м/с, /(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 (i) = R(i)w + J(i)B + G(i), (11)
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) и
М(2)
— удельные скорости
отмирания бактерий и фитопланктона соответственно, сут-1; К(3) — скорость разложения детрита, сут-1. Скорость убыли детрита из-за седиментации в водохранилище-регуляторе К8еа = К0Л0/Л, где Л — средняя глубина водохранилища в текущий момент времени. Для водотока нижнего бьефа используется постоянное значение К8е^ = 0.1, сут-1 [7]. Величина К характеризует процесс взмучивания и оценивается на этапе калибровки.
Далее представлены уравнения, используемые в имитационной модели для описания жизнедеятельности функций микроорганизмов и других важных с экологической точки зрения процессов.
Рост бактерий и фитопланктона:
и(1) к2^тв и(2)
С(1) ’ С(2)
1 + 1 +
С(3) вС(4)
Метаболические выделения бактерий и фитопланктона:
Ь(1) = пи(1), Ь(2) = г 2 и(2),
где
и(1) а1 и (2)
а4 | /1 а2 , /1 а1\
Г1 = ^-4-----------------+ (1--------), Г2 = -+ (1-).
_____+ и(1) а4 ______+ и(2) а2
а4 а2
Смертность бактерий и фитопланктона:
С(1) \ С(2)
М(1) = “2 + “3 Цй). М(2) = «1 й(5) ■
Скорость разложения детрита:
К(3)
1.2 ■ 10-4(в°'351Т - 1)
1 + 3.0 ■ 10-4(е°-351Т)'
Использованы следующие зависимости трофических функций от условий внешней среды — температуры воды и освещенности [10]:
= КТ(е-гх - е-г),
Ке = Ка + Кь ■ СЫа, СЫс
Г е-КЛ, Г
С(2)
а
I = ^ах 1 + со8 ^ , I
2 ^ 1 ^ \ I ? тах ^ 5
О = 0 3 , 3.68 ■ 10-3(е°'4°3Т - 1) о = 02 , 2.2 ■ 10-2 (е°-21Т - 1)
Отв = 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,
(l2)
Fr_
H+1 -^n/ i+l
2ui+1 Qi+1 Qi
Ai
1 _ Fr
(13)
Сначала по (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 (тп _
(l4)
Схема безусловно устойчива. В, и, ш определяются из решения предыдущей задачи. Начальное распределение температуры в потоке находится из решения стационарного варианта исходной задачи.
Результаты расчетов температурных режимов рек отличаются от натурных данных не более чем на 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
(15)
1
n
n
В случае наличия притока в каком-либо створе реки средняя по сечению температура вычисляется по формуле:
Т _ Т-^~ + Тпр^пр (16)
Q- + Q
пр
где Тпр и Qпр — температура и расход притока, знак (+), (—) означает величину параметра соответственно ниже и выше створа впадения притока.
Концентрация каждого моделируемого вещества в этом случае также пересчитывается как средневзвешенная, аналогично (16):
См _ С-У + Спр^пР (17)
Q- + Q
пр
где Спр) — концентрация к-го вещества в притоке.
6. Структура комплекса программ
Описанные выше алгоритмы реализованы в виде комплекса программ на 1ВЫ ГО. Тексты программ написаны на алгоримическом языке РОКТИА^77 с использованием технологии модульного программирования. Это позволяет компановать программы в соответствии с особенностями решаемых задач, используя готовые блоки, и дополнять их блоками для моделирования ранее не учитываемых физических процессов [9].
Программы, входящие в структуру комплекса, можно разбить на две основные группы: сервисные и функциональные [9]. К первым относятся программы предварительной обработки исходных данных. Они преобразуют исходную информацию из текстовых файлов и записывают в файлы базы данных, которая функционирует в рамках комплекса программ. Информация в текстовые файлы вводится пользователями комплекса при помощи текстовых редакторов из таблиц. К группе функциональных модулей относятся а) программы начала и продолжения счета для проведения вычислительных экспериментов, требующих больших времен счета на ЭВМ; б) программы для проведения расчетов на системе русел; в) программы краткосрочных расчетов; г) программы для проведения расчетов отдельных физических параметров.
База данных, функционирующая в рамках комплекса программ, имеет два иерархических уровня. Первый уровень представляют управляющие файлы, создаваемые и корректируемые пользователем комплекса при помощи стандартных текстовых редакторов. В эти файлы в определенной последовательности заносятся имена файлов исходных данных и файлов с выходной информацией. Файлы исходных данных представляют второй иерархический уровень базы данных; в них содержится информация о метеорологических параметрах, морфометрические характеристики русел, данные о схематизации водоемов и водотоков, информация о притоках, константы для расчета трофических функций. Другая группа файлов второго уровня представляет собой набор выходной информации, предназначенной для графической обработки или вывода в виде таблиц. Наконец, третья группа файлов этого же иерархического уровня — файлы межпрограммного интерфейса, через которые передаются данные для продолжения счета либо по временной координате, либо по постранству, при проведении вычислительного эксперимента на системе русел.
Проектнлш вариант 1
Проехтнвариант £
Обводной канал
Рис. 1. Схемы проектных вариантов Верхне-Урюпского гидроузла.
7. Вычислительный эксперимент для
Верхне-Урюпского водохранилища-охладителя (Березовская ГРЭС-11)
С помощью описанного выше комплекса программ был проведен ряд вычислительных экспериментов для проектируемого водохранилища-охладителя Березовской ГРЭС-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 г.
