## 2 Model Formulation  
The methodology described in this paper is based on a magnetoquasistatic formulation (i.e., \(\frac{\partial\mathbf{D}}{\partial t}=0\)) of the physical problem. The formulation holds if the condition \(\sigma\gg(2\pi/\tau)\epsilon\) is verified for the modeled conductors, where \(\sigma\) and \(\epsilon\) are the electrical conductivity and permittivity, respectively, and \(\tau\) the characteristic time of the considered electrical phenomena. The current densities and the magnetic vector potential were assumed to have a preferential direction, i.e., \(z\). Hence, these can be defined as \(\mathbf{J}=J_{z}\hat{\mathbf{k}}\) and \(\mathbf{A}=A_{z}\hat{\mathbf{k}}\). By neglecting the spatial derivatives \(\frac{\partial}{\partial z}\) of \(\mathbf{J}\) and \(\mathbf{A}\) along \(z\), the problem can be formulated in the two-dimensional plane perpendicular to \(z\). The governing equation is:  
\[-\mathbf{\
abla}\cdot\left(\frac{1}{\mu}\mathbf{\
abla}A_{z}\right)+\sigma\frac{ \partial A_{z}}{\partial t}=J_{0,z}, \tag{1}\]  
where \(\mu\) is the magnetic permeability. The system's forcing term \(J_{0,z}=-\sigma\Delta\sigma^{\prime}\) can be regarded as the current density that would flow in the conductors in the absence of the electromagnetic induction mechanism as a response to a per-unit length voltage \(\Delta\sigma^{\prime}\) applied along \(z\). Assuming that all materials behave linearly, a Laplace transform can be applied to Equation (1) to obtain a differential equation in the complex frequency \(s\) domain:  
\[-\mathbf{\
abla}\cdot\left(\frac{1}{\mu}\mathbf{\
abla}\tilde{A}_{z}\right)+s\sigma \tilde{A}_{z}=\tilde{J}_{0,z}+\sigma A_{z}(0^{-}). \tag{2}\]  
Equation (2) is discretized by means of an FEM approach. The unknown is approximated by means of a piecewise polynomial representation \(\tilde{A}_{z}(x,y)=\{N\}\{\tilde{A}_{z}\}\), where \(\{N\}\) is the set of shape functions and \(\{\tilde{A}_{z}\}\) is an array constituted by the nodal values of \(\tilde{A}_{z}\). Using the generic shape function \(N_{k}\) as the weighting function, the weak formulation resulting from Equation (2) is obtained from the weighted residual on the calculation domain \(\Omega\) bounded by the closed curve \(\partial\Omega\):  
\[\int_{\Omega}\mathbf{\
abla}N_{k}\cdot\left(\frac{1}{\mu}\mathbf{\
abla}\tilde{A}_{z} \right)\mathrm{d}S+s\int_{\Omega}N_{k}\sigma\tilde{A}_{z}\mathrm{d}S=\int_{ \Omega}N_{k}\tilde{J}_{0,z}\mathrm{d}S+\int_{\Omega}N_{k}\sigma A_{z}(0^{-}) \mathrm{d}S+\oint_{\partial\Omega}N_{k}\frac{1}{\mu}\frac{\partial\tilde{A}_{ z}}{\partial n}\mathrm{d}l. \tag{3}\]  
From Equation (3), a complex linear system is derived:  
\[[\mathrm{K}(s)]\{\tilde{A}_{z}\}=\big{\{}RHS\big{[}\tilde{J}_{0,z},A_{z}(0^{- })\big{]}\big{\}}. \tag{4}\]  
The coefficient matrix \([K]\) in Equation (4) is a sparse complex matrix. The values of \([\mathrm{K}]\) depend on the properties (i.e., electrical conductivity and magnetic permeability) of the materials and on the complex frequency \(s\). The right-hand-side term \(\{RHS\}\) is a function of the transformed imposed current densities \(\tilde{J}_{0,z}\), the initial distribution \(A_{z}(0^{-})\) of the magnetic vector potential, and the boundary conditions. That is, for a given set of initial and boundary conditions, Equation (4) allows the distribution of the vector potential \(\tilde{A}_{z}\) to be found at a generic complex frequency \(s\) as the response of a forcing term \(\tilde{J}_{0,z}\). Once Equation (4) has been solved, the current density in the complex frequency domain can be computed as:  
\[\tilde{J}_{z}=\tilde{J}_{0,z}-s\sigma\tilde{A}_{z}+\sigma A_{z}(0^{-}). \tag{5}\]  
The electric current flowing through a conductor is then found by integrating the current density computed by means of Equation (5) over the conductor cross-section.