Rear-surface integral method for calculating thermal diffusivity: finite pulse time correction and two-layer samples

TL;DR

This paper extends the rear-surface integral method for calculating thermal diffusivity from laser flash experiments, incorporating realistic heat pulse shapes and two-layer samples, resulting in more accurate and applicable formulas.

Contribution

It introduces new formulas for thermal diffusivity calculation that account for finite pulse shapes and layered samples, improving upon previous assumptions.

Findings

New formulas accurately estimate thermal diffusivity from realistic pulse shapes.

Method effectively handles two-layer heterogeneous samples.

Numerical experiments validate the improved accuracy of the formulas.

Abstract

We study methods for calculating the thermal diffusivity of solids from laser flash experiments. This experiment involves subjecting the front surface of a small sample of the material to a heat pulse and recording the resulting temperature rise on the opposite (rear) surface. Recently, a method was developed for calculating the thermal diffusivity from the rear-surface temperature rise, which was shown to produce improved estimates compared with the commonly used half-time approach. This so-called rear-surface integral method produced a formula for calculating the thermal diffusivity of homogeneous samples under the assumption that the heat pulse is instantaneously absorbed uniformly into a thin layer at the front surface. In this paper, we show how the rear-surface integral method can be applied to a more physically realistic heat flow model involving the actual heat pulse shape from the laser flash experiment. New thermal diffusivity formulas are derived for handling arbitrary pulse shapes for either a homogeneous sample or a heterogeneous sample comprising two layers of different materials. Presented numerical experiments confirm the accuracy of the new formulas and demonstrate how they can be applied to the kinds of experimental data arising from the laser flash experiment.

Tables

Table 1. Table 1: Heat flow parameter values used in all numerical experiments. In each case, the target value of the thermal diffusivity is displayed rounded to five significant figures.

Homogeneous (Single-layer) $k=222\text{W}{\text{m}}^{-1}{\text{K}}^{-1},\rho =2700\text{kg}{\text{m}}^{-3},c=896\text{J}{\text{kg}}^{-1}{\text{K}}^{-1},L=0.002\text{m},$ Target value: $\alpha =k/\left(\rho c\right)=9.1766\times {10}^{-5}{\text{m}}^2{\text{s}}^{-1}$. Heterogeneous (Two-layer) $k_1=222\text{W}{\text{m}}^{-1}{\text{K}}^{-1},{\rho }_1=2700\text{kg}{\text{m}}^{-3},c_1=896\text{J}{\text{kg}}^{-1}{\text{K}}^{-1},{\ell }_1=0.00176\text{m},$ $k_2=16.3\text{W}{\text{m}}^{-1}{\text{K}}^{-1},{\rho }_2=7810\text{kg}{\text{m}}^{-3},c_2=480\text{J}{\text{kg}}^{-1}{\text{K}}^{-1},{\ell }_2=0.00024\text{m},$ Target value: ${\alpha }_1=k_1/\left({\rho }_1{c}_1\right)=9.1766\times {10}^{-5}{\text{m}}^2{\text{s}}^{-1}$, Target value: ${\alpha }_2=k_2/\left({\rho }_2{c}_2\right)=4.3481\times {10}^{-6}{\text{m}}^2{\text{s}}^{-1}$.

Table 2. Table 2: Thermal diffusivity estimates, α ~ ~ 𝛼 \widetilde{\alpha} , α ~ 1 subscript ~ 𝛼 1 \widetilde{\alpha}_{1} and α ~ 2 subscript ~ 𝛼 2 \widetilde{\alpha}_{2} , and corresponding signed relative errors, ε 𝜀 \varepsilon , ε 1 subscript 𝜀 1 \varepsilon_{1} and ε 2 subscript 𝜀 2 \varepsilon_{2} , obtained from applying the formulas ( 17 )–( 19 ) and ( 44 )–( 47 ) to the parameter values in Table 1 and the synthetic rear-surface temperature rise history discussed in Section 4.1 . Results are given for the rectangular, triangular and exponential pulses given in Eqs ( 20 )–( 22 ).

	Rectangular	Triangular	Exponential
Homogeneous (single-layer)
$\tilde{\alpha }$ [${\text{m}}^2{\text{s}}^{-1}$]	$9.1766\text{×}{10}^{-05}$	$9.1766\text{×}{10}^{-05}$	$9.1766\text{×}{10}^{-05}$
$\epsilon$ [%]	$2.0077\text{×}{10}^{-04}$	$2.0078\text{×}{10}^{-04}$	$2.0078\text{×}{10}^{-04}$
Heterogeneous (two layer)
${\tilde{\alpha }}_1$ [${\text{m}}^2{\text{s}}^{-1}$]	$9.1766\text{×}{10}^{-05}$	$9.1766\text{×}{10}^{-05}$	$9.1766\text{×}{10}^{-05}$
${\epsilon }_1$ [%]	$1.5619\text{×}{10}^{-04}$	$1.8399\text{×}{10}^{-04}$	$1.8065\text{×}{10}^{-04}$
${\tilde{\alpha }}_2$ [${\text{m}}^2{\text{s}}^{-1}$]	$4.3480\text{×}{10}^{-06}$	$4.3480\text{×}{10}^{-06}$	$4.3480\text{×}{10}^{-06}$
${\epsilon }_2$ [%]	$2.0250\text{×}{10}^{-04}$	$2.3855\text{×}{10}^{-04}$	$2.3421\text{×}{10}^{-04}$

