Агрохимия, 2020, № 4, стр. 85-93

ВВЕДЕНИЕ

Кадмий (Cd) входит в группу токсичных тяжелых металлов (ТМ), особенностью которого является накопление в почве и аккумуляция в растениях. По трофической цепочке он накапливается в организмах животных и человека, что создает техногенную угрозу жизни народонаселения [1]. Интенсивные выбросы промышленности, транспорта и использование различных химических средств (в том числе и удобрений) привели к накоплению ТМ на значительных территориях, что отрицательно влияет на почву, растения и другие живые организмы, и, как следствие, возрастает фоновый уровень содержания тяжелых металлов в биосфере [2].

Поступление Cd в почву определяют ее температурный режим и режим влажности, наличие в ней неорганических и органических веществ. Например, при снижении кислотности почвы с рН 7.0 до 5.5 содержание металла в корне райграса возрастает в 4 раза. Кроме того, способность растений к поглощению Cd во многом зависит от их видовых особенностей [3, 4]. В частности, среди культурных видов злаков в гораздо большей степени его поглощают рис, рожь и ячмень, и в меньшей степени – кукуруза и пшеница [5, 6].

Для снижения содержания доступного Cd в почве и растениях имеются специальные регламенты и контролируемые государством нормы. Они определяют предельно допустимые концентрации (ПДК) элементов в почве и растениях. В России, Европейском Союзе (ЕС) и других экономически развитых странах обсуждается его допустимое содержание в фосфорных минеральных удобрениях, т.к. Cd имеет тенденцию к преимущественному накоплению в фосфоритах. Номенклатура доз изменяется в зависимости от содержания кадмия: до 20 (слабо загрязненные), от 40 до 60 (средне загрязненные) и от 60 и более (сильно загрязненные) мг Cd/кг Р₂О₅ [7, 8].

Для внедрения этих правил и норм требуются сведения о реакции системы почва–удобрение–растение в ответе на внесение удобрений с различным содержанием Cd. Интерес представляет его содержание в почве и растениях после окончания вегетационного периода [9]. Такие сведения могут оказаться доступными в 2-х случаях: 1 – когда уже имеется репрезентативный экспериментальный материал, 2 – при наличии некоторого аппарата/механизма, посредством которого можно получить интересующий прогноз состояния этой системы при соответствующем на нее воздействии [10].

Интенсивность возделывания сельскохозяйственных культур обусловливает применение разнообразных минеральных удобрений, что может привести (и приводит) к фрагментарности и недостаточности экспериментального материала (случай 1). Поэтому актуальной является задача построения аппарата, при помощи которого реализуется случай 2. Для его построения требуется математическая модель совместной динамики биомассы и содержания в ней этих элементов питания.

В практике математического моделирования обычно постулируется логистический характер распределения биомассы растения и содержания в ней элементов питания. Однако при ее использовании как аппроксимации эмпирического распределения принимается допущение: удельный рост некоторой характеристики объекта определяется ее приростом и величиной самой характеристики, ограничивающей этот прирост. Так вот, если для содержания в растении элементов его питания принимается логистическое распределение, то получается следующее: находясь уже в самом растении, они сами себя и ограничивают. Но ведь это функция непосредственно растения и ризосферы, а не элементов питания, содержащихся в самом растении. Кроме того, их источник находится в ризосфере, а не в растении.

Согласно сказанному, актуальной представляется задача построения модели динамики элементов питания растений, оценка ее коэффициентов и оценка адекватности модели. Актуальна также задача проведения на основании моделей вычислительных экспериментов для оценки периода времени, при котором содержание тяжелого металла (ТМ) Cd в растении не превышает установленной нормы ПДК.

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

МЕТОДИКА ИССЛЕДОВАНИЯ

Все расчеты и вычислительные эксперименты выполняли на основании экспериментальных распределений, которые характеризуют влияние элементов питания растений (K, P, N, Ca, Mg, Cd) на рост яровой пшеницы сорта Ленинградская-97 на дерново-подзолистой почве. Структурированные цифровые распределения были сформированы на основании работ [11–14].

