Сравнение математических моделей тепломассопереноса в почвогрунтах

ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ

УДК 519.63

В. И. Васильев, Ю. Г. Данилов, И. С. Еремеев, В. В. Попов, Г. Г. Цыпкин, Юйжуй Сун, Чжао Яндун

СРАВНЕНИЕ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ ТЕПЛОМАССОПЕРЕНОСА В ПОЧВОГРУНТАХ

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

Ключевые слова: температура, давление, пористость, проницаемость, капиллярное давление, диффузия, концентрация, термодинамическое равновесие, влажность.

V. I. Vasiliev, Yu. G. Danilov, I. S. Eremeev, V V Popov,

G. G. Tsypkin, Yurui Sun, Zhao Yandong

The comparison of mathematical models of heat and mass transfer in potting soil

Various mathematical models of processes in potting soil are observed. With a help of them capillary force, moisture diffusion and volume phase transfer can be considered. The comparison of numerical implementation results of three mathematical models

of frost penetration and thawing, containing brine with natural researches results got with a help of sensors, worked by scientists

from China Agricultural University, Beijing, China is held.

Key words: temperature, pressure, porosity, permeability, capillary pressure, diffusion, concentration, thermodynamical equilibrium, moisture.

ВАСИЛЬЕВ Василий Иванович - д. ф.-м. н., проф. кафедры прикладной математики и информатики Института математики и информатики СВФУ им. М.К. Аммосова.

E-mail: vasvasil@mail.ru

ДАНИЛОВ Юрий Георгиевич - к. г. н., доцент кафедры географии Института естественных наук СВФУ им. М.К. Аммосова.

E-mail: dan57@rambler.ru

ЕРЕМЕЕВ Иннокентий Спартакович - м. н. с. отдела прикладного программирования Центра вычислительных технологий Института математики и информатики СВФУ им. М.К. Аммосова

E-mail: innokentiy.eremeiev@chpc.ru

ПОПОВ Василий Васильевич - к. ф.-м. н., доцент кафедры прикладной математики и информатики Института математики и информатики СВФУ им. М.К. Аммосова

E-mail: pvvchoxoon@gmail.com

ЦЫПКИН Георгий Геннадьевич - д. ф.-м. н., вед. н. с. Института проблем механики им. А. Ю. Ишлинского РАН, г. Москва.

E-mail: tsypkin@ipmnet.ru

SUN Yurui - Ph. D., professor, Research Center for Precision Agriculture, China Agricultural University, Beijing, China (Исследовательский Центр точного сельского хозяйства, Китайского сельскохозяйственного университета, г. Пекин, Китай).

E-mail: pal@cau.edu.cn

YANDONG Zhao - professor, Dept. of Engineering, Beijing Forestry University, Beijing, China (инженерный отдел Пекинского лесного уинверситета, г. Пекин, Китай).

E-mail: yandongzh@bjfu.edu.cn

Первые исследования влияния поверхностных сил на движение воды в грунтах обобщены в классических трудах Шеддегера [1], академиков П. Я. Кочиной [2], Н. А. Цытовича [3]. В них, в частности, предлагались теория пленочной миграции, теория осмотических давлений и теория сил всасывания. При таком подходе задача сводится к нахождению эффективных коэффициентов для уравнения диффузии [5]. Другой подход математического моделирования процессов массопереноса в ненасыщенных грунтах основан на фундаментальных законах механики гетерогенных сред. Так, в работе Menot J. M. предложена формулировка фронтовой модели кристаллизации воды в ненасыщенной среде, но приведенная система уравнений является не замкнутой [6]. В статье Г. Г. Цыпкина была сформулирована замкнутая математическая модель промерзания ненасыщенного грунта, учитывающая сжимаемость воздуха, при этом действия капиллярных сил не учитывались [7]. В коллективной монографии [4] рассмотрена задача на протаивание мерзлого грунта, где для содержания влаги в области талого грунта получено уравнение первого порядка и показано, что влажность в ней практически постоянна. В работах Г. Г. Цыпкина и Calore C. [9, 10] учитывается действие капиллярных сил на поверхностях фазовых переходов. Математическая модель совместной фильтрации воды и воздуха с учетом капиллярных сил, действующих в зоне смеси, сформулирована в статье Г. Г. Цыпкина [11]. Было показано, что давление воздуха можно считать постоянным при отсутствии внешнего градиента давления. В таком случае из уравнений движения жидкости и газа следует нелинейное уравнение параболического типа для водонасыщенности, которое отличается от уравнения диффузии влаги, рассмотренной Hoekstra P. [5] способом определения коэффициента диффузии. На основе фундаментальных законов механики гетерогенных сред удалось вычислить коэффициент уравнения диффузии, который зависит от капиллярного давления, вязкости воды, пористости и проницаемости грунта [11]. В качестве примера была рассмотрена задача о влиянии капиллярного давления на перераспределение воздуха и воды при протаивании насыщенного мерзлого грунта, когда в начальном состоянии поровое полупространство X > 0 заполнено льдом. В этой работе в случае линейного приближения функции Леверетта было получено аналитическое решение задачи [11].

