Finite Element Modeling of Power Cables
using Coordinate Transformations

Albert Piwonski1, Julien Dular2, Rodrigo Silva Rezende1, Rolf Schuhmann1 1Theoretische Elektrotechnik, Technische Universität Berlin, Berlin, Germany 2 TE-MPE-PE, CERN, Geneva, Switzerland

Abstract

Power cables have complex geometries in order to reduce their ac resistance. Although there are many different cable designs, most have in common that their inner conductors’ cross-section is divided into several electrically insulated conductors, which are twisted over the cable’s length (helicoidal symmetry). In previous works, we presented how to exploit this symmetry by means of dimensional reduction within the $\mathbf{H}-\varphi$ formulation of the eddy current problem. Here, the dimensional reduction is based on a coordinate transformation from the Cartesian coordinate system to a helicoidal coordinate system. This contribution focuses on how this approach can be incorporated into the magnetic vector potential based $\mathbf{A}-v$ formulation.

Index Terms:

Coordinate transformations, dimensional reduction, eddy currents, finite element modeling, magnetic vector potential formulation, power cables, tree-cotree gauging.

I Introduction

Power cables are essential elements in the transmission chain of electric power from generator to consumer. Special inner conductor designs have been developed for ac operation to minimize the undesirable current displacement caused by the skin and proximity effect [1]. Although there are many different cable designs, most have in common that their inner conductors’ cross-section is devided into several conductors, which are twisted and electrically insulated from each other (see Fig. 1).

Investigating ac losses in power cables is an important but challenging task: Solving numerically eddy current boundary value problems (BVP) in 3-D, which model the cable’s electromagnetic behaviour in the magnetoquasistatic limit, leads to tremendous computational efforts, due to arising multiscale problems (e.g., thin insulation layers).

However, if one models the power cable as a symmetric BVP, computational costs can be scaled down significantly by means of dimensional reduction. Due to the conductors’ twist, however, no conventional translational symmetry is valid, as recently shown in a formal framework using Lie derivatives [3]. In particular though, when choosing proper boundary conditions, an eddy current BVP posed on a domain as in Fig. 1 has a so-called helicoidal symmetry. In contrast to applying periodic boundary conditions, here the model can be solved in 2-D. This special symmetry was exploited before to calculate hysteresis losses in twisted superconductors [4, 5], but also for solving full Maxwell’s equations in twisted waveguides [6]. Further, in previous works, we presented how to incorporate this symmetry property within the $\mathbf{H}-\varphi$ finite element formulation applied to a power cable problem [7].

This contribution focuses on an integration into the magnetic vector potential based $\mathbf{A}-v$ formulation: In Section II, we define the coordinate transformation and its inverse that allow a bidirectional transition between Cartesian and helicoidal coords. Subsequently, Section III presents how a 3-D $\mathbf{A}-v$ formulation in Cartesian coords can be reformulated equivalently into a 2-D formulation in helicoidal coords. Further, implementation details of the 2-D model are given. In Section IV, we compare the results of this contribution and of our previous work [7] with a 3-D reference model.

Refer to caption — Figure 1: Generic model of a cable’s inner conductor: Modeled using [2].

II Coordinate transformations

In computational electromagnetism, symmetries occur quite differently: E.g., when computing electrical machines, often a translational symmetry is assumed. Further, azimuthal symmetries appear frequently when modeling cylindrical waveguides. The common key idea is the use of a coordinate system in which the geometries appear the same in one direction. Then, partial derivatives w.r.t. to this particular direction of electromagnetic field quantities are either vanishing or constant, leading to drastically reduced computational efforts. Therefore, for describing twisted conductors, a helicoidal coordinate sytem is the most suitable.

In the following, we denote points represented in the Cartesian coordinate system $(x,\,y,\,z)$ as $\mathbf{p}_{xyz}\coloneqq[x,\,y,\,z]^{\top}$ , whereas points represented in the helicoidal coordinate system $(u,\,v,\,w)$ are denoted as $\mathbf{p}_{uvw}\coloneqq[u,\,v,\,w]^{\top}$ . The birectional change of coordinates is achieved by the map $\boldsymbol{\phi}:\Omega_{xyz}\rightarrow\Omega_{uvw}$ and its inverse $\boldsymbol{\phi}^{-1}:\Omega_{uvw}\rightarrow\Omega_{xyz}$ , where $\Omega_{xyz},\,\Omega_{uvw}\subset\mathbb{R}^{3}$ :

\displaystyle\boldsymbol{\phi}(\mathbf{p}_{xyz})

