⬇ Download Markdown  ·  ⬇ Download LaTeX  ·  ⬇ original  ·  view original

1

Efficient calculation of self magnetic field, self-force, and self-inductance for electromagnetic coils. II. Rectangular cross-section

Matt Landreman 1 , Siena Hurwitz 1 , and Thomas M Antonsen, Jr. 1

arXiv:2310.12087v1 [physics.app-ph] 18 Oct 2023

1 Institute for Research in Electronics and Applied Physics, University of Maryland, College Park, MD 21043 USA For designing high-field electromagnets, the Lorentz force on coils must be computed to ensure a support structure is feasible, and the inductance should be computed to evaluate the stored energy. Also, the magnetic field and its variation inside the conductor is of interest for computing stress and strain, and due to superconducting quench limits. For these force, inductance, energy, and internal field calculations, the coils cannot be naively approximated as infinitesimally thin filaments due to divergences when the source and evaluation points coincide, so more computationally demanding calculations are usually required, resolving the finite cross-section of the conductors. Here, we present a new alternative method that enables the internal magnetic field vector, self-force, and self-inductance to be computed rapidly and accurately within a 1D filament model. The method is applicable to coils for which the curve center-line can have general noncircular shape, as long as the conductor width is small compared to the radius of curvature. This paper extends a previous calculation for circular-cross-section conductors [Hurwitz et al, arXiv:2310.09313 (2023)] to consider the case of rectangular cross-section. The reduced model is derived by rigorous analysis of the singularity, regularizing the filament integrals such that they match the true high-dimensional integrals at high coil aspect ratio. The new filament model exactly recovers analytic results for a circular coil, and is shown to accurately reproduce full finite-cross-section calculations for a non-planar coil of a stellarator magnetic fusion device. Due to the efficiency of the model here, it is well suited for use inside design optimization.

Index Terms—Electromagnetic coil, self-inductance, magnetic energy, self-force, Lorentz force, stellarator, nuclear fusion, critical current

I. I NTRODUCTION

I N a configuration of electromagnets, important quantities include the Lorentz forces, inductances, and magnetic fields inside the conductors. The forces determine the necessary support structure, and the inductances determine the stored magnetic energy and dynamics. The magnetic field in the conductor determines local forces and their shear, and in superconductors the field also determines the critical current. For general steady currents, the field is determined via a threedimensional (3D) Biot-Savart integral,

B(r) = µ0 4π d3 ˜r J(˜r) × (r − ˜r) |r − ˜r|3 , (1)

where J is the current density. Similarly, the force per unit length dF/dℓ is determined by a 5D integral, given by integrating (1) over the conductor’s cross-sectional area, and the self-inductance L is given by a 6D integral:

L = µ0 4πI 2 d3 r d3 ˜r J(r) · J(˜r) |r − ˜r| , (2)

where I is the current, and the stored energy is W = (1/2)LI 2 . If one considers two thin coils A and B that are separated by more than the conductor width, the effect of one on the other can be computed by approximating them as infinitesimally thin filaments. The field and force at a point r A reduce to 1D integrals over coil B,

B(rA ) = µ0 I B 4π drB × (rA − rB ) |rA − rB |3 , (3) dF dℓ A = µ0 I A I B 4π tA × drB × (rA − rB ) |rA − rB |3 , (4)

Corresponding author: M. Landreman (email: [email protected]).

where t is the unit tangent vector. Similarly, the mutual inductance M can be computed with a 2D integral over both

inductance M can be computed with a 2D integral over both coils: M = µ0 4π drA · drB |rA − rB | . (5)

This reduction in dimensionality is extremely advantageous for simplicity and speed of computation. The self-field, selfforce, and self-inductance of a coil cannot be computed by setting B → A in (3)-(5) due to the singularity in the denominators when the source and evaluation points coincide. However, in this paper we will derive modified singularityfree 1D integrals for efficient calculation of the self-field and self-force, and a 2D integral for calculation of the selfinductance. The new formulas here are derived by rigorous analysis of the region where source and evaluation points nearly coincide, resolving the singularities in the filamentary integrals appropriately to have the correct asymptotic behavior. For one example documented in this paper, calculating the selfforce per unit length using the reduced model derived here is ∼ 18,000 times faster than using the original high-dimensional integral.

One motivating application for our work is the design of electromagnets for magnetic confinement fusion. These magnets have a large bore to enclose the plasma that they confine. Confinement of the plasma improves with field strength [1], so high fields are employed, and supporting the large Lorentz forces is a serious concern [2]. The reduction of the critical current for superconductivity due to the large field is also important. Figure 1 shows an illustrative example of a fusionrelevant magnet, one of the non-planar coils for the Helically Symmetric eXperiment (HSX) at the University of Wisconsin [3]. The shapes of fusion magnets, particularly for the stellarator concept, may be designed by optimization. It is therefore 2

Fig. 1. Coil from the HSX stellarator, shown from two angles. The cross- section dimensions a and b shown are those of the experiment. Vectors and color indicate the self-force, averaged over the cross-section. valuable to compute the Lorentz force and other quantities rapidly inside the optimization loop. The results in this paper are also applicable to the many other contexts in which shapes of electromagnets are optimized. These applications include magnetic resonance imaging [4], [5] and accelerator beam optics [6]. One might hope that the singularities in the filamentary integrals (3)-(5) with B → A could be resolved in numerical calculations by merely neglecting the singular quadrature point [7]. However, this approach leads to substantial errors in the self-field, self-force, and self-inductance. The results of such calculations diverge logarithmically as the number of grid points increases, rather than converging to any limit. Conse- quently, analytic results for geometries such as a circular coil cannot be reproduced. These problems with simply dropping the singular grid point are demonstrated in figure 3 of [8] for the force, and in figure 2 here for the field. Proper methods for calculating the self-field, self-force, and/or self-inductance using reduced models have been previ- ously proposed by several authors. For coils with rectangular cross-section composed of straight sections and circular arcs, formulas for the internal field have been derived analytically [9], [10], [11]. However reductions of these formulas to a 1D integral for the force or 2D integral for the inductance and energy do not appear to have been reported. Also, many coil shapes of interest are not composed of straight sections and circular arcs, such as the one in figure 1. In [12], the self-force was evaluated for a sheet current, which may be a reasonable approximation when there are a large number of adjacent filaments. Garren & Chen [13] evaluated the self- Fig. 2. The field strength on the inner edge of a circular coil (inset) with major radius 1 m, square cross-section of size 1 cm × 1 cm, and 100 kA current, computed three ways. Even with few grid points, our reduced model agrees with a high-fidelity calculation. A simplistic 1D filament model in which the singular point is ignored is highly inaccurate, not converging to any limit. force on current filaments of general shape, motivated by the effect of self-force on driving solar flare eruptions. In that work, a method to compute the force using only a 1D integral was proposed in which a specific interval of the integration domain was removed, avoiding the singularity. Dengler [14] proposed a similar method of removing a specific interval from 2D integrals for computing self-inductance. These last three works did not specifically consider the magnetic field inside the conducting filament. In the present paper, we derive computationally efficient re- duced integrals for the self-field, self-force and self-inductance with several advantages compared to these previous ap- proaches. All three of these quantities are computed in a unified framework. The coil geometry is allowed to be quite general, not necessarily composed of straight sections and circular arcs. We assume only that the cross-section is rect- angular, with dimensions that do not vary along the coil, and that the radius of curvature of the coil is large compared to the cross-section dimensions (disallowing sharp corners). No assumption is made about the rectangle’s orientation or its rotation along the coil. Unlike the approaches in [13], [14], our reduced integrals cover the full periodic domain along the coil. This means that the simplest possible quadrature grid, points uniformly spaced in the curve parameter, gives spectrally accurate integrals. Also, data such as the position vector and curvature from each quadrature point can naturally be re-used to compute the self-field and self-force at the other quadrature points. In contrast, the approaches in [13], [14] require custom quadrature grids for evaluating the self-field or self-force at different points along the coil. The reduced integrals derived here resemble the integrals one would get by assuming an infinitesmally thin coil (setting B → A in (3)-(5)), except that an extra term is present in the denominator, so there is no division by zero. The extra

