Chapter 10

The connected task of non-stationary heat exchange between

mine air and mining massif

B.P. Kazakov	A.W. Shalimov
Mining Institute Ural Branch of Russian Academy of Science Perm, Russia	Mining Institute Ural Branch of Russian Academy of Science Perm, Russia

ABSTRACT

The mathematical model of heat exchange between the ventilating air and the mining massif without use of factor of non-stationary heat exchange is offered. The numerical calculation of air temperature lengthwise drifts is made depending on time. The comparative analysis of the obtained results with similar results received with the help of traditional mathematical model of heat exchange is carried out.

KEYWORDS

Heat exchange, factor of non-stationary heat exchange, massif, laplace transform, bessel function, heat capacity, heat conductivity, convection, humidity, saturated vapor, factor of heat emission

A lot of works of different authors is devoted to research heat exchange processes between mine air and rocks. All the possible variety of variants of heat exchange in real potash and other mines conditions is considered with sufficient completeness. All authors are building their calculations practically on the same physical model of heat exchange developed by Sherban A.N. and Kremnev O.A. (1954), the central part of which is the so-called «factor of non-stationary heat exchange» k_t. Essence of inputing k_t - to "coordinate" the basic characteristic of heat exchange process - density of the heat flow from air into a massif j (or on the contrary) - to well known values of air temperatures T and the "not disturbed" massif T, i.e. j=k_t (T_∞ -T). Further, if k_t=k_t (t, z) - the known time function and coordinate, it is possible to say, that the task is solved for air, as it is known, how many heat is taken away or increases at air at each moment of time in every place. Thus, the task to define the air temperature is shown as a task of definition k_t.

The calculation modes of the non-stationary heat exchange factor can be parted into two types: 1) as more exact calculation is possible, in spite of the complexity, and 2) The approximate calculation with the purpose to receive the analytical formula for k_t. In general, the exact analytical formula for k_t turns out only in one case - flat task with preset, identical at any point of air temperature. As to the cylindrical task, so it has no analytical solution, even in case of constant air temperature. To the first type of calculation k_t it is necessary to refer the solution mode of the heat conductivity cylindrical equation with the preset ventilating jet temperature as a time functions with the help of the Laplace transform (Sherban, 1953). To the second type of calculation kt it is possible to refer, for example, the method offered by Voropaev A.F. (1966), in which the flat task was solved, and then the conversion was carried out to the cylindrical with the help of so-called «factor of the body form». As a result of this and a lot of other simplifications it has got rather simple analytical dependence for k_t. The formulas for k_t, received by a solution of the cylindrical heat conductivity equation with the use of the Laplace transform, represent integrals with infinite limits from combinations of Bessel functions. Besides these formulas are exact only by constant on time air temperature. Thus k_t (t, z) depends on T (t, z), and T (t, z) depends on k_t (t, z), it means that there is no task solution in an obvious way. Because, it is not obviously possible, to divide k_t and T in view of very difficult dependence.

In the work of Braitcheva, et al., (1981) it is marked, that even if to consider the transfer processes in a drift as quasistationary (in the sense, that the change of air temperature in time will be caused only by change of temperature and heat flow on an interfacial area), this task remain extremely difficult for an analytical solution. Actually, in article (Kozdoba and Tchernyk, 1977) the system of the equations of energy transfer in is shown a most general view but the solution of this task is absent.

The task is put to calculate numerically T (z, t) in a horizontal cylindrical drift. Thus it is supposed, that air temperature is given only on an input in a drift (z=0), i.e. its value is unknown, to be exact - required value. By such a task statement it is hardly meaningful to enter and to estimate such performance, as factor of non-stationary heat exchange k_t.

The task, in general, is essentially non-stationary, because the intensity of air heat exchange with walls depends on heating depth (or cooling) of a massif, and it, in turn, will be the more, than more time will pass. In an ideal case, after infinite time passing, the whole massif should have temperature of entering air.