\displaystyle=\mathbf{p}_{uvw}=\begin{bmatrix}+x\cos(z\alpha/\beta)+y\sin(z\alpha/\beta)\\ -x\sin(z\alpha/\beta)+y\cos(z\alpha/\beta)\\ +z\end{bmatrix},

(1)

\displaystyle\boldsymbol{\phi}^{-1}(\mathbf{p}_{uvw})

\displaystyle=\mathbf{p}_{xyz}=\begin{bmatrix}+u\cos(w\alpha/\beta)-v\sin(w\alpha/\beta)\\ +u\sin(w\alpha/\beta)+v\cos(w\alpha/\beta)\\ +w\end{bmatrix}.

(2)

Here, parameters $\alpha$ , $\beta$ are related to the number of turns and to the total longitudinal length of the helical object of interest, i.e., for different geometries $\boldsymbol{\phi}$ resp. $\boldsymbol{\phi}^{-1}$ are defined differently as well. The effect of the global coordinate transformation is demonstrated in Fig. 2.

It is important to note, that the maps (1) & (2) allow a unique bidirectional transition of points between both coordinate systems, s.t. the diagram in Fig. 3 commutes. This property ensures that also the electromagnetic field quantities have a unique representation in the other coordinate system, which will be needed in the following.

III $\mathbf{A}-v$ formulation in different coordinates

In previous works, we presented how to exploit helicoidal symmetries, in the context of eddy current problems, by means of dimensional reduction within the $\mathbf{H}-\varphi$ finite element formulation [7]. In this paper, on the other hand, we present how to incorporate this approach into the more widespread $\mathbf{A}-v$ formulation.

From a topological point of view, the latter formulation is simpler, since typically no topological pre-processing of the computational domain $\Omega$ is needed. In a nutshell, this can be explained by the fact that the magnetic flux density $\mathbf{B}$ in $\Omega$ is closed ( $\mathrm{div}\,\mathbf{B}=0$ ) and exact ( $\exists\,\mathbf{A}$ , s.t. $\mathbf{B}=\mathbf{curl}\,\mathbf{A}$ ) [8]. In contrast, within the $\mathbf{H}-\varphi$ formulation, we consider the magnetic field $\mathbf{H}$ in the insulating domain $\Omega_{i}$ , which is also closed ( $\mathbf{curl}\,\mathbf{H}=\mathbf{0}$ ) as there are no currents, but not exact ( $\nexists\,\varphi$ , s.t. $\mathbf{grad}\,\varphi=\mathbf{H}$ ). This means in particular, if one represents the magnetic field purely as the gradient of a scalar potential $\varphi$ , Ampere’s law would be violated in the presence of conductors, since then no net circulations of $\mathbf{H}$ along closed curves could be generated.

III-A Weak formulation in Cartesian coordinates

In the following, we present the weak $\mathbf{A}-v$ formulation in a concise manner. More details can be found in [5]. In the Cartesian domain $\Omega_{xyz}$ with boundary $\partial\Omega_{xyz}$ , conducting subdomain $\Omega_{xyz,c}$ ( $N$ connected components), and in frequency domain with frequency $\omega/2\pi$ , the $\mathbf{A}-v$ formulation consists of two sets of equations:

	$\displaystyle\left(\nu_{xyz}\mathbf{curl}\,\mathbf{A}_{xyz},\,\mathbf{curl}\,\mathbf{A}_{xyz}^{\prime}\right)_{\Omega_{xyz}}$
$\displaystyle+$	$\displaystyle\left(j\omega\sigma_{xyz}\,\mathbf{A}_{xyz},\,\mathbf{A}_{xyz}^{\prime}\right)_{\Omega_{xyz,c}}-\langle\widetilde{\mathbf{H}}_{xyz}\times\mathbf{n},\,\mathbf{A}_{xyz}^{\prime}\rangle_{\partial\Omega_{xyz}}$
$\displaystyle+$	$\displaystyle\left(\sigma_{xyz}\,\mathbf{grad}\,{v}_{xyz},\,\mathbf{A}_{xyz}^{\prime}\right)_{\Omega_{xyz,c}}=0,$	(3)
	$\displaystyle\left(j\omega\sigma_{xyz}\,\mathbf{A}_{xyz},\,\mathbf{grad}\,v_{xyz}^{\prime}\right)_{\Omega_{xyz,c}}$
