Analytical model of eddy current bobbin coil probe responses at support plates in steam generator and heat exchanger tubes

The steam generators (SGs) of CANadian Deuterium Uranium nuclear reactors require periodic inspections to ensure their safe operation. Eddy current (EC) testing is the primary method by which the SG tubes are inspected. Many conditions in the SG tubes affect the EC response, such as fretting, pitting, cracking, as well as tube expansion, and the presence of tubesheet and support structures. When two or more of these parameters overlap, the EC signals from degradation induced wall losses are frequently difficult to interpret because these wall loss signals are often distorted due to interference with background signal variations. In this paper, a novel analytical model that describes the impedance of a bobbin coil with a plate encircling an arbitrary number of cylindrical conductors is developed. The model is validated against finite element method modelling and experiment. This model can be used to simulate the main probe response, while inspecting SG tubes near the tubesheet and support structures and thereby, provides the potential to separate out smaller flaw responses.


Introduction
Conventional eddy current (EC) testing (ECT) is a common method of non-destructive evaluation (NDE) and is the main method by which steam generator (SG) tubes are inspected [1].Analytical models are useful to quickly predict the impedance of an EC probe for a specific geometry.However, conventional ECT has limitations to address.Inspection at the tubesheet and Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence.Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
support structures is challenging because of all the parameters that affect ECT measurement at that location.These include magnetite sludge, which builds up on the bottom of the SG, the expansion of the tube into the tubesheet, the tubesheet's edge, including any overlay that may be present, in addition to potential flaws, such as pitting, cracking, fretting, stress corrosion cracking (SCC) and intergranular attack [1].When two or more of these parameters overlap, the EC signal becomes difficult to interpret.The challenge lies in detecting flaw signals within larger geometric and material artefacts, such as SCC at SG tube expansion, pitting at tubesheet edges or fretting at support structures.
Analytical models are useful to quickly predict the impedance of an EC probe for a specific geometry.Dodd and Deeds [1] derived a model for the impedance of a coil carrying a constant amplitude oscillating current in the presence of nearby layered materials underneath or inside a coil.This model has been extended to the case of an arbitrary number of cylindrical conductors inside and outside a bobbin probe [2].Truncated region eigenfunction expansion (TREE) has been used extensively to derive models for bobbin coil impedance in more complicated geometries, such as near the end of a nonferromagnetic tube [3], a borehole [4], a conducting quarter space [5] or plate [6].It has also been used to model an EC probe with a ferrite core above stratified conducting regions [7].Desjardins et al [4,8] extended Dodd and Deed's method to account for the feedback effects that arise from controlled voltage in a transmit-receive (T-R) probe configuration.
The finite element method (FEM) is an approach to finding approximate solutions to partial differential equations for more complicated geometries.Mokros [6] used FEM to simulate the response of a T-R probe near a trefoil broach support for SG tubes, Klein et al [7] used FEM to simulate T-R probe response in pressure tube to calandria tube gap measurement.Additionally, FEM has been used to aid validation of several analytical models [3,7,9,10].
Analytical and computer modelling of EC signals has been useful in EC NDE.Computer modelling of EC signals has been used to improve analysis and flaw characterisation [11].A model that describes probe wobble in SG tube inspection has been derived using the second-order vector potential [12].Probe wobble has also been removed using a wavelet transform technique [13].A time-harmonic analytical solution has been extended for pulsed ECT (PECT) using Fourier superposition for inspection of cracks [10].A method to reduce lift-off in PECT using the slope of the linear curve of the peak value of the difference signal and lift-off for determination of defect depth or width has been proposed [14].
Analytical models show potential for use in model based inverse algorithms, either derivative free methods such as the Nelder-Mead simplex method [15,16], or methods that utilise derivatives, such as Newton's method [17], the Gauss-Newton method [18].Inverse models have also been developed using error minimisation algorithms [15,16].The models may also be used to generate training data for artificial neural networks [15,19].
This paper analytically derives the impedance of a coil within a bobbin probe.The bobbin probe includes a second coil that it is mutually coupled to, reflecting the actual bobbin probe configuration used in SG inspection, which records both absolute and differential EC data.The mutual coupling between the coils is important under conditions where the leading coil passes into the region of a ferromagnetic support plate, interacting with the structure, before the second coil does.The bobbin probe is considered inside an arbitrary number of cylindrical regions, with the outermost region containing a conducting ferromagnetic slab.Examples with explicit solutions for the bobbin probe in geometries of interest in SG tube inspection are presented including: the slab alone, a tube held snugly by a supporting slab, and a tube with an air gap between it and its support.Analysis of this problem applies the method used by Dodd et al [2], which relates the coefficients of the vector potential in one region to those of the next region by a matrix equation.Applying this method using TREE [20], with a different external region based on similar models by Theodoulidis and Bowler [3,9,21], leads to useful solutions for bobbin coil impedance in the vicinity of a support plate structure and tube with an air gap.
The third example presented in this paper, an SG tube within a support plate, has been analysed by Theodoulidis and Reboud [22], who used a Leontovich impedance boundary condition [23] on the support plate surface.However, key differences in the methodology and validation are presented here.In this model, the magnetic vector potential inside the support structure is considered, allowing the structure to be as thin as desired without loss of accuracy.The model incorporates mutual inductance between the two coils, allowing for more accurate modelling of the actual SG bobbin probe configuration, which consists of two inductively coupled coaxial coils where one coil is monitored in the absolute configuration.A novel eigenvalue search and calculation algorithm for EC signal modelling is presented.The eigenvalue calculation method is based on Luck and Stevens' [24], but with customisations, such as function scaling and choice of search region, that make calculations efficient and reliable.In addition to a comparison with FEM results, the model is validated against experimental data, which was absent in [22].Experimental validation is performed against five frequencies and a comparison of the analytical and FEM predictions with experimental measurements is performed.The format of the analytical model solution allows for its future incorporation within inverse models for the extraction of multiple parameters at support plate structures from bobbin coil data.

