# Study on the use of motion compensation techniques to determine heat sources. Application to large deformations on cracked rubber specimens

Thomas Pottier, M.P. Moutrille, Jean-Benoit Le Cam, Xavier Balandraud,

Michel Grediac

<figure><img src="image_1.png" /></figure>

## To cite this version:

Thomas Pottier, M.P. Moutrille, Jean-Benoit Le Cam, Xavier Balandraud, Michel Grediac. Study on the use of motion compensation techniques to determine heat sources. Application to large deformations on cracked rubber specimens. Experimental Mechanics, 2009, 49, pp.561-574. 10.1007/s11340-008-9138-0 . hal-00454679

HAL Id: hal-00454679

https://hal.science/hal-00454679v1

Submitted on 18 May 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Study on the Use of Motion Compensation Techniques to Determine Heat Sources. Application to Large Deformations on Cracked Rubber Specimens

T. Pottier · M.-P. Moutrille · J.-B. Le Cam ·

X. Balandraud · M. Grédiac

**Abstract** This paper deals with the determination of the thermal response of elastomeric materials subjected to cyclic loading. In this case, the material undergoes large deformations, so a suitable motion compensation technique has been developed to track the material points and their temperature during the test. Special attention is paid to the Narcissus effect and to the detector matrix of the infrared camera used in the study. Heat sources are then derived from the temperature maps. The thermoelastic inversion phenomenon has been experimentally evidenced during a cyclic test performed on an elastomeric notched specimen. The heat source distribution close to the crack tip has also been deduced from the temperature maps, thus highlighting the relevance of the approach.

**Keywords** Infrared thermography ·

Large displacements · Heat sources · Rubber ·

Motion compensation

## Introduction

The present work deals with the mechanical behavior of cracked elastomeric materials. Studying the crack propagation in elastomeric materials is a critical issue that has been thoroughly studied in the past [1–5] for instance. The experimental techniques used for studying the phenomena, which take place in cracked coupons, are often limited to classical tools such as scanning electron microscopy [6] or X-rays microtomography [7] for instance. To date, infrared thermography has proved to be a relevant technique for studying engineering materials such as steel, aluminium and composites. Various studies previously carried out by Chrysochoos and coworkers have shown that heat sources produced by the material itself were more relevant than temperatures when analyzing various phenomena such as Luder’s bands [8], fatigue [9] or strain localization [10]. The main reason is that the temperature field is influenced by conduction as well as heat exchanges with ambient air and grips, unlike local heat sources.

The strain level in the above cases remains small compared to that obtained in elastomeric materials. Consequently, the tools developed for processing temperature maps provided by infrared cameras are no longer suitable for large deformations, especially when heterogeneous temperature fields are processed. It seems that this issue has only seldom been addressed in

T. Pottier
Laboratoire SYMME, Polytech'Savoie, BP 80439,
74944 Annecy le Vieux Cedex, France
e-mail: Thomas.Pottier@univ-savoie.fr

M.-P. Moutrille · J.-B. Le Cam (✉)
X. Balandraud · M. Grédiac (SEM member)
Laboratoire de Mécanique et Ingénieries (LaMI) Blaise Pascal University (UBP), French Institute for Advanced Mechanics (IFMA) Campus de Clermont-Ferrand / les Cézeaux, B.P. 265 63175 Aubière cedex, France
e-mail: lecam@ifma.fr

M.-P. Moutrille
e-mail: marie-pierre.Moutrille@ifma.fr

X. Balandraud
e-mail: xavier.balandraud@ifma.fr

M. Grédiac
e-mail: Michel.Grediac@univ-bpclermont.fr
the literature. In [11] for instance, the authors propose two motion compensation techniques to account for this effect, but it seems that the strain amplitude is not significant since only thermoelastic effects are studied while calculations are carried out within the framework of large deformations.

The aim here is to account for large deformations in the data processing. A suitable technique, which enables us to track certain material points, is presented in the first section. These points are used to describe the current geometry of the coupons under test. Special attention is paid to the treatment of the biais in the temperature fields due to the non-uniformity of the camera's detectors. This leads to a suitable motion compensation technique which, in turn, enables the correct calculation of the heat source distribution in any configuration, especially near the crack tip. This procedure is described in detail and the assumptions under which it is constructed are justified. Two final examples are presented to illustrate the relevancy of the approach.

## Background

### Heat Diffusion Equation

This section briefly summarizes the thermomechanical framework that is usually used to calculate heat sources [12]. The local state axiom [13] is assumed. Any thermodynamical system out of equilibrium is considered as the sum of several homogeneous subsystems at equilibrium. The thermodynamic process is considered as a quasi-static phenomenon. The state of any material volume element is defined by N state variables: temperature T, one of the strain tensors denoted E and some internal variables $V_1, V_2, \dots, V_{N-2}$ such as plastic strain or volume fractions of certain phases. The specific free energy potential is written as $\Psi(T, E, V_k)$. Considering the first and second principles of thermodynamics and assuming that Fourier's law is used to model the heat conduction, the heat diffusion equation can be written as follows [14]:

$$ \rho C_{E,V_k} \dot{T} - \operatorname{div}(K \operatorname{grad} T) - r = \underbrace{d_1 + \rho T \frac{\partial^2 \Psi}{\partial T \partial E} \dot{E} + \rho T \frac{\partial^2 \Psi}{\partial T \partial V_k} \dot{V}_k}_{s} \quad (1) $$

where $\rho$ is the density, $C_{E,V_k}$ is the specific heat at constant values of E and $V_k$, K is the thermal conductivity tensor and r is the external heat source. The right-hand side of equation (1) represents the heat sources

s produced by the material itself. It can be split into separate terms:

- the mechanical dissipation $d_1$ (or intrinsic dissipation): this positive quantity corresponds to the heat production due to various mechanical irreversibilities such as internal friction;

