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

Computation of the Biot-Savart line integral with higher-order convergence using straight segments

Computation of the Biot-Savart line integral with higher-order convergence using straight segments

N. McGreivy, 1, 2 C. Zhu, 2 L. Gunderson, 3 and S. Hudson 2 1) Department of Astrophysical Sciences, Princeton University, Princeton NJ 08540, USA 2) Princeton Plasma Physics Laboratory, Princeton NJ 08543-0451, USA 3) Gatsby Computational Neuroscience Unit, University College London, WC1N 3AR, United Kingdom

(*Electronic mail: [email protected])

arXiv:2105.12522v3 [physics.comp-ph] 27 Jul 2021

(Dated: 28 July 2021)

One common approach to computing the magnetic field produced by a filamentary current-carrying coil is to approximate the coil as a series of straight segments. The Biot-Savart field from each straight segment is known analytically. However, if the endpoints of the straight segments are chosen to lie on the coil, then the accuracy of the Biot-Savart computation is generally only second-order in the number of endpoints. We propose a simple modification: shift each endpoint off the coil in the outwards normal direction by an amount proportional to the local curvature. With this modification, the Biot-Savart accuracy increases to fourth order and the numerical error is dramatically reduced for a given number of discretization points.

I. INTRODUCTION

Students of first-year college physics are taught the BiotSavart formula for the magnetic field produced by a static current-carrying wire loop:

B(r) = µ0 I 4π dr′ × (r − r′ ) |r − r′ |3 .

(a)

(1)

These students learn to use the Biot-Savart law to compute the magnetic field produced by a variety of symmetric configurations such as straight segments and circular coils. Meanwhile, as researchers in plasma physics, we numerically compute the Biot-Savart law in a variety of applications including field line tracing, 1 particle tracking, magnetohydrodynamic equilibrium codes, 2 and stellarator coil optimization codes. 3 But if the magnetic field cannot be calculated analytically, then equation 1 must be computed numerically. Ideally, we would have a method for computing the Biot-Savart line integral that would be simple enough to be taught to students of even firstyear physics yet accurate and efficient enough to be used in research applications.

There are three main approaches to computing the BiotSavart line integral given by equation 1. The first approach involves using a numerical quadrature and is usually highly accurate, but gives ∇ × B 6 = 0 and requires storing the coil parameterization r 0 (s) or a series of position vectors r i 0 and tangent vectors δ r i 0 . The second approach, described in the next paragraph, is inaccurate relative to the first approach but is conceptually simple and is widely used because it only requires storing a series of discrete points r i 0 and because, to machine precision, ∇ × B = 0 when the evaluation point does not lie on the curve. The third approach relies on the Fast Multipole Method 4 and is discussed in the related work section (section III). The purpose of this paper is to propose a new approach to evaluating the Biot-Savart line integral which is both highly accurate and conceptually simple while giving ∇ × B = 0, 5 thereby enjoying the benefits of both the first and second approaches.

The second approach to numerically computing equation 1 is to approximate the coil-carrying filament as a series of

(b)

FIG. 1: (a) Standard piecewise linear approach. Here the endpoints (red) lie on the filamentary curve (black). (b) Shifted piecewise linear approach. Here the endpoints are shifted in the outwards normal direction relative to the curve.

straight segments, then sum the magnetic field from each segment. As a convenient analytic expression is known for the magnetic field resulting from a straight segment 6 , a filamentary coil may be simply represented as an array of discrete endpoints. These endpoints lie on the coil and are typically equally spaced in distance or parameterization angle. This approach to computing equation 1, which we call the ‘standard piecewise linear approach’, is illustrated in figure 1a. In figure 1a, the filamentary curve is shown in black, the discrete endpoints are shown in red, and the piecewise linear segments are shown in blue.

Although the standard piecewise linear approach is convenient, simple, and flexible, it is not a particularly accurate way of computing the Biot-Savart line integral. One source of numerical inaccuracy can be understood from figure 1a: because the red endpoints lie on the curve, then the blue straight segments always lie locally inside the black curve. Because the straight segment always lies locally inside the curve, the field from the segment is biased and the result is that the accuracy of the Biot-Savart computation is only second-order in the number of endpoints.

We therefore propose a slightly different approach. Instead of choosing endpoints which lie on the filament, we choose Computation of the Biot-Savart line integral with higher-order convergence using straight segments

FIG. 2: A close-up of the geometry of a shifted linear segment. The filamentary coil is shown in black. The normal vector to the curve N is drawn in green, and the straight segment is drawn in blue. The black crosses represent the endpoints of the straight segments, and the magnitude of the shift is given by α. To lowest order, the area shaded in grey (positive) is equal to the area shaded in brown (negative) if α is given by equation 2.