Equations

$$\alpha \approx \frac{1.37L^2}{{\pi }^2{t}_{0.5}},$$

$$\displaystyle\alpha=\frac{T_{\infty}\bigl{(}L^{2}-\ell^{2}\bigr{)}}{6\int_{0}^{\infty}\left[T_{\infty}-T(L,t)\right]\,\text{d}t}.$$

$$\frac{\partial T}{\partial t}(x,t)=\alpha \frac{{\partial }^{2}T}{\partial x^2}(x,t),0<x<L,t>0,$$

$$T(x,0)=0,0<x<L,$$

$$-k\frac{\partial T}{\partial x}(0,t)=q(t),\frac{\partial T}{\partial x}(L,t)=0,$$

$$Q_{\infty }={\int }_0^{\infty }q(s)\text{d}s<\infty ,$$

$$T_{\infty }=\underset{t\to \infty }{\lim }T(x,t)=\frac{Q_{\infty }}{\rho cL},0<x<L.$$

$$\rho c{\int }_0^{L}T(x,t)\text{d}x=Q(t),$$

$$Q(t)={\int }_0^{t}q(s)\text{d}s.$$

$$u(x)={\int }_0^{\infty }\left[T_{\infty }-T(x,t)\right]\text{d}t.$$

$$\frac{{\text{d}}^{2}u}{\text{d}{x}^2}=-\frac{T_{\infty }}{\alpha },0<x<L,$$

$$u(x)=c_0+c_{1}x-\frac{T_{\infty }{x}^2}{2\alpha }.$$

$$\frac{\text{d}u}{\text{d}x}(0)$$

$$\frac{\text{d}u}{\text{d}x}(L)$$

$$\rho c{\int }_0^{L}u(x)\text{d}x={\int }_0^{\infty }\left[Q_{\infty }-Q(t)\right]\text{d}t.$$

$$u(x)=\frac{{\int }_0^{\infty }[Q_{\infty }-Q(t)]\text{d}t}{\rho cL}+\frac{T_{\infty }L}{\alpha }\left[x-\frac{L}{3}-\frac{x^2}{2L}\right].$$

$$\alpha =\frac{L^2}{6(I_T-I_q)},$$

$$I_T={\int }_0^{\infty }[1-T(L,t)/T_{\infty }]\text{d}t,$$

$$I_q={\int }_0^{\infty }[1-Q(t)/Q_{\infty }]\text{d}t.$$

$$q(t)$$

$$\alpha$$

$$q(t)$$

$$\alpha$$

$$q(t)$$

$$\alpha$$

$$\frac{\partial T_1}{\partial t}(x,t)={\alpha }_1\frac{{\partial }^2{T}_1}{\partial x^2}(x,t),0<x<{\ell }_1,t>0,$$

$$\frac{\partial T_2}{\partial t}(x,t)={\alpha }_2\frac{{\partial }^2{T}_2}{\partial x^2}(x,t),{\ell }_1<x<L,t>0,$$

$$T_1(x,0)=0,0<x<{\ell }_1,$$

$$T_2(x,0)=0,{\ell }_1<x<L,$$

$$-k_1\frac{\partial T_1}{\partial x}(0,t)=q(t),t>0,$$

$$\frac{\partial T_2}{\partial x}(L,t)=0,t>0,$$

$$T_1({\ell }_1,t)=T_2({\ell }_1,t),t>0,$$

$$k_1\frac{\partial T_1}{\partial x}({\ell }_1,t)=k_2\frac{\partial T_2}{\partial x}({\ell }_1,t),t>0,$$

$$T_{\infty }=\underset{t\to \infty }{\lim }{T}_1(x,t)=\underset{t\to \infty }{\lim }{T}_2(x,t)=\frac{Q_{\infty }}{{\rho }_1{c}_1{\ell }_1+{\rho }_2{c}_2{\ell }_2}.$$

$${\rho }_1{c}_1{\int }_0^{{\ell }_1}{T}_1(x,t)\text{d}x+{\rho }_2{c}_2{\int }_{{\ell }_1}^L{T}_2(x,t)\text{d}x=Q(t),$$

$$u_1(x)={\int }_0^{\infty }\left[T_{\infty }-T_1(x,t)\right]\text{d}t,$$