term is negligible in most of the integration domain, but it becomes significant when the distance between source and evaluation points is small, comparable to the coil thickness. Physically, this term reflects the fact that when the source and evaluation points are this close, the finite thickness of the conductor cannot be neglected in the true high-dimensional integrals. In the reduced integrals, a factor of order unity is present multiplying the regularization term, and it is calculated rigorously to match the behavior of the true high-dimensional integrals.

This paper builds on our previous work in [8] for thin coils with circular cross-section. In the present work, similar techniques are applied to coils where the cross-section is rectangular. The case of rectangular cross-section is more complicated because the orientation of the cross-section must be specified, affecting some aspects of the calculations. Many details of the calculation here differ from and are more complicated than the circular-cross-section case, although the final results have similarities in form.

One important assumption in our method is that the current density is uniform across the conductor cross-section. At a coarse-grained level, this model of the current distribution is accurate for a coil composed of multiple turns in both the transverse directions, as each turn carries the same current. This model may not be highly accurate for high temperature superconducting (HTS) tapes, in which case the current may distribute itself non-uniformly in the transverse dimension along the tape [15], [16]. The uniform-current-density model also does not account for skin currents associated with diamagnetism that may be relevant for some superconductors. Nonetheless, the model here may still be useful for superconducting coils for initial estimations or use inside design optimization. A uniform current density approximation has been used for superconducting coils by other authors coils, e.g. section 2.2 of [17].

In the following section, the problem geometry is defined and coordinates are established. Then in section III, the main results are collected together for convenience. The methodology for deriving these results is explained in more depth in section IV. Sections V-VII give specific details of the selfinductance, self-field, and self-force calculations respectively. The reduced integrals have fine-scale structure in the integrand, and in section VIII we show how to remove this structure analytically so the integrals can be evaluated accurately on coarse grids. Finally the results are discussed in section IX, where we conclude.

II. P ROBLEM GEOMETRY

The conductor is considered to have a rectangular crosssection with side dimensions a and b at all points along its length. We define coordinates to parameterize the conductor volume. A periodic angle-like coordinate ϕ ∈ [0, 2π) will be used parallel to the conductor, while coordinates u, v ∈ [−1, 1] will be used in the transverse directions. The position vector r is expressed as

r(ϕ, u, v) = rc ( ϕ) + ua 2 p(ϕ) + vb 2 q(ϕ), (6)

3 n, but it rce and ickness. source ss of the ensional unity is lculated ensional for thin , similar ction is is more rc u = 1 u = -1 v = -1 v = 1 b a p q Fig. 3. Geometry definitions. The coil center-line rc is the curve shown in red. A section of the coil conductor surrounding the center-line is shown in

gray.

where r c (ϕ) is the position vector of the coil center-line, and p, q are unit vectors in the perpendicular plane aligned with the conductor cross-section. These definitions are illustrated in figure 3. The goal of our approach will be to express the inductance, self-force, and field strength in terms of integrals over ϕ involving r c .

The unit tangent of the center-line is t = |r ′ c | −1 r ′ c , where primes denote d/dϕ. The set {t, p, q} forms an orthonormal curve frame for r c . We consider this frame to be right-handed:

curve frame for rc . We consider this frame to be right-handed: t · p × q = 1. Differentiating p · p = 1, p · q = 0, etc., one

finds [18]

d dϕ t p q = drc dϕ 0 κ1 κ2 − κ1 0 κ3 −κ2 −κ3 0 t p q , (7) for some functions {κ1 ( ϕ), κ2 ( ϕ), κ3 ( ϕ)}. The standard cur-

vature κ(ϕ) of the center-line is the magnitude

pof the derivative of t with respect to arclength, hence κ = κ 21 + κ 22 . Noting the top row of (7) is also equal to κn|r ′ c |, where n is the unit normal, then κ 1 p + κ 2 q = κn. Crossing with t gives

κ1 q − κ2 p = κb, (8)

where b = t × n is the binormal. To compute κ 1 and κ 2 , one

where b = t × n is the binormal. To compute κ1 and κ2 , one can use κ1 = κn·p = κb·q, and κ2 = κn·q = −κb·p. Note

that since the {t, p, q} curve frame is defined to align with the rectangular conductor cross-section, the frame is smooth at any points or segments where κ = 0, in contrast to the Frenet frame.

From (6)-(7), the Jacobian of the (ϕ, u, v) coordinates is

√ g = ab |r′c | 4 1 − κ1 a u 2 − κ2 b v 2 . (9)

We consider the current density to be

J = I ab t, (10)

where I is the total current, i.e. uniform in magnitude throughout the conductor. Using the formula for the divergence in curvilinear coordinate systems, it can be verified that (10)

curvilinear coordina implies ∇ · J = 0.

4 where