$\displaystyle+$	$\displaystyle\left(\sigma_{xyz}\,\mathbf{grad}\,v_{xyz},\,\mathbf{grad}\,v_{xyz}^{\prime}\right)_{\Omega_{xyz,c}}=\sum_{i=1}^{N}I_{i}V_{i}^{\prime},$	(4)

with the magnetic vector potential $\mathbf{A}_{xyz}$ and electric scalar potential $v_{xyz}$ (with test functions $\mathbf{A}_{xyz}^{\prime}$ , resp. $\mathbf{grad}\,{v}_{xyz}^{\prime}$ ), magnetic reluctivity $\nu_{xyz}$ and electric conductivity $\sigma_{xyz}$ . The term $\widetilde{\mathbf{H}}_{xyz}\times\mathbf{n}$ denotes the fixed tangential magnetic field at the boundary of the domain. Further, the notations $(\mathbf{f},\,\mathbf{g})_{\Omega_{xyz}}$ and $\langle\mathbf{f},\,\mathbf{g}\rangle_{\partial\Omega_{xyz}}$ denote $\int_{\Omega_{xyz}}\mathbf{f}\cdot\mathbf{g}\,\mathrm{d}V_{xyz}$ and $\int_{\partial\Omega_{xyz}}\mathbf{f}\cdot\mathbf{g}\,\mathrm{d}S_{xyz}$ .

The electric scalar potential $v_{xyz}$ is only defined in the conducting domain and is further used to associate global quantities to the $i$ -th conductor, namely, its total current $I_{i}$ and its voltage drop $V_{i}$ [5]. In $\Omega_{xyz,c}$ , the magnetic vector potential is directly linked to the physical electric field $\mathbf{E}_{xyz}=-j\omega\mathbf{A}_{xyz}-\mathbf{grad}\,v_{xyz}$ , s.t. here no gauging is needed (modified vector potential). However, in the insulating domain $\Omega_{xyz,i}$ , gauging is necessary to find a unique solution.

III-B Weak formulation in helicoidal coordinates

The key idea for the dimensional reduction is the following: Instead of computing the integrals in (3) & (4) on the Cartesian domain $\Omega_{xyz}$ , we transform them on the domain $\Omega_{uvw}$ , i.e., we shift their computation to a $w$ -invariant space (see Fig. 2).

This transformation requires a reexpression of the involved electromagnetic field quantities. The one-forms $\mathbf{A}_{xyz}$ , $\mathbf{grad}\,v_{xyz}$ and $\widetilde{\mathbf{H}}_{xyz}$ and the two-form $\mathbf{curl}\,\mathbf{A}_{xyz}$ (magnetic flux density $\mathbf{B}_{xyz}$ ), transform as follows [5, 8]:

	$\displaystyle\mathbf{A}_{xyz}(\mathbf{p}_{xyz})$	$\displaystyle=J_{\boldsymbol{\phi}^{-1}}^{-\top}\ \mathbf{A}_{uvw}(\boldsymbol{\phi}(\mathbf{p}_{xyz})),$		(5)
	$\displaystyle\ \mathbf{curl}\,\mathbf{A}_{xyz}(\mathbf{p}_{xyz})$	$\displaystyle=\frac{J_{\boldsymbol{\phi}^{-1}}}{\mathrm{det}(J_{\boldsymbol{\phi}^{-1}})}\ \mathbf{curl}\,\mathbf{A}_{uvw}(\boldsymbol{\phi}(\mathbf{p}_{xyz})),$		(6)

where $J_{\boldsymbol{\phi}^{-1}}$ denotes the Jacobian of $\boldsymbol{\phi}^{-1}$ evaluated at point $\mathbf{p}_{uvw}$ . Here, it is important to mention, that the $\mathbf{curl}$ operator on the right hand side of eq. (6) is not the actual differential operator in the helicoidal coordinate system. Rather, this operator is to be understood as blindly applied (applied as in Cartesian coords with a relabeling of variables and components) – its results are then directly mapped back into the Cartesian coordinate system as the equality sign in eq. (6) holds componentwise for vector fields expressed in Cartesian coordinates. Also, the qualitatively same statement applies, when inserting $\mathbf{grad}\,v_{xyz}$ into formula (5). Changing variables also introduces, rather conceptually than computationally, a factor $\mathrm{det}(J_{\boldsymbol{\phi}^{-1}})=1$ in the integrals in eq. (3) & (4).

