Модель оптимального роста биомасс растительного покрова

УДК 574.4+631.43+631.6+626.8

А. В. Воротынцев Вычислительный центр им. А. А. Дородницына ФИЦ ИУ РАН

Модель оптимального роста биомасс растительного

покрова

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

Ключевые слова: моделирование, растительный покров, распределение ассимилятов.

A. V. Vorotyntsev Dorodnicyn Computing Centre, FRC CSC RAS

Optimal growth model of plant cover biomass

We develop a biomass growth model for four plant organs coupled with the model of heat and moisture transfer in the plant cover and soil, with account taken of the energy balance. Using growth functions, we study the problem of rational distribution of photosynthesis assimilates and optimal biomass growth to maximize the final value of the biomass reproductive plant organ.

Key words: modeling, plant cover, assimilates allocation.

1. Введение

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

© Воротынцев Л. В., 2017

© Федеральное государственное автономное образовательное учреждение высшего образования «Московский физико-технический институт (государственный университет)», 2017

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

Приведем систему уравнений модели переноса, описанной в [1—5]. В слое 0 ^ х ^ Щ растительного покрова (РП) и в корнеобитаемом слое -Из ^ х ^ 0 почвы рассматриваются две подсистемы уравнений (1) - (4) для температуры воздуха Та, листьев Тг, концентрации водяного пара в межлистном воздухе да, в устьичных полостях листьев дг:

Jтa = -Срка дХ, 3да = -ка , 0 ^ X ^ Иг,

дТа дЗта , , дда д,1 ча , , /гЛ

= ^^ + Т д = ^^ + ^

¡тг + XIдг = дХ,

и две подсистемы (5) - (7) для температуры Тз и водного потенциала фз < 0 почвы:

3тз = -Срктз^^, 3-з = -кфз, -Щз ^ X ^ 0,

дх дх

дТз = _ д3тз дфз = _ д3фз .

Сз дЬ дх , Сф дЬ дх Тфз

¡фз = ЗдГРк + Б кБ к

а также подсистема (8) - (10) для концентрации углекислого газа СО2:

два д

дСа дЗса

3 = _ к — (8

а — дх ,

+ 1е1, (9

дЬ дх ¡с1 = ВСБ1С -

с краевыми условиями (11) - (14):

дЬ дх

¡с1 = ВсБг С - Са)р1,

3та = СрВа

3са — Ба (Са Са);

3та = СрБтз

3са — 0;

3та + Х3да - ¿тз = Я

Тз = Т0, фз = фЗ, х = -Из.

Транспорт воды в растениях и испарение ими влаги 3дг (транспирация) регулируются водным потенциалом листьев фг, их устьичным сопротивлением т^ и описывается нелинейными выражениями:

гщ , _ _

3д1 = ¡дг йх, 3д1 = Бк

Jo

где фз = фзРк йх, щ = фг/фт, 1/Бд = 1/Бт + Тзг, 1/Бк = Ткс + 1/(БкБк),

3-нк

Тзг = Тзз

3дг = 0 при ф3 ^ фт ^ 0.

Подсистема уравнений, описывающая распределение углеводов, преобразование поглощенного СО2 в углеводы, затем в биомассу органов и рост биомассы листьев (г = I), корней

(г = к), репродуктивных органов (г = К) и стеблей (г = в), имеет вид [6], [7]:

Ашг = рг(г)Ф Аг - Кг(Т )шгАг - ИгАшг, или

Ф' ' т-^ Шг = Рг(г)-^ - ЩШг, = 1>Рг&) ^ 0,% € (1,к,К,в),

1+ К

где Ф = БгтФо, Кг = Кг(Т)(1 + К), Лг = е^а(Тг), Фо = П = БгБ-1 + Б-1,т = Л1Г1 + тс.

Л1

(аТ^)-1 + ГсС- 1;

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

Рис. 1. Схема переноса тепла, влаги и CO2 в РП

В системе (1) - (20) функции ,1та(х,г), ^а(х,г), .1са(х,г), .1тз(х,г), '1фа(х, г), 1дг(г), Та(х,г), т(г), Т3(х,г), да(х,г), са(х,г), фг(г), фа(х,г) являются искомыми. Остальные величины - заданные функции и константы; Кц(х,Ь) — заданная поглощенная слоем (0,х) длинноволновая и коротковолновая радиация; К(0, г) — радиационный баланс у поверхности почвы; К[(г) = Еи(И1,г) — радиация, поглощенная покровом, Ка = К(0,г) + ,1та

