Determination of Soil Thermophysical Properties with Various Models Across Seasons in Alkaline-Nonalkaline Soils of Igdır, Türkiye

1. Introduction

In soil formation, climate, topography, parent material and organisms outstand as effective factors. Soil salinity or alkalinity, which is caused by the effects of climate,topography and anthropogenic factors, is one of the important problems of arid-semi-arid regions and has negative effects on soil health [1,2].

Soil temperature, in the context of climate change may affect soil functioning and quality, which is one of the important factors of climate, is effective on physical (soil thermal properties, water flow etc.), chemical (mineralization, nutrients etc.) and biological (microbial activity, microbial community etc.) properties of soils [3,4,5]. Soil thermal properties are closely related to moisture content, temperature, bulk density, organic matter, salt content and porosity [5]. Some researchers have determined that in non-sodic soils, as the salt content of the soils increases, the thermal diffusivity decreases. In sodic soils sSodicity (ESP), which refers to the presence of excessive sodium ions (Na⁺) in the soil exchange complex, can significantly alter soil structure and behavior, which in turn influences how heat or water are stored and transferred through the soil[1,6,7,8].

Noborio & McInnes (1993) [8] reported that the apparent thermal conductivity of soils decreases as the salt concentration increases. Abu-Hamdeh & Reeder (2000) [9,10] also detected that an increase in the amount of added salts, at a given moisture content, decreased thermal conductivity. Morozova [11] determined that the thermal conductivity increased with the increase of KCl. Tikhonravova [12,13] revealed that the thermaldiffusivity depends not only on the salt content but also on the moisture content.

Relevant research has shown that the field studies cost money, labor, time, etc. in determining the thermal properties of the soils. This shows that laboratory measurements under controlled conditions and in shorter periods are common [14,15,16,17,18].