In terms of the $(u,\,v,\,w)$ coordinates, introducing trial function spaces $A(\Omega_{uvw})$ & $V(\Omega_{uvw})$ and test function spaces $A_{0}(\Omega_{uvw})$ & $V_{0}(\Omega_{uvw})$ which will be defined in Section III-C, we can reformulate the weak formulation as follows: Seek $\mathbf{A}_{uvw}\in A(\Omega_{uvw})$ and $\mathbf{grad}\,v_{uvw}\in V(\Omega_{uvw,c})$ , s.t. $\forall\,\mathbf{A}_{uvw}^{\prime}\in A_{0}(\Omega_{uvw})$ and $\forall\,\mathbf{grad}\,v_{uvw}^{\prime}\in V_{0}(\Omega_{uvw,c})$ :

	$\displaystyle\left(\boldsymbol{\nu}_{uvw}\mathbf{curl}\,\mathbf{A}_{uvw},\,\mathbf{curl}\,\mathbf{A}_{uvw}^{\prime}\right)_{\Omega_{uvw}}$
$\displaystyle+$	$\displaystyle\left(j\omega\boldsymbol{\sigma}_{uvw}\,\mathbf{A}_{uvw},\,\mathbf{A}_{uvw}^{\prime}\right)_{\Omega_{uvw,c}}-\langle\widetilde{\mathbf{H}}_{xyz}\times\mathbf{n},\,\mathbf{A}_{uvw}^{\prime}\rangle_{\partial\Omega_{uvw}}$
$\displaystyle+$	$\displaystyle\left(\boldsymbol{\sigma}_{uvw}\,\mathbf{grad}\,{v}_{uvw},\,\mathbf{A}_{uvw}^{\prime}\right)_{\Omega_{uvw,c}}=0,$	(7)
	$\displaystyle\left(j\omega\boldsymbol{\sigma}_{uvw}\,\mathbf{A}_{uvw},\,\mathbf{grad}\,v_{uvw}^{\prime}\right)_{\Omega_{uvw,c}}$
$\displaystyle+$	$\displaystyle\left(\boldsymbol{\sigma}_{uvw}\,\mathbf{grad}\,v_{uvw},\,\mathbf{grad}\,v_{uvw}^{\prime}\right)_{\Omega_{uvw,c}}=\sum_{i=1}^{N}I_{i}V_{i}^{\prime},$	(8)

where the effect of the change of variables is fully contained in two anisotropic, $w$ -invariant material parameters, written as tensors:

	$\displaystyle\boldsymbol{\sigma}_{uvw}(\mathbf{p}_{uvw})$	$\displaystyle=\sigma_{xyz}(\boldsymbol{\phi}^{-1}(\mathbf{p}_{uvw}))J_{\boldsymbol{\phi}^{-1}}^{-1}J_{\boldsymbol{\phi}^{-1}}^{-\top}\mathrm{det}(J_{\boldsymbol{\phi}^{-1}}),$		(9)
	$\displaystyle\boldsymbol{\nu}_{uvw}(\mathbf{p}_{uvw})$	$\displaystyle=\nu_{xyz}(\boldsymbol{\phi}^{-1}(\mathbf{p}_{uvw}))J_{\boldsymbol{\phi}^{-1}}^{\top}J_{\boldsymbol{\phi}^{-1}}/\mathrm{det}(J_{\boldsymbol{\phi}^{-1}}).$		(10)

Now, since neither the integrands nor the computational domain $\Omega_{uvw}$ depend on the $w$ -coordinate, one may solve the problem on any $uv$ -plane for a fixed value of $w$ . For simplicity, we choose $w=0$ , because then $x=u$ and $y=v$ , see eq. (2).

III-C Space discretization and implementation details

We solve the systems of equations (7) & (8) by using the open-source finite-element software GetDP [9], which allows for flexible function space definitions, whereas the meshing process is performed by Gmsh controlled via the Julia API [10, 11]. In previous works, we derived the 2-D computational domain (symmetry cell) analytically using the theory of envelopes [7]. Now, the cross-section is computed as the intersection of 3-D helicoidal conductors with a plane using built-in CAD routines. This will allow us to investigate elliptically shaped conductor cross-sections in future work, which can be created under mechanical stress in the cable manufacturing process [1].

Although we solve the problem on a 2-D domain, we still assume, that the magnetic vector potential $\mathbf{A}_{uvw}$ has full three components, which are interpolated separately: The in-plane components ( $A_{u}$ & $A_{v}$ ) are discretized using 2-D Whitney edge functions. In the domain $\Omega_{uvw,i}\setminus\partial\Omega_{uvw,c}$ , though, only degrees of freedom (DOF) associated with the cotree of the mesh are non-zero, i.e., a tree-cotree gauge is applied to achieve a unique solution [12]. This defines the function space $A_{uv}(\Omega_{uvw})$ .

