Методологические основы математического моделирования гидротермического режима почвы в зимний период Текст научной статьи по специальности «Науки о Земле и смежные экологические науки»

УДК 333.93

Л.А. Хворова, О.А. Иванова

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

В моделировании агроэкосистем наиболее разработанной является проблема формирования водно-теплового режима почв для теплого полугодия. В данной работе рассмотрено формирование водно-теплового режима в холодный период года. За основу взята монография Е.М. Гусева [1], в которой представлен комплекс физико-математических моделей формирования гидротермического режима почвы в холодное полугодие.

Зимний режим запасов почвенной влаги определяется в основном двумя процессами:

1. Обогащением запасов почвенной влаги за счет талых вод и атмосферных осадков в период оттепелей.

2. Внутрипочвенным передвижением воды в парообразном и жидком виде в верхние слои почвы, которое обеспечивает поступление ее в значительных количествах.

Зимний период влажности почвы для периода времени от перехода среднесуточной температуры воздуха через 5 °С осенью до перехода ее через 5 °С весной.

Холодный сезон характерен тем, что в отдельные его периоды основными являются разные стороны гидротермических процессов в почвенном профиле, поэтому для разных этапов рассматриваемого сезона предлагается [1] строить свою модель, отражающую основную физическую сущность процесса именно для этого этапа.

В работе [1] предлагается выделение следующих характерных периодов:

1. Промерзание почвы под снежным покровом, сопровождающееся миграционным подтоком воды со стороны талой зоны и длящееся фактически с момента наступления отрицательных температур до схода снежного покрова.

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

3. Весеннее снеготаяние и впитывание талых вод в почву.

4. Оттаивание почвы после схода снежного покрова.

5. Зимние оттепели.

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

Зимний период характеризуется тем, что с момента наступления отрицательных температур воздуха большая часть осадков выпадает в твердом виде (снег, иней, изморозь), накапливаясь на поверхности почвы в виде снежного покрова и гололеда, а в глубь почвы распространяется зона с температурами ниже 0 еС. В результате для полевых агроэкосистем можно выделить три характерных слоя: снежный покров, мерзлая зона почвы, ее талая зона. Для каждой из указанных зон можно записать в соответствующей форме уравнения тепло-, влагопереноса, которые являются основой описания гидрологических и гидрофизических процессов формирования почвенных вод. Далее все характеристики, используемые в подобных уравнениях и связанные с соответствующим слоем, будут иметь индексы: 1 - для снежного покрова; 2 -для мерзлой почвы; 3 - для талой. Теплогидрофизические параметры талой и мерзлой почвы зависят от типа почвы, ее плотности, влажности и льдистости.

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

продвижения фронта промерзания X [1]:

Т * ЛХ _

Т ,, _ ум Чт,

М

с \т I

Т*_ 1р№ (V - ин) + -М, (1)

где чт - кондуктивный поток тепла к границе промерзания со стороны талой зоны, чм - кондуктивный отток тепла от границы промерзания в талую зону.

68

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

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

qM=

X + H

Я2

где Н _ —^к - приведенная высота снега, а температура поверхности почвы

т _ ТсХ

п X + н '

В работе рассмотрен вариант учета теплового потока со стороны талой зоны при произвольном характере движения фронта промерзания

(необязательно г1/2 , как в задаче Стефана), использующий достаточно простой алгоритм.

В основу расчета положен интегральный метод теплового баланса [1]. При этом также будем считать, что на больших глубинах температура почвы постоянна и равна некоторой среднегодо-

V т * V

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

на проникания, аппроксимируем профиль температуры параболической зависимостью от ъ, где ъ - вертикальная координата:

Т3 = Аі 2 + & + С (2)

Здесь коэффициенты А , Ви С - функции времени.

С учетом условий на границах X и и

Тз(Х) = 0; (3)

Тъ(и) = Т * (4)

и условия непрерывности потока на и

дтз

дz

_ 0 .

Уравнение (2) приведем к виду

т3 _ т - т

(z - s)

(5)

(6)

(X - а)

Соответственно поток тепла со стороны талой зоны к границе промерзания будет равен

_ я дтз

qT — Я3

213T

(7)

ді х и - X где 13 - коэффициент теплопроводности талой зоны.

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

Приведем задачу расчета динамики глубины проникания и . При этом за нулевой момент времени принимается момент перехода температуры поверхности через 0 еС. Интегрируя по ъ в пределах от X до и уравнение теплопроводности для талой зоны

с3 _-

дтз

дt

_ я

д 2T3

