# Simulation of induction thermography NDT technique using SIBC

Abdoulaye Ba, Huu Kien Bui, Gerard Berthiau, Didier Trichet, Guillaume Wasselynck

► To cite this version:

Abdoulaye Ba, Huu Kien Bui, Gerard Berthiau, Didier Trichet, Guillaume Wasselynck. Simulation of induction thermography NDT technique using SIBC. COMPEL: The International Journal for Computation and Mathematics in Electrical and Electronic Engineering, 2020, 39 (5), pp.1071-1083. 10.1108/COMPEL-01-2020-0007. hal-04086607

HAL Id: hal-04086607

https://hal.science/hal-04086607

Submitted on 11 Apr 2024

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
# Simulation of induction thermography NDT technique using SIBC

Abdoulaye Ba, Huu Kien Bui, Gerard Berthiau, Didier Trichet, Guillaume Wasselynck

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

## To cite this version:

Abdoulaye Ba, Huu Kien Bui, Gerard Berthiau, Didier Trichet, Guillaume Wasselynck. Simulation of induction thermography NDT technique using SIBC. COMPEL: The International Journal for Computation and Mathematics in Electrical and Electronic Engineering, 2020, 39 (5), pp.1071-1083. 10.1108/COMPEL-01-2020-0007. hal-04086607

HAL Id: hal-04086607

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

Submitted on 11 Apr 2024

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 Attribution 4.0 International License
# Simulation of Induction Thermography NDT
## Technique Using SIBC

Abdoulaye Ba, Huu-Kien Bui, Gérard Berthiau, Didier Trichet, Guillaume Wasselynck

University of Nantes

Lab. IREENA, 44600 Saint-Nazaire Cedex, France

huu-kien.bui@univ-nantes.fr

**Purpose** - We present a lightened 3D finite element model for coupled electromagnetic thermal simulation of the induction thermography non-destructive testing technique to reduce the computation time.

**Design/methodology/approach** - The time harmonic electromagnetic problem is expressed in $A - \phi$ formulation and lightened by using the Surface Impedance Boundary Condition (SIBC) applied to both the massive induction coil surface and the surface of conductor workpiece including open cracks. The external circuit is taken into account by using the impressed voltage or the impressed current formulation. The thermal diffusion in the workpiece is solved using surface electromagnetic power density as thermal source.

**Findings** - The accuracy and the usefulness of the method for the design of the induction thermography NDT technique have been shown with acceptable deviation compared with a full FEM model. It is also observed that at high frequency, when the ratio between the local radius of the conductor and the skin depth is high, a very good accuracy can be obtained with the SIBC methods. At lower frequency, the effect of the curvature of the surface becomes significant. In this case, the use of the Mitzner's impedance can help to correct the error.

**Originality/value** - The SIBC can be used for both massive coil and workpieces with open cracks to alleviate 3D finite element models of the coupled electrothermal model. The implementation in matrix form of the coupled electrothermal formulation is given in details. The comparisons with reference analytical solution and full 3D finite element model show the accuracy and performance of the method. In the test case presented, the computation time is 6.6 times lower than the classical model.

**Keywords:** Impedance boundary condition, massive coil, induction thermography

# 1 Introduction

Induction thermography can be applied to the Non-Destructive Testing (NDT) of various kind of conductor workpiece such as metal and carbon fiber reinforced polymers. Using the induction heating principle, the piece is heated up to a few dozen degrees above the ambient temperature. The presence of defects deviates the eddy-current and the heat flow. This leads to abnormal temperature distribution on the surface of the piece which can be revealed on thermal images measured by a thermal camera. This technique turns out to be very efficient for open crack detection.

The simulation tools for the design of the NDT technique must be able to deal with massive coils of complex shape, strong skin and proximity effect. The classical full finite element model with impressed external circuit constraints, such as current and voltage, is well adapted for the modeling of complex shape coil. However, this model is very time consuming due to the necessity of fine mesh in the skin depth layer. It is actually true since one is usually working at high frequency with high conductive materials making the coil and the workpiece [1] [2].
# Simulation of Induction Thermography NDT
Technique Using SIBC

Abdoulaye Ba, Huu-Kien Bui, Gérard Berthiau, Didier Trichet, Guillaume Wasselynck

University of Nantes

Lab. IREENA, 44600 Saint-Nazaire Cedex, France

huu-kien.bui@univ-nantes.fr

**Purpose** - We present a lightened 3D finite element model for coupled electromagnetic thermal simulation of the induction thermography non-destructive testing technique to reduce the computation time.

**Design/methodology/approach** - The time harmonic electromagnetic problem is expressed in $A - \phi$ formulation and lightened by using the Surface Impedance Boundary Condition (SIBC) applied to both the massive induction coil surface and the surface of conductor workpiece including open cracks. The external circuit is taken into account by using the impressed voltage or the impressed current formulation. The thermal diffusion in the workpiece is solved using surface electromagnetic power density as thermal source.

**Findings** - The accuracy and the usefulness of the method for the design of the induction thermography NDT technique have been shown with acceptable deviation compared with a full FEM model. It is also observed that at high frequency, when the ratio between the local radius of the conductor and the skin depth is high, a very good accuracy can be obtained with the SIBC methods. At lower frequency, the effect of the curvature of the surface becomes significant. In this case, the use of the Mitzner's impedance can help to correct the error.

**Originality/value** - The SIBC can be used for both massive coil and workpieces with open cracks to alleviate 3D finite element models of the coupled electrothermal model. The implementation in matrix form of the coupled electrothermal formulation is given in details. The comparisons with reference analytical solution and full 3D finite element model show the accuracy and performance of the method. In the test case presented, the computation time is 6.6 times lower than the classical model.

**Keywords:** Impedance boundary condition, massive coil, induction thermography

## 1 Introduction

Induction thermography can be applied to the Non-Destructive Testing (NDT) of various kind of conductor workpiece such as metal and carbon fiber reinforced polymers. Using the induction heating principle, the piece is heated up to a few dozen degrees above the ambient temperature. The presence of defects deviates the eddy-current and the heat flow. This leads to abnormal temperature distribution on the surface of the piece which can be revealed on thermal images measured by a thermal camera. This technique turns out to be very efficient for open crack detection.

The simulation tools for the design of the NDT technique must be able to deal with massive coils of complex shape, strong skin and proximity effect. The classical full finite element model with impressed external circuit constraints, such as current and voltage, is well adapted for the modeling of complex shape coil. However, this model is very time consuming due to the necessity of fine mesh in the skin depth layer. It is actually true since one is usually working at high frequency with high conductive materials making the coil and the workpiece [1] [2].
<figure><img src="image_2.png" /><figcaption>Figure 1: A studied domain comprising a massive coil and a workpiece with open crack flaw. The voltage can be impressed on the contours $Γ_1^s$ and $Γ_2^s$.</figcaption></figure>