Объект изучения представлен распределениями биомассы надземной части растения и содержания в этой части элементов минерального питания. Методы обработки этого материала основаны на приемах математического моделирования и технологии решения задач нелинейного оценивания параметров средствами системы MATLAB.

В работе использовали экспериментальные распределения биомассы и соответствующие ей распределения выносов химических соединений элементов (K, P, N, Ca, Mg, Cd). Они были сформированы при экспериментальных исследованиях на яровой пшенице [15]. Методика их проведения состояла в следующем. На окультуренной дерново-подзолистой средне суглинистой почве выращивали яровую пшеницу сорта Ленинградская-97 (вегетационный период в зависимости от районов возделывания и погодных условий меняется от 85 до 110 сут [16]). Почва имела нейтральную реакцию среды, высокое содержание обменного калия и подвижного фосфора, низкое – содержание гумуса (1.5%), характеризовалась высокой насыщенностью основаниями и не нуждалась в известковании. При закладке опыта содержание кислорастворимых соединений Cd в почве было равно 0.5 мг/кг. Для повышения валового содержания ТМ и возможности оценки влияния большей концентрации на общий ценоз замкнутой системы в нее искусственно был внесен раствор сульфата кадмия из расчета 1 мг Cd/кг дерново-подзолистой почвы. В экспериментах рассматривали отклик системы на внесение различного количества калия, который применили в виде сульфата калия (K₂SO₄, 46% д.в. в дозах от K₂₀ до K₂₂₀). Для обеспечения комплексного питания растений применяли азотные удобрения (сульфат аммония, 21% д.в. в дозе N₁₅₀) и фосфорные удобрения (Р_с, 26% д.в. в дозе Р₁₀₀).

В опытах использовали сосуды Кирсанова, где масса почвы составила 5.5 кг при ее влажности 50% (из чего следует, что сухая масса равнялась 2.25 кг) [17, 18]. Посев пшеницы производили проросшими семенами по 25 растений на сосуд. После прореживания в каждом сосуде было оставлено по 20 растений. Влажность почвы поддерживали на оптимальном уровне (70% ПВ) за счет контролируемого полива, температуры – за счет специально отапливаемых тепличных помещений. Всего было 8 вариантов, где на каждый из них приходилось 9 временных точек отбора образцов растений пшеницы на 10, 14, 19, 24, 29, 37, 47, 57, 67-е сутки после посева. Отчет времени отбора начали с момента первых всходов. В опытах определяли сухую надземную массу растений, а также вынос этой частью растений как основных элементов питания (макроэлементов), так и микро- и мезоэлементов (включая такой элемент-антагонист как кадмий) по стандартной методике [19].

Модель динамики биомассы надземной части растений и динамики выноса элементов его питания. В практике изучения агроценозов широкое применение получили регрессионные уравнения [20, 21], где зависимой переменной выступает содержание ТМ в растении, а предикторами – их подвижные соединения в почве. В частности, для Cd в [22] приведены такие уравнения (табл. 1).

Проблемы использования этих уравнений для прогноза зависимости выноса Cd от его содержания в почве вызваны масштабом времени их корректного применения. Действительно, согласно им, внесение в почву удобрений сразу приводит к изменениям его содержания в растении, что уместно только для всего вегетационного периода. В связи с этим актуальна задача построения моделей, где реализован механизм динамики выноса для более адекватной оценки процессов в биогеоценозе.

При построении модели принимали такие допущения:

1 – удельный прирост в единицу времени биомассы $\{ {{B}_{k}},k = 1{\kern 1pt} - {\kern 1pt} 8\} $ определяется ее естественным приростом, внутривидовой конкуренцией и текущим содержанием в ризосфере элементов питания растения;

2 – удельный прирост содержаний элементов питания $\{ {{C}_{{i,k}}}$ $i = 1{\kern 1pt} - {\kern 1pt} 6$, $k = 1{\kern 1pt} - {\kern 1pt} 8\} $ в растении пропорционален их текущим содержаниям в ризосфере (допускается влияние на этот прирост остальных элементов);

