The finite build coil computation would be be an implementation of these papers by Hurtwitz et. al and Landreman et. al:

• Circular cross section finite build paper
• Rectangular cross section finite build paper

An existing Julia implementation of the finite build rectangular cross section method is here.

• if we do a new _FiniteBuildCoil ABC, this will need a new optimizable param for geom which is r for circular and L W. otherwise the new subclasses need it defined
• i.e. FourierPlanarCoil -> FourierPlanarFiniteBuildCoil
  • maybe try doing FourierPlanarFiniteBuildCoil( _FiniteBuildCoil, FourierPlanarCoil) and see if it works
• basically just copy the way we have compute_magnetic_field setup but rename compute_self_force etc and define new formulas in coils.py as needed
  • need a cross_section="circular" or "rectangular" argument

@ryanwu4 @ddudt We should have a solution for the twist angle as well, so we don't implement the single thing for planar coils

Regarding twist angle, I think the approach in Singh et. al 2020 has been used in the original Julia implementation for finding a coordinate frame for this.

Specifically, Section 3.1 of that paper would define the coil centroid coordinate frame at each coil point, and we can initialize the FiniteBuildCoil object with a $\alpha_n$ vector with corresponding Fourier modes such that the rotation vector $\alpha$ at each $\phi$ represents the rotation of the winding pack cross section within that coordinate frame at the corresponding value of $\phi$ (Section 3.3).

Another thought: how useful are circular cross sections in practice? Do people use them? Would just rectangular be enough? (I would assume rectangular is closer to reality?)

I'm starting with implementing just the rectangular cross section coils, but I think I'll leave the functionality open for circular cross sections since the formulas have been derived there anyways.

Notes from meeting with Paul:

• be sure to have ability to consider the winding angle of a curve/coil without also needing the winding pack information, to allow for optimizing the filamentary coil's winding frame without needing to consider the finite build effects
• Rectangular is 99% of required cases, probably don't initially need pack XS varying along coil
• need separate, large weights on binormal/torsional strain to ensure it is respected during optimization

Regarding the first and third points, it almost sounds like the winding angle feature would be implemented separately from the finite build coil feature in that case? Should that be something else that the FiniteBuildCoil inherits from?

Yeah, I think the order would be something like

Coil(MagneticField): # base class
FramedCoil(Coil): # coil plus rotating reference frame, but no width/height - still filamentary
FiniteBuildCoil(FramedCoil): # as above but with width/height of winding pack

also here's some ideas from last week:

"""
A winding pack is defined by a width (w) and height (h), centered around a 
filamentary curve (x), carrying a total current Itot

           ^ e_2
           |
           |

  x10      w       x11
    |---------------|
    | *  *  *  *  * |
  h | *  *  x  *  * |  --> e_1
    | *  *  *  *  * |
    |---------------|
  x00              x01
  
The x coordinate of the winding pack is the centerline, parameterized by a fourier 
series or spline
the xij coordinates define the corners of the winding pack

x00 = x - h/2 e_2 - w/2 e_1
x10 = x + h/2 e_2 - w/2 e_1
x01 = x - h/2 e_2 + w/2 e_1
x11 = x + h/2 e_2 + w/2 e_1

where e_1 and e_2 are the normal and binormal basis vectors in the rotating 
reference frame about the centerline (eg frenet seret, or centroid frame 
is usually smoother). Using the corner coordinates we can easily 
discretize the finite build coil into polytopes (rectangular sections 
or tetrahedra).

For calculating magnetic fields, it is discretized into individual 
filaments (*) - nw in the width direction and nh in the height 
direction (in the example above, nw=5, nh=3), and integration over 
the cross section is done with gaussian quadrature (maybe lobatto 
if we want to include endpoints?). 

eg

y1, w1 = leggauss(nw)
y2, w2 = leggauss(nh)

B = biot_savart(x=(x + y2*h/2 e_2 + y1*w/2 e_1), 
                I=(I_tot*w1*w2))

(I feel like for uniform J in rectangular cross section there is an analytic answer for B?...)

"""