endpoints which are slightly shifted in the outwards normal direction. This alternative approach, which we call the ‘shifted piecewise linear approach’, is sketched in figure 1b. Now the blue segments sometimes lie inside and sometimes lie outside the black curve, which eliminates the source of bias mentioned in the previous paragraph. The result is that, if the proper shift is chosen, then the accuracy of the Biot-Savart computation will be fourth order in the number of endpoints and the error will be dramatically reduced.

It remains to determine by what magnitude to shift the endpoints in the outwards normal direction. Figure 2 sketches a close-up of the geometry of a shifted linear segment, where the magnitude of the shift is labeled by α. In appendix A, we compare the Biot-Savart fields from a curved segment and shifted straight segment, and show that a shift of

α = κ|δ r′ |2 12

(2)

in the outwards normal direction cancels the second-order error in the Biot-Savart computation near the center of the coil. Similarly, in appendix B we show that the same shift minimizes the mean squared deviation between the curved segment and the straight segment. Here, the curvature κ ≡ 1/R, where R is the radius of curvature, δ r 0 ≡ (dr 0 /ds)δ s, s is a coordinate parametrizing the coil, δ s is the coordinate spacing between points, and the normal direction is defined by the Frenet-Serret formulas.

Inspired by the calculations in appendices A and B, we propose choosing α using equation 2. Although we have not discovered an analytic proof for the order of convergence with this shift, we do present careful numerical convergence studies (see section II) to test our method. Based on these numerical convergence experiments, we conclude that once there are enough endpoints such that the curve between endpoints is locally well-described by a parabola, then the convergence will be fourth-order for a closed curve with continuous curvature. To ensure fourth-order convergence of a piecewise differentiable curve or non-closed curve requires a slightly modified approach, which is described in appendix C.

In summary, the shifted piecewise linear approach is a highly accurate method of computing the magnetic field produced by a twice-differentiable filamentary coil. The piecewise linear approach results in a continuous curve and a piece2 wise constant representation of the tangent vector and therefore gives ∇ × B = 0 and only requires storing a series of discrete points r i 0 . This is in contrast to the numerical quadrature method, which gives ∇ × B 6 = 0 and which requires storing a set of discrete points r i 0 along with their associated tangent

a set of discre vectors δ r′i .

II. NUMERICAL TESTS