We present, in this paper, an accurate lightened 3D finite element model for the coupled electromagnetic thermal simulation of the induction heating process using the Surface Impedance Boundary Condition (SIBC). This boundary condition is applied to both the massive induction coil and the conductor workpiece with open cracks. The finite element model takes into account the external electrical circuit by means of voltage (or current) constrained formulation. The latter allows for accurate computation of the coil impedance which is necessary for setting frequency and power factor. The implementation in $A - φ$ formulation will be given in details. To couple with thermal problem, the surface electromagnetic power density is used as thermal heat flux density.

The paper is organized as follows. In the section 2, the finite element $A - φ$ formulations with SIBC and external circuit constraints are presented. The formulations are then validated on a simple test case by an exact analytical solution. In the section 3, the coupling with thermal problem is then introduced. Finally, an application in Induction Thermography simulation is presented in the section 4.

# 2 Weak formulation of SIBC method

## 2.1 The use of surface impedance

A general configuration of an induction thermography testing is shown in the Fig. 1. When using finite element with SIBC, the volume of the conductor has to be excluded from the finite element resolution, only the surface meshes are required. If all conductors are replaced by their surface, the Maxwell’s equation to be solved is $curlH = 0$ in the remaining nonconducting domains and can be written in term of magnetic vector potential as:

$$ curl(\frac{1}{\mu}curlA) = 0 \qquad (1) $$

where the condition $A = 0|_Γ$ is applied on the boundary $Γ$ of the surrounding air box. By applying the weighted residual method to the equation (1) with the test functions $(w^{e_i} + gradw^{n_i})$ where $w^{e_i}$ and $w^{n_i}$ are respectively the edge and nodal shape function and by using Green’s formulae, we obtain the weak formulation of the electromagnetic problem:

$$ \int_{\Omega} curlw^{e_i} \frac{1}{\mu}curlA \cdot dv - \int_{\Gamma_c} (w^{e_i} + gradw^{n_i})(n \times H)ds = 0, \quad \forall w^{e_i}, \forall w^{n_i} \qquad (2) $$
<figure><img src="image_3.png" /><figcaption>Figure 1: A studied domain comprising a massive coil and a workpiece with open crack flaw. The voltage can be impressed on the contours $Γ_1^s$ and $Γ_2^s$.</figcaption></figure>

We present, in this paper, an accurate lightened 3D finite element model for the coupled electromagnetic thermal simulation of the induction heating process using the Surface Impedance Boundary Condition (SIBC). This boundary condition is applied to both the massive induction coil and the conductor workpiece with open cracks. The finite element model takes into account the external electrical circuit by means of voltage (or current) constrained formulation. The latter allows for accurate computation of the coil impedance which is necessary for setting frequency and power factor. The implementation in $A-φ$ formulation will be given in details. To couple with thermal problem, the surface electromagnetic power density is used as thermal heat flux density.

The paper is organized as follows. In the section 2, the finite element $A-φ$ formulations with SIBC and external circuit constraints are presented. The formulations are then validated on a simple test case by an exact analytical solution. In the section 3, the coupling with thermal problem is then introduced. Finally, an application in Induction Thermography simulation is presented in the section 4.

# 2 Weak formulation of SIBC method

## 2.1 The use of surface impedance

A general configuration of an induction thermography testing is shown in the Fig. 1. When using finite element with SIBC, the volume of the conductor has to be excluded from the finite element resolution, only the surface meshes are required. If all conductors are replaced by their surface, the Maxwell’s equation to be solved is $curlH = 0$ in the remaining nonconducting domains and can be written in term of magnetic vector potential as:

$$ curl(\frac{1}{\mu}curlA) = 0 \qquad (1) $$

where the condition $A = 0|_Γ$ is applied on the boundary $Γ$ of the surrounding air box. By applying the weighted residual method to the equation (1) with the test functions $(w^{e_i} + gradw^{n_i})$ where $w^{e_i}$ and $w^{n_i}$ are respectively the edge and nodal shape function and by using Green’s formulae, we obtain the weak formulation of the electromagnetic problem:

$$ \int_{\Omega} curlw^{e_i} \frac{1}{\mu}curlA \cdot dv - \int_{\Gamma_c} (w^{e_i} + gradw^{n_i})(n \times H)ds = 0, \quad \forall w^{e_i}, \forall w^{n_i} \qquad (2) $$
where $\mathbf{n}$ is the outward normal vector of $\Gamma_c = \Gamma_{c_1} \cup \Gamma_{c_2}$ as shown in the Fig. 1 and $\mathbf{A}$ discretized by edge element. On the surface of the conductors $\Gamma_c$, one applies the SIBC [3]:

$$ \mathbf{n} \times \mathbf{H}|_{\Gamma_c} = \mathbf{n} \times \frac{1}{\mu} \text{curl} \mathbf{A}|_{\Gamma_c} = \frac{1}{Z_c} (\mathbf{n} \times \mathbf{E}) \times \mathbf{n}|_{\Gamma_c} \quad (3) $$

where $Z_c$ is the surface impedance which depends on the conductor properties. In the literature, two expressions of the surface impedance $Z_c$ known respectively as Leontovich, Mitzner approximation [6] can be given as:

$$ Z_c^{\text{Leontovich}} = \frac{1+j}{\sigma\delta} \quad \text{with} \quad \delta = \sqrt{\frac{2}{\omega\sigma\mu}} \quad (4) $$

$$ Z_c^{\text{Mitzner}} = Z_c^{\text{Leontovich}} \left( 1 + \frac{1-j}{4} \delta \left( \frac{1}{r_{H_t}} - \frac{1}{r_{E_t}} \right) \right) \quad (5) $$

where $\omega = 2\pi f$ is the angular frequency, $j$ the imaginary unit, $\sigma$ the electrical conductivity and $\delta$ the skin depth. In the expression of Mitzner's impedance, $r_{H_t}$ and $r_{E_t}$ are the local radii of curvature of the surface. The Mitzner's impedance is an improvement of Leontovich's one that can be used to correct the effect of the curvature of the conductor. Another surface impedance expression known as Rytov's impedance involving partial derivative of the tangential component of $\mathbf{H}$ is not considered in this paper. By noting that $\mathbf{E} = -j\omega(\mathbf{A} + \text{grad}\phi)$ where $\phi$ is the primitive function in time of the electric scalar potential, the weak formulation (2) can be rewritten using the SIBC as:

