# Numerical simulation of the heat generation process at vibrating contact defects

K. Truyaert, S. Delrue, V. Aleshin, K. van den Abeele

► To cite this version:

K. Truyaert, S. Delrue, V. Aleshin, K. van den Abeele. Numerical simulation of the heat generation process at vibrating contact defects. 2018 Quantitative InfraRed Thermography, Jun 2018, Berlin, Germany. 10.21611/qirt.2018.088. hal-02364174

HAL Id: hal-02364174

https://hal.science/hal-02364174

Submitted on 6 Jul 2022

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Numerical simulation of the heat generation process at vibrating contact defects

by K. Truyaert*, S. Delrue*, V.Aleshin** and K. Van Den Abeele*

* Wave Propagation and Signal Processing Research Group, KU Leuven Campus Kulak Kortrijk, E. Sabbelaan 53, Kortrijk, Belgium, kevin.truyaert@kuleuven.be

**Univ. Lille, CNRS, Centrale Lille, ISEN, Univ. Valenciennes, UMR 8520 – IEMN, LIA LICS/LEMAC, F-59000 Lille, France

## Abstract

Ultrasonic thermography is a promising technique for non-destructive evaluation of vibrating defects, due to its ability to perform fast, full field material inspections. A major drawback, however, is the inefficient conversion of vibrational energy to heat. To overcome this problem a solution based on Local Defect Resonance (LDR) excitation has been proposed. To validate this new concept and to find an optimal set of operational parameters it is crucial to have realistic models that can describe the ultrasonic wave-defect interaction. An important aspect therein is the modelling of frictional heat generation at the defect interface, which will be considered here.

## 1. Introduction

Realistic modelling of vibrating contacts is important as it may support experimental techniques for ultrasonic Non-Destructive Testing and Evaluation (NDT&E) of structures, buildings, materials or components with internal damage. Indeed, dedicated comparison between experimental data and numerical results is crucial to characterise internal defects and to decide whether to repair or replace the damaged piece of material.

In this study, we propose a 2D numerical model that allows to describe the specific vibrational and thermal behaviour of ultrasonically excited contact defects. Apart from clapping phenomena (i.e. repetitive opening and closing of the contact defect), often implemented by introducing different elastic moduli for normal compression and tension, friction effects also need to be included. This of course makes the model much more complicated as phenomena such as hysteresis and memory-dependence also have to be taken into account. Due to this hysteresis, the energy of the motion at the contact points will be dissipated to its environment in the form of thermal energy. Moreover, roughness of the contact profile also plays an important role as it may have a significant influence on the contact vibrational and thermal behaviour. To take all these effects into account, we use the semi-analytic Method of Memory Diagrams (MMD) that was especially developed to calculate the hysteretic frictional response of rough surfaces and has the advantage that no discontinuous contact model needs to be used, making all calculations explicit. The method uses a simple memory function to describe the different contact properties, instead of updating complex stress and displacement fields at the contact, drastically simplifying the contact description. Moreover, the instantaneous friction-induced loss of mechanical energy is also easily calculated. The full MMD contact model is directly integrated into a finite element modelling environment enabling the study of the dynamic and thermal behaviour of structures with vibrating contact defects.

The description of the multiphysics model is given in the next sections, where the contact model and the energy dissipation calculations are discussed in more detail. The potential of the method is illustrated by means of an example where the interaction between an ultrasonic wave and an internal defect are discussed.

## 2. Description of the frictional contact model

This section includes the model that is used to describe the specific behaviour of a dynamically excited contact defect. First, the mechanical contact problem is discussed. After this, the model is extended to also consider the friction-induced heat generation.

### 2.1. Mechanical contact problem

In the model that is proposed here, a crack is split up into several mesoscopic cells (see figure 1). The size of each cell is chosen to be smaller than the characteristic scale of the macroscopic elastic field, but large enough to host a large number of microscopic features. This way, a cell includes the average normal (N) and tangential (T) load per unit area and average half normal (a) and half tangential (b) displacements, all linked by a relationship that can theoretically be derived for axisymmetric, convex contact profiles. The use of such contact profiles also allows to include roughness in the model, as every complex contact profile with asperities can be replaced by an equivalent axisymmetric contact profile [1]. Moreover, the mesoscopic cell formulation with the axisymmetric, convex contact profile as the link between the normal force and displacement, allows rough crack faces to be modelled as flat in the final numerical model. Finally, by describing the contact in each mesoscopic cell separately, it is possible to solve the full contact mechanics problem.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /><figcaption>Fig. 1: Division of a crack into different mesoscopic cells.</figcaption></figure>