- the thermomechanical couplings: they correspond to the couplings between the temperature and the other state variables. The coupling between temperature and strain $\rho T \frac{\partial^2 \Psi}{\partial T \partial E} \dot{E}$ is the thermoelastic coupling [15, 16]. In case of purely thermoelastic behavior, this leads to a temperature decrease (negative heat source) for a positive strain rate $\dot{E}$ and conversely. Coupling effects involving $V_k$ may create significant heat sources. For instance, when $V_k$ represents the volume fraction of a given phase k, the $\rho T \frac{\partial^2 \Psi}{\partial T \partial V_k} \dot{V}_k$ quantity is a latent heat production [17].

### Usual Assumptions to Calculate Heat Sources

The classical approach used to assess the heat sources from the temperature fields obtained by an infrared camera [18, 19] is detailed in this section.

**Bidimensional formulation of the heat diffusion equation** It is first assumed that the problem is bidimensional. For the sake of simplicity, heat conduction will be considered isotropic. Equation (1) can be then rewritten as follows:

$$ \rho C_{E,V_k} \dot{T} - k \Delta T - r = s \quad (2) $$

where $\Delta$ is the laplacian operator and k the conductivity coefficient that replaces the thermal conductivity tensor K when the material is isotropic. Let us now consider a flat specimen. $(x, y, z)$ is a cartesian coordinate system such that the x-direction is the main direction of the specimen and the z-direction is perpendicular to the mid-plane. The specimen is assumed to be thin, so the through-thickness temperature T is assumed to be almost constant. On the other hand, the temperature gradient $\frac{\partial T}{\partial z}$ is not constant near the front and back sides of the specimen because of heat exchanges with ambient air [14]. Averaging equation (2) through the thickness of the specimen leads to:

$$ \rho C_{E,V_k} \dot{T} - k \left( \frac{\partial^2 \bar{T}}{\partial x^2} + \frac{\partial^2 \bar{T}}{\partial y^2} \right) - \frac{1}{e} k \left[ \frac{\partial T}{\partial z} \right]_{-e/2}^{e/2} - r = s \quad (3) $$

where e is the specimen thickness and $\bar{T}$ is the mean through-thickness temperature. The material
parameters $\rho$, $C_{E,V_k}$, $k$, the external heat sources $r$ and the heat sources $s$ produced by the material itself are also assumed to be constant through the thickness of the specimen. The heat exchanges with the ambient air on the external surfaces are modeled as follows:

$$ \begin{cases} k \frac{\partial T}{\partial z} \left( \frac{e}{2} \right) = -h (\bar{T} - T_{amb}) \\ k \frac{\partial T}{\partial z} \left( \frac{-e}{2} \right) = h (\bar{T} - T_{amb}) \end{cases} \quad (4) $$

where $h$ is a convection coefficient and $T_{amb}$ is the ambient temperature. Introducing equation (4) in equation (3) leads to:

$$ \rho C_{E,V_k} \dot{\bar{T}} - k \left( \frac{\partial^2 \bar{T}}{\partial x^2} + \frac{\partial^2 \bar{T}}{\partial y^2} \right) + \frac{2h}{k} (\bar{T} - T_{amb}) - r = s \quad (5) $$

For the sake of simplicity, the overlined symbol is removed from the equation and the $2h/(k \times \rho C)$ ratio is denoted $1/\tau$ hereafter. $\tau$ can be considered as a time constant characterizing the heat exchange with ambient air. This leads to the 2D-version of the heat diffusion equation:

$$ \rho C_{E,V_k} \left( \dot{T} + \frac{T - T_{amb}}{\tau} \right) - k \Delta_{2D} T - r = s \quad (6) $$

where $\Delta_{2D}$ is the two-dimensional laplacian operator in the $(x, y)$ plane.

**Reference temperature field** It is difficult to estimate the external heat sources $r$, so the objective now is to remove these quantities from the equation using temperature variations $\theta$:

$$ \theta = T - T_0 \quad (7) $$

where $T_0$ is the reference temperature field. This reference temperature is captured by the IR camera when the material does not produce any heat source. The $T_0$ field can be measured just before loading the specimen for instance. In this case, equation (6) becomes:

$$ \rho C_{E,V_k} \left( \dot{T}_0 + \frac{T_0 - T_{amb}}{\tau} \right) - k \Delta_{2D} T_0 - r = 0 \quad (8) $$

The external heat sources $r$, the time constant $\tau$ and the material parameters $\rho$, $C_{E,V_k}$ $k$ are assumed to be constant during the test. In this case, the following equation is obtained by subtracting equations (8) and (6):

$$ \rho C_{E,V_k} \left( \dot{\theta} + \frac{\theta + T_0 - T_{amb}}{\tau} \right) - k \Delta_{2D} \theta = s \quad (9) $$

In practice, temperature fields are captured by the IR camera. To estimate the heat sources $s$ produced by the material, the procedure consists in calculating the left-hand side of equation (9) by processing the bidimensional temperature fields.

## Experimental Setup

### Material and Specimens

The material considered here is a 34 part per hundred of rubber in weight (phr) carbon black filled natural rubber (CB-NR). It is obtained from a compound that is allowed to cure for 7 minutes at 160 degree C. The chemical composition and some of the thermomechanical characteristics that will be used for the heat source calculation are given in Table 1. The specific heat and the thermal conductivity have been quoted from the literature [20]. The other quantities have been measured using suitable standard tests that are not detailed here. The specimen geometry is depicted in Fig. 1(a). This is a 4 mm thick specimen whose sides are slightly curved to ensure a higher local deformation state in the gauge section. A 2 mm long crack is initiated with a razor blade at the center of the specimen before testing. Two metallic inserts have been bonded at the top and bottom of the elastomeric material to be able to grip the specimen in the testing machine. The specimens are then prepared for temperature measurements. The surfaces are slightly polished and cleaned. They are considered sufficiently flat and plane to generate an