Further, the component $A_{w}$ is spanned by nodal functions and is essentially fixed to zero on $\partial\Omega_{uvw}$ as we model the cable’s shield as a perfect electrically conductive cylinder. Since we do not introduce an explicit gauge for $A_{w}$ , this component is implicitly Coulomb-gauged. This constitutes the function space $A_{w}(\Omega_{uvw})$ , s.t. in total $A(\Omega_{uvw})=A_{uv}(\Omega_{uvw})\oplus A_{w}(\Omega_{uvw})$ . The function $\mathbf{grad}\,v_{uvw}\in V(\Omega_{uvw})$ consists of $N$ $w$ -directed constants (one constant per connected component of $\Omega_{uvw,c}$ ).

IV Results

The results of the 2-D $\mathbf{A}-v$ formulation and of our previously presented 2-D $\mathbf{H}-\varphi$ formulation [7] are compared with a 3-D reference cable model that is implemented using the commercial software CST Studio Suite [2], referred to as CST. Each of the 13 helicoidal conductors (longitudinal length $0.2\,\text{m}$ , cross-section radius $5\,\text{mm}$ ) carries a total current of amplitude $\sqrt{2}/13\,\mathrm{A}$ at $f=50\,\text{Hz}$ . In all models, we considered annealed copper for the conductors’ material (conductivity $\approx 58.12\cdot 10^{6}\,\text{S}/\text{m}$ ) and further assumed a non-magnetic material in the whole domain ( $\nu_{0}=\mu_{0}^{-1}$ with $\mu_{0}=4\pi\cdot 10^{-7}\,\text{H}/\text{m}$ ). Since the current constraint is differently implemented in the 3-D model (current ports located at the cuboid bounding box of the domain), we compare local quantities only at the cable’s center (see Fig. 4).

In the following, we compare the 2-D models discretized by $58.64\text{k}$ triangles with a 3-D model discretized by $1.19\text{M}$ tetrahedra: We used first order finite elements resulting in $98\text{k}$ DOF in the 2-D $\mathbf{A}-v$ formulation, $57\text{k}$ DOF in the 2-D $\mathbf{H}-\varphi$ formulation and $1.21\text{M}$ DOF in the 3-D model. The DOF of the 2-D models differ so much on the same mesh, since only within the $\mathbf{H}-\varphi$ formulation we can directly incorporate and exploit the fact that the $w$ -component of the magnetic field $\mathbf{H}_{uvw}$ is a constant in the insulating domain $\Omega_{uvw,i}$ [7]. So far, no way has been found to utilize this knowledge also in the $\mathbf{A}-v$ formulation. The 2-D $\mathbf{A}-v$ model outputs the finite-element approximated vector potential $\mathbf{A}_{uvw}$ and the gradient of the electric scalar potential $\mathbf{grad}\,v_{uvw}$ .

Using transformation rules (5) & (6), all electromagnetic field quantities in the Cartesian coordinate system can be derived, e.g.:

	$\displaystyle\mathbf{J}_{xyz}$	$\displaystyle=-\sigma_{xyz}J_{\boldsymbol{\phi}^{-1}}^{-\top}\,\left(j\omega\mathbf{A}_{uvw}+\mathbf{grad}\,v_{uvw}\right),$		(11)
	$\displaystyle\mathbf{H}_{xyz}$	$\displaystyle=\mu_{xyz}^{-1}\frac{J_{\boldsymbol{\phi}^{-1}}}{\mathrm{det}(J_{\boldsymbol{\phi}^{-1}})}\ \mathbf{curl}\,\mathbf{A}_{uvw}.$		(12)

As a local comparison, $\mathbf{J}_{xyz}$ and $\mathbf{H}_{xyz}$ are evaluated along the $x$ -axis. The results depicted in Fig. 5 show an accurate agreement between all models. In the $\mathbf{A}-v$ formulation, the current density $\mathbf{J}_{xyz}$ is multiplicatively linked to the linear interpolated magnetic vector potential $\mathbf{A}_{uvw}$ (and the piecewise-constant vector field $\mathbf{grad}\,v_{uvw}$ ), see eq. (11). In contrast, in the $\mathbf{H}-\varphi$ formulation, the current density is multiplicatively linked to the $\mathbf{curl}$ operator applied to the linear interpolated magnetic field $\mathbf{H}_{uvw}$ [7]. Therefore, the current density resulting from the $\mathbf{A}-v$ formulation appears comparatively smooth, as shown in the zoom window in the upper plot of Fig. 5.