B0 = µ0 I 4πab su , sv su s v [G (b(v − sv ), a(u − su ) ) q

B0 = µ0 I 4πab su , sv su s v [G (b(v − sv ), a(u − su ) ) q (17) −G (a(u − su ) , b(v − sv )) p] , G(x, y) = y tan− 1 x y + x 2 ln 1 + y2 x2 , (18) Bκ = µ0 I 64π su , sv su s v K (u − su , v − sv ) , (19) K (U, V ) = −2U V (κ1 q − κ2 p ) ln aU 2 b + bV 2 a (20) + (κ2 q − κ1 p ) aU 2 b + bV 2 a ln aU 2 b + bV 2 a + 4aU 2 κ2 p b tan− 1 bV aU − 4bV 2 κ1 q a tan− 1 aU bV , Bb = µ0 I κb 8π (4 + 2 ln 2 + ln δ) . (21)

Fig. 4. The dimensionless quantities k and δ defined in (12)-(13) appearing in the reduced models.

III. M AIN RESULTS

For convenience, our main results are collected in this section. These expressions are derived in later sections of the paper. The self-inductance L is given by

L = µ0 4π 2π 0 dϕ 2π 0 d˜ϕ r′c · ˜r′c |rc − ˜rc | 2 + δab , (11) where ˜r′c = drc ( ˜ ϕ )/d˜ϕ, δ = exp − 25 6 + k , (12) k = 4b 3a tan− 1 a b + 4a 3b tan− 1 b a + b2 6a2 ln b a + a2 6b2 ln a b 4 2 2 4

k = 4b 3a tan− 1 a b + 4a 3b tan− 1 b a + b2 6a2 ln b a + a2 6b2 ln a b − a4 − 6a2 b2 + b4 6a2 b2 ln a b + b a . (13)

These last two quantities are plotted as a function of the aspect ratio b/a in figure 4. For a cross-section that is square, k =

ratio b/a in figure 4. For a cross-section that is square, k = 2.5565 and δ = 0.19985. For a ≫ b, k = 7/6 + ln(a/b) and

δ = a/(be 3 ). As usual, the stored energy follows immediately from W = (1/2)LI 2 . For simplicity, the number of turns N is taken to be one; in the multi-turn case, one merely needs to multiply the inductance expressions throughout this paper by N 2 .

The self-force per unit length is

dF dℓ = It × Br eg , (14)

where

re Br eg (ϕ) = µ0 I 4π 2π 0 d˜ϕ ˜r′c × (rc − ˜rc ) |rc − ˜rc | 2 + δab 3 /2 . (15)

In (17), (19), and throughout this paper, the sums are over s u = ±1 and s v = ±1. The term B 0 is the field of an infinite straight wire with the same rectangular cross-section as our original general coil.

In the formula (16) for B, observe that B 0 ∼ µ 0 I/a, whereas B κ and B b are much smaller: ∼ µ 0 Iκ. Thus, for most applications, B κ and B b can be neglected. We nonetheless retain them here because for computing the cross-sectionaveraged Lorentz force, B 0 gives a vanishing contribution to leading order, so these higher order terms are required.

As anticipated in the introduction, the reduced integrals (11) and (15) have a regularizing term in the denominator, δab. Without this term, which prevents division by zero, the integrals would reduce to the expected integrals for infinitesmally thin coils. This new term is negligible compared to the adjacent term |r c − r̃ c | 2 in most of the integration region. However, where the separation between the source and evaluation points r̃ c and r c is comparable to the dimensions of the cross-section, the regularizing term δab becomes important.

Note that the inductance (11) and force (14) are independent of the angle by which the cross-section is oriented.

Finally, we will show that equivalent formulas for computing (11) and (15) which can be more accurate numerically are

L = µ0 4π 2π 0 dϕ |r′c | ln 64 δab |r′c | 2 (22) + µ0 4π 2π 0 dϕ 2π 0 d˜ϕ r′c · ˜r′c |rc − ˜rc | 2 + δab − |r′c | 2 2 − 2 cos ˜ϕ − ϕ |r′c | 2 + δab

and

The magnetic field inside the conductor is B = Br eg + B0 + Bκ + Bb , (16) and Br eg = µ0 I κb 8π −2 + ln 64 δab |r′c | 2 (23)

5

+ µ0 I 4π 2π 0 d˜ϕ ˜r′c × (rc − ˜rc ) |rc − ˜rc | 2 + δab 3 /2 −r′c × r′ ′ c 1 − cos ˜ϕ − ϕ 2 − 2 cos ˜ϕ − ϕ |r′c | 2 + δab 3 /2 .

In both (22) and (23), the last term cancels the near-singularity in (11) and (15), so the integrands are significantly smoother. Therefore the integrals can be evaluated to reasonable accuracy with many fewer quadrature points.

We now proceed to show how the formulas of this section are derived.

IV. M ETHODOLOGY

To rigorously justify the reduced formulas (11)-(21), we begin by establishing orderings. Let R indicate the large spatial

begin by establishing orderings. Let R indicate the large spatial scale of the coil center-line, with R ∼ |r′c | ∼ κ− 1 ∼ κ− 1 1 ∼ κ− 1 2 ∼ κ− 1 3 . Ordering a ∼ b, we consider an expansion in the

small parameter a/R ≪ 1. Moreover, an intermediate scale d is introduced, satisfying

a ≪ d ≪ R.

V. S ELF - INDUCTANCE AND STORED MAGNETIC ENERGY

The self-inductance L is defined by equating (1/2)LI 2 to the magnetic energy W = d3 rB2 /(2µ0 ) = d3 rA · J/2,

where A is the vector potential satisfying B = ∇ × A. Thus,

L = 1 I 2 d3 rA · J, where A = µ0 4π d3 ˜r ˜J |r − ˜r| . (25)

Tildes indicate quantities evaluated at the source location specified by the integration variables. As explained in the previous section, the total potential is written as a sum of near and far contributions, partitioned in ϕ. Using (10),

A = An ear + Af ar , An ear = µ0 I 4πab 1 −1 d˜u 1 −1 d˜v ϕ+ϕ0 ϕ−ϕ0 d˜ϕ √ ˜g˜t |r − ˜r| , Af ar = µ0 I 4πab 1 −1 d˜u 1 −1 d˜v 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ √ ˜g˜t |r − ˜r| .

(27)

(28)

Based on this decomposition, (25) gives the self-inductance as a sum of the two corresponding terms:

(24)

Associated with d, we define an angle ϕ 0 = d/R. In the original high-dimensional integrals for the inductance and field, we separate the integrals into a sum of two parts, depending on whether the source and evaluation ϕ locations are within ϕ 0 or not. This condition reflects whether the source and evaluation points are within a distance d or not. In the “far” term, where the difference in ϕ exceeds ϕ 0 , the finite cross-section of the conductor can be neglected, resulting in a substantial simplification of the integrals. In the “near” term, where the source and evaluation ϕ values are within ϕ 0 , the shape of the coil can be Taylor-expanded. The resulting simplified integrals can then be performed analytically. Now that the total quantity (inductance or field) has been expressed as a sum of two simplified terms, a reduced 1- or 2-dimensional integral is proposed as an ansatz. The proposed integral is similarly written as a sum of near and far contributions. The near parts of the high-fidelity and reduced models are shown to be equal. Similarly, the far parts of the high-fidelity and reduced models are equal. Thus, the high-fidelity and reduced models are equivalent in the limit a/R ≪ 1.

A related method of introducing an intermediate scale has been used previously for inductance problems. See e.g. page 123 of [19], or [14].

We now proceed to apply this procedure to the inductance, force, and field in the following sections. The force can be computed either directly, by integrating the force density J × B, or indirectly by applying the principle of virtual work to the inductance. The latter turns out to be convenient because it requires only lowest-order terms, whereas the direct calculation of the force requires finding higher-order terms in the field.

(26)

L = Ln ear + Lf ar , (29) Ln ear = 1 abI 1 −1 du 1 −1 dv 2π 0 dϕ √ gAn ear · t, (30) Lf ar = 1 abI 1 −1 du 1 −1 dv 2π 0 dϕ √ gAf ar · t. (31)

In these last two lines, note that the ϕ integrals without tildes still cover the full interval [0, 2π).

For the far contribution, |r − r̃| ≫ d, so all of the nonzerothickness terms can be dropped: r − r̃ ≈ r c − r ˜ c . Also the Jacobian (9) can be approximated as

√ g ≈ ab |r′c | 4 . (32) The integrals (28) and (31) then reduce to Af ar = µ0 I 4π 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ ˜r′c |rc − ˜rc | , (33) Lf ar = µ0 4π 2π 0 dϕ 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ r′c · ˜r′c |rc − ˜rc | . (34)

To evaluate A near and L near , quantities in (27) are Taylorexpanded about ϕ ′ ≈ ϕ. Details of the evaluation of A near and L near are given in appendix A. After a lengthy calculation, the result is

Ln ear = µ0 4π 2π 0 dϕ |r′c | 2 ln 2ϕ0 | r′c | √ ab + 25 6 − k , (35)

with k given by (13).

We now propose an ansatz of the form (11), where δ remains to be determined. Calling this expression L f il , we can write it without approximation as a similar sum of two terms:

Lf il = Lf il f ar + Lf il n ear , Lf il n ear = µ0 4π 2π 0 dϕ ϕ+ϕ0 ϕ−ϕ0 d˜ϕ r′c · ˜r′c |rc − ˜rc | 2 + δab ,

(36)

(37)

6 Julia implementation of the algorithm in [21]. Code used for all the numerical computations in this paper can be obtained at [22] and [23].

Figure 5 shows a comparison of the reduced model (11) with the high-fidelity 6D integral (2) for the HSX coil. To test the dependence on a and b in (12)-(13), the aspect ratio of the cross-section b/a is varied from 0.1 to 10 for three values of a. For comparison, the length of the coil is l = 2.05 m, comparable to a circle with radius l/(2π) = 0.327 m. The figure shows that the reduced model agrees extremely well with the high-fidelity model. Differences are imperceptible on the scale of the plot except at the largest values of a and b, where the expansion is not as accurate.

If the coil consists of N turns, with current I N = I/N per turn, the current density and hence magnetic energy W is unchanged. The inductance L N in this case is defined by W =

unchanged. The inductance LN in this case is defined by W = (1/2)LN I 2 N , so LN = N 2 L, hence the inductance formulas throughout the paper are multiplied by N 2 .

In the case of a circular coil, the inductance expressions here for both the high-fidelity and reduced models simplify to a classical formula. This result is shown in appendix D.

Fig. 5. The self-inductance of the HSX coil in figure 1, computed with the 2D reduced model (11) (dashed curves), closely agrees with highfidelity calculations using the 6D integral (2) (solid curves). The agreement is excellent across a wide range of values for the conductor cross-section dimensions a and b.

Lf il f ar = µ0 4π 2π 0 dϕ 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ r′c · ˜r′c |rc − ˜rc | 2 + δab . (38)

In the far term, the δab term in the denominator is negligible

In the far term, the δab term in the denominator is negligible compared to |rc − ˜rc | 2 . Then comparing with (34), we see f il f il 2

compared to |rc − ˜rc | 2 . Then comparing with (34), we see that Lf il f ar = Lf ar . In Lf il n ear , approximating |rc − ˜rc | 2 ≈ ( ˜ ϕ − ϕ)2 |r′c | 2 , the integral can be evaluated using (69), with the

result

Lf il n ear = µ0 4π 2π 0 dϕ |r′c | 2 ln 2ϕ0 |r′c | √ ab − ln δ . (39)

Comparing this expression to (35), we see that Lf il n ear = Ln ear

if we choose δ according to (12). Thus, the high-fidelity calculation of L by (25) rigorously agrees with (11) in the

calculation of limit a ≪ R.

We can now confirm the accuracy of the reduced model (11) by numerical comparison to the high-fidelity 6D integral (2). Throughout the paper, such comparisons between high-fidelity and reduced models will be made using the HSX coil displayed in figure 1. The rectangular cross-section is chosen to be oriented along the “coil centroid frame” described in [20].

RThis frame is defined in terms of the centroid C = (1/l) dℓ

Rr c , where dℓ = |r ′ c |dϕ is the incremental arclength and l = dℓ is the total length. The frame is defined by

p = w − (w · t)t |w − (w · t)t| , q = t × p, (40)

where w = r c − C. Thus, p is a vector pointing outward from the centroid, projected to be orthogonal to the tangent. Here, and in later sections of the paper, the multi-dimensional integrals are evaluated using the HCubature.jl package, a

VI. I NTERNAL FIELD

Now we derive the model for the magnetic field inside the conductor, (16)-(21). As explained in section IV, we will proceed by showing that the near and far parts of the high fidelity model match those of the reduced model. There are two options for computing B: either compute A to the requisite order and then take the curl, or directly evaluate the Biot-Savart law (1). Here we will follow the latter approach. (To find B to the desired accuracy from A, A would need to be computed to higher order than the calculation in Appendix A.)

Similar to section (V), we begin by writing the highfidelity expression (1) as a sum of near and far contributions,

fidelity expression (1) as a sum of near and far contributions, partitioned in ϕ: B = Bn ear + Bf ar , (41) Bn ear = µ0 I 4πab 1 −1 d˜u 1 −1 d˜v ϕ+ϕ0 ϕ−ϕ0 d˜ϕ √ ˜g˜t × (r − ˜r) |r − ˜r|3 , (42) 1 1 2π+ϕ−ϕ √

4πab −1 −1 ϕ−ϕ0 |r − ˜r| (42) Bf ar = µ0 I 4πab 1 −1 d˜u 1 −1 d˜v 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ √ ˜g˜t × (r − ˜r) |r − ˜r|3 . (43)

In the far term, the source and evaluation points are at least a distance ∼ d ≫ a apart, so the thickness of the conductor can be neglected. Hence we can approximate r ≈ r c and use (32), with the result

Bf ar = µ0 I 4π 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ |˜r′c | ˜t × (rc − ˜rc ) |rc − ˜rc | 3 . (44)

To evaluate B near , quantities in (42) are Taylor-expanded about ϕ̃ ≈ ϕ. Details are relegated to Appendix B. The result obtained is

Bn ear = B0 + Bκ + µ0 I κb 4π 1 + ln 4ϕ0 | r′c | √ ab (45)

7

Fig. 6. Magnetic field inside the conductor of the HSX coil, for the cross-section highlighted in red at left. The reduced model (16), requiring only a 1D integral, agrees well with the high fidelity model - the 3D integral (1). Both the field magnitude and one component (z) are shown.

where B 0 is the field of an infinite straight wire in (17), and

where B0 is the field Bκ is given by (19).

Now that the far and near behavior of the high-fidelity model is determined, we can guess a form for a reduced model with only a 1D integral. In contrast to the inductance and force, the near internal field depends on u and v, which would not be captured by a 1D integral alone like (11) or (14). Therefore we must add this local variation. A natural ansatz is

B = Br eg + Bu v , (46)

for some B uv (ϕ, u, v) to be determined, where B reg is given

for some Bu v (ϕ, u, v) to be determined, where Br eg is given by (15). We can write Br eg = Br eg n ear + Br eg f ar where Br eg n ear = µ0 I 4π ϕ+ϕ0 ϕ−ϕ0 d˜ϕ ˜r′c × (rc − ˜rc ) |rc − ˜rc | 2 + δab 3 /2 , (47) Br eg f ar = µ0 I 4π 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ ˜r′c × (rc − ˜rc ) |rc − ˜rc | 2 + δab 3 /2 . (48)

In Br eg f ar , the δab in the denominator can be neglected com- pared to |rc − ˜rc | 2 , so Br eg f ar ≈ Bf ar . To evaluate Br eg n ear , we ′

approximate r c − r̃ c ≈ (ϕ − ϕ̃)r ′ c in the denominator, and use

c ˜r′c × (rc − ˜rc ) ≈ 1 2 ( ˜ ϕ − ϕ)2 |r′c | 3 κb (49) in the numerator. Then the integral may be evaluated using (80), with the result Br eg n ear = µ0 I κb 8π −2 + 2 ln 2ϕ0 | r′c | √ ab − ln δ . (50) Equating (41) with (46), the filament model and high-fidelity calculation agree if one chooses Bu v = B0 + Bκ + µ0 I κb 8π (4 + 2 ln 2 + ln δ) . (51)

The reduced model derived here is demonstrated in figure

6. Here, the field inside the conductor is shown for the HSX

coil. In the the left panel of figure 6, the cross-section is shown as a red rectangle. The field magnitude and one field component are displayed in the right four panels, both for the reduced model (16), and for the original 3D Biot-Savart integral (1). The current is taken to be 150 kA. The highfidelity and reduced models agree well. For coils in which the cross-section dimensions are smaller relative to the coil curvature, the agreement improves further.

VII. S ELF - FORCE

The self-force can be computed from the self-inductance using the principle of virtual work, as explained in section 3 of [13]. A perturbation ∆r to the curve shape is considered, with the current varying to keep the flux Φ = LI fixed so there is no induced voltage. The change to the energy W = (1/2)LI 2

no induced voltage. The change to the energy W = (1/2)LI is then equated to the mechanical work dℓ ∆rc ·dF/dℓ, with

the result

dℓ ∆rc · dF dℓ = I 2 2 ∆L. (52)

This same result can also be derived by considering the shape perturbation ∆r c to be at fixed current. In this case, in addition to the mechanical work, the energy varies also due to the time integral of the IV power associated with the voltage V = dΦ/dt since the flux varies.

Thus, we consider a perturbation to r c in the inductance formula (11): r c → r c + ∆r c . For the terms in which the perturbation acts on the tilde position vector, we can swap the tilde and non-tilde variables. Integrating by parts to collect factors of ∆r c , we obtain

∆L = dℓ ∆rc · t × µ0 2π 2π 0 d˜ϕ ˜r′c × (rc − ˜rc ) |rc − ˜rc | 2 + δab 3 /2 . (53)

Equation (52) with (53) must hold for any perturbation ∆r c (ℓ), implying the integrands are equal. Thus, the force formula (14) with (15) follows.

8 We can also show the same result (14) can be obtained by integrating the Lorentz force density over the conductor crosssection. In this approach we proceed as discussed in section IV: the high fidelity calculation is divided into near and far regions, and the reduced model (14) with (15) is shown to give the same near and far contributions.

RBeginning with the volume-integrated Lorentz force F = d 3 rJ × B, applying

volume-integrated Lorentz force F = d rJ × B, applying d/dℓ = |r′c | − 1 d/dϕ gives the 5D integral for the high-fidelity

force per unit length:

dF dℓ = µ0 I 2 16πab 1 −1 du 1 −1 dv 1 −1 d˜u 1 −1 d˜v 2π 0 d˜ϕ (54) × 1 − uaκ1 2 − vbκ2 2 ˜g t × ˜ t × (r − ˜r) |r − ˜r|3 .

Taking advantage of our earlier decomposition (41)-(43), we have

dF dℓ = dF dℓ near + dF dℓ f ar (55)

where

dF dℓ near = I 4 1 −1 du 1 −1 dv 1 − uaκ1 2 − vbκ2 2 t × Bn ear , (56) dF dℓ f ar = I 4 1 −1 du 1 −1 dv 1 − uaκ1 2 − vbκ2 2 t × Bf ar = µ0 I 2 4π 2π+ϕ−ϕ0 ϕ+ϕ0 d˜ϕ|˜r′c | t × ˜ t × (rc − ˜rc ) |rc − ˜rc | 3 . (57)

Fig. 7. The self-force of the HSX coil, computed with the 1D reduced model (14) (dashed curves), closely agrees with high-fidelity calculations using the 5D integral (54) (solid curves). The z component of the force per unit length is displayed, averaging over the red cross-section shown in figure 6. The agreement is excellent across a wide range of values for the conductor crosssection dimensions a and b.

The integral (56) is done in appendix C, with the result

dF dℓ near = µ0 I 2 κn 4π k 2 − 13 12 − ln 2ϕ0 | r′c | √ ab , (58)

where k is again given by (13). Next, we consider an ansatz for the force of the form (14) with (15), but not assuming any particular value for δ. Dividing up the range of the ϕ integral as usual, the 1D filament model for the force per unit length is

dF dℓ f il = dF dℓ f il near + dF dℓ f il f ar (59)

where

dF dℓ f il near = It × Br eg n ear , dF dℓ f il f ar = It × Br eg f ar , (60) and Br eg n ear and Br eg f ar are given in (47) and (48). In Br eg f ar , the

f ar f ar δab in the denominator can be neglected compared to |rc − ˜rc | 2 , so we find (dF/dℓ)f il f ar = (dF/dℓ)f ar . Using (50) in (dF/dℓ)f il n ear , we find that (dF/dℓ)f il n ear = (dF/dℓ)n ear if we

choose δ in (15) according to (12). Thus, the near part of the high-fidelity model equals the near part of the reduced model, and the far part of the high-fidelity model equals the far part of the reduced model. As discussed in appendix C, this direct approach to justifying the reduced model required high-order terms in B near , whereas the virtual work approach required only the leading-order vector potential via the inductance.

The reduced 1D expression for the force agrees quantitatively with the high fidelity 5D integral (54) in numerical calculations. One demonstration of this agreement is shown Fig. 8. Another demonstration that the self-force of the HSX coil, computed with the 1D reduced model (14) (dashed curve), closely agrees with highfidelity calculations using the 5D integral (54) (solid curve). Here the true cross-section dimensions of the HSX coil are used. The z component of the force per unit length is displayed, averaging over the cross-section. The agreement is excellent at all points along the coil.

in figure 7. Here, the two models are evaluated for the HSX coil at the red cross-section highlighted in figure 6. Again, the current is taken to be 150 kA. A particular component (z) of the self-force is plotted in figure 7 for a range of crosssection dimensions a and b. Differences between the highfidelity and reduced models are imperceptible on the scale of 9

Fig. 9. Integrand for B reg without and with the singularity subtraction method of section VIII. The coil geometry is shown in the 3D inset, with the red cross-section indicating where B reg is being evaluated. Singularity subtraction removes the fine structure from the integrand, easing the resolution required for accurate integration.

Fig. 10. Convergence of the self-force (14) at one point along the coil (ϕ = 0) with respect to the number of quadrature points, without and with the singularity subtraction method of section VIII. With singularity subtraction, comparable accuracy can be obtained with many fewer quadrature points.

ϕ, the integrand in (11) is the plot except for the thickest cross-section, for which the a ≪ R expansion is less accurate.

A second demonstration of the accuracy of the reduced model is presented in figure 8. Here, the cross-section dimensions from figure (1) are used: a = 13 cm and b = 6 cm, matching the real HSX experiment. In figure 8, a component of the self-force is plotted at all points along the coil, computed using both the high-fidelity and reduced models. Even for this thick conductor, the agreement is excellent.

VIII. E FFICIENT QUADRATURE

Here we explain the derivation of (22) and (23). Although in sections IV-VII we have reduced the integrals to two or one dimension, the integrals (11) and (15) still require high resolution in the remaining dimension(s) due to fine structure in the integrands in the nearly singular region. This structure can be seen in figure 9. Here, the integrand for one component of B reg is plotted for the HSX coil at the point highlighted in figure 6. The dimensions a and b are reduced to 2 cm since the issue becomes more pronounced for thin cross-sections. The fine structure in the integrand is apparent near ϕ̃ = 0. Further computational efficiency can therefore be gained by manipulating the reduced integrals (11) and (15) to remove the fine-scale structure. We proceed by expanding the integrands of (11) and (15) about the nearly singular point. A term with the same behavior which can be integrated analytically is then subtracted and added. This calculation is similar to the one in [8].

Beginning with the inductance, and expanding r̃ c about ϕ̃ ≈

r′c · ˜r′c |rc − ˜rc | 2 + δab ≈ |r′c | 2 ˜ϕ − ϕ 2 |r′c | 2 + δab .

To make the quantity on the right periodic, we replace

˜ϕ − ϕ 2 → 2 − 2 cos ˜ϕ − ϕ

(61)

(62)

since the right and left expressions are the same to second order near ϕ̃ ≈ ϕ. When the result is subtracted and added to (11), one of the terms is integrated analytically using

2π 0 dϕ (2 − 2 cos ϕ + ∆)1 /2 = 4 √ 4+∆ K 4 4+∆ ≈ ln 64 ∆ (63)

where K is the complete elliptic integral, and ∆ =

where K is the complete elli δab/|r′c | 2 ≪ 1. Eq (22) follows.

Next, we proceed in a similar way for B reg . Expanding r̃ c about ϕ̃ ≈ ϕ, the integrand in (15) is

˜r′c × (rc − ˜rc ) |rc − ˜rc | 2 + δab 3 /2 ≈ 1 2 ˜ϕ − ϕ 2 r′c × r′ ′ c ˜ϕ − ϕ 2 |r′c | 2 + δab 3 /2 . (64)

Again, (62) is applied to make the quantity on the right of (64) periodic. When the result is subtracted and added to (15), one of the terms is integrated analytically using

2π 0 dϕ 2 − 2 cos ϕ (2 − 2 cos ϕ + ∆)3 /2 (65) = − 4 √ 4+∆ E 4 4+∆ − K 4 4+∆ ≈ −2 + ln 64 ∆

for ∆ ≪ 1, where E and K are complete elliptic integrals.

for ∆ ≪ 1, where E and K are complete Using r′c × r′ ′ c /|r′c | 3 = κb, we obtain (23).

10 The blue curve in figure 9 shows how this singularitysubtraction method removes the fine-scale structure from the integrand for B reg . We therefore expect the number of quadrature points required for accurate integration to be substantially reduced. This expectation is borne out in figure 10, which shows the convergence of the self-force with respect to the number of quadrature points. The coil parameters are the same as in figure 9, and uniformly spaced quadrature points in ϕ̃ with uniform weights are used. It is clear that the singularitysubtraction method greatly improves the convergence. Similar behavior is observed for the inductance (not shown).

For a comparison of the computational savings that can be gained by the reduced models in this paper, we have computed the self-force per unit length for the point along the coil in figure 6, both by integrating over the full finite-thickness cross-section (eq (54) and using our reduced 1D model (eq (14)-(15)). Adaptive quadrature is used for the 5D integral, while uniformly spaced quadrature points are used for the 1D integral. Both formulas are evaluated in Julia on one core of a 2020 M1 Macbook laptop. For three significant digits in the result, the high fidelity calculation takes 1.0 s, whereas the reduced model takes 57 µs, a factor of ∼ 18,000 speed-up.

IX. C ONCLUSION

In this work we have derived reduced expressions for efficiently calculating several quantities of interest for electromagnets: the self-inductance, magnetic energy, field inside the conductor, and self-force per unit length. These quantities normally must be computed by numerical calculations that resolve the cross-section of the conductor, involving high-dimensional integrals or equivalent partial differential equations. Here, by using analytic methods to reduce the problem, the quantities can be computed using only 1D integrals (for the self-field and self-force) or 2D integrals (for the self-inductance and energy). The main results have been collected in section III. The reduced models were formulated to rigorously match the original high-fidelity formulas in the limit that the conductor cross-section dimensions were small compared to the curvature of the coil center-line. Interestingly, the inductance, energy, and force-per-unit-length are found to be independent of the angle by which the rectangular cross-section is oriented. Throughout the paper, the reduced models have been compared to high-fidelity calculations in which the conductor crosssection is fully resolved. The agreement is excellent even for complicated coil shapes such as that in figure 1.

A key assumption of this work was that the current density was uniform across the cross-section. This approximation is accurate if there are many turns of the conductor in both dimensions of the cross-section. However the uniform-currentdensity approximation might not be sufficiently accurate in superconducting coils for some quantities of interest. For example, the distribution of current across a high-temperaturesuperconducting tape is typically nonuniform [15], [16], determined by physics not included here. The critical current of a superconducting coil is a property of practical importance, and precise computation of the critical current likely requires determination of this nonuniform current distribution by a method as in [17]. However, the internal field computed by (16) might at least provide an estimate for the critical current, as follows. First the coil current I is initialized to a plausible value; I will be updated later. Given B from (16), a formula for the local critical current density j c in terms of the local field (such as eq (6) of [16]) gives the local critical current density. This current density could be integrated over the cross-section to get a critical current at each ϕ, and minimization over ϕ would then give an estimate for the overall critical current I c :

Ic ≈ min ϕ ab 4 1 −1 du 1 −1 dv jc ( B(ϕ, u, v, I)). (66)

Root-finding in one variable could then be applied to solve for the coil current I such that I c = I, consistent with (66). This approach corresponds to the first iteration of the self-consistent numerical solution in [15], [16].

One other limitation of the work here is that sharp corners in the coil center-line are disallowed. It might be possible to address this case in future research.

In this work, the force per unit length was computed by averaging the Lorentz force over the cross-section. However other quantities related to Lorentz forces may also be of interest, including local force variation within the conductor crosssection, and stresses associated with gradients of the local force density. These additional quantities can be computed also using the approach here, since the full magnetic field vector within the conductor has been derived, (16).

A potential application for the results of this paper is in design optimization of electromagnets. One might wish to optimize the shapes of magnets for reduced peak forces, reduced magnetic energy, or reduced peak field. In optimization, it is necessary to evaluate the objective function many times, so speed of evaluation is a priority. Therefore the efficient formulas here will be advantageous.

A PPENDIX A D ETAILS OF THE NEAR VECTOR POTENTIAL AND

SELF - INDUCTANCE CALCULATION Now we evaluate A near . Again (32) can be used. In (27), r̃ ′ c and t̃ can be approximated by the corresponding quantities at the non-tilde (evaluation) point. In the denominator, for most of the integration domain, ϕ − ϕ̃ ∼ ϕ 0 = d/R, in which case |r − r̃| ≈ ϕ − ϕ̃ |r ′ c |. However if this expression is used in (27), the integral diverges near ϕ ≈ ϕ 0 . Evidently, more care is required in an inner layer where ϕ − ϕ̃ ∼ a/R. In this layer, we have to leading order

|r − ˜r|2 ≈ ϕ − ˜ϕ 2 |r′c | 2 + a2 4 (u − ˜u)2 + b2 4 (v − ˜v)2 . (67)

This expression is valid also in the rest of the near integration region where ϕ − ϕ̃ ∼ ϕ 0 . Thus,

An ear = tµ0 I 16π |r′c | 1 −1 d˜u 1 −1 d˜v ϕ+ϕ0 ϕ−ϕ0 d˜ϕ (68) × ˜ϕ − ϕ 2 |r′c | 2 + a2 4 (u − ˜u)2 + b2 4 (v − ˜v)2 − 1/2 .

11 The ϕ̃ integral is done using

ϕ0 −ϕ0 dϕ 1 (x + yϕ2 )1 /2 = 1√ y ln 1 + 2ϕ0 yϕ0 + √ y x + yϕ20 x ≈ 1√ y ln 4yϕ20 x ,

where (8) has been used. The terms on the right-hand side have a relative ratio of d 2 /(aR) which is not obviously large or small. Using (32) and |r̃ ′ c | ≈ |r ′ c | for the Jacobian, we have

Bn ear ≈ B( 0) + B( 1)

= √ y ln1 + x ≈ 1√ y ln 4yϕ20 x , (69) where x/(yϕ20 ) ∼ a2 /d2 ≪ 1 has been used. The remaining ˜u and ˜v integrals can be done using 1 −1 d˜u 1 −1 d˜v ln a b (u − ˜u)2 + b a (v − ˜v)2 (70) = −12 + 4 su , sv su s v F (u − su , v − sv ) where the sum is over su , sv ∈ {−1, 1} and F (U, V ) = aU 2 4b tan− 1 bV aU + bV 2 4a tan− 1 aU bV (71) + UV 4 ln aU 2 b + bV 2 a . To derive (70, the ˜u integral is straightforward, but the ˜v inte- gral that follows requires care because terms in the integrand are non-smooth at ˜v = v. These terms are correctly integrated using 1 −1 d˜v (v − ˜v) tan− 1 c v − ˜v (72) = c + sv sv 2 c2 tan− 1 v − sv c − (v − sv )2 tan− 1 c v − sv , for sv ∈ {−1, 1} and any c. Finally we obtain An ear = tµ0 I 4π 3 + ln 16ϕ20 |r′c | 2 ab (73) − su , sv su s v F (u − su , v − sv ) . With An ear now determined, the u and v integrals in (30) (with (32)) can be done using 1 −1 du 1 −1 dv su , sv su s v F (u − su , v − sv ) = 8 ln 2 − 14 3 + 4k. (74) The result is (35). AP PENDIX B DE TAILS OF THE NEAR MAGNETIC FIELD CALCULATION Here we give details of the calculation of Bn ear . We first evaluate the leading-order contributions to (42). In the denominator, (67) can be used, as in appendix A. In the numerator of (42), we have to leading order ˜t×(r−˜r) ≈ (u− ˜ u ) a 2 q−(v−˜v) b 2 p+ κ 2 ( ˜ ϕ −ϕ)2 |r′c | 2 b, (75) where B( 0) = µ0 I |r′c | 32π 1 −1 d˜u 1 −1 d˜v ϕ+ϕ0 ϕ−ϕ0 d˜ϕ (77) × (u − ˜u)aq − (v − ˜v)bp ϕ − ˜ϕ 2 |r′c | 2 + a2 4 (u − ˜u)2 + b2 4 (v − ˜v)2 3 /2 , B( 1) = µ0 I |r′c | 3 κb 32π 1 −1 d˜u 1 −1 d˜v ϕ+ϕ0 ϕ−ϕ0 d˜ϕ (78) × ˜ϕ − ϕ 2 ϕ − ˜ϕ 2 |r′c | 2 + a2 4 (u − ˜u)2 + b2 4 (v − ˜v)2 3 /2 . In (77), the ˜ϕ integral is first evaluated using ϕ0 −ϕ0 dϕ 1 (x + yϕ2 )3 /2 = 2ϕ0 x x + yϕ20 ≈ 2 x√ y , (79) for x ≪ yϕ20 . Then doing the ˜u and ˜v integrals, we find B( 0) = B0 with B0 given by (17)-(18). In (78), the ˜ϕ integral is first evaluated using ϕ0 −ϕ0 dϕ ϕ2 (x + yϕ2 )3 /2 = 2 y3 /2 ln x 2 − ϕ0 y x + yϕ20 − ln −√ yϕ0 + x + yϕ20 ≈ 1 y3 /2 −2 + 2 ln 2 − ln x y + 2 ln ϕ0 , (80) where the last line holds for x ≪ yϕ20 . The ˜u and ˜v integrals are performed using (70), and the result is B( 1) = µ0 I κb 8π (81) × 1 + 2 ln 4ϕ0 | r′c | √ ab − su , sv su s v F (u − su , v − sv ) , where the sum is over su , sv ∈ {−1, 1}, and F is given by (71). Observe that B( 0) is ∼ µ0 I /a, whereas B( 1) is for- mally smaller, with terms of size ∼ µ0 I /R and ∼ (µ0 I /R) ln(Rϕ0 / a). Also Bf ar has terms of size ∼ µ0 I /R and ∼ (µ0 I /R) ln ϕ0 , which can be seen by evaluating Bf ar for a circular coil in (96). Let us therefore compute all of the terms in B through order ∼ µ0 I /R. Indeed, we will need all terms of this order to compute the self-force by direct integration of J × B. The terms in B( 1) are already small so no high-order corrections are needed there, but corrections do need to be kept in the larger term B( 0) . Notice that any corrections to the integrand with an odd power of (ϕ − ˜ϕ)

12

vanish. Any corrections to the integrand that are small by a factor (ϕ − ˜ϕ)2 will give a contribution to B that is smaller by (a2 /R2 ) ln(ϕ0 R /a), which can be seen by comparing (79) with (80). Hence, we only need to keep corrections to the integrand that are smaller than the leading terms by a/R. There are no such corrections to the numerator of the Biot- Savart law, ˜t × (r − ˜r), but we do need to keep the corrections to the Jacobian (9), and use the following more accurate expression in the denominator: |r − ˜r|2 ≈ 1 − κ1 a 2 (u + ˜u) − κ2 b 2 (v + ˜v) ϕ − ˜ϕ 2 |r′c | 2 + a2 4 (u − ˜u)2 + b2 4 (v − ˜v)2 . (82) Then applying (79), ϕ+ϕ0 d˜ϕ √ ˜g 3 ≈ 2ab 1 + aκ1 4 (u − ˜u) + bκ2 4 (v − ˜v) (83) by replacing u ↔ v, a ↔ b, su ↔ sv , and κ1 ↔ κ2 . The necessary integral is 1 −1 du 1 −1 dv su , sv su s v a(u − su ) u b tan− 1 b(v − sv ) a(u − su ) + b(v − sv )2 2a tan− 1 a(u − su ) b(v − sv ) + (3u − su ) (v − sv ) 4 ln a(u − su ) 2 b + b(v − sv )2 a = − 2 3 + 8 ln 2 + 4k. (87) Applying this result in (86), we arrive at (58). AP PENDIX D CI RCULAR COILS In this section, we show how the formulas in this paper simplify when the coil center-line is a circle. In this case it

factor (ϕ − ϕ) will give a contribution to B that is smaller by (a2 /R2 ) ln(ϕ0 R /a), which can be seen by comparing (79) with (80). Hence, we only need to keep corrections to the integrand that are smaller than the leading terms by a/R. There are no such corrections to the numerator of the Biot- Savart law, ˜t × (r − ˜r), but we do need to keep the corrections to the Jacobian (9), and use the following more accurate expression in the denominator: |r − ˜r|2 ≈ 1 − κ1 a 2 (u + ˜u) − κ2 b 2 (v + ˜v) ϕ − ˜ϕ 2 |r′c | 2 + a2 4 (u − ˜u)2 + b2 4 (v − ˜v)2 . (82) Then applying (79), ϕ+ϕ0 ϕ−ϕ0 d˜ϕ √ ˜g |r − ˜r|3 ≈ 2ab 1 + aκ1 4 (u − ˜u) + bκ2 4 (v − ˜v) a2 (u − ˜u)2 + b2 (v − ˜v)2 (83) necessary integral is 1 −1 du 1 −1 dv su , sv su s v a(u − su ) u b tan− 1 b(v − sv ) a(u − su ) + b(v − sv )2 2a tan− 1 a(u − su ) b(v − sv ) + (3u − su ) (v − sv ) 4 ln a(u − su ) 2 b + b(v − sv )2 a = − 2 3 + 8 ln 2 + 4k. (87) Applying this result in (86), we arrive at (58). AP PENDIX D CI RCULAR COILS In this section, we show how the formulas in this paper simplify when the coil center-line is a circle. In this case it will be shown that the sum of near and far contributions to

We thus find the complete magnetic field through order µ 0 I/R is

Bn ear ≈ B( 0) + B( 1) + B( 0) 1 (84)

where the new term is

B( 0) 1 = µ0 I 16π 1 −1 d˜u 1 −1 d˜v (u − ˜u) aq − (v − ˜v) bp a2 (u − ˜u)2 + b2 (v − ˜v)2 × [(u − ˜u) aκ1 + (v − ˜v) bκ2 ] . (85)

Evaluating the ũ and ṽ integrals, and forming the sum (84) with (81), we arrive at (45).

A PPENDIX C D ETAILS OF THE NEAR SELF - FORCE CALCULATION In this appendix some details of the derivation of (58) are presented. The starting point is (56) with (45). Note that in B near in (45), the largest term is B 0 ∼ µ 0 I/a, whereas the other terms are ∼ µ 0 I/R. If the quantity in parentheses in (56) is approximated by 1, B 0 gives a contribution that integrates to zero by symmetry. Therefore, the small uaκ 1 and vbκ 2 terms in (56) must be kept for B 0 . However the quantity in parentheses in (56) can be approximated by 1 for the other terms in B near . Thus, we have

dF dℓ near = − µ0 I 2 κn 4π 1 + ln 4ϕ0 | r′c | √ ab (86) + I 8 1 −1 du 1 −1 dv [−(uaκ1 + vbκ2 ) t × B0 + 2t × Bκ ] .

In the B κ term, the contribution from the second line of (20) vanishes due to symmetry. In the p component, the vbκ 2 term vanishes due to symmetry, as does the uaκ 1 term in the q component. In the remaining integral, it is sufficient to evaluate just one of the p or q components, since the other follows each physical quantity becomes independent of the choice of ϕ 0 , as it should. It will also be shown that the self-inductance formula reduces to a classical expression.

For a circular coil that is concentric with the z axis, typically the cross-section would be oriented such that p = n = −e R and q = b = e z , where e R and e z are cylindrical unit vectors. Then a and b are the radial and axial extent of the conductor,

Then a and b are the radial and κ2 = 0, and κ1 = κ = 1/R0 .

To evaluate the far contributions to the quantities of interest, a necessary integral is

2π−ϕ0 ϕ0 c os χ dχ √ 1 − cos χ = − 4 √ 1 − cos ϕ0 ln tan ϕ0 4 sin ϕ0 2 + sin ϕ0 ≈ 2 √ 2 (−2 + 2 ln 2 − ln ϕ0 ) , (88)

where the last line is the leading behavior for ϕ 0 ≪ 1.

Using (88), the far contribution to the vector potential (33) can be evaluated as

Af ar = tµ0 I 2π (−2 + 2 ln 2 − ln ϕ0 ) . (89)

Adding (73), the total vector potential at leading order is

An ear = tµ0 I 4π −1 + ln 256R20 ab (90) − su , sv ∈{−1,1} su s v F (u − su , v − sv ) .

As desired, the ϕ 0 terms have cancelled.

The far contribution to the inductance (34) evaluates to

Lf ar = µ0 R 0 (−2 + 2 ln 2 − ln ϕ0 ) , (91)

where (88) has been used. Adding the near contribution (35), the total inductance becomes

L = µ0 R 0 ln 8R 0 √ ab + 1 12 − k 2 , (92)

13 with k given by (13). This classical result can be found as the leading terms on page 353 of [24], page 423 of [25], or page 137 of [26]. To show the equivalence,

tan− 1 1 x = − tan− 1 (x) + π 2 sign (x) (93)

can be used.

The same result can also be obtained from the reduced model (11). A necessary integral is

2π 0 dχ cos χ √ 1 − cos χ + ϵ (94) = 4 √ 2+ ϵ (1 + ϵ) K 2 2+ ϵ − (2 + ϵ) E 2 2+ ϵ ≈ √ 2 (−4 + 5 ln 2 − ln ϵ) ,

where K and E are the complete elliptic integrals, and the last line gives the leading terms for ϵ ≪ 1. Application of this integral to (11) in the case of a circle again gives (92).

For the magnetic field in a circular coil, evaluating the highfidelity 3D integral, (44) can be evaluated using

2π−ϕ0 ϕ0 dχ √ 1 − cos χ = 2 √ 2 ln cot ϕ0 4 ≈ 2 √ 2(2 ln 2 − ln ϕ0 ) (95)

as

Bf ar = µ0 I b 4πR0 (2 ln 2 − ln ϕ0 ) . (96)

Adding (45), the total field is

B = B0 + Bκ + µ0 I b 4πR0 1 + ln 16R0 √ ab , (97)

which is manifestly independent of ϕ 0 . The same result can be obtained from the filament model (16).This can be done by using (15) with (65) to obtain

Br eg = µ0 I ez 4πR0 ln 8R 0 √ ab + 13 12 − k 2 . (98)

Plugging this expression into (16), we again find (97)

Finally, the force per unit length for a circular coil can be evaluated for the high-fidelity model by substituting (96) into (57). Adding (58), we obtain

dF dℓ = µ0 I 2 eR 4πR0 ln 8R 0 √ ab + 13 12 − k 2 . (99)

The same result can be obtained from the reduced model (14) by applying (98).

A CKNOWLEDGMENT

We thank Robert Granetz, Stuart Hudson, Thomas Kruger, Jorrit Lion, Nicolo Riva, Felix Warmer, Yuhu Zhai, and Caoxiang Zhu for useful discussions. We are also grateful to Sarah Weiss and Wrick Sengupta for locating and translating [24]. This work was supported by a grant from the Simons Foundation (No. 560651, T. A.). This work was also supported by the U.S. Department of Energy under Contract DE-FG02- 93ER54197.

R EFERENCES

[1] D. Whyte, “Small, modular and economically attractive fusion enabled

by high temperature superconductors,” Philosophical Transactions of the Royal Society A, vol. 377, no. 2141, p. 20180354, 2019. [2] A. Panin, W. Biel, P. Mertens, F. Nunio, K. Sedlak, L. Zani, and M. Cole-

man, “Approaches to analyze structural issues of the european demo toroidal field coil system at an early design stage,” IEEE transactions on applied superconductivity, vol. 26, no. 4, pp. 1–5, 2016. [3] F. S. B. Anderson, A. F. Almagri, D. T. Anderson, P. G. Matthews, J. N.

Talmadge, and J. L. Shohet, “The helically symmetric experiment, (hsx) goals, design and status,” Fusion Tech., vol. 27, p. 273, 1995. [4] S. S. Hidalgo-Tobon, “Theory of gradient coil design methods for

magnetic resonance imaging,” Concepts in Magnetic Resonance Part A, vol. 36, no. 4, pp. 223–242, 2010. [5] X. Chen and M. Steckner, “Electromagnetic computation and modeling

in MRI,” Medical Physics, vol. 44, no. 3, pp. 1186–1203, 2017. [6] S. Russenschuck, Field computation for accelerator magnets: analytical

and numerical methods for electromagnetic design and optimization. John Wiley & Sons, 2011. [7] R. Gatto, F. Bombarda, S. Gabriellini, S. Murgo, and V. Zotta, “Time-

extended inductive tokamak discharges with differentially-tilted toroidal field coils,” Nucl. Fusion, vol. 63, no. 4, p. 046008, 2023. [8] S. Hurwitz, M. Landreman, and T. M. Antonsen, Jr., “Efficient cal-

culation of the self magnetic field, self-force, and self-inductance for electromagnetic coils,” arXiv:2310.09313, 2023. [9] S. J. Sackett, “Calculation of electromagnetic fields and forces in coil

systems of arbitrary geometry,” Lawrence Livermore National Laboratory, Tech. Rep., 1975. [10] L. Urankar, “Vector potential and magnetic field of current-carrying

finite arc segment in analytical form, Part III: Exact computation for rectangular cross section,” IEEE Transactions on Magnetics, vol. 18, no. 6, pp. 1860–1867, 1982. [11] J. Lion, F. Warmer, H. Wang, C. Beidler, S. Muldrew, and R. Wolf, “A

general stellarator version of the systems code process,” Nuclear Fusion, vol. 61, no. 12, p. 126021, 2021. [12] R. Robin and F. A. Volpe, “Minimization of magnetic forces on

stellarator coils,” Nucl. Fusion, vol. 62, no. 8, p. 086041, 2022. [13] D. A. Garren and J. Chen, “Lorentz self-forces on curved current loops,”

Phys. Plasmas, vol. 1, p. 3425, 1994. [14] R. Dengler, “Self inductance of a wire loop as a curve integral,”

Advanced Electromagnetics, vol. 5, p. 1, 2016. [15] L. Rostila, J. Lehtonen, and R. Mikkonen, “Self-field reduces critical

current density in thick YBCO layers,” Physica C: Superconductivity, vol. 451, pp. 66–70, 2007. [16] F. Gömöry and B. Klinčok, “Self-field critical current of a conductor

with an elliptical cross-section,” Superconductor Science and Technology, vol. 19, no. 8, p. 732, 2006. [17] N. Riva, R. Granetz, R. Vieira, A. Hubbard, A. Pfeiffer, P. Harris,

C. Chamberlain, A. Watterson, D. Anderson, R. Volberg et al., “De-

velopment of the first non-planar REBCO stellarator coil using VIPER cable,” Superconductor Science and Technology, 2023. [18] R. L. Bishop, “There is more than one way to frame a curve,” The

American Mathematical Monthly, vol. 82, p. 246, 1975. [19] L. D. Landau, E. Lifshitz, and L. Pitaevskii, Electrodynamics of contin-

uous media, 2nd ed., 1984, vol. 8. [20] L. Singh, T. Kruger, A. Bader, C. Zhu, S. Hudson, and D. Anderson,

“Optimization of finite-build stellarator coils,” J. Plasma Phys., vol. 86, p. 905860404, 2020. [21] A. C. Genz and A. A. Malik, “Remarks on algorithm 006: An adaptive

algorithm for numerical integration over an n-dimensional rectangular region,” Journal of Computational and Applied mathematics, vol. 6, p. 295, 1980. [22] M. Landreman, https:// github.com/ landreman/ CoilForces.jl, 2023. [23] ——, https:// www.doi.org/ 10.5281/ zenodo.10013332, 2023. [24] B. Weinstein, “Zur berechnung des potentials von rollen,” Annalen der

Physik, vol. 21, p. 329, 1884. [25] T. R. Lyle, “On the self-inductance of circular coils of rectangular

section,” Philosophical Transactions of the Royal Society of London. Series A, Containing Papers of a Mathematical or Physical Character, vol. 213, p. 421, 1914. [26] E. B. Rosa and F. W. Grover, Formulas and tables for the calculation

of mutual and self-inductance. US Government Printing Office, 1948, no. 169.