© 2004 Электронный журнал "Jahrbuch fur EcoAnalytic und EcoPatologic" |
Приведено в приложении к книге "Количественная гидроэкология..."
Опубликовано: Крестин С.В., Розенберг Г.С. Двухмерная модель "цветения воды" в водохранилище равнинного типа // Изв. Самарского научн. центра РАН. 4. № 2(8). 2002. Стр. 276 - 279.
Модель трансформации азотосодержащих веществ Куйбышевского водохранилища
С.В. Крестин
Институт экологии Волжского бассейна РАН, г. Тольятти
Считается, что самым распространенным биогенным элементом в пресноводных водоемах является фосфор, на втором месте находится азот. Поэтому концентрации азота и фосфора в воде часто являются решающими при прогнозировании "цветения воды". Во многих работах (см., например, [Сердюцкая, 1984; Леонов, 1989; Леонов и др., 1991; Леонов, Цхай, 1995]) рассматриваются математические модели трансформации соединений азота, фосфора и кислорода в водной среде, а так же их взаимодействие с гидробионтами.
Ниже рассматривается еще одна оригинальная модель, опирающаяся на модели А.В. Леонова и описывающая трансформацию соединений азота в водной среде.
Модель включает в себя 17 переменных состояния водной среды Ci (i=1-17): соединения азота (биологические – в биомассах гетеротрофных бактерий BN, фитопланктона FN, простейших PRN, зоопланктона ZON; химические – в форме растворенного органического DON, аммонийного NH4, нитритного NO2, нитратного NO3, детритного ND [i=1-9, соответственно]), соединения фосфора (биологические – в биомассах гетеротрофных бактерий BP, фитопланктона FP, простейших PRP, зоопланктона ZOP; химические – в форме растворенного органического DOP и минерального DIP; детритного PD [i=10-16, соответственно]) и содержание растворенного в воде кислорода O2 (i=17).
Для оценки условий трансформации биогенных соединений в экосистеме водохранилища используется, так называемое, нульмерное приближение (или приближение полного перемешивания). При этом заведомо пренебрегают пространственной неоднородностью экологических процессов в водохранилище.
Модель экосистемы водохранилища задается уравнением общего вида:
d(CiW)/dt = WRi +
Σ Qj qji – Q Ci + Yi ,где Ri – скорость биохимической трансформации соответствующего соединения Ci, г/(м3 сут); W – объем водохранилища, м3; t – время, сут; Qj и qji – расходы реки и боковых притоков и концентрации компонентов в них, м3/сут и г/м3, соответственно; Q – расход пропуска из водохранилища, м3/сут; Yi – скорости выноса биогенных веществ в водоем со дна и с водосбора (без учета притоков; для O2 эта переменная описывает реаэрацию и поглощение донными отложениями, г/сут.).
Изменение биомассы гидробионтов характеризуется выражением
Ri = (Ui – Li – Si) Ci – Gi ,
где i = 1-4 и 10-13 для компонентов N и P, соответственно; Ui, Li, Si – удельные скорости потребления веществ, выделения продуктов метаболизма и отмирания гидробионтов, соответственно, сут-1; Gi – скорость выедания гидробионтов согласно схеме их трофических взаимодействий, сут-1 (см., например, [Уморин, 1983].
Для каждой группы гидробионтов Ci задается ассортимент Pooli взаимозаменяемых соединений N и P:
Pooli =
Σ dij Cj ,где dij – коэффициенты предпочтения в потреблении субстратов,
Σ dij = 1, 0 ≤ dij ≤ 1.Общее выражение для удельных скоростей потребления задается уравнениями:
где Ki – максимальная скорость потребления субстратов сут-1, скорректированная по температуре и освещенности.
Скорости потребления Uij гидробионтами Ci отдельных соединений Cj задаются:
Нетрудно видеть, что удельные скорости потребления Ui представляют собой суммы скоростей потребления отдельных Uij веществ:
Ui =
Σ Uij .Скорости метаболических выделений гидробионтов Li и их смертности Si составляют определенную долю удельных скоростей потребления и задаются уравнениями:
Li = Ri Ui ,
Ri
= α i1Ui / (1 + α i2Ui) + 1 – α i1/ α i2 ,Si= Vi1 + Vi2 Ci / Ui ,
где αi1, αi2, Vi1, Vi2 – внутренние параметры модели.
Выражения для скоростей трансформации остальных соединений Ci, где i=5-9 и 14-17, имеют вид:
R5 = K5 C9 – U1 C1 + C3(L3 – U3,5) + L2 C2 + L4 C4 ,
R6 = L3 C3 – U2,6 C2 – K7 C6 ,
R7 = K7 C6 – U2,7 C2 – K8 C7 , R8 = K8 C7 – U2,8 C2 ,
R9 = C3(S3 – U3,9) + S2 C2 + S1 C1 – K5 C9 + C4(S4 – U4,9) ,
R14 = K6 C16 + L11 C11 + L13 C13 – K9 C14 – U10 C10 + C12(L12 – U12,14) ,
R15 = K9 C14 – U11 C11 + L10 C10 ,
R16 = C12 (S12 – U12,16) + S11 C11 + S10 C10 – K6 C16 + C13(S13 – U13,16) ,
R17 = p1U1C1/(1 + p2U1) – H1 K7 C6 – H2 K8 C7 – H3(K1 C1 + L2 C2 + L3 C3 + L4 C4) ,
где Кi – константы трансформации веществ, Нi – стехиометрические коэффициенты, рi – константа фотосинтетического процесса.
При описании кислородного режима учитываются температурные зависимости скоростей потребления O2 при нитрификации и его выделения при фотосинтезе. В расчетах оцениваются затраты на окисление продуктов метаболизма гидробионтов и поглощение O2 донными отложениями. Использована зависимость коэффициента реаэрации водной среды от скорости ветра, предложенная для водохранилищ [Bahks, 1975, цит. по: Бреховских, 1988] .
Значения максимальной скорости потребления веществ гидробионтами Ki0 корректируются по температуре и освещенности:
Ki = Ki0 RTi RIi ,
где i – индекс рассматриваемого гидробионта (изменяется от 1 до 4); Ki0 – оптимальное значение скорости потребления веществ; RTi и RIi – коэффициенты коррекции по температуре и освещенности соответственно.
Для описания температурной коррекции максимальных скоростей потребления веществ использована экспоненциальная функция:
где T – температура воды, 0C; Aik – константы.
Температурная коррекция скорости трансформации детритных компонентов N и P, соответственно в DON и DOP, задается в модели линейной зависимостью:
K5 = K50 T , K6 = K60 T ,
где K50 и K60 – значения параметров при T = 1 0C.
Температурная коррекция скоростей трансформации неорганических химических соединений осуществлялась в соответствии с зависимостью
Ki = Ki0 Ξ (T-20) ,
где i = 7-10 при описании химической трансформации NH4 до NO2, NO2 до NO3, DOP до DIP; Ξ – температурный коэффициент, равный 1.05; Ki0 – значения параметров при 20 0C.
Зависимость коэффициента коррекции RI2 от условий освещенности имеет вид [Леонов, Цхай, 1995]:
RI2 = (e / Kl h0) [ exp(-rx) – exp(-r1) ] ,
r1 = I / I0 , rx = r1 [ exp(-Kl h0) ] ,
K1 = Ka + Kb C4 ,
I = Ia [ 1+ cos(2 π (t – tp) / f ] ,
где t – время суток, Ia – среднесуточная освещенность, f – фотопериод, tp=12 ч, I0=350 кал/(см2 сут); Ka, Kb и h0 – внутренние параметры модели.
Для бактерий, простейших и зоопланктона коррекция соответствующих коэффициентов по освещенности не проводилась и RI1 = RI3 = RI4 = 1.
Величина нагрузки Yi в общем случае может быть представлена в виде:
Yi = Ji Ω + ω Ei ,
где ω – длина (в м), Ω – площадь зеркала водохранилища (м2); Ji – поток вещества через межфазные поверхности, г/(м2 сут); Ei – путевая нагрузка, связанная с поступлением с берегов, кроме учитываемых притоков, г/(м сут).
Для описания потока оседающего детрита используется следующее выражение:
Jiout = -k Ci , i=9 и 16 ,
где k – скорость осаждения (гидравлическая крупность), м/сут.
В след за А.В. Леоновым [1989] считается, что величина потока взвешенного материала со дна в водную среду прямо пропорциональна скорости ветра V, м/с, и обратно пропорциональна глубине водохранилища H, м:
Jiin = bi V / H, i=9 и 16 ,
Значение k = 2.0, b9 = 6.4*10-3, b16 = 0.1 определены при параметризации модели.
Поток O2 через свободную поверхность водохранилища представлен в виде
Jiin = kr (C0i – Ci) , i=17 .
Здесь коэффициент переноса через границу воздух-вода, имеет вид
kr = 0.728 V1/2 – 0.317V + 0.037V2 ,
а C017 – содержание растворенного кислорода в воде при насыщении, рассчитываемое в зависимости от температуры [Леонов, 1989]:
C017 = 14.61996 – 0.4042T + 0.00842T2 – 0.00009T3.
В практике моделирования кислородного режима водных объектов процесс поглощения O2 донными отложениями обычно воспроизводится в виде реакции нулевого порядка, так как в широком диапазоне концентраций скорость поглощения не зависит от содержания O2 в воде [Бреховских, 1988; Мизандронцев, 1991]. Поэтому учитывается лишь температурная зависимость соответствующего потока O2:
Jiout = J0i (1.09)T-2 , i=17 .
Здесь J017 – поток O2 при T=20 0C, который соответствует конкретному типу донных отложений и характеру окислительно-восстановительных реакций вблизи границы вода-дно. Величина этого потока принимается неизменной и определяется при параметризации модели.
Описанная выше модель достаточно громоздка. Поэтому для иллюстрации проанализируем меньшую по объему модель трансформации только азотосодержащих соединений.
В работах [Леонов, 1989; Леонов, Цхай, 1995] изучаются органический, аммонийный, нитритный, нитратный и дитритный азот, а так же взаимодействие этих форм азота с бактериями, фитопланктоном, простейшими и зоопланктоном.
В нашей модели ради простоты ограничимся рассмотрением взаимодействий только трех видов азотосодержащих соединений (NH4, NO2, NO3) и фитопланктона. Здесь применен камерный подход, состоящий в том, что в водоеме исключается пространственная неоднородность и все компоненты считаются идеально перемешанными. Учтены следующие процессы взаимодействия. Аммоний солевой частично поглощается фитопланктоном, частично перерабатывается в нитриты. Некоторое количество нитритов преобразуется в нитраты. Нитраты поглощаются фитопланктоном. Кроме взаимодействия их друг с другом учитывается отток и приток воды и связанный с этим перенос компонентов. Все выше сказанное позволяет записать следующую систему дифференциальных уравнений:
N'1 = QПР N1ПР/W – QОТN1/W + (S + M – (1 – δ)μ)νB – K1N1
N'2 = QПР N2ПР/W – QОТN2/W + K1N1 – K2N2
N'3 = QПР N3ПР/W – QОТN3/W + K2N2 – δμνB ,
где N1, N2, N3 – концентрации, соответственно NH4, NO2, NO3; K1, K2 – коэффициенты не консервативности для аммонийного и нитритного азота; S, M, μ – скорости осаждения, смертности и роста фитопланктона, ν – стехиометрический коэффициент, характеризующий содержание азота в единице сухого веса биомассы фитопланктона; δ – доля нитратного азота в потреблении азота фитопланктоном; QПР, QОТ – количество втекающей и вытекающей из камеры воды, W – объём камеры, B – концентрация фитопланктона (см. [Рекомендации по прогнозированию.., 1984]):
B = -b/a + (B0 + b/a)exp(at/ε1) ,
где:
a = ε
1(ρ+S+M+Q/W) , b = μ/H – ε0(ρ+S+M+QОТ/W)+QПРε1ВПРФ/W ,ρ, S, M
– скорость респирации, осаждения и гибели фитопланктона, соответственно; ε0 – коэффициент поглощения света водой при отсутствии фитопланктона; ε1 – коэффициент поглощения света в результате развития фитопланктона; В0 – начальная концентрация фитопланктона; ВПРФ – концентрация фитопланктона в поступающей в камеру воде; H – средняя глубина камеры; t – время.Введем обозначения:
β = B0 – α; γ = (ρ+S+M+Q/W) .
Тогда:
b/a = - α , Β0 – b/a = B0 – α = β , a/ε1 = - γ и B = α + βexp(- γt) .Решая уравнения системы, после несложных выкладок получаем:
N1 = A1 exp(-a1t) + E1 + F1 exp(-γt) ,
где: a1 = K1 + QОТ/W ,
N01 – начальная концентрация аммонийного азота.
N2 = A2 exp(-a2 t) + E2 + D2 exp(-a1t) + F2 exp
(-γt) ,где: a2 = K2 + QОТ/W ,
N20 – начальная концентрация NO2.
N3 = A3 exp(-a3 t) + E3 + D3 exp(-a1 t) – A2 exp(-a2 t) + F3
exp(-γt) ,где: a3 = QОТ/W ,
N03 – начальная концентрация NO3.
Расчеты производились по однокамерной модели с шагом в 10 дней на общем промежутке времени в 1 год. Коэффициенты, при которых производились вычисления, приведены в следующей таблице:
Таблица 1.
Коэффициенты, при которых производились расчеты
K1 |
0.207 |
H |
30.0 |
K2 |
10.0 |
BO |
0.000305733 |
Ρ |
0.5 |
N01 |
0.001 |
M |
0.5 |
N02 |
0.00002 |
S |
0.6 |
N03 |
0.14 |
μ |
0.02 |
QОТ |
0.503 |
ν* |
0.08 |
QПР |
0.503 |
ε 1 |
9.0 |
W |
6.642 |
δ |
0.1 |
Параметры, отмеченные *, взяты из работы [Рекомендации по прогнозированию.., 1984], остальные получены экспериментально в ходе расчетов.
На вход камеры подавались концентрации, взятые из наблюдений, произведенных в 1988 г. на Приплотинном плёсе Куйбышевского водохранилища [Паутова, Номоконова, 1994]. В ходе расчетов выявилось хорошее совпадение данных, взятых из наблюдений с расчетными (см. рис. 1 - 4). На каждом рисунке оба графика нормированы на единицу, то есть максимальное значение как экспериментальных, так и теоретических данных равно единице.
1) Концентрация NH4
2) Концентрация NО2
3) Концентрация NО3
4) Концентрация фитопланктона
Рис. 1-4. Зависимость нормированных концентраций ингредиентов от времени t
(1 – результаты моделирования, 2 – наблюдения)
Графики на этих рисунках состоят из достаточно сложных сочетаний "впадин" и "поднятий", как для кривых наблюдений, так и для расчетных кривых. Сам характер совпадений поведения экспериментальных и теоретических кривых свидетельствует о высокой степени количественной адекватности рассмотренной модели. Сезонная динамика концентрации фитопланктона (рис. 4), кроме того, демонстрирует качественную картину трех пиков "цветения": весеннего (незначительный), летнего и осеннего.
Дальнейшее совершенствование модели может идти, например, в направлении "корректировки" превышения расчетных значений над экспериментальными "поднятиями" и "впадинами".