<table><caption>**Table 1** Material formulation and thermomechanical properties [6, 20]</caption><thead><tr><th>Formulation/property</th><th>Value</th></tr></thead><tbody><tr><td>Material formulation (phr)</td><td></td></tr><tr><td>Natural rubber</td><td>100</td></tr><tr><td>Zinc oxide</td><td>9.85</td></tr><tr><td>Oil</td><td>3</td></tr><tr><td>Carbon black</td><td>34</td></tr><tr><td>Sulfur</td><td>3</td></tr><tr><td>Stearic acid</td><td>3</td></tr><tr><td>Antioxidant</td><td>2</td></tr><tr><td>Accelerator</td><td>4</td></tr><tr><td>Thermomechanical properties</td><td></td></tr><tr><td>Shore A hardness</td><td>58</td></tr><tr><td>Elongation at break (%)</td><td>635</td></tr><tr><td>Stress at break (MPa)</td><td>22.9</td></tr><tr><td>Density</td><td>1.13</td></tr><tr><td>Specific heat (J/Kg.K)</td><td>2100</td></tr><tr><td>Thermal conductivity (W/m.K)</td><td>0.8</td></tr></tbody></table>
<figure><img src="image_3.png" /><figcaption>**Fig. 1** Specimen geometry (a) and domains #1 (b) and #2 (c)</figcaption></figure>

emissivity close to one. Carbon black fillers (34 phr, see Table 1) make surfaces naturally black, so no special surface preparation is required.

Loading and Measurements

Three different specimens have been tested under uni-
axial cyclic loading using a 15 kN MTS testing machine.
The signal is sinusoidal and the loading frequency fL is
set to 0.5 Hz. The stretch ratio λ is defined by the ratio
between the maximum and the initial lengths. It varies
between 1.03 and 1.67. Such a stretch ratio level is

usually applied to rubber specimens subjected to cyclic
loading [6]. In order to avoid the Mullins effect [21],
specimens are previously tested under the same loading
conditions and over 10 cycles.

Temperature measurements are performed at am-
bient temperature with a Cedip Jade III-MWIR in-
frared camera which features a local plane array of
320 × 240 pixels and detectors with a wavelength range
of 3.5 – 5 µm. The integration time is 1500 µs and
the acquisition frequency fa is 150 Hz. The data sheet
provided by the manufacturer of the camera provides
a noise equivalent temperature difference of 20 mK
for a temperature range of 5 – 40 degree C. It has
been checked through some simple experiments that
the actual thermal resolution exhibits the same order
of magnitude. In order to ensure that the internal
temperature of the camera is optimal for performing
the measurements, it is set up and switched on for one
hour before the experiment. The stabilization of the
temperature in the camera is necessary to avoid any
drift of the measurements during the test.

Since the specimen undergoes large deformations,
any given material point at its surface clearly moves
during a loading cycle while the zone captured by
the detector array of the infrared camera remains un-
changed. The smaller the distance to the fixed grip,
the greater the amplitude of this displacement. A key
issue is to establish the correspondence between the
geometry of the coupon at any stage of the loading
and the reference geometry to be able to analyze and
process correctly the temperature maps at different
loading amplitudes. It seems that constructing a suit-
able motion compensation technique in the context of
mechanical testing and temperature measurement us-
ing an infrared camera has only seldom been addressed
in the literature. Sakagami et al. [11] recently proposed
two techniques. The first one consists in installing a
brindled pattern with different infrared emissivities on
the coupon. Displacements and deformations on the
surface are obtained by processing these brindled pat-
terns with a digital image correlation software. The
second technique does not use any brindled pattern.
Visible images are taken simultaneously by a digital
camera and an infrared camera. The images captured
by the first camera are processed with a suitable digital
image correlation software to provide the displacement
field and the IR camera directly displays the tempera-
ture field. It seems however that these two techniques
have been developed for smaller strain amplitudes than
those given in the current work.

The procedure developed in the present study is
somewhat different. The idea is to track a limited
number of material points during their displacement.
These points are in fact reflective spots that are plotted on the rubber surface before testing with a PILOT G-1 0.7 gold pen. The spot density is adjusted in such a way that the spot density increases as the expected strain gradient increases. This density is a trade-off between two phenomena: the higher the density, the better motion compensation within the high-gradient zones, but the lower the surface of measurement. Indeed, the measurement performed at the spot surface is not considered in this analysis. Far from the crack, the distance between the reflective spots is about 4 mm [see Fig. 1(a)]. The spot density is higher around the crack. Due to the difference of emissivity between ink and plain elastomer, reflective spots clearly appear on the temperature map. Their area is about 2 pixels × 2 pixels. Reflective spots are considered small enough to avoid a significant loss of pixels used for processing the temperature fields.

Two main zones are investigated in this work. Their location and size are shown in Fig. 1(b) and (c). The first zone enables us to observe the surface of the specimen when the maximum stretch ratio (1.67) is reached. The second zone focuses on the area around the crack. Here, the spatial resolution is the size of the surface observed by any pixel. This depends on the magnification of the objective, and therefore on the size of the zone under investigation. For zones 1 (the whole surface) and 2 (the zone surrounding the crack), it is equal to 219 and 131 µm, respectively.

## Characterization of Time Constant τ

The time constant τ involved in equation (6) can be determined experimentally. This constant depends on the convection phenomenon that defines the heat exchange at the interface with the ambient air. In the present case of large cyclic displacements, it can be assumed that the convection coefficient h appearing in equation (9) varies along the surface of the sample. This coefficient depends on the distance from the fixed grip since the velocity in the ambient air also varies along this direction. The convection coefficient h is expected to be lower near the fixed grip than near the moving grip. The loading is cyclic, so h also depends on time.

As a first approach and in order to focus on the post-processing technique, it is assumed here that parameter τ is a constant that is measured when the specimen is static. The following procedure is used for measuring this quantity. First, the specimen is cooled while resting for thirty minutes in a refrigerator. It must be pointed out that one of the ends of the specimen is fixed to a support made of a thermal insulating material. Second, the specimen and its support are removed from the