— энергия, поглощенная верхним слоем почвы; Бг(г), Б к(г) — площадь поверхности соответственно листьев и корней РП на единицу поверхности почвы; Бгрг(х,г), БкРк(х,г) — плотность листовой и всасывающей корневой поверхности РП; Бд(фг), Бс(фг) — коэффициенты проводимости для водяного пара и С02 на границе лист-воздух. Т°, д°а, са, и Та0,

— функции времени, измеряемые на высоте метеобудки На > Нг над покровом и глубине -На почвы; Q(t) — осадки; р(Т) — насыщенная влажность воздуха при температуре Т; ср, са — теплоемкости воздуха и почвы; ка(х), кта(х) — коэффициенты турбулентной проводимости воздуха и теплопроводности почвы, кфа(фв) — влагопроводность почвы; сф = сф(фа); X — теплота парообразования; сг — концентрация СО2 в хлоропластах листьев, Бт, Бк, Ба, Бта, Бда, — постоянные коэффициенты проводимости.

Продукты Ф Аг фотосинтеза распределяются по органам растения долями рг(г) и используются для роста Ашг биомасс органов; при этом часть углеводов расходуется на дыхание поддержания Кг(Тг)шгАг и дыхание роста КгАшг. Температурный фактор а(Тг) и константы £ф, а, тс определяют влияние потока фотосинтетически активной радиации (ФАР) и концентрации сг на скорость Фо преобразования СО2 в ассимиляты единицей поверхности листа.

Дифференциальные уравнения систем связаны алгебраическими выражениями (4), (10), описывающими балансы энергии в межлистном воздухе при х > 0 и на поверхности х = 0 почвы.

Введем обозначения:

йа

аа = р(Т0) - яЧ, Р(Т) = р(Тг) + 5(Т - Тг),

к

