Hyperbolic Basis Functions for Time-Transient Analysis of Eddy Currents in Conductive and Magnetic Thin Sheets

This paper presents a new time-domain finite-element approach for modelling thin sheets with hyperbolic basis functions derived from the well-known steady-state solution of the linear flux diffusion equation. The combination of solutions at different operating frequencies permits the representation of the time-evolution of field quantities in the magnetic field formulation. This approach is here applied to solve a planar shielding problem in harmonic and time-dependent simulations for materials with either linear or nonlinear characteristics. Local and global quantities show good agreement with the reference solutions obtained by the standard finite element method on a complete and representative discretization of the region exposed to a time-varying magnetic field.


I. INTRODUCTION
T HIN sheets of high permeability and/or conductivity are often employed to mitigate stray fields produced by electric and electronic devices such as rotating machines, large power transformers, induction heating equipment, welding and forming machines [1]- [5]. In the surroundings of these devices, the field intensity needs to be at acceptable value to comply with the exposure limits for humans and for electromagnetic compatibility and interference reasons [1]. However, the shielding efficiency of thin sheets can be directly affected by its material characteristics, shape and position [4]. Consequently, the availability of models able to predict the electromagnetic (EM) behavior in such structures at an affordable computational cost is key to optimizing these devices.
In terms of shape, the high aspect ratio of thin-sheet shields presents a challenge to numerical simulation. Indeed, the direct application of a numerical method, such as the Finite Element Method (FEM), can be computationally expensive or even prohibitive due to the associated meshing difficulties [6]. On the one hand, a coarse mesh inside these sheets is unable to capture the EM phenomena and may lead to elements with high aspect ratio, which affect the FEM solution accuracy and convergence [7]. On the other hand, a high-density mesh can increase inordinately the number of unknowns in the problem and therefore the computational cost. The EM problem is even harder to solve if nonlinear materials characteristics are considered in time-transient analysis.
An efficient way to overcome these difficulties is to use the classical Thin-Shell (TS) model [8]- [11]. In this model, a reduced-dimension geometry replaces the actual volume of the thin regions, and suitable impedance boundary conditions (IBCs) account for the EM behavior within the original volume. These IBCs are defined from the analytical solution of the EM problem throughout the volume of the thin sheet. Thus, since the smallest dimension of the layer is neglected in a geometric sense, errors are avoided that would have been caused by the original anisotropic meshing of the thin structure with poor aspect ratios [2].
The TS model provides a good compromise between accuracy and computational cost [6], but its application is still mostly restricted to linear and harmonic regime analysis since the analytical solution is known a priori and the IBCs can be easily established. Currently available time-domain and nonlinear TS-FEM approaches are often derived from classical IBCs, whether using orthogonal polynomial basis functions to express the magnetic flux density through the shell thickness [12]- [17], Fast Fourier Transform (FFT) to update the residual from the harmonic solution [1], [4], or simply a linear field variation through its thickness (strip approximations) [18]- [21]. However, when considering timetransient analysis of nonlinear thin sheets, more representative models are required.
In [22], low-order surface impedance boundary conditions (SIBCs) are defined using basis functions derived from the steady-state analytical solution of semi-infinite slab problems. These SIBCs are applied in time-transient FE simulations to remove large conducting regions from the computational domain. Although the problems involving thin sheets are different, their nature is the same. In SIBCs, field quantities penetrate the surface from one boundary of a bulk domain, whereas in the TS model, the penetration occurs simultaneously from the two extended faces of a thin sheet. Moreover, the IBCs in the TS approach proposed in [9] are derived in a similar way than the SIBCs in [23]. For these reasons, a time-domain extension of the classical TS model with basis functions derived from the steady-state solution of a slab of finite thickness, equivalent to the model proposed in [22] for SIBCs, is a natural approach to pursue. This paper presents a novel time-domain extension of the classical TS model to solve 2-D shielding problems. The physics inside the thin region is captured by hyperbolic basis functions derived from the steady-state analytical solution of the linear flux diffusion equation. We demonstrate that the use of two hyperbolic basis functions leads to IBCs equivalent to the classical TS model in harmonic regime. In the time-transient analysis, the use of n pairs of hyperbolic basis functions, each representing different frequencies, coupled to FEM models outside the TS, allows computation of the time evolution of the physical quantities throughout the domain without resolving the TS. The method is here developed for a magnetic field (h-)formulation and extended to nonlinear cases. Results show good agreement with the 2-D FE reference solutions, with a greatly reduced number of degrees of freedom (DoFs) and therefore, at a lower computational cost.

