3D modeling of a superconducting dynamo-type flux pump

High temperature superconducting (HTS) dynamos are promising devices that can inject large DC currents into the winding of superconducting machines or magnets in a contactless way. Thanks to this, troublesome brushes in HTS machines or bulky currents leads with high thermal losses will be no longer required. The working mechanism of HTS dynamo has been controversial during the recent years and several explanations and models have been proposed to elucidate its performance. In this paper, we present the first three-dimensional (3D) model of an HTS flux pump, which has good agreement with experiments. This model can be beneficial to clarify the mechanism of the dynamo and pinpoint its unnoticed characteristics. Employing this model, we delved into the screening current and electric field distribution across the tape surface in several crucial time steps. This is important, since the overcritical screening current has been shown to be the reason for flux pumping. In addition, we analyzed the impact of both components of electric field and screening current on voltage generation, which was not possible in previous 2D models. We also explored the necessary distance of voltage taps at different airgaps for precise measurement of the voltage across the tape in the dynamo.

Since 1960s, when the first generation of superconducting flux pumps have been designed and operated, the flux pump technology has undergone a remarkable developement [1][2][3][4][5] . Within around 60 years, various mechanically or electrically driven superconducting flux pumps have been operated using low-temperature superconductor (LTS) materials in their previous generation and high-temperature superconductor (HTS) materials in their last generation.
Superconducting dynamos convert mechanical energy into electromagnetic energy. They use permanent magnets as the source of varying magnetic field to induce DC volatge inside a superconducting wire. In their previous generation with LTS wires, creating normal region in superconductor was necessary to generate voltage and the local magnetic field needed to exceed the critical magnetic field of the wire. Generating this high amount of magnetic field in commercial LTS wires is not always possible, which limits the application of LTS superconducting dynamos. However, since Hoffmann and et al. designed and implemented the new generation of flux pumps with HTS materials 6 , the application of flux pumps and in particular HTS dynamos has became popular and widespread. In this type of dynamos, the existence of normal regions in superconductor is no more required and the voltage can be created even with very low amount of alternating magnetic field 7 . The only necessary characteristic for a material to be used in the new generation of flux pump is having non-linear resistivity 8 . This makes the HTS dynamos a reliable and efficient candidate to be utilized in superconducting electrical machines [9][10][11][12][13][14][15][16] .
The simplified working mechanism of an HTS flux pump can be explained as follows. Any magnitude of varying magnetic field causes circulating screening currents in the tape surface. These screening currents result in electric field in the areas above the local critical current of the tape. The accumulation of this induced electric field in one complete cycle results in certain net voltage, which causes flux pumping. However, this is not the case in materials with linear resistivity, where the net voltage is always zero. In other words, having non-linear resistivity is the key and the only necessary characteristic of a material to be used in an HTS flux pump. In addition, the local critical current is subjected to alternation based on the local magnetic field density and its angle with respect to tape surface. The J c (B, θ) dependency can be measured via experiment, where θ is the magnetic field angle with respect to the tape surface.
Although the mechanism of HTS dynamo have been investigated and explained in 7,8,17 using experiments and 2D modeling, there are still many aspects that have been left untouched. These issues can hardly be disclosed by experiments and are not feasible to be studied via 2D modeling. The first one is the complete distribution of screening current on tape surface during flux pumping process. The current density and electric field on the tape www.nature.com/scientificreports/ surface have two components: J x , E x along the tape width and J y , E y along the tape length. 2D modeling assumes infinitely long tape along the length, thus there are only J y and E y components for current density and electric field and J x and E x are missing. Considering these facts and limitations of 2D modeling, the distribution of overcritical screening current, which has been shown to be the reason for flux pumping and its resultant electric field was not thoroughly pinpointed in previous literature. That is why the necessity of a 3D model to certify the previous proposed mechanisms of flux pumping and to better clarify the procedure is inevitable. By using an accurate and efficient 3D model, it is possible to pinpoint the mechanism of flux pumping and its details. In addition, 3D modeling gives the opportunity to analyze more complex geometries, regarding both tape or magnet shapes, which is also not possible via 2D modeling. In previous literature, there have been several works focusing on modeling of flux pumps in order to investigate its performance. In 18 a simplified 2D model using A-formulation was proposed, which could simulate the HTS dynamo mechanism. Although this model could describe some key aspects of flux pumping, it was unrealistic. In 8 a 2D finite element numerical model based on H-formulation was presented, which had fairly good agreement with experiment and could explain well the mechanism behind the voltage generation in an HTS dynamo in open-circuit condition. An efficient and fast model, yet with good agreement with experiment based on minimum electromagnetic entropy production method (MEMEP) was proposed to model the procedure of voltage generation in an HTS dynamo in open-circuit mode 7 . In 17 , Mataira and et al. extended their 2D model based on H-formulation to current-driven mode to further explain the underlying physics of voltage generation in an HTS dynamo. They also validated their modeling results with experiments. In 19 a benchmark of several different models for modeling of HTS dynamo in open-circuit mode was presented. This benchmark includes H-formulation with shell current method 8,17,20 , MEMEP method 7 , coupled H-A method 21 , coupled T-A formulations 22 , Segregated H-formulation 23 , integral equation 24 and volume integral equation-based equivalent circuit method 25 . All of the methods showed very good agreement with each other either qualitatively or quantitatively and the MEMEP method was the fastest in terms of computational speed.
In this article, we present the first 3D model of a dynamo-type HTS flux pump. Our modeling method, being MEMEP 3D, utilizes optimum number of degrees of freedom by solving all the variables of the problem only inside the superconducting region, which makes the model very fast and efficient. Employing the proposed model, the performance of the dynamo is investigated and analyzed. In particular, we pay special attention to the necessary distance of voltage taps at different airgaps for precise measurement of the voltage across the tape in the dynamo. The article structure is as follows. First, we explain the modeling method, including the introduction of the MEMEP 3D method, magnet modeling process and the model configuration. Afterwards, the article analyses the screening current and electric field in several critical time steps of the magnet movement over the HTS tape. Then, the impact of length on voltage generation and measurement is investigated and, at the end, the model results are compared to experiments in various airgaps.
The main results of this article were presented in a poster in the Applied Superconductivity Conference 2020, available on-line from 23rd of October 2020 26 , and a pre-print of this article was available in arXiv (arXiv:2011.00989) from 2nd of November 2020. A few days afterwards, a pre-print from L Prigozhin and V Sokolovsky was uploaded to arxiv (arXiv:2011.03524) reporting also modeling of flux pumps of finite length, which was finally published in 27 .