In this section, we perform a few numerical tests of our proposed approach, the so-called ‘shifted piecewise linear approach’. To do so, we compare the Biot-Savart accuracy of the shifted piecewise linear approach to the accuracy of the standard piecewise linear approach as a function of the number of discretization points N for three separate test problems: (1) an off-axis measurement from a circular coil, (2) an offaxis measurement from a non-planar modular stellarator-like coil, and (3) an on-axis measurement of the magnetic field from a D-shaped, tokamak-like coil. We measure the error as a function of the number of discretization points N, ε N , using the normalized ` 2 norm between the produced magnetic field for a given number of discretization points, B N (r), and the magnetic field as N → ∞, B ∞ (r):

εN = |BN (r) − B∞ ( r)| |B∞ ( r)| .

(3)

To avoid possible biases related to the particular choice of measurement point r, in each test problem we average the error ε N over a set of R = 100 measurement points {r j } Rj=1 where r j = r + δ j and δ j is a small random vector. In each test problem, the maximum deviation of ε N relative to the average is less than an order of magnitude. For test problem 1, B ∞ is known analytically and involves a computationally efficient elliptic integral. 7

Our results are shown in figure 3. In each case, ε N is smaller using the shifted piecewise linear approach by at least an order of magnitude compared to the standard piecewise linear approach for all N > 20. Additionally, while the accuracy of the standard piecewise linear approach is second-order, the accuracy of the shifted piecewise linear approach is fourth-order for each experiment. The shifted piecewise linear approach is significantly more accurate than the standard piecewise linear approach for a given number of discretization points.

We also have found that our results are qualitatively robust to the choice of measurement point r: it doesn’t matter whether r is on-axis, off-axis, near the coil or far from the coil, inside or outside the coil. In each case, the convergence of the shifted piecewise linear approach remains fourth-order.

III. RELATED WORK

Early work by Urankar has developed analytic expressions for the magnetic field produced by filamentary circular 7 and elliptic 8 arc segments, thin conic current sheets 9 , rectangular cross-section circular arcs with azimuthal 10 and arbitrary 11 Computation of the Biot-Savart line integral with higher-order convergence using straight segments 3

FIG. 3: The average error εN as a function of the number of discretization points N for (left) a circular coil (center) a stellarator-like, non-planar coil, and (right) a tokamak-like, D-shaped coil. The shifted piecewise linear approach is fourth-order in each experiment, while the standard piecewise linear approach is second-order. Code for these experiments, as well as additional experiments, can be found at https://github.com/nickmcgreivy/biot-savart-line-integral.

current densities, and polygonal cross-section circular arcs 12 . Hanson and Hirshman present new compact expressions for the analytic fields from a straight line segment. 6 Suh develops closed-form expressions for the magnetic field from volumetric elements with planar edges and constant and/or linear current density which involve line integrals over the boundary edges of the element. 13 Fontana examine a variety of methods for computing the Biot-Savart fields from a volumetric current-carrying coil, including (i) splitting the domain into tetrahedral regions and performing a quadrature on each domain, and (ii) numerically integrating a semianalytic one-dimensional integral based on Urankar’s analytic expressions. 14 Babic and Akyel generate analytic expressions for the magnetic field from trapezoidal prisms. 15 Nunes et al. approximate filamentary, surface current, and volumetric current densities as first-order or second-order finite elements, then use a quadrature to evaluate the current from each finite element and an adaptive method which reduces the required number of Gauss points. 16

Wobig et al. rewrite the Biot-Savart volume integral for a rectangular cross-section finite-build modular stellarator coil as a simpler, tractable integral over finite-build coil coordinates. 17 Although the integrals over the coil crosssection have no analytic solution, they can either be evaluated numerically or semi-numerically in a power series. McGreivy, Hudson, and Zhu have approximated the same integral using a multi-filament approximation and a quadrature-based approach to computing the Biot-Savart line integral. 18 The approach of using a numerical quadrature requires computing the tangent vector to the filamentary coil; this approach is often also significantly more accurate than the standard piecewise linear approach but is not examined in this paper to simplify the presentation.

A different approach to computing the Biot-Savart law uses the Fast Multipole Method (FMM) to improve the runtime scaling of the Biot-Savart law from O(NM) to O(N + M) where N is the number of sources and M is the number of evaluation points. This approach to computing the Biot-Savart law will significantly improve the runtime of a Biot-Savart calculation which requires a large number of evaluation points M, but a runtime comparison of these approaches as N and M vary is outside of the scope of this paper. Recently, van Nugteren and Deelen have developed an open-source library called Project Rat which implements this approach. 19

IV. CONCLUSION

We investigated the standard piecewise linear approach to computing the magnetic field from a filamentary coil. We found that this approach is only second-order accurate if the endpoints lie on the coil, because the straight segments are biased to always lie locally inside the coil. We propose a new approach, the shifted piecewise linear approach, where the endpoints are shifted in the outwards normal direction by an amount proportional to the local curvature. We find that a shift of κ|δ r 0 | 2 /12 minimizes the squared distance in real space between the filament and the straight segment, while also cancelling the second-order errors in the Biot-Savart field near the radius of curvature of the filament. We conduct simple numerical tests, and find that the shifted piecewise linear approach has higher-order convergence than the standard piecewise linear approach and is dramatically more accurate for a given number of discretization points.

There are two main conclusions. First, the shifted piecewise linear approach is superior to the standard piecewise linear approach, and should be used when possible. Second, when filamentary coils are stored in databases, the stored endpoints should not lie on the coil but rather should be shifted normally outwards.

Computation of the Biot-Savart line integral with higher-order convergence using straight segments 4 ACKNOWLEDGMENTS The magnetic field for the curved segment is given by This work was supported by the U.S. Department of Energy under Contract No. DE-AC02-09CH11466 through the Princeton Plasma Physics Laboratory.

B(r0 ) = −ˆz L/2 −L/2 R + L2 8R + x′ 2 2R x′ 2 + (R + L2 8R − x′ 2 2R )2 3 /2 dx′ .

(A5)

Again defining β ≡ x 0 /L, this can be written as DATA AVAILABILITY STATEMENT

The data that support the findings of this study are available from the corresponding author upon reasonable request. Appendix A: Derivation of Piecewise Linear Shift Here we show that by shifting each point on the fil- amentary coil in the outwards normal direction by an amount κ|δ r′ |2 /12, where κ is the local curvature and dr′ ≡ (dr′ /ds)δ s, cancels the second-order error in the Biot-Savart field near the center of the coil. To do this, we compare the Biot-Savart fields for the straight segment and curved section as shown in figure 2. With the origin defined as the midpoint between the two endpoints of the curved segment in figure 2, the straight segment can be written as r′ = x′ ˆx + α ˆy for x′ ∈ [−L/2(1 + α/R), L/2(1 + α/R)] while the curved segment can be written as r′ = x′ ˆx + (L2 /8R − x′ 2 /2R) ˆ y for x′ ∈ [−L/2, L/2]. This quadratic parametrization of the curved segment is equal, through second order, to that of any arbitrary smooth curve. We compute the magnetic field at the point r0 = (0, −R, 0) where R ≡ 1/κ. We compute the Biot-Savart fields for the straight segment and curved segment and work in the limits L/R ∼ α/L ∼ 𝒪(ε) and α/R ∼ 𝒪(ε 2 ), keeping all terms of up to second order. The magnetic field for the straight segment at r0 is given by the Biot-Savart law (equation 1) and equal to B(r0 ) = −(R + α)ˆz L/2(1+α/R) −L/2(1+α/R) dx′ x′ 2 + (R + α)2 3 /2 . (A1) Defining β ≡ x′ /L, this can be written as B(r0 ) = − L R2 (1 + α R )ˆz (1+ α R )/2 −(1+ α R )/2 dβ (1 + α R )2 + L2 R2 β 2 3 /2 . (A2) The integral can be performed analytically, it gives B(r0 ) = − L R2 ˆz L2 4R2 (1 + α R )2 + (1 + α R )2 . (A3) Expanding the denominator gives B(r0 ) = −ˆz L R2 1 − α R − L2 8R2 + 𝒪(ε 4 ). (A4) B(r0 ) = −ˆz L R2 1/2 −1/2 1 + L2 8R2 (1 + 4β 2 ) L2 R2 β 2 + (1 + L2 8R2 (1 − 4β 2 ))2 3 /2 dβ . (A6) Expanding the denominator, this gives B(r0 ) = − ˆz L R2 1/2 −1/2 1 + L2 8R2 (1 + 4β 2 ) − 3L2 2R2 β 2 − 3L2 8R2 (1 − 4β 2 ) dβ (A7) B(r0 ) = −ˆz L R2 1 − 5L2 24R2 + 𝒪(ε 4 ) (A8) Comparing equations A4 and A8, we see that set- ting α = L2 /12R gives the same second-order error in Bz at r0 = (0, −R, 0). Here, R ≡ 1/κ and to fourth order L2 = |δ r′ |2 ≡ |(dr′ /ds)δ s|2 . Appendix B: Minimization of Mean Squared Distance With the same parameterization as above, we have that the squared distance from the curved segment to the straight segment is (αx′ /R)2 + (L2 /8R − x′ 2 /2R − α)2 for x′ ∈ [−L/2, L/2]. Taking the average over the segment gives ⟨|∆r′ |2 ⟩ = 1/2 −1/2 L2 8R − α 2 + α 2 R2 − 1 R L2 8R − α L2 β 2 + 1 4R2 L4 β 4 dβ = L4 64R2 − L2 4R α + α 2 + 1 12 − L4 8R2 + L2 R α + L2 R2 α 2 + 1 80 L4 4R2 . (B1) The derivative with respect to α, d⟨|∆r′ |2 ⟩/dα = − L2 4R + 2α + L2 12R + L2 6R2 α, (B2) when set to zero, gives the optimal shift of the straight seg- ment: α = 1 12 L2 R 1 + 𝒪(ε 2 ) (B3)

Computation of the Biot-Savart line integral with higher-order convergence using straight segments Appendix C: Prescription for Piecewise Continuous Curvature 5 6 J. D. Hanson and S. P. Hirshman, “Compact expressions for the Biot-Savart fields of a filamentary segment,” Physics of Plasmas 9, 4410–4412 (2002).

If the coil curvature is only piecewise continuous, one should approximate each segment individually. Although we have found that closed curves with continuous curvature have fourth-order error convergence using the shifted piecewise linear approach, to obtain fourth-order convergence from a nonclosed segment requires a slightly modified method. We describe this method below.

Consider approximating the magnetic field coming from one such section, parameterized by s i ∈ [0, 1] for i = 0, 1, . . . , N. The procedure is the same as before, with the following modulations:

1. Set s 0 = 0 and s N = 1. These points will not be shifted. 2. For i ∈ 1, 2, . . . , N − 1, we used to place s i uniformly.

Instead, spread them out slightly, such that the spacing in s between (s 0 and s 1 ) and √ (s N−1 and s N ) are smaller than the others by a factor of 2.

3. Shift the points 1, 2, . . . , N − 1 using equation (2). For

all shifts, use the δ s between the interior points, not the shorter length of the endpoints.

1 T. S. Pedersen, M. Otte, S. Lazerson, P. Helander, S. Bozhenkov, C. Biedermann, T. Klinger, R. C. Wolf, and H. S. Bosch, “Confirmation of the topology of the Wendelstein 7-X magnetic field to better than 1:100,000,” Nature Communications 7, 13493 (2016). 2 S. Lazerson, J. Schmitt, C. Zhu, J. Breslau, and A. STELLOPT Developers, “STELLOPT,” [Computer Software] https://doi.org/10.11578/dc. 20180627.6 (2020). 3 C. Zhu, S. R. Hudson, Y. Song, and Y. Wan, “New method to design stellarator coils without the winding surface,” Nuclear Fusion 58, 016008 (2017). 4 L. Greengard and V. Rokhlin, “A new version of the fast multipole method for the laplace equation in three dimensions,” Acta Numerica 6, 229–269 (1997). 5 It is well known that the Biot-Savart integral has zero curl in regions where J = 0, as predicted by Ampere’s law. However, the integrand of the BiotSavart law has non-zero curl even where J = 0. Therefore, in contrast with many other methods (e.g., quadrature, discontinuous piecewise linear), the explicit construction of a closed filamentary curve ensures that the resulting magnetic field is curl-free to machine precision.

7 L. Urankar, “Vector potential and magnetic field of current-carrying finite arc segment in analytical form, part i: Filament approximation,” IEEE Transactions on Magnetics 16, 1283–1288 (1980). 8 L. Urankar, “Vector potential and magnetic field of current-carrying finite elliptic arc segment in analytical form,” Zeitschrift für Naturforschung A 40, 1069 – 1074 (1985). 9 L. Urankar, “Vector potential and magnetic field of current-carying finite arc segment in analytical form, part ii: Thin sheet approximation,” IEEE Transactions on Magnetics 18, 911–917 (1982). 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 18, 1860–1867 (1982). 11 L. Urankar, “Vector potential and magnetic field of current-carrying finite

arc segment in analytical form, part iv: General three-dimensional current density,” IEEE Transactions on Magnetics 20, 2145–2150 (1984). 12 L. Urankar, “Vector potential and magnetic field of current-carrying circu-

lar finite arc segment in analytical form. v. polygon cross section,” IEEE Transactions on Magnetics 26, 1171–1180 (1990). 13 J.-C. Suh, “The evaluation of the biot–savart integral,” Journal of Engineer-

ing Mathematics 37, 375–395 (2000). 14 M. Fontana, “Integration methods for the calculation of the magnetostatic

field due to coils,” Preprint 2001-7, Chalmers Finite Element Center, University of Chalmers, Göteborg (2001). 15 S. I. Babic and C. Akyel, “An improvement in the calculation

of the magnetic field for an arbitrary geometry coil with rectangular cross section,” International Journal of Numerical Modelling: Electronic Networks, Devices and Fields 18, 493–504 (2005), https://onlinelibrary.wiley.com/doi/pdf/10.1002/jnm.594. 16 A. Nunes, O. Chadebec, P. Kuo-Peng, P. Dular, and B. Cucco, “Contribu-

tion on the source field calculation through the biot-savart equation using curvilinear elements and an adaptive process,” Progress In Electromagnetics Research C 88, 145–161 (2018). 17 H. Wobig, E. Harmeyer, F. Herrnegger, and J. Kisslinger, “The magnetic

field of modular coils,” IPP Report (2009). 18 N. McGreivy, S. R. Hudson, and C. Zhu, “Optimized finite-build stel-

larator coils using automatic differentiation,” (2020), arXiv:2009.00196 [physics.plasm-ph]. 19 J. van Nugteren and N. Deelen, “Project rat,” https://gitlab.com/

Project-Rat (2020). 20 S. L. Gralnick and F. H. Tenney, “Analytic solutions for constant-

tension coil shapes,” Journal of Applied Physics 47, 2710–2715 (1976), https://doi.org/10.1063/1.322993. 21 H. Lei, L.-Z. Wang, and Z.-N. Wu, “Integral analysis of a magnetic field

for an arbitrary geometry coil with rectangular cross section,” IEEE Transactions on Magnetics 38, 3589–3593 (2002). 22 D. Wilton, S. Rao, A. Glisson, D. Schaubert, O. Al-Bundak, and C. Butler,

“Potential integrals for uniform and linear source distributions on polygonal and polyhedral domains,” IEEE Transactions on Antennas and Propagation 32, 276–281 (1984).