II. 1-D FLUX DIFFUSION PROBLEM IN A SLAB OF FINITE THICKNESS
By assuming a thin region with a sufficiently high aspect ratio, the EM problem in a sheet can be formulated as a 1-D flux diffusion problem in a slab of finite thickness. In Fig. 1, we consider a thin sheet of thickness d whose normal is parallel to the y-axis. The tangential component of the magnetic field (h x ) is in the x-direction (Fig. 1a) and the tangential component of the electric field (e z ) is in the z-direction (Fig. 1b). The slab problem can be then formulated in terms of h x or e z , i.e.
where ρ is the electric resistivity (σ = 1/ρ) and µ is the magnetic permeability (ν = 1/µ). These expressions are derived from Faraday's and Ampere's laws, respectively. We assume that ρ and µ are constants and that we have harmonic boundary conditions (BCs), i.e.
whereĥ ± x andê ± z are respectively the magnetic and electric field magnitudes, and φ ± h and φ ± e their related phase shifts. Using a phasor representation (symbols with a bar), i.e.
where {.} is the real part of the argument. The solutions to (1)-(2) with the boundary conditions (5)-(6) are given complex notation (symbols with right arrow on top) and can be written as: where ψ ± (y) = sinh ad 2 ± ay sinh ( ad) , and a = 1+ δ ,  = √ −1, δ = 2/(µσω), ω = 2πf , and f is the operating frequency. Note that the functions ψ ± (y) appear in both (7) and (8). These functions are used later in this paper to define the basis functions required for the timetransient analysis of thin sheets.

A. Hyperbolic Basis Functions and Classical TS Model
We propose the use of basis functions obtained from steadystate solutions of the 1-D flux diffusion problem governed by (1)-(4). Taking as example the problem in terms of the magnetic field (1), together with the BCs (3), the steady-state solution for h x can be written as where the analytical expressions for θ ± c and θ ± s given in Table 4.2-II of [24]. Here, θ c and θ s are obtained directly from the real ( ) and imaginary ( ) parts of (9), i.e.
(13) Note that the solution (10) can be interpreted as a least squares approximation of h ± x in Ω s using the hyperbolic functions ψ ± (y). In addition, if we consider the harmonic solution (7) with δ d, ψ ± (y) in (9) reduces to ψ ± (y) which is equivalent to the Lagrange polynomials of first order defined across the thickness d of the sheet. Indeed, with δ d, the field quantities have a linear variation through the sheet thickness, and the functions ψ ± (y) can account for this behavior. To illustrate this, the functions θ ± c and θ ± s for δ d are plotted in Fig. 2. The equivalent solution in terms of the electric field can be obtained by replacingĥ ± x and φ ± h byê ± z and φ ± e in (10). However, in this paper, we are solely interested in magnetic field quantities and the h-formulation. The solution depending on the electric field would be useful, say, to implement the proposed approach in the magnetic vector potential (a-)formulation.
Besides, it can be demonstrated that the application of the ψ ± (y) functions as basis functions in the variational form of a 1-D finite element problem in the harmonic regime leads to the same IBCs used in the classical TS model [9], which are, using a vector representation (bold symbols), with η h = − ωµ a tanh ad 2 and η e = σ a tanh ad 2 . Expression (15) connects the discontinuity of the tangential components of the magnetic field to the mean value of the tangential electric field. This discontinuity is related to the total net current flowing in the sheet [8]. Moreover, equation (16) connects the discontinuity of the tangential electric field to the mean value of the tangential magnetic field, which is related to the amount of perpendicular flux absorbed in the plane of the sheet. When δ d, the coefficients 1 a tanh ad 2 in η h and η e can be approximated by d/2 [6].
The classical TS model in the form of IBCs has been used extensively to tackle problems involving thin regions in harmonic regime simulations [4]- [6], [8]- [10]. These IBCs were originally defined from the analytical solution for the field distribution, and the integration of the analytical expressions of the electric and magnetic current densities over the thickness of the thin sheet, which gives equivalent surface currents representatives [9]. However, to the best of our knowledge, no definition of the TS model in the form of the hyperbolic basis functions has been proposed in the literature before. These functions appear naturally in the solution (10). Therefore, the proposed approach can be easily extended to time-transient analysis, as described next.