The task is simulated as follows. There is a massif half space with volumetric heat capacity
and thermal conductivity χ_m. In a massif there is a drift of radius r₀, located on an axis z. The beginning of a drift is set up at the beginning of an axis z. The temperature of air T₀, submitted to this drift depends on time (T₀=T₀ (t)) with volumetric heat capacity
with speed v=v (z), which, generally speaking, is strictly decreasing function z (taking into account the opportunity to reduce air flow due to „branches"). The task is cylindrical symmetric, coordinate r is counted from an axis of the channel. It is supposed, that at the initial moment of time t=0 all massif has identical temperature T_∞, and all air in the channel, with an exception of cross-section z=0, has the same temperature T, the temperature "of the not disturbed" massif. At the initial moment of time the air temperature in the cross-section z=0 is equal to T₀ (0), and, further, changes in time as T₀ (t).

Further the model becomes simpler, proceeding by the following reasons. The turbulent heat conductivity of air is supposed much greater, than the molecular heat conductivity of a massif, but only in a radial direction. And in the direction of an axis z air heat conductivity is considered as molecular and, so it is much smaller than the massif heat conductivity. Such an anisotropy is given by reason, that the intensive heat transfer inside an air flow is caused by presence of whirlwinds in cross section of a flow, arising as a result of a convection. In a direction of an axis z the convection force does not work. Proceeding from above-stated, heat conductivity of air in the radial direction is considered as infinite, and in the direction of an axis z - equal to zero.

One more simplification is made, that the speed of an air movement is much greater, than speeds of heat distribution in a massif, the temperature drops on an axis z and in the radial direction will be of the different order (if on z are hundreds meters, on r-tens centimeters). It allows to neglect the heat distribution in the direction of an axis z in a massif. Thus, only radial heat conductivity in a massif is considered.

For simplification of mathematical record of the task the dimensionless coordinates are entered: distance at r and at z is measured in r₀ (radius of a drift), and time in
. The temperature will be measured from T_∞. The heat conductivity equation in a massif at cylindrical coordinates looks like:

(1)

The density of a heat flow in an air in direction to boundary j_a is equal to the density of a heat flow j_m in a massif from boundary deep into (boundary r=1). Let's write down the balance equation of a heat content Q in elementary air volume V by section
and thickness z. A heat exchange surface of this volume is
. Because an air heat conductivity in the radial direction is considered as infinite, so the heat content balance is reduced to equality of the whole heat flow outside limits of this surface and heat content changes in the given air volume:

(2)

Taking into account, that
and
, where T- temperature of air and massif at boundary r=1 (air temperature at boundary and at all cross-sections of the channel is identical, because the air heat conductivity in the radial direction - infinity), we receive a condition at the boundary:

(3)

where
,
. As the condition (3) contains a derivative on z, it should be complemented by value T in a point (r=1, z=0), i.e. air temperature in the input in the channel should be preset (this is the second boundary condition):

(4)

The entry condition should be added to the equations (1), (3) and (4):

(5)

everywhere, except of a point (r=1, z=0), where condition (4) with T₀ (0) is satisfied.

The non-stationary cylindrical task (1) - (5) was solved with the help of Laplace transform. The function T=T (r, t, z) is put in conformity with its image:

,

where p- is the complex parameter. In this case a definition range  - Re (p) > 0. The equation in the partial derivatives (1) for the original T is reduced to the ordinary differential equation for the image :

(6)

where
, and T_f (z) T(r=1, t=t_f(z), z) - temperature of " air front" at the moment of time
, when it has flown a distance z, with the initial conditions

(7)

(8)

which are received from condition (3) and (4) accordingly.

Temperature of "the front" T_f=T_f (z) can be easily calculated, as the heat exchange in this case occurs between air and "the not disturbed" massif with temperature T=0 (is counted from T_). Let " the air Echelon" of radius r₀ be isolated, which flows distance z during a time t=z/v, thus its heat content change

will be

or
,

where V and S - volume and surface of heat exchange "of the isolated «air echelon", and j=T_f (t) - density of a heat flow from air into a massif (α - heat emission factor). Taking into account, that S/V=2/r0, we receive by t→0 the differential equation for T_f (t) in a dimensional kind:

, where
.

Now, going to dimensionless variables and having put t=t_f (z), we receive
, where
. In a limiting case =  - T_f (z > 0) =0, that corresponds to indefinitely fast cooling (heating) of the "air front" up to temperature of the "not disturbed" massif, as in this case heat exchange is formally infinite.

The equation (6) is reduced to the Bessel equation, which solution is:

(9)

where J₀ both N₀- are Bessel and Neumann function of the zero order, the factors f(z) and g(z) are subjects to definition. There should be only one unknown factor, as function from z. Then at substitution (9) in (7) the differential equation for this function with the entry condition (8) turns out which unequivocally defines τ. The missing information is given by a condition on infinity (r=, p, z) =0. The correlation of factors f(z) and g(z) should be those, that at r→ - →0. It means, that (at r >1 -  =0):

(10)

Proceeding from asymptotic decomposition of functions J₀ and N₀ at r→ (Duait, 1973), it is possible to conclude, that in (10) k =-i at
and k=i at
, where i- is imaginary unit.

The Laplace transform allows to divide variables z and r and, thus, to reduce a task dimension. Now, if to put r=1, there is no necessity to consider a heat distribution in a massif. There is only information left on dependence on z, and, further, after transition to the original, - on t. In further the coordinate r is reduced, that means r=1. Substituting (9) in (7), we receive the heterogeneous differential equation for f:

(11)

,

,

here J₁ both N₁- are Bessel and Neumann function of the first order. Solving (11), and, taking in the account, that (8), we receive for:

(12)

here
, where
- the dimensionless speed, and from parameters a(z) and (z) is taken "out" the dependence from z:

and

Here
.

The restoration of the original T(t, z) using the image (p, z) is given by the formula of Laplace return transform:

(13)

where the integration is conducting along any straight line with the material coordinate x, greater than the parameter of growth function T.

In the practice the speed v(z) is a stepping decreasing function z. With the every "bifurcation" the speed decreases by jump, from "bifurcation" up to "bifurcation" - does not change. If to number all "bifurcations", the integral W(z) can be shown to the sum:

,

where i - quantity of "bifurcations" up to coordinate z, j- "bifurcation" number, z_j- the coordinate of j-"bifurcation", _j- air speed on a site [z_j, z_j+1].

As a result of heat exchange of ventilating air with a massif there can be a condensation of a moisture from air, if air is cooled by a massif, and the evaporation of a moisture from brine, if air, on the contrary, is heated up. First takes place in the summer, second - in winter. In general, the process of the "evaporation "- condensation" depends on a lot of parameters, and it is difficult enough to calculate it in details. For the estimated calculation of heat influence of phase transition on a temperature structure the following model is offered. First: the relative humidity of entering air - limiting, i.e. 100 % at contact to a dry surface and 70 % at contact to saturated brine (pressure of the saturated vapor above brine makes 70 % from pressure of the saturated vapor above fresh water). Second: the relative humidity remains still constant (limiting) along the whole way of air flow, i.e. at heating the relative humidity should decrease, but due to the evaporation remains former, and similar - for cooling. The formula (2) will change as follows:

(14)

where j_f - is the density of heat outflow from air to the evaporation (< 0) or the density of heat inflow from the condensation (> 0). j_f - there is a birth of heat at the whole cross-section in unit of time to unit of a contact surface F, i.e.

,

Where r- specific heat of evaporation, and dm_f - mass allocated (> 0) or absorbed water (< 0).

,

C₀- the concentration of saturated vapor in air (100 % or 70 %), V- volume of considered "echelon" of air. Substituting all ¸ in (14), we receive the same formula (3), with only that difference, that v (z) is replaced on

.

For real parameters the numerical value in brackets is equaled approximately to three, therefore the heat estimation of the phase transition is similar to the speed increase of air in three times.

Analytical analogue of the formula (13), received by Voropaev A.F. (1966), looks in our designations like:

(15)

The formula satisfies to all limiting cases: at z=0 or at
t = - T=T0, at z = ∞ - T=0 (temperature is counted from T), but in a fact at t=0 is impossible T=0 (T=0 only in the case α=∞). However the author notices, that the formula works in the time interval from 1 hour till
50 years to within 6 %.

The calculation under the formulas (13) and (15) was made for the potash mines conditions of the Verhnekamsk deposit at the following parameter values of the task: temperature "of the not disturbed" massif +8°C, temperature of air in an input is constant and equal +15°C, the speed of air is 3 m/s, the drift diameter is 4 m, the dimensionless complexes a=1200000 and b=2692 (contain the information about the thermal conductivity of a massif and about the relation of volumetric heat capacities of a massif and air). The absorption and separation of heat by the evaporation and the condensation were not considered, the heat emission factor α was assumed equal to infinity, and the air speed - constant, not dependent on z. The results are submitted at figures 1a (formula (13)) and 1b (formula (15)) for the time interval per 10 years from a beginning of an ventilation with the interval per 1 year and at figures 2a (formula (13)) and 2b (formula (15)) for the time interval per 100 years from a beginning of a ventilation with the interval per 10 years. The comparative analysis of the profiles has shown the significant quantitative divergences of results, namely, according to calculations by the method of Voropaev, the air warms up a massif much faster, than by the calculation with (13). The warming up intensity differs approximately in 5-8 times. So, in 10 years at distance of 10 km according to (13) the temperatures should to rise on 1/4 degrees and almost to 2 degrees, according to (15), and in 100 years - on 3/4 and on 4 degrees accordingly.

The numerical calculation (13) for a case T₀=T₀ (t) was also made with the purpose to take into consideration the seasons temperature fluctuations of entering air and to define, on what distance from a shaft these fluctuations are appreciable. The reference temperature T₀(t) undertook as the periodic function with the period per one year. Thus the year Θ in dimensionless units is broken in 12 equal intervals appropriate to months of year. T0 (t) was simulated as follows. During n-th of quantity of months (from 0 up to 12) T₀(t) =T₁, and during the left 12-n of months T₀ (t) =T₂. In this case image T₀(t) looks like:

.

For the Urals climate n=7 (the heated period of year since September till March), T₁ = + 2C, T₂ = + 15C (in a dimensional kind). The calculation structures of the temperature air fluctuations in drifts within one year are submitted at a figure 3 (2 winter and 2 summer months - 2, 5-th, 8-th and 11 months, beginning from September). The calculations were carried out for different intervals of time from a beginning of a drift operation (ventilation) from 1 year till 100 years. It turned out, that the difference in the received curves in an interval from 2 till 100 years is not appreciable, that makes the task "as if " stationary -the temperature structures practically do not drive (as to 1-st year, so, the profiles differ a little from the others, because the structures T were not established). Actually, the speed of heating (cooling) is so small, that it is simply imperceptible. And the reason it that T1 < T_, and T2 > T_ and approximately to identical value.

Thus, it is possible to tell, that at a fig. 3 the temperature profiles of the established thermal mode are represented at preset values of the task parameters (with a clause, certainly, that it is only « the apparent stability »). The temperature fluctuations are appreciable on a distance up to 5 km from a shaft, and it will be well coordinated to the experimental results.

The advantage of the offered model of heat exchange of ventilating air with a mining massif, on our sight, is that a subject of calculations is the real physical characteristic - temperature, which calculation error can be estimated and made as much as small by the numerical calculation. Models, to calculate k_t, even, if the error estimation of kt is made, hardly can unequivocally define the calculation error of temperature with the help of received k_t. You see, k_t is necessary to be integrated according to time, and the mistake will be accumulating the more, than the more the time interval will be.

references

Braitcheva N.A., Chernyak W.P., Sherban A.N., 1981, The calculation methods of ventilating air temperature in underground structures, Kiev

Kozdoba L.A., Chernyak W.P., 1977, “The physical characteristics and the mathematical description of system " massif - drift" in a connection with the problem of the forecast and the regulation of a thermal mode of deep mines”. - In the book: The thermal mode of the deep coal mines and metal ore mines. - Kiev, Nauk. Dumka, pp. 40-49

Kremnjev O.A. 1954, The heat exchange between a ventilating jet and mining massifs of old mines and drifts. A. of Sc. of the Ukraine, N10, pp. 12-17

Sherban A.N., 1953, The bases of the theory and methods of thermal calculations of a mine air, Ugletehizdat, Moscow

Voropaev A.F. 1966, The theory of heat exchange of a mine air and rocks in the deep mines, 252 p

Duait G.B. 1973, The tables of integrals and the other mathematical formulas. 228 p

2

I SZKOŁA AEROLOGII GÓRNICZEJ 1999

5

68

PROCEEDINGS OF THE 7^TH INTERNATIONAL MINE VENTILATION CONGRESS

67

The connected task of non-stationary heat exchange between

---