$$u_2(x)={\int }_0^{\infty }\left[T_{\infty }-T_2(x,t)\right]\text{d}t.$$

$$\frac{{\text{d}}^2{u}_1}{\text{d}{x}^2}=-\frac{T_{\infty }}{{\alpha }_1},0<x<{\ell }_1,$$

Full text

††journal: International Journal of Heat and Mass Transfer

1 Introduction

The most popular technique for measuring the thermal diffusivity of a solid material is the laser flash method [1, 2, 3]. Originally developed by Parker et al. [4], this method involves subjecting the front surface of a small sample of the material to a heat pulse of radiant energy and recording the resulting temperature rise on the opposite (rear) surface of the sample (Figure 1). The thermal diffusivity is then calculated from the rear-surface temperature rise curve by developing a mathematical model describing the heat flow emanating from the front surface and analysing the solution of this model at the rear surface.

Parker et al. [4] assumed the sample is homogeneous and thermally insulated; the heat pulse is instantaneously and uniformly absorbed by a thin layer at the front surface; and the heat flow is one-dimensional in the direction of the thickness of the sample extending from the front surface to the rear surface (Figure 1b). These simplifying assumptions yield the now famous and widely-used formula [5, 6, 7, 8] for the thermal diffusivity

[TABLE]

where $L$ is the thickness of the sample and $t_{0.5}$ is the half-rise time, the time required for the rear-surface temperature rise to reach one half of its maximum (steady-state) value $T_{\infty}$ (Figure 1d).

Since 1961, many papers have presented advances and extensions on Parker et al.’s original method. In particular, numerous papers have focussed on the finite pulse-time effect, that is, the observation that the heat pulse occurs over a finite duration [9, 10] and is rarely well approximated by an ideal instantaneous pulse (as assumed by Parker et al. [4]). In the one-dimensional case, the boundary condition imposing thermal insulation at the front surface [4, 11] is replaced by a boundary condition specifying the heat flux applied at the front surface based on the actual shape of the heat pulse [9] (Figure 1a). Commonly, simple functions are used to represent the pulse shape such as those describing rectangular [12], triangular [9, 12] and exponential [13] pulses. For the half-rise time approach, accounting for the finite-pulse time leads to small corrections [12, 9] to the formula for the thermal diffusivity (1).

Considerable research activity has also focussed on the application of the laser flash method to layered samples consisting of two or more adjacent homogeneous layers with distinct thermal diffusivities [13, 14, 15, 16]. For a two layer sample, only one of the thermal diffusivities can be calculated from the temperature rise recorded at the rear surface [17, 15]. In this case, assuming the thermal diffusivity of one of the layers is known, extension of the half-rise time approach yields a formula for the thermal diffusivity of the other layer involving the half-rise time and the volumetric heat capacity of both layers [13].

Recently, Carr [11] presented a method for calculating the thermal diffusivity from the rear-surface temperature rise history. This so-called rear-surface integral method expresses the thermal diffusivity exactly in terms of the area between the theoretical rear-surface temperature rise curve, $T(L,t)$, and the steady-state temperature, $T_{\infty}$, over time (see Figure 1d). For a homogeneous and thermally insulated sample, where the heat pulse is instantaneously and uniformly absorbed by a thin layer of thickness $\ell$ at the front surface, the thermal diffusivity formula is given by [11]:

[TABLE]

Numerical experiments carried out on synthetic data generated by adding Gaussian noise to the theoretical rear-surface temperature rise curve, $T(L,t)$, demonstrated that the formula (2) produces estimates of the thermal diffusivity that are more accurate than the standard formula (1).

For the case of an ideal instantaneous heat pulse ($\ell=0$), the formula (2) reduces to one derived by Baba [18] using their areal heat diffusion method. The areal heat diffusion method obtains a closed-form expression for the area $\int_{0}^{\infty}\left[T_{\infty}-T(L,t)\right]\,\text{d}t=\lim_{s\rightarrow 0}[s^{-1}T_{\infty}-\widetilde{T}_{r}(s)]$ involving $\widetilde{T}_{r}(s)$ the Laplace transform of $T_{r}(t):=T(L,t)$. Our rear-surface integral method for deriving (2) avoids explicit calculation of $T_{r}(t)$ and its Laplace transform. This means that it is easier to apply to more sophisticated heat flow models as considered in the current paper.

Carr [11] essentially provided a proof-of-concept for the rear-surface integral approach applied to the classical heat flow model of the laser flash experiment [4]. In the current paper, we present two key contributions. Firstly, we extend the rear-surface integral approach to account for finite pulse time effects deriving new formulas for the thermal diffusivity for rectangular, triangular, exponential and arbitrary pulse shapes (Section 2). Secondly, we extend these new theoretical results to two layer samples, developing new formulas for the thermal diffusivity of both layers (Section 3). After deriving the new formulas, we verify their accuracy and demonstrate their application to the types of experimental data that arise from the laser flash experiment (Section 4). We then conclude with a summary of the work and a discussion regarding the limitations of the new formulas (Section 5).