1. Модель 1.

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

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

д (kfw dJ dSw ] + k dfw dSw

dt dx

V Mw dSw dx J Mw dSw dx ’

0 < x < l, t >0.

Здесь Р - давление, Т - температура, v - скорость фильтрации, £ - насыщенность, k - проницаемость,

9 - пористость, ¡л - коэффициент динамической вязкости, / - относительная фазовая проницаемость флюида. Нижние индексы соответствуют: g - газу (воздуху), w - воде. Капиллярное давление по определению равно разности давлений между газом

и водой Рс = Рт - Рё. В случае смачиваемой среды Рс < 0.

Капиллярное давление Рс зависит от водонасыщен-ности и определяется соотношением Леверетта

Р ) = х№)1/2 ] ) [12].

Параметр X характеризует смачиваемость пород, а J ) - функция Леверетта. В области перед

фронтом плавления лед сосуществует с воздухом, справедлива следующая система уравнений.

Уравнение сохранения энергии в мерзлой зоне:

СтРт дТ = дХ[Лт дХ) 0< X < < > 0 .

В талой области уравнение теплопроводности имеет вид:

с'А1г=дХ1Т) ^

На границе фазового перехода выполняется условие:

Функция Леверетта [12]:

J (SJ--

0.4916667 + 1/120 8 + 0.00195, 0 < Sw < 0.5,

(1 + Sw )8

0.5 1/256 8 + 0.00195, 0.5 < Sw < 1.

(1.5 - Sw)8

Относительная фазовая проницаемость воды определяется формулой:

^ (К) =

2. Модель 2.

Область условно разделяется на зоны с фазовыми переходами и зоны без фазовых переходов. На границах раздела можно было включить Дельтафункции для учета условий на границах двухфазной зоны, но, как показывают расчеты, они на результаты практически не влияют. В итоге уравнение теплопроводности имеет вид [13], [14]:

я Эт Л Э (. ЭтЛ в Э ( Эт Л

р трЯ, т Эг Л Эх [ Эх } + Т Эх [ * Эх /

0<х<1, г >0.

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

Для коэффициента К(Бк) будем использовать формулы Аверьянова-Будаковского [15]:

Т = -ас

( к - К ^

^ тіп

К — К

V тах тіп J

К < К < К

тіп — ^ — тах ’

где п = 3.5.

Коэффициент диффузии будем вычислять по формуле Гарднера [15]:

) = А 0

можно записать в виде, зависимом от температурного поля:

¥=IК В т °

3. Модель 3.

Для описания движения воды в грунте используется уравнение Ричардса:

дт=К +дк дЬ  0< 1 < '• '>а

дг дх ^ дх ) дх

Для определения температурного поля будем использовать модель Стефана. В талой зоне имеем уравнение для температуры:

СшРш дх=дх[Лт Ш0< х < ^% г > 0-

В мерзлой зоне распределение температуры описывается уравнением:

с,р, дт=£ [л> дх ) ^)<х< 1 ’ ? > 0-

На границе фазового перехода имеем условие Стефана:

= Я,И - Я , х = %

4. Результаты натурных измерений

Измерения температуры и влажности проводились с использованием китайских приборов, измеряющих диэлектрическую проводимость почв и грунтов [16, 17]. Калибровка приборов была проведена китайскими учеными из Пекинского сельскохозяйственного университета, дальнейшие измерения и пересчет показаний приборов в единицы измерения влажности и температуры производился аспирантами и учеными ИМИ СВФУ Измерения влажности и температуры почв и грунтов производились периодически один раз в неделю до максимальной глубины заложенных скважин через каждые 10 см, начиная с глубины

10 см. Глубина заложения скважин зависит от глубины протаивания и колеблется от 150 до 180 см. Для равноценного сравнения результатов учитывались данные только до 150 см (наименьшей максимальной глубины заложенных скважин). Скважины приборов расположены в различных ландшафтно-геоморфологических условиях. Первая - на охраняемой площадке полигона ИМЗ СО РАН. В ландшафтном отношении это средняя часть склона небольшого (60х70 м) термокарстового понижения глубиной около 2 метров, расположенного на второй надпойменной террасе р. Лены. Трансаккумулятивные условия, преобладают супесчаные почвы. Растительность естественная - мезофитный разнотравный луг. Глубина залегания грунтовых вод - около 2 м. Полигон расположен на краю города Якутска, огорожен, вблизи на расстоянии не менее 100 метров отсутствуют какие-либо строения. Ближайшие постройки преимущественно одноэтажные, поэтому существенного влияния на розу ветров и выпадение осадков они не имеют. Полигон расположен менее чем в 1 км к западу от Якутского управления Федеральной службы по гидрометеорологии и мониторингу окружающей природной среды, т. е. данные этой организации можно использовать в работе без интерпретации. Кроме того, в 2013 г. на этом же полигоне будет установлена метеоплощадка лаборатории метеорологии кафедры географии ИЕН СВФУ, что позволит вести мониторинг метеорологических элементов и их влияния на температуру и влажность почв непосредственно

на месте. Глубина заложения скважины 150 см. Вторая скважина расположена на северной окраине г. Якутска, на охраняемой территории фермерского хозяйства Г. Эма. В ландшафтном отношении это 1-я надпойменная терраса. Элювиальные условия, преобладают супесчаные почвы. Скважина расположена на краю обрабатываемого сельскохозяйственного поля. Ближайшая метеостанция расположена на территории Якутского аэропорта в 2 км южнее. Результаты данной метеостанции с высокой степенью корреляции можно интерпретировать для исследований. Глубина заложения скважины 170 см. Третья скважина расположена на территории студенческого городка СВФУ на охраняемой территории. В ландшафтном отношении это 2-я надпойменная терраса. Элювиальные условия, преобладают супесчаные почвы. Скважина расположена возле учебного корпуса КФЕН СВФУ в зоне теневого воздействия здания. В зимнее время здесь накапливаются сугробы снега, что будет влиять как на степень промерзания, так и длительность оттаивания в весенний период. До ближайшей метеостанции Якутского управления Федеральной службы по гидрометеорологии и мониторингу окружающей природной среды 1,5 км. В этой точке необходимо учитывать отепляющее и ветроломное значение здания. Глубина заложения скважины 160 см. Четвертая скважина расположена за городом на территории дачного участка. В ландшафтном отношении это 1-я надпойменная терраса. Скважина расположена на верхней трети склона старичного террасового озера в трансэлювиальных условиях, преобладают супесчаные почвы. Глубина залегания грунтовых вод более 2 метров. Нижняя часть скважины расположена в водонасыщенных песках. Глубина заложения скважины 180 см. На изменение температуры и влажности почвогрунтов в скважинах оказывают влияние температура и влажность воздуха, а также количество и интенсивность выпадающих осадков, мощность снежного покрова, его плотность и длительность залегания, интенсивность образования и таяния. Кроме того, несомненно, оказывает влияние поверхностный и подземный сток, а также механический состав почв и грунтов, влияющий на водо- и воздухопроницаемость. От водонасыщеннос-ти почвогрунтов, от их минералогического состава зависит тепло- и электропроводность. Все эти факторы, в конечном итоге, влияют на диэлектрическую проводимость, и их необходимо учитывать при проведении замеров и расчетов. Делать обобщающие выводы еще рано в связи с коротким рядом наблюдений, малым количеством осадков, выпавших в наблюдаемый период, а также отсутствием учета испаряемости. Также, чтобы уменьшить или учесть при расчетах влияние внешних факторов, необходимо

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

5. Результаты расчетов

Расчеты проводились при значениях исходных параметров:

ф = 0.25, q = 3.34-105Дж/кг, pw =1000кг/м3, p¡ =910кг/м3, ps = 2000кг/м3, Л = 2Вт/

Cw = 4200Дж/(кг• K), C, =2000Дж/(кг• K), Cs =1000Дж/(кг• K),

=0.002Па• с, х = -10-3Дж/м2, k = 10-17м2.

На рис. 1 приводится сравнение распределений влажности, соответствующие первой, второй и третьей моделям. На рис. 2 показаны распределения температуры для первой и второй моделей.

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

Численные эксперименты показали, что распределение влажности наиболее точно по характеру описывают модель Г. Г. Цыпкина [11], учитывающую действие капиллярных сил, и модель, приведенную в монографии Кучмента [15]. По значениям более точной оказалась математическая модель Г. Г. Цып-кина [11] (на рис. 1 - средняя линия).

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

Работа выполнена при финансовой поддержке РФФИ (грант 12-01-98514-р-восток-а).

Л и т е р а т у р а

1. Шейдеггер А. Физика течения жидкостей через пористые материалы. - М.: Гостоптехиздат, 1960. - 250 с.

2. Кочина П. Я. Теория движения грунтовых вод. -М.: Наука, 1977. - 664 с.

3. Цытович Н. А. Механика мерзлых грунтов. -М.: Высшая школа, 1973. - 446 с.

4. Жесткова Т. Н. Формирование криогенного строения грунтов. - М.: Наука, 1982. - 216 с.

5. Hoekstra P. Moisture movement in soils under temperature gradient with the cold side temperature below freezing // Water Resources Res. 1966. V. 2. No 2. P. 241-250.

6. Menot J. M. Equations of frost propagation in

О 0.2 0.4 0.6 0.8

х, m

Рис. 1. Распределение влажности

1 1 1 1 1 1 Модель 1 *

Модель 2 0

Измерения —

^eOGOOOOCJOOOOOOOOOOOOOOOOOUOOOOOOOOOOOOOOOOGOOO '

XX .. ... xxXxXXX XX xxXxX X

- 1 1 1 1

0 0.2 0.4 0.6 0.8

х, m

Рис. 2. Распределение температуры

1 1.2 1.4

unsaturated porous media // Eng. Geol. 1979. V. 13. No 1-4. P. 101-109.

7. Цыпкин Г. Г. Линейная задача о фазовых переходах вода-лед в ненасыщенных грунтах // Изв. АН СССР. МЖГ. - 1990. - № 3. - С. 68-73.

8. Васильев В. И., Максимов А. М., Петров Е. Е., Цыпкин Г. Г. Тепломассоперенос в промерзающих и протаивающих грунтах. - М.: Наука. Физматлит, 1996. - 224 с.

9. Tsypkin G. G., Calore C. Role of capillary forces in vapour extraction from low permeability, water saturated geothermal

reservoir // Geothermics. 2003. V. 32. P. 219-237.

10. Tsypkin G. G., Calore C. Influence of capillary forces on water injection into hot rock, saturated with superheated vapour // Int. J. Heat and Mass Transfer. 2007. V. 50. No 15, 16. P. 3195-3202.

11. Цыпкин Г. Г. Влияние капиллярных сил на распределение влагонасыщенности при протаивании мерзлого грунта // Изв. АНСССР. МЖГ. - 2010. - № 6. - С. 130-140.

12. Leverett M. C. Capillary behavior in porous media // Trans. A.I.M.E. 1941. V. 142. Р. 341-358.

13. Ентов В. М., Максимов А. М., Цыпкин Г. Г. Об образовании двухфазной зоны при кристаллизации смеси в пористой среде // Докл. АН СССР. 1986. Т. 288, № 3. С. 621- 624.

14. Васильев В. И., Попов В. В. Численное решение задачи промерзания грунта // Матем. моделирование. -2008. - Т. 20. - № 7. - С. 119-128.

15. Кучмент Л. С. Модели процессов формирования речного стока. - Л: Гидрометеоиздат, 1980. - 142 с.

16. Image-based comparison between a-ray scanner and a dual-sensor penetrometer technique for visual assessment of bale density distribution Sun, Y (College of Information and

Electrical Engineering, China Agricultural University, 100083 Beijing, China); Cheng, Q.; Meng, F.; Buescher, W.; Maack,

C.; Ross, F.; Lin, J. Source: Computers and Electronics in Agriculture, v. 82, P. 1-7, March 2012 Database: Compendex.

17. A novel method to determine the volume of sensitivity for soil moisture sensors Sun, Yurui (College of Information and Electrical Engineering, China Agricultural Univ., 100083 Beijing, China); Sheng, Wenyi; Cheng, Qiang; Chai, Jin; Yun, Yuliang; Zhao, Yandong; Xue, Xuzhang; Lammers, Peter Schulze; Damerow, Lutz; Cai, Xiang Source: Soil Science Society of America Journal, V. 76, № 6, P. 1987-1991, 2012 Database: Compendex.