БИОФИЗИКА, 2020, том 65, № 6, с. 1219-1229
БИОФИЗИКА СЛОЖНЫХ СИСТЕМ
УДК 519.876.5
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ МИНЕРАЛЬНОГО
ПИТАНИЯ РАСТЕНИЙ В СИСТЕМЕ «УДОБРЕНИЕ-ПОЧВА-РАСТЕНИЕ»
© 2020 г. В.А. Четырбоцкий*, А.Н. Четырбоцкий**, Б.В. Левин***
*Центр инновации АО «Апатит», 119333, Москва, Ленинский просп., 55/1, стр. 1
**Дальневосточный геологический институт ДВО РАН, 690022, Владивосток, просп. 100-летия Владивостока, 159
***ПАО «ФосАгро», 119333, Москва, Ленинский просп., 55/1, стр. 1
E-mail: vchetyrbotskii@gmail.com
Поступила в редакцию 02.09.2019 г.
После доработки 29.06.2020 г.
Принята к публикации 19.08.2020 г.
Разработана численная модель пространственно-временной динамики многопараметрической си-
стемы, составляющими которой выступают биомасса растений, подвижные и неподвижные формы
элементов минерального питания среды, ризосферные микроорганизмы и параметры состояния
внешней среды (температура, влагосодержание, кислотность). Проведена параметрическая иден-
тификация и проверка адекватности разработанной модели на имеющейся выборке эксперимен-
тальных данных по росту яровой пшеницы «Красноуфимская-100» на торфяной низинной почве.
Результаты представлены временными распределениями биомассы исследуемой сельскохозяй-
ственной культуры и данными по содержанию основных элементов питания внутри самого расте-
ния (азота, фосфора, калия). Дана агрономическая оценка и интерпретация полученных резуль-
татов.
Ключевые слова: элементы минерального питания, минеральная емкость растений, кривая логистиче-
ского роста, параметрическая идентификация модели.
DOI: 10.31857/S0006302920060228
Применение большинства известных динами-
водятся характерные для замкнутых морских
ческих моделей агросистем ограничено изучени-
экосистем уравнения, составляющими которых
ем лишь небольшого числа составляющих эле-
выступают морская среда, биомасса фитопланк-
ментов системы ввиду, во-первых, большого чис-
тона, минеральные вещества и отмершая органи-
ла прямых и косвенных связей между самими
ка. Минеральное питание организмов подразде-
элементами и, во-вторых, отсутствием достовер-
ляется на группы сходных веществ (на основе
ного экспериментального материала для верифи-
азота, фосфора, кремния и т.д.), где морская сре-
кации выполненных работ. Так, в работе [1] при-
да играет роль ризосферы, фитопланктон - роль
водится модель динамики роста нескольких ви-
растений, минеральные вещества - элементов
дов сельскохозяйственных культур, где основным
питания. В работе [3], аналогично работам [1] и
предикатом продуктивности агросистемы высту-
[2], не учитываются факторы изменений состоя-
пает температура. При неблагоприятном темпе-
ния внешней среды, пространственного перерас-
ратурном режиме принимается отрицательный
пределения элементов ссреды и учет внутриси-
рост биомассы, что не отражает целесообразность
стемных взаимодействий. Модель динамики тра-
применения такой модели в реальных системах.
востоя в зависимости от климатических факторов
Рост численности популяции замкнутого сооб-
без учета внутренних особенностей отдельных
щества микроорганизмов в однородной стацио-
элементов приведена в работе [4]. Распределение
нарной среде изучался в работе [2]. Модель отра-
почвенного органического вещества без отсылки
жает превращение элементов питания среды
к продуктивности самой среды изучалось в рабо-
(жертвы) и микроорганизмов (хищника) друг в
тах [5, 6].
друга, тем самым имитируя круговорот превра-
Вне рамок перечисленных исследований ока-
щения питательных веществ, однако не учитыва-
зываются важные практические случаи простран-
ет внутривидовые взаимодействия популяции и
ственной неоднородности ризосферы, приуро-
влияния параметров внешней среды на динами-
ченность сообществ микроорганизмов к местам
ческие переменные выражений. В работе [3] при-
распределений питательных субстратов, внесе-
1219
1220
ЧЕТЫРБОЦКИЙ и др.
Рис. 1. Блок-схема модели химико-микробиологических процессов ризосферы.
ние удобрений в наиболее благоприятные фазы
ее плодородие, которое поступает извне посред-
роста растений. Вопрос о распределении раство-
ством внесения удобрений; ризосферные микро-
римых и нерастворимых форм питательных ве-
организмы, которые выступают катализатором
ществ, а также их связь между собой и растением
метаболических реакции в почве и которые слу-
был рассмотрен лишь в ограниченном числе слу-
жат универсальным механизмом следования пи-
чаев и не имеет дополнительных связывающих
тательных элементов до клеток растений; сами
факторов между системой «питательный эле-
растения.
мент-почва-растение». Также стоит отметить
Процесс взаимодействия приведенных объек-
отсутствие методов и алгоритмов решения задач
тов состоит в следующем.
параметрической идентификации моделей (спо-
собов оценивания численных значений парамет-
Удобрения попадают в почву в виде минераль-
ров модели) и проверки степени адекватности
ных солей. Определенная их часть связывается с
изучаемым реальным системам.
мезо- и микроэлементами почвы, образуя соли
труднорастворимых соединений. Эти элементы
Перечисленные особенности подчеркивают
неподвижных форм, как правило, статичны и
актуальность разработки модели динамики ком-
практически не меняют своей координаты с тече-
плексной системы, динамическими переменны-
нием времени, сорбируются в почвенные твердые
ми которой выступают растения, ризосфера и со-
породы и не растворяются в подвижный пита-
общества населяющих ее микроорганизмов, эле-
тельный почвенный раствор [8]. Остальная часть
менты подвижной и неподвижной форм как
либо улетучивается в атмосферу, либо сорбирует-
минерального (минеральные удобрения), так и
ся в жидкую фазу почвы (почвенный раствор),
органического (навоз, компост, торф, солома и
которая непосредственно прилегает к корневой
др.) питания. Всюду в дальнейшем под термином
зоне растений и доставляет питательные элемен-
«минеральное питание растений» будет пони-
ты к корневым клеткам. При этом соли одних
маться совокупность процессов поглощения, пе-
элементов находятся в определенных отношени-
редвижения и усвоения клетками корней расте-
ях с солями других элементов, катализируя либо
ний химических элементов, получаемых из поч-
ингибируя их реакции с почвенным раствором.
вы в форме ионов минеральных солей [7].
Перешедшая в почвенный раствор определен-
Блок-схема модели продукционных процессов
ная часть минеральных солей потребляется поч-
ризосферы представлена на рис. 1.
венными микроорганизмами. Дополнительным
Объектами исследования выступают: основ-
источником их питания служат корневые выделе-
ное питательное вещество почвы, определяющее ния растений. Конечным продуктом их реакции
БИОФИЗИКА том 65
№ 6
2020
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ
1221
являются биологически активные вещества и
ми распределениями элементов подвижной фор-
белковые соединения, поглощаемые вместе с пи-
мы в ризосфере, содержанием в ней органики и
тательными элементами среды растением. Такое
состоянием среды (распределениями температу-
поглощение служит основой для формирования
ры, влажности, плотности ризосферы, кислотно-
урожая и органических остатков ризосферы.
сти и т.п.), периодичностью внесения удобрений,
Цель работы состоит в построении математи-
взаимосвязями между элементами системы, про-
ческой модели динамики системы «удобрение-
странственными распределениями растений и
почва-растение», составляющими компонента-
ми которой выступают совокупность элементов
рядом других причин. Принимаются внутри- и
питания почвы, сообщества микроорганизмов
межвидовые отношения между содержанием эле-
ризосферы и биомасса растений. В качестве ос-
ментов (конкуренция, симбиоз, отношения «па-
новных типовых уравнений, используемых для
разит-хозяин» и т.д.), где каждый элемент «бо-
построения выражений модели, были использо-
рется» за вакантное место в системе «почва-рас-
ваны уравнения самолимитированного роста
тение» и тем самым вытесняет другие элементы
Ферхюльста и ограничения потребления субстра-
из почвенного раствора. Примером такого пове-
та Моно [3, 7]. Для каждой из переменных выпол-
нено построение модели, учитывающей опреде-
дения может служить отношение между цинком и
ленные связи и самосогласованный характер
кадмием, а также между кадмием и кальцием:
составляющих элементов системы. Далее выпол-
цинк и кальций вытесняют кадмий из раствора
нено исследование роста яровой пшеницы на
его солей [9]. Принимаются также определенные
торфяной низинной почве с помощью программ-
взаимоотношения между содержаниями элемен-
ной реализации, адаптированной под уникаль-
тов в самом растении. Допускаются простран-
ный репрезентативный экспериментальный ма-
ственно-временные перераспределения элемен-
териал.
тов неподвижной формы (их миграция).
МОДЕЛЬ ДИНАМИКИ ЭЛЕМЕНТОВ
Согласно принятым допущениям запись урав-
МИНЕРАЛЬНОГО ПИТАНИЯ
нений модели динамики элементов питания по-
Динамика содержания в растениях элементов
движных и неподвижных форм принимает следу-
минерального питания определяется начальны- ющий в
nC
C
i
=
f
С
(T,W
)[q
i,0
R
i
+
q
i, j
C
j
]B
t
j=1
nK
(aw)
(
a)
(s)
)
(wa w)
R
=
1
a
(C
C
)
+
a
p
δ
t
t
+
a
S +a
C
i
(
i
)
i,0
i
i
ik
(
i
,k
)
i
i
i
,
(1)
k=
1
(w)
n
K
C
i
(aw)
(w)
)
(wa w)
(w)
2
(w)
=
a
(
C
C
)
+
a
p
δ(t
t
)
a
C
+
D
C
i
i,0
i
i
ik
i
,k
i
i
i
i
t
k=1
(w
C
(0,X )
=
0 и
C
)(
0,
X
)
=
0, где
X
(x,y,z)
i
i
где {Ci, i = 1 ÷ nC} - выраженное в единицах мас-
циенты {qij, i,j = 1 ÷ nC} - интенсивность влияния
сы содержание в растении элементов его мине-
элементов питания на их прирост в самом расте-
рального питания в подвижной форме, nC - чис-
нии. Наличие такой связи обусловлено химиче-
скими свойствами элементов питания: в одних
ло видов элементов питания; Ci(w) - выраженное
случаях они могут выступать ингибиторами дина-
в единицах массы содержание этих же элементов
мики, а в других - активаторами.
неподвижной формы в ризосфере. Функция
Набор {Ci,0, i = 1 ÷ nC} определяет начальное
fC(T,W) характеризует влияние температурного
распределение элементов питания подвижной
(T) и влажностного (W) режимов. Понятно, что ее
формы в ризосфере (Ri вводится для удобства за-
набор аргументов не ограничивается только эти-
ми переменными: сюда следует также добавить
писей). Коэффициенты ai(a) и ai(w) характеризуют
показатель pH, плотность среды ρ и т. д. Для ком-
те доли внесенных удобрений, которые приходят-
пактности записей ограничимся только темпера-
ся на элементы подвижной и неподвижной форм.
турой и влажностью. Далее, коэффициенты про-
Их сумма не превышает единицу (и может ока-
порциональности {qi,0, i = 1 ÷ nC} характеризуют
заться меньше ее), поскольку удобрения могут
интенсивности поступления элементов подвиж-
частично улетучиваться либо не попадать в рас-
ной формы из ризосферы в растение, а коэффи-
твор. Набор {pik, i = 1 ÷ nC, k = 1 ÷ nK} характери-
БИОФИЗИКА том 65
№ 6
2020
1222
ЧЕТЫРБОЦКИЙ и др.
зует внесенные в ризосферу в момент времени tik
естественным приростом и убылью, элементами
удобрения, где δ(t - tik) - дельта-функция Дира-
подвижной формы (азот, фосфор и т.д.), внутри-
и межвидовой конкуренцией за пищевой суб-
ка; ai(aw) и ai(wa) характеризуют взаимные перехо-
страт (он продуцируется самим растением, стеб-
ды элементов между подвижной и неподвижной
ли, листья и корни которого выделяют углеводы и
формами; ai(s) - доля опада, которая приходится
органические кислоты) и пространственный ре-
на i-й элемент подвижной формы; nK - порядко-
сурс (они расселяются на поверхности и внутри
растений). Эти выделения (корневые экссудаты)
вое число вносимых в почву видов удобрений.
Коэффициенты (1) подлежат оцениванию на ос-
состоят преимущественно из различных органи-
новании вычислительных экспериментов.
ческих соединений. Они активно выделяются
корнями растений в окружающую среду и обеспе-
Согласно уравнению (1), динамика элементов
чивают питательным субстратом микроорганиз-
минерального питания подвижной формы в рас-
тении следует положениям, принятым в модели
мы ризосферы, создавая тем самым благоприят-
«ресурс-потребитель» [10]. При расчете элемен-
ные условия для их существования в прикорне-
тов подвижной формы в ризосфере следует учи-
вых зонах [11].
тывать, что их количество Ci перешло в само рас-
В выделениях корней содержатся органиче-
тение, а количество ai(aw)(Ci,0 - Ci) перешло в не-
ские соединения большой физиологической ак-
подвижную форму. Добавлены же сюда элементы
тивности - витамины, ростовые вещества, ино-
удобрений, опада и перешедшие в подвижную
гда алкалоиды и т. д. Многие из них выделяются и
форму элементы неподвижной формы.
надземными органами растений. Поэтому на
Динамика элементов неподвижной формы
корнях и надземных органах растений размножа-
Сi(w) определяется суммой перешедших из по-
ется обильная сапрофитная микрофлора [9]. Пи-
тание микроорганизмов в ризосфере обеспечива-
движной формы элементов (первый член в пра-
ет коровые ризодепозиты, включая высокомоле-
вой части) и той части массы удобрений, которая
кулярные слизи полисахаридной и белковой
приходится на элементы неподвижной формы.
Допускается их пространственно-временное пе-
природы, утраченные части растений (корневой
рераспределение (миграция), которую тут выпол-
чехлик, корневые волоски, отмершие части кор-
ня) [12, 13]. Ряд видов микроорганизмов находят-
няет диффузионный механизм (Di(w) - коэффи-
ся в антагонистических отношениях с другими
циент диффузии и 2 - оператор Лапласа).
видами, что ведет к регрессу системы. Например,
среди эпифитных бактерий встречаются такие,
которые образуют действующие отрицательно на
МОДЕЛЬ ДИНАМИКИ
фитопатогенные бактерии вещества и тем самым
МИКРООРГАНИЗМОВ
предохраняют растения от заболеваний [11].
Микроорганизмы выступают агентами пере-
носа питательных для растений элементов и су-
Согласно принятым допущениям, запись
щественно стимулируют их рост. Динамика мик-
уравнений модель динамики сообществ микро-
робиоты определяется состоянием ризосферы,
организмов принимает вид:
n
Y
∂Y
μ(S,B)
k
=
f
Y
(T,W
)
b
k,j
Y
j
Y
k
+∇(A
k
∇μ(S,B))
(μ)
t
b
+
b
μ(S,B)
k
k
j=1
(s)
(b)
μ(S,B)
=
b
k
S +b
k
B
,
(2)
(
Y x,y,
z)
=
Y
(
x y
z)
k
k,0
Y
k
|
=
0
n Γ
где fY(T,W) характеризует влияние температурно-
низмов их питательный субстрат; bk,j - коэффи-
го T и водного W режимов ризосферы на динами-
циенты, посредством которых учитывается внут-
ку микроорганизмов (надо также учитывать кис-
ри- и межвидовая конкуренции микроорганиз-
лотность среды pH, ее плотность ρ и т. д. ); Yk -
мов. Учет разномасштабного вклада опада и
растений в питательный субстрат микроорганиз-
численность/биомасса микроорганизмов k-го
вида, k = 1 ÷ nY; S и B - масса органики и масса
мов фиксируется с помощью коэффициентов
растений, которые представляют для микроорга-
bk(s) и bk(s). Ak - функция эффективности поиска
БИОФИЗИКА том 65
№ 6
2020
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ
1223
питательного субстрата [14]; слагаемое вне квад-
мирование органики. При построении моделей
ратных скобок определяет механизм простран-
динамики биомассы и органики были приняты
ственного перераспределения микроорганизмов,
следующие допущения:
согласно которому они следуют дивергенции гра-
1. Динамика биомассы растений определяется
диента массы их питательного субстрата (хемотак-
ее естественной скоростью роста, содержанием в
сис [10,14]). Последнее соотношение в выражениях
ризосфере элементов минерального питания,
(2) указывает отсутствие потоков на границах обла-
внутри- и межвидовых отношений между расте-
ниями (конкуренция, симбиоз и т. д.).
сти рассмотрения модели; ∇ =
,
,
)
- опера-
x y z
2. Характер потребления растением элементов
тор градиента, а ∇ - дивергенция выражения. Ко-
минерального питания определяется положени-
эффициенты (2) подлежат оцениванию на осно-
ем модели «ресурс-потребитель», где ресурсом
вании вычислительных экспериментов.Согласно
выступает элементы минерального питания, а по-
первому уравнению
(2), локальная динамика
требителем - растение.
микроорганизмов (выражение в квадратных
3. Отдельные виды микроорганизмов вызыва-
скобках) следует принятым для микроорганизмов
ют стимулирующий рост растений; в частности,
положениям естественного прироста в форме
продуцируемые ризобактериями антибиотики и
Моно [14]. Питательным субстратом для них вы-
токсины воздействуют на фитопатогены, приво-
ступает органика и растения. Микроорганизмы
дя к ингибированию их роста и гибели [17].
следуют распределению своего питательного суб-
страта, что определяется здесь дивергенцией от
4. Учитываются отношения между элементами
градиента массы субстрата. Такой механизм про-
питания (конкуренция, симбиоз и т.д.), которые
странственного распределения микроорганизмов
уже употребились растением.
подобен механизму режима с обострением [15,
5. Попавшие в ризосферу отмершие части рас-
16], когда происходит взрывной характер роста
тений трансформируются в элементы органики,
биомассы микроорганизмов за конечное время.
где интенсивность их накопления определяется
И чем выше масса субстрата (в том числе и био-
функцией возраста растения (чем оно старше, тем
масса растений), тем выше масса микроорганиз-
выше интенсивность отмирания его частей) и те-
мов.
кущего временного сезона года (весна, лето или
осень).
МОДЕЛЬ ДИНАМИКИ БИОМАССЫ
6. Динамика органики определяется отмерши-
РАСТЕНИЙ И ОРГАНИКИ
ми частями растений и ее потреблением микро-
организмами, для которых она выступает их суб-
Динамика биомассы растений определяется
стратом.
элементами его питания, конкуренцией за ре-
сурс, результатами жизнедеятельности микроор-
Следуя этим допущениям, запись уравнений
ганизмов и интенсивностью отмирания расте-
модели динамики биомассы растений и органики
ния, что существенным образом влияет на фор-
принимает следующий вид:
˙
n
C
n
Y
∂
(с)
(y)
B = f
(T,W
)[d
φ
t)
+
d
C
+
d
Y -d
B]B
B
0
i
i
k k
1
t
i=1
k=1
(
B x,y,z)
=
B
(x y
z
)
min
,
(3)
n
Y
S
(s)
(s)
2
t)B
(d
+
d
Y
)S +D
S
0
k k
S
t
k=1
(
S x,y,
z
)
=
S
(x y
z)
0
где fB(T,W) характеризует влияние температурно-
которая существенным образом зависит от воз-
го T и влажностного W режимов ризосферы на
раста растения и сезона года (например, с наступ-
динамику биомассы (как и в предыдущих случа-
(c)
лением осени деревья теряют свою листву); di
ях, необходимо учитывать кислотность среды рН,
характеризует отнесенную к единице биомассы
ее плотность ρ и т. д.); d0 - коэффициент есте-
интенсивность потребления элементов подвиж-
ственного прироста биомассы; ϕ(t) характеризует
ной формы, i = 1 ÷ nC; посредством dk(y) учитыва-
интенсивность формирования опада растения
(стволы, сухие вершины, ветви, сухие стебли), ется влияние микроорганизмов на динамику ро-
БИОФИЗИКА том 65
№ 6
2020
1224
ЧЕТЫРБОЦКИЙ и др.
ста растения; последний член в квадратных скоб-
ТЕМПЕРАТУРНЫЕ И ВЛАЖНЫЕ РЕЖИМЫ
ках характеризует самолимитирование растения
(влияние самого растения на свой рост).
При построении модели полагаем, что
Последние два уравнения определяют дина-
воздействие среды ризосферы на динамические
мику органики, где d0(s) - коэффициент есте-
переменные системы задается функциями
fC(T,W,pH), fY(T,W,pH,ρ), fB(T,W,pH,ρ) (как
ственного разрушения органики или ее есте-
ственных потерь (например, выветривание орга-
правило, кусочно-непрерывными функциями [1,
ники), а dk(s) - коэффициенты потребления
18, 19]). Представляется естественным, что их
органики микроорганизмами.
запись принимает вид:
(C)
(C)
(C)
(C)
(C)
(C)
(C)
(C)
f
T,W,
pH
=
C
ψ
T
,T
,T
ψ
W
,W
,W
ψ
pH
,
pH
,
pH
ψ
ρ
C
(
)
C
(
1
2
)
(
1
2
)
(
1
2
)
(
1
2
)
(Y )
(Y )
(Y )
(Y )
(Y )
(Y
)
(Y )
(Y )
f
Y
(T,W,
pH)
=
C
Y
ψ
(
T
1
,T
2
,T
)
ψ
(
W
1
,W
2
,W
)
ψ
(
pH
1
,
pH
2
,
pH
)
ψ
(
ρ
1
2
)
,
(4)
(B)
(B)
(B)
(B)
(B)
(B)
(B)
(B)
f
(T,W,
pH)
=
C
ψ
T
,T
,T
ψ
W
,W
,W
ψ
pH
,
pH
,
pH
ψ
ρ
(
)
(
)
(
)
(
)
1
2
1
2
1
2
1
2
B
B
где CC, CY, CB - эмпирические константы, а
(X - X
min
)(X
max
X
)
при
X
min
<
X < X
max
ψ(X
,X
, X )
min max
=
0
При записи (4) полагается, что первые два ар-
где K(W) - коэффициент фильтрации в ненасы-
гумента функций в правых частях характеризуют
щенной зоне, η - высота всасывающего давле-
оптимальный для жизненного цикла объектов
ния, высота над плоскостью сравнения..
диапазон соответствующих величин. Выбор
именно такого представления обусловлен тем,
ЧИСЛЕННАЯ РЕАЛИЗАЦИЯ
что за пределом оптимального для диапазона тем-
пературы, влажности и т. д. имеет место его адек-
Для изучения динамики роста яровой пшени-
цы «Красноуфимская-100» на торфяной низин-
ватная реакция. Кроме того, такое представление
ной почве была выполнена адаптация модели, где
удобно для программной реализации. При небла-
учитывается особенность доступного для иссле-
гоприятном для объектов состоянии ризосферы
дований экспериментального материала. Здесь
правая часть первого уравнения (3) становится
он представлен временными распределениями
отрицательной, что приводит к падению биомас-
биомассы растения и выраженных в единицах
сы. Когда потери биомассы ниже допустимого
массы содержания в растении его основных эле-
предела, растение погибает.
ментов питания (калий, фосфор и азот).
Для нахождения динамики температуры в ри-
Схема опытов и условия закладки экспери-
зосфере допустимым представляется использова-
ментов состоят в следующем [22]. Перед посе-
ние уравнения теплопроводности в виде
вом специально отобранную торфяную низин-
ную почву обработали азотными, фосфорными
T
и калийными удобрениями (сульфат аммония,
ρc
p
=∇D
T
T
,
(5)
простой суперфосфат и сульфат калия соответ-
t
ственно). Эксперименты проводили на десяти
где cp - удельная теплоемкость вещества ризосфе-
типах объектов (варианты торфяной низинной
почвы), отличающихся различием содержания
ры, DT - коэффициент его теплопроводности.
K2O. Сам объект представлял собой сосуд с поч-
Поскольку моделью ризосферы выступает ка-
вой, масса которой составляла 3.5 кг. Числен-
пиллярно-пористое тело, то здесь выполняется
ные значения измеряли в масштабных единицах
закон Дарси [20]: поток воды направлен в сторону
«г/сосуд». Набор {24,28,34,39,46,55} определяет
уменьшения давления почвенной влаги и про-
даты отбора растительных проб с момента появ-
порционален ее градиенту. В ризосфере динами-
ления первых всходов, которые соответствуют
ка распределения влаги по вертикальному про-
наиболее критическим фазам роста растений
филю, как правило, задается уравнением Ричард-
(единица времени - сутки). Поэтому начальный
са [21]:
момент приходится на 24-е сутки проведения
экспериментов. Характер жизненного цикла та-
W

