Abstract
An improved method which can analyze the eddy current density in conductor materials using finite volume method is proposed on the basis of Maxwell equations and T-formulation. The algorithm is applied to solve 3D electromagnetic nondestructive evaluation (E’NDE) benchmark problems. The computing code is applied to study an Inconel 600 work piece with holes or cracks. The impedance change due to the presence of the crack is evaluated and compared with the experimental data of benchmark problems No. 1 and No. 2. The results show a good agreement between both calculated and measured data.
1. Introduction
The finite volume method (FVM) is a discretization method which is well suited for the various kinds of the conservation laws of numerical simulations. Originally, the FVM was extensively used in solving the dynamic fluid problems [1, 2]. Since the beginning of the eighties, it has been introduced for the solution of Maxwell’s equations, precisely in wave propagation domain [3, 4]. In 2001, the FVM was presented as a new method to investigate the eddy current distribution problems from the teamwork problem 21 [5]. After that, FVM was used to build a 2D model to solve magnetic problems with a thin conductor compared with its skin depth [6]. Recently, the FVM has been introduced for the solution of 3D practical nonlinear magnetostatic problems [7]. Similarly, a model of 3D eddy current nondestructive testing (NDT) problems by FVM method is presented for benchmark problem JSAEM No. 6 and furthermore, the movement of the probe coil is taken into account thanks to the nonconforming mesh generation [4]. However, the problem domain should be remeshed when the position of the coil is changed.
The FVM is quite different from the finite difference method (FDM) or finite element method (FEM). Roughly speaking, the FDM is the oldest numeric method to approximate the differential equations through the use of Taylor’s expansions. The FDM becomes difficult to use when the coefficients involved in the equation are discontinuous (in the case of heterogeneous media) which is not a problem for FVM if discontinuities of coefficients occur corresponding with the boundary of “control volume” [2]. There is also a potential bottleneck of the FDM when handling the multidimension geometry [8]. These shortcomings drive to replace the basic of differential equation by integral formulation and subsequently the development of FVM and the FEM. From the industrial point of view, the FVM is known as a robust and cheap method for the discretization of conservation laws [2].(i)“Robust” means that the FVM behaves well for different formulations and focuses on the physical meaning; the obtained solution does not appear without physical meaning.(ii)“Cheap” is obtained thanks to the simple and reliable computational coding for the complex problem.
The FVM maintains the FEM’s advantages, such as the facility to deal with the complex geometry [9] and boundary conditions, but the mathematical expressions is easier than FEM. For simple geometry, the FVM uses less CPU time and memories than the FEM to solve the IEEJ problem [10] and to solve composite materials’ problems [11].
This paper proposes an investigation work on -iterative technique based on -Omega formulation associated with the FVM method. With the algorithm we have developed, only the electrical vector potential is calculated. Moreover, thanks to the application of Biot-Savart’s law [11–13], the movements of the coil are taken into account without re-meshing the problem domain and it permits to calculate directly the magnetic field generated by the eddy current in the conductor. This paper discusses also how to deal with the crack or hole using the FVM.
To test the proposed technique, the simulated results are compared with the experienced ones which are obtained from the E’NDE benchmark problems. The aim is to evaluate the impedance change due to the presence of a crack in the center of the isotropic conductor, taking the probe movements into account. The test system is showed in Figure 1. The pancake coil is located in the center over a plate conductor. The conductor has a thickness-through hole in the center. The coil moves along the direction of the crack length (-axis).
2. Formulation and Algorithm
A-V formulation and -Omega formulation are used a lot to solve electromagnetic problems. Compared with A-V formulation, -Omega formulation has the advantage of diminishing the number of unknown elements being solved in the system [14]. In -Omega formulation, the electric vector potential exists only in the conductor and the magnetic scalar potential Omega does not exist only in the conductor but also in the air. However, solely the conductor is interesting in our problem. Therefore, we want to use a simple method similar with -Omega formulation, considering only the electric potential . Some authors have already proposed -iterative method where the Omega is not included [13]. To begin with the Maxwell equations and the supplementary equations we can obtain the final formulation as where is the electric vector potential, is the magnetic field, and are electrical conductivity and magnetic permeability; respectively, is the electric field, and is the magnetic flux density.
To explain the formulation to calculate the eddy current in the conductor, the scheme is shown in Figure 2.
At the beginning, the magnetic field created by the coil is applied directly to the conductor by using Biot-Savart’s law [15]: where is the applied source magnetic field, is the source current density, and is the distance between the source and the conductor. Obviously, the eddy current excited by this magnetic field in the conductor is obtained by using (5). According to the equation of electric current conservation, the eddy current density can be derived as Meanwhile, the eddy currents in the conductor generate another magnetic field which is opposite to the source magnetic field in the conductor. This field can be calculated by reusing the Biot-Savart law: where represents the distance between the volumes in the conductor. Afterwards, it is added to the source magnetic field to replace the former source field. Thus, the cycle continues to form an iterative process, while the vector and its normal component on the surface are concerned to calculate by using computer code “ELMES” [13, 16].
To carry out such a modeling, the association of FVM to this iterative process is used. Additionally, it is necessary to control the change of the variables in the iterative process of solution. A relaxation is introduced in the following equation which is based on relaxation factor . The iteration continues until the convergence is reached. The convergence is assigned to . It can accelerate the convergence and reduce the number of iterations, consequently the computation time. For each iteration, the correction added to new estimation of the solution vector is given by where and present the th and th iteration of electric vector potential, respectively. The algorithm of this iterative process is shown in Figure 3.
3. Finite Volume Method
As introduced in [1, 4], the principle of the FVM method is to divide the whole domain into discrete control volume. The surfaces of control volume are positioned midway between adjacent points. Thus, each point is surrounded by a control volume or cell. A general point is indentified by and its neighbors in three-dimensional geometry, the nodes to the west and east along axis are identified by and , the nodes to the west and east along axis are and , the nodes to the top and bottom along axis are, respectively, and (see Figure 4). The integration of differential forms of this formulation for each volume belonging to the conductor can be written as where represents the volume of element . By developing (10), we find that the inverses of the surface conductivities are used on the left side. The conductivity is constant everywhere in the conductor for the materials without crack. But if there is a crack, the conductivity is considered as zero which is a different material connected to the conductor. It brings a problem when calculating the inverse of it. To deal with this issue, the conductivity is set at a value which is several thousand times smaller than the one of the conductor [17]. Moreover, the control volumes are divided in order to coincide with the boundary of two materials; thus the method with which the conductivity of the interface is treated. As shown in Figure 4, the physical interface represents a surface in the y-z plane. Inside, the elements and are two connected control volumes in two different materials with a common interface .
In fact, they belong to different volumes but should have the same conductivity. A simple approximation for the conditions of passage between two media is to make the geometric mean of the two physical volumes at the interface [7, 18]. Thus, for the interface , the inverse of conductivity can be written as where and are the inverse of the conductivities of the volumes and , respectively. Equally, the inverses of the conductivities on the other interfaces keep their values.
4. Validation and Application
There were a lot of discussions on the E’NDE benchmark problems. Bowler et al. [19] describe theoretical predictions of eddy current probe responses for NDT and compared with the benchmark problem. A boundary integral vector potential formulation has been developed to evaluate the eddy current interactions [20]; the formulas of the electric dyadic Green’s function for a limited thickness plate are derived [21]. Some commercial software is used to simulate the eddy current distribution in the plate of this benchmark problem [22, 23]. The proposed algorithm associated with FVM method is applied to the E’NDE benchmark problems No. 1 and No. 2 (as shown in Figure 1). These problems deal with a pancake coil, placed above a conductive square plate in the centre with a lift-off of 0.5 mm. The thickness of the plate is 1.25 mm. The coil has ax-symmetric shape and is made of 140 turns, supplied with an alternative current of 0.007 A and 150 kHz or 300 kHz. Its inner and outer diameters are mm and mm. The height is 0.8 mm. The parameters of the coil are shown in Figure 5. The plate has an electric conductivity and a relative magnetic permeability .
For problem No. 1, the dimension of the plate is 40 mm × 40 mm. The coil is fixed over the centre of the plate. In this problem, the solution requires the determination of the impedance change of the plate. When we use Biot-Savart’s law to calculate the magnetic field in the conductor created by eddy current, it risks running “out of memory” in Matlab as the matrix is full. To reduce the full matrix to a “sparse” one, we suppose that the eddy currents far away from the coil create a small magnetic field that can be neglected. It means that the maximum distance for the computation is limited to several times the outer radius of the coil. The impedance change of the plate and the difference between the measured and simulated resulted are shown in Figure 6 at the frequency of 150 kHz.
(a)
(b)
From these figures, we find that the impedance change approaches the measurement from and it begins to remain stable. An approximation of this fact is illustrated in Figure 7.
The impedance changes which are going be evaluated by the proposed method are compared with those obtained from measurements at the frequencies of 150 kHz and 300 kHz (see Table 1). , and represent the impedance change, real part, and imaginary part of impedance change, respectively.
From the table, there is, respectively, 4.72% and 6.96% differences for the two frequencies between the simulated and measured results. The fact that we use the same meshes in both cases explains why in case of 300 kHz the error is higher than that in case of 150 kHz. To solve the system, a preconditioning method developed in [24] is used, and the stabilized biconjugate gradient solver (BICGSTAB) is applied. The algorithm used 237 iterations, and BICGSTAB used 27 iterations to inverse the matrix on a 2.53 GHz Intel Core 2 PC, where the plate has been divided into 14400 volumes.
For problem No. 2, the size of the plate is ; the rectangular thickness-through crack has the dimension of 10 mm length and 0.2 mm width. The coil moves parallel to the axis, along the crack length direction from the centre to , with a step of . The solution of this problem requires calculating the impedance change due to the presence of the hole. It can be regarded as the difference between the impedances of the plate without hole and with hole. The impedance change in the function of the movements of the coil is represented in Figure 8. These figures show a comparison between the calculated values obtained with the proposed method and the measured results.
(a)
(b)
We can see that the impedance change becomes smaller when the coil is far from the centre of the hole. When the coil is at about 10 mm from the centre of the hole, the impedance remains the same as if there was not a hole in the plate. For the frequency of 150 kHz, the skin depth is a little greater than the thickness of the plate while it is smaller than the thickness for the frequency of 300 kHz. Due to the issue of skin depth, BICGSTAB used 500 iterations to inverse the matrix and the algorithm used 1500 iterations for every position on the same PC, where the plate has been divided into 34000 volumes at the frequency of 150 kHz. However, 1200 iterations to inverse the matrix, 2100 iterations for every position, and 45000 volumes are needed for 300 kHz.
5. Conclusion
The E’NDE benchmark problems No. 1 and No. 2 have been solved by the -iteration technique associated with the FVM method. Thanks to Biot-Savart’s law, the mesh of the plate remains the same when the coil moves and the magnetic fields created by the eddy current are evaluated. A limitation of reducing the full matrix is applied. The approximation of interface conductivity calculation for the FVM has shown its efficiency to deal with cracks or holes in the conductor. The simulation results are in agreement with the benchmark problems.