<figure><img src="image_4.png" /><figcaption>Fig. 2 Temperature vs. time. Experimental curve and exponential fitting</figcaption></figure>

refrigerator and temperature fields are captured during the return to thermal equilibrium at ambient temperature. Conduction between specimen and support is reduced thanks to the nature of the constitutive material of the support. Figure 2 presents a typical temperature evolution for a given point of one of the specimens. As shown in this figure, this curve can be reasonably fitted by the following exponential function:

$$T = T_{\text{amb}} + (T_{\text{init}} - T_{\text{amb}}) e^{-t/\tau} \quad (10)$$

where $T_{amb}$ is the ambient temperature and $T_{init}$ is the initial temperature (for $t = 0$). $τ$ is identified during the last stage of this heating phase (1 degree C range) to obtain a quantity that is as representative as possible of the actual testing conditions of the specimen when subjected to cyclic loading. This quantity has been measured at four different points. The scatter is low and the mean value is $τ = 345$ s. The tests performed on the other two specimens led to similar values.

## Data Post-processing for Large Displacement Cases

### Motion Compensation Technique

The main objective of the post-processing technique described in this section is to calculate the left-hand side of equation (9) by processing the experimental data.

At this stage, two configurations must be clearly distinguished: the reference geometry and the current one. In our case, the reference is a situation for which slight preloading is applied to the specimen to avoid any local buckling. Calculations are carried out in the current configuration, results are displayed in the reference one.
As explained above, the large displacements undergone by the specimen do not allow direct calculation of heat sources in the reference geometry since a given material point plotted in the reference geometry does not match its counterpart in the deformed configuration when observed through the infrared camera, especially away from the fixed grip. Consequently, a given pixel of the IR detector matrix does not correspond to the same material point plotted during the deformation of the specimen surface. The objective is then to track the material points before calculating the left-hand part of equation (9). A suitable algorithm is developed using the Matlab package and a 'reshaping' operation of the temperature maps is performed. For this purpose, the procedure is illustrated by some schematic views in Fig. 3. It relies on the following steps:

**Step 1:** As explained above, the reference geometry is the undeformed geometry. As a result, all temperature maps will be reshaped and displayed in this given geometry. The spots plotted on the specimen surface are visually localized. From a practical point of view, the co-ordinates of the center of each spot are considered as the co-ordinates of the spot itself. In the following, the line and column indexes of a spot center in the 320 × 240 matrix are denoted *ispot*<sub>low</sub> and *jspot*<sub>low</sub> respectively [see Fig. 3(a)].

**Step 2:** The spots are then tracked during the cyclic loading. For the sake of simplicity, their locations are interpolated between the lowest and highest imposed displacements.

**Substep 2.1:** The spots at the highest imposed displacement are localized by hand on the screen and recorded. The line and column indexes of this spot are denoted (*ispot*<sub>high</sub>, *jspot*<sub>high</sub>) [see Fig. 3(b)].

**Substep 2.2:** At a given time $t_k$, the displacement of the spot from the reference geometry is obtained using a sine function as a first approximation since the applied loading is a sinusoidal function of time. Thus

