Метод вычисления коэффициента эффективной температуропроводности по данным измерений температуры почвы

Управление тепловым режимом почв, разработка различных почвенных конструкций, анализ термических эффектов, связанных с рыхлением или уплотнением поверхностных слоев почвы, применением мульчи из сыпучих материалов и многие другие явления и процессы связаны с использованием прогнозных математических моделей теплопереноса в почвах. Для разработки и функционирования таких моделей необходимым экспериментальным обеспечением является функция температуропроводности, получение которой возможно с помощью методов, основанных на решении прямых и обратных задач теплопереноса. Кроме того, знание функции температуропроводности как функции температуры среды позволит изучить пампинг-эффект, возникающий при распространении температурных волн внутри почвы и имеющий и другие приложения в других средах и рассмотренный в достаточно простых случаях [1,2]. В работе [3] предложен ряд формул для определения коэффициента температуропроводности в различных ситуациях. При этом предполагается, что коэффициент постоянен внутри слоя, что является иногда достаточно грубым предположением, так как известно, что в пределах от -0.5 до 0 коэффициент температуропроводности резко меняется в разы, связанный с фазовыми превращениями воды и льда, содержащимися в почве [4]. Кроме того, во многих работах [5-8] для определения коэффициента температуропроводности используется упрощенное, линейное уравнение теплопроводности, что не позволяет учесть реально нелинейные процессы тепло-переноса. Формулы для температуропроводности в случае одной суточной гармоники предложены в работе [9], но опять-таки в предположении о постоянстве этого коэффициента внутри слоя. В некоторых работах предлагаются формулы для расчета коэффициента температуропроводности, в которых учитываются такие характеристики, как влажность, плотность, температура, степень дисперсности и содержание органического вещества [3]. В данной работе на основе обратного метода решения нелинейного уравнения теплопроводности предлагается другая формула для расчета данного коэффициента, которая приводит к сопоставимым результатам.

Постановка задачи и метод решения

Есть данные измерений температуры почвы на различных глубинах с одинаковой частотой измерений на протяжении достаточно длинного периода времени (порядка года). На практике встречаются измерения с частотой от 15 мин. до 4 ч., это зависит от объекта, на котором проводятся мониторинг температуры почвы. Требуется оценить значения коэффициента эффективной температуропроводности. Численные методы далеко не всегда подходят для оценки этого коэффициента, поэтому желательно получить аналитическую формулу.

Будем рассматривать суточные и полусуточные колебания температуры, т.е. температуру почвы представим в виде:

T = T 0 + A о e ^{e 0} z sin( ш о t + ф о)+

+ A i e ' ¹ z sin( ш 1 t + ф i) , (1)

где T 0 – среднесуточное значение температуры (является функцией глубины на масштабах времен порядка сутки и к тому же и времени при рассмотрении масштабов времени более суток), A ₀ – амплитуда суточных колебаний температуры на поверхности почвы, β 0 – коэффициент затухания колебаний с глубиной, A ₁ – амплитуда полусуточных колебаний на поверхности почвы, β ₁ – коэффициент затухания колебаний с глубиной, ω 0 – частота суточных колебаний, выраженная в часах, ω ₁ – частота полусуточных колебаний ш 1 = 2 ш ₀ . При этом функции ф _о ( z ) и Ф 1( z ) можно представить в виде:

Ф о( z )= a + bz, ф i( z ) = a 1 + b 1 z. (2)

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

Для определения коэффициента эффективной температуропроводности будем использовать нелинейное уравнение теплопроводности, записанное в виде:

= d (KtT) ’ (3) ∂t ∂z ∂z где ось z направлена вертикально вниз, z = 0 -поверхность почвы. Естественно потребовать от решения выполнения следующего условия: ∂T /∂z → 0 при x ^ го. Учитывая (1) и (2), можем многие слагаемые уравнения (3) вычислить аналитически.

Рассмотрим временной интервал, равный суткам. Предположим,что можно пренебречь изменчивостью среднесуточного значения во времени. Среднесуточную температуру представим в виде:

T о = Ce ^{- Yz} + d. (4)

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

времени и глубине могут быть записаны в аналитическом виде:

∂T

= A о e ∂t

cos( ш 0 t + ф о) ш о+

+ A i e ^{в 1} z cos(2 ш о t + ф i)2 ш о .        (5)

— = -CYe ^Yz - в о A о e ^{в 0} z sin( ш о t + ф о) +

∂z

+ A о e 0° cos( ш о t + ф о) b—

• — в i A i e ^{в 1} z sin(2 ш о t + ф i) +

+ A _i e ^{- e 1} z cos(2 ш _о t + ф _i ) b _i .              (6)

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

• ^z ∂T ^z ∂ ∂T

J. »^ds = Lд-AKT ⁽⁷⁾

Исходим из определения несобственного интеграла