2 Homogeneous sample

2.1 Heat flow model

Consider a thermally insulated homogeneous sample of uniform thickness $L$, thermal conductivity $k$, volumetric heat capacity $\rho c$ and thermal diffusivity $\alpha=k/(\rho c)$ (Figure 1b). Let $\mathcal{T}(x,t)$ be the temperature of the sample at location $0\leq x\leq L$ and time $t\geq 0$ and assume the heat pulse applied at the front surface, $x=0$, is initiated at $t=0$. We assume the temperature rise, $T(x,t)=\mathcal{T}(x,t)-T_{0}$, above the initial uniform temperature of the sample, $T_{0}$, resulting from the heat pulse satisfies the following heat flow model [9, 12, 2]:

[TABLE]

where $q(t)$ is the heat flux function representing the laser pulse applied at the front surface of the sample (Figure 1a). Provided $\lim_{t\rightarrow\infty}q(t)=0$ and the total amount of heat absorbed at the front surface

[TABLE]

the steady-state temperature rise is given by [9]:

[TABLE]

This result is derived by first applying the thermally insulated boundary conditions (5) that arise in the limit $t\rightarrow\infty$ to the linear steady-state solution of the heat equation (3) yielding $\lim_{t\rightarrow\infty}T(x,t)=T_{\infty}$, where $T_{\infty}$ is a constant. The value of $T_{\infty}$ is then identified by noting that the change in the amount of heat in the sample must be balanced by the amount of heat entering through the front surface:

[TABLE]

where

[TABLE]

At steady state, Eq (8) reduces to $\rho cLT_{\infty}=Q_{\infty}$ yielding the stated result (7).

2.2 Thermal diffusivity

We now extend the rear-surface integral method for calculating thermal diffusivity [11] to non-instantaneous heat pulse durations. The aim of this method is to derive a formula for the thermal diffusivity, $\alpha$, in terms of an integral involving $T(L,t)$, the rear-surface temperature rise extracted from the heat flow model (3)–(5). The starting point for the derivation is the function:

[TABLE]

Differentiating $u(x)$ twice with respect to $x$, noting $T_{\infty}$ (7) is a constant and using the heat equation (3) and initial condition (4) yields the differential equation [11]:

[TABLE]

which identifies $u(x)$ as a quadratic function:

[TABLE]

To identify the integration constants, $c_{0}$ and $c_{1}$, we combine the form of $\text{d}u/\text{d}x$ with the boundary conditions (5) [11] yielding:

[TABLE]

Imposing both of these boundary conditions on $u(x)$ (12) yields the same result, namely $c_{1}=T_{\infty}L/\alpha$, after noting the relationship between $T_{\infty}$ and $Q_{\infty}$ (7). Therefore, to uniquely identify $u(x)$ we require an auxiliary condition analogous to the one required for the steady-state solution (8). Integrating $u(x)$ (10) from $x=0$ to $x=L$, multiplying by the volumetric heat capacity $\rho c$ and reversing the order of integration yields the appropriate condition:

[TABLE]

Imposing this auxiliary condition on $u(x)$ (12) uniquely identifies $c_{0}$ and hence $u(x)$:

[TABLE]

In summary, we have derived a closed-form expression for the integral (10) by formulating and solving the boundary value problem consisting of Eqs (11) and (13)–(15). Equating both expressions for $u(x)$, Eqs (10) and (16), at the rear surface, $x=L$, and rearranging yields a formula for the thermal diffusivity that can be expressed as:

[TABLE]

where

[TABLE]

Eqs (17)–(19) express the thermal diffusivity explicitly in terms of $L$, the thickness of the sample, $T(L,t)$, the rear-surface temperature rise history, $T_{\infty}$, the steady-state temperature rise value (7), $Q_{\infty}$, the total amount of heat absorbed through the front surface (6), and $Q(t)$, the amount of heat that has been absorbed through the front surface at time $t$ (9). For the types of functions commonly used to represent the shape of the heat pulse, which specify the finite pulse time $\tau$, and/or peak of the pulse $\beta$ (Figure 1a), the formula for the thermal diffusivity (17)–(19) simplifies significantly:

Rectangular pulse [12]:

[TABLE]

Triangular pulse [12, 9]:

[TABLE]

Exponential pulse [13]:

[TABLE]

In the limit $\tau\rightarrow 0$, the rectangular pulse tends to the instantaneous pulse $q(t)=Q_{\infty}\delta(t)$, where $\delta(t)$ is the Dirac delta function, yielding $\alpha=L^{2}/(6I_{T})$ as derived previously [11]. Hence, the above formulas provide simple finite pulse time corrections to the thermal diffusivity involving the duration and peak of the pulse. Application of the thermal diffusivity formula (17)–(19) to arbitrary pulse shapes is addressed later in Section 4.2.

An interesting question is the following: if one assumes the heat pulse occurs instantaneously but actually occurs over a finite duration, what is the effect on the error in the thermal diffusivity estimate? The answer to this question is evident from the formulas for $\alpha$ in Eqs (20)–(22). For example, if the true pulse was triangular, the thermal diffusivity would be estimated by $\alpha_{I}=L^{2}/(6I_{T})$ instead of $\alpha_{T}=L^{2}/[6(I_{T}-\frac{\tau+\beta}{3})]$ leading to an underestimate since $\tau$ and $\beta$ are both positive numbers. Going further, the relative error in the calculation can be expressed as $\left|\alpha_{I}-\alpha_{T}\right|/\alpha_{T}=2(\tau+\beta)\alpha_{I}/L^{2}$. For the parameter values considered in this paper (Table 1), this result leads to a large relative error of approximately $27\%$, which highlights the importance of incorporating finite pulse time effects in the calculations.

3 Two layer sample

3.1 Heat flow model

We now consider a two-layer sample comprised of two homogeneous materials with constant but different thermodynamic properties. The first layer has thickness $\ell_{1}$ and extends from the front surface of the sample at $x=0$ to the interface between the two layers at $x=\ell_{1}$ while the second layer has thickness $\ell_{2}$ and extends from the interface at $x=\ell_{1}$ to the rear surface of the sample at $x=\ell_{1}+\ell_{2}=:L$ (Figure 1c). The first layer is assumed to have thermal conductivity $k_{1}$, volumetric heat capacity $\rho_{1}c_{1}$ and thermal diffusivity $\alpha_{1}=k_{1}/(\rho_{1}c_{1})$ while the second layer has thermal conductivity $k_{2}$, volumetric heat capacity $\rho_{2}c_{2}$ and thermal diffusivity $\alpha_{2}=k_{2}/(\rho_{2}c_{2})$. Let $\mathcal{T}_{1}(x,t)$ be the temperature of the sample in the first layer at location $0\leq x\leq\ell_{1}$ and time $t\geq 0$ and $\mathcal{T}_{2}(x,t)$ be the temperature of the sample in the second layer at location $\ell_{1}\leq x\leq L$ and time $t\geq 0$. We assume the temperature rise in the first and second layers, $T_{1}(x,t)=\mathcal{T}_{1}(x,t)-T_{0}$ and $T_{2}(x,t)=\mathcal{T}_{2}(x,t)-T_{0}$, after initiation of the heat pulse at $t=0$, satisfies the following heat flow model [13, 15]:

[TABLE]

where the internal boundary conditions (29)–(30) specify continuity of temperature and heat flux at the interface between the two layers.

In a similar manner to the homogeneous case, provided $\lim_{t\rightarrow\infty}q(t)=0$ and the total amount of heat $Q_{\infty}$ absorbed at the front surface (6) is finite, the steady-state solution of the two-layer heat flow model (23)–(30) is given by [13, 19]:

[TABLE]

This result is derived by applying the interface conditions (29)–(30) and the thermally insulated boundary conditions (27)–(28) that arise in the limit $t\rightarrow\infty$ to the linear steady-state solutions of the heat equations (23)–(24) yielding $\lim_{t\rightarrow\infty}T_{1}(x,t)=\lim_{t\rightarrow\infty}T_{2}(x,t)=T_{\infty}$, where $T_{\infty}$ is a constant. As with the homogeneous case, the value of $T_{\infty}$ can be identified by noting that the change in the amount of heat in the sample must be balanced by the amount of heat entering the sample through the front surface:

[TABLE]

where $Q(t)$ is as defined previously (9). At steady state this heat conservation statement reduces to $T_{\infty}\left[\rho_{1}c_{1}\ell_{1}+\rho_{2}c_{2}\ell_{2}\right]=Q_{\infty}$ yielding the stated result (31).

3.2 Thermal diffusivity

We now extend the analysis of Section 2.2 to two-layer samples. The aim is again to derive a formula for the thermal diffusivity in terms of an integral involving the theoretical rear-surface temperature rise history, which for the two-layer heat flow model (23)–(30) is represented by $T_{2}(L,t)$. For the two-layer case, we have separate definitions of Eq (10) in each layer:

[TABLE]

Using an identical procedure to the one carried out for the single-layer case yields the following differential equations and boundary conditions:

[TABLE]

The solution of Eqs (35)–(37) is

[TABLE]

where $a_{0}$ and $b_{0}$ are arbitrary constants. The appropriate interface conditions are given by

[TABLE]

which are derived by combining the interface conditions of the heat flow model (29)–(30) with the definitions of $u_{1}(x)$ and $u_{2}(x)$ (33)–(34). For example, the second interface condition (41) is derived as follows:

[TABLE]