[щ 1 Ьх 1 1 + 1 [щ ас 1 5х

Бак Уо ка(СУ Б'ак срБак БНк Бак Б'ак Л к(С )' Бт а сРБа

± = ± + __L= 5х 6х ± = = ^ +

Б а Б а БТа' БТа срБтв Бт срБт Бг, Бд Бя Бг, Б\ Б'я Бт, 1 1 1 , ^Б' „_л Бв

^-г- ,Ьг = О -г- ,Ьв = О-1—в Бе = Ба + А + Б3

1—1 г ь _о— 1 Б

Бв Бдз ' Б^'" ^ Бт ,в = БТз

Б е = Ба + Бт + Бв.

Отметим, что введенные величины положительны, кроме, возможно, Ог и Оа. Очевидно,

Б[/Бт < 1, Б'а/Б'Та < 1.

3. Приближенная модель переноса влаги и тепла

Для модели (1) - (11) получим приближенную модель при условиях

Бе/Бт < 1, Б/Б'Та < 1, ЛБе/Бнк < 1.

Здесь запись а < 1 математически эквивалентна а ^ 0; она используется из-за своей большей физической содержательности. Введем обозначения:

в =

й =

1

З-а = Ъ^а, 4а = Х'^О + Б'Х;

4а = ЪТ +

Я- = Яа + З'/Б а, Я+ = 4° + Ъ-1

Т- = Т0 + 3'а/

В [5] доказано, что в условиях (22) справедлива система

дв-/дЬ = -д,1ф8/дх - ЗдРк, зф = з+а/

-^|*=0 + Здо = Я(Ь), Зф = въз+а/(1 + в), X = 0; (31) ^ = ф8,х = -Н3;

где

Зфз = -кф3дфз/дх, Сф = дв-/дфэ < 0,

в = -ё--1-=--+ Акс^~л-—; Акс = -;-р.-,

Э 1 - Ц8 Эк 1 - И8 -фтБк

при 0 < —8 < 1, в = ж при 1 < —8 .

Проинтегрируем (30) на отрезке [-Нк, 0] с учетом Зда = Зф + Зд0:

Г -

+ Зда - Я

Чтобы функционально в (34) связать и> и ф8 и сделать нижеследующие оценки, примем два допущения. Во-первых, допустим, что для производной Сф = двш/дф8 основной гидрофизической характеристики (ОГХ) вш = вш (ф8) с достаточной здесь точностью можно принять Сф = двш/дф8 для зависимости средней влажности вш от среднего потенциала ф8 по слою [-Нк, 0]. Во-вторых, допустим простейшую зависимость ф8 = Ае-а™в™. Тогда (34), (33) преобразуется к уравнению вида

(

1 + Ъв

1 + в''

Зда(—8) = З+а -, + о ' при 0 < —8 < 1' Зда(—8) = ЪЗ+а' при 1 < —8 (36

Переменную —8 = ф8/фт, пропорциональную среднему водному потенциалу слоя [-Нк, 0], будем называть приведенным водным потенциалом корнеобитаемого слоя почвы.

4. Модель фотосинтеза растительного покрова

Получим приближенную модель фотосинтеза растительного покрова, предположив, что в системе уравнений (8) - (10): а) концентрации са, С\, фотосинтетически активная радиация незначительно изменяются внутри растительного покрова; б) покров находится в состоянии квазистационарного равновесия в отношении С\, т.е. дС\/дЬ ~ 0.

С учетом сделанных предположений проинтегрируем (9) - (10) по высоте покрова на отрезке 0 ^ х ^ Н1:

г Н

Зса 1х=Н + Зса |х=0 + / БсБг (С1 - Са)'[П(х = 0.

Jо

С учетом (10) - (12) получим —Ба(са — + БсБг(сг — са) = 0. Выразив отсюда са — сг и подставив в (10), получим выражение для потока СО2, поглощаемого листьями в единице объема РП:

БаБсБг о

1с1 = — тл , п ? (са — 01 )р1'

Па + БсБг

Интегрируя это выражение на отрезке 0 ^ х ^ Нг, т.е. по высоте РП, получим величину БгФ поглощенного СО2 столбом РП с единичным основанием за единицу времени. Отсюда величина СО2, поглощенного единицей поверхности листа, будет равна Фо = (С°а — С1)/т\. Подставляя сг = са — т1Ф0 в (19), получим для Ф0 уравнение

Ф

о =

Г1Ф0 — (с°а + та!? )Фо + Лса а!? = 0,

А1 Г /-

са + га!? — у (с0а — та!?)2 + 4т са!?с

2

Л1с0а V = О?.

1+ V + л/(1 — V )2 + 4утс/т с

Легко видеть, что Ф0 & Л1соа/т при V ^ то, например, при I? или т ^ то. Однако в дальнейшем для получения явных оценок в условиях дефицита почвенной влаги будем использовать следующие более сильные условия:

1 =

д 1 — Г,

(1 — Лг Г,) +

Лкс 'да Б1

т30 Б'а + П Б к В самом деле, (38) прямо следует из (33), (23), (22):

в = (Па + Б,)(1/Б'д + 1/Пт) & Па + Б,)/(БдБг). Пусть ттез < о-1 < тспг, Бг Б-1 < Б-1. Из (16) следует

Б-1 = т* + П-1 & тяЬ, Б-1 = 1/(т- + т-1) + тте, & Б-1,

т1 = Бг Б-1 + Б-1 & Б-1 & Б-1.

Таким образом, при выполнении приведенных выше неравенств

с°а/П & с°аБд

Введем обозначения

Ф0 & с°а/п & с°аБд.

тк

А = —, Бг = вгт, Бк = вк тк, (41) тг

В 1 Л Г Я Лкс 'да в Ф в Тс°а Б Б

В1 = 1 — ЛгГs, в2 =-п/ . ^тт Ф1 = ——¡с Б = Б

т,0 Б а + Б, вк 1 + Б

0(А, Г,) = Б1Аа+ В2 (1 — Г,), ^ (А,^) = Уточним выражение (40):

т1 & Б-1 = т,0

Ф/ Бгтс°а

1 + Б (1 + кг)п

= тгФ1С(А, г,^).

Для приближения (40) модели фотосинтеза выпишем модель роста четырех биомасс mi(t) растительного покрова и сопряженную с ней модель динамики приведенного водного потенциала jjs корнеобитаемого слоя почвы:

d I

dfrni = Pi

d 1 + Ьв

dt js = awijs[Jqa{jjs) - Q

E Ф>=1, Pi

mi(0) = т0, js(0) = j0

5. Рациональное распределение ассимилятов

Рост биомасс Шг(Ь) в (43) зависит от значений ростовых функций р¿(¿), моделирующих сложные механизмы распределения и транспорта ассимилятов фотосинтеза в растении. Обычно функции рг(Ь) определяются натурными измерениями.

Варьируя рг(Ь) в (43), можно в значительной степени регулировать рост всех биомасс Шг(£) и, следовательно, величину биомассы ше{Ь) репродуктивного органа в конце жизненного цикла растения [0,Т]. Представляет интерес определить такие рг(£), при которых ше(¿) достигает максимума. Разумно предполагать, что такая структура распределения ассимилятов может возникать эволюционно у некоторых видов растений, например злаковых, подвергаемых отбору на максимальную урожайность репродуктивных органов.

Рассмотрим задачу выбора функций рг(£), максимизирующих для системы (43) - (45) на заданном отрезке [0, Т] функционал

ше

при ограничениях (45).

Обозначим через ф(С) = (ф1^),ф2^),фз^)) сопряженные с фазовыми переменными ш^) = (ш1^),ш2^),шз^)),г = 1, 2, 3, функции; пусть ф4(С) — сопряженная функция к Гамильтониан для системы (43) - (45) имеет вид