Two non-flat contact profiles will result in a mechanical contact that will be in one of the following three contact states: contact loss, full sliding and partial slip. The latter only occurs due to the presence of asperities at the contact interface. Surface points with different heights will be compressed differently upon a normal loading scenario and in combination with a tangential loading, these points will experience a different tangential behaviour (stick or slip), based on their compression.

The link between the load properties (N and T') and the displacement properties (a and b) is determined by the Method of Memory Diagrams (MMD) [2, 3]. This method uses a semi-analytical internal memory function that describes the hysteretic effects and is based on the Reduced Elastic Friction Principle (REFP, [4, 5]). This principle allows the tangential variables of an axisymmetric contact profile to be expressed as function of the normal variables in the following way:

$$ \left\{ \begin{array}{l} b = \theta\mu(a - q), \\ T = \mu(N(a) - N(q)), \end{array} \right. \qquad (1) $$

where $\mu$ is the friction coefficient between the two surfaces in contact, $\theta$ is a material parameter that is based on Poisson's ratio $\nu$:

$$ \theta = \frac{2-\nu}{2(1-\nu)}, \qquad (2) $$

$N(a)$ is the relation between the normal force and normal displacement which is considered to be known either by experimental observations or from theoretical calculations using contact mechanics formulations, and $q$ is a normal indentation parameter that can be used to separate the stick from the slip region. In the MMD language, Equation (1) is translated into:

$$ \left\{ \begin{array}{l} b = \theta\mu \int_{0}^{a} D(\alpha) d\alpha, \\ T = \mu \int_{0}^{a} D(\alpha) \left. \frac{dN}{da} \right|_{a=\alpha} d\alpha. \end{array} \right. \qquad (3) $$

The memory function, $D(\alpha)$, contains all information of the current and previous contact situation. Instead of redistributing a complex traction distribution in the contact area, this memory function is updated. The rules about how to update the memory function based on the history of its drive parameter, the tangential displacement ($b$) in this work, and how to calculate the output, the tangential load ($T'$) in this work, can be found in previous works about the MMD [2, 3]. The shape of the memory function, $D(\alpha)$, depends on the loading history.

With this description, the partial slip phase can be described. This phase is active as long as the surfaces are in contact and Coulomb's frictional description holds: in force terms $|T'| < \mu N$, or in displacement terms $|b| < \theta\mu a$. If the threshold value is reached, the partial slip phase switches to the full sliding contact phase. In case of total sliding, the tangential displacement ($b$) is split into two components: $b = b_0 + \tilde{b}$. The first component, $b_0$, is the total slipped distance. This distance is defined as the current tangential distance between two points on opposite contact faces, that have the same tangential coordinate in the unperturbed geometry. The second parameter, $\tilde{b}$, takes into account the fact that asperities recede under tangential loading. It is this parameter that will be used to update the memory function in the MMD during the loading process. The memory function for total sliding equals $\pm 1$ over the whole interval $0 \le \alpha \le a$, depending on the direction of the tangential load in the 2D model.
<table><caption>**Table 1.** The three possible contact states in the model of contact defects with rough surfaces, with for each state the conditions under which the contact state occurs, the solutions for the tangential displacement components $b_0$ and $\tilde{b}$, the solutions for the normal and tangential forces $N$ and $T$, and the memory diagram $D(\alpha)$. (Taken from [3])</caption><thead><tr><th>Crack state</th><th>If</th><th>Then</th><th>Memory diagram</th></tr></thead><tbody><tr><td>(i) Contact loss</td><td>$a(t) < 0$</td><td>$\tilde{b}(t) = 0$<br>$b_0(t) = b(t)$<br>$T(t) = N(t) = 0$</td><td>1↑<br>-1↑<br>D(α)<br>→ α</td></tr><tr><td>(ii) Total sliding</td><td>$a(t) \ge 0$<br>$|\tilde{b}(t - \Delta t)| \ge \theta \mu a(t)$</td><td>$\tilde{b}(t) = \pm \theta \mu a(t)$<br>$b_0(t) = b(t) - \tilde{b}(t)$<br>$T(t) = \pm \mu N(t)$</td><td>1↑<br>-1↑<br>|D(α)|<br>→ α</td></tr><tr><td>(iii) Partial slip</td><td>$a(t) \ge 0$<br>$|\tilde{b}(t - \Delta t)| < \theta \mu a(t)$</td><td>$b_0(t) = b_0(t - \Delta t)$<br>$\tilde{b}(t) = b(t) - b_0(t)$<br>$T(t) = MMD(\tilde{b}(t))$</td><td>1↑<br>-1↑<br>D(α)<br>→ α</td></tr></tbody></table>

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

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

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

The description for the contact loss phase is simple. In this case ($a < 0$) no contact happens and hence, $N = T = 0$. The tangential displacement component $\tilde{b}$ is zero as well, as the asperities remain unstrained. The memory function is zero for all $\alpha$.

Table 1 contains a summary of the three possible contact states and the different criteria that depict the differentiation between them. An illustrative example of a memory function $D(\alpha)$ is given for all three cases.

## 2.2. Energy dissipation induced through friction

A contact experiencing a dynamic loading, will experience friction between its two contact surfaces during the partial slip and full sliding phases. As a result of this friction, energy will be dissipated in the form of heat. In the case of a periodic normal and tangential loading, a constant hysteretic loop is run through each cycle. It has been shown that the amount of dissipated energy during such a periodic loading corresponds to the area of the hysteretic curve [6]. This property can be used to calculate the energy released by friction during each periodic loading cycle [7].

However, in the more general case of an arbitrary loading procedure, the hysteretic loop is not closed and non-periodic. Therefore, a different method capable of determining the instantaneous energy dissipation, $\Delta W$, needs to be developed. As the MMD formalism is used to describe the relation between the tangential properties of the contact, it is practical to write the instantaneous friction-induced energy dissipation $\Delta W$ in terms of parameters that are intrinsic to the MMD as follows:

$$ \Delta W = -2\mu(|\Delta b| - \theta\mu\Delta a) \cdot [N(a) - N(q) + (q - a)\frac{dN}{da}]_{a=\alpha}, \quad (4) $$

where all variables are already present in the above described contact model. A detailed derivation and analysis of this general instantaneous energy dissipation equation will be published soon [8]. Applying this equation to the problem of a contact between two flat surfaces, the known solution for work done by friction is obtained: $\Delta W = 2|\Delta b|T$. The other terms have a specific meaning that is related to the axisymmetric contact profile, where energy is stored in the deformation of this profile.

## 3. Coupled numerical finite element model

The simulation of ultrasonic thermography experiments requires a coupled thermo-mechanic model consisting of three components in order to simulate the behaviour of an ultrasonically excited contact where heat is generated at the contact interface: a crack model, a solid mechanics unit and a thermal unit. The crack model takes into account the behaviour of the normal and tangential parameters at the contact interface. The solid mechanics unit is needed to solve the elastic wave equations for the ultrasonic wave passing through the specimen and the thermal unit is required to define the contact interface as a heat source and to solve the thermal diffusion equation. In this work, the commercially available software package COMSOL Multiphysics®, which uses finite elements to solve the differential equations, is used. More particularly, the above described MMD is used to solve the contact problem, through a direct link between the FEM-model and an external procedure implemented in Matlab®. The Structural Mechanics module of COMSOL® [9] is used to solve the elastic wave propagation problem and the Heat Transfer module [10] is used to solve the heat diffusion equation.
<figure><img src="image_7.png" /><figcaption>Fig. 2: Representation of the studied specimen with an inclined internal defect.</figcaption></figure>

The three components of the numerical model are linked to each other to perform a time-dependent simulation as follows:

(a) *Update of displacement and temperature distribution*

Using the normal ($\sigma$) and tangential ($\tau$) stresses, defined by appropriate internal and external boundary conditions, the Structural Mechanics module calculates the displacements in the computational domain. This is done for each time step $t$ during the procedure. The normal ($u_n$) and tangential ($u_t$) displacements at the contact interface are stored as output of the COMSOL® module and are used as an input for the MMD in the external procedure. The Heat Transfer module is responsible for solving the heat diffusion equation to retrieve the temperature distribution in the computation region. This module is also responsible for the calculation of the temperature increase due to the dissipated energy at the defect's faces.

(b) *Update of contact stresses and energy dissipation*

In a second step, the displacement driven MMD is used to calculate the new contact stresses at the integration points of the contact faces using the displacements from step (a). Apart from the contact stresses, the dissipated energy, induced by the frictional behaviour of the material faces in contact, is also calculated. Both the calculated contact stresses and dissipated energy are returned to COMSOL® to update the internal mechanical and thermal boundary condition in, respectively, the Structural Mechanics and Heat Transfer module. It is assumed that all dissipated energy is converged to thermal energy.

(c) *Update time step*

Steps (a) and (b) are repeated for the next time step, until the total calculation time is reached.

This approach is an extension of the implementation described in Delrue et al. [3], where the mechanical contact problem was considered without the account for thermal phenomena.

**4. Numerical example of friction-induced heat generation at an ultrasonically excited crack**

An example of a propagating shear wave in a 2D aluminium sample having an internal crack is discussed here. Apart from nonlinear vibrations [3], the internal crack will generate a thermal signature induced through friction on the surface faces in contact.

**4.1. Model specifications**

A rectangular sample with a width of 5 cm and a height of 10 cm is studied (see figure 2). The material of the sample is aluminium, with a density $ρ = 2700 \text{ kg m}^{-3}$, a Young's modulus $E = 70 \text{ GPA}$ and a Poisson's ratio $ν = 0.33$. The sample contains a crack, the center of which is located at a depth of 5 cm. The crack has a length of 2 cm and is inclined over 20 degrees with respect to the horizontal. The frictional coefficient $μ$ is set equal to one. At the top horizontal boundary, a continuous shear wave is excited at a frequency $f = 100 \text{ kHz}$. The bottom horizontal boundary, as well as the two vertical boundaries are defined as low reflecting boundaries. This eliminates wave reflections on these boundaries, meaning that
<figure><img src="image_8.png" /><figcaption>**Fig. 3:** Tangential stress-displacement curve for the central point on the crack. The sample is excited with a shear wave excitation at 100 kHz, with an amplitude of 100 nm. The hysteretic behaviour of the tangential properties indicates the dissipation of internal energy, which is transformed into thermal energy. The hysteretic behaviour is non-periodic and the instantaneous energy has to be calculated using equation (4).</figcaption></figure>

the model can be seen as a part of a larger component, extending to the left, right and below the considered computational geometry. Additionally, all four (exterior) boundaries are defined as open thermal boundaries. At the internal crack interface, thin elastic layer and heat source boundary conditions are defined to mimic the mechanical and thermal crack behaviour using the MMD formalism described earlier. The normal reaction curve for contact between two aluminium samples is taken from literature [11].

## 4.2. Generation and diffusion of heat

The excited ultrasonic shear wave propagates through the aluminium sample and interacts with the defect, introducing dynamic clapping and friction effects at the contact faces. This is illustrated by the tangential reaction curve of the center point on the crack in figure 3 where a switching between different contact regimes is observed. The hysteretic loops that are created while in contact prove that energy will be released during the frictional process. Apart from the thermal crack behavior, the defect will also introduce nonlinear wave distortions, as discussed in Delrue et al. [3].

Figure 4 shows the temperature evolution at the center of the defect. Heat is generated during each moment of contact (i.e. $N \neq 0, T \neq 0$) and is diffused when the contact is lost (i.e. $N = T = 0$). This causes the alternations in temperature observed in the inset of figure 4. When a non-zero tangential stress is experienced, friction is present and the

<figure><img src="image_9.png" /><figcaption>**Fig. 4:** Temperature evolution at the center of the defect. During the ultrasonic shear wave excitation, an increase in temperature is observed while in contact (i.e. $N \neq 0, T \neq 0$), whereas during contact loss (i.e. $N = T = 0$), thermal energy is diffused, causing a drop in temperature. This is depicted in the inset figures. After 1 ms, the excitation is deactivated and the thermal energy is dissipated.</figcaption></figure>

dissipated energy is transferred into thermal energy, causing the temperature to increase locally. During contact loss, the generated heat will diffuse away from the considered point, causing a drop in temperature. Note that not all energy is diffused during contact loss. Hence, the temperature is globally increasing until a stationary situation is reached. This occurs around 1 ms, when the ultrasonic source is turned off. After 1 ms, the simulation runs further and the thermal energy is diffused through the body. No new contact interactions occur, which means that no extra heat will be generated at the contact interface. This can be seen in figure 4, where the temperature at the middle of the delamination drops exponentially. Snapshots of the temperature distribution inside the aluminium sample at a time of 1 ms and 1 s are shown
in figure 5. The boundaries of the specimen are defined as a thermal sink at a temperature of 293.15 K, the original temperature of the specimen, which clarifies the temperature distribution at 1s.

<figure><img src="image_10.png" /><figcaption>**Fig. 5:** Snapshots of temperature distribution in the region of the crack in the aluminium sample. The first snapshot shows the temperature distribution at the end of the 1 ms excitation. The second figure shows how the generated heat is distributed around the crack after 1 s, where the boundaries of the sample have boundary conditions representing a heat sink at the original temperature of the sample. This figure nicely illustrates the heat generation at the crack interface and the heat diffusion in the aluminium sample.</figcaption></figure>

**5. Conclusion**

The development of realistic thermo-mechanical models for elastic wave propagation in structures with internal defects is important for the field of thermographic non-destructive testing and evaluation as it allows internal defects not only to be detected, but also to be localised and characterised by comparing numerical and experimental data. In the present paper, a 2D thermo-mechanical finite element model, implemented in the commercially available software package Comsol Multiphysics®, has been proposed in which the particular mechanical and thermal behaviour of vibrating contacts is described using the Method of Memory Diagrams.

A particular example of an ultrasonic shear wave propagating in an aluminium sample with internal contact defect has been considered. It was shown that during different contact states (full sliding and partial slip), friction-induced heat is generated at the defect interface. The obtained results thereby illustrate the potential of the code to assist ultrasonic and thermosonic NDT methods for defect characterization.

In the near future, we plan to extend the model by also including thermoelastic effects, where deformation of the material results in a local temperature increase. However, for contact defects, we foresee that friction will be the main source of temperature increase. The final goal is to exploit the model to validate the recent concept of LDR-thermosonics [12, 13], and to perform an in-depth parametric study in view of determining the most important operational parameters for LDR-thermosonics to be used for the quality assessment of different types of materials and structures.

**6. Acknowledgements**

The research leading to these results has gratefully received funding from Internal Funds KU Leuven (C24/15/021).

**REFERENCES**

[1] V. L. Popov et al., *Method of Dimensionality Reduction in Contact Mechanics and Friction*. Berlin, Heidelberg: Springer Berlin Heidelberg, 2015.

[2] V. Aleshin et al., "Two dimensional modeling of elastic wave propagation in solids containing cracks with rough surfaces and friction - Part I: Theoretical background," *Ultrasonics*, vol. 82, pp. 11–18, 2018.

[3] S. Delrue et al., "Two dimensional modeling of elastic wave propagation in solids containing cracks with rough surfaces and friction - Part II: Numerical implementation," *Ultrasonics*, vol. 82, pp. 19–30, 2018.

[4] J. Jäger, "Axi-symmetric bodies of equal material in contact under torsion or shift," *Arch. Appl. Mech.*, vol. 65, no. 7, pp. 478–487, 1995.

[5] M. Ciavarella, "The generalized Cattaneo partial slip plane contact problem. I—Theory," *Int. J. Solids Struct.*, vol. 35, no. 18, pp. 2349–2362, 1998.

[6] J. R. Barber et al., "Frictional elastic contact with periodic loading," *Int. J. Solids Struct.*, vol. 48, no. 13, pp. 2041–2047, Jun. 2011.

[7] C. Putignano et al., "Frictional energy dissipation in contact of nominally flat rough surfaces under harmonically varying loads," *J. Mech. Phys. Solids*, vol. 59, no. 12, pp. 2442–2454, Dec. 2011.

[8] K. Truyaert et al., "Theoretical calculation of the instantaneous friction-induced energy losses in arbitrarily excited axisymmetric mechanical contact systems," *Preprint*.

[9] COMSOL et al., "Structural Mechanics Module User's Guide."
[10] COMSOL et al., "Heat Transfer Module User's Guide."

[11] S. Biwa et al., "On the acoustic nonlinearity of solid-solid contact with pressure-dependent interface stiffness," *Trans. Soc. Mech. Eng. J. Appl. Mech.*, vol. 71, no. 4, pp. 508–515, 2004.

[12] I. Solodov et al., "Highly-efficient and noncontact vibro-thermography via local defect resonance," *Quant. Infrared Thermogr. J.*, vol. 12, no. 1, pp. 98–111, Jan. 2015.

[13] M. Rahammer et al., "Fourier-transform vibrothermography with frequency sweep excitation utilizing local defect resonances," *NDT E Int.*, vol. 86, no. March, pp. 83–88, 2017.