Both $u_{1}(x)$ and $u_{2}(x)$ (38)–(39) already satisfy flux continuity at the interface (41) and therefore this interface condition yields no information about $a_{0}$ and $b_{0}$. As with the single-layer case, an additional condition is required to completely identify the solution. The analogue of the auxiliary condition (15) for the two-layer case is:

[TABLE]

which is derived by substituting the definitions of $u_{1}(x)$ and $u_{2}(x)$ (33)–(34) into the left-hand side of Eq (42), reversing the order of integration, combining the integrals and making use of the expression for $T_{\infty}$ (31) and the heat conservation statement (32). Applying the auxiliary equation (42) and the first interface condition (40), produces a pair of linear equations that allow $a_{0}$ and $b_{0}$ to be identified. This determines the final forms for $u_{1}(x)$ and $u_{2}(x)$, which together comprise the solution to the boundary value problem consisting of Eqs (35)–(37) and (40)–(42). Since we are only interested in the rear-surface behaviour, we give only the final expression for $u_{2}(x)$:

[TABLE]

Equating both expressions for $u_{2}(x)$, Eqs (34) and (43), at the rear surface, $x=L$, and rearranging yields formulas for the thermal diffusivity in either the first or second layers:

[TABLE]

where

[TABLE]

Both formulas (44)–(45) assume the layer thicknesses, $\ell_{1}$ and $\ell_{2}$, and the volumetric heat capacities, $\rho_{1}c_{1}$ and $\rho_{2}c_{2}$, are known, expressing the thermal diffusivity explicitly in terms of $T_{2}(L,t)$, the rear-surface temperature rise history, $T_{\infty}$ the steady-state temperature rise value (31), $Q_{\infty}$ the total amount of heat absorbed through the front surface (6), and $Q(t)$ the amount of heat that has been absorbed through the front surface at time $t$ (9). Calculating the thermal diffusivity of either layer also requires the thermal diffusivity of the other layer to be known, which is consistent with the well-known result that only one thermal diffusivity of a two-layer sample can be calculated using the rear-surface temperature [17, 15]. Since the formula for $I_{q}$ is identical for the homogeneous and two-layer case, the two-layer thermal diffusivity formulas for the rectangular, triangular and exponential pulses, analogous to Eqs (20)–(22), are easily obtained by substituting the expressions for $I_{q}$ from Eqs (20)–(22) into Eqs (44)–(45). Finally, we confirm that removing the heterogeneity by setting $\alpha_{1}=\alpha_{2}=\alpha$, $\rho_{1}c_{1}=\rho_{2}c_{2}=\rho c$ and $T_{2}(L,t)=T(L,t)$ in both formulas (44)–(45) and rearranging for $\alpha$ yields the single-layer formula (17) as expected.

4 Numerical experiments

4.1 Verification of formulas

We now verify the thermal diffusivity formulas derived in Sections 2.2 and 3.2 for the homogeneous and two-layer samples. In practice, the rear-surface temperature rise history takes the form of discrete values: $\widetilde{T}_{0},\ldots,\widetilde{T}_{N}$, where $\widetilde{T}_{i}$ is the sampled value at $t_{i}=i\Delta t$, $\Delta t=t_{N}/N$ and $t_{N}$ is the final sampled time. In the discrete case, provided $t_{N}$ is large enough (see, e.g. [20, 21], for how to estimate steady state times), the integral $I_{T}$ appearing in the thermal diffusivity formulas, Eqs (17)–(19) and (44)–(47), can be approximated using the trapezoidal rule [11]. For example, for the homogeneous sample:

[TABLE]

Verification is carried out using synthetic rear-surface temperature rise data generated from solving the heat flow model using a known set of parameter values (Table 1) given previously by Czél et al. [2]. Both the homogeneous model (3)–(5) and two-layer model (23)–(30) are solved numerically using finite volume schemes consisting of $n$ nodes, as outlined in A and B, respectively. These numerical solutions produce sampled values of the rear-surface temperature rise at equally-spaced discrete times such that the $i$th sampled value $\widetilde{T}_{i}$ is given by $T_{n}(t_{i})$, where $T_{n}$ represents the temperature rise value in the finite volume scheme at the node located at the rear surface ($x=L$, $n$th node). Each estimate of the thermal diffusivity is compared to the target value in Table 1 by calculating the signed relative errors $\varepsilon=(\alpha-\widetilde{\alpha})/\alpha$ (homogeneous) and $\varepsilon_{k}=(\alpha_{k}-\widetilde{\alpha}_{k})/\alpha_{k}$ for $k=1,2$ (two-layer) with the tilde used to indicate the estimated value produced from the rear-surface integral method.

