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

УДК 622.272: 536.24

А.П. Неустроев, Ю.А. Хохолов

УЧЕТ ИНФИЛЬТРАЦИИ АТМОСФЕРНЫХ ОСАДКОВ И ФИЛЬТРАЦИИ ВЛАЖНОГО РУДНИЧНОГО ВОЗДУХА ПРИ НАКОПЛЕНИИ ЛЬДА В СЛОЕ ПОРОДНОЙ ПРЕДОХРАНИТЕЛЬНОЙ ПОДУШКИ В УСЛОВИЯХ КРИОЛИТОЗОНЫ

Аннотация. Для отработки подкарьерных запасов трубки «Удачная» принята система этажного обрушения с одностадийной выемкой под породной предохранительной подушкой. Для безопасного выпуска руды подушка должна сохранять сыпучие свойства в течение всего периода существования рудника. Для изучения процесса льдонакопления разработана общая математическая модель тепло- и массопереноса, учитывающая три процесса, происходящие в породной подушке: инфильтрация жидкой влаги за счет атмосферных осадков, фильтрация влажного воздуха и теплоперенос с фазовыми переходами влаги. Рассмотрено воздействие атмосферных осадков и рудничного воздуха на льдонакопление в предохранительной подушке. Математическая модель позволяет рассчитать образование льда в пустотах материала подушки за счет инфильтрации воды и конденсации влаги из влажного рудничного воздуха. Рассмотрена породная подушка толщиной 30 м. Численные расчеты проведены для хорошо проницаемых и слабопроницаемых горных пород. Разработанные математические модели тепло- и массообмена и программы расчета температурно-влажностного состояния предохранительной подушки при фильтрации влаги могут служить основой для прогноза интенсивности смерзания и для управления процессом безопасного равномерного опускания подушки.

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

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

DOI: 10.25018/0236-1493-2018-12-0-39-47

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

ISSN 0236-1493. Горный информационно-аналитический бюллетень. 2018. № 12. С. 39-47. © А.П. Неустроев, Ю.А. Хохолов. 2018.

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

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

• температура влаги и скелета породы одинакова;

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

• перенос жидкой влаги осуществляется только в жидкой фазе под действием силы тяжести и потенциала влаги;

• фильтрации влаги и воздуха происходят по закону Дарси;

• рассматривается одномерная задача, причем начало отсчета совпадает с поверхностью подушки, ось направлена вертикально вниз.

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

д_ дх

кф дР

Рв -й дх 1 дР

-К

+ 1

вертикальная координата с положительным направлением вниз от поверхности подушки, м; V — скорость фильтрации жидкой влаги, м/с.

Фильтрация воздуха через пористую среду описывается уравнением Л.С. Лей-бензона [13]:

т

д£= кп -рд2р2 дг цв дх2

(3)

где р — давление воздуха, Па; т — пористость с учетом наличия влаги, доли ед.; кп — коэффициент проницаемости породной подушки для расчета давления воздуха, м2; цв — динамическая вязкость воздуха, Па • с.

При заполнении пор пород жидкой влагой и льдом происходит уменьшение коэффициента проницаемости породной подушки, который определяется из соотношения [14]:

К

184аф (1- т)

(4)

(1)

(2)

где аф — коэффициент формы зерен породы; б — эквивалентный диаметр зерен, м.

Таким образом пористость и коэффициент проницаемости для расчета фильтрации воздуха являются переменными по времени величинами.

Скорость фильтрации воздуха (м/с) в соответствии с законом Дарси рассчитывается по формуле:

Ч--^ (5)

цв дх

Уравнениетеплопереноса в паровоздушной смеси имеет вид:

Рв -ё дх

где № — объемная влажность, доли. ед; t — время; с; Р — потенциал влаги, Па; кф — коэффициент влагопроводности, м/с; рв — плотность воды; м3/с; ё — ускорение свободного падения, м/с2; х —

(Л,Рп+ А,Р.)[

(6)

= аДт-е)+12-]

где 9 — температура воздуха, °С; сп, рп, св, рв — удельные теплоемкости, плотности соответственно пара и сухого воздуха; Дж/(кг • К), кг/м3, а — объемный коэф-

