Математическое и компьютерное моделирование продуктивности растений в зависимости от динамики влажности почвы

Математическое и компьютерное моделирование продуктивности растений в зависимости от динамики влажности почвы

В. М. Казиев, С. К. Кирьязева, Д. А. Кирьязев

При разработке различных систем автоматизированного прогнозирования урожайности, при расчете максимальных урожаев и их агротехническом, экономическом, экологическом обеспечении важное место занимают модели роста и развития растений. Растение - сложная стохастическая система, содержащая множество параметров состояния, количественные изменения которых ведут к количественному и качественному изменениям всей системы в целом. Математическая модель роста и развития растений должна описывать основные процессы, на которые влияет управляющее воздействие. В первом приближении (достаточном для моделирования ростовых функций) система “растение - среда обитания” может быть интерпретирована как динамическая система с распределенными параметрами, а математические модели системы могут быть описаны с помощью дифференциальных уравнений. При построении таких моделей необходимо принимать во внимание те значительные трудности, которые возникают при идентификации моделей, а также невозможность точно и полно описать такую сложную динамическую систему как “растение - среда обитания”. В связи с этим целесообразно создание достаточно простых моделей процесса роста (банка таких моделей), с небольшим числом неизвестных параметров – параметров агроэкосистемы, без которых растение не может существовать, не может функционировать как система. При таком подходе выигрыш может быть достигнут за счет использования более тонких и точных математических методов идентификации и прогноза, более интеллектуального, эффективного и гибкого математического и программного обеспечения, эффективных критериев адекватности и устойчивости моделей, а также технологии моделирования.

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

Описание математической модели и процедуры ее идентификации

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

W (t) = q(t)P(t) + P1(t) — E(t) — H(t), (1)

где P(t) - величина осадков; q(t) - коэффициент использования осадков (определяется, например, экспертно или по формуле Харченко С.И. [6], через Wmin - наименьшую влагоёмкость почвы и Wz - влажность завядания); P1(t) - подпитывание (приток) из грунтовых вод; E(t) - суммарное испарение с корнеобитаемого слоя; H(t) - уровень (сток) грунтовых вод, W(t) - средняя по слою влажность почвы (с учётом поливов или на межполивной период).

Оценим и учтём влияние накопившейся к некоторому моменту времени биомассы растений на экологически обоснованную величину суммарного испарения в каждый момент времени.

Величину суммарного испарения с корнеобитаемой зоны растений представим в виде суммы интенсивности транспирации растениями E1(t) и интенсивности испарения с поверхности почвы E0(t):

E(t) = E0(t) + E1(t) (2)

Известно [1], что динамика прироста биомассы в предположении, что прирост биомассы хорошо коррелирует с интенсивностью транспирации растительного покрова, описывается уравнением:

x (t) = a(t)E1(t) — b(t)x(t), (3)

где x(t) - биомасса культуры; a(t) - эффективность транспирации; b(t) - коэффициент расхода на дыхание.

Для определения динамики накопления биомассы может быть использован банк различных моделей, из которых подбирается по тем или иным критериям адекватности наилучшая модель (по результатам идентификации).

В рассматриваемой нами процедуре моделирования будем использовать

достаточно простую для идентификации модель Ферхюльста - Вольтерра [2]:

x (t) = [ (t) —  (t)x(t)]x(t), (4)

где  - коэффициент роста (автоприроста),  - коэффициент сопротивления среды (нехватки воды).

Известно также, что динамика прироста биомассы хорошо описывается уравнением Давидсона - Филиппа [3]:

x (t)= e0(t)(F(t) — R(t)), (5)

где e0 - коэффициент перехода от массы усвоенной СО2 к сухой фитомассе; F - суммарный фотосинтез растений; R - суммарное дыхание растений.

Интенсивность дыхания за сутки зависит от величины накопившейся биомассы. Экспериментально получено [4], что

R(t) = b(t)x(t) + e1F(t), (6)

где e1 — коэффициент затрат на рост биомассы растений.

Коэффициенты  0 ,  1 - экспериментально определяемые; для большинства культур можно полагать  0 = 0,68,  1 = 0,27.

С учетом уравнений (1) — (4) имеем следующую модель расчета влажности почвы с учетом динамики накапливаемой биомассы:

W (t) = q(t)P(t) + P1(t) — E(t) — H(t),

(7)

E1 (t) =[ (t) —  (t)x(t) + b(t)]x(t)/a(t),

Из (4), (5) и (6) имеем:

b(t) = (1-e1)F(t)/x(t) — ( (t) —  (t)x(t))/e0, (8)

Для нахождения влажности почвы нам необходимо идентифицировать переменные  и  . Эта задача достаточно сложна из-за сложности и дороговизны проведения экспериментальных исследований (мониторинга). Мы продемонстрируем имитационную процедуру её решения для случая постоянных параметров модели (4); случай кусочно-постоянных параметров - аналогичен и влияет только на размерность задачи, а случай произвольных функции сводим к проблеме аппроксимации их некоторой системой базисных функций.

Решение уравнения (4), как легко проверить, имеет вид:

x(t)= /( + Ce —  t). (9)

Теперь для того, чтобы найти  и  нужно, согласно метода наименьших квадратов, решить задачу минимизации квадратичного функционала вида:

n

f( ,  , c) =  (xi0 — xi)2  min, (10)

i=1

где i - номер фазы вегетации растения (i=1,2,...,n); n - число фаз вегетации; xi0 - экспериментальные величины урожайности культуры за репрезентативный период времени; xi - теоретические величины урожайности сельхозкультур, определяемые по формуле (9).

Для нахождения решения задачи (10) необходимо решить нелинейную систему уравнений:

df / d =0, df / d = 0, df / dc = 0. (11)

Решаем эту систему численно (например, методом Зейделя), с требуемой точностью  и критерием адекватности вида:

( i+1 —  i)2 + ( i+1 — i)2 + (сi+1 — сi)2

Больше работ по теме: