# Semi-analytical modeling of eddy current inspection of stratified and homogenized composite materials

Houssem Chebbi, Denis Prémel

► To cite this version:

Houssem Chebbi, Denis Prémel. Semi-analytical modeling of eddy current inspection of stratified and homogenized composite materials. 10th International Symposium on NDT in Aerospace, Oct 2018, Dresden, Germany (AeroNDT 2018), Oct 2019, Dresden, Germany. cea-04442312

HAL Id: cea-04442312

https://cea.hal.science/cea-04442312

Submitted on 6 Feb 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
# Semi-analytical modeling of eddy current inspection of stratified and homogenized composite materials

Houssem.CHEBBI¹, Denis.PRÉMEL¹

¹ Commissariat à l'Energie Atomique et aux énergies alternatives (CEA) LIST, Gif-sur-Yvette cedex, France

Contact e-mail: houssem.chebbi@cea.fr

## 1 abstract

In the present paper, we deal with the development of a semi-analytical model for the fast computation of a quasi-static field induced by a 3D eddy current probe in a flawless composite material. After a prior phase of homogenization, the workpiece is assumed to be homogenized and stratified. Therefore, each layer of the material is characterized by a conductivity tensor with axial or biaxial form according to the main orientation of the fibers. The approach is based on a modal decomposition of the Electromagnetic (EM) field in each layer. To take into account the wave propagation through the multilayered structure and boundary conditions at each interface, a recursive scattering matrix (S-matrix) algorithm is implemented. Numerical results provided by the model are compared to other simulated data obtained by a Finite Element commercial software.

Keywords:

computational electromagnetics, semi-analytic model, composite material, scattering matrix, homogenization.

## 2 introduction

Composite materials, such as Carbon Fiber Reinforced Polymer (CFRP), are increasingly used in the aeronautic industry to reduce the weight of aeronautical structures. Therefore, the inspection of these complex materials for characterization or the detection of potential defects is a major challenge due to their complex physical behavior, to ensure their good performance in different operations.

In one hand, Eddy Current testing technique (ECT) is sensitive to the variations of the conductivity of the material, it can be applied for inspection (crack detection, thickness measurements...) and characterization (bulk conductivity, fiber orientation...). In the other hand, the interactions between electromagnetic waves and composite materials are highly depending on the microstructure of the material itself. This gives rise to challenging modeling issues. Many homogenization approaches were proposed [1, 2] to overcome high anisotropy, random positioning of the fibers. Therefore, each layer of the stratified media is characterized by an equivalent conductivity tensor and equivalent permittivity tensor. In the same framework, many research works [3, 4, 5] present a variety of EC techniques providing interesting results in characterizing homogenized composite material.
# Semi-analytical modeling of eddy current inspection of stratified and homogenized composite materials

Houssem.CHEBBI¹, Denis.PRÉMEL¹

¹ Commissariat à l'Energie Atomique et aux énergies alternatives (CEA) LIST, Gif-sur-Yvette cedex, France

Contact e-mail: houssem.chebbi@cea.fr

## 1 abstract

In the present paper, we deal with the development of a semi-analytical model for the fast computation of a quasi-static field induced by a 3D eddy current probe in a flawless composite material. After a prior phase of homogenization, the workpiece is assumed to be homogenized and stratified. Therefore, each layer of the material is characterized by a conductivity tensor with axial or biaxial form according to the main orientation of the fibers. The approach is based on a modal decomposition of the Electromagnetic (EM) field in each layer. To take into account the wave propagation through the multilayered structure and boundary conditions at each interface, a recursive scattering matrix (S-matrix) algorithm is implemented. Numerical results provided by the model are compared to other simulated data obtained by a Finite Element commercial software.

Keywords:

computational electromagnetics, semi-analytic model, composite material, scattering matrix, homogenization.

## 2 introduction

Composite materials, such as Carbon Fiber Reinforced Polymer (CFRP), are increasingly used in the aeronautic industry to reduce the weight of aeronautical structures. Therefore, the inspection of these complex materials for characterization or the detection of potential defects is a major challenge due to their complex physical behavior, to ensure their good performance in different operations.

In one hand, Eddy Current testing technique (ECT) is sensitive to the variations of the conductivity of the material, it can be applied for inspection (crack detection, thickness measurements...) and characterization (bulk conductivity, fiber orientation...). On the other hand, the interactions between electromagnetic waves and composite materials are highly depending on the microstructure of the material itself. This gives rise to challenging modeling issues. Many homogenization approaches were proposed [1, 2] to overcome high anisotropy, random positioning of the fibers. Therefore, each layer of the stratified media is characterized by an equivalent conductivity tensor and equivalent permittivity tensor. In the same framework, many research works [3, 4, 5] present a variety of EC techniques providing interesting results in characterizing homogenized composite material.
In general, analytical models are restricted to cases with simple geometries, although providing exact solutions. In contrast, purely numerical methods are so general and used for complex problems, however, considering the computing time and memory storage, they are quite expensive. In our modeling approach, we consider a compromise that we call semi-analytical (SA) model, to insure accuracy, efficiency and rapidity. Our semi-analytic modeling framework consists in writing Maxwell's equations in each layer, in the Fourier domain in order to obtain an algebraic form. Moreover, the EM field is expanded as sum of Eigen functions. Applying boundary conditions (BC) leads us to estimate the coefficients of the expansion and to reconstruct the EM field in any region. A stable recursive algorithm like the S-matrix one, is implemented in order to take into account the wave propagation through the entire structure and BC at each interface.

the paper is organized as follows, first we describe the formalism used for the SA model, then we illustrate the keystone of the approach which is the modal decomposition of the EM field in both isotropic and anisotropic media. Next we introduce the S-matrix algorithm used for stratified media. Finally, the numerical results of the model are compared to those of a FE commercial software and results from literature.

# 3 semi-analytical model framework

## 3.1 summary

Our semi-analytical modeling framework is based on previous research works [6, 7, 8]. This work constitutes the preliminary task for the development of a complete numerical model, which has the capability to address stratified and anisotropic layered media with geometrical complexity. As a start point, we adopt the case of homogenized and anisotropic semi-infinite media excited by a 3D eddy current probe placed in air. Then we move to the case of one layered media and multilayered media with different orientation of the fibers. We attribute to each layer an equivalent conductivity tensor taking into account the preferred anisotropy orientation.

## 3.2 material characteristics

<figure><img src="image_1.png" /><figcaption>Figure 1: Structure of the considered CFRP</figcaption></figure>

The structure of CFRP is usually described as in Fig1. These materials are made of several unidirectional plies of about 125 µm thickness, stacked with different fiber orientations. Each unidirectional ply is constituted of carbon fibers embedded in an electrically non conductive
In general, analytical models are restricted to cases with simple geometries, although providing exact solutions. In contrast, purely numerical methods are so general and used for complex problems, however, considering the computing time and memory storage, they are quite expensive. In our modeling approach, we consider a compromise that we call semi-analytical (SA) model, to insure accuracy, efficiency and rapidity. Our semi-analytic modeling framework consists in writing Maxwell's equations in each layer, in the Fourier domain in order to obtain an algebraic form. Moreover, the EM field is expanded as sum of Eigen functions. Applying boundary conditions (BC) leads us to estimate the coefficients of the expansion and to reconstruct the EM field in any region. A stable recursive algorithm like the S-matrix one, is implemented in order to take into account the wave propagation through the entire structure and BC at each interface.

the paper is organized as follows, first we describe the formalism used for the SA model, then we illustrate the keystone of the approach which is the modal decomposition of the EM field in both isotropic and anisotropic media. Next we introduce the S-matrix algorithm used for stratified media. Finally, the numerical results of the model are compared to those of a FE commercial software and results from literature.

# 3 semi-analytical model framework

## 3.1 summary

Our semi-analytical modeling framework is based on previous research works [6, 7, 8]. This work constitutes the preliminary task for the development of a complete numerical model, which has the capability to address stratified and anisotropic layered media with geometrical complexity. As a start point, we adopt the case of homogenized and anisotropic semi-infinite media excited by a 3D eddy current probe placed in air. Then we move to the case of one layered media and multilayered media with different orientation of the fibers. We attribute to each layer an equivalent conductivity tensor taking into account the preferred anisotropy orientation.

## 3.2 material characteristics

<figure><img src="image_2.png" /><figcaption>Figure 1: Structure of the considered CFRP</figcaption></figure>

The structure of CFRP is usually described as in Fig1. These materials are made of several unidirectional plies of about 125 µm thickness, stacked with different fiber orientations. Each unidirectional ply is constituted of carbon fibers embedded in an electrically non conductive
matrix. The conductivity in the direction transverse to the fibers is not null since there are contacts between fibers not perfectly aligned. The conductivity tensor of each ply depends on the type, the volume fraction, and the orientation of the fibers in the ply; it varies between 5.10³ and 5.10⁴ S/m in the direction of the fibers, and between 10 and 5.10² S/m in the transverse direction [9]. When the plies are assembled, there are also contacts between the fibers of adjacent plies, giving rise to a cross ply conductivity which is of lower value being typically half the transverse value [10]. Conductivity effects are therefore introduced in the conductivity tensor by involving the longitudinal and transversal conductivities ($\sigma_l, \sigma_t$) of the plies, the orientation $\theta$ of each ply, and the cross ply conductivity $\sigma_n$ [11]:

$$ \bar{\bar{\sigma}} = \begin{bmatrix} \sigma_l \cos^2(\theta) + \sigma_t \sin^2(\theta) & \frac{\sigma_l - \sigma_t}{2} \sin(2\theta) & 0 \\ \frac{\sigma_l - \sigma_t}{2} \sin(2\theta) & \sigma_t \cos^2(\theta) + \sigma_l \sin^2(\theta) & 0 \\ 0 & 0 & \sigma_n \end{bmatrix} \quad (1) $$

Thus, we define the permittivity tensor as:

$$ \bar{\bar{\epsilon}} = \begin{bmatrix} \epsilon_{xx} & \epsilon_{xy} & 0 \\ \epsilon_{xx} & \epsilon_{xy} & 0 \\ 0 & 0 & \epsilon_{zz} \end{bmatrix} = \epsilon_r \bar{\bar{\mathbf{I}}} - \frac{i \bar{\bar{\sigma}}}{\omega \epsilon_0} \quad (2) $$

where $\epsilon_0$ is the permittivity in air. Capacity effects may be introduced by different values of the relative permittivity $\epsilon_r$. In our case we choose $\epsilon_r = 1$ but it is not required in general.

### 3.3 formulation

Assuming the absence of external electric source, working in the harmonic regime and with a non magnetic material, we focus on Maxwell-Faraday law for induction and Ampère's circulation law equations. we introduce the permeability in the air $\mu_0$, the impedance coil in the air $Z_0 = \sqrt{\frac{\mu_0}{\epsilon_0}}$ and the wave number $k = Z_0 \omega \epsilon_0$, thus:

$$ \nabla \times \mathbf{E} = -ik\mu_r Z_0 \mathbf{H} \quad (2) $$

$$ Z_0(\nabla \times \mathbf{H}) = ik\bar{\bar{\mathbf{E}}} \quad (3) $$

## 4 model decomposition (M-D)

In this section, we illustrate the way we decompose the EM field into sum of eigenmodes in both sides of a separating interface. we start with the isotropic media which is the air in our case then the anisotropic media which is the slab.

### 4.1 M-D for homogenized isotropic material

In isotropic material, we consider a Cartesian coordinate system and planar interfaces, thus we adopt Rayleigh's decomposition. We divide the field into transversal (t index) and longitudinal (l index) component ($\mathbf{E}_t, \mathbf{H}_t$) and ($\mathbf{E}_z, \mathbf{H}_z$), so that we obtain a new expression of $\Psi_z$ ($\Psi$ could be $E$ or $H$) in function of transversal components. At this stage, we use the 2D discrete Fourier transform ($TF(f) = \tilde{f}$) and we define the truncate numbers $M_u$ and $M_v$ as the numbers of modes along the 2 dimensions X and Y. In what follows and to avoid ambiguity, $f = \tilde{f}$.
matrix. The conductivity in the direction transverse to the fibers is not null since there are contacts between fibers not perfectly aligned. The conductivity tensor of each ply depends on the type, the volume fraction, and the orientation of the fibers in the ply; it varies between 5.10³ and 5.10⁴ S/m in the direction of the fibers, and between 10 and 5.10² S/m in the transverse direction [9]. When the plies are assembled, there are also contacts between the fibers of adjacent plies, giving rise to a cross ply conductivity which is of lower value being typically half the transverse value [10]. Conductivity effects are therefore introduced in the conductivity tensor by involving the longitudinal and transversal conductivities ($\sigma_l, \sigma_t$) of the plies, the orientation $\theta$ of each ply, and the cross ply conductivity $\sigma_n$ [11]:

$$ \bar{\bar{\sigma}} = \begin{bmatrix} \sigma_l \cos^2(\theta) + \sigma_t \sin^2(\theta) & \frac{\sigma_l - \sigma_t}{2} \sin(2\theta) & 0 \\ \frac{\sigma_l - \sigma_t}{2} \sin(2\theta) & \sigma_t \cos^2(\theta) + \sigma_l \sin^2(\theta) & 0 \\ 0 & 0 & \sigma_n \end{bmatrix} \quad (1) $$

Thus, we define the permittivity tensor as :

$$ \bar{\bar{\epsilon}} = \begin{bmatrix} \epsilon_{xx} & \epsilon_{xy} & 0 \\ \epsilon_{xx} & \epsilon_{xy} & 0 \\ 0 & 0 & \epsilon_{zz} \end{bmatrix} = \epsilon_r \bar{\bar{\mathbf{I}}} - \frac{i \bar{\bar{\sigma}}}{\omega \epsilon_0} \quad (2) $$

where $\epsilon_0$ is the permittivity in air. Capacity effects may be introduced by different values of the relative permittivity $\epsilon_r$. In our case we choose $\epsilon_r = 1$ but it is not required in general.

### 3.3 formulation

Assuming the absence of external electric source, working in the harmonic regime and with a non magnetic material, we focus on Maxwell-Faraday law for induction and Ampère's circulation law equations. we introduce the permeability in the air $\mu_0$, the impedance coil in the air $Z_0 = \sqrt{\frac{\mu_0}{\epsilon_0}}$ and the wave number $k = Z_0 \omega \epsilon_0$, thus:

$$ \nabla \times \mathbf{E} = -ik\mu_r Z_0 \mathbf{H} \quad (2) $$

$$ Z_0(\nabla \times \mathbf{H}) = ik\bar{\bar{\epsilon}}\mathbf{E} \quad (3) $$

## 4 model decomposition (M-D)

In this section, we illustrate the way we decompose the EM field into sum of eigenmodes in both sides of a separating interface. we start with the isotropic media which is the air in our case then the anisotropic media which is the slab.

### 4.1 M-D for homogenized isotropic material

In isotropic material, we consider a Cartesian coordinate system and planar interfaces, thus we adopt Rayleigh's decomposition. We divide the field into transversal (t index) and longitudinal (l index) component ($\mathbf{E}_t, \mathbf{H}_t$) and ($\mathbf{E}_z, \mathbf{H}_z$), so that we obtain a new expression of $\Psi_z$ ($\Psi$ could be $E$ or $H$) in function of transversal components. At this stage, we use the 2D discrete Fourier transform ($TF(f) = \tilde{f}$) and we define the truncate numbers $M_u$ and $M_v$ as the numbers of modes along the 2 dimensions X and Y. In what follows and to avoid ambiguity, $f = \tilde{f}$.
Conventionally, $\partial_i \equiv \frac{\partial}{\partial i}$ and in Fourier domain the partial derivatives become $\partial_x \equiv -i\alpha$ and $\partial_y \equiv -i\beta$. We take also the assumption of variable separation: $f(x, y, z) = f(x, y)e^{-i\gamma}$ so $\partial_z \equiv -i\gamma$. Therefore we notice that the z components verify Helmholtz's equation:

$$[\partial_z^2 + k_c^2 - \alpha^2 - \beta^2]\Psi_z = 0 \qquad (4)$$

Where $k_c^2 = k^2\mu_r\epsilon_r$. therefore, we have two scalar potentials:

$$E_z = a^{TM-}(\alpha, \beta)e^{-i\gamma z} + a^{TM+}(\alpha, \beta)e^{i\gamma z} \qquad (5)$$

$$H_z = a^{TE-}(\alpha, \beta)e^{-i\gamma z} + a^{TE+}(\alpha, \beta)e^{i\gamma z} \qquad (6)$$

Coming back to Maxwell equations, we use a $TE_z(H_z = 0)$, $TM_z(E_z = 0)$ decomposition to rewrite the transversal components of the EM field in function of modes.

$$\begin{bmatrix} E_y \\ Z_0 H_x \\ Z_0 H_y \\ E_x \end{bmatrix} = a^{TE\pm}\Phi^{TE\pm} + a^{TM\pm}\Phi^{TM\pm} \qquad (7)$$

Such:

$$\Phi^{\pm} = \frac{1}{k_c^2 - \gamma^2} \begin{bmatrix} \pm \gamma \alpha E_z + \mu_r \alpha Z_0 H_z \\ \varepsilon \beta E_z \pm \gamma \alpha Z_0 H_z \\ -\varepsilon \alpha E_z \pm \gamma \beta Z_0 H_z \\ \pm \gamma \beta E_z - \mu_r Z_0 H_z \end{bmatrix} \qquad (8)$$

## 4.2 M-D for homogenized anisotropic material

In the general case of biaxial anisotropy form of the permittivity tensor, Maxwell's equations can be written as two sets of three equations each:

$$E_{sys} = \begin{cases} \partial_y E_z - \partial_z E_y = -ik\mu_r Z_0 H_x \\ \partial_z E_x - \partial_x E_z = -ik\mu_r Z_0 H_y \\ \partial_x E_y - \partial_y E_x = -ik\mu_r Z_0 H_z \end{cases} \qquad (9)$$

$$H_{sys} = \begin{cases} \partial_y H_z - \partial_z H_y = \frac{ik}{Z_0}(\varepsilon_r^{xx} E_x + \varepsilon_r^{xy} E_y) \\ \partial_z H_x - \partial_x H_z = \frac{ik}{Z_0}(\varepsilon_r^{yy} E_y + \varepsilon_r^{yx} E_x) \\ \partial_x H_y - \partial_y H_x = \frac{ik}{Z_0} \varepsilon_r^{zz} E_z \end{cases} \qquad (10)$$

Likewise, we write $E_z = \text{fct}(H_y, H_x)$ and $H_z = \text{fct}(E_y, E_x)$. Thus we determine two decoupled matrix systems of dimension (2x2) to be solved, verified by the transversal components of the electric field in one hand and the magnetic field in the other hand presented in (11) and (13). In the same framework, Sabbagh and M. Robert in [12, 13, 14, 15] present a similar approach but do not take into account the decoupled variables, they resolve a (4x4) matrix system to describe the EM field in the layers.
Conventionally, $\partial_i \equiv \frac{\partial}{\partial i}$ and in Fourier domain the partial derivatives become $\partial_x \equiv -i\alpha$ and $\partial_y \equiv -i\beta$. We take also the assumption of variable separation: $f(x, y, z) = f(x, y)e^{-i\gamma}$ so $\partial_z \equiv -i\gamma$. Therefore we notice that the z components verify Helmholtz's equation:

$$[\partial_z^2 + k_c^2 - \alpha^2 - \beta^2]\Psi_z = 0 \qquad (4)$$

Where $k_c^2 = k^2\mu_r\epsilon_r$. therefore, we have two scalar potentials:

$$E_z = a^{TM-}(\alpha, \beta)e^{-i\gamma z} + a^{TM+}(\alpha, \beta)e^{i\gamma z} \qquad (5)$$

$$H_z = a^{TE-}(\alpha, \beta)e^{-i\gamma z} + a^{TE+}(\alpha, \beta)e^{i\gamma z} \qquad (6)$$

Coming back to Maxwell equations, we use a $TE_z(H_z = 0)$, $TM_z(E_z = 0)$ decomposition to rewrite the transversal components of the EM field in function of modes.

$$\begin{bmatrix} E_y \\ Z_0 H_x \\ Z_0 H_y \\ E_x \end{bmatrix} = a^{TE\pm}\Phi^{TE\pm} + a^{TM\pm}\Phi^{TM\pm} \qquad (7)$$

Such:

$$\Phi^{\pm} = \frac{1}{k_c^2 - \gamma^2} \begin{bmatrix} \pm\gamma\alpha E_z + \mu_r\alpha Z_0 H_z \\ \varepsilon\beta E_z \pm \gamma\alpha Z_0 H_z \\ -\varepsilon\alpha E_z \pm \gamma\beta Z_0 H_z \\ \pm\gamma\beta E_z - \mu_r Z_0 H_z \end{bmatrix} \qquad (8)$$

## 4.2 M-D for homogenized anisotropic material

In the general case of biaxial anisotropy form of the permittivity tensor, Maxwell's equations can be written as two sets of three equations each:

$$E_{sys} = \begin{cases} \partial_y E_z - \partial_z E_y = -ik\mu_r Z_0 H_x \\ \partial_z E_x - \partial_x E_z = -ik\mu_r Z_0 H_y \\ \partial_x E_y - \partial_y E_x = -ik\mu_r Z_0 H_z \end{cases} \qquad (9)$$

$$H_{sys} = \begin{cases} \partial_y H_z - \partial_z H_y = \frac{ik}{Z_0}(\varepsilon_r^{xx} E_x + \varepsilon_r^{xy} E_y) \\ \partial_z H_x - \partial_x H_z = \frac{ik}{Z_0}(\varepsilon_r^{yy} E_y + \varepsilon_r^{yx} E_x) \\ \partial_x H_y - \partial_y H_x = \frac{ik}{Z_0} \varepsilon_r^{zz} E_z \end{cases} \qquad (10)$$

Likewise, we write $E_z = \text{fct}(H_y, H_x)$ and $H_z = \text{fct}(E_y, E_x)$. Thus we determine two decoupled matrix systems of dimension (2x2) to be solved, verified by the transversal components of the electric field in one hand and the magnetic field in the other hand presented in (11) and (13). In the same framework, Sabbagh and M. Robert in [12, 13, 14, 15] present a similar approach but do not take into account the decoupled variables, they resolve a (4x4) matrix system to describe the EM field in the layers.
$$ \partial_z \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \mathcal{L}_{EG} \begin{bmatrix} H_x \\ H_y \end{bmatrix} \qquad (11) $$

where

$$ \mathcal{L}_{EG} = \begin{bmatrix} -\frac{Z_0 \partial_x \partial_y}{ik \varepsilon_r^{zz}} & \frac{Z_0 \partial_x^2}{ik \varepsilon_r^{zz}} - ik \mu_r \\ \frac{Z_0 \partial_y^2}{ik \varepsilon_r^{zz}} + ik \mu_r & \frac{Z_0 \partial_y \partial_x}{ik \varepsilon_r^{zz}} \end{bmatrix} \qquad (12) $$

$$ \partial_z \begin{bmatrix} H_x \\ H_y \end{bmatrix} = \mathcal{L}_{GE} \begin{bmatrix} E_x \\ E_y \end{bmatrix} \qquad (13) $$

where

$$ \mathcal{L}_{GE} = \begin{bmatrix} \frac{\partial_x \partial_y}{iZ_0 k \mu_r} + ik\varepsilon_r^{vx} & -\frac{\partial_x^2}{iZ_0 k \mu_r} + ik\varepsilon_r^{yy} \\ \frac{\partial_y^2}{iZ_0 k \mu_r} - ik\varepsilon_r^{xx} & -\frac{\partial_x \partial_y}{iZ_0 k \mu_r} - ik\varepsilon_r^{xy} \end{bmatrix} \qquad (14) $$

By applying $\partial_z$ to the first matrix equation verified by the transversal components of the Electric field and the TF, we obtain the eigenvalues system to be resolved in the Fourier domain.

$$ -\gamma^2 \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \mathcal{L}_{EG} \mathcal{L}_{GE} \begin{bmatrix} E_x \\ E_y \end{bmatrix} \qquad (15) $$

This system provides two eigenvalues $\pm\sqrt{\lambda_1}$ et $\pm\sqrt{\lambda_2}$ corresponding to the direction of the propagation of the wave, and two associated eigenvectors $\vec{v}_1$ et $\vec{v}_2$. So $E_x$ and $E_y$ are a linear combination of eigenvectors.

$$ \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \begin{bmatrix} \nu_x^- \\ \nu_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} + \begin{bmatrix} \nu_x^+ \\ \nu_y^+ \end{bmatrix} \begin{bmatrix} a_x^+ \\ a_y^+ \end{bmatrix} \qquad (16) $$

The eigenvalues are used to calculate the wave attenuation factors $\psi_i^{\pm}(z)$ for several depths. From (13) we reconstruct the $\mathbf{H}_t$ components:

$$ \begin{bmatrix} H_x \\ H_y \end{bmatrix} = \frac{-i}{\gamma} \mathcal{L}_{GE}^{-1} \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \begin{bmatrix} w_x^- \\ w_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} + \begin{bmatrix} w_x^+ \\ w_y^+ \end{bmatrix} \begin{bmatrix} a_x^+ \\ a_y^+ \end{bmatrix} \qquad (17) $$

Finally the modal decomposition in the anisotropic layer is defined by:

$$ \begin{bmatrix} \mathbf{E}_t \\ Z_0 \mathbf{H}_t \end{bmatrix} = a_x^{\pm} \mathbf{\Phi}_x^{\pm} + a_y^{\pm} \mathbf{\Phi}_y^{\pm} \qquad (18) $$
$$ \partial_z \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \mathcal{L}_{EG} \begin{bmatrix} H_x \\ H_y \end{bmatrix} \qquad (11) $$

where

$$ \mathcal{L}_{EG} = \begin{bmatrix} -\frac{Z_0 \partial_x \partial_y}{ik \varepsilon_r^{zz}} & \frac{Z_0 \partial_x^2}{ik \varepsilon_r^{zz}} - ik \mu_r \\ -\frac{Z_0 \partial_y^2}{ik \varepsilon_r^{zz}} + ik \mu_r & \frac{Z_0 \partial_y \partial_x}{ik \varepsilon_r^{zz}} \end{bmatrix} \qquad (12) $$

$$ \partial_z \begin{bmatrix} H_x \\ H_y \end{bmatrix} = \mathcal{L}_{GE} \begin{bmatrix} E_x \\ E_y \end{bmatrix} \qquad (13) $$

where

$$ \mathcal{L}_{GE} = \begin{bmatrix} \frac{\partial_x \partial_y}{iZ_0 k \mu_r} + ik\varepsilon_r^{vx} & -\frac{\partial_x^2}{iZ_0 k \mu_r} + ik\varepsilon_r^{yy} \\ \frac{\partial_y^2}{iZ_0 k \mu_r} - ik\varepsilon_r^{xx} & -\frac{\partial_x \partial_y}{iZ_0 k \mu_r} - ik\varepsilon_r^{xy} \end{bmatrix} \qquad (14) $$

By applying $\partial_z$ to the first matrix equation verified by the transversal components of the Electric field and the TF, we obtain the eigenvalues system to be resolved in the Fourier domain.

$$ -\gamma^2 \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \mathcal{L}_{EG} \mathcal{L}_{GE} \begin{bmatrix} E_x \\ E_y \end{bmatrix} \qquad (15) $$

This system provides two eigenvalues $\pm\sqrt{\lambda_1}$ et $\pm\sqrt{\lambda_2}$ corresponding to the direction of the propagation of the wave, and two associated eigenvectors $\vec{v}_1$ et $\vec{v}_2$. So $E_x$ and $E_y$ are a linear combination of eigenvectors.

$$ \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \begin{bmatrix} v_x^- \\ v_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} + \begin{bmatrix} v_x^+ \\ v_y^+ \end{bmatrix} \begin{bmatrix} a_x^+ \\ a_y^+ \end{bmatrix} \qquad (16) $$

The eigenvalues are used to calculate the wave attenuation factors $\psi_i^{\pm}(z)$ for several depths. From (13) we reconstruct the $\mathbf{H}_t$ components:

$$ \begin{bmatrix} H_x \\ H_y \end{bmatrix} = \frac{-i}{\gamma} \mathcal{L}_{GE}^{-1} \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \begin{bmatrix} w_x^- \\ w_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} + \begin{bmatrix} w_x^+ \\ w_y^+ \end{bmatrix} \begin{bmatrix} a_x^+ \\ a_y^+ \end{bmatrix} \qquad (17) $$

Finally the modal decomposition in the anisotropic layer is defined by:

$$ \begin{bmatrix} \mathbf{E}_t \\ Z_0 \mathbf{H}_t \end{bmatrix} = a_x^{\pm} \mathbf{\Phi}_x^{\pm} + a_y^{\pm} \mathbf{\Phi}_y^{\pm} \qquad (18) $$
$$ \Phi_x^{\pm} = \begin{bmatrix} v_x^{\pm} \\ v_y^{\pm} \\ w_x^{\pm} \\ w_y^{\pm} \end{bmatrix} \qquad (19) $$

## 4.3 applying Boundary Conditions (air-plate interface)

In the absence of any surface current ($J_s = 0$), the BC of an EM field at a separating interface are defined by the continuity of the tangential components of the EM field:

$$ \mathbf{E}_1 \times \mathbf{n} = \mathbf{E}_2 \times \mathbf{n} \qquad (20) $$

$$ \mathbf{H}_1 \times \mathbf{n} = \mathbf{H}_2 \times \mathbf{n} \qquad (21) $$

At this stage, we introduce a reference state field $\Phi_{ref}$ which a field du a source term in air:

$$ \Phi_{ref} = \begin{bmatrix} E_{0y} \\ Z_0 H_{0x} \\ Z_0 H_{0y} \\ E_{0x} \end{bmatrix} \qquad (22) $$

And we aim to evaluate a perturbation field so that :

$$ \Phi_i = \Phi_p + \Phi_{ref} \qquad (23) $$

Fig. 2 illustrates the incidence, the reflected and the transmitted field with the unknown coefficients to be calculated by applying BC.

<figure><img src="image_2.png" /><figcaption>Figure 2: BC illustration</figcaption></figure>

Then we write:

### 4.3.1 in the air

$$ \begin{bmatrix} \mathbf{E}_t \\ \mathbf{H}_t \end{bmatrix} = \begin{bmatrix} \Phi_0^{TE+} & \Phi_0^{TM+} \end{bmatrix} \begin{bmatrix} a_0^{TE+} \\ a_0^{TM+} \end{bmatrix} + \Phi_{ref} \qquad (24) $$
Such

$$ \Phi_x^{\pm} = \begin{bmatrix} v_x^{\pm} \\ v_y^{\pm} \\ w_x^{\pm} \\ w_y^{\pm} \end{bmatrix} \qquad (19) $$

## 4.3 applying Boundary Conditions (air-plate interface)

In the absence of any surface current ($J_s = 0$), the BC of an EM field at a separating interface are defined by the continuity of the tangential components of the EM field:

$$ \mathbf{E}_1 \times \mathbf{n} = \mathbf{E}_2 \times \mathbf{n} \qquad (20) $$

$$ \mathbf{H}_1 \times \mathbf{n} = \mathbf{H}_2 \times \mathbf{n} \qquad (21) $$

At this stage, we introduce a reference state field $\Phi_{ref}$ which a field du a source term in air:

$$ \Phi_{ref} = \begin{bmatrix} E_{0y} \\ Z_0 H_{0x} \\ Z_0 H_{0y} \\ E_{0x} \end{bmatrix} \qquad (22) $$

And we aim to evaluate a perturbation field so that :

$$ \Phi_i = \Phi_p + \Phi_{ref} \qquad (23) $$

Fig. 2 illustrates the incidence, the reflected and the transmitted field with the unknown coefficients to be calculated by applying BC.

<figure><img src="image_3.png" /><figcaption>Figure 2: BC illustration</figcaption></figure>

Then we write:

### 4.3.1 in the air

$$ \begin{bmatrix} \mathbf{E}_t \\ \mathbf{H}_t \end{bmatrix} = \begin{bmatrix} \Phi_0^{TE+} & \Phi_0^{TM+} \end{bmatrix} \begin{bmatrix} a_0^{TE+} \\ a_0^{TM+} \end{bmatrix} + \Phi_{ref} \qquad (24) $$
### 4.3.2 in the plate

$$ \begin{bmatrix} \mathbf{E}_t \\ \mathbf{H}_t \end{bmatrix} = \begin{bmatrix} \Phi_x^- & \Phi_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} \qquad (25) $$

Numerically, the reference field is calculated by some external modules already implemented into CIVA. the definition of $\Psi_{ref}$ leads to assume that $a_x^- = 0$ and $a_y^- = 0$. From (20) and (21) applied at the first separating surface, we have:

$$ \begin{bmatrix} \Phi_0^{TE+} & \Phi_0^{TM+} \end{bmatrix} \begin{bmatrix} a_0^{TE+} \\ a_0^{TM+} \end{bmatrix} + \Phi_{ref} = \begin{bmatrix} \Phi_x^- & \Phi_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} \qquad (26) $$

The constants are determined by resolving the system :

$$ \begin{bmatrix} a_0^{TE+} \\ a_0^{TM+} \\ a_x^- \\ a_y^- \end{bmatrix} = - \begin{bmatrix} \Phi_0^{TE+} & \Phi_0^{TM+} & \Phi_x^- & \Phi_y^- \end{bmatrix}^{-1} \Phi_{ref} \qquad (27) $$

finally, we have the total field in the whole region but in the Fourier domain, therefore, we apply the inverse Fourier transform to reconstruct the perturbation field in the spatial domain and from (23) we obtain the total field.

## 4.4 scattering matrix approach

the main idea of the S-matrix algorithm is to represent the different interfaces as four-port black boxes with two inputs and two outputs (respectively ($a_p^-$, $a_{p+1}^+$) and ($a_{p+1}^-$, $a_p^+$) for the ($p+1$)th interface). the inputs-outputs relationship of such a box is then entirely described by a (4x4) matrix, called intermediate S-matrix and denoted by S(p+1,p) for the ($p+1$)th interface). The so-called scattering matrices algorithm is a recursive routine that aims at linking all the intermediate S-matrices in order to compute the global S-matrix, which relates directly ($a_0^-$, $a_N^+$) and ($a_N^-$, $a_0^+$), thus leading to a black box representation of the global mock-up, including all the internal reflections.

<figure><img src="image_3.png" /><figcaption>Figure 3: Stratified media</figcaption></figure>

We refer the reader to [16] for more detailed set up of a general S-matrix algorithm as numerical scheme and implementation.
### 4.3.2 in the plate

$$ \begin{bmatrix} \mathbf{E}_t \\ \mathbf{H}_t \end{bmatrix} = \begin{bmatrix} \Phi_x^- & \Phi_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} \qquad (25) $$

Numerically, the reference field is calculated by some external modules already implemented into CIVA. the definition of $\Psi_{ref}$ leads to assume that $a_x^- = 0$ and $a_y^- = 0$. From (20) and (21) applied at the first separating surface, we have:

$$ \begin{bmatrix} \Phi_0^{TE+} & \Phi_0^{TM+} \end{bmatrix} \begin{bmatrix} a_0^{TE+} \\ a_0^{TM+} \end{bmatrix} + \Phi_{ref} = \begin{bmatrix} \Phi_x^- & \Phi_y^- \end{bmatrix} \begin{bmatrix} a_x^- \\ a_y^- \end{bmatrix} \qquad (26) $$

The constants are determined by resolving the system :

$$ \begin{bmatrix} a_0^{TE+} \\ a_0^{TM+} \\ a_x^- \\ a_y^- \end{bmatrix} = - \begin{bmatrix} \Phi_0^{TE+} & \Phi_0^{TM+} & \Phi_x^- & \Phi_y^- \end{bmatrix}^{-1} \Phi_{ref} \qquad (27) $$

finally, we have the total field in the whole region but in the Fourier domain, therefore, we apply the inverse Fourier transform to reconstruct the perturbation field in the spatial domain and from (23) we obtain the total field.

## 4.4 scattering matrix approach

the main idea of the S-matrix algorithm is to represent the different interfaces as four-port black boxes with two inputs and two outputs (respectively ($a_p^-$, $a_{p+1}^+$) and ($a_{p+1}^-$, $a_p^+$) for the ($p+1$)th interface). the inputs-outputs relationship of such a box is then entirely described by a (4x4) matrix, called intermediate S-matrix and denoted by S(p+1,p) for the ($p+1$)th interface). The so-called scattering matrices algorithm is a recursive routine that aims at linking all the intermediate S-matrices in order to compute the global S-matrix, which relates directly ($a_0^-$, $a_N^+$) and ($a_N^-$, $a_0^+$), thus leading to a black box representation of the global mock-up, including all the internal reflections.

<figure><img src="image_4.png" /><figcaption>Figure 3: Stratified media</figcaption></figure>

We refer the reader to [16] for more detailed set up of a general S-matrix algorithm as numerical scheme and implementation.
# 5 results and discussion

In this section, we illustrate the modal validation through comparison with FEM results. Two test-cases was considered:

* one layered media with $\theta = 0^\circ$ and $45^\circ$

* 4 layered slab with different orientations:

- $[0^\circ, 0^\circ, 90^\circ, 90^\circ]$

- $[0^\circ, 90^\circ, 0^\circ, 90^\circ]$

- $[0^\circ, -45^\circ, 90^\circ, 45^\circ]$

## 5.1 semi-analytic (SA) model features:

<figure><img src="image_4.png" /><figcaption>Figure 4: CIVA configuration for reference field calculation</figcaption></figure>

The model is developed on Matlab software. The test case configurations and parameters we adopted here are:

* Circular coil: its characteristics are described in Tab. 1.

* Anisotropic plate: its characteristics are described in Tab. 2.

* reference field is calculated on a surface delimited by the red outline.

* $M_u=25, M_v=25$.

* Spatial domain: [-23, 32]x[-32,32] mm².

<table><caption>Table 1: Coil characteristics</caption><thead><tr><th>R<sub>int</sub></th><th>R<sub>ext</sub></th><th>Height</th><th>Frequency</th><th>Intensity</th><th>Lift-off</th></tr></thead><tbody><tr><td>4 [mm]</td><td>6 [mm]</td><td>1 [mm]</td><td>10<sup>4</sup> [Hz]</td><td>0.4 [A]</td><td>0.6 [mm]</td></tr></tbody></table>

<table><caption>Table 2: specimen characteristics</caption><thead><tr><th>Thickness</th><th>(σ<sub>xx</sub>, σ<sub>yy</sub>, σ<sub>zz</sub>)</th><th>depths</th><th>θ</th></tr></thead><tbody><tr><td>1[mm]</td><td>(10<sup>4</sup>, 10<sup>3</sup>, 10<sup>3</sup>) [S/m]</td><td>7</td><td>0°</td></tr></tbody></table>
# 5 results and discussion

In this section, we illustrate the modal validation through comparison with FEM results. Two test-cases was considered:

* one layered media with $\theta = 0^\circ$ and $45^\circ$

* 4 layered slab with different orientations:

- $[0^\circ, 0^\circ, 90^\circ, 90^\circ]$

- $[0^\circ, 90^\circ, 0^\circ, 90^\circ]$

- $[0^\circ, -45^\circ, 90^\circ, 45^\circ]$

## 5.1 semi-analytic (SA) model features:

<figure><img src="image_5.png" /><figcaption>Figure 4: CIVA configuration for reference field calculation</figcaption></figure>

The model is developed on Matlab software. The test case configurations and parameters we adopted here are:

* Circular coil: its characteristics are described in Tab. 1.

* Anisotropic plate: its characteristics are described in Tab. 2.

* reference field is calculated on a surface delimited by the red outline.

* $M_u=25, M_v=25$.

* Spatial domain: [-23, 32]x[-32,32] mm².

<table><caption>Table 1: Coil characteristics</caption><thead><tr><td>$R_{int}$</td><td>$R_{ext}$</td><td>Height</td><td>Frequency</td><td>Intensity</td><td>Lift-off</td></tr></thead><tbody><tr><td>4 [mm]</td><td>6 [mm]</td><td>1 [mm]</td><td>10<sup>4</sup> [Hz]</td><td>0.4 [A]</td><td>0.6 [mm]</td></tr></tbody></table>

<table><caption>Table 2: specimen characteristics</caption><thead><tr><td>Thickness</td><td>($\sigma_{xx}, \sigma_{yy}, \sigma_{zz}$)</td><td>depths</td><td>$\theta$</td></tr></thead><tbody><tr><td>1[mm]</td><td>(10<sup>4</sup>, 10<sup>3</sup>, 10<sup>3</sup>) [S/m]</td><td>7</td><td>0°</td></tr></tbody></table>
## 5.2 FEM-(SA) Model validation

In the representation of the model results, we are limited to the $H_z$ components which is the most complicated one to reconstruct since it depends on the two transversal components of the electric field $E_x$ and $E_y$, so if there is a good agreement between data, it will confirm that the other components also agree with FEM results. in figure (5) and (6) we illustrate the 3D view of the real part $\Re(\cdot)$ and imaginary part $\Im(\cdot)$ of $H_z$. It is obvious the effect of the anisotropy with the orientation. This allows us to be able to detect the fiber orientation is each layer by plotting the current density. In figure(7) we show a slice view of $H_z$ at 7 different depths in the slab.

<figure><img src="image_5.png" /><figcaption>Figure 5: 3D view of $H_z$ for one layer case with $\theta = 0^\circ$</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 6: 3D view of $H_z$ for one layer case with $\theta = 45^\circ$</figcaption></figure>
## 5.2 FEM-(SA) Model validation

In the representation of the model results, we are limited to the $H_z$ components which is the most complicated one to reconstruct since it depends on the two transversal components of the electric field $E_x$ and $E_y$, so if there is a good agreement between data, it will confirm that the other components also agree with FEM results. in figure (5) and (6) we illustrate the 3D view of the real part $\Re(\cdot)$ and imaginary part $\Im(\cdot)$ of $H_z$. It is obvious the effect of the anisotropy with the orientation. This allows us to be able to detect the fiber orientation is each layer by plotting the current density. In figure(7) we show a slice view of $H_z$ at 7 different depths in the slab.

<figure><img src="image_6.png" /><figcaption>Figure 5: 3D view of $H_z$ for one layer case with $\theta = 0^\circ$</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6: 3D view of $H_z$ for one layer case with $\theta = 45^\circ$</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Figure 7: 3D view of $H_z$ for one layer case with the sequence [0°, -45°, 90°, 45°]</figcaption></figure>

In Tab. 3 we illustrate the efficiency of the semi-analytic model developed through the calculation time. It's obvious that the model present a exceptional rapidity.

<table><caption>Table 3: Comparison of calculation time</caption><thead><tr><th>FEM</th><th>SA Model</th></tr></thead><tbody><tr><td>30 min to 2h depending on the mesh</td><td>55 [s]</td></tr></tbody></table>

# 6 conclusion

The work presented in this paper constitutes a contribution to a three-dimensional numerical modeling for the simulation of the interaction of an Eddy Current 3-D probe with a anisotropic and stratified planar media, as an homogenized composite material, for simulating eddy current non destructive testing applications. We have developed a semi-analytic model that solves steady-state Maxwell's equations for materials with biaxial conductivity based on modal decomposition of the field. We started with the semi infinite case and then we use the recursive S-matrix algorithm to resolve for the multi-layered case. Finally, the model is numerically validated by using FE Data. The EM field generated in a plate with a 3D air-core cylindrical coil is clearly reproduced at different depths. By using all previous developments already implemented into the CIVA platform, we will be able soon to tackle any kind of EC probe (rectangular, racetrack, D-coil, spiral circular ...).

# acknowledgment

This paper is supported by European Union's Horizon 2020 research and European Training Network NDToN AIR in NDT and SHM of aircraft structures H2020-MSCA-ITN-2016 under grant agreement No 722134.
<figure><img src="image_8.png" /><figcaption>Figure 7: 3D view of $H_z$ for one layer case with the sequence [0°, -45°, 90°, 45°]</figcaption></figure>

In Tab. 3 we illustrate the efficiency of the semi-analytic model developed through the calculation time. It's obvious that the model present a exceptional rapidity.

<table><caption>Table 3: Comparison of calculation time</caption><thead><tr><th>FEM</th><th>SA Model</th></tr></thead><tbody><tr><td>30 min to 2h depending on the mesh</td><td>55 [s]</td></tr></tbody></table>

# 6 conclusion

The work presented in this paper constitutes a contribution to a three-dimensional numerical modeling for the simulation of the interaction of an Eddy Current 3-D probe with a anisotropic and stratified planar media, as an homogenized composite material, for simulating eddy current non destructive testing applications. We have developed a semi-analytic model that solves steady-state Maxwell's equations for materials with biaxial conductivity based on modal decomposition of the field. We started with the semi infinite case and then we use the recursive S-matrix algorithm to resolve for the multi-layered case. Finally, the model is numerically validated by using FE Data. The EM field generated in a plate with a 3D air-core cylindrical coil is clearly reproduced at different depths. By using all previous developments already implemented into the CIVA platform, we will be able soon to tackle any kind of EC probe (rectangular, racetrack, D-coil, spiral circular ...).

# acknowledgment

This paper is supported by European Union's Horizon 2020 research and European Training Network NDToN AIR in NDT and SHM of aircraft structures H2020-MSCA-ITN-2016 under grant agreement No 722134.
References

[1] O. Ouchetto, S. Zouhdi, A. Bossavit, G. Griso, B. Miara and A. Razek, "Homogenization of structured electromagnetic materials and metamaterials". Journal of Materials Processing Technology 181 (2007) 225-229.

[2] G. Wasselynck, Didier Trichet and Javad Fouladgar, "Determination of the Electrical Conductivity Tensor of a CFRP composite Using a 3-D Percolation Model". IEEE Transaction on Magnetics, Vol. 49, No. 5, May 2013.

[3] H. Menana and Mouloud Féliachi, "3-D Eddy Current Computation in Carbon Fiber Reinforced Composites". IEEE TRANSACTIONS ON MAGNETICS, VOL. 45, NO. 3, MARCH 2009.

[4] W. Yin, Philip J. Withers, Umesh Sharma and Anthony J. Peyton, "Noncontact Characterization of Carbon Fiber Reinforced Plastics Using Multifrequency Eddy Current Sensors". IEEE Transactions on instrumentation and measurement, Vol. 58, No. 3, MARCH 2009.

[5] T. M. Roberts, Harold A. Sabbagh and L. David Sabbagh, "Electromagnetic scattering for a class of anisotropic layered media". Journal of mathematical physics, 29, 2675 (1988); doi: 10.1063/1.528060.

[6] L. David Sabbath and Harold A. Sabbagh, "Eddy-Current Modeling and Flaw Reconstruction". Journal of nondestructive evaluation, Vol. 7, No. 1/2, 1988.

[7] H.A. Sabbagh, L. David Sabbagh and Thomas M. Roberts, "An Eddy Current Model For Three-Dimensional Nondestructive Evaluation Of advanced Composites". Naval Surface Weapons Center Research and Technology department, 28 July 1985.

[8] T.M. Roberts, Harold A. Sabbagh and L. David Sabbagh, "Electromagnetic Interactions with an anisotropic Slab". IEEE TRANSACTIONS ON MAGNETICS, VOL. 24, NO. 6, NOVEMBER 1988.

[9] T.M. Roberts, "Explicit Eigenmodes for Anisotropic Media". IEEE Transactions on Magnetics, Vol. 26, No. 6, Nov. 1990.

[10] F. Caire, D. Prémel and G. Granet, "Fast computation of the field diffracted by a multilayered conductor with non-parallel rough interfaces. Application to Eddy Current non Destructive Testing simulation". IEEE Transactions on Magnetics, Vol.51, No.3, March 2015.

[11] F. Caire, D. Prémel and G. Granet, "Semi-Analytic computation of quasi-static field induced by a 3-D eddy current probe scanning a 2D layered conductor with parallel rough interfaces". International Journal of Numerical Modeling: Electronic Networks, Devices and fields int. J. Numer. Model. 2013; 00:1-13.

[12] D. Prémel and G. Granet, "Fast computation of the response of any 3D Eddy Current probe scaling a 3D stratified conductor characterized by a set of arbitrary interfaces of complex shape". International Journal of Applied Electromagnetics and Mechanics 41 (2014) 1-0.
References

[1] O. Ouchetto, S. Zouhdi, A. Bossavit, G. Griso, B. Miara and A. Razek, "Homogenization of structured electromagnetic materials and metamaterials". Journal of Materials Processing Technology 181 (2007) 225-229.

[2] G. Wasselynck, Didier Trichet and Javad Fouladgar, "Determination of the Electrical Conductivity Tensor of a CFRP composite Using a 3-D Percolation Model". IEEE Transaction on Magnetics, Vol. 49, No. 5, May 2013.

[3] H. Menana and Mouloud Féliachi, "3-D Eddy Current Computation in Carbon Fiber Reinforced Composites". IEEE TRANSACTIONS ON MAGNETICS, VOL. 45, NO. 3, MARCH 2009.