Dynamo modeling
MEMEP 3D method. In this article, we model the HTS dynamo by means of the Minimum Electro-Magnetic Entropy Production method in 3D (MEMEP 3D). This variational method is based on T-formulation and is able to exploit several strategies to speed up calculation such as parallel computing, dividing into sectors and symmetry. The method works based on minimizing a functional containing the magnetic vector potential A and the current density J . It is proven that the minimum of the 3D functional for any given time step is the unique solution of Maxwell differential equations [28][29][30][31] .
For magnetization problems, the method uses the effective magnetization T as state variable, defined as where J is the current density. Configurations with transport current are possible after adding an additional term to the equation above 29 . The general relation between the electric field and vector and scalar potentials is and the current conservation equation is where ∂ t A is the change of vector potential with respect to time at a given location r and ϕ is the scalar potential. Equation (3) is always satisfied because ∇ · (∇ × T) = 0 . Thus, we only need to solve Eq. (2). In addition, the Coulomb's gauge, ∇ · A = 0 , has been assumed to solve the 3D problem, and hence ϕ becomes the electrostatic potential 32 .
Since the current density J and T only exist inside the material, meshing is only needed in this region, increasing the calculation speed significantly 29 . The vector potential A in the functional has two contributions including A a and A J standing for the vector potential due to the applied field and the vector potential due to the current density in the superconductor, respectively. The A a component can be replaced by the vector potential due to the magnet in a dynamo, while A J is calculated with the following volume integral of current density where U is dissipation factor, which can include any E − J relation of superconductor To solve the problem in time dependent mode, the functional is minimized in discrete time steps. Assuming the functional variables in a certain time step t 0 as T 0 , A J0 and A a0 , the upcoming time step will be t = t 0 + t and the variables in the time step t will be T = T 0 + T , A J = A J 0 + A J and A a = A a0 + A a , respectively, where T , A T and A a are the variables difference between two time steps and t is the time difference between the two time steps, which does not need to be constant.
The optimum number of mesh for a tape with dimensions 12 × 48 × 0.001 in mm in each direction is chosen as 50 × 60 × 1 along the x, y and z axis, respectively inside the superconductor. The results remain almost unchanged for higher number of meshes. However, for having better resolution on the tape surface, instead of 60 mesh, 200 meshes has been chosen along the tape length. The thin film approximation (one mesh along the tape thickness) has been considered since it has no effect on the results due to low thickness of the tape. However, the model enables to take several elements in the thickness into account 33 . 3D cylindrical magnet modeling. As explained in section "MEMEP 3D method", for our MEMEP 3D method, we need to calculate the vector potential due to the magnet as A a in the functional. For calculations using J c (B, θ) , the numerical method also requires the magnetic field contribution from the magnet. Therefore, all we need is to calculate the vector potential and magnetic field of the magnet in certain observation points inside the tape, which are the central locations of each cell. Since MEMEP 3D only needs to mesh the superconductor but not the air, there are no rotating mesh issues due to moving magnets. Instead, there is just a time-varying external vector potential and magnetic field created by the magnet, A M and B M , in the superconductor volume. These quantities are calculated only once at the beginning of each time step, representing a small fraction of the total computing time. The numerical method how to calculate these quantities is detailed in the Appendix.
The modeled permanent magnet selected for this study is a cylinderical magnet with the diameter and height equal to 10 mm 34 . The magnet type is a N42 (Nd-Fe-B magnet), possessing the remanence magnetic field of around 1.3 T.