3 дz2

(8)

где С3 - объемная теплоемкость талой почвы, получим

дт3

дт3

X

(9)

С учетом равенства (6) последнее уравнение приведем к виду

d(s -X) + 3^ _ 6-

dt

dt

s - X

(10)

s

c

3

я

3

s

a

3

69

Я3

где a3 _— - температуропроводность талой Где

A _ ЯзГ_ я, т.

. Заметим, что равенство (14), запи-

зоны.

Введем обозначение g _ а - X . Тогда окончательно формулировку задачи динамики глубины проникания запишем в виде

+ 6gdX _ 12а3 (11)

с начальными условиями

Х(0) _ 0, g(0) _ g0. (12)

Интегрируя уравнение (11), получим

2 2 Х(г)

g (г) - go + 6 |gdX _ 12азг. (13)

0

Используя теорему о среднем, оценим интеграл в (13):

X(t) і

}gdX »2g(t)X(t) , О2

что приводит к

g 2(t) + 3g (t)X _ gО2 + 12aзt. И окончательно

qT _

g(t) _ - — X +J4X + 12a3t + g0

2Язт *

(14)

(15)

(16)

—x2 + l2aзt+gQ -2X

L* dX _ - 12TC dt

213T *

X + H j— X2 + ^t + gо2 -fX • (17)

9

^ 2

Приняв H ° 0 , TC = Tn = const , а также

g0 = 0 , придем к задаче Стефана, точное решение которой имеет вид

| = bt1/2, (18)

где b есть решение некоторого трансцендент-

ного уравнения, содержащего ТП , Т * и теплофизические параметры почвы. Можно показать, что приближенное решение задачи Стефана, основанное на уравнении (17), также имеет вид (18), причем

санное как приближенное, в этом случае является точным.

Приняв гипотетические, но типичные по порядку величин следующие значения параметров:

Т* _ -Тс _ 10°С , 12 _ 13 _ 10-3 кал/(см град с), с2 _ с3 _1 кал/(см3град), LрV (V - иН) _ 15 кал/(см3), получим в случае точного решения задачи Стефана и ее приближенного решения, основанного на интегральном методе теплового баланса, соответственно Р _ 0,022 см/с1/2 и Р _ 0,0217 см/с1/2. Таким образом, ошибка при упрощенном расчете глубины промерзания с использованием равенства (17) при типичных по порядку величин значениях параметров в примере задачи Стефана не превышает 2%, что вполне приемлемо в гидрологических расчетах.

Рассмотрим определение величины начальной глубины проникания g0 . В первом приближении можно считать, что поверхностный слой начинает охлаждаться и соответственно развивается глубина проникновения с того момента времени, когда температура поверхности почвы становится меньше Т* В рамках выбранного приближения примем, что температура поверхности с того момента времени, когда она равна Т', до момента наступления нулевой температуры

изменяется линейно, т.е. Тп » Т*| 1--I, (20)

где т - интервал времени между указанными моментами.

Профиль температуры в пределах глубины проникания с учетом условия на границах г = 0 и г = и зададим в виде

т _ т*-i~*

(т *- тп )

(z - s)

П> S 2

(21)

Интегральная форма уравнения теплопроводности выглядит в этом случае следующим образом:

d

ds

дт3

dt Jтзdz - T3(s) = a3^t

dt

дz

-a3

дтз

дz

(22)

A

4 6

2 Язт*

і - і я1+1 Li+ Гі _ і a)2 - і La.

46

(19)

и с учетом уравнений (20) и (21) дает

d(s2) + 2 s2 і2

--------+ 2--_ !2a3

dt t 3 при начальном условии

(2З)

3

т

s

О

s

і/2

a

і/2

70

Рис. 1. Зависимость глубины промерзания почвы от высоты снежного покрова в январе 1989 г.

ст(0) = 0 . (24)

Решением задачи (23), (24) является