In Table 2, results for both the homogeneous and two-layer samples are presented for the rectangular (20), triangular (21) and exponential (22) pulses with $\tau=0.005\,\text{s}$, $\beta=0.001\,\text{s}$, $Q_{\infty}=7000\,\text{J}\,\text{m}^{-2}$, $N=1000$ temperature rise values, $\Delta t=10^{-4}\,\text{s}$, $t_{N}=0.1\,\text{s}$ and $n=500$ nodes. In all cases, the estimate agrees with the target values of the thermal diffusivity in Table 1 to 4–5 significant digits with a relative error of approximately $10^{-4}$. These small non-zero values are explained by the error present in the numerical experiment incurred from truncating the infinite upper limit of the integral $I_{T}$ at the finite value $t_{N}$ and approximating the resulting integral numerically (48) as well as numerically solving the heat flow models to generate the sampled values. In summary, the results in Table 2 confirm the correctness of the analysis in Sections 2.2 and 3.2.

4.2 Application to noisy data

We now investigate the accuracy of the thermal diffusivity formulas when applied to noisy rear-surface temperature data. To mimic the data arising from the laser flash experiment, we follow Carr [11] and add Gaussian noise to the sampled values from the previous section: $\widetilde{T}_{i}=T_{n}(t_{i})+z_{i}$, where $z_{i}$ is a random number generated from a normal distribution with mean zero and standard deviation $\sigma$. Repeatedly generating synthetic rear-surface temperature rise datasets in this manner and calculating the corresponding thermal diffusivities allows distributions of the signed relative error to be constructed.

In Figure 2, results are reported for 10,000 realisations, low ($\sigma=0.005\,^{\circ}\mathrm{C}$), moderate ($\sigma=0.02\,^{\circ}\mathrm{C}$) and high ($\sigma=0.05\,^{\circ}\mathrm{C}$) levels of noise, and the exponential pulse (22) with $\beta=0.001\,\mathrm{s}$ (as considered in the previous section). Example rear-surface temperature rise histories for each level of noise are given in Figures 2(a)–(c). The relative error distributions in Figures 2(d)–(f) are smoothed using MATLAB’s inbuilt kernel smoothing function ksdensity. These figures also include the $0.5\%$ and $99.5\%$ quantiles for the data; for example, 99% of the estimates of $\alpha_{1}$ have a relative error of between $-1.62\%$ and $+1.51\%$ for $\sigma=0.02\,^{\circ}\mathrm{C}$ (Figure 2e). Similar results are obtained for the rectangular and triangular pulses and are therefore not reported. Comparing Figures 2(e)–(f), we see that the thermal diffusivity estimates are slightly less accurate for the second layer than the first layer. This is because for this particular problem $\alpha_{1}>\alpha_{2}$ (Table 1) so the small error incurred through approximating $I_{T}$ via the trapezoidal rule (48) is amplified by a greater amount in Eq (45), where it is multiplied by $\alpha_{1}$, than in Eq (44), where it is multiplied by $\alpha_{2}$. Overall, the results in Figures 2(d)–(f) are comparable in accuracy to those reported previously for the rear-surface integral method [11].

Our numerical experiments so far have assumed that the laser pulse shape can be described exactly by either the rectangular (20), triangular (21) or exponential (22) pulse. We have also calculated the steady-state temperature, $T_{\infty}$, assuming known values for the amount of heat absorbed at the front surface, $Q_{\infty}$, and the volumetric heat capacity(s), $\rho c$, $\rho_{1}c_{1}$ and $\rho_{2}c_{2}$ (i.e., by evaluating Eq (7) and Eq (31) directly using the parameter values in Table 1). We now explore calculating the thermal diffusivity directly using the discrete pulse shape measurements that arise in a laser flash experiment and a value of $T_{\infty}$ inferred from the rear-surface temperature rise history ($\widetilde{T}_{i}$ where $i=0,\ldots,N$). Given discrete values of the heat pulse over time, $q(\hat{t}_{j})$ where $\hat{t}_{j}=j\delta t/M$ and $j=0,\ldots,M$ (note that the time between pulse measurements $\delta t$ differs from the time between rear-surface temperature rise measurements $\Delta t$), we approximate $Q_{\infty}$ using the trapezoidal rule:

[TABLE]

and $I_{q}$ using two applications of the trapezoidal rule

[TABLE]

where $N_{i}=t_{i}/\delta t$ is the number of intervals of length $\delta t$ comprising the interval $[0,t_{i}]$ (i.e., the approximation of $I_{q}$ assumes that $\Delta t/\delta t$ is an integer). To approximate $T_{\infty}$ from the rear-surface temperature rise history, we use the height of the horizontal line of best fit through the last $p$ measurements (Figure 3b), which yields the mean value

[TABLE]

