- Email: [email protected]
Development of a three-dimensional finite-element model for high-temperature superconductors based on the H-formulation
Cryogenics 53 (2013) 142–147
Contents lists available at SciVerse ScienceDirect
Cryogenics journal homepage: www.elsevier.com/locate/cryogenics
Development of a three-dimensional finite-element model for high-temperature superconductors based on the H-formulation Francesco Grilli a,⇑, Roberto Brambilla b, Frédéric Sirois c, Antti Stenvall d, Steeve Memiaghe c a
Karlsruhe Institute of Technology, Karlsruhe, Germany Ricerca Sistema Energetico S.p.A., Milano, Italy c Ecole Polytechnique de Montréal, Montréal, Canada d Tampere University of Technology, Tampere, Finland b
a r t i c l e
i n f o
Article history: Available online 8 May 2012 Keywords: Finite-element models High-temperature superconductors Ac losses 3-D simulations
a b s t r a c t Finite-element models are a powerful and widely used tool for evaluating the ac losses of HTS tapes and wires as well as of assemblies such as cables and coils. The H-formulation, which uses the magnetic field components as state variables, has proved to be an efficient implementation to solve 2-D problems, involving infinitely long or axially-symmetric geometries; an excellent agreement with experimental data has been found in many cases. However, the simulation of certain applications requires a full 3-D model. In this paper we report on the development of a 3-D model based on the H-formulation. We describe the implementation of Maxwell equations, the imposition of current constraints and we discuss the issues related to meshing 3-D volumes. The model is validated by comparing the results with those obtained with 2-D models in cases that can be investigated in 2-D; then, it is used to simulate cases that can be handled only in 3-D. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction In recent years finite-element simulations have become a powerful tool to investigate in detail the current density and magnetic field distributions inside high-temperature superconductors (HTS) and to compute their ac losses. In particular, the H-formulation [1,2] has increased in popularity both in the academic and the industrial world and it has been successfully used to study the ac loss performance not only of individual wires, but also of more complex systems like cables and coils [3–5]. The model has been adapted to take into account the presence of ferromagnetic materials as well [6,7]. Until now, the model has been only used to simulate 2-D problems, i.e. geometries representing conductors either infinitely long or characterized by rotational symmetries (coils), and excellent agreement with experimental data has been observed, both for individual and interacting tapes [8]. However, many problems require the simulation of 3-D effects, such as coupling currents, local non-uniformities of the physical or geometrical properties, and twisted or transposed geometries. An extension of the model to 3-D is therefore necessary. In this paper we present our first results of such extension. The paper is organized as follows. At the beginning we give a brief description of the model implementation. Then we describe and discuss the results of our simulations: first we validate the model ⇑ Corresponding author. E-mail address: [email protected] (F. Grilli). 0011-2275/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.cryogenics.2012.03.007
by comparing its results with those that are obtained in cases that can be simulated in 2-D; then we provide results for cases that can be analyzed only with the new model, like the current sharing between superconductors caused by the presence of localized defects and a superconductor with varying cross-section along its length. Finally, we draw the conclusions, with focus on the open issues and the main challenges of 3-D numerical modeling of superconductors. 2. Model description The basic equation of the model, which is implemented in the Partial Differential Equation (PDE) module of the finite-element software package Comsol Multiphysics [9], is Faraday’s equation
l
@H þ r E ¼ 0; @t
ð1Þ
where the components of the magnetic field H = (Hx, Hy, Hz) are the state variables. For non-magnetic materials, the relation B = l0H is assumed. The superconductor material is modeled by a non-linear resistivity
q¼
n1 Ec J ; J c J c
ð2Þ
where Ec is the critical electric field (usually equal to 1 lV/cm), Jc is the critical current density, J is the absolute value of the current density, and n is the power-law exponent. This latter typically
F. Grilli et al. / Cryogenics 53 (2013) 142–147
assumes values between 15 and 50 for high-temperature superconductors. For the simulations contained in this paper we used n = 35. The conductors are usually surrounded by vacuum, which is modeled as a resistive material with very high resistivity (1 X m). For all the simulations of this paper, the frequency of the magnetic field and current sources is 50 Hz. It is very important to remark that the use of an isotropic power-law resistivity such as (2) is not a very realistic assumption, both from the physics and modeling point of view. 2-D models have proven to be in excellent agreement with measurement when the 2-D symmetry is also realized in the measurements – see for example [8,10]. In those cases, the magnetic field lines are always perpendicular to the current, and pinning and flux flow in the microscopic scale can be used to describe what happens inside the superconductor [11]. However, in real 3-D cases the current is not perpendicular to applied field and a better determination of the microscopic model is needed. Additionally, the microscopic description needs to give input for the resistivity tensor used in the eddy current type 3-D models. Furthermore, the resistivity tensor is not necessarily only a function of the anisotropic properties in the conductor, like different Jc along c-axis or in ab-plane of coated conductors, but may depend also on the direction of the magnetic field. The so-called double critical state model was proposed by Clem and Perez-Gonzalez [12]. They utilized different Jc values in perpendicular direction and along the flux lines. There is also a lot of recent interest in the topic [13–15], but these works do not provide the solution for resistivity in eddy current type solvers yet, nor methods to easily measure just a few parameters to characterize the behavior of Jc in all possible field orientations with respect to all possible current directions in a sample. This being said, since this work is the first extension of the Hformulation from 2-D to 3-D, the simplifying assumptions of (2) are used in this paper, also because they allow a direct validation of the new model against the 2-D one and an easy verification that the electromagnetic quantities are computed correctly. Similar to the 2-D model [2], an external magnetic field can be easily applied by appropriately setting the Dirichlet boundary conditions on the domain’s boundary. A transport current I(t) can be imposed by means of a current constraint on the desired face (domain boundary), as follows. First one defines the variable describing the current flowing through a face S as
IS ¼
Z
ðJ x nx þ J y ny þ J z nz ÞdS;
ð3Þ
S
where (nx, ny, nz) is the unit vector defining the normal direction of the face S and (Jx, Jy, Jz) are the three coordinate components of the current density; then, one imposes I(t) IS = 0 as an algebraic constraint. This latter forces the current flowing through the face S to be equal to I(t) at each time step of the simulation. In this way, when different conductors are present, multiple independent current constraints can be imposed. The choice of the mesh is very important in order to have sufficiently precise results, yet keeping the total number of degrees of freedom of the problem at a reasonable level. In this regard, the use of the so-called swept mesh can be very useful. A swept mesh is created by starting from a mesh of the 2-D cross-section of the conductor region; the 2-D mesh (which can be quite fine) is replicated along the length of the conductor. An example is represented in Fig. 1, where one can also notice that the mesh density along the conductor’s length (x-axis) is not uniform, but is refined in the proximity of the region where important changes of the current density and magnetic field distributions occur – see Section 3.3 for more details. This type of mesh is usually more convenient than the standard free tetrahedral mesh, even if its use might not be possible or straightforward with certain geometries. For particularly complex geometries, the use of specifically dedicated meshing
143
Fig. 1. Example of swept mesh. The mesh is first generated on a 2-D cross-section (yz plane) and then extruded and replicated along the third dimension (x direction). In this example the mesh is refined near mid-length along the x-axis, because it is used to study the case of a localized defect there – see Section 3.3 for details.
programs can be helpful. As far as the discretization of the elements is concerned, first-order edge elements, which guarantee the continuity of the tangential field component across adjacent elements, are used. The problems shown in this paper contain approximately between 15,000 and 200,000 mesh elements. As a consequence, the computation time is very variable. Even for the same problem, the computation time can change dramatically depending on the value of the applied field or current. Sometimes, in order to obtain convergence, a very small time step has to be set for the solver, which inevitably results in longer computation times. As an approximate figure, we can say that the simulation of one cycle at 50 Hz takes between a few hours and 1–2 days on a 8-core workstation equipped with 8 GB of RAM. 3. Results In this section we present the results for different 3-D geometries. The first two cases serve also as validation of the model, since they can be solved by 2-D models (at least to a certain extent). For the case of helicoidal superconductors (Section 3.2), a comparison with a 2-D model utilizing the combination of translation and rotation symmetry is performed. This eddy current model is formed from the magnetostatic model of helicoidal conductors [16] coupled to Faraday’s law, similarly to what has been done in [17].
Fig. 2. Schematic illustration of the coupling currents induced by a varying magnetic field in two thin superconductors separated by a normal metal barrier. The field is applied perpendicular to the flat face of the superconductors (z-axis). In the illustrated case part of the current induced in the xy-plane flows from one superconductor to the other through the barrier (partial coupling).
144
F. Grilli et al. / Cryogenics 53 (2013) 142–147
(a)
(b)
Fig. 3. Magnetic flux (a) and current density (b) profiles across the width (y-axis in Fig. 2) of two thin superconductors separated by a resistive barrier for different values of the barrier resistivity. The profiles computed with a 2-D model for the fully uncoupled and fully coupled situations are also shown. The profiles are taken at the peak of a field of 50 mT of amplitude at 50 Hz, applied perpendicular to the superconducting strips.
Fig. 4. View of the simulated helicoidal conductors (one pitch length). The plane at half length represents the 2-D cross-section, where the magnetic flux density distributions of Fig. 5 are taken.
The other simulated problems are examples of applications of the 3-D model in cases of practical interest, such as the current sharing due to local non-uniformities of the superconductor and the variation of the superconductor’s cross-section. 3.1. Coupling between superconductors in external ac magnetic field We first considered the case of two superconducting thin strips separated by a resistive barrier under the action of a perpendicular ac field. The magnetic field induces current loops. According to the considered conditions and geometry (variation rate of the magnetic field, length of the conductors, resistance of the barrier), the induced currents can flow separately in each superconducting strip (uncoupled situation) or from one strip to the other (coupled situation) [18]. The cases of full uncoupling and full coupling can be simulated in 2-D by appropriately setting the current constraints,
and they can be used to validate the 3-D model. The case of partial coupling (when only a fraction of the induced current crosses the resistive barrier) is schematically represented in Fig. 2 and can be simulated only with a 3-D model. We simulated two superconducting strips, 2 mm wide, 5 lm thick, each with Ic = 200 A, separated by a normal metal barrier 0.5 mm wide. We applied a field of 50 mT and observed the magnetic field and current density profiles along the y-axis of Fig. 2. The resistivity of the barrier was varied in order to consider different degrees of coupling between the strips. Fig. 3 shows the magnetic field and current density profiles at the peak of the applied field for different values of the resistivity of the metal barrier, along with the results obtained with 2-D simulations by imposing perfect uncoupling and coupling between the superconductors, respectively. It can be seen that for a sufficiently high value of the resistivity the 3-D results agree very well with the 2-D uncoupled case; and for a sufficiently low value of the resistivity, the 3-D results agree very well with the 2-D coupled case. This shows that the 3-D model works well and the electromagnetic quantities are computed correctly. For an intermediate value of the resistivity, one has a situation of partial coupling: in both superconductors the induced current density has both positive and negative components, but in one superconductor it is predominantly positive, in the other one it is predominantly negative. This means that part of (but not all) the induced current flows from one superconductor to the other – see also Fig. 2. This is a typical case that cannot be handled with 2-D simulations. A similar change of the degree of coupling can be obtained by varying the frequency of the applied magnetic field or the conductors’ length.
3.2. Helicoidally-wound superconductors As a second example, we simulated three helicoidally wound round superconductors carrying transport current. The wires are considered to be in parallel, so that the transport current is injected by using only one current constraint. We simulated one pitch length, as schematically represented in Fig. 4. The figure also shows the 2-D cross-section at half-length used to evaluate the current density distribution. The radius of each wire is 2.5 mm and the pitch length is 10 cm. The critical current value is 108 A/m2. The amplitude of the applied current is 0.7Ic. Fig. 5 shows the x, y, z components of the magnetic flux density taken on the plane situated at mid-length at the peak current instant. The magnetic flux distributions have the following features:
F. Grilli et al. / Cryogenics 53 (2013) 142–147
145
Fig. 5. Distribution of the different components of the magnetic flux density taken on the 2-D cross-section at half length represented in Fig. 4: the x, y, z components are shown in subfigures (a)–(c), respectively, whereas the absolute value is shown in subfigure (d). The values on the color scale are in teslas. The arrows indicate the direction of the magnetic flux vector.
Outside the conductor region the magnetic flux lies on the xyplane with circular patterns, similarly to the flux generated by a straight conductor carrying current. The flux partially penetrates the superconductor material. The flux generated by the three wires interacts, such that in the central region it is oriented mostly along the z direction. The computation results have been compared with those obtained with a 2-D helicoidal model based on the use of Riemannian manifolds [16]. The results for the current density distributions are very similar and so are the losses computed for different values of the transport current – see Fig. 6. The slight observed discrepancies are probably due to the coarser mesh utilized in the 3-D model.
Fig. 7 shows the current density distribution for a transport current of 350 A. The current is initially distributed between the two strips, but in the vicinity of the defect (characterized by a much lower Jc) it crosses over to the ‘‘good’’ superconductor through the resistive barrier. One can notice that far from the defect the ‘‘good’’ superconductor is not saturated with current, whereas in the vicinity of x = 0 it is, due to the extra current it must carry. This type of simulation can be useful to study current sharing phenomena, for example between filaments in striated coated conductors or between strands in Roebel-assembled coated conductor cables. The model can also be coupled to a thermal model, in order to study the propagation of the normal zone during quench phenomena, similarly to what has already been done in 2-D [19].
3.3. Current sharing caused by localized defect
3.4. Superconductor with varying cross-section
As a further case we considered the same geometry of Section 3.1, but with transport current instead and in the presence of a localized defect in one of the superconducting strips. The transport current is injected at one end of the conductors (which are therefore considered in parallel). The defect is simulated by using a Jc varying with position:
As a final example we simulated a superconductor with varying cross-section along its length carrying transport current, as represented in Fig. 8. Due to the axial symmetry, the displayed case is actually a 2-D problem, but we solve it in 3-D to show the possibility of handling conductors with varying cross-section. More realistic geometries can be built from images of real superconducting samples. Since Jc is assumed to be a constant parameter, the critical current of the displayed conductor varies along the length, according to the transversal cross-section. In the displayed case, we inject a transport current that is lower than the highest Ic and higher than
2
J c ðxÞ ¼ J c0 ½1 0:6eðx=x0 Þ ;
ð4Þ
with x0 = 0.5 mm. In other words, at x = 0 Jc drops down to 40% of its original value, with a characteristic length of 0.5 mm.
146
F. Grilli et al. / Cryogenics 53 (2013) 142–147
(a)
(b)
Fig. 8. Current density distribution on different cross-sections of a conductor with varying cross-section along its length. The transport current is applied along the zaxis. The local critical current limit is exceeded in the bottleneck.
so-called filament sausaging) of practical superconductors, such as MgB2 wires.
(c) Fig. 6. Distribution of the modulus of the current density (normalized by Jc) computed with the 3-D model (a) and the 2-D model described in [16,17] (b). Comparison of the ac losses as a function of the transport current (normalized to the critical current of each conductor), computed with the two models (c).
4. Conclusion and outlook We successfully implemented a 3-D model for time-dependent simulations of high-temperature superconductors using the H-formulation. The model has been validated by comparing the results with cases that can be solved in 2-D and has been used to show different examples of application on 3-D geometries. The model can be used to study very interesting cases for practical applications, such as the current sharing in filamentarized structures and in composite cables, the variation of the superconductor’s cross-section and the presence of localized defects. This latter is particularly interesting if thermal effects are also taken into account, as it has already been done in 2-D. The choice of an appropriate mesh is of paramount importance in order to keep the number of degrees of freedom of the problem at an acceptable level. In this regard, one can take advantage of the swept mesh, which allows to propagate a 2-D mesh in the third dimension without an excessive increase of the number of mesh elements. Acknowledgments
Fig. 7. Distribution of the absolute value of the current density (normalized by Jc0) generated by the transport current (injected along the x direction) in two thin superconductors separated by a resistive barrier. The strip on the left (y < 0) has uniform critical current density Jc = Jc0; the strip on the right (y > 0) has a Jc varying with position according to Eq. (4). Due to the presence of the defect in the right filament, part of the current crosses through the resistive barrier to the left filament, saturating it near x = 0, where one can observe that the Jc value is exceeded in all the transversal cross-section of the left tape.
The authors acknowledge the following sources of funding: Helmholtz-University Young Investigator Grant VH-NG-617 (F. Grilli); Research Fund for the Italian Electrical System under the Contract Agreement between RSE and the Ministry of Economic Development – General Directorate for Nuclear Energy, Renewable Energy and Energy Efficiency stipulated on July 29, 2009 in compliance with the Decree of March 19, 2009 (R. Brambilla); FQRNT and MITACS (F. Sirois and S. Memiaghe); Academy of Finland, project number 131577, Foundation for Technology promotion in Finland, Emil Aaltonen foundation (A. Stenvall). References
the lowest Ic. As expected, this results in the full saturation of the conductor in the ‘‘bottleneck’’ region; on the other hand, where the conductor is wider, the transport current represents only a fraction of the local Ic and the current density distribution assumes the typical Bean-like profile, with a circular band where J ’ Jc near the surface and lower values inside. This simple example shows that our 3-D model can be used to study the influence of the variation of the cross-section (the
[1] Hong Z, Campbell AM, Coombs TA. Numerical solution of critical state in superconductivity by finite element software. Supercond Sci Technol 2006;19:1246–52. [2] Brambilla R, Grilli F, Martini L. Development of an edge-element model for AC loss computation of high-temperature superconductors. Supercond Sci Technol 2007;20:16–24. [3] Thakur KP, Staines MP, Lakshmi LS, Long NJ. Numerical computation of ac losses and flux profiles in high-aspect-ratio superconducting strips in perpendicular ac magnetic field. IEEE Trans Appl Supercond 2009;19:3770–8.
F. Grilli et al. / Cryogenics 53 (2013) 142–147 [4] Grilli F, Sirois F, Brault S, Brambilla R, Martini L, Nguyen DN, et al. Edge and top/ bottom losses in non-inductive coated conductor coils with small separation between tapes. Supercond Sci Technol 2010;23:034017. [5] Rodriguez-Zermeno VM, Mijatovic N, Træholt C, Zirngibl T, Seiler E, Abrahamsen AB, et al. Towards faster FEM simulation of thin film superconductors: a multiscale approach. IEEE Trans Appl Supercond 2011;21:3273–6. [6] Nguyen DN, Ashworth SP, Willis JO, Sirois F, Grilli F. A new finite-element method simulation model for computing ac loss in roll assisted biaxially textured substrate ybco tapes. Supercond Sci Technol 2010;23:025001. [7] Ainslie MD, Rodriguez-Zermeno VM, Hong Z, Yuan W, Flack TJ, Coombs TA. An improved FEM model for computing transport AC loss in coils made of RABiTS YBCO coated conductors for electric machines. Supercond Sci Technol 2011;24:045005. [8] Nguyen DN, Grilli F, Ashworth SP, Willis JO. Ac loss study of antiparallel connected YBCO coated conductors. Supercond Sci Technol 2009;22:055014. [9] Finite-Element Software Package Comsol Multiphysics, Version 4.2., 2011.. [10] E. Pardo, J. Souc, J. Kovac, AC Loss in ReBCO Pancake Coils and Stacks of them: Modelling and Measurement, arxiv:1109.2526v1 [cond-mat.supr-con] (2011). [11] Campbell AM. An introduction to numerical methods in superconductors. J Supercond Novel Magn 2011;24:27–33.
147
[12] Clem JR, Perez-Gonzalez A. Flux-line-cutting and flux-pinning losses in type-II superconductors in rotating magnetic fields. Phys Rev B 1984;30:5041–7. [13] Badía A, López C. Vector magnetic hysteresis of hard superconductors. Phys Rev B, 2002;65:104514. [14] Badı´a-Majós A, López C, Ruiz HS. General critical states in type-II superconductors. Phys Rev B 2009;80:144509. [15] Ruiz H, Badı´a-Majós A. Smooth double critical state theory for type-II superconductors. Supercond Sci Technol 2010;23:105007. [16] Stenvall A, Tarhasaari T, Grilli F, Raumonen P, Vojencˇiak M, Pellikka M. Manifolds in electromagnetism and superconductor modelling: using their properties to model critical current of twisted conductors in self-field with 2-D model. Cryogenics, submitted for publication. [17] Stenvall A, Tarhasaari T. An eddy current vector potential formulation for estimating hysteresis losses of superconductors with FEM. Supercond Sci Technol 2010;23:125013. [18] Wilson MN. Superconducting magnets. Oxford: Clarendon Press; 1983. [19] Roy F, Dutoit B, Grilli F, Sirois F. Magneto-thermal modeling of secondgeneration HTS for resistive fault current limiter design purposes. IEEE Trans Appl Supercond 2008;18:29–35.
Contents lists available at SciVerse ScienceDirect
Cryogenics journal homepage: www.elsevier.com/locate/cryogenics
Development of a three-dimensional finite-element model for high-temperature superconductors based on the H-formulation Francesco Grilli a,⇑, Roberto Brambilla b, Frédéric Sirois c, Antti Stenvall d, Steeve Memiaghe c a
Karlsruhe Institute of Technology, Karlsruhe, Germany Ricerca Sistema Energetico S.p.A., Milano, Italy c Ecole Polytechnique de Montréal, Montréal, Canada d Tampere University of Technology, Tampere, Finland b
a r t i c l e
i n f o
Article history: Available online 8 May 2012 Keywords: Finite-element models High-temperature superconductors Ac losses 3-D simulations
a b s t r a c t Finite-element models are a powerful and widely used tool for evaluating the ac losses of HTS tapes and wires as well as of assemblies such as cables and coils. The H-formulation, which uses the magnetic field components as state variables, has proved to be an efficient implementation to solve 2-D problems, involving infinitely long or axially-symmetric geometries; an excellent agreement with experimental data has been found in many cases. However, the simulation of certain applications requires a full 3-D model. In this paper we report on the development of a 3-D model based on the H-formulation. We describe the implementation of Maxwell equations, the imposition of current constraints and we discuss the issues related to meshing 3-D volumes. The model is validated by comparing the results with those obtained with 2-D models in cases that can be investigated in 2-D; then, it is used to simulate cases that can be handled only in 3-D. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction In recent years finite-element simulations have become a powerful tool to investigate in detail the current density and magnetic field distributions inside high-temperature superconductors (HTS) and to compute their ac losses. In particular, the H-formulation [1,2] has increased in popularity both in the academic and the industrial world and it has been successfully used to study the ac loss performance not only of individual wires, but also of more complex systems like cables and coils [3–5]. The model has been adapted to take into account the presence of ferromagnetic materials as well [6,7]. Until now, the model has been only used to simulate 2-D problems, i.e. geometries representing conductors either infinitely long or characterized by rotational symmetries (coils), and excellent agreement with experimental data has been observed, both for individual and interacting tapes [8]. However, many problems require the simulation of 3-D effects, such as coupling currents, local non-uniformities of the physical or geometrical properties, and twisted or transposed geometries. An extension of the model to 3-D is therefore necessary. In this paper we present our first results of such extension. The paper is organized as follows. At the beginning we give a brief description of the model implementation. Then we describe and discuss the results of our simulations: first we validate the model ⇑ Corresponding author. E-mail address: [email protected] (F. Grilli). 0011-2275/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.cryogenics.2012.03.007
by comparing its results with those that are obtained in cases that can be simulated in 2-D; then we provide results for cases that can be analyzed only with the new model, like the current sharing between superconductors caused by the presence of localized defects and a superconductor with varying cross-section along its length. Finally, we draw the conclusions, with focus on the open issues and the main challenges of 3-D numerical modeling of superconductors. 2. Model description The basic equation of the model, which is implemented in the Partial Differential Equation (PDE) module of the finite-element software package Comsol Multiphysics [9], is Faraday’s equation
l
@H þ r E ¼ 0; @t
ð1Þ
where the components of the magnetic field H = (Hx, Hy, Hz) are the state variables. For non-magnetic materials, the relation B = l0H is assumed. The superconductor material is modeled by a non-linear resistivity
q¼
n1 Ec J ; J c J c
ð2Þ
where Ec is the critical electric field (usually equal to 1 lV/cm), Jc is the critical current density, J is the absolute value of the current density, and n is the power-law exponent. This latter typically
F. Grilli et al. / Cryogenics 53 (2013) 142–147
assumes values between 15 and 50 for high-temperature superconductors. For the simulations contained in this paper we used n = 35. The conductors are usually surrounded by vacuum, which is modeled as a resistive material with very high resistivity (1 X m). For all the simulations of this paper, the frequency of the magnetic field and current sources is 50 Hz. It is very important to remark that the use of an isotropic power-law resistivity such as (2) is not a very realistic assumption, both from the physics and modeling point of view. 2-D models have proven to be in excellent agreement with measurement when the 2-D symmetry is also realized in the measurements – see for example [8,10]. In those cases, the magnetic field lines are always perpendicular to the current, and pinning and flux flow in the microscopic scale can be used to describe what happens inside the superconductor [11]. However, in real 3-D cases the current is not perpendicular to applied field and a better determination of the microscopic model is needed. Additionally, the microscopic description needs to give input for the resistivity tensor used in the eddy current type 3-D models. Furthermore, the resistivity tensor is not necessarily only a function of the anisotropic properties in the conductor, like different Jc along c-axis or in ab-plane of coated conductors, but may depend also on the direction of the magnetic field. The so-called double critical state model was proposed by Clem and Perez-Gonzalez [12]. They utilized different Jc values in perpendicular direction and along the flux lines. There is also a lot of recent interest in the topic [13–15], but these works do not provide the solution for resistivity in eddy current type solvers yet, nor methods to easily measure just a few parameters to characterize the behavior of Jc in all possible field orientations with respect to all possible current directions in a sample. This being said, since this work is the first extension of the Hformulation from 2-D to 3-D, the simplifying assumptions of (2) are used in this paper, also because they allow a direct validation of the new model against the 2-D one and an easy verification that the electromagnetic quantities are computed correctly. Similar to the 2-D model [2], an external magnetic field can be easily applied by appropriately setting the Dirichlet boundary conditions on the domain’s boundary. A transport current I(t) can be imposed by means of a current constraint on the desired face (domain boundary), as follows. First one defines the variable describing the current flowing through a face S as
IS ¼
Z
ðJ x nx þ J y ny þ J z nz ÞdS;
ð3Þ
S
where (nx, ny, nz) is the unit vector defining the normal direction of the face S and (Jx, Jy, Jz) are the three coordinate components of the current density; then, one imposes I(t) IS = 0 as an algebraic constraint. This latter forces the current flowing through the face S to be equal to I(t) at each time step of the simulation. In this way, when different conductors are present, multiple independent current constraints can be imposed. The choice of the mesh is very important in order to have sufficiently precise results, yet keeping the total number of degrees of freedom of the problem at a reasonable level. In this regard, the use of the so-called swept mesh can be very useful. A swept mesh is created by starting from a mesh of the 2-D cross-section of the conductor region; the 2-D mesh (which can be quite fine) is replicated along the length of the conductor. An example is represented in Fig. 1, where one can also notice that the mesh density along the conductor’s length (x-axis) is not uniform, but is refined in the proximity of the region where important changes of the current density and magnetic field distributions occur – see Section 3.3 for more details. This type of mesh is usually more convenient than the standard free tetrahedral mesh, even if its use might not be possible or straightforward with certain geometries. For particularly complex geometries, the use of specifically dedicated meshing
143
Fig. 1. Example of swept mesh. The mesh is first generated on a 2-D cross-section (yz plane) and then extruded and replicated along the third dimension (x direction). In this example the mesh is refined near mid-length along the x-axis, because it is used to study the case of a localized defect there – see Section 3.3 for details.
programs can be helpful. As far as the discretization of the elements is concerned, first-order edge elements, which guarantee the continuity of the tangential field component across adjacent elements, are used. The problems shown in this paper contain approximately between 15,000 and 200,000 mesh elements. As a consequence, the computation time is very variable. Even for the same problem, the computation time can change dramatically depending on the value of the applied field or current. Sometimes, in order to obtain convergence, a very small time step has to be set for the solver, which inevitably results in longer computation times. As an approximate figure, we can say that the simulation of one cycle at 50 Hz takes between a few hours and 1–2 days on a 8-core workstation equipped with 8 GB of RAM. 3. Results In this section we present the results for different 3-D geometries. The first two cases serve also as validation of the model, since they can be solved by 2-D models (at least to a certain extent). For the case of helicoidal superconductors (Section 3.2), a comparison with a 2-D model utilizing the combination of translation and rotation symmetry is performed. This eddy current model is formed from the magnetostatic model of helicoidal conductors [16] coupled to Faraday’s law, similarly to what has been done in [17].
Fig. 2. Schematic illustration of the coupling currents induced by a varying magnetic field in two thin superconductors separated by a normal metal barrier. The field is applied perpendicular to the flat face of the superconductors (z-axis). In the illustrated case part of the current induced in the xy-plane flows from one superconductor to the other through the barrier (partial coupling).
144
F. Grilli et al. / Cryogenics 53 (2013) 142–147
(a)
(b)
Fig. 3. Magnetic flux (a) and current density (b) profiles across the width (y-axis in Fig. 2) of two thin superconductors separated by a resistive barrier for different values of the barrier resistivity. The profiles computed with a 2-D model for the fully uncoupled and fully coupled situations are also shown. The profiles are taken at the peak of a field of 50 mT of amplitude at 50 Hz, applied perpendicular to the superconducting strips.
Fig. 4. View of the simulated helicoidal conductors (one pitch length). The plane at half length represents the 2-D cross-section, where the magnetic flux density distributions of Fig. 5 are taken.
The other simulated problems are examples of applications of the 3-D model in cases of practical interest, such as the current sharing due to local non-uniformities of the superconductor and the variation of the superconductor’s cross-section. 3.1. Coupling between superconductors in external ac magnetic field We first considered the case of two superconducting thin strips separated by a resistive barrier under the action of a perpendicular ac field. The magnetic field induces current loops. According to the considered conditions and geometry (variation rate of the magnetic field, length of the conductors, resistance of the barrier), the induced currents can flow separately in each superconducting strip (uncoupled situation) or from one strip to the other (coupled situation) [18]. The cases of full uncoupling and full coupling can be simulated in 2-D by appropriately setting the current constraints,
and they can be used to validate the 3-D model. The case of partial coupling (when only a fraction of the induced current crosses the resistive barrier) is schematically represented in Fig. 2 and can be simulated only with a 3-D model. We simulated two superconducting strips, 2 mm wide, 5 lm thick, each with Ic = 200 A, separated by a normal metal barrier 0.5 mm wide. We applied a field of 50 mT and observed the magnetic field and current density profiles along the y-axis of Fig. 2. The resistivity of the barrier was varied in order to consider different degrees of coupling between the strips. Fig. 3 shows the magnetic field and current density profiles at the peak of the applied field for different values of the resistivity of the metal barrier, along with the results obtained with 2-D simulations by imposing perfect uncoupling and coupling between the superconductors, respectively. It can be seen that for a sufficiently high value of the resistivity the 3-D results agree very well with the 2-D uncoupled case; and for a sufficiently low value of the resistivity, the 3-D results agree very well with the 2-D coupled case. This shows that the 3-D model works well and the electromagnetic quantities are computed correctly. For an intermediate value of the resistivity, one has a situation of partial coupling: in both superconductors the induced current density has both positive and negative components, but in one superconductor it is predominantly positive, in the other one it is predominantly negative. This means that part of (but not all) the induced current flows from one superconductor to the other – see also Fig. 2. This is a typical case that cannot be handled with 2-D simulations. A similar change of the degree of coupling can be obtained by varying the frequency of the applied magnetic field or the conductors’ length.
3.2. Helicoidally-wound superconductors As a second example, we simulated three helicoidally wound round superconductors carrying transport current. The wires are considered to be in parallel, so that the transport current is injected by using only one current constraint. We simulated one pitch length, as schematically represented in Fig. 4. The figure also shows the 2-D cross-section at half-length used to evaluate the current density distribution. The radius of each wire is 2.5 mm and the pitch length is 10 cm. The critical current value is 108 A/m2. The amplitude of the applied current is 0.7Ic. Fig. 5 shows the x, y, z components of the magnetic flux density taken on the plane situated at mid-length at the peak current instant. The magnetic flux distributions have the following features:
F. Grilli et al. / Cryogenics 53 (2013) 142–147
145
Fig. 5. Distribution of the different components of the magnetic flux density taken on the 2-D cross-section at half length represented in Fig. 4: the x, y, z components are shown in subfigures (a)–(c), respectively, whereas the absolute value is shown in subfigure (d). The values on the color scale are in teslas. The arrows indicate the direction of the magnetic flux vector.
Outside the conductor region the magnetic flux lies on the xyplane with circular patterns, similarly to the flux generated by a straight conductor carrying current. The flux partially penetrates the superconductor material. The flux generated by the three wires interacts, such that in the central region it is oriented mostly along the z direction. The computation results have been compared with those obtained with a 2-D helicoidal model based on the use of Riemannian manifolds [16]. The results for the current density distributions are very similar and so are the losses computed for different values of the transport current – see Fig. 6. The slight observed discrepancies are probably due to the coarser mesh utilized in the 3-D model.
Fig. 7 shows the current density distribution for a transport current of 350 A. The current is initially distributed between the two strips, but in the vicinity of the defect (characterized by a much lower Jc) it crosses over to the ‘‘good’’ superconductor through the resistive barrier. One can notice that far from the defect the ‘‘good’’ superconductor is not saturated with current, whereas in the vicinity of x = 0 it is, due to the extra current it must carry. This type of simulation can be useful to study current sharing phenomena, for example between filaments in striated coated conductors or between strands in Roebel-assembled coated conductor cables. The model can also be coupled to a thermal model, in order to study the propagation of the normal zone during quench phenomena, similarly to what has already been done in 2-D [19].
3.3. Current sharing caused by localized defect
3.4. Superconductor with varying cross-section
As a further case we considered the same geometry of Section 3.1, but with transport current instead and in the presence of a localized defect in one of the superconducting strips. The transport current is injected at one end of the conductors (which are therefore considered in parallel). The defect is simulated by using a Jc varying with position:
As a final example we simulated a superconductor with varying cross-section along its length carrying transport current, as represented in Fig. 8. Due to the axial symmetry, the displayed case is actually a 2-D problem, but we solve it in 3-D to show the possibility of handling conductors with varying cross-section. More realistic geometries can be built from images of real superconducting samples. Since Jc is assumed to be a constant parameter, the critical current of the displayed conductor varies along the length, according to the transversal cross-section. In the displayed case, we inject a transport current that is lower than the highest Ic and higher than
2
J c ðxÞ ¼ J c0 ½1 0:6eðx=x0 Þ ;
ð4Þ
with x0 = 0.5 mm. In other words, at x = 0 Jc drops down to 40% of its original value, with a characteristic length of 0.5 mm.
146
F. Grilli et al. / Cryogenics 53 (2013) 142–147
(a)
(b)
Fig. 8. Current density distribution on different cross-sections of a conductor with varying cross-section along its length. The transport current is applied along the zaxis. The local critical current limit is exceeded in the bottleneck.
so-called filament sausaging) of practical superconductors, such as MgB2 wires.
(c) Fig. 6. Distribution of the modulus of the current density (normalized by Jc) computed with the 3-D model (a) and the 2-D model described in [16,17] (b). Comparison of the ac losses as a function of the transport current (normalized to the critical current of each conductor), computed with the two models (c).
4. Conclusion and outlook We successfully implemented a 3-D model for time-dependent simulations of high-temperature superconductors using the H-formulation. The model has been validated by comparing the results with cases that can be solved in 2-D and has been used to show different examples of application on 3-D geometries. The model can be used to study very interesting cases for practical applications, such as the current sharing in filamentarized structures and in composite cables, the variation of the superconductor’s cross-section and the presence of localized defects. This latter is particularly interesting if thermal effects are also taken into account, as it has already been done in 2-D. The choice of an appropriate mesh is of paramount importance in order to keep the number of degrees of freedom of the problem at an acceptable level. In this regard, one can take advantage of the swept mesh, which allows to propagate a 2-D mesh in the third dimension without an excessive increase of the number of mesh elements. Acknowledgments
Fig. 7. Distribution of the absolute value of the current density (normalized by Jc0) generated by the transport current (injected along the x direction) in two thin superconductors separated by a resistive barrier. The strip on the left (y < 0) has uniform critical current density Jc = Jc0; the strip on the right (y > 0) has a Jc varying with position according to Eq. (4). Due to the presence of the defect in the right filament, part of the current crosses through the resistive barrier to the left filament, saturating it near x = 0, where one can observe that the Jc value is exceeded in all the transversal cross-section of the left tape.
The authors acknowledge the following sources of funding: Helmholtz-University Young Investigator Grant VH-NG-617 (F. Grilli); Research Fund for the Italian Electrical System under the Contract Agreement between RSE and the Ministry of Economic Development – General Directorate for Nuclear Energy, Renewable Energy and Energy Efficiency stipulated on July 29, 2009 in compliance with the Decree of March 19, 2009 (R. Brambilla); FQRNT and MITACS (F. Sirois and S. Memiaghe); Academy of Finland, project number 131577, Foundation for Technology promotion in Finland, Emil Aaltonen foundation (A. Stenvall). References
the lowest Ic. As expected, this results in the full saturation of the conductor in the ‘‘bottleneck’’ region; on the other hand, where the conductor is wider, the transport current represents only a fraction of the local Ic and the current density distribution assumes the typical Bean-like profile, with a circular band where J ’ Jc near the surface and lower values inside. This simple example shows that our 3-D model can be used to study the influence of the variation of the cross-section (the
[1] Hong Z, Campbell AM, Coombs TA. Numerical solution of critical state in superconductivity by finite element software. Supercond Sci Technol 2006;19:1246–52. [2] Brambilla R, Grilli F, Martini L. Development of an edge-element model for AC loss computation of high-temperature superconductors. Supercond Sci Technol 2007;20:16–24. [3] Thakur KP, Staines MP, Lakshmi LS, Long NJ. Numerical computation of ac losses and flux profiles in high-aspect-ratio superconducting strips in perpendicular ac magnetic field. IEEE Trans Appl Supercond 2009;19:3770–8.
F. Grilli et al. / Cryogenics 53 (2013) 142–147 [4] Grilli F, Sirois F, Brault S, Brambilla R, Martini L, Nguyen DN, et al. Edge and top/ bottom losses in non-inductive coated conductor coils with small separation between tapes. Supercond Sci Technol 2010;23:034017. [5] Rodriguez-Zermeno VM, Mijatovic N, Træholt C, Zirngibl T, Seiler E, Abrahamsen AB, et al. Towards faster FEM simulation of thin film superconductors: a multiscale approach. IEEE Trans Appl Supercond 2011;21:3273–6. [6] Nguyen DN, Ashworth SP, Willis JO, Sirois F, Grilli F. A new finite-element method simulation model for computing ac loss in roll assisted biaxially textured substrate ybco tapes. Supercond Sci Technol 2010;23:025001. [7] Ainslie MD, Rodriguez-Zermeno VM, Hong Z, Yuan W, Flack TJ, Coombs TA. An improved FEM model for computing transport AC loss in coils made of RABiTS YBCO coated conductors for electric machines. Supercond Sci Technol 2011;24:045005. [8] Nguyen DN, Grilli F, Ashworth SP, Willis JO. Ac loss study of antiparallel connected YBCO coated conductors. Supercond Sci Technol 2009;22:055014. [9] Finite-Element Software Package Comsol Multiphysics, Version 4.2., 2011.
147
[12] Clem JR, Perez-Gonzalez A. Flux-line-cutting and flux-pinning losses in type-II superconductors in rotating magnetic fields. Phys Rev B 1984;30:5041–7. [13] Badía A, López C. Vector magnetic hysteresis of hard superconductors. Phys Rev B, 2002;65:104514. [14] Badı´a-Majós A, López C, Ruiz HS. General critical states in type-II superconductors. Phys Rev B 2009;80:144509. [15] Ruiz H, Badı´a-Majós A. Smooth double critical state theory for type-II superconductors. Supercond Sci Technol 2010;23:105007. [16] Stenvall A, Tarhasaari T, Grilli F, Raumonen P, Vojencˇiak M, Pellikka M. Manifolds in electromagnetism and superconductor modelling: using their properties to model critical current of twisted conductors in self-field with 2-D model. Cryogenics, submitted for publication. [17] Stenvall A, Tarhasaari T. An eddy current vector potential formulation for estimating hysteresis losses of superconductors with FEM. Supercond Sci Technol 2010;23:125013. [18] Wilson MN. Superconducting magnets. Oxford: Clarendon Press; 1983. [19] Roy F, Dutoit B, Grilli F, Sirois F. Magneto-thermal modeling of secondgeneration HTS for resistive fault current limiter design purposes. IEEE Trans Appl Supercond 2008;18:29–35.