[4] W. Yin, Philip J. Withers, Umesh Sharma and Anthony J. Peyton, "Noncontact Characterization of Carbon Fiber Reinforced Plastics Using Multifrequency Eddy Current Sensors". IEEE Transactions on instrumentation and measurement, Vol. 58, No. 3, MARCH 2009.

[5] T. M. Roberts, Harold A. Sabbagh and L. David Sabbagh, "Electromagnetic scattering for a class of anisotropic layered media". Journal of mathematical physics, 29, 2675 (1988); doi: 10.1063/1.528060.

[6] L. David Sabbath and Harold A. Sabbagh, "Eddy-Current Modeling and Flaw Reconstruction". Journal of nondestructive evaluation, Vol. 7, No. 1/2, 1988.

[7] H.A. Sabbagh, L. David Sabbagh and Thomas M. Roberts, "An Eddy Current Model For Three-Dimensional Nondestructive Evaluation Of advanced Composites". Naval Surface Weapons Center Research and Technology department, 28 July 1985.

[8] T.M. Roberts, Harold A. Sabbagh and L. David Sabbagh, "Electromagnetic Interactions with an anisotropic Slab". IEEE TRANSACTIONS ON MAGNETICS, VOL. 24, NO. 6, NOVEMBER 1988.

[9] T.M. Roberts, "Explicit Eigenmodes for Anisotropic Media". IEEE Transactions on Magnetics, Vol. 26, No. 6, Nov. 1990.

[10] F. Caire, D. Prémel and G. Granet, "Fast computation of the field diffracted by a multilayered conductor with non-parallel rough interfaces. Application to Eddy Current non Destructive Testing simulation". IEEE Transactions on Magnetics, Vol.51, No.3, March 2015.

[11] F. Caire, D. Prémel and G. Granet, "Semi-Analytic computation of quasi-static field induced by a 3-D eddy current probe scanning a 2D layered conductor with parallel rough interfaces". International Journal of Numerical Modeling: Electronic Networks, Devices and fields int. J. Numer. Model. 2013; 00:1-13.

[12] D. Prémel and G. Granet, "Fast computation of the response of any 3D Eddy Current probe scaling a 3D stratified conductor characterized by a set of arbitrary interfaces of complex shape". International Journal of Applied Electromagnetics and Mechanics 41 (2014) 1-0.
[13] D. Prémel and G. Granet, "Development of the curvilinear coordinate method for the computation of quasi-static fields induced by an eddy current probe scanning a 3D conductor of complex shape characterized by an arbitrary 2D surface". International J. Numer Model. 2017; e2219.

[14] R. C. Rumpf, "Improved Formulation of Scattering Matrices for Semi-analytical Methods that is consistent with convention". Progress In Electromagnetics Research B, Vol. 35, 241-261, 2011.

[15] E.l. Tan, "Note on formulation of the enhanced scattering-(transmittance-) matrix approach". J. Opt. Soc. Am. A, Vol. 19, No. 6, JUNE 2002.

[16] L. Li, "Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings". IJ. Opt. Soc. Am. A, Vol. 13, No. 5, JUNE 1996.

[17] L. Li, "Note on the S-matrix propagation algorithm". IJ. Opt. Soc. Am. A, VOL. 20, NO. 5, April 2003.

[18] Shaoni Jiao, Jian Li, Fei Du, Lei Sun and Zhiwei Zeng, "Characteristics of Eddy Current Distribution in carbon fiber Reinforced Polymer". Journal of Sensors, Volume 2016, Article ID 4292134, 8 pages.

[19] Wuliang Yin, Philip J. Withers, Umesh Sharma and Anthony J. Peyton, "Noncontact characterization of carbon fiber reinforced plastics using multifrequency eddy current sensors". IEEE Transactions on instrumentation and measurement, Vol. 58, No. 3, MARCH 2009.

[20] J. B. Park, T. K. Hwang, H. G. Kim, and Y. D. Doh, "Experimental and numerical study of the electrical anisotropy in unidirectional carbonfiber-reinforced polymer composites," Smart Mater. Struct., vol. 16, pp. 57-66, 200.

[21] S. B. Pratap and W. F. Weldon, "Eddy current in anisotropic composites applied to pulsed machinery," IEEE Trans. Magn., vol. 32, no. 2, pp. 437-444, Mar. 1996.

[22] G. Granet, "Analysis of diffraction by surface-relief crossed gratings with use of the Chandezon method: application to multilayer crossed gratings". Journal of the Optical Society America A 1998; 15:1121-1131.

[23] N P K. Cotter et al. "Scatteringmatrix approach to multilayer diffraction". Journal of the Optical Society America 1995; 12: 1097-1103.

[24] L. Li, "Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings". Journal of the Optical Society America A 1996; 13: 1024-1035.
[13] D. Prémel and G. Granet, "Development of the curvilinear coordinate method for the computation of quasi-static fields induced by an eddy current probe scanning a 3D conductor of complex shape characterized by an arbitrary 2D surface". International J. Numer Model. 2017; e2219.

[14] R. C. Rumpf, "Improved Formulation of Scattering Matrices for Semi-analytical Methods that is consistent with convention". Progress In Electromagnetics Research B, Vol. 35, 241-261, 2011.

[15] E.l. Tan, "Note on formulation of the enhanced scattering-(transmittance-) matrix approach". J. Opt. Soc. Am. A, Vol. 19, No. 6, JUNE 2002.

[16] L. Li, "Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings". IJ. Opt. Soc. Am. A, Vol. 13, No. 5, JUNE 1996.

[17] L. Li, "Note on the S-matrix propagation algorithm". IJ. Opt. Soc. Am. A, VOL. 20, NO. 5, April 2003.

[18] Shaoni Jiao, Jian Li, Fei Du, Lei Sun and Zhiwei Zeng, "Characteristics of Eddy Current Distribution in carbon fiber Reinforced Polymer". Journal of Sensors, Volume 2016, Article ID 4292134, 8 pages.

[19] Wuliang Yin, Philip J. Withers, Umesh Sharma and Anthony J. Peyton, "Noncontact characterization of carbon fiber reinforced plastics using multifrequency eddy current sensors". IEEE Transactions on instrumentation and measurement, Vol. 58, No. 3, MARCH 2009.

[20] J. B. Park, T. K. Hwang, H. G. Kim, and Y. D. Doh, "Experimental and numerical study of the electrical anisotropy in unidirectional carbonfiber-reinforced polymer composites," Smart Mater. Struct., vol. 16, pp. 57-66, 200.

[21] S. B. Pratap and W. F. Weldon, "Eddy current in anisotropic composites applied to pulsed machinery," IEEE Trans. Magn., vol. 32, no. 2, pp. 437-444, Mar. 1996.

[22] G. Granet, "Analysis of diffraction by surface-relief crossed gratings with use of the Chandezon method: application to multilayer crossed gratings". Journal of the Optical Society America A 1998; 15:1121-1131.

[23] N P K. Cotter et al. "Scatteringmatrix approach to multilayer diffraction". Journal of the Optical Society America 1995; 12: 1097-1103.

[24] L. Li, "Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings". Journal of the Optical Society America A 1996; 13: 1024-1035.