• ^z ∂T             ^z ∂T

ds = li^m / ds •      (8)

∞ ∂t      L→∞ L ∂t

Вычислим сначала левую часть равенства (7)

z dT^-

∞ ∂s

lim

L→∞

UL ( A

0 e

β 0 z

cos( ш о t + ф о) ш о+

+ A i e ^{в 1} z cos(2 ш о t + ф i)2 ш о) ds ,      (9)

где s – переменная интегрирования. Интегралы в правой части выражения (9) берутся аналитически дважды интегрированием по частям. Еще учтем, что φ является линейной функцией по z, представимой в виде (2). Выпишем представление для первого слагаемого в правой части (9) без учета постоянного множителя. Заметим при этом, что второе слагаемое без учета множителя такого же типа, как и первое.

(i+;) г

e ^es cos( шt + ф ) ds

= ——e ^e z cos( a + bz + шt )+

+— e ^eL cos( a + bL + шt )+

b

⁺ в 2 e

-βz

sin( a + bz + шt ) —

-

^b -βL β 2 ^e

sin( a + bL + шt ) .

Учитывая теоремы о бесконечно малых и их пределах, получим, что при L → ∞ два слагаемых стре-

мятся к нулю и поэтому

z

(1+;)/ β_∞

^es cos( ш ₀ 1 + ф ) ds —

^- β^e

^ez cos( a + bz + ш ₀ 1 )+

b

⁺ e 2 e

^ez sin( a + bz + ш ₀ 1 ) .

Аналогично для второго слагаемого

(1+ b 2)Л β _{1    ∞}

^{e 1} s cos(2 ш ₀ 1 + ф 1) ds —

-

1 β 1 ^e

b 1

⁺ в 2

^{e 1 z} cos( a i + b i z + 2 ш 0 1 ) +

^{e 1} z sin( a 1 + b 1 z + 2 ш ₀ 1 ) .

После некоторых сокращений получим, что часть выражения (7) будет иметь вид:

левая

∞

∂T ∂t

— ^Ш 0 A 0 e

^e 0 z ( Я2 b К 2 sin( a + bz + Ш 0 1 ) — в 0 ² + b ²

-

β в 2 + b ²

cos( a + bz + ш ₀ 1 )^ +

+2 ш о A i e

β 1 z

b 1

в 2 + ь 2

sin( a 1 + b 1 z + 2 ш 0 1 ) —

-

β 1

в 2 + b 1

cos( a i + b i z + 2 ш 0 1 )^ .     (13)

Вычислим теперь интеграл в правой части выражения (7). Из свойств интеграла следует, что

Г ^КтdT^dS — K t dT Г

_∞ ∂s ∂s          ∂z ∞

.

Учитывая асимптотические свойства температуры, получаем выражение для коэффициента эффективной температуропроводности k _T :

к т ( t,z ) —

z ∂T _ds

∞ ∂t

∂T dz 1 x t,z

.

Исходя из всех предыдущих расчетов, получаем следующую аналитическую формулу для расчета k _T :

K T —

K о + K i

—YCe ^Yz + B 0 + B i

где

K 0 ^{— Ш} 0 A 0 e

^в 0 z ( № \ к 2 sin( a + bz + ш 0 1 ) — в ₀ ² + b ²

-

-r^⁰— cos( a + bz + ш 0 1 ) ) , в 0² + b ²

K i — 2 ш 0 A i e

β 1 z

b 1

в 1 ² + b ² 1

sin( a i + b i z + 2 ш 0 1 ) —

-

β 1

в 2 + b 1

cos( a 1 + b 1 z + 2 ш ₀ 1 )^ ,

B 0 — — в 0 A 0 e ^{в 0} z sin( ш 0 1 + a + bz )+

+ bA ₀ e ^{- e} 0 z cos( ш ₀ 1 + a + bz ) ,

B 1 — —в i A 1 e ^{- e 1} z sin(2 ш ₀ 1 + a 1 + b 1 z )+ + b 1 A 1 e 1 z cos(2 ш 0 1 + a 1 + b 1 z ) .

Формула (16) представляет собой аналитическую формулу для вычисления коэффициента эффективной температуропроводности, в котором отражены все виды теплообмена, происходящие внутри почвы. Заметим, что если невозможно пренебречь изменчивостью во времени, то среднесуточная температура представима как сложная функция времени и пространства, и формула (16) должна быть преобразована с учетом этого факта. Кроме того, в формуле (16) учтена линейная зависимость сдвига фаз по глубине. В отдельные моменты времени это может быть не так, учет этого требует развития предложенного подхода. Расчеты по формуле (16) для одного дня наблюдения на точке Бакчарского болота в Томской области показали вполне приемлемые результаты. Формула (16) легко может быть преобразована при использовании годовых колебаний температуры.

Выводы

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