3 – коэффициенты уравнений не зависят от опытов,

4 – модель охватывает весь вегетационный период растения, где временной масштаб составляет 1 сут.

Согласно этим допущениям, запись модели совместной динамики биомассы растения и элементов его питания принимает вид:

где i = 1–6 и k = 1–8, ${{a}_{{0,1}}}$ – коэффициент естественного прироста, ${{a}_{{0.2}}}$ характеризует внутривидовую конкуренцию, $\{ {{C}_{{i,k,0}}}$, $i = 1{\kern 1pt} - {\kern 1pt} 6$, $k = 1{\kern 1pt} - {\kern 1pt} 8\} $ – начальное содержание элементов питания в ризосфере, $\{ {{a}_{{i,1}}}$, ${{a}_{{i,2}}}$, $i = 1{\kern 1pt} - {\kern 1pt} 6\} $ характеризуют доступность элементов ризосферы и их динамику перехода из ризосферы в растение, $\{ {{b}_{{i,j}}}$, ${{d}_{{i,j}}}$, $i,j = 1{\kern 1pt} - {\kern 1pt} 6\} $ характеризует влияние элементов ризосферы на потребление растением элементов его питания, $\{ {{B}_{{k,0}}}$, $k = 1{\kern 1pt} - {\kern 1pt} 8\} $, $\{ C_{{i,k}}^{{\left( 0 \right)}}$, $i = 1{\kern 1pt} - {\kern 1pt} 6$, $k = 1{\kern 1pt} - {\kern 1pt} 8\} $ – начальные показатели биомассы и содержаний в растении элементов питания.

Когда в ризосфере содержание элементов питания достаточно низкое или оно уже полностью израсходовано (выражение в круглых скобках равно нулю), кривая биомассы выходит на свой предельный уровень ${{a}_{{0,1}}}{\text{/}}{{a}_{{0,2}}}$.

Оценка параметров модели. Разностные схемы системы (1) линейны по коэффициентам. Поэтому для их нахождения можно применить метод наименьших квадратов (МНК). Подходящие для его применения разностные схемы принимают вид

где t – текущие сутки и суточный шаг по времени.

Если в (2) ограничиться только первыми двумя членами его правой части, то вычисленные МНК доверительные интервалы коэффициентов принимают величины: a_0, 1 = (1.546 ± 0.101) × 10^–1 сут^–1 и a_0, 2 = (–0.033 ± 0.002) × 10^–1 (сут г)^–1.

Биомасса следует логистическому уравнению

Тогда к окончанию вегетационного периода (при больших t) биомасса растения близка к своей максимальной величине $B{\text{*}} = {{a}_{{0.1}}}{\text{/}}{{a}_{{0,2}}}$.

На рис. 1 представлены распределения: а – сглаженное выборочное распределение, б – модельное распределение для всего вегетационного периода роста пшеницы Ленинградская-97 на дерново-подзолистой почве (как правило, 90 сут [23]). Актуальность построения модели распределения для всего вегетационного периода обусловлена важностью задач прогноза урожая. Для этого случая искомыми являются начальная величина биомассы ${{B}_{{k,0}}}$, ${{a}_{{0,1}}}$ и ${{a}_{{0,2}}}$. Их оценка следует решению задачи поиска минимума функционала, который выступает мерой отклонений выборочных и модельных распределений. Ее решение определяет такие доверительные интервалы: B_k, 0 = 0.265 ± ± 0.034 г, a_0, 1 = 0.132 ± 0.024 сут^–1 и a_0, 2 = (–0.054 ± ± 0.008) × 10^–1 (сут г)^–1.