а = 2а3( , (25)

откуда

И о = а(**) = 2 ^ а3т* . (26)

Подставив уравнение (26) в (17), получим окончательно

qT _

213T

L* _ - 12TC

213T

dt

X + H

9 2 3

4x + i2a3(t+10)- 2 X

(27)

(28)

X T (z) 0 dz T (z)

I _J dz J Сф (T ')dT '_ J —dT J сэф (T ')dT'. dT

(30)

Здесь сэф - эффективная теплоемкость мерзлой почвы, являющаяся функцией температуры почвы,

ды н (Т)

эф

_ с0 + lPw

dT

(31)

где с0 - аддитивная объемная теплоемкость

в отсутствие фазовых переходов, ын (Т) - зависимость от температуры количества незамерзающей воды. Учитывая, что

dz

dT

Tn

(32)

где 'о = у •

Умножив уравнение (28) на ({ + Н) и проинтегрировав полученное равенство по ( от (1 до

(2 = (1 + Д( , получим следующий вариант рекуррентной формулы пошагового расчета глубины промерзания:

{((2) = -Н ((2) - +

2 212TC Dt

[X(ti)+н (t 2 )]2 —+ L (t2)

qT (ti )Dt L (t 2)

. (29)

Общее теплосодержание мерзлой зоны

dI 1 Tn T

L* _ — _- — J dT J сэф (T)dT. T

(33)

^{ 1 п о о

Интегрирование уравнения (33) с учетом (31) дает

Таблица 1 Фактическая и расчетная глубина промерзания почвы (январь 1989 г.)

Число месяца Фактическая глубина промерзания, см Расчетная глубина промерзания, см Отклонение

1 129 127 2

2 130 130 0

3 131 132 1

4 132 133 1

5 133 134 1

6 135 136 1

7 137 137 0

8 138 138 0

9 139 139 0

10 141 140 1

П

X

2

71

Тп

| ы Н (Т)ёТ

Ь = ЬршШ - Ьрш --------------

1П

-^ = Ьршш - Ьр,гыН’П -^, <34)

Тп

|ыН(Т )ёТ

-ТП о

где ыН П =----------- - среднее в диапазоне тем-

ТП

ператур (о...ТП) количество незамерзшей воды. При интегрировании учтена непрерывность по температуре зависимости ыН (ТП ) , т.е.

Ын (о) = Ш .

Таким образом, в общем случае коэффициент Ь определяется выражением (34). Однако во многих случаях для описания зависимости

-Т п {Т \

ыН (Т П) можно использовать ее линейную аппроксимацию, т.е. представить

ыН[П = а + ЪТП, (35)

где а и Ъ - аппроксимационные коэффициенты. В результате с учетом равенства (35) имеем

Ь = ЬрШ(Ш - а) - (2ЬрШЪ + со)Т- (36)

и в качестве необходимых для расчета промерзания параметров ыН и 2 соответственно полагаем

ыН = а ; (37)

с2 = со + 2ЬрШЪ . (38)

Теплопроводность снега рассчитывалась по

формуле Янсона, связывающей коэффициент теплопроводности снега в кал/(с см град) с его плот-

НОСТЬЮр1 в г/см3: 1 = 1о-3(о,о5 + 1,9р1 + бр^).

Для апробации модели глубины промерзания почвы использовались материалы наблюдений Каменской метеорологической станции.

При изучении модели глубины промерзания почвы проведены численные эксперименты по выявлению зависимости глубины промерзания

Таблица 2 Фактическая и расчетная глубина промерзания почвы (апрель 1989 г.)

Число месяца Фактическая глубина промерзания, см Расчетная глубина промерзания, см Отклонение

1 212 213 1

2 212 212 0

3 212 211 1

4 211 211 0

5 211 210 1

6 210 209 1

7 210 209 1

8 208 209 1

9 207 208 1

10 207 208 1

почвы от высоты снежного покрова. На рисунке 1 представлена зависимость глубины промерзания почвы от высоты снежного покрова по материалам наблюдений Каменской метеорологической станции.

На данном этапе рассмотрен упрощенный вариант расчета глубины промерзания почвы по

формуле { = а'1/2 + Ъ , где коэффициенты а и Ъ

определяются методом наименьших квадратов. Для первой декады января 1989 г., когда идет промерзание почвы, значения коэффициентов равняются а = 5,647449, Ъ = 121,8112, т. е. уравнение имеет вид { = 5,647449/1/2 +121,8112. Результаты расчетов приведены в таблице 1.

Для первой декады апреля 1989 г., когда идет оттаивание почвы, значения коэффициентов равняются а = -2,57Ю5, Ъ = 215,7767, уравнение

имеет вид { =-2,57Ю5(1/2 + 215,7767. Результаты расчетов приведены в таблице 2.

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

Литература

1. Гусев Е.М. Формирование режима и ресурсов почвенных вод в зимне-весенний период. - М., 1993.

2. Хворова Л.А. Задача типизации агрометеорологических условий вегетационного периода озимых

культур с использованием алгоритмов распознавания образов / Л.А. Хворова, О.А. Иванова // Материалы IX краевой конференции по математике. - Барнаул, 2006.