Theory
Derivation of the analytical model uses magnetic vector potential A under time-harmonic A (r, ϕ , z, t) = A 0 (r, ϕ , z) e jωt and axisymmetric conditions (A = A ϕ φ ) in the Coulomb gauge (∇ • A = 0).Under these conditions, the vector potential and the differential equation describing can be written where ω is the driving angular frequency, µ and σ are, respectively, the permeability and conductivity of the medium.A ϕ in equation ( 1) can be separated into a product of radial and axial components, i.e.A ϕ (r, z) = R (r) Z (z) .Applying the vector Laplacian allows each component to be separated.One part can be set to a separation variable k 2 and the other can be set to −k 2 , so that both parts will sum to zero.Which part is set to what sign of the number and which part gets the jωµσ term is a matter of convenience.Here, the radial part will be set to the positive variable and the axial part will be set to the negative.Each part has two cases, one in which it inherits the jωµσ term, and one where it does not: Radial part, air : Radial part, material : Axial part, air : Axial part, material : Equations ( 2) and ( 3) can be rewritten respectively as These are the modified Bessel equations, and have the solutions [25] where I 1 (x) , K 1 (x) are the modified Bessel functions of the first and second kind, respectively, with order 1.Equations ( 4) and ( 5) have solutions respectively.In ( 8)- (11), are functions that are constant with respect to position and depend only on the separation variable.R (r) and Z (z) in the forms of ( 8)- (11), represent the general solution to equation (1).The particular solution is a summation over all values that k is allowed to have such that the boundary conditions are satisfied.

Boundary conditions
Boundary conditions on the vector potential state that it is continuous at the boundary between two regions [25] and its derivative is derived from boundary conditions on the magnetic field [25,26].All boundaries in the considered geometries will have either constant radius or constant height.The general case of two regions, 1 and 2, that share a boundary with a normal unit vector n is shown in figure 1.Consider the case in which the boundary has constant height, z = z b .The unit vector normal to this surface is z.The vector potential is continuous across this boundary, for all r in the boundary.The condition on the derivative simplifies to In the case when these two regions share a radial boundary at r = r b .The normal unit vector is r.The vector potential is still continuous, for all z in the boundary.The condition on the radial derivative is written If there happen to be any surface currents in the circumferential direction, they are added to the left-hand side of ( 13) or (15).Modified Bessel functions have a convenient identity involving derivatives of the form in equation ( 15) [27,28]: The asymptotic conditions on the vector potential are straightforward.In a region that contains the z-axis (r = 0), the vector potential remains bounded.Since K 1 (x) diverges for x → 0, its coefficient in that region is necessarily 0. Likewise, the vector potential remains bounded as r → ∞, for external regions (regions where r has a lower bound, but no upper bound), P (k) = 0, because I 1 (x) increases without bound for x → ∞.