HTS tape modeling.
The isotropic E-J power law is implemented in the functional within the dissipation factor of (6), assuming the following non-linear characteristic of the superconductor where E c = 10 −4 V/m is the critical electric field, J c is the critical current density and n is the n-value of the HTS tape.
The parameters of the modeled Superpower SF12050CF wire including J c (B, θ) data are derived from 8 . The tape has 12 mm width, 1 µ m thickness, 48 mm length, n-value of 20, and the critical current at self field of 281 A. The modeled tape has similar (but not identical) characteristics to the HTS tape used in 34 , where we aim to compare our modeling results with the experimental results of that article. Figure 1a shows the experimental J c (B, θ) dependence, where θ is the angle of the applied magnetic field with respect to the normal vector of the tape surface (Fig. 1b). Since the measurements were not performed in the whole 360 • range, we only use the measured critical current in the range between 0° and 180° and then we assume that it is symmetrical in the range between − 180° and 0° for the model. In contrast to 2D modeling, in 3D modeling there is another angle, showed as φ in Fig. 1b, which is the angle of applied magnetic field with respect to the x axis. Although our modeling method allows this kind of J c (B, θ , φ) dependence 31 , the J c measurements were done with φ = π/2 only. Indeed, complete J c (B, θ , φ) measurements are scarce for any type of sample due to experimental complexity. Therefore, our model only considers the angle θ , while angle φ is neglected with the assumption that B x plays the same role as B y . In other words, we assume that J c (B, θ , φ) = J c (B, θ , π/2) for any φ . As seen for other types of experiments, such as cross-field demagnetization, this assumption does not have a severe influence on the electromagnetic behavior 33 .
Model configuration. Figure 2a,b shows the configuration of the 3D model in 2 different views. The magnet rotates in the xz-plane on top of the tape in counter clockwise direction, where the rotor external radius is 35 mm. The magnet magnetization is pointed out to the outside of rotation circle. θ M is defined as the magnet angle, where the magnet starting position ( θ M = 0 ) is set when its center is aligned with the positive direction of z-axis. The airgap is defined as the minimum distance between the magnet outer surface and tape upper surface happening at ( θ M = 180 ). The rotating frequency is set to 12.3 Hz to be comparable with measurements conducted in 34 . About the voltage in a flux pump. In this article, we take into account that the voltage taps are placed at the ends of the tape, they are made of a conductor of negligible resistance, and they are very far away from the magnet at any time (Fig. 2c). With these assumptions, the output instantaneous open-circuit voltage of the flux pump is calculated using Eq. (2) as where l is the tape length and w is the tape width, and A y and E y are y components of total vector potential and electric field, respectively. The DC value of the instantaneous open-circuit voltage is where f is the rotation frequency of the magnet. V DC only depends on the electric field generated by the non-linear resistivity of tape, since the total vector potential A is periodic in a full cycle 7 .
Another important parameter is V since it can be obtained directly from measurements, which is defined as output voltage difference of the tape at 77 K (superconducting state) and at 300 K (normal state)  where V oc denotes the open-circuit voltage, A J,av is the volume-average vector potential due to superconducting screening current, and l is the tape effective length.