$$ 
\left\{ 
\begin{aligned}
\Delta \text{ispot} &= \frac{\text{ispot}_{\text{high}} - \text{ispot}_{\text{low}}}{2} + \frac{\text{ispot}_{\text{high}} - \text{ispot}_{\text{low}}}{2} \sin \left( 2\pi f_L \frac{t_k - \varphi}{f_a} \right) \\
\Delta \text{jspot} &= \frac{\text{jspot}_{\text{high}} - \text{jspot}_{\text{low}}}{2} + \frac{\text{jspot}_{\text{high}} - \text{jspot}_{\text{low}}}{2} \sin \left( 2\pi f_L \frac{t_k - \varphi}{f_a} \right)
\end{aligned}
\right. 
\quad (11) 
$$

<figure><img src="image_5.png" /><figcaption>(a) step 1: choice of the reference geometry</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>(b) step 2: choice of a current geometry at time $t_k$</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(c) step 3: determination of the displacement field</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>(d) step 4: temperature evolution map plotted in the reference geometry</figcaption></figure>

<caption>**Fig. 3** Motion compensation processing. (a)–(d) Steps 1–4</caption>
where $\varphi$ is the time shift of the reference image, $f_L$ is the loading frequency and $f_a$ is the acquisition frequency of the IR camera. For a given value $t_k$, the location of the spot centers is obtained as follows:

$$ \left\{ \begin{aligned} ispot &= ispot_{low} + \Delta ispot \\ jspot &= jspot_{low} + \Delta jspot \end{aligned} \right. \qquad (12) $$

**Step 3:** For a given time $t_k$, the displacement at any pixel $(i, j)$ of the IR detector matrix is sought. This displacement is obtained by interpolating the displacements of the spots using a finite element approach where the nodes are the spot centers. The specimen is then meshed with triangular T3 elements using the Delaunay procedure [22]. The spots recorded from the picture are the nodes. The $(\Delta i, \Delta j)$ displacement is deduced at each pixel $(i, j)$ from the displacement $(\Delta ispot, \Delta jspot)$ of the nodes using linear shape functions [22] [see Fig. 3(c)].

**Step 4:** At any time $t_k$, the temperature value measured at any pixel $(i, j)$ is shifted to the pixel $(i_{ref}, j_{ref})$ of the reference geometry. Thus:

$$ \left\{ \begin{aligned} i_{ref} &= i - \Delta i \\ j_{ref} &= j - \Delta j \end{aligned} \right. \qquad (13) $$

The procedure described above provides the temperature evolutions in the reference geometry [see Fig. 3(d)]. It must be emphasized that it has an important consequence in terms of spatial resolution of the temperature fields. Since the surface of each finite element changes during a loading cycle, the number of pixels included in each element also evolves. Two cases are considered to account for this feature in the calculations. First, several pixels $(i, j)$ in the deformed configuration may correspond to the same pixel in the reference geometry $(i_{ref}, j_{ref})$. In this case, the different values are simply averaged to find the final value allocated to a given point defined in the reference geometry. Second, certain pixels $(i_{ref}, j_{ref})$ of the reference geometry may not have any allocated temperature, thus leading to some voids in the temperature map. This appears for instance when the area of a given finite element is smaller in the current geometry than in the reference geometry. It has been presently decided to fill in these 'blanks' with the average of the surrounding non-empty pixels. The eight adjacent pixels are used for this purpose.

The calculation of heat sources relies on the temperature variations $\theta$ [see equation (9)]. A specific

procedure when subtracting the reference temperature $T_0$ [see equation (7)] must be used because the IR detectors of any matrix array camera feature a slight non-uniformity. This procedure is described in detail in the following section.

## Influence of the IR Detector Matrix Non-uniformity

It is well known that a matrix of IR detectors is characterized by a non-uniformity in terms of gain and offset. The impact of this non-uniformity can be limited by performing a so-called Non-Uniformity Correction (NUC). A so-called offset correction is performed for this purpose. It consists in allocating the same encoded temperature to the pixels by modifying the offsets while capturing a thermal scene that is assumed to be homogeneous. Whatever the result of the NUC operation, the temperature maps are 'altered' by an additional non-uniform temperature distribution $H$ due to a halo. As shown in Fig. 4, its amplitude is about a few tenths of a degree. It has been verified that if the camera used for performing the experiment is running for more than one hour before starting the test, this field remains approximately constant: $H(i, j, t_k) = H(i, j)$. This property is very useful for small displacements of the material points since subtracting the reference temperature field $T_0$ does not affect the temperature variations $\theta$. Let us now examine the consequence of the non-uniformity of the IR detector matrix in the present case of large displacements.

At a given material point $M$ located at pixel denoted $(i_{ref}, j_{ref})$ in the initial temperature field, the temperature is the sum of the actual temperature $T_0$ $(i_{ref}, j_{ref})$ and some small additional value $H$ $(i_{ref}, j_{ref})$ due to non-uniformity. Thus, one can write:

$$ T_{0\text{camera}}(M) = T_0(i_{ref}, j_{ref}) + H(i_{ref}, j_{ref}) \qquad (14) $$

<figure><img src="image_9.png" /><figcaption>Fig. 4 Temperature distribution due to $H$</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Fig. 5** Temperature variation in a given measurement zone. Curve a: without motion compensation processing; curve b: with a motion compensation where the H field is not subtracted; curve c: with a motion compensation where the H field is subtracted</figcaption></figure>

The actual difference between $\theta(M)$ and $T_{\text{camera}}(M) - T_{\text{0 camera}}(M)$ is equal to a few tenths of a degree. This is a problem if actual thermal scenes of the same order of magnitude are processed. The calculation of heat sources relies on the measurement of the actual temperature variations $\theta$ at each time. It is therefore proposed to subtract field $H$ from reference field $T_{\text{0 camera}}$ and from any other temperature map $T_{\text{camera}}$ before launching the motion compensation technique described above. In practice, field $H$ is captured before performing the mechanical test on a black body at ambient temperature, with the IR camera running for more than one hour before capturing the temperature fields. More precisely, a one-second sequence at acquisition frequency $f_a = 150$ Hz is captured and the 150 images are averaged to obtain a filtered field.

At time $t_k$, the same material point is located at a new pixel $(i, j)$. The motion compensation technique described above gives the link between $(i, j)$ and $(i_{\text{ref}}, j_{\text{ref}})$. Its temperature provided by the camera at time $t_k$ can also be split in two parts denoted $T$ and $H$:

$$T_{\text{camera}}(M) = T(i, j) + H(i, j) \quad (15)$$

Subtracting equation (15) and equation (14) to obtain the temperature variation from the camera data leads to:

$$T_{\text{camera}}(M) - T_{\text{0 camera}}(M) = \theta(M) + H(i, j) - H(i_{\text{ref}}, j_{\text{ref}}) \quad (16)$$

where $\theta(M) = T(i, j) - T_{\text{0 camera}}(i_{\text{ref}}, j_{\text{ref}})$ is the actual temperature variation of a given material point $M$. So equation (16) proves that the following inequality is verified in the general case:

$$T_{\text{camera}}(M) - T_{\text{0 camera}}(M) \neq \theta(M) \quad (17)$$

<figure><img src="image_11.png" /><figcaption>**Fig. 6** Temperature variations versus time for 7 zones along the specimen</figcaption></figure>

## Application to a Rigid-body Motion

A first test is performed on a specimen subjected to rigid-body displacement. One side of the specimen is fixed in the moving grip and the other side remains free. The moving grip describes a cyclic linear translation whose amplitude is 29.45 mm. A temperature gradient of about 6 degrees C is generated within the specimen using a frozen steel block placed on the free side. The temperature field is measured by the IR camera. It is assumed that the temperature field in the specimen does not change during a short acquisition time of a few seconds. The motion compensation technique is applied to obtain the temperature evolution at each material point. No heat is produced by the material during its rigid-body motion. Thus, no temperature variation is expected once the post-processing method is applied. The results obtained are shown in Fig. 5. In this figure, temperature variations $\theta$ are plotted for three different cases. In each case, the considered temperature of one
material point is obtained by averaging temperatures over a zone of 9 pixels. The first curve (curve a) shows the temperature variation obtained without applying the motion compensation technique. At a given time, the temperature of one given material point is measured in the zone of 9 pixels. The material point moves but the measurement zone remains fixed. Thus, the temperature measured in this zone is the temperature of a set of material points whose temperature is different. The second curve (curve b) represents the temperature variation of the zone discussed above when the motion compensation technique is applied without subtracting the H field. The temperature variation calculated by the program is different to zero, so this quantity does not correspond to the actual temperature variation of the material point. As shown in the last case (curve c), applying the motion compensation technique that takes into account the influence of the H field leads to temperature variations lower than 0.04 degree C. This quantity is considered as small enough to validate the processing.

## Numerical Implementation and Calculation of the Heat Sources

This section describes the numerical calculation of the left-hand side terms in equation (9) once the procedure described above is applied.

### Absorption term $\rho C_{E,V_k} \frac{\partial \theta}{\partial t}$
This involves the calculation of a first-order derivate with respect to time. It is approximated at each pixel ($i_{ref}$, $j_{ref}$) using the following finite difference scheme:

$$ \frac{\partial \theta}{\partial t} (i_{ref}, j_{ref}, t_k) \approx \frac{\theta (i_{ref}, j_{ref}, t_k + p) - \theta (i_{ref}, j_{ref}, t_k - p)}{2p \times \Delta t} \quad (18) $$

where $p$ is the half-width of the temporal window and $\Delta t$ is the time step. The temporal window is defined by the number of successive pictures used to perform the temporal derivation. Because of noise, parameter $p$ is set at 4 in the present case to smoothen the result. Preliminary numerical tests, which are not given here, have shown that this value is optimal.

### Heat exchange with ambient air $\rho C_{E,V_k} \frac{\theta + T_0 - T_{amb}}{\tau}$
This term is obtained at each pixel ($i_{ref}$, $j_{ref}$). It is simply equal to:

$$ \rho C_{E,V_k} \frac{\theta (i_{ref}, j_{ref}, t_k) + T_0 (i_{ref}, j_{ref}, t_k) - T_{amb}}{\tau} $$

If the amplitude of temperature variation is considered to be equal to some tenths of a degree, the order of magnitude of $\frac{\partial \theta}{\partial t}$ is $\frac{0.1}{\text{second}}$ for a loading frequency equal to 0.5 Hz. The order of magnitude of $\frac{\theta + T_0 - T_{amb}}{\tau}$ is roughly equal to $\frac{0.1C}{345 \text{ seconds}}$. As a result, it can be

<figure><img src="image_12.png" /><figcaption>Fig. 7 Stretch ratio, temperature variations and heat sources vs. time (a)–(c). Zone 3</figcaption></figure>
easily checked that the heat exchange with ambient air is negligible compared to the absorption term.

In-plane heat conduction k Δ₂Dθ A finite difference scheme is used for the two terms of the laplacian operator:

$$ \Delta_{2D}\theta(i_{\text{ref}}, j_{\text{ref}}, t_k) = \frac{\partial^2\theta}{\partial x^2}(i_{\text{ref}}, j_{\text{ref}}, t_k) + \frac{\partial^2\theta}{\partial y^2}(i_{\text{ref}}, j_{\text{ref}}, t_k) \quad (19) $$

The present approach considers the temperature variations θ in the current geometry. Since this quantity has been transferred in the reference geometry thanks to the above data processing method, the usual three-point discretization with constant spatial increments cannot be used here. Once the motion compensation procedure is applied, the Δx (as well as Δy) increment is not constant in the reference geometry. The following equations are preferred in this case:

$$ \left\{ \begin{aligned} \frac{\partial^2\theta}{\partial x^2}(i_{\text{ref}}, j_{\text{ref}}, t_k) &= \frac{2\theta(i_{\text{ref}}, -m, j_{\text{ref}}, t_k)}{\Delta x^-(\Delta x^- + \Delta x^+)} - \frac{2\theta(i_{\text{ref}}, j_{\text{ref}}, t_k)}{\Delta x^- \Delta x^+} \\ &\quad + \frac{2\theta(i_{\text{ref}}, +m, j_{\text{ref}}, t_k)}{\Delta x^+(\Delta x^- + \Delta x^+)} \\ \frac{\partial^2\theta}{\partial y^2}(i_{\text{ref}}, j_{\text{ref}}, t_k) &= \frac{2\theta(i_{\text{ref}}, j_{\text{ref}}, -n, t_k)}{\Delta y^-(\Delta y^- + \Delta y^+)} - \frac{2\theta(i_{\text{ref}}, j_{\text{ref}}, t_k)}{\Delta y^- \Delta y^+} \\ &\quad + \frac{2\theta(i_{\text{ref}}, j_{\text{ref}}, +n, t_k)}{\Delta y^+(\Delta y^- + \Delta y^+)} \end{aligned} \right. \quad (20) $$

where m and n are the half-widths (in pixel) of the spatial windows along the x- and y-directions respectively. The spatial window is defined by the number of pixels used to calculate the spatial derivation. The distances from the central pixel (i ref, j ref) of the window in the current geometry are denoted Δx-, Δx+, Δy- and Δy+.

<figure><img src="image_13.png" /><figcaption>**Fig. 8** Temperature variations vs. time for 7 zones located around the crack tip</figcaption></figure>

The numerical procedure to obtain the above quantities is not detailed here.

The weight of the in-plane heat conduction on the total sources s is estimated at each pixel by calculating the following ratio denoted R:

$$ R = \frac{k |\max(\Delta_{2D}\theta)|}{\rho C |\max(\frac{\partial\theta}{\partial t}) - \min(\frac{\partial\theta}{\partial t})|} \times 100 \quad (21) $$

The following section presents the results obtained in domains #1 and #2 defined in Fig. 1. In both cases, it has been observed that ratio R is very small. The in-plane conduction term is negligible compared to the absorption term: less than 1% for most of the specimen and less than 1.5% around the crack. Consequently, this term is ignored in the following sections.

## Results

The procedure described above is now applied to the determination of the heat sources. Representative results from one of a series of three specimens are presented here. In this section, results obtained in domains #1 and #2 are presented [see Fig. 1(b) and (c)]. The first one enables us to observe the whole specimen up to its maximum stretching.

### Domain #1: Whole Surface

The post-processing method is used to track the material points and their temperature variation during cyclic loading. Then, variations of temperature can be plotted in the reference geometry. Here, the spot size is small, so the emissivity difference between spots and substrate does not disturb the temperature measurement of the substrate itself. Seven zones are considered along the surface of the specimen. They are defined in practice by 3 × 3 = 9 pixels. Figure 6 shows that each zone is equidistant from the others and from the metallic
inserts. The temperature reported here is the average
temperature of the seven zones.

Three main conclusions can be drawn from this
figure:

- The thermal response of the material is not strictly sinusoidal, especially for very small strain amplitudes that correspond to the smallest temperature variations in Fig. 6. This is due to a well-known phenomenon in elastomeric material for very small strain amplitudes [23, 24]. In fact, the thermal response is caused by two types of coupling, namely the thermoelastic and the isentropic couplings. Thermoelastic and isentropic couplings are either first or second order phenomena, depending on the deformation level. For example, a small increase of the deformation from the undeformed state induces a temperature decrease. This is a consequence of the first order thermoelastic effects in this range of deformation. However, if the deformation is large enough, the isentropic coupling becomes greater than the thermoelastic coupling. It induces a temperature increase. A minimum is obtained in the thermal response curve when the two coupling phenomena exhibit the same order of magnitude. This phenomenon may be referred to as thermoelastic inversion [24].

- The maximum temperature variation in zones 3, 4 and 5 is greater than the temperature variation in zones 1, 2, 6 and 7 at any time. The maximum temperature is obtained for zone 4. This is simply due to the fact that the gauge section is the most stretched one, thus leading to the highest temperature variation in this zone.

- For the same stretch ratio level, the temperature of zones 1 and 7 on the one hand, and the temperature of zones 2 and 6 on the other hand, are different. This is due to the moving grip whose temperature increases during fatigue tests.

The latter result illustrates that the temperature is
not actually a relevant indicator of the thermomechan-
ical response of materials because of the parasitic phe-
nomena such as conduction between grips and coupon.
Heat sources are more suitable quantities, as already
mentioned by Chrysochoos [18]. To illustrate this, heat
sources are calculated at each material point and for
each different time.

The results obtained are presented in Fig. 7. For
the sake of simplicity, only the results obtained for
zone 3 are plotted. Other zones exhibit a similar re-
sponse. In this figure, the stretch ratio, the temperature

variations and the heat sources are plotted versus time.
It is to be noted that the stretch ratio considered here
is a macroscopic quantity (λ<sub>macro</sub>), i.e. it is calculated
from the displacement of the moving grip. Calculation
performed by the Finite Element Method has shown
that the microscopic stretch ratio (λ<sub>micro</sub>) reaches 1.14

<figure><img src="image_14.png" /><figcaption>**Fig. 9** Temporal plot of stretch ratio, temperature variations and heat sources in zones 1 and 5</figcaption></figure>
and 5 at the crack tip while $\lambda_{macro}$ is equal to 1.03 and 1.67, respectively. The following conclusions can be drawn:

- as expected, the heat sources are almost equal to zero when the strain rate defined by $\frac{d\lambda}{dt}$ is zero.

- as explained above, temperature variations in Fig. 7(b) are not strictly sinusoidal. For small strain amplitudes ($\lambda < 1.1$), a short plateau may be observed. It is related to the thermoelastic inversion that takes place when both the thermoelastic and the isentropic phenomena exhibit the same order of magnitude. The thermoelastic inversion does not clearly appear in the curve in Fig. 7(b), but a significant change in slope is observed in Fig. 7(c), at about $\lambda_{macro} = 1.1$.

The results presented in this section show that the method developed here allows to track material points and to calculate heat sources produced at each of these material points and at any time. Moreover, when the whole specimen surface is considered, the spatial resolution around the crack is not sufficient to perform the calculation. Consequently, a new analysis of certain zones surrounding the crack (domain #2) has been performed using a better spatial resolution. The results obtained are presented below.

## Domain #2: Zone Surrounding the Crack

The same method as above is applied here. The purpose is to show that the proposed motion compensation technique allows us to assess the heat sources from the temperature field measurement around the crack which features more significant strain gradients.

Similarly to the previous analysis, the reflective spots located around the crack are tracked during the material deformation. The temperature of any material points are then deduced. Seven zones are chosen to present the results in the vicinity of the crack. Figure 8 shows the location of the seven zones and their temperature evolution. Processing the temperature fields using the procedure described above shows that temperature evolutions in zones 3, 4 and 5 are quite different from their counterparts obtained in zones 1, 2, 6 and 7. Zone 3, 4 and 5 are located near the crack tip, i.e. in a high-stretched area, and zones 1, 2, 6 and 7 are located in low-stretched zones. The temperature variations in zones 1 and 7 as well as zones 2 and 6 are quite similar. This can be explained by the symmetry of the crack and by the fact that the zones are close to each other. Moreover, the temperature evolution in zones 1 and 7 (i.e. the less stretched zones), clearly highlights the thermoelastic effects discussed above. This phenomenon does not appear in zones 3, 4 and 5 in which the minimum strain level is higher. As explained above, as $\lambda_{macro}$ varies between 1.03 and 1.67 at the macroscopic scale, $\lambda_{micro}$ varies between 1.14 and 5 at the crack tip. It is worth noting that the temperature amplitude in zone 4 is lower in Fig. 6 than in Fig. 8. This is due to the higher strain level in this zone in the latter case. On the other hand, temperature variations are different between the two cases (here, lower in the second case). This is simply due to the fact that the two measurements have not been performed at the same time during the test.

Figure 9 is the temporal plot of the macroscopic stretch ratio $\lambda_{macro}$, the temperature variations and the heat sources in zones 1 and 5. As expected, the maximum value of the heat sources is obtained in the crack tip zone (illustrated by zone 5 here). This observation is in agreement with the high stretch level in this zone.

<figure><img src="image_15.png" /><figcaption>**Fig. 10** Heat source map in the zone surrounding the crack tip</figcaption></figure>
Due to the high strain rate, the evolution of heat source in this zone is similar to the evolution of zone 3 in the previous analysis. The fact that the thermoelastic inversion is correctly described in zone 1 is due to:

i- the very low strain rate that allows a correct sampling of the signal;

ii- the minimum stretch ratio reached during one cycle that is close to one, i.e. below the value of the stretch ratio at which the thermoelastic inversion occurs. This value is referred to as $\lambda_{ti}$ in the following. It is equal to 1.06 according to common values found in the literature (see for example [25]).

On the other hand, in zone 5, FE simulations have shown that the minimum stretch ratio during one cycle does not reach one, but a quantity that remains slightly greater than $\lambda_{ti}$. As a result, the slope of the heat sources decreases, but does not reach zero.

Finally, Fig. 10 shows an example of the heat sources map in the zone surrounding the crack tip ($\lambda_{macro} = 1.3$, during unloading). As may be seen, the crack is slightly open because of the preloading. In this figure, it clearly appears that the maximum value in terms of heat absorption is obtained at the crack tip and that a high gradient of heat sources is observed when going from the crack tip to the bulk material. Moreover, heat sources are close to zero in the low-stretched areas on both sides of the crack.

## Conclusion

A motion compensation technique has been developed in this study. It enables us to track the temperature variation of any material point subjected to large displacements. A first analysis, carried out on elastomeric material at a large scale corresponding to the whole sample surface, highlights the competition between thermoelastic and isentropic couplings. Both of them contribute to the thermal response of this type of material. A second analysis is performed with a refined spatial resolution in the domain surrounding the crack. A map of heat sources around the crack is obtained using the procedure. It is observed that the heat source density is highest near the crack tip while no heat sources are detected along the crack lips, close to the free boundary. The main perspective is to go further in terms of spatial resolution in order to investigate more closely the crack tip and to track crack propagation during the fatigue test.

**Acknowledgements** The support of this research by the "Agence Nationale pour la Recherche" is gratefully acknowledged (PHOTOFIT project). Thanks are also due to the French Laboratory of the Trelleborg company for providing the material tested in this study.

## References

1. Lindley PB, Gent AN, Lake GJ (1964) Cut growth and fatigue of rubbers. I. The relationship between cut growth and fatigue. J Appl Polym Sci 8:455-466

2. Gent AN, Hindi M (1990) Effect of oxygen on the tear strength of elastomers. Rubber Chem Technol 63:123-134

3. Greensmith HW (1956) Rupture of rubber. IV. Tear properties of vulcanizates containing carbon black. J Polym Sci 21:175-187

4. Rivlin RS, Thomas AG (1953) Rupture of rubber. I. Characteristic energy for tearing. J Polym Sci 3:291-318

5. Thomas AG (1958) Rupture of rubber. V. Cut growth in natural rubber vucanizates. J Polym Sci 31:467-480

6. Le Cam JB, Huneau B, Verron E, Gornet L (2004) Mechanism of fatigue crack growth in carbon black filled natural rubber. Macromolecules 37:5011-5017

7. Le Gorju Jago K (2007) Fatigue life of rubber components: 3D damage evolution from X-ray computed tomography. In: Boukamel A, Laiarinandrasana L, Méo S, Verron E (eds) Constitutive models for rubber. V. Balkema, London

8. Chrysochoos A, Louche H (2001) Thermal and dissipative effects accompanying luders band propagation. Mater Sci Eng A Struct 307:15-22

9. Berthel B, Chrysochoos A, Wattrisse B, Galtier A (2007) Infrared image processing for the calorimetric analysis of fatigue phenomena. Exp Mech (in press)

10. Wattrisse B, Chrysochoos A, Muracciole JM, Némoz-Gaillard M (2001) Analysis of strain localization during tensile tests by digital image correlation. Exp Mech 41(1):29-39

11. Sakagami T, Nishimura T, Yamaguchi T, Kubo N (2006) Development of a new motion compensation technique in infrared stress measurement based on digital image correlation method. Nihon Kikai Gakkai Ronbunshu A / Trans Jpn Soc Mec Eng 72:1853-1859

12. Nguyen QS, Germain P, Suquet P (1983) Continuum thermodynamics. J Appl Sci 50:1010-1020

13. Boccara N (1968) Les principes de la thermodynamique classique. In: PUF coll. SUP

14. Chrysochoos A, Louche H (2000) An infrared image processing to analyse the calorific effects accompanying strain localisation. Int J Eng Sci 38:1759-1788

15. Boulanger T, Chrysochoos A, Mabru C, Galtier A (2004) Calorimetric analysis of dissipative and thermoelastic effects associated with the fatigue behavior of steels. Int J Fatigue 26:221-229

16. Berthel B, Wattrisse B, Chrysochoos A, Galtier A (2007) Thermoelastic analysis of fatigue dissipation properties of steel sheets. Strain 43(3):273-279

17. Balandraud X, Chrysochoos A, Leclercq S, Peyroux R (2001) Influence of the thermomechanical coupling on the propagation of a phase change front. C R Acad Sci Ser IIB Mec 329(8):621-626
18. Chrysochoos A, Maisonneuve O, Martin G, Caumon H, Chezeau JC (1989) Plastic and dissipated work and stored energy. Nucl Eng Des 114:323-333

19. Chrysochoos A, Belmhajoub F (1992) Thermographic analysis of thermomechanical couplings. Arch Mech 44(1):55-68

20. Granta (2006) CES selector. http://www.grantadesign.com/products/ces/

21. Mullins L (1948) Effect of stretching on the properties of rubber. Rubber Chem Technol 21:281-300

22. Zienkiewicz OC (1977) The finite element method in engineering science. McGraw-Hill, London

23. Gough J (1985) Proc Lit and Phil Soc Manchester, 2nd, ser. 1:288

24. Joule JP (1857) On some thermodynamic properties of solids. Philos Mag 4th 14:227

25. Anthony RL, Caston RH, Guth E (1942) Equations of state for naturals and synthetic rubber like materials: unaccelerated natural soft rubber. J Phys Chem 46:826