TREE of coil coefficients
Before analysing the full system, the vector potential due to the coil alone shall be presented.This is because the analysis of the full system depends on the expansion coefficients of only the coil being known.The derivation follows that which is presented in [5] most closely, but for a different geometry.
Consider a delta coil with radius r 0 and axial position z 0 that carries current i 0 .Here, the region is truncated at z = ±h.The domain can be broken into two regions: region + where r > r 0 and region −, where r < r 0 .
As per equation ( 10), the z-dependence of the vector potential can be written as a combination of trigonometric functions with the assertion that the vector potential is 0 on the truncation boundaries, which allows it to be written as a discrete sum of sine terms alone.Applying this condition, (18) where k i = iπ 2h , C +,i are the coefficients in region r > r 0 , C −,i are the coefficients in region r < r 0 .For compactness and clarity, new notation is introduced.Notice that equations ( 17) and ( 18) may be written as a product of a row vector, a diagonal matrix and a column vector.For any column vectors u, v, w, whose ith components are u i , v i and w i , respectively, where λ is a scalar, superscript T indicates matrix transpose and the vector in angled brackets indicates a diagonal matrix such that each element of the vector is the corresponding diagonal element of the matrix, i.e. (⟨v⟩) ii = (v) i .This is equivalent to the diag() function in MATLAB [28].Since the order of vectors on the right-hand side of ( 19) can be changed, the relations may be obtained.Additionally, it is defined that vectors inside angled brackets multiply element-wise, i.e.
With this notation in mind, equations ( 17) and ( 18) may be rewritten (23) where (k) i = k i , likewise for C + and C − .It is implied that functions with vector input are evaluated elementwise.For example, K 1 (kr) = (K 1 (k 1 r) , K 1 (k 2 r) , . ..) T .Applying boundary conditions ( 14) and ( 15) with surface current i 0 δ (z − z 0 ), The orthogonality of trigonometric functions can be exploited with the relation [29], where δ ij is the Kronecker delta symbol.By multiplying both sides of ( 24) and ( 25) by column vector sin (k (z + h)) and integrating z from −h to h, Isolating 26 for C − and substituting into 27, a formula for C + is derived: Inverting a diagonal matrix is straightforward, each diagonal element is replaced with its reciprocal [30], i.e. ⟨v⟩ −1 = 1 v .Using the Wronskian of K 1 (x) , I 1 (x) [27], (28) simplifies to A similar procedure may be performed to find C − .Now consider a coil with rectangular cross section, extending from r = r 1 to r 2 and z = z 1 to z 2 as shown in figure 2. Because it is sufficient to know the vector potential, due to the coil in the region where r > r 2 , to determine the vector potential in the material and beyond, attention shall be focused there.As before, the total vector potential in region r > r 2 can be written as an integral of the vector potential due to a delta coil over current density Therefore, the total vector potential can be written as where the ith coil coefficient C c is given by where i 0 is the current through the coil and

