Determination of Thermal Diffusivity of Material by the Numerical-Analytical Model of a Semi-Bounded Body

A mathematical description of a material’s thermal diffusivity ae in a semi-bounded body is proposed with a relatively simple algorithm for its numerical and analytical calculation by solving the inverse problem of thermal conductivity. To solve the problem, it is necessary to obtain the temperature values of the unbounded plate as a result of a thermophysical experiment. A plate can be conditionally considered as a semi-bounded body as long as Fourier number Fo ≤ Foe (Foe = 0.04–0.06). It is assumed that the temperature distribution over a cross-section of the heated layer of the plate with thickness R is sufficiently described by a power-like function whose exponent depends linearly on the Fourier number. A simple algebraic expression is obtained for calculating ahc in time interval Δτ from the dynamics of temperature change T(Rp, τ) of a plate surface with thickness Rp heated at boundary conditions of the second kind. Temperature T(0, τ) of the second surface of the plate is used only to determine end time τe of the experiment. Moment of time τe, at which the temperature perturbation reaches adiabatic surface x = 0, can be set by the condition T(Rp, τe) – T(0, τ = 0) = 0.1 K. An approximate method of calculating the dynamics of changes in depth of heated layer R by values of Rp, τe, and τ is proposed. The calculation of ahc for time interval Δτ is reduced to an iterative solution of a system of three algebraic equations by matching the Fourier number, for example, using a standard Microsoft Excel procedure. Estimation of the accuracy of calculation of ahc at radiation-convective heating was performed using the initial temperature field of the refractory plate with thickness Rp = 0.05 m, calculated by the finite difference method under initial condition T(x, τ = 0) = 300 (0 ≤ x ≤ Rp). The heating time was 260 s. Calculation of ahc, i has been performed for ten time moments τi + 1 = τi + Δτ, Δτ = 26 s. Average mass temperature of the heated layer for the entire time was τeT = 302 K. The arithmetic-mean absolute deviation of ae (T = 302 K) from the initial value at the same temperature was 2.8%. Application of the method will simplify conducting and processing experiments to determine thermal diffusivity of materials.