Although the effects of salts on thermal properties of soils have been studied, these studies were carried out under controlled conditions in the laboratory and may not reflect the field conditions, particularly for sodic or alkline soils[]`. In field conditions, it is important to calculate thermal properties with long-term temperature measurements ot the main thermal properties of soils such as volumetric heat capacity (C_v), thermal diffusivity (κ), thermal conductivity (λ), damping depth (d), thermal effusivity (e), heat wave velocity (υ), heat wave length (Λ) and heat flow (q).

The determination of these parameters is important in shaping the thermal regimes of soils. In this work, classical (amplitude, arctangent, logarithmic, and phase) methods are first used to determine the thermal conductivity coefficient [19,20,21,22,23,24,25,26,27,28,29,30,31]. Subsequently, a new "point" method was developed by our team [32,33,34,35,36,37,38], which, unlike classical approaches, also incorporates the amplitude of daily surface soil temperature (Ta)—a critical parameter that significantly influences thermal diffusivity in soils. Moreover, the proposed method accounts for the depth of the soil profile (i.e., z = L). Notably, in this approach, the thermal diffusivity parameter can be calculated at any investigation point (zi) within the soil profile. Selecting the most suitable predictive model is essential for the effective monitoring and management of soil thermal dynamics. Accurate models provide deeper insight into soil temperature fluctuations, which directly influence key processes such as irrigation scheduling, crop development, nutrient cycling, and microbial activity. Reliable predictions of soil thermal properties are particularly valuable for informed decision-making in areas such as precision agriculture, land reclamation, and environmental monitoring. Consequently, evaluating and selecting models based on their accuracy and performance under real field conditions is critical to promoting long-term sustainability and efficient resource use. [3-5, 9, 10, fm32-38].

he objectives of this study were: (i) to calculate thermal diffusivity (κ) using both classical methods and the new "point" method developed by our team, and to identify the most accurate method for the study area; (ii) to determine key thermal properties of the soil—namely volumetric heat capacity (Cv), thermal diffusivity (κ), thermal conductivity (λ), damping depth (d), thermal effusivity (e), heat wave velocity (υ), heat wavelength (Λ), and heat flow (Q)—across different seasons; and (iii) to investigate the effect of soil alkalinity on these thermal properties.

2. Mathematical formulation of the problem

The propagation of heat waves in soil represents one of the earliest applications of Fourier's theory of heat conduction [19] to the study of natural phenomena. Numerous theoretical and experimental studies have since been conducted—both in the field and under laboratory conditions—using measured soil temperature data and various heat transfer models to determine the thermal properties of soils. In most of these studies, analytical solutions were derived based on the unrealistic boundary condition T(z→∞,t)=T0T(z \rightarrow \infty, t) = T_0T(z→∞,t)=T0​, assuming a constant temperature at an infinite soil depth. Algorithms developed from these solutions have been widely used to estimate the thermal diffusivity parameter.

In contrast, the present study employs analytical solutions that are based on more realistic boundary conditions: this study uses analytical solutions obtained for the boundary conditions ∂T(z=L, t)/∂z=0 at finite (i.e.: z=L) depth and ∂T(∞, t)/∂z=0 at indefinite (i.e.: z→∞) depth, which reflects reality [35,36,37,38].

2.1. The equation of heat movement in soil by conduction

The differential equation that describes the heat movement occurring through conduction in a three-dimensional heterogeneous-isotropic porous medium is expressed as follows [23,24,26,27]:

$$\left({\rho }_s{c}_s\right){\displaystyle \frac{\partial {\rm T}}{\partial t}}={\displaystyle \frac{\partial }{\partial x}}\left({\lambda }_x{\displaystyle \frac{\partial {\rm T}}{\partial x}}\right)+{\displaystyle \frac{\partial }{\partial y}}\left({\lambda }_y{\displaystyle \frac{\partial {\rm T}}{\partial y}}\right)+{\displaystyle \frac{\partial }{\partial z}}\left({\lambda }_z{\displaystyle \frac{\partial {\rm T}}{\partial z}}\right)+r,$$

(1)

where T(x,y,z, t) - is the soil temperature (K or °C) at point z (m) at time moment t (sec); ∂T/∂x, ∂T/∂y, ∂T/∂z - is the temperature change in the direction of the ox, oy and, oz axis; ∂T/∂t - is the temperature change in unit time; ρ_s - is the soil bulk density (kg m^-3); c_s is the specific heat capacity (J kg^-1 °C^-1); λ_x λ_y λ_z – is the component of the thermal conductivity of the soil in the x, y and z directions, respectively (Wm^-1 °C^-1); r - the source of heat inside the soil is (°C).

Experiments have shown that ∂ρ_s/∂T~0, ∂C_s/∂T~0, ∂λ/∂T~0 are correct for the important soil properties: density, volumetric heat capacity and thermal conductivity parameters when the temperature ranges from -50 ≤ T ≤ + 50 ºC [26, p. 13].

The general one-dimensional heat transfer equation for a homogeneous isotropic porous medium can be written as follows:

$$\left({\rho }_s{c}_s\right){\displaystyle \frac{\partial {\rm T}}{\partial t}}={\displaystyle \frac{\partial }{\partial z}}\left({\lambda }_z{\displaystyle \frac{\partial {\rm T}}{\partial z}}\right)+r,\left(z,t\right)\in \left(0,t\right)\times \left(0<z<L\mathrm{or}\infty \right)$$

(2)

where T(z, t) - is the soil temperature (°C) at point z (m) at time moment t (sec); ∂T/∂z - is the temperature change in the direction of the oz axis; ∂T/∂t - is the temperature change in unit time; ρ_s - is the soil bulk density (kg m^-3); c_s is the specific heat capacity (J kg^-1 °C^-1); λ_z – is the component of the thermal conductivity of the soil in the z direction (Wm^-1 °C^-1); r - the source of heat inside the soil is (°C).

It is known that when modeling soil processes occurring in the soil-plant-atmosphere system, the main stages are identification (selection of equations describing these processes and boundary conditions) and realization implementation (solution of direct and inverse problems) of the model

Therefore, a systematic approach to analyzing the thermal regime of soil (depending on laboratory or field conditions) should sometimes involve simplification of the heat transfer equation and boundary conditions [36].

Essential simplification of equation (1) is possible for some problems if the coefficients of heat capacity, thermal conductivity and thermal diffusivity of soil are set constant

It is known that one-dimensional heat transfer in soil is described by the classical thermal conductivity equation, which has the form [19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41]:

$${\displaystyle \frac{\partial T}{\partial t}}=\mathsf{\kappa }{\displaystyle \frac{{\partial }^{2}T}{\partial z^2}}\left(\kappa ={\displaystyle \frac{\lambda }{C_v}}\right),\left(z,t\right)\in \left(0,t\right)\times \left(0<z<L\mathrm{or}\infty \right)$$

(3)

where T(z, t) is the soil temperature (^°C) at point z(m) at time t – time (hour, day, year); κ – thermal diffusivity of the soil (m² s^-1); λ – is the thermal conductivity of the soil (Wm^-1 °C^-1); C_v = ρ_s c_s – volumetric heat capacity of the soil (J m^{-3 °}C^-1); ρ_s - is the soil bulk density (kg m^-3); and c_s - is the specific heat capacity (J kg^-1 °C^-1); L is soil depth (m) starting from which T(z, t) = const or ∂T(L, t)/∂z = 0.

There are a number of works [19,20,21,23,24,25] in which solutions of the heat conduction equation (2.3) under different boundary conditions are considered. In this work we will also consider equation (2.3). In order to determine the change in soil temperature under the influence of various environmental factors depending on a time and depth, equation (1) must be solved analytically or numerically. For this purpose, Equation (1) must be defined and supplemented with initial and boundary conditions that account for the various factors influencing temperature variations in soils.

2.3. Identification boundary-value conditions of the soil heat conduction equation

To solve Equation (2.3) analytically, it is necessary to define initial and boundary conditions that accurately reflect heat transfer within the soil. Depending on the nature of these conditions, three types of initial-boundary value problems can be distinguished.

2.3.1. Identification ınitial condition

It should be noted that, regardless of dimensionality, only one initial condition is required for the thermal conductivity problem because Equation (2.3) is first-order in time (it involves the first derivative of temperature with respect to time). To study the influence of boundary conditions on heat transfer modeling in soil, both initial and boundary conditions must be specified. The temperature distribution of the soil profile at the initial time is given by the following equation, which serves as the initial condition:

$${\rm T}\left(z,t=0\right)=\phi \left(z\right),0\le z\le L\mathrm{or}0\le z<\infty$$

(4)

Here, φ(z) is the function expressing the distribution of temperature values in the [0, L] or [0,∞) region of the soil.

For modeling, the analytical expression of this function needs to be determined in advance. Some researchers define the function φ(z) as follows:

$${\rm T}\left(z,t=0\right)={\mathrm{T}}_0,const,0\le z\le L\mathrm{or}0\le z<\infty$$

(5)

For convenience of solving the boundary value problem, the function φ(z) is sometimes approximated by a polynomial of degree n. Then the initial condition (2.4) will have the form:

$${\rm T}\left(z,\text{}t=0\right)={\displaystyle {\sum }_{i=0}^{i=n}{a}_i{z}^i},0\le z\le L\mathrm{or}0\le z<\infty$$

(6)

Here a_n are the coefficients of the polynomial.

It is easy to see that the special case (2.6) for n = 0 corresponds to the uniform initial condition (2.5).

Since the temperature in the soil decreases with depth, the function φ(z) is sometimes approximated by a decreasing step function:

$${\rm T}\left(z,t=0\right)=\left\{\begin{array}{l}{{\rm T}}_i,z_i<z\le z_{i+1}\left(i=0,1,2,\mathrm{...},n\right)\\ {{\rm T}}_i=const,x_0=0;x_{n+1}=L\&\infty \end{array}\right.$$

(7)

It is known that the influence of the initial condition practically does not affect the soil temperature distribution at the moment of observation.

If the moment of time we are interested in is sufficiently distant from the initial one, it makes sense to neglect the initial conditions, since their influence on the process with the passage of time weakens even more, especially in the case of periodic formulation of the problem. In the periodic formulation of the problem, the initial condition does not need to be given (the so-called “tasks without an initial conditions”) [23,24,25,26,35,36].

2.3.2. Identification of boundary conditions

The condition known as Boundary Condition 1, which expresses the temperature change at the soil surface, is written as a function of time in the following form:

$$T\left(z=0,t\right)=\phi (t)$$

(8)

where φ(t) is a function of the soil surface temperature; the analytical form of the function must be determined in advance.

Various empirical functions can be determined by the least squares method (LSM) to analytically express φ(t) based on daily, monthly, or annual soil surface temperature measurements.

As is known, the temperature of the earth has a distinct periodic character, both daily and annual. Sometimes there may be deviations from this periodic distribution.

For example, if the weather is cloudy, the periodicity of the measurement values of the thermal sensors used to determine the temperature values will be disturbed.

In these and similar meteorological (extreme wind, rain, etc.) conditions, periodicity will be disturbed in the distribution of temperature values measured by sensors on the soil surface at some moments of time.

To overcome this deficiency, it is appropriate to use a trigonometric polynomial to determine the analytical expression of the distribution of soil surface temperature values. That is, the first type of boundary condition on the soil surface will be as follows ​​[23,24]:

$$T(z=0,t)=T_0+{\displaystyle {\sum }_{j=1}^m{T}_j\cdot \cos \left(j\cdot \mathsf{\omega }t+{\epsilon }_{\mathrm{j}}\right)},t\ge 0$$

(9)

where T₀ is the average daily (or annual) temperature of the active soil surface (^oC); T_j is the amplitude of the wave at the surface level for the j^th harmonic (^oC or K); j is the index of the harmonic in the series; m is the harmonic number; ω = 2π/τ₀ is the angular daily (or annual) frequency; τ_o is the temperature wave period (days or years), for τ₀=24 hours: ω = 7,27∙10^-5 (rad/s); t – time; ε_j is the phase angle of the wave at the surface level for the j^th harmonic (radians).

Sometimes condition (2.9) can be written in a more practical form:

$$T(z=0,t)=T_0+{\displaystyle {\sum }_{j=1}^m\left[A_j\cos \left(j\cdot \omega t\right)+B_j\cos \left(j\cdot \omega t\right)\right]},t\ge 0$$

(10)

Here, the parameters T₀, T_a, A_j and B_j are calculated by the following formulas:

$$T_0={\displaystyle \frac{1}{N}}{\displaystyle {\sum }_{i=1}^{N}T\left(0,t_i\right)},A_j={\displaystyle \frac{2}{N}}{\displaystyle {\sum }_{i=1}^{N}T\left(0,t_i\right)\cos \left(j\omega t_i\right)},B_j={\displaystyle \frac{2}{N}}{\displaystyle {\sum }_{i=1}^{N}T\left(0,t_i\right)\sin \left(j\omega t_i\right)}$$

(11)

Where; T(0, t_i) - is the temperature of the soil surface at time ti (^oC); ti - is the time of observations, N - is the number of observations. The phase shifts (phase angle) ε_j and amplitude terms (T_j) are calculated with the coefficients Aj and B_j as follows:

$$T_j=\sqrt{A_j^2+B_j^2},{\epsilon }_j=\left\{\begin{array}{l}\arctan \left(-B_j/A_j\right),A_j>0\\ \pi /2,A_j=0\\ \pi -\arctan \left(B_j/A_j\right),A_j<0\end{array}\right.,\left(j=1,2,\mathrm{...},m\right)$$

(12)

The boundary conditions can be set at any depth of the profile, similar to the boundary condition at the soil surface.

The oscillation of the soil temperature wave becomes increasingly relaxed as depth increases, and after a certain depth z > L the temperature does not change.

Therefore, the boundary conditions are written in two forms that can express the situation where there is no heat source and transfer at the depth of the soil:

In most studies [19,20,21,22,23,24,25,26,27,28,29,30,31], the soil is usually considered as a semi-infinite region, and the lower boundary condition is set by considering the fact that the soil temperature at great depth is constant and equal to T₀. In this case, a boundary condition of the 1st kind is set:

$$T\left(z\to \infty ,t\right)=T_0\mathrm{or}T\left(z\to \infty ,t\right)=T_{\infty }=const$$

(13)

It should be especially emphasized that at infinity the temperature of the Earth: T(z→∞,t) cannot be known (since the radius of the Earth's sphere R≠∞ is not infinite, since R < 6.371 km).

Therefore, it is important to emphasize that in practical calculations it is not reasonable or possible to give the boundary condition at z→∞.

Moreover, the study shows that temperature fluctuations in the soil attenuate already at a depth of 0.5—1 m of soil [23].

Consequently, at a certain soil depth, one of the following conditions, known as boundary condition 2.type, more plausibly reflects reality:

$$\partial T\left(z=L,t\right)/\partial z=0$$

(14)

or

$$\partial T\left(z\to \infty ,t\right)/\partial z=0$$

(15)

2.4. Solution of the direct problem of thermal conductivity in soil

Given the initial and boundary conditions, different solutions can be proposed for the direct problem of heat transfer, i.e., the determination of the temperature dynamics at a specific depth.

In particular, of great practical interest are the problems of thermal conductivity in soils with periodically varying temper atures on the surface of soil.

Problem-1. First, let us consider the heat transfer problem for a semi-infinite homogeneous soil region, i.e. [0 ≤ z < ∞):

$${\displaystyle \frac{\partial T}{\partial t}}=\kappa {\displaystyle \frac{{\partial }^{2}T}{\partial z^2}}\left(\kappa =\lambda /C_v\right),\left(z,t\right)\in \left(0,t\right)\times \left(0<z<\infty \right),$$

(16)

$$T(z,0)=T_0$$

(17)

$$T\left(z=0,t\right)=T_0+{\displaystyle {\sum }_{j=1}^{j=m}{T}_j\cdot \cos \left(j\omega t+{\epsilon }_j\right)},\left(t\ge 0\right)$$

(18)

$$\partial T\left(z\to \infty ,t\right)/\partial z=0,\left(t\ge 0\right)$$

(19)

The solution of equation (2.16), with the initial (2.17) and boundary conditions (2.18)-(2.19), in dimensionless variables and parameters is as follows [33]:

$$T\left(y,\tau \right)=T_0+{\displaystyle {\sum }_{j=1}^{j=m}{\Phi }_j\left(b_j,y\right)\cdot \cos \left[j\overline{\omega }\tau +{\alpha }_j\left(b_j,y\right)\right]}$$

(20)

where, $y={\displaystyle \frac{z}{L}},\tau ={\displaystyle \frac{\kappa }{L^2}}t,b_j=\sqrt{j{\displaystyle \frac{\overline{\omega }}{2}}}=L\sqrt{j{\displaystyle \frac{\pi }{{\tau }_0\kappa }}},\overline{\omega }=\omega {\displaystyle \frac{L^2}{\kappa }}$

$${\Phi }_j\left(b_1,y\right)=T_j{e}^{-b_{j}y},{\alpha }_j\left(b_j,y\right)={\epsilon }_j-b_{j}y$$

(21)

Ф_j (b_j, у) is the amplitude of temperature fluctuations at the dimensionless depth y; α_j (b_j, у) is the phase angle of the wave at the surface level for the j^th harmonic at dimensionless depth y (radians).

This problem was studied by Fourier; it was first used by Kelvin to determine the course of temperature in the soil of Edinburgh [23]:

As correctly noted in [33,35,38], when performing practical calculations, it is impossible to set as input data the values of soil temperature at infinity, since they are unknown. And it is impossible to measure it. Therefore, instead of T(z→∞, t) it is necessary to set the temperature at some depth L, starting from which at z > L the value T(z, t) = const or ∂T(L, t) / ∂z = 0. Thus, condition (2.15) (∂T(L, t) / ∂z = 0) is more consistent with the real conditions than boundary condition (2.14) (T(z→∞, t) = T₀=const) [34]. Therefore, the following boundary value problem for equation (2.16) should also be considered.

Problem-2. Next, consider the boundary value problem for the thermal conductivity equation in a finite [0 ≤z ≤ L] homogeneous soil region in the following form:

$${\displaystyle \frac{\partial T}{\partial t}}=\kappa {\displaystyle \frac{{\partial }^{2}T}{\partial z^2}},\left(z,t\right)\in \left(0,t\right)\times \left(0<z<L\right)$$

(16)

$$T\left(0,t\right)=T_0+{\displaystyle {\sum }_{j=1}^{j=m}{T}_j\cdot \cos \left(j\omega t+{\epsilon }_j\right)},\left(t\ge 0\right)$$

(18)

$$\partial T\left(z=L,t\right)/\partial z=0,\left(t\ge 0\right).$$

(22)

The solution of the problem (2.16), (2.18) и (2.22) with dimensionless variables and parameters is as follows [33,34,36,38]:

$${\rm T}\left(y,\tau \right)={{\rm T}}_0+{\displaystyle {\sum }_{j=1}^{j=m}{\Phi }_j\left(b_j,y\right)\cdot \cos \left[j\overline{\omega }\tau +{\alpha }_j\left(b_j,y\right)\right]},$$

(23)

Where $y=z/L,\tau =\kappa t/L^2$ and

$$\begin{array}{l}{\Phi }_j\left(b_j,y\right)=T_j\cdot \sqrt{{\displaystyle \frac{\mathrm{ch}\left(d_j\right)+\cos \left(d_j\right)}{\mathrm{ch}\left(2b_j\right)+\cos \left(2b_j\right)}}},d_j=2b_j\left(1-y\right),b_j=L\sqrt{j{\displaystyle \frac{\mathsf{\omega }}{2\kappa }}}\\ {\alpha }_j\left(b_j,y\right)={\epsilon }_j-\arctan \left[{\displaystyle \frac{\mathrm{sh}\left(q_j\right)\sin \left(b_{j}y\right)+\mathrm{sh}\left(b_{j}y\right)\sin \left(q_j\right)}{\mathrm{ch}\left(q_j\right)\cos \left(b_{j}y\right)+\mathrm{ch}\left(b_{j}y\right)\cos \left(q_j\right)}}\right],q_j=b_j\left(2-y\right)\end{array}$$

(24)

$\mathrm{ch}\left(z\right)=\left(e^z+e^{-z}\right)/2,\mathrm{sh}\left(z\right)=\left(e^z-e^{-z}\right)/2$- are the hyperbolic cosine and sine, respectively.

In the analytical methods developed for the determination of the parameters of the heat and mass transport models, it has been shown that the data at any depth z = z_i of the soil profile are determined with statistically greater variation than the average data in the [0, L] layer of the soil. In other words, the average temperature of a certain soil layer (e.g. 0-20 or 0-40 cm layer), varies less than the temperature at a certain depth z = z_i [39,40].

Therefore, the average values of soil layer temperature should be used to determine the thermal diffusivity coefficient from the data of field and laboratory experiments.

The use of average temperature values ​​in determining the heat emission parameter is called the average integral method.

For this purpose, to determine the average temperature in the layer 0 ≤ z ≤ L, let us calculate the definite integral of the dimensionless solution (2.20) of equation (2.16) in the interval 0 ≤ y ≤ 1 with respect to the variable y.

$$\overline{{\rm T}}\left(b_1,\mathrm{...},b_j,t\right)=T_0+{\displaystyle \sum _{j=1}^m{{\rm M}}_j\left(b_j\right)\cdot \cos \left[j\omega t+{\alpha }_j\left(b_j\right)\right]},b_j=L\sqrt{j\pi /{\tau }_0\kappa }$$

(25)

$${{\rm M}}_j\left(b_j\right)=T_j\sqrt{{\displaystyle \frac{\cosh \left(b_j\right)-\cos \left(b_j\right)}{b_j^2{e}^{b_j}}}},{\alpha }_j\left(b_j\right)={\epsilon }_j-\arctan \left\{{\displaystyle \frac{1-e^{-b_j}\left(\sin b_j+\cos b_j\right)}{1+e^{-b_j}\left(\sin b_j-\cos b_j\right)}}\right\}$$

(26)

Similarly, for the solution (2.23) with arbitrary harmonic m, the average temperature in the layer 0 ≤ z ≤ L can be determined as follows:

$$\overline{{\rm T}}\left(b_1,\mathrm{...},b_j,t\right)={{\rm T}}_0+{\displaystyle \sum _{j=1}^m{\Phi }_j\left(b_j\right)\cdot \cos \left[j\omega t+{\gamma }_j\left(b_j\right)\right]}$$

(27)

$$\begin{array}{l}{\Phi }_j\left(b_j\right)=T_j\left[{\displaystyle \frac{\sqrt{{\sinh }^2\left(2b_j\right)+{\sin }^2\left(2b_j\right)}}{\sqrt{2}{b}_j\Delta (b_j)}}\right],\Delta (b_j)=\cosh \left(2b_j\right)+\cos \left(2b_j\right)\\ {\gamma }_j\left(b_j\right)={\epsilon }_j-\arctan \left[{\displaystyle \frac{\sinh \left(2b_j\right)-\sin \left(2b_j\right)}{\sinh \left(2b_j\right)+\sin \left(2b_j\right)}}\right],b_j=L\sqrt{j{\displaystyle \frac{\pi }{{\tau }_0\kappa }}}.\end{array}$$

(28)

These solutions are widely used in soil science for two purposes: estimating temperature values in the soil profile and calculating thermal conductivity (k).

2.5. Solution of the inverse problem of thermal conductivity in soil

Most methods for determining the thermal diffusivity of soil are based on solving inverse problems of thermal conductivity, which are obtained under given boundary conditions. Depending on the boundary conditions, different methods are obtained.

A more detailed review of studies on the determination of the diffusion parameter in soils is given in many works [19,20,21,22,23,24,25,27,30,31,32,33,34,35,36,41,42,43,44,45,46,47,48,49].

This paper details methods for determining the thermal diffusivity using analytical solutions (2.20) and (2.23) found by second-order boundary conditions ∂T(z→∞, t)/∂z=0 and ∂T(z=L, t)/∂z=0.

Depending on the use of measured temperature values of the soil profile, these methods are divided into 3 classes: layered, point, and numerical (harmonic) methods.

These methods are based on using 1) temperature values at two different depths of the soil profile, e.g. z=z_i and z=z_i+1, 2) temperature values at a certain depth z=zi and 3) the value of the average temperature in the [0, L] layer of the soil.

2.5.1. a. layered methods

These methods are based on temperature measurements at two different soil profile depths z=z_i and z=z_i+1 and are called layerd methods. In addition, these layered (classical) methods will be denoted as M1, M2, M3, M4 algorithms.

2.5.1. b. Point methods

These methods are based on the use of measured temperature values at a certain depth z = z_i of the soil profile and the amplitude (Ta) of the soil surface temperature wave and are called point methods [32,33,34,35,36,38]. In addition, these point (proposed) methods will be denoted as algorithms M5, M6, M7, M8.

2.5.1. c. Numerical (Harmonic) methods

The value of heat diffusion (κ) is selected iteratively in a way that minimizes the sum of the differences between the measured temperature values ​​T_mes(z_i, t_i) at the desired depth zi and the temperature values ​​T_cal (z_i, t_i) estimated according to solutions (2.20) or (2.23) [27,28,29,30,33]:

$$\underset{i,j}{\min }{\displaystyle {\sum }_{j=1}^{j=r}{\left[T_{mes}\left(z_i,t_j\right)-T_{cal}\left(z_i,t_j\right)\right]}^2}$$

(29)

For example, the iteration can be continued until the condition is met. The sum of squares will be composed of 24 whole (r=0,1,2,…,24) or half-hour (r=0.5,1.5,…,23.5) differences between measured and calculated temperature.

The existing and proposed methods for calculating the thermal diffusivity parameter, developed respectively for the boundary conditions ∂T(∞, t)/∂z=0 and ∂T(z=L, t)/∂z=0, are given below.

The existing (Classical) methods for calculating the thermal diffusivity parameter derived from the boundary condition ∂T(z→∞, t)/∂z=0 are as follows

Classical methods. Below we present classical algorithms based on the solution of inverse problems of the heat transfer equation obtained for a semi-infinite homogeneous soil region under the boundary conditions ∂T(∞, t)/∂z=0 when the diurnal temperature variation at the soil surface is represented by one (M1, M4) and two (M2, M3) harmonics;

М1. Amplitude Method. [19,21,23]:

$${\mathsf{\kappa }}_i={\displaystyle \frac{\pi }{{\tau }_0}}\cdot {\displaystyle \frac{{\left(z_{i+1}-z_i\right)}^2}{{\ln }^2\left[A\left(z_{i+1}\right)/A\left(z_i\right)\right]}},A\left(z_i\right)={\displaystyle \frac{{{\rm T}}_{\max }\left(z_i\right)-{{\rm T}}_{\min }\left(z_i\right)}{2}}$$

(30)

M2. Arctangent Method [22]:

$${\kappa }_i={\displaystyle \frac{\pi }{{\tau }_0}}{\displaystyle \frac{{\left(z_{i+1}-z_i\right)}^2}{{\arctan }^2\frac{\left[T_1(z_i)-T_3(z_i)\right]\left[T_2(z_{i+1})-T_4(z_{i+1})\right]-\left[T_2(z_i)-T_4(z_i)\right]\left[T_1(z_{i+1})-T_3(z_{i+1})\right]}{\left[T_1(z_i)-T_3(z_i)\right]\left[T_1(z_{i+1})-T_3(z_{i+1})\right]+\left[T_2(z_i)-T_4(z_i)\right]\left[T_2(z_{i+1})-T_4(z_{i+1})\right]}}}$$

(31)

M3. Logarithmic Method [20]:

$${\kappa }_i={\displaystyle \frac{4\pi {\left(z_{i+1}-z_i\right)}^2}{{\tau }_0}}{\ln }^{-2}\left\{{\displaystyle \frac{{\left[T_1(z_i)-T_3(z_i)\right]}^2+{\left[T_2(z_i)-T_4(z_i)\right]}^2}{{\left[T_1(z_{i+1})-T_3(z_{i+1})\right]}^2+{\left[T_2(z_{i+1})-T_4(z_{i+1})\right]}^2}}\right\}$$

(32)

M4. Phase Method [24]:

$${\kappa }_i={\displaystyle \frac{\pi }{{\tau }_0}}{\left({\displaystyle \frac{z_{i+1}-z_i}{{\epsilon }_{i+1}-{\epsilon }_i}}\right)}^2={\displaystyle \frac{\pi }{{\tau }_0}}{\left({\displaystyle \frac{\Delta z_i}{\Delta {\epsilon }_i}}\right)}^2,\left\{\begin{array}{l}\Delta z_i=z_{i+1}-z_i\\ \Delta {\epsilon }_i={\epsilon }_{i+1}-{\epsilon }_i\end{array}\right.$$

(33)

where τ₀ in algorithms (2.30)– (2.33) is the period of the heat wave (for example, 24 hours for daily observations); T_min (z) and T_max (z) are minimum and maximum tem-perature values during measurements at depths z = z_i and z = z_i+1, respectively; T_i(z) and T_i+1(z) in algorithms (2.31) and (2.32) are the values of soil temperature at depths z = z_i and z = z_i+1, respectively, at time t_i = iτ₀/4 (i=1,2,3,4) (for this example, τ₀ = 24 h and t₁ = 6, t₂ = 12, t₃ = 18 and t₄ = 24 h); εi is the initial phase of the soil temperature at depth z_i in algorithm (2.33).

Improved methods. In contrast to the above algorithms, we continued to develop new methods that take into account the real process of heat exchange in the soil. For this purpose, we used the solution of the heat transfer equation, under the boundary conditions ∂T(z→∞, t)/∂z=0 and ∂T(z=L, t)/∂z=0, obtained for limited and semi-limited soil thickness. These methods are based on the use of measured temperature values at a certain depth z = z_i of the soil profile and the amplitude (T_a) of the soil surface temperature wave and are called point methods [32,33,34,35,36,37,38].

Tthese point (proposed) methods will be denoted as algorithms M5, M6, M7, M8. The proposed methods for calculating the thermal diffusivity of soil are presented below.

The following methods (M5 and M6) are developed on the basis of a boundary condition of the second kind in the infinity z→∞, i.e. ∂T(z→∞, t)/∂z=0, when the diurnal temperature fluctuations at the soil surface is represented by one (m=1) (M5) and two (m=2) (M6) harmonics.

These methods, unlike the classical ones (М1-М4), include the amplitude of fluctuations (T_a) of the soil surface temperature

M5. Proposed Point method.

Algorithm including the amplitude of fluctuations (T_a) of the soil surface temperature and the logarithm (m=1)[36]:

$${\kappa }_i={\displaystyle \frac{\pi }{{\tau }_0}}\cdot {\displaystyle \frac{{\left(2z_i\right)}^2}{{\ln }^2\left\{{\displaystyle {\sum }_{j=1}^2{\left[T\left(z_i,t_j\right)-T\left(z_i,t_{j+2}\right)\right]}^2}/4{{\rm T}}_a^2\right\}}}$$

(34)

where; T(z_i, t_j) in the algorithm (2.34) are the soil temperature values at the depth z=z_i, at the time at the time moment t_j=j·τ₀ / 4 (j=1,2,3,4) (for τ₀ = 24 h and t₁=6, t₂=12, t₃=18 and t₄=24 h).

M6. Proposed Point method.

Algorithm including the amplitude of fluctuations (T_a) of the soil surface temperature and the logarithm (m=2) [36,38]:

$${\kappa }_i={\displaystyle \frac{\pi }{{\tau }_0}}\cdot {\displaystyle \frac{{\left(2z_i\right)}^2}{{\ln }^2\left\{{\displaystyle {\sum }_{j=1}^4{\left[T\left(z_i,t_j\right)-T\left(z_i,t_{j+4}\right)\right]}^2}/8{{\rm T}}_a^2\right\}}}$$

(35)

where; T(z_i, t_j) in the algorithm (2.35) are the soil temperature values at the depth z=z_i, at the time at the time moment t_j=j·τ₀ / 8 (j=1,2,…,8) (for τ₀ = 24 h and t₁=3, t₂=6, t₃=9,…, t₈=24 h).

The following methods (M7 and M8) are developed on the basis of the boundary condition of the second kind at a finite depth z=L, i.e. at ∂T(z=L, t)/∂z=0, when the daily temperature fluctuation at the soil surface is represented by one (m=1) (M7) and two (m=2) (M8) harmonics.

These methods, in addition to the amplitude of fluctuations (T_a) of the soil surface temperature, also take into account the depth of the soil profile (L).

M7. Proposed Point method.

Algorithm including the amplitude of fluctuations (T_a) of the soil surface temperature and also take into account the depth of the soil profile (L) (m=1) [36]:

$$M_{mes}^{m=1}\left(b;y_i\right)={\displaystyle \frac{{\displaystyle \sum _{j=1}^2{\left[T_j\left(y_i,t_j\right)-T_{j+2}\left(y_i,t_{j+2}\right)\right]}^2}}{4{{\rm T}}_a^2}}={\displaystyle \frac{\cosh \left[2b\left(1-y_i\right)\right]+\cos \left[2b\left(1-y_i\right)\right]}{\cosh \left(2b\right)+\cos \left(2b\right)}}=M_{cal}^{m=1}\left(b;y_i\right)$$

(36)

where; T(z_i, t_j) in the algorithm (2.36) are the soil temperature values at the depth z=z_i, at the time moment t_j=j·τ₀ / 4 (j=1,2,3,4) (for τ₀ = 24 h and t₁=6, t₂=12, t₃=18 and t₄=24 h).

Using the amplitude of fluctuations (T_a) of the soil surface temperature and the temperature values T(z_i,t_j) measured at different moments of time (t_i = 1, 2, 3 and 4) at a given depth z = z_i (in the dimensionless notation 0 ≤ y_i ≤ 1) of the soil layer [0, L], we first find the value of the left side of equation (2.36), i.e. $M_{mes}^{m=1}(b;y_i)$. Then, the values ​​of the expression $M_{cal}^{m=1}\left(b;y_i\right)$ on the right-hand side of equation (2.36) are calculated using different values ​​of b_i > 0 for each y_i=z_i/L.

Then the value of b* satisfying the condition $\left|M_{mes}^{m=1}(b^{\ast };y_i)-M_{cal}^{m=1}(b^{\ast };y_i)\right|\le {10}^{-4}$ is determined. Finally, using the found value of b*, the thermal diffusivity (κ_i) of soil at a given depth is calculated z=z_i according to the following equation:

$$b_i^{\mathrm{*}}=L\sqrt{\pi /{\tau }_0{\kappa }_i}\Rightarrow {\mathsf{\kappa }}_i=\left(\pi /{\tau }_0\right)\cdot {\left(L/b_i^{\mathrm{*}}\right)}^2$$

(37)

M8. Proposed Point method.

Algorithm including the amplitude of fluctuations (T_a) of the soil surface temperature and also take into account the depth of the soil profile (L) (m=2) [32,33,34]:

$$M_{mes}^{m=2}\left(b;z_i\right)={\displaystyle \frac{{\displaystyle \sum _{i=1}^4{\left[T\left(z_i,t_i\right)-T\left(z_i,t_{i+4}\right)\right]}^2}}{8{Т}_a^2}}={\displaystyle \frac{\cosh \left[2b\left(1-y_i\right)\right]+\cos \left[2b\left(1-y_i\right)\right]}{\cosh \left(2b\right)+\cos \left(2b\right)}}=M_{cal}^{m=2}\left(b;z_i\right)$$

(38)

This method differs from M7 in that the number of harmonics on the soil surface is m=2. Also, after calculating the values of the left M_mes^m=2 (b) and right M_cal^m=2(b) sides of equation (2.38), the value of b*_i is determined which satisfies the condition $\left|M_{mes}^{m=2}\left(b^{\ast };z_i\right)-M_{cal}^{m=2}\left(b^{\ast };z_i\right)\right|\le {10}^{-4}$.

Finally, using the b* value, the value of the thermal diffusivity (κ_i) is calculated using the equation (2.37).

2.5.2. Calculation of Thermal Properties of Soil

2.5.2. a. Volumetric heat capacity

The volumetric heat capacity of the soil was calculated using the generally accepted formula [43,51,52]:

$$C_v=C_{m,s}\cdot {\rho }_b+C_{v,w}\cdot \theta =\left[C_{m.org}\cdot {\displaystyle \frac{m_{org}}{m}}+C_{m.\min }\cdot \left(1-{\displaystyle \frac{m_{org}}{m}}\right)\right]\cdot {\rho }_b+C_{v,w}\cdot {\rho }_w\cdot \theta$$

(39)

where, the C_m,s is the specific heat of the soil's solid part, J/(kg·^oC); ρ_b is the soil bulk density, kg/m³; C_v,w is the volumetric heat capacity of the soil moisture equal to 4186.6 kJ/( m^{3 o}C); C_w is the specific heat of water, J/(kg ^oC); ρ_w is the water density, kg/m³; θ is the volumetric moisture content (m³/m³); C_{m, org} and C_{m, min} are the specific heats of the organic and mineral components of the soil solid phase respectively, J/(kg ^oC); m_org is the mass of soil organic matter, kg; m is the soil mass (kg); m_org /m is the content of organic substance in soil, %.

2.5.2. b. Calculation of other thermal properties

Other thermophysical parameters such as thermal conductivity (λ), damping depth (d), thermal effusion (e), heat wave velocity (υ) and heat wave length (Λ), will be calculated using the following formulas, respectively [23,24,27,53,54]:

$$\lambda =\kappa C_v,d=\sqrt{{\displaystyle \frac{{\tau }_0\kappa }{\pi }}},e=C_v\sqrt{\mathsf{\kappa }},\vartheta =2\sqrt{{\displaystyle \frac{\pi \kappa }{{\tau }_0}}},\Lambda =2\sqrt{\pi {\tau }_0\kappa }=2\pi d$$

(40)

2.5.2. c. Calculation of heat flux in soil

Calculation of heat flux in soil is based on the use of values of known thermal properties. If these properties (volumetric heat capacity - Cv, diffusivity - κ and average daily soil surface temperature - T₀, amplitudes - T_i and phase shifts - ε_i) are determined, it is possible to calculate the amount of heat and heat flux (J) passing through the soil surface and causing a certain temperature change.

The soil heat flux in a soil (J) at depth z and time t is given by the Fourier law of heat conduction can be expressed as [19,21,23,24,50]:

$$J\left(z,t\right)=-\mathsf{\lambda }{\displaystyle \frac{\partial T\left(z,t\right)}{\partial z}}$$

(41)

where J is the soil heat flux (Wm⁻²); λ is the thermal conductivity of the soil (Wm⁻¹ K⁻¹); T is the soil temperature (^oC or K) and z is the depth from the soil surface (m); ∂T/∂z is the gradient of temperature. When referring to J at the soil surface (z=0), it is denoted G.

Existing methods for calculating heat flow based on the found solution for the boundary condition ∂T(z→∞, t)/∂z=0.

The heat flux into the soil q (z, t) is calculated by substituting the solution of equation (2.3) into equation (2.41). To calculate the heat flux from the soil surface to the depth, we first use the general solution (2.20) obtained for the boundary condition of the second kind, i.e. ∂T(z→∞, t)/∂z=0 in equation (2.41). It can be easily shown that the heat flux into the soil at any depth z=h and time t is of the form [23,24]:

$$J\left(z,t\right)={\displaystyle {\sum }_{j=1}^m{T}_{\mathrm{j}}{C}_v\sqrt{j\mathsf{\omega }\kappa }\cdot \exp \left(-h\sqrt{j\frac{\mathsf{\omega }}{2k}}\right)\cdot \cos \left[\frac{\pi }{4}+\left(j\mathsf{\omega }\mathrm{t}+{\mathsf{\epsilon }}_j-h\sqrt{j\frac{\mathsf{\omega }}{2k}}\right)\right]}$$

(42)

Proposed methods for calculating heat flow based on the found solution for the boundary condition ∂T(z=L, t)/∂z=0.

Similarly, using the general solution (2.23) obtained for the boundary condition of the second kind, i.e., ∂T(z=L, t)/∂z=0 into equation (2.41), we have [32,33]:

$$\begin{array}{l}J\left(z,t\right)=-\mathsf{\lambda }{\left.{\displaystyle \frac{\partial T\left(z,t\right)}{\partial z}}\right|}_{z=h}=\\ =-{\displaystyle \frac{\lambda }{L}}{\displaystyle \sum _{j=1}^m\frac{T_j{b}_j}{\Delta \left(b_j\right)}\left\{\left[{\Pi }_{1j}\left(b_j,y_h\right)\cos \left(j\mathsf{\omega }t+{\mathsf{\epsilon }}_j\right)+{\Pi }_{2j}\left(b_j,y_h\right)\sin \left(j\mathsf{\omega }t+{\mathsf{\epsilon }}_j\right)\right]\right\}}\end{array}$$

(43)

where,

$$\begin{array}{l}{\Pi }_{1j}\left(b_j,y_h\right)=\mathrm{sh}\left(b_j{y}_h\right)\cos \left[b_j\left(2-y_h\right)\right]-\mathrm{sh}\left[b_j\left(2-y_h\right)\right]\cos \left(b_j{y}_h\right)+\\ +\mathrm{ch}\left(b{y}_h\right)\sin \left[b_j\left(2-y_h\right)\right]-\mathrm{ch}\left[b_j\left(2-y_h\right)\right]\sin \left(b{y}_h\right)\\ {\Pi }_{2j}\left(b_j,y_h\right)=\mathrm{sh}\left[b_j\left(2-y_h\right)\right]\cos \left(b_j{y}_h\right)-\mathrm{ch}\left[b_j\left(2-y_h\right)\right]\cdot \sin (b_j{y}_h)+\\ +\mathrm{ch}\left(b_j{y}_h\right)\cdot \sin \left[b_j\left(2-y_h\right)\right]-\mathrm{sh}\left(b_j{y}_h\right)\cos \left[b_j\left(2-y_h\right)\right]\\ \Delta \left(b_j\right)=\mathrm{ch}\left(2b_j\right)+\cos \left(2b_j\right),b_j=L\sqrt{j{\displaystyle \frac{\mathsf{\omega }}{2\kappa }}},y_h={\displaystyle \frac{h}{L}},0\le h\le L\end{array}$$

(44)

2.5.3. Comparison of Methods, Models assessment

To determine the values of the parameters of empirical models (type (2.10)), as well as to compare the measured T_mes(z_i,t_j) and calculated T_cal(z_i,t_j) temperature values of the studied soils using formulas (2.20) and (2.23) for one and two harmonics, various criteria are used, which are given below [54,55,56]:

1. Pearson’s Correlation Coefficient (r):

$$r={\displaystyle \frac{{\displaystyle {\sum }_{i=1}^{i=n}\left(T_i-{\overline{T}}_i\right)\left({\tilde{T}}_i-{\overline{\tilde{T}}}_i\right)}}{\sqrt{{\displaystyle {\sum }_{i=1}^{i=n}{\left(T_i-{\overline{T}}_i\right)}^2}\cdot {\displaystyle {\sum }_{i=1}^{i=n}{\left({\tilde{T}}_i-\overline{\tilde{T}}\right)}^2}}}}$$

(45)

2. Nash–Sutcliffe efficiency coefficient (NSE) [56].

$$NSE=1-{\displaystyle \frac{{\displaystyle {\sum }_{i=1}^n{\left(T_i-{\tilde{T}}_i\right)}^2}}{{\displaystyle {\sum }_{i=1}^n{\left(T_i-\overline{T}\right)}^2}}}$$

(46)

This formula can be applied for linear regression and original data on any model.

3. Adjusted coefficient of determination or adjusted R-squared (R²_adj)

$$R_{adj}^2=1-\left(1-R^2\right)\cdot {\displaystyle \frac{n-1}{n-p}}$$

(47)

4. Root mean squared error (RMSE):

$$R{\rm M}S{\rm E}={\sigma }_{{\rm T}/t}=\sqrt{{\displaystyle \frac{1}{q}}{\displaystyle {\sum }_{i=1}^{i=n}{\left({{\rm T}}_i-{\widetilde{{\rm T}}}_i\right)}^2}},q=\left\{\begin{array}{l}n-p-1,n\le 30\\ n-p,n>30\end{array}\right.$$

(48)

5. Normalized root mean square error (NRMSE),

$$NRMSE={\displaystyle \frac{1}{{{\rm T}}_{i,\max }-{{\rm T}}_{i,\min }}}\sqrt{{\displaystyle \frac{1}{n}}{\displaystyle {\sum }_{i=1}^n{\left({{\rm T}}_i-{\widetilde{{\rm T}}}_i\right)}^2}}$$

(49)

6. Mean absolute percentage error (MAPE):

$$MAPE=\%100{\displaystyle \frac{1}{n}}{\displaystyle {\sum }_{i=1}^n\left|\frac{T_i-{\tilde{T}}_i}{T_i}\right|}$$

(50)

7. Mean absolute error (MAE),

$$MAE={\displaystyle \frac{1}{n}}{\displaystyle {\sum }_{i=1}^n\left|T_i-{\tilde{T}}_i\right|}$$

(51)

8. The bias is used as a measure for a systematic underlying mismatch between the observed and predicted temperature values of soils

$$bias={\displaystyle \frac{1}{n}}{\displaystyle {\sum }_{i=1}^n\left(T_i-{\tilde{T}}_i\right)}$$

(52)

9. Agreement index or Willmott’s index of agreement (D):

$$D=1-{\displaystyle \frac{{\displaystyle {\sum }_{i=1}^n{\left(T_i-{\tilde{T}}_i\right)}^2}}{{\displaystyle {\sum }_{i=1}^n{\left(\left|T_i-{\overline{T}}_i\right|+\left|{\tilde{T}}_i-{\overline{T}}_i\right|\right)}^2}}}$$

(53)

10. Theil's forecast accuracy coefficient:

$$UI={\displaystyle \frac{\sqrt{ESS}}{\sqrt{{\displaystyle {\sum }_{i=1}^n{T}_i^2}}+\sqrt{{\displaystyle {\sum }_{i=1}^n{\tilde{T}}_i^2}}}},UII={\displaystyle \frac{\sqrt{ESS}}{\sqrt{{\displaystyle {\sum }_{i=1}^n{T}_i^2}}}},UIII={\displaystyle \frac{\sqrt{ESS}}{\sqrt{{\displaystyle {\sum }_{i=1}^n{\tilde{T}}_i^2}}}}$$

(54)

11. Akaike information criterion:

$$AIC=\left\{\begin{array}{l}\ln \left({\displaystyle \frac{ESS}{n}}\right)+{\displaystyle \frac{2p}{n}},\left(n/p\ge 40\right)\\ \ln \left({\displaystyle \frac{ESS}{n}}\right)+{\displaystyle \frac{2p(p+1)}{n-\left(p+1\right)}},\left(n/p<40\right)\end{array}\right.$$

(55)

In equations (2.46)-(2.53), the following notation is adopted: T_i represents the observed values of the dependent variable; ${\tilde{T}}_i$ represents the estimated values of the dependent variable; $\overline{T}={\displaystyle {\sum }_{i=1}^{i=n}{T}_i}/n-$ represents the average of the observed values; n represents the number of data points; p represents the number of estimable parameters in the approximating model, including the intercept term, where p<n; R²represents the coefficient of determination; ESS represents the sum of squared errors and is given by $ESS={\displaystyle {\sum }_{i=1}^{i=n}{\left(T_i-{\tilde{T}}_i\right)}^2}$.

In our study, the computation of model parameters is carried out using the implementation of the Levenberg-Marquardt method in the STATISTICA-7 software package.

3. Materials and Methods

3.1. Study area

The experiment was established at Iğdır University Agricultural Research and Application Center. The climate is hot in summer and mild in winter in the region. According to Koppen-Geiger climate classification, the climate is Bsh (Hot semi-arid climates).

The mean, maximum, and minimum temperature values from 1941 to 2021 were 12.2, 42, and -30.3 °C, respectively.

The highest precipitation is in May and the lowest is in August. While the annual average precipitation is 254.2 mm, the evaporation is 1,094.9 mm [57]. The soil of the region were calcareous alluvial material which were formed as a result of floods of the Aras River. The horizons of the soil profile include Ap (0–35 cm), AB (35–55 cm), B (55–96 cm), BC (96–170 cm) in the region. According to the WRB system, the soil is a Calcaric Pantofluvic Fluvisol (Loamic, Aric, and Densic) [58]. The soil physical and chemical properties are shown in Table 1.

Temperature sensors (Elitech RC-4) were used in the study. The sensors were placed at different depths (0, 5, 10, 15, 20, 40 cm) for both alkaline and non-alkaline soil. They were programmed to collect temperature data hourly and data were received via computer. The sensors collected data for a year to determine the seasonal changes in soil thermal properties of both non-alkaline and alkaline soils.

A basic diagram of the sensor’s installation position is shown in Figure 1.

3.2. Soil analysis

For determining soil properties, disturbed and undisturbed soil samples were collected from both non alkaline and alkaline soils with three replicates. Besides, disturbed soil samples were collected to determine soil moisture content for each season from both non alkaline and alkaline soils with three replicates. Soil texture, soil organic carbon, pH and electrical conductivity were determined in disturbed soil samples while soil bulk density was determined in undisturbed soil samples [50,59,60,61,62,63]. Exchangeable sodium percentage (ESP) was calculated according to Shahid et al. [64].

3.3. Calculation of thermal properties of soil

Thermal properties of the soils: volumetric heat capacity (C_v), thermal diffusivity (κ), thermal conductivity (λ), attenuation depth (d), thermal effusion (e), heat wave velocity (u), and length heat wave (Λ) were calculated using the above-described formulas (2.30)-(2.40).

3.3.1. Volumetric heat capacity

For determining volumetric heat capacity of the soil; specific heat capacity, bulk density and soil moisture content were calculated with the Eq. (2.39):

$$C_v=C_{m,s}\cdot {\rho }_b+C_{v,w}\cdot \theta ,\mathrm{where}{C}_{m,s}=C_{m.org}\cdot {\displaystyle \frac{m_{org}}{m}}+C_{m.\min }\cdot \left(1-{\displaystyle \frac{m_{org}}{m}}\right)$$

(3.1)

where C_m,s is the specific heat capacity of the soil solid phase, J/(kg °C); ρ_b is the soil bulk density, kg/m³; C_v,w is the volumetric heat capacity of soil moisture, equals to 4186.6 kJ/(m³ С); θ is the volumetric moisture content, m³/m³; C_m.org and C_m.min are the specific heat capacities of organic and mineral components of the soil solid phase, respectively, J/(kg °C); m_org is the mass of soil organic matter, kg; m is the soil mass, kg; and m_org/m is the organic matter content in the soil, %.

3.3.2. Calculation of thermal diffusivity of soil

The determination of thermal diffusivity coefficient of the soil has been discussed in many experimental and theoretical works [19,21,23,24,26,27,40].

Thermal diffusivity will be determined first by existing (classical) layerd methods using formulas (2.30)-(2.33) and then by the proposed point methods using formulas (2.34)-(2.38). These methods are developed based on the boundary condition ∂T(z→∞, t)/∂z=0 and ∂T(x=L, t)/∂z=0.

3.3.3. Calculation of other thermal properties

Other thermophysical parameters such as thermal conductivity (λ), damping depth (d), thermal effusion (e), heat wave velocity (υ) and heat wave length (Λ) will be calculated using using formulas (2.40).

3.3.4. Heat Flux

Using the formula (2.42) found for the boundary condition ∂T(z→∞, t)/∂z=0, for m=1 and m=2, to calculate the heat flux q₁(z=0,t) on the soil surface (z=0) at time t, we have [23,24]:

$$q_1\left(z=0,t,m=1\right)=q_{1,z=0}^{m=1}\left(t\right)=B_1\cdot \cos \left[{\displaystyle \frac{\pi }{4}}+\left(\omega t+{\epsilon }_1\right)\right]$$

(3.2)

$$q_1\left(z=0,t,m=2\right)=q_{1,z=0}^{m=1}\left(t\right)=B_1\cos \left[{\displaystyle \frac{\pi }{4}}+\left(\omega t+{\epsilon }_1\right)\right]+B_2\cos \left[{\displaystyle \frac{\pi }{4}}+\left(2\omega t+{\epsilon }_2\right)\right]$$

(3.3)

Similarly, using the formula (2.43)-(2.44) found for the boundary condition ∂T(z=L, t)/∂z=0 for m=1 and m=2, to calculate the heat flux q₂(z=0,t) on the soil surface (z=0) at time t, we have:

$$q_2\left(z=0,t,m=1\right)=q_{2,z=0}^{m=1}\left(t\right)={\Gamma }_1\left\{\sin \mathrm{h}\left(2b\right)\cos \left[{\displaystyle \frac{\pi }{4}}+\left(\omega t+{\epsilon }_1\right)\right]-\sin \left(2b\right)\cos \left[{\displaystyle \frac{\pi }{4}}-\left(\omega t+{\epsilon }_1\right)\right]\right\}$$

(3.4)

$$\begin{array}{l}{q}_2\left(z=0,t,m=2\right)=q_{2,z=0}^{m=2}\left(t\right)={\Gamma }_1\left\{\sin \mathrm{h}\left(2b_1\right)\cos \left[{\displaystyle \frac{\mathsf{\pi }}{4}}+\left(\omega t+{\epsilon }_1\right)\right]-\sin \left(2b_1\right)\cos \left[{\displaystyle \frac{\mathsf{\pi }}{4}}-\left(\omega t+{\epsilon }_1\right)\right]\right\}+\\ +{\Gamma }_2\left\{\sin \mathrm{h}\left(2b_2\right)\cos \left[{\displaystyle \frac{\mathsf{\pi }}{4}}+\left(2\omega t+{\epsilon }_2\right)\right]-\sin \left(2b_2\right)\cos \left[{\displaystyle \frac{\mathsf{\pi }}{4}}-\left(2\omega t+{\epsilon }_2\right)\right]\right\}\end{array}$$

(3.5)

$$\mathrm{Where}{{\rm B}}_1={{\rm T}}_1{C}_v\sqrt{\mathsf{\omega }\mathsf{\kappa }},{{\rm B}}_2={{\rm T}}_2{C}_v\sqrt{2\omega \kappa }$$

(3.6)

$$\begin{array}{l}{\Gamma }_1={\displaystyle \frac{{{\rm T}}_1}{\mathsf{\Delta }\left(b_1\right)}}\left(\lambda \sqrt{{\displaystyle \frac{\omega }{\kappa }}}\right),{\Gamma }_2={\displaystyle \frac{{{\rm T}}_2}{\mathsf{\Delta }\left(b_2\right)}}\left(\lambda \sqrt{2{\displaystyle \frac{\omega }{\kappa }}}\right),b_1=L\sqrt{{\displaystyle \frac{\omega }{2\kappa }}},b_2=L\sqrt{2{\displaystyle \frac{\omega }{2\kappa }}},\\ \mathsf{\Delta }\left(b_1\right)=\cos \mathrm{h}\left(2b_1\right)+\cos \left(2b_1\right),\mathsf{\Delta }\left(b_2\right)=\cos \mathrm{h}\left(2b_2\right)+\cos \left(2b_2\right),\omega =2\mathsf{\pi }/{\tau }_0\end{array}$$

(3.7)

3.3.5. Calculation of parameters of the soil surface

Daily average temperature of the soil surface (T₀), amplitude (T_a), phase (ε) were calculated from the Eqs (2.11) and (2.12).

3.4. Comparison of Methods

The performance of the classical-layer (M1-M4) and proposed (M5-M8) methods will be evaluated using the most common model selection criteria described above by formulas (2.45)-(2.52). These criteria will also be used in determining the soil surface parameters, i.e. T₀, T₁, T₂, ε₁ and ε₂.

4. Results and Discussion

The different soil temperature values were observed in non-alkaline and alkaline soil. Figure 2a – 2b show the mean temperature value of the seasons. Figures show that, the highest soil temperatures in the investigated soil depths and seasons were determined in alkaline soils. While it was measured as 17.59 °C in alkaline soils at 5 cm soil depth in spring, it is 15.91 °C in non-alkaline soils. This can be explained with the higher bulk density and higher volumetric heat capacity in alkaline soils.

As a matter of fact, soil temperature increases with the increase in particle contact as a result of the decrease in air-filled pores between soil particles due to the increased bulk density and the combination of solid particles. In addition, soil alkalinity causes soil aggregates to break down and micro aggregates to increase [65].

4.1. Calculations of the parameters of the soil surface zone

To determine the parameters of the active soil surface (T₀, T_i and ε_i) and their statistical approximation indices, we adopted one and two harmonics in condition (2.9).

Based on the results of temperature measurements for z = 0, i.e. T (z = 0, t_i), the parameters of the distribution of the surface temperature of the studied soils were determined using the least squares method.

At the same time, six metrics were used by formulas (2.46)-(2.48), (2.53)-(2.55) i.e: R², R²_adj, D, UI, σ, AIC_c.

Table 2 and Table 3 gives results of calculation of the parameters (T₀, T_i, and ε_i), and also statistical characteristics of approximation between T_mes (0, t_i)– the measurement initial data, and T_cal (0, t_i) – the calculatin data computed from formula (2.9) at m = 1 and m = 2 for non-alkaline and alkaline soils.

As can be seen from Table 2 and Table 3 the introduction of the second harmonic makes it possible to determine with high accuracy the parameters of the temperature distribution on the soil surface. It was determined that the average temperature (T₀), amplitude (T_a), phase (ε) of the soil surface were higher in alkaline soils in all seasons (Table 3). Daily average temperature of the soil surface (T₀), amplitude (T_a), phase (ε) were calculated from the Eqs (2.11) and (2.12). The results of these calculations are given in Table 4.

	Nonalkaline soils					Alkaline soils
	Winter	Spring	Summer	Autumn		Winter	Spring	Summer	Autumn
To	1.8643	12.9127	24.5478	15.0455	To	2,4052	14,8979	28,9658	17,3838
Ta1	1.0327	3.4375	5.4270	3.9363	Ta2	4,4773	7,9065	14,6342	13,1626
ε1	1.8468	2.1760	2.3607	2.2011	ε2	2,3250	2,5422	2,5776	2,6118
Cv	1878.49	2126.30	1596.04	1746.15	Cv	2014.11	2349.05	1637.30	1846.64

4.2. Soil Thermal Properties

Thermal properties of the soil (volumetric heat capacity, thermal diffusivity, thermal conductivity, damping depth, thermal effusivity, heat wave velocity, and heat wave length) changing depending on the season, soil and calculation models (М1-М8), are given in Table 5 and Table 6.

The calculations of these parameters are described in more detail below.

4.2.1. Volumetric heat capacity

The calculated values ​​of the volumetric heat capacity of soils are given in Table 5 and Table 6. As can be seen, the volumetric heat capacity of soils varied depending on the soil and the time of year.

For example, the values of specifi heats of the organic and mineral components of the soils under study are equal to C_m,org = 1925.928 J kg^-1 °C^-1 or 0.46 cal g^{-1 o}C and Cm,min = 753.624 J kg^-1 °C^-1 or 0.18 cal g^{-1 o}C respectively.

Based on results of the analyses of soil properties at the experimental site for layer 0.1–0.2 m, the mean bulk density of (nonalkaline) soil ρ_b was 1023.3 kg m^-3, the content of organic substance, to m_org/m = 1.62 % and volumetric moisture content, to θ = 0.1521 m³m^-3(tabl.1)

$$\begin{array}{l}{C}_{ms}=C_{morg}\cdot {\displaystyle \frac{m_{org}}{m}}+C_{m\min }\cdot \left(1-{\displaystyle \frac{m_{org}}{m}}\right)=1925.928{\displaystyle \frac{J}{kg°C}}\cdot 0.0162+753.624{\displaystyle \frac{J}{kg°C}}\cdot \left(1-0.0162\right)=\\ =31,2000{\displaystyle \frac{\mathrm{J}}{\mathsf{к}\mathsf{г}\cdot {}_{\mathrm{o}}\text{С}}}+741,4152912{\displaystyle \frac{\mathrm{J}}{\mathrm{kg}\cdot {}_{\mathrm{o}}\text{С}}}=772,6153{\displaystyle \frac{\mathrm{J}}{\mathrm{kg}\cdot {}_{\mathrm{o}}\text{С}}}=0.772615{\displaystyle \frac{\mathrm{kJ}}{\mathrm{kg}\cdot {}_{\mathrm{o}}\text{С}}}\end{array}$$

Next, using the formulas for calculating C_v and considering that for the 0.1-0.2 m layer of our soil θ = 0.1521 m³m^-3, C_ms = 0.772615 kJ/(kg ^oC), and the specific volumetric heat capacity of water C_v,w = 4186.8 kJ/(m ^oC), we calculate the value of volumetric heat capacity for the 0.1-0.2 m layer of non-alkaline soil as follows:

$$\begin{array}{l}{C}_v=C_{ms}\cdot {\mathsf{\rho }}_b+C_{vw}\cdot \mathsf{\theta }=0.772615{\displaystyle \frac{\mathrm{kJ}}{\mathrm{kg}\cdot {}_{\mathrm{o}}\text{С}}}\cdot 1023.3{\displaystyle \frac{kg}{m^3}}+4186.8{\displaystyle \frac{\mathrm{J}}{m^3\cdot {}_{\mathrm{o}}\text{С}}}\cdot 0.1521{\displaystyle \frac{m^3}{m^3}}=\\ =790.627162{\displaystyle \frac{kJ}{m^3{}_{\mathrm{o}}\text{С}}}+636.81228{\displaystyle \frac{kJ}{m^3{}_{\mathrm{o}}\text{С}}}=1427,44{\displaystyle \frac{kJ}{m^3{}_{\mathrm{o}}\text{С}}}\end{array}$$

Similarly, the values of volumetric heat capacity of non-alkaline soil can be easily calculated for all soil layers and for all seasons.

For alkaline soil, ρ_b was 1125.6 kg m^-3, organic matter content was m_org/m = 1,24%, and volumetric moisture content was θ = 0.1628 m³ m^-3.(tabl.1)

To give an example, the value of specific heat capacity (C_ms) for a 0,1–0,2 m layer of solid part of non-alkaline soil calculated by formula (2) would be as follows:

$$\begin{array}{l}{C}_{ms}\left(10-20cm\right)=C_{morg}\cdot {\displaystyle \frac{m_{org}}{m}}+C_{m\min }\cdot \left(1-{\displaystyle \frac{m_{org}}{m}}\right)=1925,928{\displaystyle \frac{\mathrm{J}}{\mathrm{kg}\cdot {}_0\text{С}}}\cdot 0.0124+753,624{\displaystyle \frac{\mathrm{J}}{kg\cdot {}_0\text{С}}}\left(1-0.0124\right)=\\ =23,8815072{\displaystyle \frac{\mathrm{J}}{\mathrm{kg}\cdot {}_0\text{С}}}+744,2790624{\displaystyle \frac{\mathrm{J}}{\mathrm{kg}\cdot {}_0\text{С}}}=768,16057{\displaystyle \frac{\mathrm{J}}{\mathrm{kg}\cdot {}_0\text{С}}}=0,768161{\displaystyle \frac{\mathrm{kJ}}{\mathrm{kg}\cdot {}_0\text{С}}}\end{array}$$

Next, using the formulas for calculating C_v and considering that for the 0.1-0.2 m layer of our soil θ = 0.1628 m³m^-3, C_ms = 0.768161 kJ/(kg ^oC), and the specific volumetric heat capacity of water C_v,w = 4186.8 kJ/(m ^oC), we calculate the value of volumetric heat capacity for the 0.1-0.2 m layer of non-alkaline soil as follows:

$$\begin{array}{l}{C}_v\left(10-20cm\right)=C_{ms}\cdot {\mathsf{\rho }}_b+C_{vw}\cdot \mathsf{\theta }=768,16057{\displaystyle \frac{J}{kg\cdot {}_0\text{С}}}\cdot 1.125,6{\displaystyle \frac{kg}{m^3}}+4186,8{\displaystyle \frac{kJ}{m^3\cdot {}_0\text{С}}}\cdot 0,1628{\displaystyle \frac{m^3}{m^3}}=\\ =864,664582{\displaystyle \frac{\mathrm{kJ}}{m^3\cdot {}_0\text{С}}}+681,5165021{\displaystyle \frac{\mathrm{kJ}}{m^3\cdot {}_0\text{С}}}=1546,181084{\displaystyle \frac{kJ}{m^3\cdot {}_0\text{С}}}\end{array}$$

Calculations of C_v values for non-alkaline and alkaline soils for all seasons and for all soil layers, are given in Table 5 and 6.

4.2.2. Thermal Diffusivity (κ)

The values of the thermal diffusivity parameter of the studied soils, for all season were calculated using the formulas (2.30)-(2.38) and are given in Table 5 and Table 6.

The estimated values of thermal diffusion calculated using classical (М1-М4) and proposed (М5-М8) methods show that there are significant deviations between them.

All the considered thermal properties varied depending on the time of year, the models used and the soils (Table 3-4). Tong et al. [44] also reported that the thermal diffusivity differed according to seasons.

Considering the differences between soils, higher thermal diffusivity was observed in non-alkaline soils, with the exception of the values ​​calculated by methods M1-M4.

These values ​​may be due to low values ​​of the average daily soil surface temperature (T_o), the amplitude of fluctuations (T_a) of the soil surface temperature, the phase shift (_ε) and the volumetric heat capacity in non-alkaline soils (Table 2 and Table 3) or Table 4.

When the thermal diffusivity results of the soils are examined, there are quite high values ​​in the classical models (Table 5).

However, more moderate results were obtained in the proposed models

For example, the thermal diffusivity values calculated by the amplitude method (М1) were 7.2882х10^-6 m s^-1 in winter, 92.768х10^-6 m s^-1 in spring, 209.4553х10^-6 m s^-1 in summer and 35.4889х10^-6 m s^-1 in autumn in alkaline soil.

Whereas, the thermal diffusivity values calculated by the (M8) method were 0.2371х10^-6 m s^-1 in winter, 0.3214x10^-6 m s^-1 in spring, 0.1935x10^-6 m s^-1 in summer and 0.1413x10^-6 m s^-1 in autumn in alkaline soil.

This deviation shows that classical models are incorrect due to theoretical assumptions in setting boundary conditions at depth, i.e. x=L.

Consequently, the parameters of the soil’s active surface (T₀, T_a, and ε), also affects soil thermal diffusion.

The results of the thermal diffusivity values calculated by methods 1-4 differ greatly depending on the time of year. This indicates the mathematical incompatibility of these methods. In contrast, the κ-values calculated by the improved method do not differ significantly.

The results of our studies (Table 5 and Table 6) showed that the thermal properties (temperature, thermal conductivity, etc.) of soils decreased with increasing salt content at a given water content are consistent with the findings [8].

Also, unlike the others, the improved method in addition to the amplitude of fluctuations (Ta) of the soil surface temper-ature, also take into account the depth of the soil profile (L), which are important parameters.

4.2.3. Thermal conductivity (λ)

The values of the thermal conductivity parameter of the studied soils were calculated using the formula λ=κC_v and varied depending on the model, season and soil (Table 5 and Table 6).

As in the case of thermal diffusivity, calculations by classical methods (M1-M4) give a scatter of values of the thermal conductivity parameter λ

When assessing the difference between soils, higher values ​​of thermal conductivity calculated using the proposed methods (M5-M8) were observed in non-alkaline soils for all seasons.

Guo et al [66] also reported that the thermal conductivity of the upper soil layers did not increase with a higher solid ratio in the solonchak and solonetz, as the salts in the upper layers contributed to a lower thermal conductivity. Besides some researchers determined that the effects of salts in soil on thermal conductivity varied with different water contents [67,68]. At medium water content, thermal conductivity decreased to varying degrees in different saline soils.

4.2.4. Damping depth (d)

The values of the damping depth of the soil temperature were calculated using the formula (2.40), i.e. d=√τ₀κ/π and varied depending on the models (M1-M8) used to calculate the κ parameter, as well as the season and soil (Table 5-6).

When analyzing the results of temperature measurements in both soils, it is evident that at depths z ≤ 15 cm, temperature waves attenuate (Fig. 1 and Fig. 2).

The estimated values ​​of damping depth calculated using classical (М1-М4) and proposed (М5-М8) methods show that there are significant deviations between them.

For example, the damping depth values calculated by the amplitude method (М1) were 44.77 cm in winter, 159.73 cm in spring, 240.01 cm in summer and 31.52 cm in autumn in alkaline soil.

This situation is similar for the other methods (M2, M3 and M4). Excessive deviations between the obtained values reduce the reliability of these methods. These results are not consistent with measured temperatures in alkaline soils (Figure 2)

The damping depth (d) calculated using the values of k obtained by formulas (2.30-2.33) turned out to be greater than 10 cm, i.e. z >10

However, the damping depths (d) according calculated using the thermal diffusivity (κ) values ​​obtained from formulas (2.34)–(2.38), was calculated as ~ 15 cm on average all seasons for non-alkaline soil (Table 5).

For alkaline soil, the results of calculations using to the point methods (formulas (2.34)–(2.38)) showed that the damping depths were approximately equal to 8 cm on average for all seasons (Table 6).

For alkaline soil, the results of calculations using to the point methods (formulas (2.34)–(2.38)) showed that the damping depths were approximately equal to 8 cm on av-erage for all seasons (Table 6).

These findings agree with the observed data. In line with these findings, it is seen that the most adaquate model for this study is the point methods.

4.2.5. Thermal effusivity (e)

The values of the thermal effusivity parameter of the studied soils were calculated using the formula e=C_v√κ and varied depending on the models (M1-M8) used to calculate the κ parameter, season and soil (Table 5-6).

The estimated values ​​of thermal effusion calculated using classical (М1-М4) and proposed (М5-М8) methods show that there are significant deviations between them.

For example, the thermal effusivity values calculated by the amplitude method (М1) were 5437.41 W·h^0.5∙m^-2∙ ^oC^-1 in winter, 22625.17 W·h^0.5∙m^-2∙ ^oC^-1 in spring, 23695.88 W·h^0.5∙m^-2∙ ^oC^-1 in summer and 11000.88 W·h^0.5∙m^-2∙ ^oC^-1 in autumn in alkaline soil.

This situation is similar for the other methods (M2, M3 and M4). Excessive deviations between the obtained values reduce the reliability of these methods.

In the point method (M8), according to the seasons; 980.67 W·h^0.5∙m^-2∙ ^oC^-1 in winter, 1331.69 W·h^0.5∙m^-2∙ ^oC^-1 in spring, 720.16 W·h^0.5∙m^-2∙ ^oC^-1 in summer and 698.66 W·h^0.5∙m^-2∙ ^oC^-1 in autumn in alkaline soil (Table 4). A similar situation is observed in non-alkaline soils; 1630.40 W·h^0.5∙m^-2∙ ^oC^-1 in winter, 2096.83 W·h^0.5∙m^-2∙ ^oC^-1 in spring, 1676.63 W·h^0.5∙m^-2∙ ^oC^-1 in summer and 1547.45 W·h^0.5∙m^-2∙ ^oC^-1 in autumn.

The highest thermal effusivity values were observed in non-alkaline soils. These results can be explained by the fact that the thermal diffusivity, which is a factor of the heat absorption, is higher in these models.

4.2.6. Heat wave velocity (ϑ)

The values of the heat wave propagation velocity were calculated using the formula (2.40), i.e. ϑ=2√πκ/τ₀ [23,52] and varied depending on the models (M1-M8) used to calcu-late the κ parameter, as well as the season and soil (Table 5-6).

Using the proposed method (model M8), the heat wave velocity values it was found as 1.0467 m s^-1 in winter, 1.1893 m s^-1 in spring, 1.2669 m s^-1 in summer and 1.0688 m s^-1 in autumn in non-alkaline soil; 0.5872 m s^-1 in winter, 0.6837 m s^-1 in spring, 0.5305 m s^-1in summer and 0.4563 m s^-1 in autumn in alkaline soil.

As can be understood from the text the heat wave velocity was higher in non-alkaline soil for all seasons. The maximum heat wave velocity was in summer at non-alkaline soil

4.2.7. Heat wave lenght (Λ )

The value of the heat wave length is the distance traveled by the heat wave per unit of time.

The values of the heat wave lenght calculated using the formula Ʌ=2πd and varied depending on the models (M1-M8) used to calcu-late the κ parameter, as well as the season and soil (Table 5-6) [53].

The heat wave lenght was detected as 0.9 m in winter, 1.03 m in spring, 1.09 m in summer and 0.92 m in autumn in non-alkaline soil; 0.5 m in winter, 0.59 m in spring, 0.45 m in summer and 0.39 m in autumn in alkaline soil.

It was concluded that the lenght between the two waves in non alkaline was much more than in alkaline soil.

The difference can be explained by the higher thermal diffusivity of non-alkaline soils.

The difference can be explained by the higher thermal diffusivity and the lower soil temperature attenuation depths of non-alkaline soils.

Soil thermal properties (κ, λ, d, e, ….) different in each model used.

4.3. Assessment of soil thermal diffusivity models

The measured (T_mes) and predicted temperature values (T_cal) (using solutions (2.20) and (2.23) for one and two harmonics) were compared to evaluate the effectiveness.eight methods (M1-M8).

To do this, we set up a linear regression equation, i.o. Times (zi,tj) =a+bt call(zi,tj), between measured and estimated value of studied soil temperature for all four season.

The effectiveness of the eight methods (M1-M8) was assessed using on the basis of two criteria: the Pearson correlation coefficient (r) and the The Root Mean Square Error (RMSE, σ), are described by the equations (2.45) and (2.48). The calculation results are presented in Table 7 and Table 8.

The analysis of the calculation results (see Table 7 and Table 8) showed that basically the proposed point method (M8) (formula 2.38) is the most efficient model as it gives more accurate predictions in soil layers x= 5, 10 and 15 cm for T(z,t) for both soils in all seasons than the other algorithms.

4.4. Heat flux (q)

According to tables 7 and 8, it was determined that the proposed method M8 is a more adequate model.

Therefore, the values of thermal diffusivity found according to the proposed method M8 were taken into account when calculating the heat flux from the soil surface.

Using the thermal diffusivity values (Table 5 and Table 6) calculated using the proposed point method (M8), the heat flux q on the surface (z = 0) of the soil at time t in both (alkaline and non-alkaline) soils was calculated.

These values ​​varied depending on the season and soils (Figure 5). In non-alkaline soils, the highest values of heat flux was 21.29 W·m^-2 at 14 o'clock in winter, 80.12 W m^-2 at 12 o'clock in spring, 106.86 W m^-2 at 12 o'clock in summer, and 76.08 W m^-2 at 13 in autumn has been determined. The highest heat flux in non-alkaline soils was observed in summer (Figure 5).

Maximum heat flux was 58.48 W·m^-2 at 12 o'clock in winter, 128.90 W·m^-2 at 11 o'clock in spring, 130.35 W·m^-2 at 11 o'clock in summer, and 121.09 W·m^-2 at 11 o'clock in autumn. The highest heat flux was observed in summer (Figure 5).

The highest values ​​of heat flow were observed in alkaline soils for all seasons. These values ​​can be explained by the higher amplitude (T_a), phase shift (ε) and volumetric heat capacity of soils (C_v), which are included in the corresponding formulas (3.2) - (3.7) and are the main factors influencing the heat flux of soils (tables 2 and 3 or table 4).

5. Conclusions

In the study, the thermal properties of different soils (alkaline and non-alkaline) were investigated using various classical (amplitude, arctangent, logarithmic and phase) and proposed point methods in different seasons from 01.11.2019 to 30.10.2020.

All the considered thermal properties varied depending on the time of year, the models used and the soils (Table 3-4).

Considering the differences between soils, higher thermal diffusivity was observed in non-alkaline soils, with the exception of the values ​​calculated by methods M1-M4.

The M8 model proved to be the best among the eight models in estimating soil temperature values due to the strong correlation between calculated and measured values.

For example, the correlation value (r) at a depth of z=5 cm in non-alkaline soil was, respectively, r=0.9838 in winter, 0.9956 in spring, 0.9944 in summer, and 0.9952 in autumn.

Whereas for alkaline soil at a depth of z = 5 cm, respectively, we have: r=0.9726 in winter, 0.9540 in spring, 0.9903 in summer and 0.9954 in autumn

Accordingly, the following values were obtained for the root mean square error (RMSE, σ) at a depth of z=5 cm: 0.2494 ^oС in winter, 0.2488 ^oС in spring, 0.5228 ^oС in summer and 0.3279 ^oС in autumn in non-alkaline soil.

Whereas for alkaline soil at a depth of z = 5 cm, respectively, we have: σ=1.1284 in winter, 1.7239 in spring, 5.2323 in summer and 2.9487 in autumn in alkaline soil.

The calculation results showed that the proposed point method (M8) (formula 36) is the most effective model, because it gives more accurate predictions for soil temperature values T(z,t) than other algorithms for both soils in all seasons.

The proposed point method M8 is more informative and adequately reflects heat transfer in the soil, since it includes, in addition to the values ​​of the soil profile temperature, also the amplitude of fluctuations (Ta) of the soil surface temperature as well as the depth of the soil profile (L), where there is no temperature gradient.

When comparing the attenuation depths that best reflect the attenuation of heat waves in soils, it was found that the best of the eight models used to calculate the soil thermal conductivity coefficient is the proposed M8 point method for both soils for the entire season.

The results of our studies (Table 4 and Table 5) showed that the thermal properties (thermal diffusivity, thermal conductivity, düşürüm depth, thermal effusivity, heat wave velocity, and heat wave length) of soils decreased with increasing salt content at a given water content are consistent with the findings of Noborio and McInnes (1993).