фициент теплообмена между компонентами, температура которых Т и 9, Вт/(м3 • К); 12 — теплота фазовых переходов пар-вода, Дж/кг; } — скорость фазового перехода пар-вода, кг/(м3 • с).

В соответствие с законом Дальтона для газовых смесей общее давление влажного атмосферного воздуха составляет [15]:

р = рв + рп, Па (7)

где рв и рп — соответственно парциальные давления сухого воздуха и водяного пара, Па.

Для расчета плотности пара и сухого воздуха воспользуемся уравнением состояния для идеальных газов [15]:

Рп = Кп(9 + 273,15) (8)

" - Рв (9)

Яв(9 + 273,15)

где Яп — газовая постоянная пара, равная 461,66 Дж/(кг • К); Яв — газовая постоянная сухого воздуха, равная 287,04 Дж/(кг • К).

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

рп = 600,36-ехр(0,0777-6) , Па. (10)

Формула для расчета скорости конденсации влаги имеет вид:

600,36 ■ ехр(0,0777-0)-[0,0777-(0 + 273,15)-1]

]= —

ае эе

——

18ху

(0 + 273,15)2 600,36- ехр(0,0777-0)5^ (0 + 273,15) дх

(11)

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

а = Эф • а (12)

где Эф — удельная поверхность фильтрации, м2/м3, а — коэффициент теплоотдачи, Вт/(м2 • К).

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

Ыи = 0,021 • Яв08 • Рг°'43(Рг/Ргст)0'25, (13)

где Ыи = а • б/Хв — критерий Нуссельта; б — канала, м; Рг = ц • св /А,в — критерий Прандтля; Ргст — критерий Прандтля у стенок канала; Яв —критерий Рейнольдса; Хв — коэффициент теплопроводности воздуха, Вт/(м • К).

Диаметр канала можно определить по формуле И.И. Вахромеева [14]:

б = 0,64 т0 бср + 0,38 • 10-3 м, (14)

где бср — средний диаметр кусков, м; т0 — коэффициент пустотности, доли ед.

Для расчета удельной поверхности фильтрующих пор применяется формула Козени-Кармана [13]:

m

S=m- —

'<k„

(15)

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

[C(T)

_d_ dx

+ Li •®-p-5(T-T*)]d- =

МЛ -\pwv+ av (0-T) dx dx

,(16)

C =

x =

cL p1; T < T • c2p2; T > T •

X^ T < T * X- T > T *

(17)

(18)

где Т — температура, °С; с1, с2 — удельные теплоемкости мерзлых и талых пород, Дж/(кг • К); Х^ Х2 — коэффициенты теплопроводности мерзлых и талых пород, Вт/(м • К); р — плотность пород, кг/м3;

— скрытая теплота фазовых переходов вода-лед, Дж/кг; ю — весовая влажность пород, долях ед.; Т* — температура фазовых переходов, °С; 8(Т — Т*) — дельта-функция Дирака, 1/К; сж, рж, сп, рп, св, рв — удельные теплоемкости, плотности соответственно воды, пара и сухого воздуха, Дж/(кг • К), кг/м3; 12 — теплота фазовых переходов пар-вода, Дж/кг; } — скорость фазового перехода пар-вода, кг/(м3 • с).

Система (1)—(22) замыкается путем задания зависимостей определяющих гидрофизические характеристики материала породной подушки,

Р = (19)

к = £,(№, Т). (20)

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

W(x, 0) = W0, 0 < x < L; (21) p(x, 0) = p0, 0 < x < L; (22) T(x, 0) = T0, 0 < x < L; (23)

где L — толщина подушки, м.

На верхней границе задается условие, которое учитывает интенсивность атмосферных осадков:

^ lT = N-k*' f> 0; (24)

p2-g dx

где N — скорость дождя (объем воды, выпавшей на единицу поверхности х = 0 за единицу времени), м/с.

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

Р = Ра^ t > 0 (25)

На поверхности подушки происходит конвективный теплообмен по закону Нью-тона-Рихмана. Тогда здесь задается граничное условие III рода: dT

X— = a(T - T.(t}), t > 0; (26) dx