Furthermore, the comparison of the ohmic losses, representing a global quantity, shows a good match: both 2-D models output a length-related power loss of $21.9\,\mu\text{Wm}^{-1}$ , whereas the 3-D model has a total loss of $4.34\,\mu\text{W}$ . Scaling the length-related losses up to the cable’s longitudinal length results into a loss of $4.38\,\mu\text{W}$ , which deviates 0.9% from the 3-D result. We suspect that this discrepancy is mainly due to the different excitation types.

V Conclusion & outlook

In this work, we have shown how a coordinate transformation can be used for the dimensional reduction of eddy current problems modeling helicoidal symmetric power cables. In particular, we have focused on integrating the symmetry approach into the vector potential based $\mathbf{A}-v$ formulation that can be solved in 2-D (Section III). The results of this 2-D model are in excellent agreement with our previously presented 2-D model [7] and with the 3-D reference model (Section IV). This reinforces the confidence in the symmetry approach for reducing computational costs when analyzing power cables’ electromagnetic behaviour numerically. Future works include the handling of also non-ideal symmetries, see e.g., [5, 3].

Acknowledgment

Special thanks to Christophe Geuzaine (University of Liège) who provided helpful advice on implementations in GetDP.

References

[1] R. Suchantke, “Alternating Current Loss Measurement of Power Cable Conductors with Large Cross Sections Using Electrical Methods,” Ph.D. dissertation, Tech. Univ. Berlin, Berlin, Germany, 2018.
[2] CST Studio Suite. (2021), Dassault Systèmes, Accessed: Mar. 18, 2022. [Online]. Available:https://www.3ds.com/products-services/simulia/products/cst-studio-suite/.
[3] A. Marjamäki, T. Tarhasaari and P. Rasilo, ”Utilizing Helicoidal and Translational Symmetries Together in 2-D Models of Twisted Litz Wire Strand Bundles,” IEEE Transactions on Magnetics, vol. 59, no. 5, pp. 1-4, May 2023, Art no. 7400504, doi: 10.1109/TMAG.2023.3237767.
[4] A. Stenvall, F. Grilli and M. Lyly, “Current-Penetration Patterns in Twisted Superconductors in Self-Field,” IEEE Transactions on Applied Superconductivity, 2013, vol. 23, no. 3, pp. 8200105-8200105, Art no. 8200105, doi: 10.1109/TASC.2012.2228733.
[5] J. Dular, “Standard and Mixed Finite Element Formulations for Systems with Type-II Superconductors,” Ph.D. dissertation, Univ. of Liège, Liege, Belgium, 2023.
[6] A. Nicolet and F. Zolla, “Finite element analysis of helicoidal waveguides,” IET Science, Measurement & Technology, 2007, vol. 1, pp. 67–70, doi: 10.1049/iet-smt:20060042.
[7] A. Piwonski, J. Dular, R. S. Rezende and R. Schuhmann, ”2-D Eddy Current Boundary Value Problems for Power Cables With Helicoidal Symmetry,” IEEE Transactions on Magnetics, vol. 59, no. 5, pp. 1-4, May 2023, Art no. 6300204, doi: 10.1109/TMAG.2022.3231054.
[8] I.V. Lindell, “Differential forms in electromagnetics,” John Wiley & Sons, 2004.
[9] P. Dular, C. Geuzaine, F. Henrotte and W. Legros, “A general environment for the treatment of discrete problems and its application to the finite element method,” IEEE Transactions on Magnetics, vol. 34, no. 5, pp. 3395–3398, 1998.
[10] C. Geuzaine and J. Remacle, “Gmsh: a three-dimensional finite element mesh generator with built-in pre- and post-processing facilities,” International Journal for Numerical Methods in Engineering, 2009, vol. 79, no. 11, pp. 1309-1331.
[11] J. Bezanson, A. Edelman, S. Karpinski and V. Shah, “Julia: A Fresh Approach to Numerical Computing,” SIAM Review, 2017, vol. 59, no. 1, pp. 65-98, doi: 10.1137/141000671.
[12] E. Creusé, P. Dular and S. Nicaise, “About the gauge conditions arising in finite element magnetostatic problems,” Computers & Mathematics with Applications, vol. 77, no. 6, pp. 1563-1582, 2019.

Finite Element Modeling of Power Cables using Coordinate Transformations