з

Н

г=1

где Н0 = рхфх + Р2Ф2 + рзфз, Г = (р = (р1,р2, рз) : р1 + р2 + рз = 1,рг > 0).

Согласно принципу максимума Понтрягина оптимальное управление

р*(^) = (р\(¿),р2&),рз&)) системой (43) - (45) необходимо подчиняется условию максимума

Н

р&

где Шг, р.8 удовлетворяют (43) - (45), а фг — соответствующей сопряженной системе:

ф = -дН/дшг, ф 4 = -дН/др,8, г =

фз

Для частного случая, когда параметры внешней среды постоянны, в частности постоянны приведенный водный потенциал р,8 почвы, температура листьев Т, заданы Q(t), ,7ф8ик-, а также Rl = Rk = RR, в [6] доказана теорема о магистрали для фазовых переменных из (43).

Опишем содержательную интерпретацию результатов формальной теоремы, представленных на рис. 2 и рис. 3, где т1 = тг, т2 = тк, т3 = ти. Заметим, что оптимальные р*(Ь) = 0, ибо из-за отсутствия в модели (43) - (45) учета архитектоники растения затраты ассимилятов на рост биомассы т, стеблей оказываются непродуктивными. Также нетрудно видеть, что затраты ассимилятов на рост ти продуктивны лишь на заключительном репродуктивном периоде [т2,Т]. Для начального вегетационного периода Т1 ^ Ь ^ Т2 рассмотрим функцию Е(А, удельной фотосинтетической мощности покрова, приходящейся на единицу суммарной биомассы:

Е(А, =

тг

тг + тк

С(А,Г,,Ф

1.

Рис. 2. Оптимальные ростовые функции модели

Рис. 3. Оптимальная траектория биомасс модели (43) - (46) Функция Е имеет единственный максимум при А* = л/В2/В1 — решение уравнения

А „ ^ _ ч —В1А2 + В2

= Ф.(1-Г )й_

дА 1

= Ф1

(А + 1)2 (В1А + В2)2

0,

равный ¥* = $1(1 — /л3)(л/Б1 + л/В2) 2. Из модели (43) для вегетационного периода т1 ^ t ^ т2

т = Р1Ш10Ф1 — Ешг, т к = рк т ОФ1 — Ешк,

при КI = Я'к = К для переменной Л(^ = тк(^)/тг(¿) получим уравнение

\ =

Уравнение интересно тем, что вместе с уравнением (44) для приведенного водного потенциала Д„ представляют замкнутую систему уравнений, описывающую динамику биомасс РП на отрезке Т1 ^ t ^ Т2. Для рассматриваемых в этом разделе постоянных условий среды постоянным является Л* = \/В2/В1, максимизирующее ¥. Полагая в (53) А = 0, получим систему

р*к — Л* р* = 0, рк + р* = 1,

из которой находим ростовые функции р**(и ркрешающие задачу (46), а также максимизирующие удельную фотосинтетическую мощность ¥: 1 Л*

р* = , рк = ГГА* , Р*к = 0,Р* = 0, t ^ [Т1,Т2)'

Сложив уравнения (52), с учетом данных в (52) значений Л*, р*, рк получим систему

выражений, определяющих переменные т0), тк({), М({), Д„Л*({) :

-

- 1п М = ¥ (Л*, Д,) — К,

тг =

в 1 + ьв

— Д8 = Дз [Зда

где

Л* В/В о (В'а + Р'з>„0

Выражения (56), заменившие на t е [т1, т2) систему (43), имеют замечательное свойство: они представляют биомассы интегралом, суммирующим по времени малые случайные вклады весьма переменчивых погодных условий, влияющих на рост через фотосинтетическую мощность ¥ (Л*, Д„,{).

Система (56) - (58) предполагает наличие в растениях некоего механизма адаптации, максимизирующего удельную фотосинтетическую мощность ¥ (Л, Д„ оптимальным распределением ассимилятов фотосинтеза между листьями и корнями растения, т.е. реализующим ростовые функции р**(Ь), рк

Итак, период [0, Т} жизненного цикла растения разбивается на три фазы: начальную [0, т1), вегетативную [т1, т2) и репродуктивную [т2,Т}. На начальной фазе система стремится достичь сбалансированного отношения биомасс корней и листьев Л*, направляя все ассими-ляты на рост органа с недостаточной биомассой. Например, начальное значение А(0) > Л* на рис. 3 характеризуется недостаточной биомассой листьев т[, поэтому все ассимиляты направляются на рост листьев, то есть

р* = 1,р*к = 0, р*К = 0, р* = 0, 0 < г < Т1.

В момент Т1 отношение Л(^ становится равным Л*, система достигает магистрали тк(^ = Л*тг (^ и переходит в состояние сбалансированного роста листьев и корней при t е [Т1,Т2). Это — вегетативная фаза, в течение которой выполняется (52). Скорость роста зависит от условий среды, в частности от приведенного водного потенциала Д„.

С момента Т2 начинается третья репродуктивная фаза, при которой все ассимиляты направляются на рост репродуктивной биомассы шд:

pi = 0,pk = 0, pR = 1, p*s =0, Т2 < t < T.

В момент времени T система заканчивает рост с максимальным значением шд (T).

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

Таким образом, оптимальное сбалансированное распределение ассимилятов соответствует максимизации удельной мощности фотосинтеза F(X, fls,t). Траектория магистрали обладает эффектом «притяжения» оптимальных траекторий с другими значениями начальных условий. Можно ожидать, что использование магистрали позволит лучше структурировать модель растительного покрова.

Литература

1. Полуэктов Р.А. Динамические модели агроэкосистем. Л.: Гидрометеоиздат, 1991.

2. Воротынцев А.В. Модель переноса тепла и влаги в системе почва-растение с учетом баланса энергии // Вестник БГУ. 2012. Вып. 1. С. 107-111.

3. Воротынцев А.В. Исследование модели водно-теплового режима системы почва-растение-приземный слой воздуха. М.: Вычислительный центр РАН. Сообщения по прикладной математике, 1991.

4. Воротынцев А.В. Исследование модели переноса тепла и влаги в системе почва-растение // Математическая биология и биоинформатика. Спецвыпуск по итогам конференции ЭкоМатМод-2011. 2012. Т. 7, № 1. С. 45-53.

5. Воротынцев А.В. Приближенная модель переноса тепла и влаги в системе почва-растение с учетом баланса энергии // Вестник БГУ. Математика, Информатика. 2014. Вып. 9, № 1. С. 22-26.

6. Воротынцев А.В., Рогожина Н.Л., Соколова И.А., Степченкова С.О. Две задачи оптимального управления для модели системы почва-растение. М.: Вычислительный центр АН СССР. Сообщения по прикладной математике, 1988.

7. Бихеле З.Н., Молдау Х.А., Росс Ю.К. Математическое моделирование транспирации и фотосинтеза растений при недостатке почвенной влаги. СПб.: Гидрометеоиздат, 1980.

References

1. Poluektov Е.А. Dynamic models of agroecosystems. L.: Gidrometeoizdat, 1991. (in Russian).

2. Vorotyntsev А.У. The model of heat and moisture transfer in soil-plant system with the energy balance. Vestnik BSU. 2012. N 1. P. 107-111. (in Russian).

3. Vorotyntsev А.У. Research of the model moisture-heat regime of soil-plant-air system. М.: Computing center of RAS. Reports on applied mathe-matics, 1991. (in Russian).

4. Vorotyntsev А.У. Research of the model of heat and water transfer in soil-plant system. Mathematical biology and bioinformatics. Special issue on the conference EcoMathMod-2011. 2012. V. 7, N 1. P. 45-53. (in Russian).

5. Vorotyntsev A.V. The approcximate model of moisture-heat transfer in soil-plant system subject to the energy balance. Vestnik BSU. Mathe-matics, Informatics. 2014. I. 9, N 1. P. 22-26. (in Russian).

6. Vorotyntsev A.V., Rogojina N.L., Sokolova I.A., Stepchenkova S.O. Two optimal control problem for a model of the soil-plant system. M.: Computing center of RAS. Reports on applied mathematics, 1988. (in Russian).

7. Bihele Z.N., Moldau C.A., Ross J.K. Mathematical modeling of transpiration and photosynthesis at the lack of soil moisture. SP.: Gidrometeoizdat, 1980. (in Russian).

Поступила в редакцию 19.07.2017