где а — коэффициент теплообмена на поверхности подушки, Вт/(м2 • К); где T1(t) — температура атмосферного воздуха, °С.

На нижней границе задаются влажность пород и давление рудничного воздуха:

W = W2, t > 0 (27)

p = p3, t > 0 (28) где W2 — влажность пород на нижней границе, доли ед.; р3 — давление рудничного воздуха, Па.

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

X- = -aL (Т - T2(t)), t > 0 (29) дх

где а1 — коэффициент теплообмена на нижней границе, Вт/(м2 • К), где T2(t) — температура рудничного воздуха, °С.

При решении уравнения влагопро-водности (1) производную dP/dx можно выразить следующим образом:

dP dP dW ndW

— =--= M- (30)

dx dW dx dx

где M = dP/dW.

Таким образом вместо уравнений (1) и (2) получаем уравнения влагопро-водности относительно влажности W:

dW

~dt

д_ дх

\ •ё 8х

■кл

1

1

g dx

(31)

(32)

Коэффициент влагопроводности рассчитывается через коэффициент проницаемости среды по формуле [13]:

fk'

(33)

где р — плотность воды, кг/м3, g — ускорение свободного падения, м/с2; ц — динамическая вязкость воды, Па • с.

В работе [11] для расчета коэффициента влагопроводности льдистых горных пород была предложена формула:

k =

—- (34)

(1 + С щ )2

где кт — коэффициент влагопроводности для талых грунтов, м/с; Мм — объемная льдистость, доли ед.; С1 — параметр, характеризующий соотношение между размерами кристаллов льда и частиц пород, среднее значение которого равно 8.

Зависимость коэффициента влагопроводности талых пород от влажности задана по формуле [12]:

к=кА

N3,5

ym-W0;

(35)

Льдистость, доли ед. О 0,05 0,1 0,15 0,2

s о-

го

X S

3"

с;

£ 5

10-

15

20

25

30

Т

ИЮЛЬ, 1 год январь, 1 год июль, 2 год январь, 2 год

б) Льдистость, доли ед.

0 0,05 0,1 0,15 0,2

s 0

го

X S

с;

£ 5

10

15

20

25

30

ИЮЛЬ, 1 год январь, 1 год июль, 2 год январь, 2 год

■1

Рис. 1. Распределения льдистости по глубине породной подушки в разные моменты времени: хорошо проницаемые породы (а), слабопроницаемые породы (б)

Fig. 1. Ice content versus depth of rock cushion at different time: readily permeable rocks (a); weakly permeable rocks (b)

в

4 8 HiVW 12 1 2 3 4 / 5 j 1' 12 1

_Л NX MS'_ЛЧХ

О 20

о

га" % 10

го о.

ё 0 5

CD

н -10

-20

-30

-40

-50

Время, мес.

Рис. 2. Динамика температур атмосферного воздуха и породной подушки возведенной из хорошо проницаемых пород на разных глубинах

Fig. 2. Temperatures in mine air and readily permeable rock cushion at different depths

T

где кф — коэффициент фильтрации, м/с, т — пористость, доли ед; №0 — влагоем-кость прекращения движения влаги в жидком виде, доли ед.

В качестве примера рассмотрим породную подушку толщиной 30 м. Расчеты проводились при следующих исходных данных: начальная температура породной подушки —4 °С; пористость 0,2; температура наружного воздуха рассчитывалась по формулам (40, 41); температура рудничного воздуха +2 °С;

20

депрессия 100 Па. Численные расчеты проведены для хорошо проницаемых и слабопроницаемых горных пород. На рис. 3 приведены распределения льди-стости в породной подушке летом и зимой по истечении 1 и 2 лет. Как видно из графиков, при поступлении воды атмосферных осадков в предохранительную подушку происходит постепенное накопление льда в порах горных пород за счет аккумулированного холода. Так для хорошо проницаемых пород летом

о

о

го" о.

н

го

Q.

0)

10 О -10 -20 -30 -40 -50

Л

4 75^6^7 8 9 10 11 ' 12, 1 2 Г 'Л \ \ 3 4, 9\10 11 1? /X- w \\

\\\

\ 1 - 2- — T

3- ~^~на гп.0,5 м