$$ \int_{\Omega} \text{curl} \mathbf{w}^{e_i} \frac{1}{\mu} \text{curl} \mathbf{A} \mathbf{d}v + \int_{\Gamma_c} (\mathbf{w}^{e_i} + \text{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw(\mathbf{A} + \text{grad} \phi)) \times \mathbf{n} \mathbf{d}s = 0 \quad (6) $$

## 2.2 Impressing voltage V of external circuit

A voltage $V$ can be impressed between the two contours $\Gamma_1^s$ and $\Gamma_2^s$ of the coil end (cf Fig. 1) on which the potential is considered constant. The function $\phi$ can be approximated as [4]:

$$ \phi = \sum_{\Gamma_c \setminus (\Gamma_1^s \cup \Gamma_2^s)} w^{n_j} \phi_j^f + \sum_{\Gamma_1^s \cup \Gamma_2^s} w^{n_j} \phi_j^i \quad (7) $$

where "f" and "i" stand for free and impressed, $\phi_j^i = 0$ on $\Gamma_1^s$ and $\phi_j^i = V/jw$ on $\Gamma_2^s$. To improve the convergence, one can decompose $\phi$ in dynamic term and static term as below:

$$ \phi = \sum_{\Gamma_c \setminus (\Gamma_1^s \cup \Gamma_2^s)} w^{n_j} \phi_j^{\text{dynamic}} + \sum_{\Gamma_{c_2} \setminus (\Gamma_1^s \cup \Gamma_2^s)} w^{n_j} \phi_j^{\text{static}} + \sum_{\Gamma_1^s \cup \Gamma_2^s} w^{n_j} \phi_j^i \quad (8) $$

where $\phi_j^{\text{dynamic}}$ are obtained by solving an electrokinetic problem on the coil domain with impressed potentials $\phi_j^i$. The potentials $\phi_j^{\text{static}}$ is to be solved together with $\mathbf{A}$ in the final dynamic eddy current problem. One can write (8) as:

$$ \phi = \phi^{\text{dynamic}} + \phi^{\text{static}} = \phi^{\text{dynamic}} + \alpha \phi^i \quad (9) $$

with $\alpha$ computed for all nodes of the coil surface. The function $\alpha$ is equal to 1 on $\Gamma_2^s$ and equal to 0 on $\Gamma_1^s$. By applying vector calculus identity $\mathbf{Y}(\mathbf{X} \times \mathbf{n}) = \mathbf{X}(\mathbf{n} \times \mathbf{Y})$ to (6), the final weak formulation of SIBC method reads:

$$ \int_{\Omega} \text{curl} \mathbf{w}^{e_i} \frac{1}{\mu} \text{curl} \mathbf{A} \mathbf{d}v + \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathbf{A}) \mathbf{d}s + \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \text{grad} \phi^{\text{dynamic}}) \mathbf{d}s \\ = - \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \text{grad} \alpha \phi^i) \mathbf{d}s \quad (10) $$