B. Hyperbolic Basis Functions in Time-Transient Analysis
In time-transient analysis, we define n pairs of ψ ± k (y), where k is the harmonic rank relative to a fundamental frequency f 1 chosen in accordance with the problem to model, and 1 ≤ k ≤ n. The number of basis functions n is defined according to the frequency content of h x and the desired accuracy. Therefore, ψ ± k (y) is still defined by (9), with ω = 2πf k (which affects the values of δ and a). Then, according to (13), each ψ ± k (y) generates the even θ ± ck (y) and odd θ ± sk (y) functions, which we write as The first two cosines satisfy θ ± c1 (y = ±d/2) = 1 while the remaining functions θ ± ck (y) and θ ± sk (y) in (18)(19) vanish at the boundaries of the thin region (Γ ± ). This allows us to connect the 1-D equations to the exterior FE global system of equations; see Section III. Examples of the proposed basis functions with n = 3 are presented in Fig. 3 for δ ≤ d. For cases with δ d, the first two cosines functions are enough to represent the profile of h x in Ω s , since it has a linear variation throughout the thickness of the sheet (Fig. 2). Furthermore, the sine basis vanish everywhere.
The expansion of h x (y, t) in terms of (17)- (19) can be written in matrix form as where h ± ck (t) and h ± sk (t) are unknowns of the problem to model. Additionally, the 1-D variational form of the partial differential equation (1), disregarding homogeneous BCs, is where h x is the test function assumed to vanish at Γ ± s . The FE discretization of (23) by means of N = 4n basis functions θ p (y) and θ q (y), with p, q ∈ [1, N ], for h x and h x respectively, and assuming isotropic linear materials, leads to a system of equations, expressed in matrix form as where the elements of [S] and [M] are calculated as which can be evaluated numerically for each pair of basis functions θ p and θ q . Then, considering the implicit Euler scheme for the time-discretization of (24), coupled to (3), one obtains a system of algebraic equations to be solved at each time-step of the simulation.
The instantaneous loss L(t) in Joule is calculated as [25] where S pq is given by (25).
In the nonlinear case, the still isotropic resistivity ρ (or the magnetic permeability µ) in the variational form (23) can depend on the magnetic field intensity h x or its derivative ∂ y h x . The resulting nonlinear system of equations is solved by the Newton-Raphson (NR) iterative method, as presented in [22], but with integral terms evaluated over the thickness of the thin region, i.e., −d/2 ≤ y ≤ d/2. These integrals are solved numerically using the Legendre-Gauss quadrature at every iteration of the NR method.