Анализ выборочного распределения биомассы (рис. 1а) указывает независимость биомассы от схемы опытов. Действительно, для них практически совпадают первые и конечные уровни. Форма кривых отвечает логистической или сигмоидальной форме (4), особенность которых состоит в наличии точки перегиба и выходе этих кривых на плато. Распределения показывают высокое подобие модельных и выборочных кривых. Корреляция случаев (а) и (б) составляет 0.995, корреляция (в) и (г) – 0.994. При ее расчете использовали только те строки модельного массива, которые приходятся на сутки проведения выборочных экспериментов.

Если в (3) ограничиться членами только в одной круглой скобке правой части (которые отвечают только Cd), то оценка МНК приводит к таким доверительным интервалам: b_6, 6 = (6.358 ± ± 0.532) × 10^–5 (сут мг)^–1 и d_6, 6 = (–1.027 ± 0.039) × × 10^–3 (сут мг)^–1. Тогда корреляция между выборочным и модельным распределениями составляет 0.942. Для этого случая решение второго уравнения (1) принимает вид

где первое слагаемое отражает влияние начального распределения, а второе – механизм динамики выноса элементов. Соотношение (5) допускает предметную интерпретацию. Действительно, поскольку биомасса растения – монотонно возрастающая функция, то при большом t показатель экспоненты стремится к –∞, а сама экспонента – малое число. Тогда вынос элементов определяется их содержанием в ризосфере ${{C}_{{6,k,0}}}~$и соотношением между показателями доступности элементов питания из ризосферы и темпа их утилизации.

Для модели распределения Cd в течение его вегетационного периода следует определить коэффициенты ${{C}_{{6,k,0}}}$, ${{b}_{{6,6}}}$ и ${{d}_{{6,6}}}$. Их оценки следуют решению задачи поиска минимума меры отклонений экспериментальных и модельных распределений. Ее решение определяет доверительные интервалы: $C_{6}^{{\left( 0 \right)}}$ = 0.02 ± 0.002 мг, b_6, 6 = (2.319 ± ± 0.251) × 10^–5 (сут мг)^–1 и d_6, 6 = (–4.149 ± 0.052) × × 10^–3 (сут мг)^–1. Анализ выборочных и модельных распределений показывает их высокую степень подобия. Так, корреляция случаев (а) и (б) составляет 0.995, корреляция (в) и (г) – 0.994.

Вычислительные эксперименты. Постановка вычислительных экспериментов ставится следующим образом. При оптимальной для роста растения температуре и влажности сосуд радиуса 10 см заполнен окультуренной дерново-подзолистой среднесуглинистой почвой массой 5.5 кг [12, 15]. В его почве нет вредных примесей ТМ, т.е. при $~t = 0$ ${{C}_{{i,k,0}}} = 0$. Принимается, что каждый год доза внесения ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}~$ в посев составляет до 120 кг/га [18, 24]. Тогда для одного сосуда эта норма в среднем составляет 377 мг. Нормы загрязнения фосфорных удобрений определяются как слабо загрязненные (20 мг Cd/кг ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}$), средне загрязненные (40 мг Cd/кг ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}$) и сильно загрязненные (60 мг Cd/кг ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}$). В этом случае они равны 0.008, 0.016 мг и 0.024 мг Cd соответственно. Посев и сбор урожая выполняли раз в год без технических обработок (мелиорации, обработки удобрениями, классической, минимальной, нулевой обработками и т.п.). В почве ТМ закрепляются в верхнем почвенном слое гумусовыми веществами, глинистыми минералами, карбонатами и минералами с переменным зарядом [25]. В почве Cd сразу связывается в нерастворимые соли, либо адсорбируется на корнях растений. Выветривание, вымывание, а также испарение с поверхности почвы для ТМ практически отсутствуют. Такая имитация достигается наличием сосуда, где потери Cd имитируют природные процессы. Поэтому эта постановка эксперимента обеспечивает наиболее достоверную оценку содержания и действия кадмия (имитация изолированного лабораторного эксперимента) на биогеоценоз.

Вычислительный эксперимент поставлен таким образом: организуют цикл вычислений для определения года, когда содержание Cd в ризосфере окажется равным или превысит заданный уровень ПДК. Согласно модели (1), растение в течение отдельного года потребляет часть внесенного в почву Cd. Максимальное число итераций задано равным 600 годам.