Modeling results and discussion
Current and electric field analysis. In this section, all the modeling results have been calculated and presented using constant J c , since it eliminates the complications of J c variations due to inhomogeneous magnetic field from the magnet and simplifies the performance study of the dynamo. Later, in section "Comparison to experiments", we take the J c (B, θ) dependence into account.
To examine better the flux pump behavior, we defined several key positions of the magnet while traversing the HTS tape. Figure 3a shows the schematic of these positions along with their magnet rotation angle θ M . These positions belong to the second cycle of the magnet rotation when the tape is already fully saturated with the screening current and has passed the transient state. Figure 3b shows the graph of frequency-normalized V versus magnet angle θ M . As it can be seen, in position A, while the magnet projection touches the tape, pumping voltage has already started. This is because the magnetic field of the magnet already appeared inside the tape before this time, starting from around θ M =154 • . At position C, although the magnet is in the middle of the tape, the curve minimum happens shortly afterwards, which can be explained by the hysteresis nature of the tape screening currents. Figure 4b displays the current density (maps of the modulus and current lines) for key magnet positions (see Fig. 4a).
Step 0 belongs to the beginning of second cycle when the magnet is still far away from the tape and the tape is fully saturated with screening currents remained from the first cycle. In this step, no overcritical currents in the tape are observed (see also Fig.4c related to J c -normalized J y profiles in the middle of the tape length). At step A, the magnet projection touches the tape and starts to enter it from the left side. Since the magnet already creates high applied magnetic fields on the tape, overcritical currents can be observed in the left side of the tape. These circulating currents are induced in such a way that their resultant magnetic field can oppose the magnetic field of magnet and eventually repel it (based on Lenz's law). At step B, half of the magnet has entered the tape and the magnet magnetic field has occupied the whole tape. Thus, the overcritical currents can be observed throughout the tape. At step C, the magnet has reached the middle of the tape while the tape is still fully occupied with the overcritical screening currents. The screening currents direction has been altered compared to step B because of changing the tendecy of the magnet magnetic field inside the tape from increasing to decreasing mode. At step D, the situation is similar to step B but with opposite direction of the overcritical screening currents. At step E, while the magnet is leaving the tape, the tape is still fully saturated with the overcritical screening currents and the direction of these currents are the same as step 0. While the magnet moves far enough from the (10)  www.nature.com/scientificreports/ magnet, the overcritical currents in the tape vanish completely and the remanent screening currents in the tape become similar to step 0. Far away from the magnet and the main screening currents, currents of lower magnitude (just below J c ) appear on the tape edges on the sides and, with even lower magnitude, at the ends. These currents have higher values along the left edge of the tape, where the magnet enters the tape and become maximum at step D. As it will be shown in Fig. 5, these low values of screening current have insignificant contribution in generating electric field in the tape. For example at step D, the electric field due to these low values of screening currents are 3 orders of magnitudes smaller than the maximum electric field in the tape at their highest values. Figure 5b depicts electric field maps related to the key magnet positions of Fig. 5a and Fig. 5c shows the profiles of the y component of the electric field E y , in the middle of the tape length (the E x component vanishes there). The main area in the tape responsible for generating voltage in the tape lies under and along the magnet cross section (area from 19 to 29 mm along y axis), where the magnet is traversing. At steps B and D, the largest yellow areas are observed, showing the high values of electric field induced in these areas (around 0.4 V/m), which justifies the occurrence of the minimum peaks in Fig. 3 at these steps.
In Fig. 6, we study the x and y components of the current density and electric field in more detail in a certain time step, being when the magnet is just on the top and concentric to the tape (step C). In Fig. 6a it is obvious that J x stays just below the critical current while J y reaches up to 1.5 times of J c [see Fig. 6b]. As a result, E y is almost 20 times higher than E x [see Fig. 6c,d]. In addition, since J x is symmetric along the length, E x is also symmetric [see Fig. 6a,c)]. Thanks to this symmetry, the average E x in the whole sample vanishes, and thence E x does not contribute to the voltage.
In addition, it is worth mentioning that the maximum value of J x occurs at steps D and E with maximum J c -normalized values of 1.23 and 1.20, respectively. The resultant E x due to these currents is still insignificant compared to E y generated by J y . This E x can be recognized as vague blue areas in Fig.5b at steps B and D outside the magnet cross section area (above and below the area between 19 and 29 mm, where the magnet traverses).