TREE of the vector potential
Figure 2 shows a cross section of the target geometry, which accommodates an SG tube separated by a support structure by an air gap.It shows a geometry of a coil within M cylindrical layers with a conducting slab (grey rectangle) present.
As shown in figure 2, the cylindrical regions are labelled 1 to M and the regions above, in, and below the slab are regions 1 ′ , 2 ′ and 3 ′ , respectively.Each region m has outer radius ρ m and inner radius ρ m−1 , conductivity σ m and relative permeability µ r,m , the radius of the hole in the slab is the outer radius of the Mth region, ρ M .The slab has conductivity σ s , relative permeability µ r,s and thickness c.
From equations ( 8)-( 11), the vector potential in each region can be obtained as where P m and Q m are the coefficients in region m.Equation (34) can also apply for m = 0, but P 0 = C ed , the coefficients for the vector potential due to ECs, Q 0 = C c for the coil's coefficients, and ρ −1 = r 2 .q and p = q 2 − jωµ 0 µ r,s σ s are the eigenvalues in the primed regions.
Applying boundary conditions ( 12) and ( 13) at the horizontal boundaries of the slab let Evaluating the vector potential at z = −h, and asserting the truncation condition leads to the equation for the eigenvalues: The set of solutions q, where each q = (q) i satisfies sin (q (c − h)) ( cos (pc) cos (qh) − p µr,sq sin (pc) sin (qh) ( cos (pc) sin (qh) + qµr,s p sin (pc) cos (qh) Boundary conditions ( 14) and ( 15) at r = ρ M say that the vector potential at this boundary is a piecewise function of z.It is still possible to exploit the orthogonality of the sine function.By multiplying both sides by the column vector sin (k (z + h)) from the left, integrating z from −h to h, and substituting the coefficient vectors with their expressions in terms of C 1 ′ the two piecewise relations collapse to a pair of z-independent equations where + sin (k i (c − h)) q j cos (p j c) cos (q j h) − p j µ r,s q j sin (p j c) sin (q j h) + sin (k i h) q j cos (q j h) The following is presented for clarity when multiplying block matrices.Suppose two large matrices of equal size can be broken into a 2 × 2 array of submatrices.The product of the large matrices is the block-wise product of the submatrices [31,32]: Following this, equations ( 44)-(48) may be packaged together as block matrix equations: To simplify, let which represents the contribution of region l at radial boundary m, and which is a diagonal matrix of ones and µ r,m 's whose inverse accounts for the contribution to the derivative condition by the relative permeability of material m.The coefficients in regions 1 to M − 1 can be written in terms of the next region's coefficients: with and The goal is to find an expression for C ed in terms of C c .To do this both terms are expressed as a linear transformation of C 1 ′ , then one is substituted into the other where ↔ U is the product of all the matrices that transform the vector.Because every matrix that contributes to where Solving C 1 ′ in terms of C c , and substituting, With the vector potential due to ECs now known, the change in coil impedance may be calculated.

Coil impedance
The change in impedance with respect to its value in air is expressed in terms of the voltage induced by ECs.This voltage, V ed , is found by integrating their vector potential A ed over the coil's cross section.
This allows change in coil impedance to be written as where Change in mutual inductance can be found by integrating the vector potential due to the ECs caused by one coil over the cross section of the other.If the 'driver' placed between z = z 1 , z 2 and the 'receiver' placed between z = z 3 , z 4 , the resulting expression is where C 1 is given in (65), C 2 is also given by ( 65), but with the axial positions, radii and turn count replaced by the equivalent parameter from the receiver.Because the only variable in equations ( 63) and ( 64) that changes with the number of cylindrical regions (other than the index of ρ M ) is ↔ U, the following examples shall only present its explicit forms, with the implication that the full solution to coil impedance is (62) with ↔ A and ↔ B given by ( 63) and (64) respectively.

Example 1: coil inside a support structure only.
If there are no cylindrical regions (M = 0), then only equation ( 55) applies.This represents the probe's response to a support structure alone.Now, the matrix ↔ R lm has a convenient analytic inverse: Knowing the Wronskian of modified Bessel functions, this is not difficult to verify by direct evaluation.From (55), (69) ↔ µ 1 contains the relative permeability of the tube µ r,t , and ↔ U will have ρ 0 and ρ 1 as scalar coefficients (i.e. ↔ U = ρ 0 ρ 1 (matrices)), which means ρ 0 and ρ 1 will also be scalar coefficients of ↔ U's submatrices.However, because of the form of (62) (a matrix-inverse matrix product), any scalar coefficients that equally affect ↔ A and ↔ B cancel.Therefore, when impedance change is desired, ρ m in (68) may be ignored without loss of accuracy.

Example 3: tube-air gap-support.
Now consider a tube with inner radius ρ 0 , outer radius ρ 1 , centred in a hole of radius ρ 2 in a slab of conducting material.There are two cylindrical regions (M = 2), the first is the tube, the second is an air gap between the tube and its support.↔ U may be derived in this case by starting with (55), applying (54) once, then capping it with (56).Doing this, while noting that the second material is air, gives A three-dimensional (3D) view of this geometry is shown in figure 3.

Eigenvalues
The novel method for eddy current signal modelling by which the eigenvalue equation is solved will be presented here.The eigenvalue equation ( 43) is a transcendental equation with complex valued solutions.Its expression becomes much simpler in the case of a thick support plate by setting the truncation where H = h δs , δ s is the skin depth of the support structure, Q = qδ s and P = Q 2 − 2j.
The algorithm for solving equation (72) comes from Luck and Stevens [24].There, a complex valued function f (z) has a singularity, or pole, at a point z 0 within a closed contour C on the complex plane.The value for z 0 can be given by Now, let C be a circle of radius R on the complex plane centred at a point Q n , i.e. z = Q n + Re jθ and dz = jRe jθ dθ.This allows Equation (73) to be written as A very efficient method to evaluate the integrals in equation ( 74) is the fast Fourier transform (FFT); the integrals in equation (74) correspond to the first and second coefficients in a Fourier series (indexed from zero).
For the implementation of the algorithm in this paper, the function was set as the reciprocal of left-hand side of equation (72).By graphical inspection, it was found that solutions for (72) with µ r,s ̸ = 1 form two lines on the complex plane, one on the bottom, just above the real axis, and the other on top, just below the line y = 1 x .The circles were chosen for maximum coverage of the plane and exactly one solution within each.For the bottom solutions, The top solutions were then transformed with the function √ x 2 + 2i, which maps values x from the real axis to the line y = 1 x .The FFT length was initialised to 32.The eigenvalues were then refined using the standard Newton-Raphson method [30].
In figure 4 below, the red and blue lines illustrate where the real and imaginary parts of the left-hand side of equation (72) evaluate to zero.A solution is located wherever the lines cross.The bottom search circles are shown in green.In certain cases, the search circle for n = 1 includes the origin, the trivial solution to equation (72).The smaller blue circle in figure 4 represents the search circle for this case.Note that dimensions have been normalised with respect to the skin depth of the support structure.

Experiment
To help validate these models, experiments were designed and carried out to measure the impedance of a coil.
The design of the coil is based on the Zetec ® A-480-ULC, the eddy current probe used to inspect CANadian Deuterium Uranium (CANDU ® ) SGs.It consists of two identical, coaxial coils with a gap of one coil length between them.The probe diameter is 12.0 mm (0.48") as available from Zetec ® [32].The probe was wound in a 3D printed former and connected to a BNC connector.The probe's coils were connected in a pseudo-T-R configuration, where the monitored leading coil was connected to the BNC connector and the other was shorted to itself, to simulate the actual bobbin coil configuration, where Z p is the total impedance of the probe, Z 1 is the total impedance of the drive coil, Z 2 is the total impedance of the second coil, M is the total mutual inductance between the two.Total impedance here includes the DC resistance of the coil, its self-inductance and the lossy inductance due to surrounding structures.Similarly, total mutual inductance is its value in air plus the lossy contribution due to surrounding structures.An apparatus, schematically shown in figure 5, was constructed to validate example 3 of the model and investigate the EC signal under multiparameter conditions.A similar apparatus was also constructed to validate example 1 of the model.A Nortec ® 600D was used to measure the change of probe impedance.This was done by measuring the analogue voltage of the output pins on the back of the instrument, using a National Instruments USB-6210 16 bit analogue-todigital converter.This voltage is always between −5 and +5 V, regardless of the instrument's settings [33], so the output was assumed to be proportional to the probe's impedance change, i.e.V out = c∆Z p for some constant c.The instrument was calibrated for impedance change measurements, and nulled (output set to 0) when the probe was in air for example 1 and inside the tube and far from a flaw for example 3. A 410 stainless steel collar was constructed to represent a support structure.Data was collected with the probe at 1 mm intervals as it passed through a collar and stainless steel 304 (non-magnetic SS 304) tube.The physical parameters of the collar and tube are shown in table 1.The conductivity of the SS 304 tube was measured using a four-point potential drop method.The test frequencies for example 1 presented here are 30, 60, 120, 240 and 180 kHz.The chosen test frequencies for example 3 were f 90 of the tube and multiples or fractions of it: f90 8 , f90 4 , f90 2 , f 90 , 2f 90 .For the tube presented here, these were 18, 37, 74, 150, 290 kHz.

FEM
To validate the models against numerical methods, each experiment was numerically simulated in FEM magnetics (FEMM) Version 4.2 [35], whose output was then compared to that of the model.Figure 6 shows the geometry and mesh in FEMM with the magnetic field output.

Results
To validate the model, data was collected using the method outlined in section 3. The results, plotted in figure 7, compare FEM and analytical calculations to measurements of   The results plotted in figure 8 show the median of three sets of data where the coil was entering and exiting the collar, in a stainless steel tube, representing example 3 of the model.
Below is a comparison of the root mean square difference (RMSD) between experimental values of the model example 3 versus analytical and FEM values.The values that were compared correspond to 46 positions.The analytical and FEM models were scaled to fit measurement data in a least-squares sense at each frequency separately.Table 2 shows a comparison of the RMSDs between analytical and experimental, FEM and experimental, and analytical and FEM at the various measured frequencies.From these results, FEM works better for lower frequencies, while the analytical model is more accurate at higher frequencies and overall.The analytical model assumes a slab extending to infinity in x-y horizontal directions, whereas the FEM model considers the actual experimental collar geometry.This effect is more important at low frequencies, where the fields extend further [36] and therefore, are more affected by the edge of the collar.Interestingly, the real SGs geometry more closely approximates an infinite slab.At higher frequencies the FEM model does not as easily accommodate the skin effect associated with the high permeability steel collar, this is also where differences in the two models are the greatest as shown by the last column in table 2. The nominal values for conductivity and permeability for the collar used in the analytical calculations produce impedance changes that agree with measurement, however, agreement may be improved with explicit measurements of these particular parameters.
Using the framework for derivation of analytical models for bobbin probe impedance in cylindrical geometries solutions that account for additional parameters could be developed.Tubesheets of CANDU ® nuclear reactors often have an overlay on the surface of the tubesheet.Modifying the external region of the model to include a thin conducting region on the surface of the slab could account for this.It may also be possible to model tubes with circumferential grooves.Building on the model, a cylindrical region could include a void between axial positions to model circumferential grooves, which could be developed using the TREE method [5].

Conclusion
Analytical models for the impedance of a bobbin probe that consisted of two coils, simulating actual probes used in SG inspection, passing through a conducting steel slab supporting layers of an arbitrary number of cylindrical conductors were derived.The method of derivation considered a TREE [5] of the magnetic vector potential in an axisymmetric geometry, representing the boundary conditions at each radial boundary with a matrix equation.Three examples were provided describing the response of (1) a coil due to a slab alone, (2) a supported tube, and (3) a supported tube with an airgap.The first and third examples were validated against experiment and FEM modelling.The second example is considered a special case of the third, so it was inferred that validation of the third example implied that the second was also valid.The analytical model provided rapid calculations that were in good agreement with FEM and experimental values, paving the way for a future implementation in an inverse algorithm.

Figure 1 .
Figure 1.Illustration of a boundary between two regions with, permeabilities µ 1 and µ 2 , respectively, defined with a unit normal vector, n.

Figure 2 .
Figure 2. Geometry of an arbitrary number of cylindrical conductors in a truncated domain.

Figure 3 .
Figure 3. 3D view of the tube-airgap-support geometry with two coils inside a supported tube with a thin gap between the tube and support.

Figure 4 .
Figure 4. Diagram of equation (72) and the search circles, shown in green.The black line is y = 1x , and the orange line is y = 1 2x .Dimensions have been normalised with respect to the skin depth of the support structure.

Figure 5 .
Figure 5. Cross section of the experimental apparatus for validation of the model used for example 3. The apparatus for example 1 omits the tube.

Figure 6 .
Figure 6.Geometry and mesh in FEMM for a probe in a tube with a collar (left) with magnetic field output (right).

Figure 7 .
Figure 7. Comparing example 1 of the model with experiment and FEM.Resistive (top left) and reactive change (top right) as the probe passes through a collar, with a Lissajous figure (bottom), with a blow-up of the lower frequencies to the right.

Figure 8 .
Figure 8. Comparing the model, example 3, with experiment and FEM.Resistive (top left) and reactive change (top right) as the probe passes through a collar outside of a tube, with a Lissajous figure (bottom).All plots show the probe exiting the collar.Measured data is the average of three sets.Lift-off is not in the horizontal direction. ) ↔U is (56) for M = 1 substituted into (55): ↔

Table 2 .
Comparison of root mean square differences (RMSD) between measurements of the model Example 3 compared to analytical (an) and FEM (fem) computations scaled to fit the measurements.canandc fem are the respective scaling factors.Frequency (kHz) RMSD (Vmeas − canZan) (V) RMSD (Vmeas − c fem Z fem ) (V) RMSD (canZan − c fem Z fem ) (V)The analytical model shows good agreement with FEM and experimental values.At high frequencies the discrepancies between analytical and FEM are likely due to the way in which FEMM accounts for high frequency effects.Within the scope of this paper, it is unclear whether these corrections are accurate.At low frequencies, the FEM model appears more accurate, while the analytical model is more accurate overall.The analytical model was computed 1000 times faster than the FEM model.A desktop computer with an AMD ® Threadripper™ 3970X 32-core processor at 4.15 GHz computed example 3 of the model at five frequencies and 46 positions in 21 min with FEM and in 1.2 s, analytically.