Identification of the Thermal Conductance of a hidden barrier from outer thermal data

Hidden defects affecting the interface in a composite slab are evaluated from thermal data collected on the upper side of the specimen. First we restrict the problem to the upper component of the object. Then we investigate heat transfer through the inaccessible interface by means of Thin Plate Approximation. Finally, Fast Fourier Transform is used to filter data. In this way we obtain a reliable reconstruction of simulated flaws in thermal contact conductance corresponding to appreciable defects of the interface.


Introduction
Consider a composite body made up of two slabs of different materials in close thermal contact. Since the contact surfaces are rough on a microscopic level, thermal contact is always imperfect (see for example [1]-Section 3.1 and [2]). With reference to figure 1, taking into account that the true contact area is a small portion of the apparent contact area, the slabs are separated by an interface of average width d, filled up with air. The present work deals with Figure 1: Contact between two real solid surfaces the nondestructive evaluation of deviations of the interface width from a given average value. Since κair d defines the thermal conductance between the slabs, a thermal model of the composite material can be implemented to evaluate local variations of the width by applying a controlled heat flux and collecting a sequence of temperature maps on the top side of the body (Active Thermography A recent effective approach to the solution of this problem, based on Reciprocity Functional technique, is described in [4]. Here, a perturbative tool like Thin Plate Approximation (TPA)(see for example [5]) is used in alternative (see section 5). The idea of using TPA comes from the general assumption (see [3] Sect 9.2.1) that thermography is effective in detection of subsurface anomalies. Here, TPA is expected to work since the Biot number of the upper slab is ≤ .1 so that it can be regarded as a thermally thin domain (see [1] Section 5.2).
We stress that Ω + and Ω − are made of different materials each characterized by density ρ ± , specific heat c ± and thermal conductivity κ ± . The third, irregular, thin slab is and ρ a , c a and κ a are its physical parameters. Let γ = a L a dimensionless parameter which represents the geometrical "thinness" of the slab Ω + .
The domain A has variable thickness (g + (x, y) − g − (x, y)) and it is assumed to be filled up with air, whose conductivity κ a = 0.002587 Wm −1 K −1 is much lower than κ ± .
It means that A opposes to heat transfer from Ω + to Ω − and the corresponding Thermal Conductance is 3 Modeling the solid interface A by means of Robin boundary conditions on the two sides of a plane. Imperfect contact.
The temperature of Ω A depends on the physical characteristics of A. We account for the different conductivities κ − , κ a , κ + in Ω A by imposing continuity of temperature and heat flux for z = g + (x, y) and z = − g − (x, y) as transmission conditions for the heat conduction equation in Ω A . The functions g + and g − are very irregular at a microscopic scale but, if their values are normally distributed in a small neighborhood of the mean values g − and g + , the set A can be successfully approximated by the parallelepiped The task of solving the heat equation in Ω A becomes much simpler if we considerĀ instead of A. The interface behaves like a controlled heat exchanger between Ω + and Ω − of conductance H = κa (g + −g − ) . A further simplification consists in squeezingĀ to the plane z = 0 and assigning the Robin boundary conditions (imperfect contact [6]) The presence of an anomaly in the interface, corresponds to a change in the heat transfer from Ω + to Ω − . Hence, in analogy with (3), we have where These boundary conditions describe a non constant imperfect contact (see [7] [9]). The pointwise evaluation of the unknown function H(x, y) is the main goal of the present work.
Remark We are assuming that variations of H in time take place in a much longer interval than (0, t max ). For this reason, H does not depend on t.

The direct model
Let u + be the solution of Boundary data are u + x (−L, y, z, t) = u + x (L, y, z, t) = u + y (x, −L, z, t) = u + y (x, L, z, t) = 0 (8) and i.e. a source of power Φ is ON for t S seconds (t S < t max while the index S is for "source"). Let u − be the solution of holds and can be written in terms of boundary data on the two sides of the interface z = 0. More precisely We stress that solutions u + and u − depend on the dimensionless parameter γ = a L introduced in section 2.

From transmission conditions to Robin Boundary Conditions
It is well known that, when κ + = κ − , the midpoint of the thermal jump due to imperfect contact is the value of the background temperature . It means that transmission conditions can be changed in ordinary boundary conditions on the two sides of the interface (see for example [10]). Here, it is for all H, x, y and t. When H = 0 (insulating interface) and H = ∞ (perfect contact) we have R = 0. R > 0 elsewhere. There is numerical evidence that R is negligible when parameters are in the range we are dealing with. Figure 2 shows the temperature gap across the interface for a heat exchange coefficient changing from 0 (blue curve) to ∞ (red curve), compared to that relative to H = 1000 (green curve). Moreover, in an analogous simpler problem for composite regions ([8] sect 14.6) in which the temperature is known explicitly, straightforward calculations lead to the following formula for "large" H: where G 0 (x, y, 0, t; 0, 0, a, 0) is the Green function of a Continuous Plane Source and the constant S is Assumed the smallness of R (confirmed by numerical computations and by the analytical estimate (13), interface boundary conditions can be written approximately as The IBVP (6)(7)(8)(9)(15) is now the direct model underlying the inverse problem of evaluating H. Observe that the problem has been restricted to the upper slab Ω + .

Dimensionless variables and Thin Plate Approximation
The dimensionless parameter γ = a L gives a measure of how much Ω + is "geometrically thin". As observed in the introduction, Ω + is considered "thermally thin" when aH κ < .1. We introduce the set of dimensionless variables ζ = 1− z γL , ξ = x L , η = y L and τ = t α+ L 2 (we recall that the numbers α ± = κ± ρ±c± are the diffusivities of upper and lower slabs respectively) and define so that, the heat equation for ζ > 0 (former upper slab) becomes Observe that in (17) we use a notation that points out the explicit dependence of v +,γ on γ. Variables ξ, η, ζ and τ are slightly different from usual ones (see for example [12] Sect 4.1) because we need to normalize ζ by means of the parameter γ which, consequently, is present in (18).
Thanks to the substitution of (11) with (15) , we can limit ourselves to the upper slab ζ ∈ (0, 1) . Furthermore, we scale suitably the heat flux

Formal expansion of v +,γ
Consider the formal expansions and It is remarkable that, though h is independent on τ (at least in the time scale at hand), its partial sums in the expansion above are dependent on τ by construction.
We observe that We claim that the limit in (25) exists as we can see by generalizing the one dimensional case with φ(x, y) = φ 0 and constant h. In this special case, it is for t ∈ (0, t S ).

Plugging the expansion of v +,γ in (18), (19) and (20) we have
with v 0ζ = v 1ζ = 0 i.e. ∂(v0+γv1) ∂ζ = 0. As suggested by (25) and (26), we have v 0 (ξ, η, ζ, τ ) ≡ v 0 (ξ, η, τ ) ≈ u + (x, y, a, t) (the symbol ≈ means that u + (x, y, a, t) is taken in practice on the top side of a thin plate of thiockness a > 0) and v 1 (ξ, η, ζ, τ ) = lim The order zero approximation of h can be obtained from second order relations. More precisely, Choosing a time valuet which corresponds to a good approximation of the "ideal" value of H (i.e. its value in absence of conduction flaws) and coming back to variables (x, y, z, t), we have where T = u(x, y, a,t) is the temperature of the accessible side z = a at timet. Recall that, in experimental real life situations, our knowledge of T come from a sequence of measurements taken by means of an infrared camera. In numerical simulations T is computed solving the direct model with a finite elements code.
Remark. When the unknown parameter H is a function of two variables, the computation of higher order term in expansion (24) is very hard. A strong smoothing procedure would be required in order to perform reliable computation of fourth or sixth order partial derivative of data T . At present we have not yet solved this problem.

Numerical computations
Simulation of experimental data collection requires the numerical solution of the direct problem described in section 4. Numerical values of parameters (in MKS units) are:

Upper slab
• κ + = 54 • ρ + = 7870 • c + = 486 Lower slab Geometrical parameters: • d 0 = 10 −5 m (average contact thickness) • a = 10 −2 m (thickness of a single slab) The upper surface is uniformly illuminated with a source (spotlight) having a power density • φ(x, y) = 10 5 Wm −2 constant in time. The contact resistance between the slabs is due to a variable thickness d(x, y). Given the thermal conductivity k a of the material between the two sheets (air) • k a = 25.87 × 10 −3 Wm −1 K −1 the unknown heat exchange coefficient is given by For the simulations we assume: where: Figure 3 shows the graph of H(x, y). Figure 4 represents a filled contour plot of H(x, y). Simulation has been coded in COMSOL Multiphysics. The model is three dimensional, with the imperfect contact between the two slabs modeled as a resistive thin layer. The continuous uniform heating is supplied for 100 seconds. Figure 5 shows the temperature distribution "recorded" after 50 seconds. A gaussian noise with σ = 0.1 o C has been added to the simulated temperature map.
6.1 Recovering H from thermal data collected on the top face of the specimen Figures 6 and 7 show the reconstruction of H(x, y) obtained by means of TPA expression (34) from exact and noisy data, respectively. In order to compute the spatial derivatives appearing in (34), a cubic smoothing spline is applied to data when computing the first-order derivatives. Moreover, the same kind of smoothing is applied to temperature maps when noise is present. The weight factors applied to obtain the result shown in figure 7 are somewhat arbitrary. Since H is actually not known, the choice of such factors does not appear to be a straightforward task.  An alternative to using cubic smoothing splines is to filter noisy data by means of Fast Fourier Transform (FFT). The following procedure is applied to a tem-perature map u i taken at time t i : 2. The spectrum is shifted to the center of the data, and scaled to have a maximum value 0: SP i = log(P i ) − max (log(P i )).
3. Looking at the graph of the scaled power spectrum (see, for example, figure 9), a suitable threshold S t is chosen to cut the image clutter. With such a threshold a mask M is computed (having 0 for SP i < S t ) and multiplied by F i .
4. The filtered temperature map at time t i is obtained by an inverse Fourier transform: The effective heat exchange coefficient computed on the Fourier-filtered map gives the result shown in figure 10.

Choice of the filter threshold
In order to obtain a smooth temperature map from the noisy data, a lowpass filter (LPF) is applied in the frequency domain. The LPF consists in applying a mask to the spectrum based on a suitable threshold θ. Denoting by U = F {u} the Fourier transform of the temperature u, and by S U = |U | 2 max(|U | 2 ) the normalized power spectrum, the mask M is defined by: A "denoised" map u D is obtained by inverse-transforming the product U · M : u D = F −1 (U · M ).The choice of the threshold is not a straightforward task: a small threshold leaves data unmodified, while a too large one smooths data too much. We adopt the following procedure, which allows an automatic computation of the threshold. If we choose a random value of the threshold, between min(S U ) and 0, an high-pass filter (HPF) is realized by the mask M c complementary to M . We shall call "noise" temperature that obtained by: Choosing θ = 0 the HPF gives a perfect replica of the map u, while the LPF gives a null power spectrum. Conversely, the LPF gives a perfect (unsmoothed) replica of the temperature data if θ = min(S U ) while the HPF gives a null spectrum in such conditions. A small value of θ > min(S U ) produces masks M and M c such to give an almost random temperature. Therefore, intuitively we should increase the value of θ, starting from the lower value of the power spectrum, until the noise temperature u N vaguely starts to look like the original map u. In those conditions, the map u D is over-smoothed.
The previous reasoning can be translated into formulae, by considering what "to look like" means in this context. If we compute the correlation among the noise (n) temperature map u N and the original map (signal, s) u, c n,s corr(u N , u), for values of θ ∈ (min(S U ), 0), the plot of c n,s versus θ is like that shown in figure 8, where a third-order polynomial fitting is superimposed to the computed correlation. The threshold θ 0 corresponding to the first minimum of the fitting curve is such to give the best low-pass filtered data, without "loosing" information (i.e. with a noise map u N actually representing noise).
Looking at figure 8 it is clear that such a minimum is not the absolute minimum of the correlation. Clearly, a zero noise temperature is the absolutely less correlated to any non-zero map! Figure 9 shows the spectrum corresponding to figure 8, actually a cut of the 3d-spectrum on a symmetry plane.

Resolution limits
Gaussian blurring is intrinsic in heat diffusion processes, so we can expect that different contact resistance profiles conduct to similar H shapes. That's indeed the case. Figure 11 shows a H profile having two deviations from the base value, respectively with rectangular and gaussian shapes.

Conclusions
Consider a composite body made up of two slabs of different materials in close imperfect thermal contact. Large deviations of the interface width from a known average value d 0 must be detected. The body is heated from above (with a spotlight) while temperature maps of the upper side of the body are collected (with an infrared camera). Since transmission conditions through the interface are transformed into Robin Boundary Conditions, we reduce the evaluation of flaws in the interface to an inverse problem for the heat equation in the upper slab only (section 4.1). The unknown is the heat transfer coefficient 2H(x, y) in the Robin condition where H is the thermal conductance of the interface. A reliable evaluation of H is obtained in section 6 by means of Thin Plate Approximation and Fast Fourier Transform with simulated data corresponding to a realistic physical situation.