Математическое моделирование экосистемы рыбоводного пруда
К р у г о в о р о т в е щ е с т в. При построении математи-
ческой модели любой экосистемы первоочередной задачей явля-
ется выбор переменных. Это наиболее субъективная и наименее
формализованная часть работы, и определяется она следующими
факторами: цели и задачи исследования; доступная эксперимен-
тальная и теоретическая информация; «обозримость модели».
В нашем случае (Богданов и др., 1990) для оценки рыбохо-
зяйственных возможностей пруда требуется подробное описание
как рационов, так и трофических взаимодействий между различ-
ными видами рыб и между рыбами и другими компонентами эко-
системы. С другой стороны, для описания процессов, происхо-
дящих в экосистеме водоема, требуется достаточно полное пред-
ставление о протекающих в рыбоводном пруду гидробиологиче-
ских процессах.
В пруду выращивается карп, белый толстолобик и буффало.
Последние два вида в модель пока не включены в силу их отно-
сительной независимости от других компонентов экосистемы и
их незначительной биомассы в общем количестве выращиваемых
рыб.
Исходя из кормовой базы карпа – С и толстолобика – Т, в
модель включены следующие переменные: F – фитопланктон, E –
бентос, Z – зоопланктон и B – бактерии.
Для описания круговорота биогенных веществ, способных
лимитировать продукционный процесс, в модель включены: Р –
растворенный минеральный фосфор и N – растворенный мине-
ральный азот.
66
Циклы биогенных элементов замыкаются детритом – D, ко-
торый, кроме того, иногда входит в рацион толстолобика.
Таким образом, мы получили 9 фазовых переменных для
модели. Входными функциями являются климатические факто-
ры: температура воды (ТТ) и интенсивность солнечной радиации
на поверхности водоема (Io). Кроме того, включены 4 управляю-
щие функции, характеризующие внесение искусственного корма
(СО – комбикорм, CV – куколки тутового шелкопряда) и мине-
ральных удобрений (SV – суперфосфат, SE – аммиачная селитра).
Круговорот веществ в модели экосистемы пруда осуществ-
ляется согласно диаграмме потоков, показанной на рис. 14.
Рисунок 14 – Концептуальная модель высокопродуктивного рыбоводного пруда
Предполагается, что такая схема достаточно полно отражает
процессы трансформации веществ в пруду. Поскольку углерод не
лимитирует продукционный процесс, но при этом является ос-
новной по весу составляющей биомассы, то его потребление мы
учитываем косвенно.
Для описания мелководного пруда глубиной порядка 1 м,
эффектами пространственного распределения организмов и ве-
ществ можно пренебречь, поэтому мы строим точечную модель.
Все переменные рассматриваются в виде концентраций (единица
измерения – г/м3 или г/м2). Под концентрацией живых объектов
понимается отношение их суммарной живой биомассы к общему
объему водоема. Моделируется нагульный пруд в течение 5 ме-
сяцев (с 1 мая по 30 сентября).
О п и с а н и е о с н о в н ы х ф у н к ц и о н а л ь н ы х
з а в и с и м о с т е й.
Основные процессы трансформации вещества качественно
примерно одинаковы для большинства пресноводных экосистем.
Поэтому за основу при построении данной модели принимались
модели экосистем водохранилищ (Воинов, Комилов, 1985; 1984;
1986), озерных экосистем (Лукьянов, 1979; Jrgensen, 1983) и ры-
боводных прудов (Svirezhtv, Krysanova, Voinov, 1984).
Переход вещества с одного трофического уровня на другой
описывается S-образными функциями (Воинов, Комилов, 1984).
Скорости процессов потребления и роста определяются количе-
ством доступного субстрата и такими физическими условиями
среды, как температура и освещенность. Предполагается, что ли-
митирование светом и температурой можно задать мультиплика-
тивными членами в общей функции потока вещества:
Qij = Q (t, i, j, TT, I) = ƒj (TT) · gj · Ф(i, j) · (1 – М Вj),
где Qij – функция, определяющая поток вещества из i-й в j-й пе-
ременную (например, из F в Z), t – время, ƒj (TT) – функция, опи-
сывающая лимитирование j-го организма температурой (j = F, Z,
B, E, T или C), gj (I) – функция, описывающая лимитирование j-го
организма светом (j = F), TT – температура воды, I – интенсив-
ность солнечной радиации [0 ≤ ƒ (TT), g (I) ≤ 1], Ф (i, j) – функция,
определяющая выедания i-го субстрата j-м организмом (напри-
мер, биогенов продуцентами или продуцентов консументами
и т. д.), М Вj – коэффициент, определяющий потери j-го организ-
ма на метаболизм. Для некоторых случаев Ф(i, j) = r(i, j), где
r(i, j) – s-образная функция.
Зависимость роста организмов от температуры воды описы-
вается модифицированной функцией Лемана, а лимитирование
роста фитопланктона светом – функцией Стила (Svirezhev,
Krysanova, Voinov, 1984). Потребление биогенных элементов –
углерода, азота и фосфора – в экосистеме синхронизировано в
68
соответствии с законами стехиометрии и принципом лимити-
рующих факторов.
Например, потребление фосфора фитопланктоном будет оп-
ределяться функцией
Ф (P, F) = min [(1/FP) · r (P, F), 1/FP · r (P, F)] (1)
Тогда потребление азота будет пропорционально
Ф (N, F) = FN · Ф (P, F).
Таким образом, для синхронизации потоков азота и фосфора
используется следующая процедура:
1. Вычисление потенциально возможных потоков биоген-
ных веществ по s-образным формулам r(i, j).
2. Определение лимитирующего биогенного вещества по
формуле (1).
3. Пересчет потоков. Потребление нелимитирующего био-
гена происходит со скоростью, определяемой потреблением
лимитирующего биогена.
Учитывая вышеизложенное, можно сказать, что потребле-
ние углерода будет пропорционально
Ф (СС, F) = FC · Ф (P, F),
где FP, FN, FC – стехиометрические коэффициенты фосфора,
азота и углерода для фитопланктона соответственно. Соответст-
вующая доля FC углерода выпадает при разложении детрита. По
разным литературным источникам, FP =1, FN равно от 5 до 16, а
FC берется около 100.
Итак, в соответствии с принципом лимитирующих факторов
потребление биогенных элементов фитопланктоном определяется
формулой:
Ф (i, F) = Ф (P, F) + Ф (N, F) + Ф (CC, F) = Ф (P, F) +
+ FN · Ф (P, F) = (1+FN + FC) · Ф (P, F).
Потребление биогенных элементов бактериями непосредст-
венно из воды описывается аналогично, но с другими стехиомет-
рическими коэффициентами (BP, BN и BC).
Процессы бактериальной деструкции аппроксимируются
уравнениями химической кинетики 1-го порядка, причем предпо-
лагается, что детрит вновь распадается на азот и фосфорсодер-
69
жащие питательные вещества в соответствии с их стехиометри-
ческим содержанием. Здесь уместны аналогии с теорией фермен-
тативных реакций. Некий «субстрат» – детрит переходит в «про-
дукт» – биогенные вещества под действием своего рода «фермен-
та» – бактерий. Часть продукта при этом используется бактерия-
ми для поддержания собственной жизнедеятельности. В итоге
протекание процесса бактериального разложения будет опис ы-
ваться формулой:
r (D, B) = μDB · DS · B/(KS DB + DS),
где μDB – максимальная скорость разложения детрита бактериями
(1/сутки), KDB – соответствующая константа полунасыщения
(г/м2).
В данном случае s (s=2) можно интерпретировать как коли-
чество молекул, «обрабатываемых» одной молекулой фермента
(одним микроорганизмом) одновременно.
На скорости процессов деструкции существенно влияет
температура. Поскольку процессы ингибируются температурами,
практически недостижимыми в условиях естественных водоемов
(+50, +60оС), для описания температурной зависимости в этом
случае используется формула Вант-Гоффа (Воинов, Комилов,
1984).
Таким образом, разложение детрита фосфора определяется
формулой:
QDP = ƒD (TT) · DES ·KPKP+KN+KC+KS· r (D, B),
где ƒD (TT) – температурная функция для процесса разложения,
DES – коэффициент деструкции (0< DES<1), KP, KN, KC – сте-
хиометрические коэффициенты для фосфора, азота и углерода
соответственно, KS – доля детрита, выходящего из круговорота в
виде трудноразлагаемых фракций и захороняемого в седиментах:
QDS = KS · QDP
Разложение детрита до азота и углерода соответственно бу-
дет:
QDN = KN · QDP, QDCC = KC · QDP.
Теперь остановимся более подробно на питании рыб. Из ли-
тературных данных следует, что некоторым рыбам свойственна
70
элективность питания, т. е. питание с переключением. Так, на-
пример, толстолобик предпочитает использовать в качестве кор-
ма фитопланктон. Но когда фитопланктона мало, или не подхо-
дит видовой состав, он переходит на питание заменяющим кор-
мом – бактериопланктоном. Если и этой пищи становится недос-
таточно, то он может переходить на вынужденный корм,
т. е. детрит. Для описания питания с переключением используем
практику моделирования по водохранилищу (Воинов, Комилов,
1986) и рыбоводным прудам (Svirezhev, Krysanova, Voinov, 1984),
где были «использованы выводы, сделанные В. С. Теном, кото-
рый проанализировал опытные данные В. С. Ивлева в зависимо-
сти питания от концентрации кормовых объектов на примере пи-
тания карпа. Интерпретация опытных данных затруднена тем, что
при модельном исследовании принято оперировать такими поня-
тиями, которые не позволяют усматривать сущность процесса.
Выбор их в значительной мере произволен. А экспериментатор
ограничен лишь эмпирически доступными величинами»
(Svirezhev, Krysanova, Voinov, 1984).
Пусть имеющиеся корма – qi (i = 1,2,3…) упорядочены по
предпочтению q1, q2, q3 …, причем q1 – концентрация излюблен-
ного корма. Тену удалось найти зависимость вероятности по-
требления объекта i-го типа ηi от суммарной концентрации объ-
ектов питания
Тогда вероятность η1 при любом наборе q1, q2, …, qn, а веро-
ятность перехода на менее предпочтительный корм ηi, есть функ-
ция концентрации qi – 1, т. е. ηi = (qi –1), i = 2, 3, 4,…
Кривые ηi (q), i = 2, 3, … имеют перевернутый s-образный
вид (рис. 2), причем ηj < ηi при j > i и lim ηi (q)= 1, lim ηi (q) = 0.
q→0 q→∞
Это означает, что переходный режим не четко выражен, а
охватывает полосу значений q, хотя зависимость ηi от q при
больших q кажется сомнительной. Видимо, ближе к реальности
зависимость
η2 (q1), η3 (q1, q2), η4 (q1, q2, q3) и т. д.
71
Рисунок 15 – Питание с переключением
Так, например, для описания переключения с питания из-
любленным кормом – фитопланктоном – на питание бактериями
в модели используется функция
,
где λF – параметр, характеризующий крутизну кривой, mF – ко-
эффициент полуобеспеченности фитопланктоном, т. е. значение
F, при котором η(F) = 0,5.
Тогда потребление бактерий толстолобиком в случае недос-
татка фитопланктона запишется в виде
QBT = ƒT (TT) · min [(r (Fcr, T) – r (F, T)], r (B, T) · η (F)] · (1 – MBT),
где Fcr – критическое значение фитопланктона, при котором тол-
столобик переходит на питание бактериями. Разность двух тро-
фических функций под знаком минимума обеспечивает ограни-
чение скорости роста рыб так, чтобы она была выше той, которая
достигнута на излюбленном корме (при F = Fcr). При нехватке
бактерий толстолобик переходит на питание детритом и т. д. В
этом случае переключение – двухступенчатое, трехступенчатое
и т. д.
72
Потребление толстолобиком фитопланктона представляется
в виде
QFT = ƒT (TT) · r (F, T) · (1 – MBT).
Остались неописанными еще два процесса: метаболизм и
смертность живых организмов экосистемы. Выделение живыми
объектами водоема продуктов метаболизма можно в первом при-
ближении считать пропорциональным суммарному потреблению
кормовых объектов. Выделение продуктов метаболизма и пере-
ход их в детрит в модели представлен следующим образом (на
примере толстолобика):
S1 (T, D) = (Q FT + QBT + Q DT) · MBT / (1 – MBT),
где MBT – коэффициент метаболизма толстолобика, Q FT, QBT,
Q DT – потоки, описывающие потребление толстолобиком фито-
планктона, бактерий и детрита соответственно.
Кроме того, для рыб учитывается зависимость усвоения
пищи от величины рациона. Так, например, при обильном пита-
нии толстолобика пища заглатывается непрерывно и проходит
через кишечник с такой быстротой, что лишь 30-40 % ее усваива-
ется, тогда как при умеренном питании усваивается в 2 раза
больше (Svirezhev, Krysanova, Voinov, 1984).
Таким образом,
S2 (T, D) = RT2 · MBBT / RTmax,
где RT = Q FT + QBT + Q DT, RTmax – максимальный рацион для тол-
столобика, MBBT – суммарный метаболический параметр для
толстолобика.
Смертность живых организмов, как известно из литератур-
ных источников (Воинов, Комилов, 1986; Svirezhev, Krysanova,
Voinov, 1984) зависит от наличия растворенного кислорода в во-
де. В нашей модели предполагается, что наличие кислорода в во-
де не лимитирует процессов жизнедеятельности и находится на
уровне насыщенности. Поэтому смертность организмов и их пе-
реход в детрит задается линейной функцией биомассы или кон-
центрации живых объектов. Естественная смертность рыб в усло-
виях оптимального режима в нагульном пруду должна отсутство-
вать. Но в связи с вероятностью паразитологических и прочих
инфекционных заболеваний рыбных популяций пруда, а также
73
других внешних воздействий, в модели учтена смертность рыб.
Итак, смертность организмов описывается следующим образом:
S (i, D) = CMORi · i,
где CMORi – коэффициент смертности i-х организмов, i = F, Z, B,
E, T или C. Образующийся при этом детрит состоит из углерода,
азота и фосфора, причем лишь азот и фосфор детрита включены в
число переменных модели. Например, QTD поток характеризуется
так:
QTD = S1 (T, D) + S2 (T, D) + S (T, D)
Таким образом, нами представлены основные потоки веще-
ства, необходимые для описания экосистемы рыбоводного пруда.
На основе диаграммы потоков вещества (рис. 15) и с учетом вы-
шесказанных соображений теперь можно написать уравнение
модели. Полностью все потоковые уравнения модели представ-
лены в табл. 17. В табл. 18 приведены формулы для всех потоков
и функциональных зависимостей, использованных в математиче-
ской модели пруда.
Таблица 17 – Потоковые дифференциальные уравнения
модели
Здесь α, β – доля комбикорма и куколки, идущие на потреб-
ление карпа, соответственно.
74
Таблица 18 – Основные потоки вещества и функциональные
зависимости
QDP = fD (TT) · DES · KP · r (D, B) / (KP + KN + KC + KS)
QDN = KN · QDP, QDCC = KC · QDP, QDS = KS · QDP
QD = f (TT) · gF (I) · Ф (P, F) · (1 – MBF)
QNF = FN · QPF, QCCF = FC · QPF
QPB = fB (TT) · Ф (P, B) · (1 – MBB)
QNB = BN · QPB, QCCB = BC · QPB
QFT = fT (TT) · r (F, T) · (1 – MBT)
QEC = fC (TT) · r (E, C) · (1 – MBC)
QBT = fT (TT) · min {Y1, Y2} · (1 – MBT)
Y1 = r (Fcr, T) – r (F, T), Y2 = r (B, T) · η(F)
QZC = fC (TT) · min {Y3, Y4} · (1 – MBC)
Y3 = r (B, C) – r (E, C), Y4 = r (Z, C) · η(E)
QDT = fT (TT) · min [min (Y1, Y5) – min (Y1, Y2)], Y 6 · (1 – M)
Y5 = r (B cr, T) ·η (F), Y6 = r (D, T) · η (F) · η (B)
QFZ = fZ (TT) · r (F, Z) · (1 – MBZ)
QBZ = fZ (TT) · r (B, Z) · (1 – MBZ)
QBE = fE (TT) · r (B, E) · (1 – MBE)
QFD = MBF · (QPF + QNF + QCCF) / (1 – MBF) + CMORF · F
QZD = MBZ · (QFZ + QBZ) / (1 – MBZ) + CMORZ · Z
QBD = MBB · (QPB + QNB+ QCCB) / (1 – MBB) + CMORB · B
QED = MBE · QBE / (1 – MBE) + CMORE · E
QTD = (MBT + MBBT · RT / RTmax) · RT / (1 – MBT) + CMORT · T
QCD = (MBC + MBBC · RC / RCmax) · RC / (1 – MBC) + CMORC · C
RT = QFT + QBT + QDT, RC = QEC + QZC
