Issue 
Mechanics & Industry
Volume 20, Number 6, 2019



Article Number  603  
Number of page(s)  8  
DOI  https://doi.org/10.1051/meca/2019026  
Published online  22 August 2019 
Regular Article
Numerical identification of the thermal conductivity tensor and the heat capacity per unit volume of an anisotropic material
^{1}
Université de Lomé, Faculté des Sciences, Département de Physique,
BP 1515,
Lomé, Togo
^{2}
Institut Pprime, UPR 3346 CNRS, Université de Poitiers, ENSMA SP2MI,
BP 30179,
86962
Futuroscope Chasseneuil Cedex,
France
^{*} email: katchonouglo@univlome.tg
Received:
19
March
2018
Accepted:
22
March
2019
In this paper, an inverse approach based on gradient conjugate method for thermal conductivity tensor and heat capacity per unit volume measurement is summarized. A suitable analysis is done for the mesh in finite element method and for the time steps for the time integration. For a composite material, it is shown the importance to identify the thermal conductivity tensor components in the principal axes.
Key words: Heat transfer / finite elements / thermal conductivity tensor / heat capacity / inverse problem
© AFM, EDP Sciences 2019
1 Introduction
When a material is subjected to mechanical solicitations, its temperature undergoes a few variation. The deformation processes induce thermal and energy effects which must be correctly predicted to properly model the mechanical behavior. For quasistatic tests, the temperature changes induced by deformation are often weak. The neglect of these variations eliminates any possibility of a complete energy balance because these small temperature variations induce large changes in energy. For more information on the mechanical response, it is necessary to perform an energy balance based on the thermodynamics of irreversible processes. This energy balance requires the knowledge of the thermal conductivity tensor k and heat capacity per unit volume c involved in the heat conduction [1–3].
Many numerical and experimental methods have been developed to determine the thermal conductivity tensor and sometimes the thermal conductivity tensor and the heat capacity per unit volume simultaneously. Inverse estimation of thermal parameters from temperature measurement has been studied in recent years [4–6]. In this paper, we propose a method using the same device than the mechanical study i.e. infrared camera and simple setup. Several materials have directiondependent thermal conductivities like for example composite materials [7,8]. Such kind of materials is denoted anisotropic, as in opposition to isotropic materials in which the thermal conductivity does not vary with the direction [9]. We are looking for the anisotropic homogeneous materials.
The evaluation of the thermal parameters as thermal conductivity tensor k and heat capacitance ρc is necessary in the characterization of new materials. They appear in the characterization of new polymers and in the analysis of thermal or heat insulation of materials [10–12].
For the estimation of the thermal properties of orthotropic materials, inverse analysis techniques based on Levenberg–Marquardt method or conjugate gradient method have been devoted [13–15]. Andersson and Ross [16] have used the transient hotwire method to measure simultaneously the thermal conductivity and heat capacity per unit volume; de Carvalho and Neto [17] combined the hotwire method and Levenberg–Marquardt method for polymer thermal properties estimation.
In this work, an identification process is developed to determine heat capacity and thermal conductivity tensor. The Fourier’s partial differential equation governing heat conduction as a differential equation is established by using the finite element method (FEM). That is well adapted to solve the problem by numerical computation. Although the FEM is well known, some terms will be developed to be included within the inverse approach. The identification of the unknown parameters k and ρc is performed from the minimization of a least square functional.
This study is a continuation of our work developed in the onedimensional (1D) case [18] where we have demonstrated the efficiency of this approach on experimental tests. The twodimensional (2D) extension is presented in this paper using a similar experimental device. The accuracy and the robustness of the inverse approach are assessed by the analysis of simulated temperature fields obtained from the 2D direct heat conduction problem. For that, different cases have been studied with several noise levels added to the temperature fields.
2 Direct problem
2.1 Modeling the thermal heat conduction in 2D
The physicalproblem consists here to solve the linear heat conduction in the case of an orthotropic material, with k a 2D tensor of the thermal conductivity and ρc a scalar value of the heat capacity per unit volume. The domain Ω is a solid initially at zero temperature. For times t > 0, two uniform heat fluxes q_{1} and q_{2} are respectively supplied at the surfaces x = 0 and y = 0 while the others are supposed insulated Figure 1. We assume that the plate is a composite with a small thickness and the length is very close or equal to the width, we can legitimately consider that the problem is 2D, where the temperature distribution within the plate is given by the function T(x, y, t).
The mathematical formulation of such physical problem is given by the Fourier’s equation (1)
with the boundary and initial conditions:
where n_{x} and n_{y} are respectively the normal at the surfaces x = 0 and y = 0 and k:
In the direct problem associated with the physical problem described above, the three thermal conductivity components k_{11}, k_{22} and k_{12} = k_{21}, the heat capacity per unit volume, as well as initial and boundary conditions, are known. The objective of the direct problem is then to determine the transient temperaturefield T(x, y, t) on the plate.
Several wellestablishment approaches could be used to solve the direct problem. We have solved this problem here using both a numerical FEM and analytic integral transform technique.
Fig. 1 Geometry of the plate. 
2.2 Finite element modeling of heat conduction in 2D
The temperature fields are computed simultaneously using a finite element formulation [19]. The region of interest is divided into an n^{e} number of threenoded elements Ω^{e}, with boundary ∂Ω^{e} with shape functions N^{e}(x, y) associated with each element e. In the isoparametric elements, the shape functions are used to transform the coordinates [19,20], and this enables better representation of any curved boundaries which may be present in the problem.
Particularly the FEM consists in a trivial choice of the decomposition functions as linear piecewise functions associated to the node i of the element e. In this paper, we choose a triangular element with linear transformation functions:
Let n_{t} and n^{e} be respectively the number of nodes and the number of elements that contains the mesh. We discretize the heat conduction equations within each element as (3)
where θ^{e} is the nodal temperature vector of the element.
The temperature function T(x, y, t) is approximated in the solution domain at any time by the vector θ(t) with (4)
where A^{e} is n^{e} × n_{t} matrix with each term is equal to 1 if the node numbered (globally) j on the structure coincides with the node numbered (locally) i on the element e, and 0 otherwise and where θ_{i} (t) is the timedependent value at node i. For each element domain Ω^{e}, the substitution of the equation (3) into the equation (1) and the application of the Galerkin method [20] results in the equation
where K^{e}, C^{e} and f^{e} are respectively the conductance matrix, the capacitance matrix and the heat load vector per element:
and the exponent T denotes the transpose of a matrix or a vector.
Let us denote by B^{e} the elementary gradient matrix:
Then the conductance elementary matrix becomes
The assemblage system equation takes the matrix form (5)
where the conductance matrix, the capacitance matrix and the heat load vector are, respectively:
The initial condition for solving the differential equation (5) follows from the equation (2d) (6)
3 Inverse problem
Let us consider the ordinary differential equation (7)
Concerning the inverse problem, the unknowns to be determined are the thermal conductivity parameters k_{11}, k_{22} and k_{12}, and the heat capacity per unit volume ρc. The other values are known.
Least squares functional
Let us assume that the temperature fields are recorded during the time t_{f}.
Consider also the following sum as squares functional (8)
where ∥.∥ is the norm associated to the scalar product in .
Our inverse problem aims the minimization of the functional under the constraint k is symmetric positive definite.
Minimization of the least squares norm (8) is done by using a Polak–Ribire conjugate gradient algorithm [21]. The main principle of the technique consists in calculating at each iteration a size and a direction of the step minimizing the objective function. The step direction is a linear combination of the descent directions of the gradient vector provided by the previous iterations.
The conjugate gradient method with an appropriate stopping criterion belongs to the class of iterative regularization techniques, in which the number of iterations is chosen so that stable solutions are obtained for the inverse problem [22–24].
For this method, we need the gradient expression . Due to the analytical expression of the solution of the direct problem given by the equation (7), we can also obtain the analytic gradient expressions for the functional .
To calculate these gradients, we calculate the variations of the functional respect to each variable k and ρc.
The variation respect to k reads:
where .
The variation of the functional respect to ρc leads to:
Since δk is symmetric, the gradient of the function versus k is equal to the matrix zeros implies:
The thermal conductivity symmetric tensor k and heat capacity per unit volume ρc are evaluated by solving the linear system:
4 Results and discussions
4.1 Temperature field simulation
Two programs have been developed for the direct and the inverse problems.
Let us consider a composite material constituted by oriented epoxy fibers in methyl methacrylate matrix. The preferential orientation of the fibers leads to an orthotropic material given two conductivity values. If fiber orientations are along the axis of the analysis, the thermal conductivity tensor is diagonal.
The dimensions of the sample are: l = L = 36.6 mm.
The diagonal thermal conductivity tensor components are:
The heat capacity per unit volume is
The initial temperature of the medium is defined by T_{0} =20 °C and the applied heat is ux q = 396.99 W/m^{2}.
We integrate according to the time the ordinary differential equation (5) with the initial condition equation (6) to simulate the temperature field for the estimation of thermal parameters.
A random white noise ωσ was added on the temperature field to mimic experimental conditions:
where θ_{i} is determined from the solution of the direct problem. σ gives the noise range and the randomnumber ω is generated in the interval [−1, 1].
4.2 Identification of the thermal parameters
In this section, the estimation of the thermal conductivity tensor k and the heat capacity per unit volume ρc is carried out using the proposed numerical approach given in equations (9).
The accuracy and robustness of the proposed identification method are assessed by the help of numerical simulations using the direct problem (5) and (6). The relative error ɛ is then obtained by comparing the calculated result α_{est} est with the imposed α_{imp}:
The estimation of the accuracy is done according to the spatial steps of the temperature fields:
where and are the number of elements according to respectively x and y directions and nb is the number of images recorded during the test.
4.2.1 Analysis with a fiber orientation along the axis analysis
The composite material and its dimensions are defined in previous section. First, we consider an analyze in the principal axis of the material then the thermal conductivity tensor k components are still:
and the cross terms of the thermal conductivity tensor k_{12} and k_{21} are null.
We show the relative errors for given dimensions of sample on (Fig. 2) and one experience duration with several spatial steps and time steps.
Without noise on the simulated temperature, we have a good agreement on the estimated parameters about 10^{−4} (Fig. 2) and accuracy seems to be insensitive to spatial and temporal parameters except for large number of elements. Indeed, in this case, the relative errors are mainly due to the numerical calculation, the same equations and same parameters (time and spatial steps) are used in direct and inverse problems. Then, this procedure can be used only to validate the implementations and to give the ultimate error of the inverse algorithm versus time and spatial steps.
The Figure 3, obtained for a noise level of σ = 0.01 °C, shows that the noise involves a direct effect on the results. Error values are multiplied by 10. We observe an evolution of the accuracy versus the parameters. Results can be improved by increasing either the spatial steps or the time steps. We can also see that the thermal conductivitytensor components depend mainly on the spatial steps while the heat capacity per volume unit depends on the time steps.
Let us analyze now the influence of the noise level on the thermal estimated parameters for these steps. Three noise levels have been studied (σ = 0.00, σ =0.01, σ =0.03 and σ =0.06) which represent the common noise levels obtained with current infrared cameras (Fig. 4). With noisy images, a filtering of the thermal field is necessary to improve the accuracy. We have tested several procedures, the best one is a spatial median filter (3 × 3 pixels) coupled with an average temporal filter (on 3 images). The efficiency of the spatial median filter is linked to a high spatial resolution that is why we have studied the noise effects for dx = 0.032 mm. This choice involves a time step of dt = 4 s to have a good accuracy.
We can observe the good precision on the estimated parameters on the filtered temperature field. For the high standard deviation σ = 0.06, the relative error on ρc is 0.054, k_{11} is 0.018 and on k_{22} is 0.082.
Fig. 2 Relative errors on the estimated parameters without noise on temperature fields. 
Fig. 3 Relative errors on the estimated parameters with noise on temperature fields: σ = 0.01 °C. 
Fig. 4 Relative errors on the estimated parameters versus the noise standard deviation. 
4.2.2 Analysis with a fiber orientation out of the axis analysis
We choose again the same composite material epoxymethyl methacrylate. The fiber orientation versus xaxis of the analysis is 30°. Then the thermal conductivity tensor components become:
The error of diagonal components is similar to the one obtained in the previous section (Fig. 5). Nevertheless, we note a low uncertainty for the thermal conductivity tensor component k_{12}.
To avoid this phenomenon, we propose to determine the thermal parameters along the principal axes of the composite material. This is possible because for this type of composite material, the angle between each fiber can be measured or defined during its elaboration.
Here, Ox and Oy are the axes of the Cartesian coordinate system, Oξ and Oη are the principal axes of the thermal conductivity tensor and ϕ is the angle between Ox and Oξ.
We denote the thermal conductivity tensor k′ in the principal axes. The relation between k′ and k is:
This expression is included in equation (8) to solve the inverse problem. Expression (8) is used again to obtain k′. The relative errors given by this new process are shown in Figure 6.
This new approach slightly disrupts the diagonal components of the thermal conductivity tensor and the heat capacity per unit volume, but it substantially improves the components k_{12} = k_{21}.
Fig. 5 Relative errors on the thermal conductivity coefficients versus the noise standard deviation. 
Fig. 6 Relative errors on the thermal conductivity coefficients versus the noise standard deviation. 
5 Conclusion
The FEM has been applied to reduce the partial differential equations of heat conduction to an ordinary differential equation system. Numerical integration of this equation system led to the simulation of the temperature field in two dimensions.
The identification algorithm based on the conjugate gradient projected allowed us to identify with a very good accuracy the coefficient of specific heat per unit volume and the components of the tensor of thermal conductivity of homogeneous orthotropic materials. This algorithm provides a simple, rapid and systematic identification of thermal parameters of materials from temperaturefield measurements.
To identify thermal parameters on orthotropic material, the uncertainties are better if the procedure is used according to principal axis of the material.
We have studied a wellknown epoxycarbon material, other types of more or less conductive material can be studied by the same procedure. This work allows us to study the thermomechanical behavior of composite to establish an energy balance.
Acknowledgements
This work was supported in part through the Agence Universitaire de la Francophonie (BAO201106U52210FT108) and the PoitouCharentes region (RPCR022).
References
 A. Chrysochoos, Infrared thermography, a potential tool for analyzing the material behaviour, Mcanique & Industries 3, 3–14 (2002) [CrossRef] [Google Scholar]
 A. Germaneau, J.C. Dupré, Thermal exchanges and thermomechanical couplings in amorphous polymers, Polym. Polym. Compos. 16, 9–17 (2008) [Google Scholar]
 A. Benaarbia, A. Chrysochoos, G. Robert, Thermomechanical analysis of the onset of strain concentration zones in wet polyamide 6.6 subjected to cyclic loading, Mech. Mater. 99, 9–25 (2016) [CrossRef] [Google Scholar]
 T. Telejko, Z. Malinowski, Application of an inverse solution to the thermal conductivity identification using the finite element method, J. Mater. Process. Technol. 146, 145–155 (2004) [CrossRef] [Google Scholar]
 C. Balaji, S.P. Venkateshan, A. Ambirajan, V.Ramakrishnan, Simultaneous estimation of principal thermal conductivities of an anisotropic composite medium, an inverse analysis, J. Heat Transfer 135, 021301 (2014) [Google Scholar]
 X. Wei, W. Chen, B. Chen, L. Suna, Singular boundary method for heat conduction problems with certain spatially varying conductivity, J. Comput. Math. Appl. 69, 206–222 (2015) [CrossRef] [Google Scholar]
 Ma.H. Zeng, J. Realff, M.L. Kumar, S.D.A. Schiraldi, Processing, structure, and properties of fibers from polyester/carbon nanofiber composites, Compos. Sci. Technol. 63, 1617–1628 (2003) [CrossRef] [Google Scholar]
 D. Hull, T.W. Clyne, An introduction to composite materials, Cambridge University Press, Cambridge, 1996 [CrossRef] [Google Scholar]
 M.M. Mejias, H.R.B. Orlande, M.N. Ozisik, A comparison of different parameter estimation technique for the identification of thermal conductivity components of orthotropic solids, in: 3rd International Conference on Inverse Problems in Engineering, Theory and Practice, Port Ludlow, WA, USA, 1999 [Google Scholar]
 V. Plana, P. Reulet, P. Millan, Experimental characterization of the thermophysical properties of composite materials by an inverse heat conduction method, J. Compos. Mater. 40, 1247–1258 (2006) [CrossRef] [Google Scholar]
 R.D. Sweeting, X.L. Liu, Measurement of thermal conductivity for fibrereinforced composites, Composites A 35, 933–938 (2004) [CrossRef] [Google Scholar]
 L. Czajkowski, W. Olek, J. Weres, R. Guzenda, Thermal properties of woodbased panels: thermal conductivity identification with inverse modelling, Eur. J. Wood Wood Prod. 74, 577–584 (2016) [CrossRef] [Google Scholar]
 B.W. Sawaf, M.N. Ozisik, Determining the constant thermal conductivities of orthotropic materials, Int. Commun. Heat Mass Transfer 22, 201–211 (1995) [CrossRef] [Google Scholar]
 B.W. Sawaf, M.N. Ozisik, Y. Jarny, An inverse analysis to estimate linearly temperature dependent thermal conductivity components and heat capacity of an orthotropic medium, Int. J. Heat Mass Transfer 38, 3005–3010 (1995) [CrossRef] [Google Scholar]
 R. Taktak, J.V. Beck, E.P. Scott, Optimal experimental design for estimating thermal properties of composite materials, Int. J. Heat Mass Transfer 36, 2977–2986 (1993) [CrossRef] [Google Scholar]
 S.P. Andersson, R.G. Ross, Thermal conductivity and heat capacity per unit volume of poly(methyl methacrylate) under high pressure, Int. J. Thermophys. 15 949–946 (1994) [CrossRef] [Google Scholar]
 G. de Carvalho, A.J.S. Neto, An inverse analysis for polymers thermal properties estimation, in: 3rd International Conference on Inverse Problems in Engineering, Port Ludlow, 1999 [Google Scholar]
 K. Atchonouglo, M. Banna, C. Vallée, J.C. Dupré, Inverse transient heat conduction problems and identification of thermal parameters, Heat Mass Transfer 45, 23–29 (2008) [CrossRef] [Google Scholar]
 J.M. Bergheau, R. Fortunier, Simulation numérique des transferts thermiques par éléments finis, Lavoisier, 2004 [Google Scholar]
 R.W. Lewis, K. Morgan, H.R. Thomas, K.N. Seetharamu, The finite element mdethods in heat transfer analysis, John Wiley, New York, 1996 [Google Scholar]
 E. Polak, G. Ribiere, Note sur la convergence de Méthodes de directions conjuguées, Rev. Franç. Informat. Rech. Opérationnelle 16, 35–43 (1969) [Google Scholar]
 O. Huajiang, Criteria for the finite element algorithm of generalized heat conduction equation, Appl. Math. Mech. 13, 587–596 (1992) [CrossRef] [Google Scholar]
 O.M. Alifanov, Inverse heat transfer problems, SpringerVerlag, Berlin, 1994 [CrossRef] [Google Scholar]
 M.N. Ozisik, H.R.B. Orlande, Inverse heat transfer: fundamentals and applications, Taylor & Francis, Pennsylvania, 1999 [Google Scholar]
 S.H. Carslaw, J.C. Jaeger, Conduction of heat in solids, 2nd edn., Oxford University Press, Oxford, 1959, p. 112 [Google Scholar]
Cite this article as: K. Atchonouglo, J.C. Dupré, A. Germaneau and C. Vallée Numerical identification of the thermal conductivity tensor and the heat capacity per unit volume of an anisotropic material, Mechanics & Industry 20, 603 (2019)
All Figures
Fig. 1 Geometry of the plate. 

In the text 
Fig. 2 Relative errors on the estimated parameters without noise on temperature fields. 

In the text 
Fig. 3 Relative errors on the estimated parameters with noise on temperature fields: σ = 0.01 °C. 

In the text 
Fig. 4 Relative errors on the estimated parameters versus the noise standard deviation. 

In the text 
Fig. 5 Relative errors on the thermal conductivity coefficients versus the noise standard deviation. 

In the text 
Fig. 6 Relative errors on the thermal conductivity coefficients versus the noise standard deviation. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.