1. Introduction
Superconducting devices are essential to speed up the world’s energy transition. Such challenges require high-quality, robust, and reliable computational models to predict and design electric equipment that applies superconductors. The literature has presented several formulations to simulate superconducting devices on a large scale. One of the most popular formulations is the H-formulation, with the magnetic field as the work variable. This formulation has several real-world applications [1], [2], [3], [4]. For example, in [2] the authors model an HTS coil in a saturated iron core using the H-formulation, the homogenization approximation is implemented, and the instantaneous power dissipation is computed. In [5], the authors simulate with the H-formulation a superconducting fault current limiter in the power electronic circuit. In this case, the H formulation provides the power dissipation profile to the thermal model, which is responsible for computing the temperature in each 2G wire layer. In [6], W Song et al. utilized the H-formulation to study the power losses in coated conductors due to nonlinear currents. The power losses due to different Total Harmonic Distortions (THD) were analysed. The result shows that harmonic components have an influence on AC losses in the bifilar stack. The H formulation solves the problem of the non-linearity in the E-J curve. However, the need to simulate the superconducting tapes leads to the mesh conditioning problem. However, the H-formulation presents limitations when dealing with nonlinear ferromagnetic materials.
Another popular technique used to simulate superconducting materials is the T-A formulation [7], [8], [9], [10], [11], [12], [2], [13], [14]. At first, the T-A formulation was developed to work with thin sheet approximation [7], reducing the simulation consumption time. In [15], Berrospe-Juarez et all. presented a multiscale and homogeneous model where, using both techniques, a superconducting coil with a large number of turns can be modelled. Moreover, some kind of application for this formulation can be found in fault current limiters [13], [10], [8], superconducting cables [16], motors [17], [18] and high field magnets [19]. In [18], F. Grilli et al. present a way to simulate the T-A formulation without the thin sheet approximation. In this paper, the authors considered both the critical state and the power law for modelling the electrical properties of a superconducting tape in an electrical machine and the AC losses in the superconducting tapes were investigated. In [20], B. Edgar et al. presented a simulation of a 32 T superconducting magnet with 20.000 turns of REBCO conductors connected in series. The T-A formulation solves the problems of the mesh in the superconducting tapes and the computation of both E-J and B-H nonlinearities. However, the T-A formulation requires Dirichlet or Neumann conditions, which sometimes are unknown. In [21] the authors present a superconductor simulation model using the magnetic vector potential as a state variable and the current density is represented as a function of the electric field. In [22] the authors present a methodology to simulate the superconductors with the A-V formulation. The methodology was used to simulate a benchmark problem and the magnetization of two cylindrical superconductors, the results show good agreement with experimental data.
In this work, we propose a new approach to simulate superconducting devices based on the A-formulation coupled with the J-formulation. Both formulations are solved using the Finite Element Method in the commercial software COMSOL Multiphysics. The A-formulation is active in all domains, except the superconducting domain, where the J-formulation is active and computes the current density directly. The proposed approach has two work variables: current density, J, and the magnetic vector potential, A. The way to couple both formulations is furtherly discussed. This approach is investigated in different scenarios for distinct shape function orders. J-A and T-A formulations were compared, and their results presented minor differences. The paper is divided as follows: Section 2 presents the J-A approach, the coupling method between formulations, and the studied cases. Section 3 discusses the results of all studied cases, where the impact of different shape function orders are analysed. The agreement of the results with the T-A formulation is discussed. Finally, Section 4 presents the conclusions of the work.
2. Methodology
This section discusses the proposed approach and the studied cases used to validate the model. Five cases are explored. The first and second cases are about the transport current and inducted current losses, respectively. The third case is a simulation of the benchmark presented in [34], where the T-A formulation is used for comparison. The voltage drop between the tape and the current density distribution are discussed. In the fourth case, the induced current in a stack of superconducting tapes is computed. The distribution of the current density obtained by the proposed model is compared with the one simulated in [35]. At last, the model of current sharing in the racetrack coil is analysed and compared with the T-A formulation implemented in [9].
2.1. J-A approach
The J-A approach has the current density and the magnetic vector potential as its working variables. Using the same principle of the T-A formulation [7], the A-formulation is active for all domains and the J-formulation is active only in the superconducting regions. The J-formulation acts as a source in the A-formulation, and the reciprocal is true as well. The A-formulation is stated as Eq. (1),
$\nabla \times\left(\frac{1}{\mu} \nabla \times \mathbf{A}\right)=\mathbf{J}$
where A,μ,J are the magnetic vector potential, the magnetic permeability, the current density. On the right side of Eq. (1), an external current from the superconducting layer is imposed with J, so that the J-formulation, which computes the current density in the superconducting domain, acts as a source in the A-formulation. With Ohm’s law, we compute the current density in the superconducting domain according to Eq. (2),
$\left(\frac{E_{c}}{J_{c}}\left|\frac{\mathbf{J}}{J_{c}}\right|^{n-1}\right) \mathbf{J}=-\frac{\partial \mathbf{A}}{\partial t}-\nabla V.$
where, $E_{c}, J_{c}(\mathbf{B}), \mathbf{B}$ and V are the critical electric field, critical current density, magnetic flux density and the voltage, respectively. In this case, ∇V is a source term that can be used to couple the superconducting device in an external electrical circuit or to impose directly in the FEM a gradient of the electric potential.
Observing the right side of (2), it is clear that the A-formulation acts as a source in the J-formulation, closing the loop of the hybrid formulation. Moreover, observing the left side of Eq. (2), it is worth mentioning that the equation is not a partial differential equation for space variables. It depends solely on a time derivative. As a result, the computational time consumption to solve the problem is reduced. If compared to the T-A formulation, less computational time is expected, as the current density is directly computed (it is the work variable). Analyzing the right side of the Eq. (2), two sources are possible, the magnetic vector potential and the electric scalar potential. This last variable may be used to enforce an external voltage. In the case that no external voltage is applied, the gradient of electrical scalar potential is zero. Therefore, Eq. (2) becomes,
$\left(\frac{E_{c}}{J_{c}}\left|\frac{\mathbf{J}}{J_{c}}\right|^{n-1}\right) \mathbf{J}=-\frac{\partial \mathbf{A}}{\partial t}.$
In the 2D case, considering an infinite Cartesian plane with a z-direction as its normal vector and without the external voltage, Eqs. (1), (3) become
$\nabla \times \frac{1}{\mu} \nabla \times\left[\begin{array}{c} 0 \\ 0 \\ A_{z} \end{array}\right]=J_{z}-\sigma \frac{\partial A_{z}}{\partial t}$
$\left(\frac{E_{c}}{J_{c}}\left|\frac{J_{z}}{J_{c}}\right|^{n-1}\right) J_{z}=-\frac{\partial A_{z}}{\partial t}$
Integral restrictions are used to enforce total currents on the superconducting domains, according to the Eq. (6),
$\iint J_{z} d S-I_{a p p}=0$
where $I_{app}$ is the applied current. In the case of a superconducting coated conductor, with the thin-film approximation in the J-A formulation, the enforcement of current is straightforward: the current computed in the superconducting line domains in J-formulation acts as an external surface current density in A-formulation, according to 7,
$J_{z}=\mathbf{J}_{\mathbf{H T S}} d_{H T S},$
where $d_{H T S}$ is the thickness of the superconducting layer and the $\mathbf{J}_{\mathrm{HTS}}$ is the current density of the superconductor domain. Considering the thin-sheet approximation, a boundary condition is used to enforce the current density in the A-formulation,
$\mathbf{n} \times\left(\mathbf{H}_{1}-\mathbf{H}_{\mathbf{2}}\right)=\mathbf{J}_{\mathbf{H T S}} d_{H T S}.$
where the vectors H1 and H2 are the magnetic field upward and downward of the superconducting tape and n is the normal unit vector. It is important to emphasize that in all cases studied in this article the critical current density is considered a constant because all reference models used this approximation. COMSOL Multiphysics software was used to model all presented cases in this work. In all simulated cases, the superconducting devices are not coupled to an external electrical circuit. Due to that, the ∇V is equal to zero in all of implemented cases.
Another advantage of the J-A formulation is the possibility of using the same order of elements in both formulations. On the other hand, in the T-A formulation, the magnetic vector potential must have a high order element than the current vector potential. Due to the J-A formulation’s ability to use the same order in both formulations, J and A, the total computation time can be decreased. This fact is an attractive asset for 3D simulations.
2.2. Cases 1 and 2: single coated conductor attached to a current source and single coated conductor with induced current
The first two cases investigate a single tape with with a 12 mm width in a 2D approximation in an x-y cut plane. The thin-sheet approximation is used for the coated conductor.
For the first case, the coated conductor is immersed in an air domain and a transport current is applied on it, as presented in Fig. 1.
Fig. 1. Illustration of case 1 and 2: Simulation of a single superconducting tape with transport or induced current. |
A sine waveform with its maximum value equal to 80% of the critical current: $i(t)=400 \cdot \sin (2 \cdot \pi \cdot 50 \cdot t)$t). Eq. (6) enforces the applied current in the superconducting domain in the J-A formulation, and appropriate Dirichlet boundary conditions are used to apply current in the T-A formulation [14]. For the second case, a coated conductor immersed in an air domain has a current induced on it by a homogenised external alternating magnetic field in the x direction, also with a sine waveform with $\mathbf{B}(t)=0.12 \cdot \sin (2 \cdot \pi \cdot 50 \cdot t)$ T as Fig. 1.
The J-A approach results are compared to T-A formulation results. The same configurations and parameters were used for both formulations, and the simulations were run on the same computer. Table 1 summarizes the parameters and configurations considered for both cases.
Table 1. Parameters and configurations for cases 1 and 2. |
| Critical current ($I_c$) | 500 | A |
| Transition index (n) | 30 | - |
| Critical Electrical Field ($E_c$) | 10-4 | V/m |
| Air conductivity ($σ_{air}$) | 0 | S/m |
| Air permeability ($μ_0$) | 4π10-7 | H/m |
2.3. Case 3: HTS dynamo benchmark
The benchmark case is about a high-temperature superconducting (HTS) dynamo. The geometry of the model is composed by a permanent magnet (PM) that rotates in anti-clockwise and passes near a single HTS-coated conductor with its terminals open. The problem is modelled in a 2D Cartesian plan assumed infinitely long. The authors simulate the problem using several formulations, and here we compare the J-A results with the T-A formulation explained and presented in [34]. Using this model, the benchmark can be evaluated with moving mesh and the E-J nonlinearity. In both formulations (J-A and T-A) the thin sheet approximation is considered. Due to the open terminals, the total current in the tape must be zero. In J-A formulation, to impose a total current equal to zero, Eq. (6) becomes
$\int J_{z} d l=0$
Table 2. HTS dynamo parameters. |
| Permanent Magnet | Width | 6 mm |
| Height | 12 mm | |
| Active length | 12.7 mm | |
| Remanent Flux | 1.25 T | |
| HTS tape | Width | 6 mm |
| thickness | 1 μm | |
| Critical current | 282 A | |
| Temperature | 77 K | |
| Transition index | 20 |
Fig. 2. Illustration of case 3: with the permanent magnet, the coated conductor and the boundary between moving and stationary meshes. |
2.4. Case 4: current induction in stack of coated conductors
This section explains the modelling of case 4: current induction on stacked coated conductors. The model was obtained from the HTS modelling Workgroup website and was implemented in [35]. It assumes 2D infinitely long xy-plan, the thin-sheet approximation is used, and the superconducting tape is oriented in the x-direction. A magnetic flux density in y-direction is enforced by a Dirichlet boundary condition. The stack has 15 tapes in open-circuit and as a result, the total current must be zero in all tapes. Eq. (9) models this condition in the J-A formulation. In the T-A formulation, to keep the same scenario, a Dirichlet boundary condition is configured at the edges of the superconducting tapes. The parameters of the model can be found in [35] and Fig. 3 depicts the studied case. For inducing the current, a homogenized magnetic field in y direction is applied as indicating in the Fig. 3.
Fig. 3. Illustration of case 4: a stack of coated conductors under the effect of an external magnetic field in the y-direction. |
2.5. Case 5: current sharing between tapes in a racetrack coil
This section explains a case study of current sharing in a racetrack coil [9], where the current distribution between two parallel superconducting tapes wound in a racetrack shape is studied. Its cross-section is simulated, as illustrated in [9]. The problem is assumed in 2D xy plane and Fig. 4 depicts the simulated case. It is a two-turn coil made of a 2 parallel conductors. Thus, tapes 1, 3, 6 and 8 are connected in series and compounding the first-turn coil. For the second-turn coil, the tapes 2, 4, 5 and 7 are in series connected. In this case study, a sine waveform current is enforced in the coil, where the current sharing between parallel tapes are evaluated.
Fig. 4. (a) Illustration of case 5: current sharing between tapes wound as a racetrack coil. Tapes 1, 3, 6 and 8 are connected in series, as well as 2, 4, 5 and 7. (b) Connections of the tapes for the studied case. |
A circuit coupling method is used to represent the electrical connections between tape segments [8], [13], [12], [10], [11]. On the other hand, in J-A formulation, this connection is implemented based on the same principle presented in the previous cases: integral restrictions, the Eqs. (10), (11), (12), (13), (14), (15), (16) model the series connections and the direction of the current. For example, tapes 1 and 8 are series connected and have opposite current directions. To model this behavior, we integrate the current density in the superconducting tapes domains and enforce a restriction that their sum is zero. To impose the current in the racetrack coil and consider the parallel tapes, we consider the integration of the current density in the superconducting domains of tapes 1 and 2 and enforce with a restriction that their sum is equal to the applied current.
$\int_{\text {tape } 1} J_{z} d l+\int_{\text {tape } 8} J_{z} d l=0$
$\int_{\text {tape } 2} J_{z} d l+\int_{\text {tape } 7} J_{z} d l=0 \text {. }$
$\int_{\text {tape } 1} J_{z} d l-\int_{\text {tape } 3} J_{z} d l=0$
$\int_{\text {tape } 2} J_{z} d l-\int_{\text {tape } 4} J_{z} d l=0$
$\int_{\text {tape } 3} J_{z} d l+\int_{\text {tape } 6} J_{z} d l=0$
$\int_{\text {tape } 4} J_{z} d l+\int_{\text {tape } 5} J_{z} d l=0$
$\int_{\text {tape } 1} J_{z} d l+\int_{\text {tape } 2} J_{z} d l-I_{a p p}=0$
3. Results
In this section, the results from the five studied cases are presented and analysed. First, we discuss the transport and induction current problems. For cases 1 and 2, the impact of element order of both J-A and T-A formulations are analysed in two scenarios, the normalised current density and the power losses. In case 3, the HTS dynamo, the voltage and the normalised current density are observed. The total computation time for the models (J-A and T-A) are compared. For case 4, induction in stacked coated conductors, magnetic flux density, and induced current results of J-A and T-A formulations are compared. The total computation time is also analysed for this case. The results of case 5, about current sharing between tapes, are the total current in each tape segment over time and the total computational time.
3.1. Results and discussions of cases 1 and 2: single coated conductor
Fig. 5 shows the normalised current density along the tape width. This problem was simulated with the elements of the finite discretisation using Lagrangian shape functions with different orders to assess the impact of this parameter in the J-A and T-A formulations. It was possible to observe that in contrast to the T-A formulation, the J-A formulation runs well using elements with the same shape functions order in J and A formulations. Moreover, this figure shows the results of the J-A formulation with second-order shape functions (in red), first-order shape functions for J with A modelled with second-order shape functions (in green), and first-order shape functions (in blue). The use of quadratic order to model A generates a smooth current density along the line width. These results show that J-A works well with all combinations presented. However, in our experience so far, using a shape function order for J higher than the order chosen for A leads to divergence problems during the simulation.
Fig. 5. Comparison of different shape function order of J-A and T-A formulation for normalised current density at 1 ms as function of the width of the tape for (a) applied current case and (b) induced current case. |
For the T-A formulation, it is necessary to have A formulation with the element order higher and T formulation. In 2D systems with a thin-sheet approximation, the T-A formulation is commonly used with second-order elements applied to A and first-order elements applied to T. This means that the resulting current density computed by the T part is a constant within each element. To check this effect, consider the example shown in Fig. 6. There is a coated conductor approximated by a line segment, defined along the x coordinate. The mesh that represents this line segment is composed of N nodes and M elements. If first-order elements are chosen for the T part, inside element 1, for example, T has the format
$T(x)=a x+b$
with $a,b$ as real constants. The Eq. (17) can be written in terms of the current potential of each node. For element 1, we consider $T_1$ as the current potential of node 1 and $T_2$ as the current potential of node 2. It can be proved that, for this case, the potential inside element 1 is given by Eq.(17)
$T(x)=\frac{T_{2}-T_{1}}{x_{2}-x_{1}} x+\frac{T_{1} x_{2}-T_{2} x_{1}}{x_{2}-x_{1}},$
with $x_1$ and $x_2$ as the coordinates of nodes 1 and 2, respectively. As the current density is computed by
$\mathbf{J}=\nabla \times \mathbf{T}, \quad J_{z}=\frac{\partial T n_{y}}{\partial x}$
it becomes
$J_{z}=\frac{T_{2}-T_{1}}{x_{2}-x_{1}}$
Fig. 6. Generic finite element division in the superconducting tape approximated by a line. |
So, the current density is directly proportional to the potential difference and inversely proportional to the distance between nodes. Interpolation is necessary to show a final continuous function for $J_z$. For the J-A approach, however, the approximation for the current density is never a constant within each element. Because the current density is directly assessed and in the combination of order for J and A, the J-formulation can have the same order as the A-formulation. It implies that the current density may be modelled by shape functions of the same order desired of the magnetic vector potential. For example, the first-order equation for J for the first element of the system in Fig. 6 is
$J(x)=\frac{J_{2}-J_{1}}{x_{2}-x_{1}} x+\frac{J_{1} x_{2}-J_{2} x_{1}}{x_{2}-x_{1}}$
For a second-order approximation, the element would have to be composed by three nodes, not two. Defining, for the first element, the overall equation for the second order case is given by the Eq. (17)
$J(x)=c x^{2}+d x+e$
With
$c=\frac{1}{D}\left[x_{2}-x_{3}, x_{3}-x_{1}, x_{1}-x_{2}\right]\left[\begin{array}{l} J_{1} \\ J_{2} \\ J_{3} \end{array}\right] $
$d=\frac{1}{D}\left[x_{3}^{2}-x_{2}^{2}, x_{1}^{2}-x_{3}^{2},-x_{1}^{2}+x_{2}^{2}\right]\left[\begin{array}{l} J_{1} \\ J_{2} \\ J_{3} \end{array}\right]$
And
$e=\frac{1}{D}\left[\begin{array}{c} x_{2}^{2} x_{3}-x_{2} x_{3}^{2} \\ -x_{1}^{2} x_{3}+x_{1} x_{3}^{2} \\ x_{1}^{2} x_{2}-x_{1} x_{2}^{2} \end{array}\right]^{T}\left[\begin{array}{l} J_{1} \\ J_{2} \\ J_{3} \end{array}\right]$
where D is the determinant of the linear system solved to find the values of constants c, d and e. This determinant is given by Eq. (17)
$\begin{array}{r} D=x_{1}^{2} x_{2}-x_{1} x_{2}^{2}-x_{1}^{2} x_{3}+x_{2}^{2} x_{3} \\ +x_{1} x_{3}^{2}-x_{2} x_{3} \end{array}$
Fig. 7 presents the power losses in the superconducting tape due to the transport current, which is computed according to Eq. (27) for both J-A and T-A formulations,
$P_{\text {losses }}=\int\left(\frac{E_{c}}{J_{c}}\left|\frac{\mathbf{J}}{J_{c}}\right|^{n-1}\right) J_{H T S}^{2} d l$
Fig. 7. Power losses as function of time for different shape functions. |
As it is possible to see in the highlight of the second peak of the curve, the J-A formulation reduces the harmonic order peaks in the power losses curves. Observing Eq. (27), power losses are computed with the square of the current density, therefore an error in the current density is amplified in the power losses’ calculation, and it generates a large peak of harmonics in power losses curve. As the J-A formulation computes directly the current density, the error in this variable is reduced and, consequently, the peak of the harmonic is reduced, as can be observed in Fig. 7.
3.2. Results and discussion of cases 3, 4 and 5: HTS dynamo, stack and racetrack coil
First, we discuss the results of case 3, the HTS dynamo. Fig. 8 presents the total voltage in the tape computed by both formulations. The agreement between the curves are evident, however spurious values are observed in the J-A formulation results at around 0.118 s. This condition is avoided if more elements are added to the superconducting tape mesh or if a smaller time-step is configured.
Fig. 8. Voltage across the superconducting tape for the dynamo studied cases. |
Fig. 9 shows the normalised current density at 0.347 s as function of the width (12 mm) of the superconducting tape for both J-A and T-A formulations. In this case, as the net current must be zero, due to the superconducting tape being an open circuit, the current density assumes positive and negatives values (currents in and out of the plane) indicating the circulation of the current within the superconducting tape. In the Fig. 9, the J-A formulation is configured with quadratic elements in both J and A formulation, while a linear T-formulation and quadratic A-formulation are used in the T-A simulation.
Fig. 9. Normalised current density at 0.347 s as function of the width of the superconducting tape. |
The instantaneous power losses for both formulations are presented in Fig. 10. The agreement between the T-A and J-A results are evident, except for one time-step close to 0.1 s. This non-agreement at this point can be fixed by reducing the simulation time-step or improving the mesh at the superconducting tape.
Fig. 10. Power losses in superconducting tape for the benchmark case. |
Table 3 presents the total computation time for three combinations of shape function orders of the J-A approach and for the T-A formulation. As expected, a faster simulation is obtained by using the linear function order for both J and A formulation, having a total simulation time equal to one minute and 28 s. The second faster simulation appears when we use the linear order in J and quadratic order in A, having a simulation time equal to one minute and thirty-three seconds. In both previous combinations, the J-A formulation is faster than the T-A formulation that was simulated using the linear order for the T and quadratic order function for A formulation. In this case, the simulation time obtained was equal to six minutes. The computation time observed at J-A formulation with both J and A with quadratic order elements is equal to two minutes and nine seconds min. In all cases, the J-A formulation is faster than the T-A formulation because: i) the current density is directly computed by the J-A formulation, whereas the T-A formulation computes the current density based on the current vector potential. ii) the J-A formulation allows the possibility of running both J and A-formulations with first order elements and the computation time is sped up, consequently. On the contrary, the T-A formulation must have the T formulation element one order lower than the A formulation. This fact causes longer computation time when compared to the J-A formulation.
Table 3. Total computation time, J-A and T-A formulation in the third case. Number of elements 1888. |
| Formulation | Simulation time (min) | Element Orders |
|---|---|---|
| J-A | 1:28 | J-1 A-1 |
| J-A | 1:33 | J-1 A-2 |
| J-A | 2:09 | J-2 A-2 |
| T-A | 6:00 | T-1 A-2 |
Next, we move to the results and analysis of case 4, the stack of coated conductors with current induced by an alternating magnetic field. Fig. 11 presents the comparison of magnetic flux density between the J-A and T-A formulations at the same time instant. In both, the maximum value of the magnetic field is 0.18 T and within the superconducting region the magnetic field is close to zero. To avoid the penetration of the magnetic field in the central part of the tape, a current density is induced in a close-loop. Fig. 12 shows the comparison of normalised current density between both formulations. Each colors represent a tape simulated. In this case, both J and T formulation use linear shape functions and the A-formulation uses quadratic shape functions. The results agree, and the error is less than 1%, according to Eq. (28).
$e r=\frac{J_{T-A}^{e}-J_{J-A}^{e}}{J_{T-A}^{e}}$
where, $J_{T-A}^{e}$ and $J_{J-A}^{e}$ are the current density computed in each element of the superconducting tape for the T-A and J-A formulations.
Fig. 11. Magnetic Flux Density comparison between the J-A and T-A formulations. |
Fig. 12. The normalised current density in the tapes 1 up to 7 for both J-A and T-A formulation at 3.2 ms. |
Table 4 presents the companions of the simulation time between the J-A and T-A formulation.
Table 4. Total computation time, J-A and T-A formulations in the fourth case. Number of elements: 20492. |
| Formulation | Simulation time (min) | Element Orders |
|---|---|---|
| J-A | 0:44 | J-1 A-1 |
| J-A | 1:32 | J-1 A-2 |
| J-A | 4:52 | J-2 A-2 |
| T-A | 4:39 | T-1 A-2 |
Moving to the final tested case, the current sharing between tapes in a racetrack coil. Fig. 13 presents the comparison between the J-A and T-A formulations. In this comparison, an error less than 1% was observed when the cases are compared. It shows that the J-A formulation can also be applied to coupling tapes at a 2D problem. However, in the J-A formulation, the implementation of coupling between tapes or tape turns is straightforward and can be done using an integral restriction.
Table 5 presents the comparison results of total computation time of J-A and T-A formulations for case 5. Again, simulations are fastest when J and A are in linear order. However, comparing the J-A and T-A formulation with the same order of elements in J, T and A, in this case, the T-A formulation is faster than J-A for 24.25% (that represents a difference between the formulations of 49 s). The authors consider that the reason for that is the coupling circuit methodology [8], [9] used to simulate the connections of the superconducting tapes. In somehow, the coupling circuit methodology can generate a faster convergence in the T-A formulation. However, it is worth mentioning, that this hypothesis should be more investigated in coming future.
Table 5. Total computation time, J-A and T-A formulations in the fifth case. Number of elements 10454. |
| Formulation | Simulation time (min) | Element Orders |
|---|---|---|
| J-A | 1:01 | J-1 A-1 |
| J-A | 3:22 | J-1 A-2 |
| J-A | 2:56 | J-2 A-2 |
| T-A | 2:33 | T-1 A-2 |
4. Conclusion
This work proposes a new approach to simulate superconducting devices in finite element method. The J-A formulation proposed here is validated for several cases studies, and was compared to T-A formulation. The new formulation can run with the J and A formulation using a linear order function approximation. Such possibility speeds up the simulation time and could be useful in 3D problems, which are yet to be studied. Other function approximation were analysed as well. In most of the cases, when using the J with a linear order and the A with a quadratic order, the J-A formulation was faster than the T-A formulation using the same respective order for T and A.
The benchmark dynamo example from HTS Modelling Workgroup was simulated in J-A formulation and the comparison with the T-A formulation reveals that results agree in for the voltage across the superconducting tape and for the current density, and the results agree for the instantaneous losses. The results on this benchmark confirm that the J-A formulation is a suitable approach to simulate superconducting devices in a moving mesh problem.
As for the stacked tapes case, the magnetic field, the current density and the simulation time were compared. The same results were observed, where using J and A with a linear function approximation gives the fastest solution and using J and A with linear and quadratic order, respectively, is faster than the T-A formulation with the T part using the linear order and A part using the quadratic order. The magnetic field and the current density agree for all tapes. In this way, we can conclude that J-A formulation could be an advantageous approach to simulate scenarios with high number of superconducting tapes.
In the last case studied, the current sharing between two parallel tapes was evaluated for both J-A and T-A formulations. The current results of all tapes agree with an error less than 1% for all the simulation time steps. Such result shows that the J-A formulation is also an appropriate approach to simulate coupling between superconducting tapes in 2D cases. The comparison between simulation time have shown that the fastest one is when J and A have the shape function order equal to one. Nevertheless, when the T and A formulation have a linear and quadratic function order, it shows the faster results than the J-A formulation.
Finally, all results in this study leads to believe that the J-A formulation can be a useful approach to simulate superconducting devices due to its simple implementation, it requires less computational effort, given by its lesser simulation time when using the linear order for both formulations (J and A), and it shows good agreement with the benchmark cases. In future studies, more cases should be investigated, including 3D modelling, but the good agreement with the actual benchmarks are already a promising indicative.
Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
Appendix A. J-A approach implementation in COMSOL Multiphysics Software
In order to implement the J-A formulation in the COMSOL Multiphysics Software, two packages should be used, the magnetic field package, named mf, and the coefficient form PDE package, a generic Partial Differential Equation (PDE) solver, that has the form equal to Eq. (A.1),
$\begin{array}{r} e_{a} \frac{\partial^{2} u}{\partial t^{2}}+d_{a} \frac{\partial u}{\partial t}+\nabla \cdot(-c \nabla u-\alpha u+\gamma) \\ +\beta \cdot \nabla u+a u=f \end{array}$
Where $e_{a,}, d_{a}, c, \alpha, \gamma, \beta$ and a are constants. To implement the J-formulation, on the left side, all constants are equal to zero, except for constant a, which is equal to the resistivity of the superconductor (a=ρ). On the right side, f, the forcing term, is equal to the variation of the magnetic vector potential with respect to the time ($f=-\partial A / \partial t$), in COMSOL software this equation is implemented by f =-d(A,t). In this case, the work variable (u) is the current density (J) and Eq. (A.1) becomes,
$a u=f$
where, in this case, no space derivative term is used. This fact speeds up the computation time. For imposing current in the superconducting domain, an integral constraint is used, according to Eq. (6). In COMSOL Multiphysics Software, the aforementioned constraint is implemented using a Global constraint, according to Eq. (6)
$\operatorname{intop}_{H T S}\left(J_{H T S}\right)-I_{a p p}=0$
where intopHTS represents the integral on the superconducting domain, which can be an area or a line (thin sheet approximation). In the case of induced current, $I_{app}$ is equal to zero. The current density computed in the J-formulation acts as a source in the A-formulation based on an external current density. Considering the thin sheet approximation, an external current is enforced in the mf using the surface current density operator, based on the following equation,
$\mathbf{n} \times\left(\mathbf{H}_{\mathbf{1}}-\mathbf{H}_{\mathbf{2}}\right)=J_{H T S} d_{H T S}$
where the n,H1 and H2 are normal unit vector, the external and the internal magnetic fields, respectively. This concludes the aspects related to the implementation of J-A in COMSOL.
Appendix B. Supplementary material
Supplementary data associated with this article can be found, in the online version, at https://doi.org/10.1016/j.supcon.2023.100049.