III. FEM IMPLEMENTATION
We study the problem of a thin region Ω s embedded in a domain Ω = Ω c ∪ Ω C c , where Ω c and Ω C c denote respectively the conducting and non-conducting parts of Ω. As depicted in Fig. 4a, the exterior boundary of Ω (∂Ω = Γ) is composed of two complementary parts Γ h and Γ e (i.e. Γ = Γ h ∪ Γ e and Γ h ∩ Γ e = ∅) that may be necessary for symmetry or physical purposes such as connecting different subproblems via their common boundaries [26]. The thin region Ω s belongs to the conductive subdomain (Ω s ⊂ Ω c ) and its interior and exterior boundaries are Γ − s and Γ + s , respectively. When coupling the TS model with the FEM, Ω s is geometrically replaced by a surface located halfway between the original boundaries (Ω s → Γ s in Fig. 4b). In the variational form, by assuming distinct BCs on both sides of Γ s , we obtain interface integrals to couple with the TS model within Ω s . Thus, the weak form of the h-formulation, obtained from the weak form of Faraday's law, is defined as follows: ∀ h ∈ H 0 (curl, Ω), where h are test functions with n × h = 0 along Γ h , n is the outward unit normal vector on Γ, and (·, ·) Ω and ·, · Γ denote respectively the volume integral over Ω and the surface integral over Γ of the scalar product of their two arguments. Note that we assume that h ∈ H(curl, Ω) already satisfies the BC on n × h along Γ h .
The two last terms of (28) express the discontinuity of the tangential components of the electric field along the surface representing the thin region and we considered that n s = −n + Γs = n − Γs (see Fig. 4a). Besides these interface terms, the weak form (28) requires the duplication of the DoFs related to the surface of the thin region. In [6], the authors propose the decomposition of the field quantities into its continuous and discontinuous parts in order to avoid nodes and edges duplication. This decomposition is also applied in [12]- [14].
Here, nodes and edges of the thin surface are duplicated, but except for the nodes located at its extremities (e.g. points p 1 and p 2 in Fig. 4b). This creates a crack in the topological structure, and the non-conducting region becomes non-simply connected. The interfaces Γ + s and Γ − s share the nodes at their extremities, such that Γ s = Γ + s ∪ Γ − s , and the tangential components of the magnetic fields on these surfaces are connected by an 1-D FE problem in the thin direction of the sheet.
In order to include the 1-D problem in the weak form (28), we express the surface integral terms on Γ s in (28) by using the variational formulation of Faraday's law over Ω s , namely The volume integrals terms in this expression have opposite signs than those presented in (28), since they are on the right side of (29). In fact, here we should consider the outward unit normal vector of the boundary of Ω s , i.e. n ± s in Fig. 4a, but for the sake of simplicity, we denoted n + s = −n − s = n s (Fig. 4b). Inside the sheet, we assume that the local magnetic field is written as h x (x, y, z, t) = h x (x, z, t)ζ(y) and the test function as h x (x, y, z) = h x (x, z)ζ (y), with h x (x, z, t) and h x (x, z, t) tangential to Γ s , and ζ(y) and ζ (y) differentiable in the interval −d/2 ≤ y ≤ d/2. The volume integrals terms in (29) are then reduced to surface integrals terms as follows where Γ s = Γ + s ∪ Γ − s . Note that the expansion in (30) has been reduced to the 2-D case, so that h x is independent of z, i.e., h x (x, t), and the terms ∇ × h x and ∇ × h x vanish. In 3-D, these terms should be taken into account.
By choosing ζ and ζ as θ p and θ q , respectively, the integral terms in (30) (25) and (26). Finally, taking the weak form of the h-formulation (28) with the interface terms rewritten (29) using (30) and (31), then estimating for each degree of freedom in Γ s (x, z) the variation in y using (20) with a system of the form (24), we obtain a coupled system of equations for the magnetic field inside and outside the TS. The IBCs in the proposed TS model are obtained from (29)-(31). Note that, with a single pair of hyperbolic basis functions, these equations become equivalent to the IBCs (15) and (16) of the classical TS model.