4 — гл.1 м

Время, мес.

Рис. 3. Динамика температур атмосферного воздуха и породной подушки возведенной возведенной из слабопроницаемых пород на разных глубинах

Fig. 3. Temperatures in mine air and weakly permeable rock cushion at different depths

первого года максимальная льдистость составляет 0,1 на глубине 1,5 м; а зимой она равна 0,15 на глубине 2,2 м. Летом второго года максимальная льдистость наблюдается также на глубине 2,2 м и происходит постепенное повышение льдистости на этой глубине, что приводит к закупорке пор. Для слабопроницаемых пород льдонакопление происходит на глубине 1,9 м.

В нижней части подушки происходит постепенное протаивание пород за счет фильтрации воздуха из рудника независимо от времени года. При этом так же происходит накопление льда в подушке за счет конденсации влаги из влажного рудничного воздуха. Кривая льдосодер-жания также имеет характерный «горб», который с течением времени перемещается к верхней границе подушки. Для слабопроницаемых пород процесс пересписок ЛИТЕРАТУРЫ

движения характерного горба менее интенсивный.

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

1. Соколов И. В., Смирнов А.А., Антипин Ю. Г., Никитин И. В., Тишков М. В. Обоснование толщины предохранительной подушки при отработке подкарьерных запасов трубки «Удачная» системами с обрушением // Физико-технические проблемы разработки полезных ископаемых. - 2018. - № 2. - С. 52-62.

2. Заровняев Б. Н. Шубин Г. В., Курилко А. С., Хохолов Ю.А. Прогноз температурно-влаж-ностного состояния предохранительной подушки при отработке подкарьерных запасов руды в условиях криолитозоны // Горный журнал. — 2016. — № 9. — С. 33—36.

3. Каймонов М.В., Хохолов Ю.А. Математическое моделирование температурно-влаж-ностного режима блока отбитой руды рудников Севера // Наука и образование. — 2010. — № 1. — С. 27—32.

4. Hansson K., Simunek J., Mizoguchi M., Lundin L.-C., van Genuchten M. Th. Water flow and Heat transport in frozen soil : Numerical solution and freeze-thaw applications // Vadose zone J. 2003. Vol. 3. P. 693—704.

5. Karra S., PainterS. L., Lichtner P. C. Three-phase numerical model for subsurface hydrology in permafrostaffected regions (PFL0TRAN-ICEv1.0). The Cryosphere. 2014. Vol. 8, P. 1935—1950.

6. Painter S. Three-phase numerical model of water migration in partially frozen geological media: model formulation, validation, and applications // Computational Geosciences. 2011. Vol. 15(1). P. 69—85.

7. Zhou Y., Zhou G. Numerical simulation of coupled heat-fluid transport in freezing soils using finite volume method // Heat and Mass Transfer. 2010. Vol. 46. P. 989—998.

8. Rawlins M.A., Nicolsky D. J., McDonald K. C., Romanovsky V. E. Simulating soil freeze thaw dynamics with animprovedpan-Arctic water balance model // Journal of Advances in Modeling Earth Systems. 2013. Vol. 5. 1—17. DOI: 10.1002/jame.20045.

9. Wicky J., Hauck C. Numerical modelling of convective heat transport by air flow in permafrost talus slopes // The Cryosphere. 2017. Vol. 11. P. 1311—1325. https://doi.org/10.5194/ tc-11-1311-2017.

10. Scott L. Painter Three-phase numerical model of water migration in partially frozen geological media: model formulation, validation and applications // Computational Geosciences. 2011. Vol. 15. P. 69—85. DOI: 10.1007/s10596-010-9197-z.

11. Мотовилов Ю.Г. Численное моделирование процесса инфильтрации воды в мерзлые почвы // Метеорология и гидрология. — 1977. — № 9. — С. 67—75.

12. Рабочев И. С., Рекс Л. М., Пягай Э. Т., Якиревич А. М. Применение модели тепло-вла-гопереноса в почвогрунтах для расчета суммарного водопотребления сельскохозяйственных культур // Почвоведение. — 1981. — № 1. — С. 50—59.