ких растений предопределяет его переход от фа-
=
K
(W
)∂η
+
1
,
(6)
зы всхода к стадии кущения и прорастанию пер-
t
z
∂z
БИОФИЗИКА том 65
№ 6
2020
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ
1225
Рис. 2. Экспериментальные распределения динамических переменных модели.
вичных (зародышевых) корней на глубину 50-
Характер распределений экспериментальных
55 см. Вегетационный период яровой пшеницы
поверхностей на рис. 2 указывает на их относи-
в зависимости от сорта, районов возделывания и
тельно плавные и согласованные между собой из-
погодных условий колеблется от 85 до 110 суток
менения, причем отмеченная тенденция является
[23, 24]. В экспериментах обеспечивали посто-
общей для всех представленных случаев. Доста-
янный температурный и водный режимы ризо-
точно отметить, что если мерой согласованности
сферы, где в граммах на сосуд учитывали только
распределений выступает коэффициент корреля-
ции между ними, то он превышает значения
биомассу надпочвенной части растения и коли-
0.918 (корреляция между наблюдениями на 39-е и
чественное содержание в ней элементов мине-
55-е сутки). Для оценки корреляций была сфор-
рального питания. Экспериментальные распре-
мирована матрица, столбцами которой выступа-
деления динамических переменных модели
ют значения биомассы и концентрации мине-
представлены на рис. 2.
ральных элементов за определенные сутки отбора
На рис. 2а представлено распределение био-
проб.
массы растений в момент отбора проб (граммы),
Форма поверхностей несущественно отлича-
а в остальных случаях - пересчитанные на окис-
ется от поверхностей, которые обычно определя-
лы (кроме азота) содержания элементов в расте-
ются принятыми в физиологии растений модели
ниях (в мг). Ось «K2O (удобрение), г» показывает
Ферхюльста или так называемыми логистически-
содержание этого элемента в ризосфере при вы-
ми уравнениями [10, 14]. Можно заметить, что
полнении экспериментов по влиянию удобрений
здесь имеет место более длительный индукцион-
на рост биомассы. Независимо от внесения ка-
ный период (от момента начала роста до момента
лийных удобрений почва также характеризуется
его быстрого роста или точки перегиба) по срав-
наличием азота и фосфора. Всего, как уже сказа-
нению с логистической кривой.
но выше, для выполнения экспериментов было
Отсутствие требуемых для оценки модели си-
подготовлено десять вариантов торфяной низин-
стемы (1)-(3) соответствующих пространствен-
ной почвы с различным содержанием K2O.
но-временных выборочных распределений обу-
БИОФИЗИКА том 65
№ 6
2020
1226
ЧЕТЫРБОЦКИЙ и др.
словили рассмотрение здесь динамики биомассы
тения. Запись редуцированных уравнений
и содержание основных элементов ее питания
систем (1)-(3) для такого случая принимает сле-
(K2O, P2O5 и N) только в надземной части рас-
дующий вид:
3
С
i
=
q
C
q
C
+
q
C
B
i,0
i,0
i,i
i
i, j
j
t
ji
(7)
B
,
,
i
=1÷3
=
(d
0
+
d
K
C
1
+
d
P
C
2
+
d
N
C
3
d
B
B)B
∂t
C
(0,
K)
=
C
(
K
)
и
B(0,
K)
B
=
B
(K)
i
i,00
0
где B0 и Ci,00 - измеренные в экспериментах за-
выполнении экспериментов они были постоян-
данные начальные распределения переменных;
ными.
аргумент K соответствует распределению участ-
Для оценки применимости модели и провер-
вующих в экспериментах типов почв с различным
ки ее адекватности выборочным распределениям
содержанием K2O. Для компактности записей в
следует выполнить процедуру параметрической
уравнениях (7) по сравнению с уравнениями (1)
идентификации, т.е. найти значения параметров
принято qi,0 = 1 - ai(aw). Коэффициенты модели
модели и установить степень соответствия моде-
оцениваются на основании серии вычислитель-
ли экспериментальному материалу. Здесь значе-
ных экспериментов. Маркировка нижних индек-
ния параметров являются решением задачи поис-
сов у коэффициентов для динамики биомассы со-
ка экстремума функционала, при построении ко-
отнесена с наименованием соответствующего хи-
торого надо учитывать различие масштабов
мического элемента (калия, фосфора и азота).
переменных (7) (величины [B] выражены в г, а
Здесь отсутствует разделение на элементы по-
[C] - в мг). Поэтому переменные этого функцио-
движной и неподвижной форм, поскольку изме-
нала должны быть безразмерными.
рения последних не были выполнены. В уравне-
ниях (7) также не учитывается влияние темпера-
Здесь задача поиска экстремума формализует-
турного и влажностного режимов, т.к. при ся следующим образом:
45
10
3
2
2
(
d)
d)
Φ(p)
=
1−
B
/
B
+
1−
(
C
)
/
(C
)(
min,
(8)

(
t,k
t,
k
)
(
i t,k
i t,k
)
t
=2
k
=1
i
=1
(d)
(d)
(d)
лось процедурой lsqnonlin программной среды
где массивы
B
=
B
и
C
=
{(
C
)
}
-
{
t,k
}
i
i t,
k
MATLAB [26]. Начальные значения параметров
выборочные распределения, а B = {Bt,k} и Ci =
оценивали согласно такой схеме.
{(Ci)t,k} следуют решению системы
(7). При
записи (8) учитывается тот факт, что всего было
Поскольку уравнения (7) линейны по своим
45 временных точек и десять типов почв с
параметрам, то их начальные значения можно
различным содержанием K2O.
оценить стандартным методов наименьших
квадратов. Для каждого i набор {qij, j = 1 ÷ 3}
Решение задачи (8) следует методологии нели-
нейного оценивания параметров [25] и выполня- определяется вектором зависимых переменных:
(
d
(d
C
t
+
(Δ
t
)
,k
C
t
,
k
i
)(
j
j
)
i
)(
j
)
,
j
=
n
t
,
k
=
10,
(d
B
)(
t
j
,
k
)(
Δ
t
)
j
где верхний индекс d указывает на выборочные
ременными выступают наборы {Ci,0(tj,k), j = 1÷ nt,
данные, tj - моменты времени проведения съе-
k = 1 ÷ 10} и {Ci(tj,k), j = 1 ÷ nt, k = 1 ÷ 10} при ко-
мок, D(t)j - промежуток между съемками в сут-
эффициентах правой части (7). При оценке коэф-
ках, nt - число съемок, k - число типов почв в
фициентов модели биомассы вектор зависимых
экспериментах. Для каждого i независимыми пе- переменных определяется выражением
БИОФИЗИКА том 65
№ 6
2020
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ
1227
(d
(d
B
t
+
(Δt)
,k
B
)(
t
,k
)
)(
j
j
)
j
,
j
=
n
t
,
k
=
1÷10.
(d
B
)(
t
,k
)(
Δ
t)
j
j
Численные оценки параметров модели и их
ментальным и модельным распределениями для
доверительные интервалы приведены в таблице,
биомассы равен 0.887.
где нижние индексы заглавий строчек и столбцов
Следуя модели (7) и используемому здесь экс-
соотнесены с наименованием соответствующих
периментальному материалу, необходимо отме-
химических элементов (K, P и N).
тить менее значимое прямое влияние калия и
фосфора на динамику биомассы растения.
Анализ элементов таблицы показывает полу-
Адекватность моделей минерального питания
торное превышение показателя интенсивности
поступления N из ризосферы в растение q0 по
их натурным экспериментальным распределени-
ям следует высоким значениям коэффициентов
сравнению с такими показателями для K2O и
корреляции между ними.
P2O5. Необходимо также отметить существенное
Период вегетации яровой пшеницы «Красноу-
превышение (более чем на два порядка) dN зна-
фимская-100» составляет в среднем 95 суток [27].
чений dK и dP, каждое из которых характеризует
На рис. 3 приведены модельные распределения
интенсивность/темп влияния N на динамику
биомассы (а) и элементов минерального питания
биомассы. Такая ситуация полностью согласует-
(б-г).
ся с натурными исследованиями его важности в
минеральном питании растений [7]. Действи-
Анализ распределений показывает, что во всех
тельно, продолжительность его поступления в
случаях они следуют близким к логистическому
растения совпадает с началом вегетации и про-
росту кривым. При этом рост значений динами-
должается до фазы молочной спелости [27]. Оп-
ческих переменных коррелирует с ростом содер-
тимальное азотное питание обусловливает усиле-
жания в ризосфере калийных удобрений.
ние синтеза пластических веществ и ускоряет
Для динамики биомассы растения степень
рост растений. При этом, в отличие от других эле-
важности его элементов питания можно оценить
ментов, N характеризуется высокой мобильно-
исходя из таких положений. Согласно уравнени-
стью и большим разнообразием форм, что и обу-
ям (7), конечные содержания элементов в над-
словливает важность его высокого положения
земной части растений (стационарное состояние)
среди остальных элементов минерального пита-
и соответствующие значения биомассы опреде-
ния. Коэффициент корреляции между экспери-
ляются соотношениями
3
q
i,
j
*
C
j
=
C
i,0
q
j=1
i
i
=
3,
*
d
d
*
d
*
d
*
0
K
P
N
B
=
+
C
+
C
+
C
K
P
N
d
d
d
d
B
B
B
B
Оценка коэффициентов параметров модели
q0
qK
qP
qN
R
CK
0.004 ± 0.001
0.029 ± 0.002
0.001 ± 0.003
0.001 ± 0.003
0.889
CP
0.004 ± 0.001
0.041 ± 0.005
0.001 ± 0.003
0.007 ± 0.004
0.901
CN
0.006 ± 0.001
0.025 ± 0.004
0.011 ± 0.001
0.001 ± 0.005
0.907
d0
dK
dP
dN
dB
B
0.128 ± 0.017
0.037 ± 0.137
0.0001 ± 0.662
5.457 ± 1.320
0.129 ± 0.032
Примечание. Число после знака ± указывает размах для доверительного интервала; в графе R для каждого вида минерального
питания указаны коэффициенты корреляции между экспериментальными и модельными распределениями.
БИОФИЗИКА том 65
№ 6
2020
1228
ЧЕТЫРБОЦКИЙ и др.
Рис. 3. Распределения динамических переменных модели (г/сосуд) на весь период вегетации растения (ось «K2
характеризует его содержание в ризосфере до посева).
где для удобства записи и последующей интер-
Результат применения адаптированной под
претации приводится смешанное использование
динамику роста яровой пшеницы модели пока-
нижних индексов. Из второго соотношения сле-
зал, что на почвах с повышенным содержанием
K2O имеет место повышенное значение биомас-
d
дует, что естественный прирост биомассы 0 со-
сы. При этом азот оказывает самое большое вли-
dB
яние на темп роста биомассы, следующим оказы-
ставляет порядка единицы, а диапазон воздей-
вается калий, а фосфор практически не оказывает
ствия азота для почв с различным содержанием
влияния. Превалирующее здесь влияние азота
удобрений - от 8.5 до 12.7.
обусловлено его ролью при росте и урожайности
растительных культур, поскольку он является
строительным материалом для их новых клеток
ВЫВОДЫ
[22-24]. Результаты моделирования как раз и от-
ражают этот факт, что подтверждает адекватность
Разработана математическая модель динами-
разработанной модели.
ки системы «удобрение-почва-растение», где
динамическими переменными выступают: содер-
жание в растении элементов ее минерального пи-
КОНФЛИКТ ИНТЕРЕСОВ
тания в подвижной форме, ризосферные микро-
Авторы заявляют об отсутствии конфликта
организмы и биомасса растений. В ризосфере
интересов.
учитываются также взаимные переходы между
элементами подвижной и неподвижной форм. На
СОБЛЮДЕНИЕ ЭТИЧЕСКИХ СТАНДАРТОВ
основании положений модели «ресурс-потреби-
тель» разработана модель динамики микроорга-
Настоящая работа не содержит описания ис-
низмов, где ресурсом выступает биомасса и от-
следований с использованием людей и животных
мершие части растений.
в качестве объектов.
БИОФИЗИКА том 65
№ 6
2020
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ
1229
СПИСОК ЛИТЕРАТУРЫ
14. Г. Ю. Ризниченко и А. В. Рубин, Математические
модели биологических продукционных процессов
1. Е. П. Колпак и М. В. Столбовая, Журн. научных
(Изд-во МГУ, М., 1993).
публикаций аспирантов и докторантов
(2013).
15. А. А. Самарский, В. А. Галактионов, С. П. Курдю-
URL: http://jurnal.org/articles/2013/mat7.html
мов и А. П. Михайлов, Режимы c обострением в за-
2. Н. С. Стригуль, Автореф. дис
канд. биол. наук
дачах для квазилинейных параболических уравнений
(Санкт-Петербург, 2000).
(Наука, М., 1987).
3. А. И. Абакумов, в кн. Математическое моделиро-
16. J. Murray, Mathematical Biology: I. An Introduction
вание в экологии (ИФХиБПП РАН, Пущино, 2011),
(Springer., New-York, 2002).
сс. 18-19.
17. A. Beneduzi, A. Ambrosini, and L. M. Passaglia, Gen-
4. И. М. Михайленко, С.-х. биология, № 1, 103
et. Mol. Biol. 35 (Suppl. 4), 1044 (2012).
(2007).
18. T. Roose, S. D. Keyes, K. R. Daly, et al., Plant Soil,
5. А. С. Комаров, Ю. С. Хораськина, С. С. Быховец и
No. 407, 9 (2016).
др., Мат. биология и биоинформатика 7 (1), 162
19. Л. А. Хворова, А. Г. Топаж, А. В. Абрамова,
(2012.
К. Г. Неупакоева, Изв. Алтайского государствен-
6. M.-C. Braakhekke, T. Wutzler, C. Beer et al., Biogeo-
ного ун-та, № 1-1 (85), 192 (2015).
sciences, No 10, 399 (2013).
20. Н. С. Бахвалов и Г. П. Панасенко, Осреднение про-
7. В. Г. Минеев, Агрохимия (Изд-во МГУ, «КолосС»,
цессов в периодических средах (Наука М., 1984).
М., 2004).
21. M. A. Celia, E. T. Bouloutas, and R. L. Zabra, Water
8. Д. С. Орлов, Химия почв (Изд-во МГУ, М., 1992).
Resource Res. 26 (7), 1483 (1990).
9. А. Кабата-Пендиас и Х. Пендиас, Микроэлементы
22. Н. А. Сладкова, Дис
канд. биол. наук (Санкт-
в почвах и растениях (Мир, М., 1989).
Петербург, 2016).
10. Ю. М. Свирежев, Нелинейные волны, диссипатив-
23. Растениеводство, под ред. Л. А. Посыпанова («Ко-
ные структуры и катастрофы в экологии (Наука,
лосС», М. 1997).
М., 1987).
24. Растениеводство с основами селекции и семеновод-
ства, под ред. Г. В. Коренева (Агропромиздат, М.,
11. Н. В. Феоктистова, А. М. Марданова, Г. М. Хадие-
1990).
ва и М. Р. Шарипова, Уч. записки Казанского ун-
та. Сер. Естественные науки 158, 207 (2016).
25. Й. Бард, Нелинейное оценивание параметров (Ста-
тистика, М., 1979).
12. О. Е. Беззубенкова, М. Н. Юхимова и Н. И. Несте-
рова, Естественные и технические науки 4, 99
26. В. П. Дьяконов, MATLAB. Полный самоучитель
(2012).
(ДМК Пресс, М., 2012).
13. В. П. Иванов, Корневые выделения и их значение в
27. Б. А. Доспехов, Методика полевого опыта (Агро-
жизни фитоценоза (Наука, М., 1973).
промиздат, М., 1985).
Mathematical Modeling of the Dynamics of Plant Mineral Nutrition
in the Fertilizer-Soil-Plant System
V.A. Chetyrbotskiy*, A.N. Chetyrbotskiy**, and B.V. Levin***
*Joint-Stock Company «Apatit», Leninskii prosp. 55/1, build. 1, Moscow, 119333 Russia
**Far East Geological Institute, Far Eastern Branch of the Russian Academy of Sciences,
prosp. 100-letiya Vladivostoka 159, Vladivostok, 690022 Russia
***Public Joint-Stock Company «PhosAgro», Leninskii prosp. 55/1, build. 1, Moscow, 119333 Russia
A numerical simulation of the spatial-temporal dynamics of a multi-parameter system is developed. The
components of this system are plant biomass, mobile and stationary forms of mineral nutrition elements, rhi-
zosphere microorganisms and environmental parameters (temperature, humidity, acidity). Parametric iden-
tification and verification of the adequacy of the model were carried out based on the experimental data on
the growth of spring wheat «Krasnoufimskaya-100» on peat lowland soil. The results are represented by tem-
poral distributions of biomass from agricultural crop under study and the findings on the content of main nu-
trition elements within the plant (nitrogen, phosphorus, potassium). An agronomic assessment and interpre-
tation of the obtained results are given.
Keywords: mineral nutrition elements, plant mineral capacity, logistic growth curve, model parametric identifica-
tion of the model
БИОФИЗИКА том 65
№ 6
2020