Impact of tape length on voltage generation.
Studying the impact of the tape length on voltage generation assists us to realize what is the minimum efficient tape length (along the y axis) that can be employed in a flux pump. From another point of view, this study help us to select the proper distance between voltage taps while measuring the voltage signals in a flux pump. This distance is of course a variable of magnet dimension, especially along the length. In this section, we study this effect with the same flux pump configuration for various airgaps. Now, we use J c (B, θ) dependency in order to obtain a realistic description that is closer to experiments. The tape length, changing between 5 mm to 48 mm, and the airgap, varying from 1 to 10 mm, are chosen as variables while the frequency is kept constant as 12.3 Hz. Figure 7a shows the result of V normalized by the frequency against the magnet angle θ M for different tape lengths and for 3.3 mm airgap. It is clear that above 20 mm length, the results of V remain almost the same  www.nature.com/scientificreports/ and the curves coincide to each other. Also, Fig. 7b shows the DC open-circuit voltage value versus tape length for various airgaps, from 1 mm to 10 mm. Up to the airgap of 6 mm, the DC voltage changes only a few percent with increasing the tape length from 24 to 48 mm. This can be justified by the distribution of current modulus on the tape surface for different airgaps, shown in Fig. 7c. In this figure, the dashed green rectangle show the area on the tape surface containing the flow of y component of current density J y , that is responsible for creation of voltage in the tape. Almost all of the created voltage in the tape is generated in this area. As the airgap increases, the rectangular becomes longer accordingly. This means that current distributes over a broader area across the tape, and thus the voltage is generated over a wider area. This fact also can be recognized in Fig. 7b by the decrease in the slope of the curves from 1 to 10 mm airgap, showing that in smaller airgaps, the induced voltage is more concentrated in the center of the tape, where the magnet traverses. For voltage signal measurements, the voltage taps should be placed in any area of the tape outside the two sides of the rectangle in order to measure the voltage signal accurately.
As a conclusion, up to airgap of 6 mm, the voltage taps distance should be at least 2.5 times larger than the magnet length. However, higher than this airgap, this distance is not sufficient anymore and for example at the airgap of 10 mm, the voltage taps distance should be around 3.5 times larger than the magnet length to be able to capture the voltage signal with highest accuracy.
Comparison to experiments. In this section, we compare 3D modeling results of the modeled dynamotype HTS flux pump with experiments published in 34 and with 2D MEMEP method presented in 7 . The calculations have been conducted with J c (B, θ) dependency and with frequency of 12.3 Hz to be comparable with experimental results. Figure 8 shows the results of the calculated and measured DC open-circuit voltage values for 5 different airgaps from 2.4 to 10.4 mm. Looking at the figure, it is obvious that the MEMEP 3D presents very good agreement with experiments and even better agreement than the MEMEP 2D method. In the case of MEMEP 2D, the decay of the DC voltage with the airgap is less pronounced than the measurements, with substantial differences at high airgaps. The reason is that the magnetic field of the magnet at 2D decays with the distance, r, as 1/r 2 www.nature.com/scientificreports/ (infinitely long dipole), while for the real cylindrical magnet the magnetic field decays as 1/r 3 (point dipole). This is the reason why 3D modeling agrees better with experiments. It is worth mentioning that, as discussed also in 7 , the airgap measurement error is not negligible due to contraction and mechanical instability of the flux pump in liquid nitrogen bath (at least ±0.5 mm airgap). This error is mostly noticeable in the lowest airgap of 2.4 mm, where a small difference in the airgap causes a high difference in the magnetic field from the magnet, and hence DC voltage. Should the gap be more accurately measured, we expect better agreement between model and experiments at low gaps. Moreover, another source of discrepancy can be that the J c (B, θ) in the model does not exactly correspond to that of the measurements, since the original tape is no longer available. As well, we needed to assume that J c (B, θ, φ) = J c (B, θ , π/2) for any φ because of unavailable measurements of this characteristics (see section "HTS tape modeling").