13. Лейбензон Л. С. Движение природных жидкостей и газов в пористой среде. — М.: Государственное издательство технико-теоретической литературы, 1947. — 244 с.

14. Хямяляйнен В. А., Понасенко Л. П., Бурков Ю. В., Франкевич Г. С., Жеребцов В. А. Тампонаж обрушенных пород. — Кемерово: Кузбас. гос. техн. ун-т, 2000. — 107 с.

15. Лариков Н. Н. Теплотехника. — М.: Стройиздат, 1985. — 432 с.

16. Кутателадзе С. С. Основы теории теплообмена. — М.: Атомиздат, 1979. — 416 с. итш

коротко ОБ АВТОРАХ

Неустроев Аркадий Петрович1 — аспирант, e-mail: [email protected], Хохолов Юрий Аркадьевич1 — доктор технических наук, ведущий научный сотрудник, e-mail: [email protected], 1 Институт горного дела Севера им. Н.В. Черского СО РАН.

ISSN 0236-1493. Gornyy informatsionno-analiticheskiy byulleten'. 2018. No. 12, pp. 39-47. Ice accumulation in rock safety cushion,

considering atmospheric precipitation seepage and moist mine air leak in the conditions of permafrost region

Neustroev A.P.1, Graduate Student, e-mail: [email protected], Khokholov Yu.A.1, Doctor of Technical Sciences, Leading Researcher, e-mail: [email protected],

1 Chersky Mining Institute of the North, Siberian Branch, Russian Academy of Sciences, 677018, Yakutsk, Republic of Sakha (Yakutia), Russia.

Abstract. For the extraction of mineral reserves of Udachnaya Pipe under the open pit bottom, the system of one-stage caving under a safety cushion is accepted. To ensure the mine safety, the cushion should preserve looseness within the whole mine life period. Developed for the analysis of ice accumulation, the general mathematical model of heat and mass transfer takes into account three interfering processes running in the safety rock cushion, namely, seepage of moisture from atmospheric precipitation, leak of moist air and heat transfer during phase transformations of moisture. It is analyzed how atmospheric precipitations and mine air affect the process of ice accumulation in the safety cushion. The mathematical model allows calculating ice formation in voids of the cushion owing water seepage and condensation of moist mine air. The model efficiency is discussed in terms of the calculations for a rock safety cushion 30 m thick. The calculations are performed for readily and weakly permeable rocks. For the readily permeable rocks, the maximum ice content is 0.1 at a depth of 1.5 m in the first summer and 0.15 at a depth of 2.2 m in the first winter. In the second summer, the maximum ice content is also observed at the depth of 2.2 m; moreover the ice content grows, which blocks pores. In case of the weakly permeable rocks, ice accumulates at a depth of 1.9 m. At the bottom of the cushion, rocks gradually thaw under leakage of air from the mine irrespectively of the season. In addition, ice accumulates in the cushion owing to condensation of moisture from the mine air. The ice content curve contains a characteristic "hump" which drifts toward the upper boundary of the cushion with time. The hump drifting is less intensive in the weakly permeable rocks. The developed mathematical models of heat and mass transfer and the programs for calculating temperature and moisture content of the safety cushion under moisture seepage can be used to predict rate of adfreesing in rocks and to control safe and uniform subsidence of the cushion. The calculation results will help develop guidelines on making safety cushions, maintaining their year-round flowability and permeability, as well as improving safety of underground extraction of diamond pipe reserves.

Key words: rock safety cushion, hybrid open pit/underground mining, permafrost region, phase transformations, moisture seepage and air leak, temperature, mathematic modeling.

DOI: 10.25018/0236-1493-2018-12-0-39-47

REFERENCES

1. Sokolov I. V., Smirnov A. A., Antipin Yu. G., Nikitin I. V., Tishkov M. V. Obosnovanie tolshchiny predokhranitel'noy podushki pri otrabotke podkar'ernykh zapasov trubki «Udachnaya» sistemami s obrush-eniem [Substantiation of safety cushion thickness for mining Udachnaya Pipe reserves under the open pit bottom by systems with caving]. Fiziko-tekhnicheskiye problemy razrabotki poleznykh iskopayemykh. 2018, no 2, pp. 52-62. [In Russ].