$$ \int_{\Gamma_c} (\mathbf{n} \times \text{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathbf{A}) \mathbf{d}s + \int_{\Gamma_c} (\mathbf{n} \times \text{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \text{grad} \phi^{\text{dynamic}}) \mathbf{d}s \\ = - \int_{\Gamma_c} (\mathbf{n} \times \text{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \text{grad} \alpha \phi^i) \mathbf{d}s \quad (11) $$
where $n$ is the outward normal vector of $\Gamma_c = \Gamma_{c_1} \cup \Gamma_{c_2}$ as shown in the Fig. 1 and $\mathbf{A}$ discretized by edge element. On the surface of the conductors $\Gamma_c$, one applies the SIBC [3]:

$$ \mathbf{n} \times \mathbf{H}|_{\Gamma_c} = \mathbf{n} \times \frac{1}{\mu} \mathrm{curl} \mathbf{A}|_{\Gamma_c} = \frac{1}{Z_c} (\mathbf{n} \times \mathbf{E}) \times \mathbf{n}|_{\Gamma_c} \quad (3) $$

where $Z_c$ is the surface impedance which depends on the conductor properties. In the literature, two expressions of the surface impedance $Z_c$ known respectively as Leontovich, Mitzner approximation [6] can be given as:

$$ Z_c^{\text{Leontovich}} = \frac{1+j}{\sigma\delta} \quad \text{with} \quad \delta = \sqrt{\frac{2}{\omega\sigma\mu}} \quad (4) $$

$$ Z_c^{\text{Mitzner}} = Z_c^{\text{Leontovich}} \left( 1 + \frac{1-j}{4} \delta \left( \frac{1}{r_{H_t}} - \frac{1}{r_{E_t}} \right) \right) \quad (5) $$

where $\omega = 2\pi f$ is the angular frequency, $j$ the imaginary unit, $\sigma$ the electrical conductivity and $\delta$ the skin depth. In the expression of Mitzner's impedance, $r_{H_t}$ and $r_{E_t}$ are the local radii of curvature of the surface. The Mitzner's impedance is an improvement of Leontovich's one that can be used to correct the effect of the curvature of the conductor. Another surface impedance expression known as Rytov's impedance involving partial derivative of the tangential component of $\mathbf{H}$ is not considered in this paper. By noting that $\mathbf{E} = -j\omega(\mathbf{A} + \mathrm{grad}\phi)$ where $\phi$ is the primitive function in time of the electric scalar potential, the weak formulation (2) can be rewritten using the SIBC as:

$$ \int_{\Omega} \mathrm{curl} \mathbf{w}^{e_i} \frac{1}{\mu} \mathrm{curl} \mathbf{A} \mathrm{Adv} + \int_{\Gamma_c} (\mathbf{w}^{e_i} + \mathrm{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw(\mathbf{A} + \mathrm{grad} \phi)) \times \mathbf{n} \mathrm{d}s = 0 \quad (6) $$

## 2.2 Impressing voltage V of external circuit

A voltage $V$ can be impressed between the two contours $\Gamma_1^s$ and $\Gamma_2^s$ of the coil end (cf Fig.1) on which the potential is considered constant. The function $\phi$ can be approximated as [4]:

$$ \phi = \sum_{\Gamma_c \setminus (\Gamma_1^s \cup \Gamma_2^s)} w^{n_j} \phi_j^f + \sum_{\Gamma_1^s \cup \Gamma_2^s} w^{n_j} \phi_j^i \quad (7) $$

where "f" and "i" stand for free and impressed, $\phi_j^i = 0$ on $\Gamma_1^s$ and $\phi_j^i = V/jw$ on $\Gamma_2^s$. To improve the convergence, one can decompose $\phi$ in dynamic term and static term as below:

$$ \phi = \sum_{\Gamma_c \setminus (\Gamma_1^s \cup \Gamma_2^s)} w^{n_j} \phi_j^{\text{dynamic}} + \sum_{\Gamma_{c_2} \setminus (\Gamma_1^s \cup \Gamma_2^s)} w^{n_j} \phi_j^{\text{static}} + \sum_{\Gamma_1^s \cup \Gamma_2^s} w^{n_j} \phi_j^i \quad (8) $$

where $\phi_j^{\text{dynamic}}$ are obtained by solving an electrokinetic problem on the coil domain with impressed potentials $\phi_j^i$. The potentials $\phi_j^{\text{static}}$ is to be solved together with $\mathbf{A}$ in the final dynamic eddy current problem. One can write (8) as:

$$ \phi = \phi^{\text{dynamic}} + \phi^{\text{static}} = \phi^{\text{dynamic}} + \alpha \phi^i \quad (9) $$

with $\alpha$ computed for all nodes of the coil surface. The function $\alpha$ is equal to 1 on $\Gamma_2^s$ and equal to 0 on $\Gamma_1^s$. By applying vector calculus identity $\mathbf{Y}(\mathbf{X} \times \mathbf{n}) = \mathbf{X}(\mathbf{n} \times \mathbf{Y})$ to (6), the final weak formulation of SIBC method reads:

$$ \int_{\Omega} \mathrm{curl} \mathbf{w}^{e_i} \frac{1}{\mu} \mathrm{curl} \mathbf{A} \mathrm{Adv} + \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathbf{A}) \mathrm{d}s + \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathrm{grad} \phi^{\text{dynamic}}) \mathrm{d}s \\ = - \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathrm{grad} \alpha \phi^i) \mathrm{d}s \quad (10) $$

$$ \int_{\Gamma_c} (\mathbf{n} \times \mathrm{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathbf{A}) \mathrm{d}s + \int_{\Gamma_c} (\mathbf{n} \times \mathrm{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathrm{grad} \phi^{\text{dynamic}}) \mathrm{d}s \\ = - \int_{\Gamma_c} (\mathbf{n} \times \mathrm{grad} \mathbf{w}^{n_i}) \frac{1}{Z_c} (\mathbf{n} \times jw \mathrm{grad} \alpha \phi^i) \mathrm{d}s \quad (11) $$
## 2.3 Impressing current I of external circuit

For the classical full finite element model, the current passing through the coil can be given by [5]:

$$I = \int_{\Omega_{\text{coil}}} \mathbf{grad} \alpha \mathbf{J} dv \qquad (12)$$

where $\Omega_{\text{coil}}$ is the coil volumic domain. In the SIBC method, the current can be computed as:

$$I = \int_{\Gamma_{c2}} \mathbf{grad} \alpha \mathbf{K} ds \qquad (13)$$

where $\mathbf{K} = \mathbf{n} \times \mathbf{H}$ is surface current. By using (3), one has:

$$\begin{aligned} I &= \int_{\Gamma_{c2}} \mathbf{grad} \alpha \frac{1}{Z_c} (\mathbf{n} \times jw(\mathbf{A} + \mathbf{grad} \phi)) \times \mathbf{n} ds \\ &= \int_{\Gamma_{c2}} (\mathbf{n} \times \mathbf{grad} \alpha) \frac{1}{Z_c} (\mathbf{n} \times jw(\mathbf{A} + \mathbf{grad} \phi)) ds \end{aligned} \qquad (14)$$

where $\phi$ is defined as in (9). In the impressed I formulation, $\phi^i$ becomes unknown and the equation (14) has to be added in the equation system (10) and (11).

## 2.4 Matrix formulation

The matrix form of the discrete $\mathbf{A}-\phi$ formulation with a voltage $V$ impressed as external circuit condition reads:

$$(\mathbf{R}^t \mathbf{M}_{ff}^{1/\mu} \mathbf{R} + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{A}) + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \phi_S = -\frac{jw}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i \qquad (15)$$

$$\frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{A} + \frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \phi_S = -\frac{jw}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i \qquad (16)$$

where $\mathbf{R}$, $\mathbf{G}$ are respectively the discrete counterparts of curl, and grad operators and $t$ denotes the transpose operator. The subscript $S$ means "defined on the surface mesh of the conductors". The matrix form of the discrete $\mathbf{A}-\phi$ formulation with a voltage $I$ impressed reads:

$$(\mathbf{R}^t \mathbf{M}_{ff}^{1/\mu} \mathbf{R} + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{A}) + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \phi_S + \frac{jw}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i = 0 \qquad (17)$$

$$\frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{A} + \frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \phi_S + \frac{jw}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i = 0 \qquad (18)$$

$$\frac{jw}{Z_c} \alpha_S^t \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{A} + \frac{jw}{Z_c} \alpha_S^t \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \phi_S + \frac{jw}{Z_c} \alpha_S^t \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i = I \qquad (19)$$

The coefficients, denoted $m_{ff}^{1/\mu}$ and $m_{ees}^i$, of the matrix $\mathbf{M}_{ff}^{1/\mu}$ and $\mathbf{M}_{ees}^i$ are computed as:

$$m_{ff}^{1/\mu} = \int_{\Omega} \frac{1}{\mu} \mathbf{w}^{f_i} \mathbf{w}^{f_j} dv \qquad (20)$$

$$m_{ees} = \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i})(\mathbf{n} \times \mathbf{w}^{e_j}) ds = \int_{\Gamma_c} \mathbf{w}_S^{e_i} \mathbf{w}_S^{e_j} ds \qquad (21)$$

$\mathbf{w}^f$ is the facet shape function and $\mathbf{w}^e$ the edge shape function.

## 2.5 A validation of electromagnetic model

The numerical models presented in the previous sections are validated by means of a simple test case for which exact analytical solution is available in [7]. We consider the case of a long cylindrical conductor carrying a alternative current and compute the current density distribution in its body. The description of this test case is given in the Fig. 2. Three models are compared :

* Exact analytical solution (Reference)
## 2.3 Impressing current I of external circuit

For the classical full finite element model, the current passing through the coil can be given by [5]:

$$I = \int_{\Omega_{\text{coil}}} \mathbf{grad} \alpha \mathbf{J} dv \qquad (12)$$

where $\Omega_{\text{coil}}$ is the coil volumic domain. In the SIBC method, the current can be computed as:

$$I = \int_{\Gamma_{c2}} \mathbf{grad} \alpha \mathbf{K} ds \qquad (13)$$

where $\mathbf{K} = \mathbf{n} \times \mathbf{H}$ is surface current. By using (3), one has:

$$\begin{aligned} I &= \int_{\Gamma_{c2}} \mathbf{grad} \alpha \frac{1}{Z_c} (\mathbf{n} \times jw(\mathbf{A} + \mathbf{grad} \phi)) \times \mathbf{n} ds \\ &= \int_{\Gamma_{c2}} (\mathbf{n} \times \mathbf{grad} \alpha) \frac{1}{Z_c} (\mathbf{n} \times jw(\mathbf{A} + \mathbf{grad} \phi)) ds \end{aligned} \qquad (14)$$

where $\phi$ is defined as in (9). In the impressed I formulation, $\phi^i$ becomes unknown and the equation (14) has to be added in the equation system (10) and (11).

## 2.4 Matrix formulation

The matrix form of the discrete $\mathbf{A}-\phi$ formulation with a voltage $V$ impressed as external circuit condition reads:

$$(\mathbf{R}^t \mathbf{M}_{ff}^{1/\mu} \mathbf{R} + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{A}) + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \phi_S = -\frac{jw}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i \qquad (15)$$

$$\frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{A} + \frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \phi_S = -\frac{jw}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i \qquad (16)$$

where $\mathbf{R}$, $\mathbf{G}$ are respectively the discrete counterparts of curl, and grad operators and $t$ denotes the transpose operator. The subscript $S$ means "defined on the surface mesh of the conductors". The matrix form of the discrete $\mathbf{A}-\phi$ formulation with a voltage $I$ impressed reads:

$$(\mathbf{R}^t \mathbf{M}_{ff}^{1/\mu} \mathbf{R} + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{A}) + \frac{j\omega}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \phi_S + \frac{jw}{Z_c} \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i = 0 \qquad (17)$$

$$\frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{A} + \frac{j\omega}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \phi_S + \frac{jw}{Z_c} \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i = 0 \qquad (18)$$

$$\frac{jw}{Z_c} \alpha_S^t \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{A} + \frac{jw}{Z_c} \alpha_S^t \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \phi_S + \frac{jw}{Z_c} \alpha_S^t \mathbf{G}_S^t \mathbf{M}_{ees} \mathbf{G}_S \alpha_S \phi^i = I \qquad (19)$$

The coefficients, denoted $m_{ff}^{1/\mu}$ and $m_{ees}^i$, of the matrix $\mathbf{M}_{ff}^{1/\mu}$ and $\mathbf{M}_{ees}^i$ are computed as:

$$m_{ff}^{1/\mu} = \int_{\Omega} \frac{1}{\mu} \mathbf{w}^{f_i} \mathbf{w}^{f_j} dv \qquad (20)$$

$$m_{ees} = \int_{\Gamma_c} (\mathbf{n} \times \mathbf{w}^{e_i})(\mathbf{n} \times \mathbf{w}^{e_j}) ds = \int_{\Gamma_c} \mathbf{w}_S^{e_i} \mathbf{w}_S^{e_j} ds \qquad (21)$$

$\mathbf{w}^f$ is the facet shape function and $\mathbf{w}^e$ the edge shape function.

## 2.5 A validation of electromagnetic model

The numerical models presented in the previous sections are validated by means of a simple test case for which exact analytical solution is available in [7]. We consider the case of a long cylindrical conductor carrying a alternative current and compute the current density distribution in its body. The description of this test case is given in the Fig. 2. Three models are compared :

* Exact analytical solution (Reference)
<figure><img src="image_3.png" /><figcaption>Figure 2: Cylindrical conductor carrying an alternative current. The SIBC can be applied on its surface as shown in the figure. Numerical entries are $r_c = 5$mm, $\sigma = 60MS/m$ (copper), $I = 100A$.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3: Comparison of the real and imaginary parts of current in the conductor in function of the radius. Case with $fr = 50kHz$. The ratio skin depth over the radius $r_c/\delta = 17.2$.</figcaption></figure>

* SIBC finite element with impressed current with Leontovich's impedance (SIBC-$Z_c^{Leontovich}$)

* SIBC finite element with impressed current with Mitzner's impedance (SIBC-$Z_c^{Mitzner}$)

Note that for methods using SIBC in which only the current density on the surface of the conductor is computed, the current distribution in the body is computed by $J(r) = J_S e^{-(1+j)\frac{r-r_c}{\delta}}$ where $J_S$ is the current density on the surface and $r_c$ the conductor radius. The Figs. 3 and 4 give the comparisons at various frequencies. The results given by the impressed voltage formulation are identical as the impressed current formulation and are not shown here for better visibility of the figures.

It can be observed that at high frequency (when the ratio between the radius of the conductor and the skin depth is high) very good accuracy can be obtained with the SIBC methods. At lower frequency, the effect of the curvature of the surface seems become significant. In this case, the use of the Mitzner's impedance can help to correct the error.
<figure><img src="image_4.png" /><figcaption>Figure 2: Cylindrical conductor carrying an alternative current. The SIBC can be applied on its surface as shown in the figure. Numerical entries are $r_c = 5$mm, $\sigma = 60MS/m$ (copper), $I = 100A$.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 3: Comparison of the real and imaginary parts of current in the conductor in function of the radius. Case with $fr = 50kHz$. The ratio skin depth over the radius $r_c/\delta = 17.2$.</figcaption></figure>

* SIBC finite element with impressed current with Leontovich's impedance (SIBC-$Z_c^{Leontovich}$)

* SIBC finite element with impressed current with Mitzner's impedance (SIBC-$Z_c^{Mitzner}$)

Note that for methods using SIBC in which only the current density on the surface of the conductor is computed, the current distribution in the body is computed by $J(r) = J_S e^{-(1+j)\frac{r-r_c}{\delta}}$ where $J_S$ is the current density on the surface and $r_c$ the conductor radius. The Figs. 3 and 4 give the comparisons at various frequencies. The results given by the impressed voltage formulation are identical as the impressed current formulation and are not shown here for better visibility of the figures.

It can be observed that at high frequency (when the ratio between the radius of the conductor and the skin depth is high) very good accuracy can be obtained with the SIBC methods. At lower frequency, the effect of the curvature of the surface seems become significant. In this case, the use of the Mitzner's impedance can help to correct the error.
<figure><img src="image_5.png" /><figcaption>Figure 4: Comparison of the real and imaginary parts of current in the conductor in function of the radius. Case with $fr = 5kHz$. The ratio the radius over the skin depth $r_c/\delta = 5.4$.</figcaption></figure>

## 3 Thermal Problem

In the thermal problem, only the workpiece domain $\Omega_c$ is considered and on its boundary $\Gamma_{c1}$, a Neumann type boundary condition is applied. The thermal problem is defined by the heat transfer equation:

$$ \rho C_p \frac{\partial T}{\partial t} + \text{div}(-\lambda \mathbf{grad} T) = P_S \quad (22) $$

with the convection type boundary condition:

$$ -\lambda \frac{\partial T}{\partial \mathbf{n}} = h(T|_{\Gamma_{c1}} - T_{\infty}) \quad (23) $$

where $T_{\infty}$ is the ambient temperature. The parameter $h$ takes into account the natural convection and the radiation on the surface of the workpiece for small variation of temperature [2]. $P_S$ is the surface power density ($W/m^2$), $C_p$ the specific heat, $\rho$ the specific mass of the material and $\lambda$ the thermal conductivity. The surface power density (flux type thermal source) is computed as:

$$ P_S = \frac{\delta}{2} \sigma \mathbf{E}_S \hat{\mathbf{E}}_S \quad (24) $$

where $\hat{\mathbf{E}}_S$ is complex conjugate of the electric field $\mathbf{E}_S$ computed on the surface of the workpiece. The time derivative of the temperature is approximated for each time step $k$ by $(T_{k+1} - T_k)/\Delta t$ with $\Delta t$ small enough ($\Delta t$ is generally chosen far lower than the thermal time constant $\Delta t \ll \tau = \rho C_p / h$) and $T_k$ is the temperature rise at the time step $k$ with respect to the initial temperature. The weak formulation of the thermal problem reads:

$$ \int_{\Omega_c} w^{n_i} \rho C_p \frac{T_{k+1}}{\Delta t} dv + \int_{\Omega_c} \mathbf{grad} w^{n_i} \lambda \mathbf{grad} T_{k+1} dv + \int_{\Gamma_{c1}} w_S^{n_i} h T_{k+1} ds = \int_{\Gamma_{c1}} w_S^{n_i} P_S ds + \int_{\Omega_c} w^{n_i} \rho C_p \frac{T_k}{\Delta t} dv \quad (25) $$

which is written in matrix form as:

$$ \left( \mathbf{M}_{nn}^{[\rho C_p / \Delta t]} + \mathbf{G}^t \mathbf{M}_{ee}^{[\lambda]} \mathbf{G} + \mathbf{M}_{nn_S}^{[h]} \right) \mathbf{T}_{k+1} = \mathbf{M}_{n_S}^{[P_S]} + \mathbf{M}_{nn}^{[\rho C_p / \Delta t]} \mathbf{T}_k \quad (26) $$
<figure><img src="image_6.png" /><figcaption>Figure 4: Comparison of the real and imaginary parts of current in the conductor in function of the radius. Case with $fr = 5kHz$. The ratio the radius over the skin depth $r_c/\delta = 5.4$.</figcaption></figure>

## 3 Thermal Problem

In the thermal problem, only the workpiece domain $\Omega_c$ is considered and on its boundary $\Gamma_{c1}$, a Neumann type boundary condition is applied. The thermal problem is defined by the heat transfer equation:

$$ \rho C_p \frac{\partial T}{\partial t} + \text{div}(-\lambda \mathbf{grad} T) = P_S \quad (22) $$

with the convection type boundary condition:

$$ -\lambda \frac{\partial T}{\partial \mathbf{n}} = h(T|_{\Gamma_{c1}} - T_{\infty}) \quad (23) $$

where $T_{\infty}$ is the ambient temperature. The parameter $h$ takes into account the natural convection and the radiation on the surface of the workpiece for small variation of temperature [2]. $P_S$ is the surface power density ($W/m^2$), $C_p$ the specific heat, $\rho$ the specific mass of the material and $\lambda$ the thermal conductivity. The surface power density (flux type thermal source) is computed as:

$$ P_S = \frac{\delta}{2} \sigma \mathbf{E}_S \hat{\mathbf{E}}_S \quad (24) $$

where $\hat{\mathbf{E}}_S$ is complex conjugate of the electric field $\mathbf{E}_S$ computed on the surface of the workpiece. The time derivative of the temperature is approximated for each time step $k$ by $(T_{k+1} - T_k)/\Delta t$ with $\Delta t$ small enough ($\Delta t$ is generally chosen far lower than the thermal time constant $\Delta t \ll \tau = \rho C_p / h$) and $T_k$ is the temperature rise at the time step $k$ with respect to the initial temperature. The weak formulation of the thermal problem reads:

$$ \int_{\Omega_c} w^{n_i} \rho C_p \frac{T_{k+1}}{\Delta t} dv + \int_{\Omega_c} \mathbf{grad} w^{n_i} \lambda \mathbf{grad} T_{k+1} dv + \int_{\Gamma_{c1}} w_S^{n_i} h T_{k+1} ds = \int_{\Gamma_{c1}} w_S^{n_i} P_S ds + \int_{\Omega_c} w^{n_i} \rho C_p \frac{T_k}{\Delta t} dv \quad (25) $$

which is written in matrix form as:

$$ \left( \mathbf{M}_{nn}^{[\rho C_p / \Delta t]} + \mathbf{G}^t \mathbf{M}_{ee}^{[\lambda]} \mathbf{G} + \mathbf{M}_{nn_S}^{[h]} \right) \mathbf{T}_{k+1} = \mathbf{M}_{n_S}^{[P_S]} + \mathbf{M}_{nn}^{[\rho C_p / \Delta t]} \mathbf{T}_k \quad (26) $$
<figure><img src="image_6.png" /><figcaption>Figure 5: Induction thermography NDT test case.</figcaption></figure>

The coefficients of the matrix $M_{nn}^{[\rho C_p/\Delta t]}$, $M_{ee}^{[\lambda]}$, $M_{nn_S}^{[h]}$, $M_{n_S}^{[P_S]}$ and $M_{nn}^{[\rho C_p/\Delta t]}$ can be calculated as below:

$$m_{nn}^{[\rho C_p/\Delta t]} = \int_{\Omega_c} \frac{\rho C_p}{\Delta t} w^{n_i} w^{n_j} dv \qquad (27)$$

$$m_{ee}^{[\lambda]} = \int_{\Omega_c} \mathbf{w}^{e_i} \lambda \mathbf{w}^{e_j} dv \qquad (28)$$

$$m_{nn_S}^{[h]} = \int_{\Gamma_{c1}} h w_S^{n_i} w_S^{n_j} ds \qquad (29)$$

$$m_{n_S}^{[P_S]} = \int_{\Gamma_{c1}} w_S^{n_i} P_S ds \qquad (30)$$

# 4 Application to induction thermography simulation

In this section, the electromagnetic thermal model developed in the previous sections is applied to simulate an induction thermography NDT test case. The configuration is shown in the Fig. 5. The workpiece has four thin open cracks. An AC 10V voltage of 10kHz is impressed between the ends of the coil. The workpiece is heated during 1 second. The temperature evolution is observed on the upper surface of the workpiece. Numerical entry data are given in the Table. 1. Two following methods are compared:

* The classical full finite element method (Full-FEM)

* The SIBC finite element with impressed voltage with Leontovich's impedance (SIBC-FEM)

The ratio thickness over skindepth in the workpiece and the coil is equal to 28 and 7.7, respectively. In Fig. 6 the distribution of induced power density is compared between full-FEM and SIBC-FEM models. The evolution of eddy-current over the x-line and the y-line on the surface of the workpiece (cf. Fig. 5) is also compared as shown in Fig. 7. The x-line is chosen to pass the cracks' tips where a high variation of eddy-current can be observed. The module of the magnetic flux density in the middle of the air-gap right above the x- and y-lines is compared in Fig. 8. These comparisons show good accordance of local quantities.

It is nothing that the presented model is limited to linear cases. In the literature, to deal with nonlinear behavior of the magnetic permeability in full FEM while keeping the use of time harmonic model in order to reduce computation time, one can use the effective permeability or equivalent reluctivity approaches [8]. The same idea can also be applied to determine the dependency of surface impedance on the
<figure><img src="image_7.png" /><figcaption>Figure 5: Induction thermography NDT test case.</figcaption></figure>

The coefficients of the matrix $M_{nn}^{[\rho C_p/\Delta t]}$, $M_{ee}^{[\lambda]}$, $M_{nn_s}^{[h]}$, $M_{n_s}^{[P_s]}$ and $M_{nn}^{[\rho C_p/\Delta t]}$ can be calculated as below:

$$m_{nn}^{[\rho C_p/\Delta t]} = \int_{\Omega_c} \frac{\rho C_p}{\Delta t} w^{n_i} w^{n_j} dv \qquad (27)$$

$$m_{ee}^{[\lambda]} = \int_{\Omega_c} \mathbf{w}^{e_i} \lambda \mathbf{w}^{e_j} dv \qquad (28)$$

$$m_{nn_s}^{[h]} = \int_{\Gamma_{c1}} h w_S^{n_i} w_S^{n_j} ds \qquad (29)$$

$$m_{n_s}^{[P_s]} = \int_{\Gamma_{c1}} w_S^{n_i} P_s ds \qquad (30)$$

# 4 Application to induction thermography simulation

In this section, the electromagnetic thermal model developed in the previous sections is applied to simulate an induction thermography NDT test case. The configuration is shown in the Fig. 5. The workpiece has four thin open cracks. An AC 10V voltage of 10kHz is impressed between the ends of the coil. The workpiece is heated during 1 second. The temperature evolution is observed on the upper surface of the workpiece. Numerical entry data are given in the Table. 1. Two following methods are compared:

* The classical full finite element method (Full-FEM)

* The SIBC finite element with impressed voltage with Leontovich's impedance (SIBC-FEM)

The ratio thickness over skindepth in the workpiece and the coil is equal to 28 and 7.7, respectively. In Fig. 6 the distribution of induced power density is compared between full-FEM and SIBC-FEM models. The evolution of eddy-current over the x-line and the y-line on the surface of the workpiece (cf. Fig. 5) is also compared as shown in Fig. 7. The x-line is chosen to pass the cracks' tips where a high variation of eddy-current can be observed. The module of the magnetic flux density in the middle of the air-gap right above the x- and y-lines is compared in Fig. 8. These comparisons show good accordance of local quantities.

It is nothing that the presented model is limited to linear cases. In the literature, to deal with nonlinear behavior of the magnetic permeability in full FEM while keeping the use of time harmonic model in order to reduce computation time, one can use the effective permeability or equivalent reluctivity approaches [8]. The same idea can also be applied to determine the dependency of surface impedance on the
<table><caption>Table 1: Numerical data for simulation.</caption><thead><tr><th colspan="3">Coil (Copper)</th></tr></thead><tbody><tr><td>Electrical conductivity</td><td>60</td><td>MS/m</td></tr><tr><td>Relative permeability</td><td>1</td><td>—</td></tr><tr><td>Impressed voltage</td><td>10</td><td>V</td></tr><tr><td>Frequency</td><td>10</td><td>kHz</td></tr><tr><td>Cross-section</td><td>5 × 5</td><td>mm × mm</td></tr><tr><td>Skin-depth (at 10kHz)</td><td>0.65</td><td>mm</td></tr><tr><th colspan="3">Workpiece (Structural steel)</th></tr><tr><td>Electrical conductivity</td><td>4</td><td>MS/m</td></tr><tr><td>Relative permeability</td><td>200</td><td>—</td></tr><tr><td>Thickness</td><td>5</td><td>mm</td></tr><tr><td>Skin-depth (at 10kHz)</td><td>0.18</td><td>mm</td></tr><tr><td>Specific heat capacity</td><td>475</td><td>J.kg<sup>-1</sup>.K<sup>-1</sup></td></tr><tr><td>Thermal conductivity</td><td>45</td><td>W.m<sup>-1</sup>.K<sup>-1</sup></td></tr><tr><td>Specific mass</td><td>7800</td><td>kg.m<sup>-3</sup></td></tr><tr><th colspan="3">Cracks</th></tr><tr><td>Size</td><td>0.5 × 7</td><td>mm × mm</td></tr><tr><td>Thickness</td><td>5</td><td>mm</td></tr></tbody></table>

<table><caption>Table 2: Comparison between Full-FEM and SIBC-FEM system.</caption><thead><tr><th></th><th>Full-FEM</th><th>SIBC-FEM</th></tr></thead><tbody><tr><td>Total number of elements</td><td>768238</td><td>472818</td></tr><tr><td>Total number of elements in the workpiece and the coil</td><td>295420</td><td>0</td></tr><tr><td>Total number of edge unknowns</td><td>1505300</td><td>964135</td></tr><tr><td>Total number of node unknowns</td><td>166443</td><td>36487</td></tr><tr><td>Matrix assembly time (s)</td><td>86</td><td>31</td></tr><tr><td>Matrix resolution time (s)</td><td>3374</td><td>508</td></tr></tbody></table>

peak value of magnetic field [9, 10]. In both cases, although high harmonic components are neglected, calculation results for induced power were satisfactory. The approaches can thus be applied to NDT computation.

The Table. 2 gives comparisons between the Full-FEM system and the SIBC-FEM system. Note that the resolution time for thermal problem for the two model is practically identical. The resolution time in the Table. 2 is given for the electromagnetic problem.

The Figs. 9 and 10 show the temperature distribution at the end of the heating phase ($t = 1s$) obtained with Full-FEM and SIBC-FEM models, respectively. A good accordance can be observed. The relative deviation is given in the Fig. 11 which is computed as:

$$ \delta T = \frac{T_{t=1s}^{SIBC-FEM} - T_{t=1s}^{Full-FEM}}{\max(T_{t=1s}^{Full-FEM})} \times 100 \quad (31) $$

In the present test case, the deviation is not higher than 7%.

# 5 Conclusion

The SIBC can be used to alleviate 3D finite element models for problem with presence of massive coils and workpieces with open cracks. In the present test case, the computation time is 6.6 times lower than the classical model. The accuracy and the usefulness of the method for the design of the induction thermography NDT technique have been shown with acceptable deviation compared with a full FEM model.
<table><caption>Table 1: Numerical data for simulation.</caption><tbody><tr><th colspan="3">Coil (Copper)</th></tr><tr><td>Electrical conductivity</td><td>60</td><td>MS/m</td></tr><tr><td>Relative permeability</td><td>1</td><td>—</td></tr><tr><td>Impressed voltage</td><td>10</td><td>V</td></tr><tr><td>Frequency</td><td>10</td><td>kHz</td></tr><tr><td>Cross-section</td><td>5 × 5</td><td>mm × mm</td></tr><tr><td>Skin-depth (at 10kHz)</td><td>0.65</td><td>mm</td></tr><tr><th colspan="3">Workpiece (Structural steel)</th></tr><tr><td>Electrical conductivity</td><td>4</td><td>MS/m</td></tr><tr><td>Relative permeability</td><td>200</td><td>—</td></tr><tr><td>Thickness</td><td>5</td><td>mm</td></tr><tr><td>Skin-depth (at 10kHz)</td><td>0.18</td><td>mm</td></tr><tr><td>Specific heat capacity</td><td>475</td><td>J.kg<sup>-1</sup>.K<sup>-1</sup></td></tr><tr><td>Thermal conductivity</td><td>45</td><td>W.m<sup>-1</sup>.K<sup>-1</sup></td></tr><tr><td>Specific mass</td><td>7800</td><td>kg.m<sup>-3</sup></td></tr><tr><th colspan="3">Cracks</th></tr><tr><td>Size</td><td>0.5 × 7</td><td>mm × mm</td></tr><tr><td>Thickness</td><td>5</td><td>mm</td></tr></tbody></table>

<table><caption>Table 2: Comparison between Full-FEM and SIBC-FEM system.</caption><thead><tr><th></th><th>Full-FEM</th><th>SIBC-FEM</th></tr></thead><tbody><tr><td>Total number of elements</td><td>768238</td><td>472818</td></tr><tr><td>Total number of elements in the workpiece and the coil</td><td>295420</td><td>0</td></tr><tr><td>Total number of edge unknowns</td><td>1505300</td><td>964135</td></tr><tr><td>Total number of node unknowns</td><td>166443</td><td>36487</td></tr><tr><td>Matrix assembly time (s)</td><td>86</td><td>31</td></tr><tr><td>Matrix resolution time (s)</td><td>3374</td><td>508</td></tr></tbody></table>

peak value of magnetic field [9, 10]. In both cases, although high harmonic components are neglected, calculation results for induced power were satisfactory. The approaches can thus be applied to NDT computation.

The Table. 2 gives comparisons between the Full-FEM system and the SIBC-FEM system. Note that the resolution time for thermal problem for the two model is practically identical. The resolution time in the Table. 2 is given for the electromagnetic problem.

The Figs. 9 and 10 show the temperature distribution at the end of the heating phase ($t = 1s$) obtained with Full-FEM and SIBC-FEM models, respectively. A good accordance can be observed. The relative deviation is given in the Fig. 11 which is computed as:

$$ \delta T = \frac{T_{t=1s}^{SIBC-FEM} - T_{t=1s}^{Full-FEM}}{\max(T_{t=1s}^{Full-FEM})} \times 100 \quad (31) $$

In the present test case, the deviation is not higher than 7%.

## 5 Conclusion

The SIBC can be used to alleviate 3D finite element models for problem with presence of massive coils and workpieces with open cracks. In the present test case, the computation time is 6.6 times lower than the classical model. The accuracy and the usefulness of the method for the design of the induction thermography NDT technique have been shown with acceptable deviation compared with a full FEM model.
<figure><img src="image_7.png" /><figcaption>Figure 6: Comparison of induced power density (Full-FEM: left, SIBC-FEM: right).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 7: Eddy-current over x- and y- lines.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 6: Comparison of induced power density (Full-FEM: left, SIBC-FEM: right).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 7: Eddy-current over x- and y- lines.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>Figure 8: Magnetic flux density over x- and y- lines.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 9: Temperature distribution on the surface of the workpiece at $t = 1s$ computed with Full-FEM model.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Figure 10: Temperature distribution on the surface of the workpiece at $t = 1s$ computed with SIBC-FEM model.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>Figure 8: Magnetic flux density over x- and y- lines.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Figure 9: Temperature distribution on the surface of the workpiece at $t = 1s$ computed with Full-FEM model.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 10: Temperature distribution on the surface of the workpiece at $t = 1s$ computed with SIBC-FEM model.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>Figure 11: Relative deviation (in %) between the two models.</figcaption></figure>

## Acknowledgment

This project has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 722134-NDTonAIR.

## References

[1] M. He, L. Zhang, W. Zheng, and Y. Feng, "Investigation on a new inducer of pulsed eddy current thermography," AIP Advances, 2016, vol. 6, no 9, p. 095221.

[2] H. K. Bui, G. Wasselynck, D. Trichet and G. Berthiau, "Performance assessment of induction thermography technique applied to carbon-fiber-reinforced polymer material," IEEE Transactions on Magnetics, vol. 51, 2015.

[3] S. M. Mimoue, J. Fouladgar, A. Chentouf and G. Develey, "A 3D impedance calculation for an induction heating system for materials with poor conductivity," IEEE Transactions on Magnetics, vol. 32, 1996.

[4] A. Desmoort, Z. De Grève, P. Dular, C. Geuzaine and O. Deblecker, "Surface Impedance Boundary Condition With Circuit Coupling for the 3-D Finite-Element Modeling of Wireless Power Transfer," IEEE Transactions on Magnetics, vol. 53, 2017.

[5] T. Henneron, S. Clénet, and F. Piriou, "Calculation of extra copper losses with imposed current magnetodynamic formulations," IEEE Transactions on Magnetics, vol. 42, no 4, p. 767-770, 2006.

[6] S. Yuferev, and N. Ida, "Selection of the surface impedance boundary conditions for a given problem," IEEE Transactions on Magnetics, vol. 35, 1999.

[7] Ramo, Simon, John R. Whinnery, and Theodore Van Duzer, "Fields and waves in communication electronics," John Wiley & Sons, p. 180-186, 1994.

[8] E. Vassent, G. Meunier and J. C. Sabonnadier, "Simulation of induction machine operation using complex magnetodynamic finite elements," IEEE Transactions on Magnetics, vol. 25, no. 4, 1989.

[9] C. Guerin, G. Meunier and G. Tanneau, "Surface impedance for 3D nonlinear eddy current problems-application to loss computation in transformers," IEEE Transactions on Magnetics, vol. 32, no 3, 1996.

[10] J. Nerg, and J. Partanen, "A simplified FEM based calculation model for 3-D induction heating problems using surface impedance formulations," IEEE Transactions on Magnetics, vol. 37, no 5, 2001.
<figure><img src="image_13.png" /><figcaption>Figure 11: Relative deviation (in %) between the two models.</figcaption></figure>

## Acknowledgment

This project has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 722134-NDTonAIR.

## References

[1] M. He, L. Zhang, W. Zheng, and Y. Feng, "Investigation on a new inducer of pulsed eddy current thermography," AIP Advances, 2016, vol. 6, no 9, p. 095221.

[2] H. K. Bui, G. Wasselynck, D. Trichet and G. Berthiau, "Performance assessment of induction thermography technique applied to carbon-fiber-reinforced polymer material," IEEE Transactions on Magnetics, vol. 51, 2015.

[3] S. M. Mimoue, J. Fouladgar, A. Chentouf and G. Develey, "A 3D impedance calculation for an induction heating system for materials with poor conductivity," IEEE Transactions on Magnetics, vol. 32, 1996.

[4] A. Desmoort, Z. De Grève, P. Dular, C. Geuzaine and O. Deblecker, "Surface Impedance Boundary Condition With Circuit Coupling for the 3-D Finite-Element Modeling of Wireless Power Transfer," IEEE Transactions on Magnetics, vol. 53, 2017.

[5] T. Henneron, S. Clénet, and F. Piriou, "Calculation of extra copper losses with imposed current magnetodynamic formulations," IEEE Transactions on Magnetics, vol. 42, no 4, p. 767-770, 2006.

[6] S. Yuferev, and N. Ida, "Selection of the surface impedance boundary conditions for a given problem," IEEE Transactions on Magnetics, vol. 35, 1999.

[7] Ramo, Simon, John R. Whinnery, and Theodore Van Duzer, "Fields and waves in communication electronics," John Wiley & Sons, p. 180-186, 1994.

[8] E. Vassent, G. Meunier and J. C. Sabonnadier, "Simulation of induction machine operation using complex magnetodynamic finite elements," IEEE Transactions on Magnetics, vol. 25, no. 4, 1989.

[9] C. Guerin, G. Meunier and G. Tanneau, "Surface impedance for 3D nonlinear eddy current problems-application to loss computation in transformers," IEEE Transactions on Magnetics, vol. 32, no 3, 1996.

[10] J. Nerg, and J. Partanen, "A simplified FEM based calculation model for 3-D induction heating problems using surface impedance formulations," IEEE Transactions on Magnetics, vol. 37, no 5, 2001.