IV. VALIDATION AND APPLICATION
We consider a 2-D planar shield (width l = 1 m and thickness d = 1 mm) placed over a pair of wires carrying a current ±I (Fig. 5). The conductors are 2x2 cm 2 separated by a distance of l 1 = 30 cm, and the distance between the conductors and the shield is l 2 = 10 cm. The free-space region is 4x4 m 2 . The coordinate system xyz is defined at the center of the shield geometry and we evaluate the magnetic field distribution along the lines AA (x = 0, y), BB (x, y = 10 cm) and CC (x = l/2 − l/100, y), and at points P 1 (x = 0,y = 10 cm), P 2 (x = l/4,0) and P 3 (x = l/2 − l/100,0).
The application of standard FEM using the h-formulation with a full 2-D representation of the shield gives the reference solution. The solutions obtained with the application of the TS model are here compared with the reference solution in terms of local field distributions, Joule losses and mesh simplification.
The relative difference (R) between the solutions is calculated as where TSS and FES are the TS and the reference FE solutions, respectively, and . 2 denote the Euclidean norm of the argument.
In terms of mesh parameters, we defined a structured rectangular mesh in Ω s with 12 elements across the shield thickness in the FE model ( Fig. 6a-left). The shield surface was discretized in 1 mm wide elements, and 100 mm wide elements were considered on the external boundary. Moreover, first-order edge elements were used in Ω.
With the described mesh configuration, the application of the TS model (Fig. 6a-right) with n = 1 represents a reduction in the total number of DoFs by 14.3% in comparison with the FE model. However, a more significant reduction can be achieved with the TS model since a coarser mesh can be   employed while maintaining a sufficiently high mesh quality and solution accuracy.
According to [7], the quality of the triangular mesh can be evaluated by computing the aspect ratio of the inscribed radius to the circumscribed radius of every triangle. For the meshes presented in Fig. 6a, the smallest aspect ratios are 0.408 and 0.513 for the 2-D FE and the TS models, respectively. If a coarser mesh is considered, e.g., elements of size 10 mm in the shield surface (Fig. 6b), these aspect ratios become 0.085 and 0.621, respectively. In this case, the low quality of the mesh in the standard FE may reduce the solution accuracy and its convergence. Therefore, with elements of 10 mm, the TS approach is preferable. It allows reducing the number of DoFs by 80.9% compared to mesh in the FE reference model while maintaining the initial mesh quality in the 2-D domain. Thus, 10 mm wide elements were used in the TS model.
All the models described in this paper were implemented in the open-source code Gmsh [27] and the solver GetDP [28]. Simulations were conducted on a personal computer with an Intel i7 2400 processor and 16 Gb of memory. Harmonic and time-transient simulations for different types of shields were performed.

A. Time-Harmonic Regime
Simulations were first performed in the harmonic regime for the sake of validation of the proposed TS model. The current I was set to 6 kA (current density of 15 A/mm 2 ) at an operating frequency of f = 50 Hz, and two shield configurations were studied: shield 1: µ r = 1 and σ = 1 MS/m ⇒ δ = 71.2 mm, shield 2: µ r = 1000 and σ = 10 MS/m ⇒ δ = 0.712 mm, where µ r = µ/µ 0 is the relative magnetic permeability with respect to µ 0 , the magnetic permeability of air.
The hyperbolic basis functions ψ ± in (9) were defined accordingly. With the first shield configuration, the basis functions are equivalent to those presented in Fig. (2), since δ d. With the second configuration, the hyperbolic basis functions were defined with δ = 0.712 mm (or δ = 1.40d, since d = 1 mm). For harmonic regime simulations, only one pair of hyperbolic basis functions was considered (n = 1).
Real and imaginary components of h x at points P 2 and P 3 inside the shields are presented in Fig. 7 and in Fig. 8, respectively. The field profiles from the TS model were obtained by evaluating the field intensities in the shield with the proposed hyperbolic basis functions. Excellent agreement with reference solutions in terms of real and imaginary components of h x were observed for both shield configurations. The excellent agreement at point P 3 shows that the proposed model can provide accurate solutions also near the extremities of the shield. Therefore, even though no special consideration has been made at its extreme points, edge effects are correctly represented.
In Fig. 9, we present a shaded plot of the magnetic flux density. With the first shield configuration (Fig. 9a), the replacement of the original 2-D region by a thin sheet has no noticeable impact on the magnetic flux density distribution. Indeed, from Fig. 7 and Fig. 8, we observe that h x is almost constant inside shield 1. The tangential components of the field are continuous on the surface representation in the TS model. However, with the second shield configuration (Fig. 9b), the physics inside the plate produces a discontinuity of h x that  deforms the flux lines in the air surrounding the edge. Thus, in both shield configurations, the solution from the TS model agrees with the FE solution in terms of field distributions inside and outside the shield. The relative differences for the local magnetic field along lines AA , BB and CC , and at points P 2 and P 3 were computed with (32). The R-values are summarized in Table I. The maximum difference is 2.46% and occurs on line CC and point P 3 with shield 2. This difference may be related to the geometrical difference between the TS and the FE models. Since the thickness of the shield is not represented in the TS model, it is expected to observe at least a slight difference near its extremities.