2. Zarovnyaev B. N. Shubin G. V., Kurilko A. S., Khokholov Yu. A. Prognoz temperaturno-vlazhnostnogo sostoyaniya predokhranitel'noy podushki pri otrabotke podkar'ernykh zapasov rudy v usloviyakh kriolitozony [Prediction and temperature of safety cushion under open pit bottom in the conditions of permafrost], Gornyy zhurnal. 2016, no 9, pp. 33-36. [In Russ].

3. Kaymonov M. V., Khokholov Yu. A. Matematicheskoe modelirovanie temperaturno-vlazhnostnogo rezhima bloka otbitoy rudy rudnikov Severa [Mathematical modeling of temperature and moisture content of broken ore blocks in mines]. Nauka i obrazovanie. 2010, no 1, pp. 27—32. [In Russ].

4. Hansson K., Simunek J., Mizoguchi M., Lundin L.-C., van Genuchten M. Th. Water flow and Heat transport in frozen soil: Numerical solution and freeze-thaw applications. Vadose zone J. 2003. Vol. 3, pp. 693—704.

5. Karra S., Painter S. L., Lichtner P. C. Three-phase numerical model for subsurface hydrology in perma-frostaffected regions (PFL0TRAN-ICEv1.0). The Cryosphere. 2014. Vol. 8, pp. 1935—1950.

6. Painter S. Three-phase numerical model of water migration in partially frozen geological media: model formulation, validation, and applications. Computational Geosciences. 2011. Vol. 15(1), pp. 69—85.

7. Zhou Y., Zhou G. Numerical simulation of coupled heat-fluid transport in freezing soils using finite volume method. Heat and Mass Transfer. 2010. Vol. 46, pp. 989—998.

8. Rawlins M. A., Nicolsky D. J., McDonald K. C., Romanovsky V. E. Simulating soil freeze thaw dynamics with animprovedpan-Arctic water balance model. Journal of Advances in Modeling Earth Systems. 2013. Vol. 5. 1—17. DOI: 10.1002/jame.20045.

9. Wicky J., Hauck C. Numerical modelling of convective heat transport by air flow in permafrost talus slopes. The Cryosphere. 2017. Vol. 11, pp. 1311—1325. https://doi.org/10.5194/tc-11-1311-2017.

10. Scott L. Painter Three-phase numerical model of water migration in partially frozen geological media: model formulation, validation and applications. Computational Geosciences. 2011. Vol. 15, pp. 69—85. DOI: 10.1007/s10596-010-9197-z.

11. Motovilov Yu. G. Chislennoe modelirovanie protsessa infil'tratsii vody v merzlye pochvy [Numerical modeling of water seepage in frozen soil], Meteorologiya igidrologiya. 1977, no 9, pp. 67—75. [In Russ].

12. Rabochev I. S., Reks L. M., Pyagay Eh. T., Yakirevich A. M. Primenenie modeli teplo-vlagoperenosa v pochvogruntakh dlya rascheta summarnogo vodopotrebleniya sel'skokhozyaystvennykh kul'tur [Application of model of heat and moisture transfer in soil for calculating total water consumption for agricultural crops]. Pochvovedenie. 1981, no 1, pp. 50—59. [In Russ].

13. Leybenzon L. S. Dvizhenieprirodnykh zhidkostey igazov vporistoy srede [Flow of natural fluids and gases in porous medium], Moscow, Gosudarstvennoe izdatel'stvo tekhniko-teoreticheskoy literatury, 1947, 244 p.

14. Khyamyalyaynen V. A., Ponasenko L. P., Burkov Yu. V., Frankevich G. S., Zherebtsov V. A. Tamponazh obrushennykh porod [Caved rock plugging], Kemerovo, KuzGTU, 2000, 107 p.

15. Larikov N. N. Teplotekhnika [Thermotechnics], Moscow, Stroyizdat, 1985, 432 p.

16. Kutateladze S. S. Osnovy teorii teploobmena [Theoretical framework of heat exchange], Moscow, Atomizdat, 1979, 416 p.

Ha 14-M ropHO-reoflomqecKOM ^opyMe «MAHHEKC POCCHA 2018»