INTRODUCTION
Thermophysical characteristics of materials are used to simulate temperature fields in elements of equipment, materials that are heat treated, or used in furnace enclosures. The reliability of mathematical modeling of temperature fields depends to a great extent on the accuracy of the values of thermophysical characteristics of materials.
Determining thermophysical characteristics of materials does not lose relevance, since new materials appear, the number of which grows continuously, and the properties of known materials vary during their operation.
Methods for determining thermophysical characteristics of materials are mainly based on solving inverse problems of heat conduction in terms of tem-perature field parameters obtained as a result of a thermophysical experiment.
Review and analysis of methods of determining thermal diffusivity by nonstationary temperature fields is given in the monographs by V.M. Fokin and V.N. Chernyshova [1] and N.P. Zhukov and N.F. Mainikov [2]. The method of complex determination of thermophysical characteristics of building materials and products by a method of nondestructive testing is presented in [1]. Extensive reviews of methods for solving inverse problems of thermal conductivity are presented in [3,4]. In [5], solutions to inverse problems are given, including an estimate of the temporal change in the interfacial conductivity between the casting and mold during aluminum solidification. Monographs [6,7] describe a technique to determine thermophysical properties of materials of metallurgi-cal production [6] and highly porous building materials [7].
Analytical solutions to the inverse problem of heat сonduction for a semi-bounded body are given in [8,9].
In [21][22][23][24][25], in order to determine thermophysical characteristics from a known temperature field, the use of a rather simple (engineering) method of numerical and analytical modeling of heat conduction processes described in [23][24][25] is proposed. In this method, analytical solutions of the differential equation of heat conduction in the form of algebraic expressions obtained for calculated time interval are used.
In this study, we propose a relatively simple method for determining thermal diffusivity a hc from the dynamics of temperatures T(R p , τ) and T(0, τ) of two surfaces of a plate with thickness R p . This approach will significantly simplify the preparation and experimental implementation of heating the body.
The values of a hc are determined by solving the inverse problem of heat conductivity according to the mathematical model of a semi-bounded body in the form of an unbounded plate, which can be treated as a semi-bounded body in the initial stage of heating when the Fourier number is less than Fo ≈ 0.06.
To obtain adiabatic conditions on the surface of the working plate, T(x = 0, τ), it is proposed to use a sec-Δτ ond adjacent plate with a greater thermal mass, for example, a plate made of the same material, but with a greater thickness. In the conjugation plane, in which the second thermocouple is installed (Fig. 1), the temperature will change only at moment τ ε , when the temperature disturbance in the working plate reaches conjugation plane x = 0. The temperature perturbation from the second plate will reach the conjugation plane later, at τ > τ e .
Measurements of temperatures T(R p , τ) and T(x = 0, τ) are performed while . Only one value a hc (T ≈ T b ) can be calculated in the case of а single heating under initial condition T(x p , τ = 0) = T b (0 ≤ x p ≤ R p ). To obtain tabulated dependence a hc (T), several experiments with different initial conditions are required.
Let us consider a mathematical description of the temperature field of an unlimited plate with thickness R p , heated only from the side x p = R p . We will assume that at moments t i , at which surface temperature T(R p , τ i ) = T 1,i is measured, the temperature distribution in the heated layer of thickness R is described by power-like function (1) We write the differential equation of thermal conductivity only for the heated layer of plate (2) with initial (3) and boundary conditions (4) Here c is the bulk heat capacity, J/(M 3 K); x is the coordinate, m; X is the relative coordinate, X = ; λ is the heat conductivity coefficient, W/(m K); q is the specific heat flux onto the surface X = 1 (x = R, Note that within an unheated layer 0 < x p < R p , while τ < τ e , the temperature will be equal to the initial one, therefore Figure 2 plots the distribution of temperatures over the thickness of the plate for several points in time, calculated using Eq. (1) for the calculation example given at the end of this study.
In general, the heat capacity, thermal conductivity, and specific heat flux can vary with time, depending The layout of the plate thicknes R p , which measures the temperature of the faces ( and ) and the plate to create adiabatic conditions on the surface of the x p = 0.
and heat transfer conditions.
Let the temperatures of the surfaces of the plate be known from the data of the thermophysical experiment: = = = and = = = , correspondingly, for X = 1 and X = 0 at moments τ i (i = 0, 1, 2, …).
When describing the temperature distribution over the thickness of the heated layer at the end of calculation time interval = using function (1) [21][22][23][24], boundary condition (5) can be written in form (6) From this we obtain the expression for determining q (7) Let us write down the equation for the heat balance of the heated layer of the plate (0 ≤ x ≤ R, 0 ≤ X ≤ 1) for ith time interval Δτ: (8) where T av,b and T av are the average-mass temperatures of heated layer R at the beginning (T av,b ) and at the end (T av ) of calculation time interval Δτ (hereinafter, the index of the time moment number is omitted for brevity).
For temperature distribution (1), T av at the end of interval Δτ is found by integration: (9) We transform heat balance equation (8) by substituting into it expressions (7) and (9) for q and T av : (10) Considering that thermal diffusivity a hc is equal to the ratio of λ and c , from Eq. (10) we obtain an expression for determining thermal diffusivity Coefficients a 0 and a 1 in formula (11) are easy to find according to temperatures and , known from the experiment, solving the system of two equations (1) written for X = 0 and X = 1 T X dX a al X dX a n λΔτ + − = +

a T a T T
As long as the plate can be treated as a semibounded body, temperature = = T b remains unchanged (τ < τ e ), since it is located at the boundary between the heated and unheated layers of the plate. Then, coefficients a 0 and a 1 are given by (13) where τ e is the moment of time when temperature T 0 (τ e ) exceeds T b by ΔT.
The value of τ e is determined by condition (14) where ΔT is the sensitivity or measurement error = T 0 . In the course of the experiment, one should only find the moment when = T 0 increases by an increment exceeding ΔT (see Fig. 2).
Thus, a relatively simple expression (11) was obtained for determining a hc for calculation time interval Δτ.
The formula for calculating a hc still contains unknown parameters: -n (exponent in Eqs. (1) and (11)); -R: the thickness of a heated layer depending on Fo, which, in turn, depends on a hc , R, and τ: The main objective of the study was to establish the relationship between parameters n, R, and Fo and function T 1 (τ) known from the experiment.
To obtain an acceptable solution to the inverse problem, various ways of determining n, R, and Fo were considered and tested.
. a R Possible methods for calculating n from the dynamics of the temperature field were considered in [23,25]. The values of n can be determined as follows [23,25]: -using function n(Fo); -using three temperatures T(X, τ), X= 0, Z, and 1 (0 < Z < 1); -using the solution to the differential equation of heat conduction for the plane of the plate, for X = 1; -using two solutions to the differential equation of heat conduction for X = 1 and 0 ≤ X ≤ 1 [23,25].
For the accepted conditions, function n(Fo) approximated by two expressions turned out to be more suitable [23,25]: Function (16) was determined by two points n(Fo 1 = 0.025) and n(Fo 2 = 0.05). Determination of regularity of the change in R(τ i ) is possible, assuming that number Fo and thermal diffusivity a hc are constant for all times τ i . Then, using the relations we obtain (18) Since plate thickness R p is known, and moment of time τ e is determined by condition -T b = ΔT (ΔT ≈ 0.1 K), then, using formula (18), it is possible to calculate the thickness of heated layer R i for each moment of time τ i in which temperature T 1 (τ i ) was measured.
Thermal conductivity coefficient λ hc can be calculated using transformed formula (7) (19) but only for those experiments in which specific heat flux q will be measured.
We consider the procedure for determining thermal diffusivity a hc using the example of the calculation of a hc from the surface temperatures of the plate, R p = 0.05 m, heated from T b = 300 K, if τ i = 0, 26, 52,… s, T 1 (τ i ) = 300, 307.1, 309.5, … K.
(2) We calculate thickness R 1 of the heated layer for τ 1 = 26 s using Eq.   (4) The system of three nonlinear equations is solved by iterative selection of number Fo Dependence n(Fo) was described by formula (16), since the iterative process using (17) was more complicated.
To find Fo with an accuracy of 0.0001 at each time interval, the "Service"/"Parameter selection" function of Microsoft Excel was used. Solving the system of equations, we obtained Fo = 0.0516, and then we find parameters n = 3.934 and a hc = 4.954 × 10 -7 m 2 /s. (5) We calculate the average-mass temperature of heated layer R = 0.0158 using Eq. (19) the average-mass temperature of the entire plate the thickness of the heated layer at the end of the second interval and the average-mass temperature at the start of the second interval of time Parameter T av,b is optional, it can be used to control the correctness of the calculation, since at τ = τ e , equality T av,b = T av,e must be fulfilled.

T R T R R T R
To assess the accuracy of the calculation method of thermal diffusivity a hc by solving the inverse problem of thermal conductivity, the method was tested. For the procedure of testing, we used a precalculated temperature field of a refractory material with specified thermophysical properties a i , c, and λ i : Calculation of the temperature field of an unbounded plate, R p = 0.05 m, was performed at T b = 300 K, gas temperature T g = 350 K, and coefficients of radiative and convective heat transfer σ = 4 × 10 -8 W/(m 2 K 4 ) and α = 30 W/(m 2 K) using the TRT dialog program [23] for modeling temperature fields.
It was necessary to determine a hc by solving the inverse problem of thermal conductivity, then comparing it with the initial values a i (23), at which the initial (test) temperature field was calculated.
The calculation results T(X = 1) = T 1 (τ i ), T(X = 0) = T(x p = 0) = T 0 (τ i ), and are given in Table 1. The results of determining a hc by the proposed method are shown in the rows of Table 1.

Temperature
= changed by 0.1 K at the moment of time τ e = 260 s (see Table 1). Coefficients a 0 and a 1 were determined by Eq. whereas the values of R i were calculated using the value τ e = 260 s and τ i from (18). The value of ΔFo j was determined as the absolute difference between the jth Fo j and the subsequent approximation of Fo j + 1 . Initial and calculated by (20) values of T av,p were almost the same, hence function (1) describes the temperature distribution fairly accurately along the thickness of the heated layer of the plate (see Fig. 2).
An arithmetic mean of ten values of a hc was attributed to the mass-average temperature of the heated layer for 260 s: T a = (300 + 303.9)/2 = 302 K. The corresponding value is a hc (T = 302) = 5.25 × 10 -7 . The deviation of δ from initial value a i (T = 302) = 5.395 × 10 -7 , calculated using Eq. (23) is δ a,av = 2.8%. To obtain a hc at other temperatures, experiments are required under different initial and boundary conditions.
Since the method presented is approximate, it is recommended to work out the experimental technique (boundary conditions) on materials with a known thermal diffusivity.

CONCLUSIONS
We proposed a mathematical description of the relationship between the thermal diffusivity of material Δτ on heating time interval a hc with surface temperatures and of plates known from a physical experiment, and a stepwise numericalanalytical algorithm for calculating a hc .
To assess the accuracy of the method, calculations of a hc were performed for ten points in time according to  the previously calculated (test) temperature field of an unbounded plate as a semi-bounded body (Fo < 0.6) at T b = 300 K, gas temperature T g = 350 K, under the condition of radiative-convective heat transfer. It is shown that for the accepted conditions of heat transfer, the average value of thermal diffusivity a hc (T = 302) was determined with an error of 2.8%. The application of the method significantly simplified conducting and processing experiments to determine thermophysical characteristics of heat-insulating and refractory materials.