Conclusion
In this article, the first 3D model of an HTS flux pump has been presented and verified against experiment. This model can help to clarify the complicated mechanism of HTS dynamo-type flux pumps and verify the previously proposed mechanism by 2D models. The unique features of the presented 3D model, based on Minimum Electro-Magnetic Entropy Production (MEMEP) enable to perform computations in conveniently short times (solving the Maxwell equations only inside the superconductor, fast solver by division into sectors and iteration, and use of parallel computing). As found for cross-sectional 2D methods 19 , MEMEP 3D could likely be faster than conventional commercial finite-element methods, since it avoids the problem of rotating mesh in the air. However, a fair comparison will require a dedicated benchmark with several methods. Employing MEMEP 3D, we studied the screening current distribution, which in its overcritical form is the reason for flux pumping phenomenon, across the tape surface and in several key positions of magnet movement. We investigate the role of the x and y components of electric field and screening current in generating the voltage in the flux pump. It was found that, while the magnet traverses the tape, the maximum value of the y component of the electric field, E y , is around 20 times larger than E x , which highlights the role of J y to generate the voltage throughout the tape surface. In addition, the average E x across the sample vanishes because of symmetry. We also studied the effect of tape length on voltage generation in the flux pump. We found that the distance between voltage taps for precise measuring of the voltage depends on not only the diameter of the magnet but also the distance between magnet and tape surface. This distance should be at least 2.5 times and 3.5 times larger than the magnet diameter for airgaps less than 0.6 and around 1 times the diameter, respectively. We expect that the minimum tape length will keep increasing with the separation. Finally, the calculated modeling results were compared against experiments for several airgaps, which showed very good agreement.

Appendix
In this appendix, we outline how to calculate A M and B M in the cylindrical coordinate system centered on the magnet and local z axis, z ′ , pointing outwards from the center of the rotor of Numerical evaluation of A and B from the magnet. First, we take into account that the magnetic flux density and vector potential in Coulomb's gauge, generated by a uniformly magnetized body is the same as those generated by the effective surface magnetization current density, K M , as www.nature.com/scientificreports/ where M is the magnetization and e n is the unit vector perpendicular to the surface pointing outwards. For the cylindrical magnet, we assume that M is uniform and parallel to the axial direction. Then, in the cylindrical coordinate system with origin at the magnet center, M = Me z ′ , where e z ′ is the unit vector in the z ′ axis. As a result, K M follows the angular direction, K M = Me ϕ ′ , being e ϕ ′ the unit vector in that direction. We numerically calculate B and A from the magnet by superposition of many current loops distributed uniformly along the magnet height 35,36 with current I l = K M z , where z is the separation between loops (see Fig. 9). The number of current loops will be increased until it reaches to a point that with further increasing, the value of vector potential (or magnetic field) in the observation points remains unchanged. This process leads to a high accuracy for calculation of magnet vector potential and it is valid for any magnet shape with rotational symmetry. Next, we detail the calculation of the vector potential of the magnet, A M , in the magnet cylindrical coordinate system, where r ′ and z ′ are its radial and axial coordinates. Thanks to cylindrical symmetry, A M (r) = A M (r ′ , z ′ )e ′ ϕ , where e ′ ϕ is the unit vector in the angular direction. By dividing the cylinder into N thin loops, the vector potential follows A M (r ′ , z ′ ) = N i A l (r ′ , z ′ , r ′ i , z ′ i ) , where r i and z i are the radial and axial coordinates ofthe center of loop i and A l is the vector potential generated by a loop, being 35,37 Above, µ 0 is the permeability of vaccum, I l is the loop current, r ′ l and z ′ l are the current loop radius and axial position, respectively, k 2 = 4r ′ l r ′ /[(r ′ l + r ′ ) 2 + (z ′ − z ′ l ) 2 ] and K(k) and E(k) are complete elliptic integrals of the first and second kinds, respectively 38 .
The magnetic field generated by the magnet, B M , has both radial and axial components, B r and B z , respectively, in the cylindrical coordinate system of the magnet. Again, the calculation of B M in the magnet coordinate system is calculated by superposition of N thin loops, B r (r ′ , z ′ ) = N i B lr (r ′ , z ′ , r ′ i , z ′ i ) and B z (r ′ , z ′ ) = N i B lz (r ′ , z ′ , r ′ i , z ′ i ) where B lr and B lz are the components of the magnetic field generated by a thin loop where z ′′ = z ′ − z ′ l and k is defined as in (13).