РЕЗУЛЬТАТЫ И ИХ ОБСУЖДЕНИЕ

Результаты экспериментов приведены на рис. 2. Случаи (а), (г) и (ж) характеризуют динамику аккумуляции Cd в ризосфере за период проведения экспериментов (600 лет); (б), (д) и (з) – динамику выноса Cd в конце вегетационного периода растений для каждого года эксперимента; (в), (е) и (и) – распределение Cd за вегетационный период в последний год эксперимента (за сезон 599–600 гг.). Форма кривых на рис. 2 показывает, что аккумуляция Cd в растении и ризосфере имеет нелинейный характер, который насчитывает 3 этапа: его линейная зависимость, постепенное снижение скорости накопления Cd в форме изгиба и стремление к постоянному уровню (плато), которое можно объяснить эффектом “насыщения”. Его интерпретация состоит в наличии такого момента, когда в растении “емкость” содержания элементов достигает предела и тем самым уравновешивает поступление элементов извне (внесение удобрений).

Таким образом, рост содержания Cd в почве и растениях прекращается и держится на постоянном уровне. Чем выше концентрация исходного элемента, тем быстрее наступает данный эффект. Такой результат подтверждается рядом авторитетных источников [3, 5, 20, 26].

Согласно расчетам для слабо загрязненных удобрений, содержание Cd в ризосфере достигнет ПДК России и Казахстана за 271 год, а для средне загрязненных и сильно загрязненных – за 89 и 55 лет соответственно. Токсичность растений в последнем варианте уже на первых итерациях на несколько порядков превосходит остальные. Такие изменения при применении удобрений с повышенным содержанием Cd будут губительным образом оказывать влияние для человека и животных, употребляющих данный продукт.

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

Следуя проведенным экспериментам, далее рассмотрен случай прогноза превышения ПДК в отдельных экономически развитых странах. Согласно табл. 2, для слабо загрязненных удобрений не наблюдается превышение допустимых норм в рамках заданных итераций (600 лет). Вариация элементов для средне загрязненных удобрений весьма значительна, а для сильно загрязненных удобрений не отмечено превышения только у исключительного ряда стран с высокими уровнями ПДК. Следует отметить, что в некоторых странах, таких как Австрия, Великобритания, Германия, ПДК значительно превышены по сравнению с другими странами, тогда как самые низкие пределы отмечены в Казахстане и России.

Понятно, что полученные результаты позволяют достаточно грубо оценить сценарии безопасной эксплуатации земель сельскохозяйственного направления. Если речь идет о проводимых экспериментах в сосудах, то вероятна ошибка измерений, связанная с большим различием процессов, которые происходят в реальных условиях в полях и вегетационных опытах в сосудах. К тому же испытания проводили всего за один период вегетации. В таком случае достоверность результатов может быть получена только на основании большой совокупности данных, полученных в полевых условиях с динамикой роста и развития растений, выноса элементов как минимум в несколько десятков лет.

На основании проведения подобных вычислительных экспериментов становится возможным оценить антропогенное влияние на среду в зависимости от определенного промежутка времени, при существенном снижении затрат на проведения полевых экспериментов. Стоимость таких работ варьирует от 500 тыс. до нескольких млн руб. за один полевой опыт.

ЗАКЛЮЧЕНИЕ

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

В рамках модели показана достоверность положения о том, что количественно динамика выноса Cd в надземную часть растения определяется биомассой этой части и его текущим распределением в ризосфере. Динамика биомассы этой части растений следует логистическому уравнению.

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

Эта оценка показала, что удобрения с низким содержанием Cd (до 20 мг Cd/кг ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}$) практически не влияли на состояние почвы и растений, тогда как средне (от 40 до 60 мг Cd/кг ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}$) и сильно загрязненные (>60 мг Cd/кг ${{{\text{P}}}_{2}}{{{\text{O}}}_{5}}$) оказывали существенное влияние на биогеоценоз.