B. Time-Transient Regime (Linear Case)
Time-transient simulations of the planar shield were also performed. In this study, we applied a pulsed current source in the wires whose waveform produces the magnetic field shown in Fig. 10 and 11 (black dashed lines) in absence of the shield. The amplitude of the current I is the same as in the harmonic regime cases (|I|= 6 kA), and the rise time (t r )  was set to 20 µs. The simulation time was t max = 50 µs with a time-step of t max /120. The implicit Euler scheme was used.
In terms of material composition, two new shielding configurations were studied: shield 3: µ r = 1000 and σ = 1 MS/m, shield 4: µ r = 100 and σ = 10 MS/m. Shields 3 and 4 are both ferromagnetic, but shield 4 is more conductive than shield 3. The simulation time was chosen to define the fundamental frequency as f = 1/(4t max ), i.e. f = 5 kHz. Therefore, these shields have same penetration depth δ = d/4.44. Consequently, the same hyperbolic basis functions can be used to tackle these problems. The first set of basis functions was defined with frequency equal to the fundamental (f 1 = f ). Additional basis functions were then calculated using odd harmonic frequencies of f 1 , i.e., f k /f 1 = 2k − 1, with 1 ≤ k ≤ n.
In Fig. 10 and 11, the time-evolution of h y at point P 1 was compared with the reference solution for shield 3 and 4, respectively. The number of basis functions in the TS model was varied from 1 to 3, and the solutions approached the reference solution as n increased.  The relative differences of the profiles of h y in the TS model to the reference solution as a function of the number of basis functions n are presented in Fig. 12. The maximum relative difference R decreases from 245.8% with n = 1 to 2.95% with n = 2, and to less than 2% for n ≥ 3. Simulations with n > 3 show little improvement in terms of solution accuracy. This is mainly due to the geometrical differences between the FE and the TS models, as discussed in the time-harmonic regime case. Despite this, the application of the TS model shows a good compromise between computational cost and solution accuracy.
The number of DoFs, the computation time and the Joule losses in the reference and TS models are summarized in Table II. Since a coarser mesh was applied in the TS model, simulations are more than five times faster with this approach than with standard FE. The number of DoFs is nearly independent of n, and is reduced by more than 70%, even with n = 5. Despite this, the relative error in losses estimation is less than 2% with n ≥ 3 in shields 3 and 4.