Results in Figures 3(c)–(e) are given for the discrete heat pulse depicted in Figure 3(a), consisting of 201 discrete values between $t=0$ and $t=0.01\,\text{s}$ (inclusive) giving $M=200$ and $\delta t=5\times 10^{-5}\,\text{s}$ in the approximations of $Q_{\infty}$ (49) and $I_{q}$ (51). The discrete values are randomly generated by adding Gaussian noise of mean zero and standard deviation $40000\,\text{J}\text{m}^{-2}$ to the exponential pulse (22) with $\beta=0.001\,\mathrm{s}$. Figure 3(b) provides an example realisation of the rear-surface temperature rise history for the homogeneous sample obtained by adding Gaussian noise to the solution of the heat flow model (3)–(5) with $q(t)$ defined as in Figure 3(a). The horizontal red dashed line in Figure 3(b) depicts the utilised approximation of $T_{\infty}$ (52) with $p=200$. Comparing Figures 2(d)–(f) and Figures 3(c)–(e) we see that the approximations (49) and (51) roughly double the length of the [0.5%, 99.5%] quantile interval.

5 Conclusion

This paper has focussed on methods for calculating thermal diffusivity from laser flash experimental data. Specifically, we have shown how to account for real heat pulse shapes and two-layer samples using the recently-proposed rear-surface integral method [11]. The main contributions of the paper are new formulas for calculating the thermal diffusivity that correct for arbitrary heat pulse shapes:

1. (i)

Eqs (17)–(19) for calculating the thermal diffusivity of a homogeneous sample; and 2. (ii)

Eqs (44)–(47) for calculating the thermal diffusivity of one layer of a two-layer sample given the thermal diffusivity of the other layer.

Computational results confirmed the accuracy of the derived formulas and demonstrated how the formulas can be applied to the kinds of experimental data arising from the laser flash experiment.

The formulas derived in this paper are limited to thermally insulated samples, in which case the heat flow model reduces to a one-dimensional model in the direction of the thickness of the sample. Incorporating heat losses from the front and rear surfaces of the sample can be carried out by combining the analysis in this paper with that presented previously by Carr [11, Section 5]. Possible directions for future research also include extension of the rear-surface integral method with finite-pulse time effects to three-layer samples or two or three dimensional heat flow models; the latter of which could allow for heat losses from the other surfaces of the sample (Figure 1b–c) to be accounted for in the formulas.

The analysis presented in this paper can also be extended to other applications involving diffusion in homogeneous or layered materials in a straightforward manner. For example, by appropriately modifying the boundary and interface conditions (27)–(30), the rear-surface integral method can help parameterise a mathematical model of heat transfer through layered skin formulated by McInerney et al. [22].

Acknowledgements

This research was partially supported by the Australian Mathematical Sciences Institute (AMSI) who provided CJW with a 2018-2019 Vacation Research Scholarship.

Appendix A Homogeneous finite volume scheme

In this appendix, we outline the finite volume scheme used to obtain a numerical solution to the homogeneous heat flow model (3)–(5). We discretise the interval $[0,L]$ using a uniform grid consisting of $n$ nodes such that the $k$th node is located at $x=kh=:x_{k}$, where $k=1,\ldots,n$ and $h=L/(n-1)$. Let $T_{k}(t)$ denote the numerical approximation to $T(x,t)$ at $x=x_{k}$. The finite volume scheme takes the form:

[TABLE]

with initial condition $T_{k}=0$ at $t=0$ (4) for $k=1,\ldots,n$ and right-hand side defined by

[TABLE]

This initial value problem is solved using MATLAB 2017B’s ode15s solver on the interval $[0,t_{N}]$ with specified absolute and relative tolerances (AbsTol and RelTol) of $10^{-12}$ and the interval of integration (tspan) specified to return the solution at $N+1$ equally-spaced discrete times: $t_{i}=i\Delta t$, where $i=0,\ldots,N$ and $\Delta t=t_{N}/N$. Using this solution, the rear-surface temperature rise history is given by the discrete values $T_{n}(t_{0}),\ldots,T_{n}(t_{N})$.

Appendix B Two-layer finite volume scheme

In this appendix, we describe the finite volume scheme used to solve the two-layer heat flow model (23)–(30). We discretise the interval $[0,L]$ using a grid consisting of $n=n_{1}+n_{2}-1$ nodes such that the $k$th node is located at $x=x_{k}$ with $x_{k}=(k-1)h_{1}$ for $k=1,\ldots,n_{1}$ and $x_{k}=\ell_{1}+(k-n_{1})h_{2}$ for $k=n_{1}+1,\ldots,n_{1}+n_{2}-1$, where $h_{1}=\ell_{1}/(n_{1}-1)$ and $h_{2}=\ell_{2}/(n_{2}-1)$. Let $T_{k}(t)$ approximate $T(x_{k},t)$ as in A. The finite volume scheme takes the form:

[TABLE]

with initial condition $T_{k}=0$ at $t=0$ (25)–(26) for $k=1,\ldots,n_{1}+n_{2}-1$ and right-hand side defined by

[TABLE]

This initial value problem is then solved to obtain the rear-surface temperature rise history in an identical manner to that described for the homogeneous case at the end of A.