C. Time-Transient Regime (Nonlinear Case)
In an attempt to extend the proposed TS model to nonlinear analysis, the shielding problem was also studied for shields with nonlinear magnetic permeability (µ = µ(h)). A sinusoidal supply current of amplitude |I| = 6 kA at an operating frequency f = 1 kHz was applied to the wires. The effects of the saturation and the influence of the number basis functions in the proposed TS model were analyzed in terms of solution accuracy. The material properties were modelled with an isotropic saturation law expressing the magnetic permeability as a function of the magnetic field intensity, i.e., where µ r0 is the relative permeability at origin and m 0 the saturation magnetic field in A/m. The differential permeability required for the application of the NR-scheme was defined as in [29]. We carried out simulations with µ 0 m 0 = 1.31 and µ r0 = 12500. The B-H curve is presented in Fig. 13. Furthermore, the electrical conductivity of the shield was fixed at σ = 1 MS/m. One time period was simulated, i.e., t max = 1/T , with T = 1/f . Moreover, the time-step was set to ∆t = T /120, and the maximum number of iterations for the NR-scheme was set to 12 in both the reference and the proposed TS models. Besides, the number of points used in the Legendre-Gauss quadrature was 20 points. This number of points is considered sufficiently high to avoid errors related to the numerical integration of the hyperbolic functions across whole the thickness of the shield. Depending on the penetration depths of the basis functions, less integration points could be used, but we kept 20 at all times to remains on the safe side.
The first set of hyperbolic basis functions was defined by taking f 1 = f and higher order basis functions that are odd multiples of f 1 . Furthermore, the magnetic permeability used in the definition of the basis functions was taken from the B-H curves corresponding to a flux density b = 1.2 T, i.e., µ r = 1000 for µ r0 = 12500. A similar approach was used in [22] to parametrize basis functions for nonlinear SIBCs. Fig. 14 shows the h x profile throughout the thickness of the shield for three specific simulation times (T /8, T /4 and T /2). Results are presented for n = 1 to 3 and compared with the 2-D FE solution. Note that the accuracy of the proposed TS model clearly improves with n. Since the penetration depth used to define the hyperbolic basis functions is at the same time inversely proportional to the square root of f and µ r , the additional frequency components can be interpreted as a way to consider an increase of the magnetic permeability. For  this reason, the saturation effects observed at t = T /8 and t = T /2 are well represented with the proposed TS model when considering higher harmonic components. The relative difference between the solutions at t = T /8 is reduced to less than 1% with n = 3 (Fig. 15).
The time-evolution of h y at point P 1 is shown in Fig. 16. The solution of the 2-D FE problem without the shield gives the field at this point, which has the same waveform as the current I applied to the wires. The solution for a linear shield problem with µ r = 12500 is presented for the sake of comparison with the nonlinear solution. Finally, the solution for the shield modeled with the TS model shows improvement as n increases, while high accuracy is observed when compared to the nonlinear reference solution. The saturation effect is clearly observed.
The relative difference of the TS model to the FE solution at P 1 at each time step is presented in Fig. 17. It was reduced from more than 5% with n = 1 to less than 1% with n 2. With a suitable choice of the basis functions, the proposed TS model can certainly be an option for simulating nonlinear thin sheets.

V. CONCLUSION
In this paper, a time-domain extension of the classical TS model for thin sheets was elaborated and discussed using the h-formulation. In our approach, the addition of n sets of basis functions derived from the steady-state solutions for the problem of a slab of finite thickness permits the representation of the time evolution of the field quantities inside the thin region and its surroundings. We did apply this method to analyze the shielding efficiency of conducting and ferromagnetic planar sheets in harmonic and time-transient regimes for linear and nonlinear shield configurations.
In the harmonic regime, the proposed model gives IBCs equivalents to those used in the classical TS model, since both methods are based on the solutions (7) and (8) and are directly included in the FE formulation. Our model, however, can also be used in time-transient FE analysis. We showed that, by adding a small number of hyperbolic basis functions, high precision can be achieved. In all cases, the relative differences were < 3% with the reference solutions with n 3, including the critical region near the extremities of the shield. Furthermore, the proposed model can achieve comparable errors with less degrees of freedom, and hence at a lower computational cost, while also avoiding meshes with poor aspect ratios.
The TS model presented in this paper is still applicationdependent since the set of hyperbolic basis functions must be defined according to the frequency content of the magnetic field inside the thin region, as well as the material composition of the sheet. Although we did not yet find a general rule to select the basis functions, the latter are easy to derive since they originate from the analytic solutions of the 1-D linear flux diffusion problem in harmonic regime. Therefore, as long as the set of basis functions is rich enough to represent a diversity of penetration depths (which varies dynamically with local magnetic saturation), nonlinear solutions can be well approximated with this approach. Finally, the presented methodology can be easily extended to other FE formulations, such as the a-formulation, as well as 3-D shielding problems.