In this article, we present a method to design a coil producing an arbitrarily shaped magnetic field by restricting the path of the coil's wires to a regular grid. The solution is then found by a simple least squares minimum. We discuss practical applications, in particular, in the active magnetic field stabilization system of the neutron electric dipole moment experiment at the Paul Scherrer Institute in Villigen, Switzerland. We also publish the software implementation of the method.

## I. INTRODUCTION

How to design a coil, or more generally, an arrangement of coils producing a desired magnetic field? In its simplest form, this is a textbook problem (e.g., ex. 6.55 and 6.62 in Ref. 1). Yet, in a general setting, it is surprisingly hard and the solutions, how the wires making up the coils should be laid, are complicated. The most widespread application of high-performance coils is Magnetic Resonance Imaging (MRI), where gradient coils give the possibility to produce spatial images. As early as the 1980's elaborate methods of MRI coil design had been developed. They range from optimizing positions of discrete windings, where use is made of symmetries specific to MRI, to analytical methods yielding surface current density, which is then discretized. A general overview can be found in Ref. 2. Another field known for complex, precise coils is plasma confinement, in particular, stellarators.^{3} There analytical solutions for the surface current density also are of use.

Here, we present a new method that may not be competitive in terms of precision, but is distinct in its simplicity, also when it comes to construction of its designs. It relies on an algebraic representation of the problem, where coil design is simplified to a simple linear least squares problem. In our method, the coils are restricted to a user-defined mesh, making it easy to deal with spatial constraints.

We base our discussion on textbook linear algebra techniques, notably solving an over-determined system of linear equations, thoroughly discussed, e.g., in Ref. 4. The main physics problem, calculating the magnetic field of coils composed of straight wire segments, is briefly discussed here. More in-depth discussions can be found, for example, in Ref. 5. Furthermore, we publish an easily accessible software implementation of the problems we discuss, including examples.^{6}

The method has originally been developed to design coils of the active magnetic field stabilization system of the neutron electric dipole moment (nEDM) experiment at the Paul Scherrer Institute in Villigen, Switzerland. This application, described in detail in Ref. 7, and discussed below, requires large (6–8 m side length) coils. In the presented method, the coil system is designed on a predefined grid. This makes the construction of complicated coils feasible, despite the size.

We begin with a description of our model for a restricted 2-dimensional case and generalize it to three dimensions. We then show how the model is used to design a coil, based on an example. In addition, we discuss possibilities for simplifying the solution. Another section is devoted to practical considerations, significant for the eventual construction. Finally, we analyze the design method in the particular case of the magnetic field stabilization system of the nEDM experiment at the Paul Scherrer Institute.

## II. COILS AS A LINEAR SPACE

Consider all possible coils that can be constructed by laying a wire on a surface of a square. The possibilities are endless. Speaking more precisely, since the wires may be shifted by arbitrarily small distances as they overlap and cross, the problem has inherently an infinite number of degrees of freedom. We present an algebraic representation that reduces the number of degrees of freedom to just a few.

We start with a straight, finite wire segment spanned between points $x1$ and $x2$ (represented by vectors in an arbitrary coordinate system) and carrying current *I*, as depicted in Fig. 1. To calculate the magnetic field it produces in the point ** p**, we use the Biot–Savart law. We consider the vector normal to the wire through the point

*p*where ** n** is a unit vector in the direction $x2\u2212x1$. The magnitude of the magnetic field in point

**is then**

*p*^{5}

where the angles *α _{i}* are not directed

and *s* is +1 if $\rho $ points onto the wire segment (between points $x1$ and $x2$) and −1 otherwise

The direction of the field is given by the right-hand rule

This formulation is independent of the coordinate system (coordinate-system dependent solutions can be found, e.g., in Ref. 8).

Let us imagine a four wire segments making up a square loop—a coil. It produces a certain magnetic field in the entire space **B**(**x**), given, according to the superposition principle, by a sum of the fields produced by each segment of the coil. By changing the current in the coil, we alter only one parameter of the magnetic field—the magnitude, but not its shape. It can therefore be said that one coil spans a one-dimensional space of magnetic fields it can produce. Adding a second, different coil creates a system spanning a two–dimensional space of fields, as the magnetic field is additive. Going a step further, four square coils tiled to form a larger square form a four-dimensional space, as shown in Fig. 2. Any coil restricted to the 2 × 2 grid can be represented in the basis of the four tile-coils.

The range of magnetic field reachable by coils restricted to a grid is a subset of all possible fields that can be created with coils constructed on the square's surface. The size of the subset is controlled by *N*, the number of tile–coils forming the grid. In this system, a coil is fully described by a vector of *N* currents, one in each of the tile–coils, denoted by $I$. The problem of coil design is thereby simplified to finding a vector $I$.

Generalization to a cube is simple, a cube being made up of six square faces. Interestingly, for the assembly in the three-dimensional space one degree of freedom is lost. Figure 3 illustrates, in the simplest case *N* = 6, a configuration in which finite currents in all six coils cancel and no magnetic field is produced. Such a combination of currents can be added to any solution with no effect on the produced field. Effectively, the space of the fields they can produce has dimension five (i.e., *N* – 1). In other words, the mapping of $I$ onto fields **B**(**x**) has in this case a one-dimensional kernel. This fact is of importance when it comes to numerically solving the system.

This is the foundation of the method. We restrict our consideration to a grid on a cuboid, but in return we can fully describe all coils in the restricted space by a vector of *N* numbers.

## III. COIL DESIGN

In the problem of coil design, one wants to create a coil, or an arrangement of coils, which best approximates a given field in a certain volume, which we will call *the volume of interest*. Rather than considering the whole volume, we pick an ensemble of *m* points of interest on its surface (the surface is sufficient because ∇⋅ **B** = 0). Hence, we look at the magnetic field **B**(**x**) only at these points and gather the values **B**(**x**_{i}) for *i* = 1,…,*m* into a vector of dimension 3*m* (*B _{x}*,

*B*, and

_{y}*B*in each point), which we shall denote $B$.

_{z}As mentioned before, the magnetic field produced by a coil at any given point in space is proportional to the current in this coil. With many coils present it is a linear combination of the currents of all coils in the system. In absence of an external magnetic field, the system of *N* tiles and *m* points of interest is thus described by a simple linear equation

where $M\u2208\mathbb{R}3m\xd7\mathbb{R}N$ is a matrix of proportionality constants. For example, the element $M(5,2)$ is the proportionality constant between the current in the second of *N* coils and the magnetic field in the *y* direction in the second of *m* points of interest, *B _{y}*(

**x**

_{2}). The matrix $M$ can be calculated analytically using the Biot-Savart law.

Equation (6), for 3*m* > *N* – 1, is an over-determined system of linear equations, $I$ being the vector of unknowns. We look for the optimal least-squares solution $I0$ to produce a **B**_{0}(**x**) in the volume of interest

The optimal solution can be calculated with the normal equation^{4}

but the problem is typically solved numerically.^{9} The majority of numerical software packages use the QR decomposition (a product of an orthogonal and upper-triangular matrix) of the matrix $M$, which is more numerically stable when compared to the normal equation.

Depending on the properties of $M$, the optimum may be multidimensional. In particular, as already mentioned, an arrangement of coils on a cube has a one-dimensional kernel, which will always cause the optimum to be at least one-dimensional. In these cases, we will call $I0$ the unique least-norm solution, which minimizes the total current in the system. $I0$ is the vector of the optimal currents in the tile arrangement of coils for approximating **B**_{0}(**x**) in the volume of interest.

Let us look at an example of a coil design on a unit cube with the number of tiles *N* = 6 × (3 × 3) (see Fig. 4). As the volume of interest we pick a cube, centered with the unit one, with side length 0.75 (with a regular mesh of 10 × 10 points on each face, a total of *m* = 488 points of interest). For the sake of simplicity, we design a coil for a homogeneous field along an axis of the cube. The solution of Eq. (7), $I0$, directly gives the currents in each tile, which are graphically depicted in Fig. 4. Note that many currents almost cancel each other, in particular, those along horizontal edges. The magnetic field produced by the solution is shown in Fig. 5(a), as a horizontal cut along the central plane. Contours show the relative deviation from the homogeneous field. Inside the volume of interest, depicted by a dashed line, the design goal of a homogeneous field is reproduced with a few per cent accuracy. The solution, and thus the contours too, depend on the choice of the volume of interest. In general, the further away the volume of interest is from the coils, the better the accuracy. If the side length of the volume of interest is decreased to 0.5, the accuracy improves to 1%, as shown in Fig. 5(b). Note that the optimal solution, and thereby the shape of the precision contours, change. Naturally, the accuracy of the field reproduction can also be improved by increasing the number of tiles.

## IV. SIMPLIFICATION OF THE TILE SYSTEM

The tile system may find an interesting practical application. Once independently controllable tiles have been built, it can be used to produce an arbitrary field. However, building many independently driven coils is a high price to pay if one wants to produce only a simple field. Additionally, note that each edge is shared between two tiles, and the effective current is the sum of two. They may add either constructively or destructively. If the given solution is dominated by subtraction of large currents, a lot of power is unnecessarily dissipated in the system. It turns out that both problems can be solved by simplifying the tile solution.

One starts by adding the currents of the adjacent tiles and assigning the sum to each common edge. The result is a complicated net of currents [upper left corner of Fig. 6(a)]. Still, each node fulfills Kirchhoff's laws. The net can then be decomposed into simple current loops by following the algorithm: First find in the net the loop with the highest current. In the example, it is either of the “597” loops on the front and back faces. This loop will make the first one in the simplified solution [depicted in Fig. 6(d), together with the next three loops]. Then subtract from the net the current of the first loop along its edges [the net that remains after subtracting the first four loops is depicted in Fig. 6(c)]. Finally, continue to find the loop with the highest current in the modified net, which will give the next loop and repeat until the current net is empty. The net remaining after eight loops are found is depicted in the bottom row of Fig. 6(e), next to the first eight loops [Fig. 6(f)]. The final simplified solution is shown in Fig. 7. The currents in the simplified coil system are much smaller, the highest being 597 instead of 1000 and they always add constructively. Also the number of separate loops is decreased from 42 to 10. Still, the total current along each edge of a tile is exactly the same as in the tile configuration.

Here we conclude here our method of coil design. The simplified arrangement of coils is the optimal one, given the grid restriction, for approximating the magnetic field in the volume of interest. We turn now to a consideration of practical aspects relevant for constructing designs with our method.

## V. PRACTICAL CONSIDERATIONS

The primary practical advantage of our coil design method is that the coils are constrained to a predefined grid. This is contrary to other methods of coil design, where the position of the wires is the output of the procedure.^{2,3} Our approach may prove useful in applications with spatial constraints. Typically, coils need to be incorporated into a setup in which other components penetrate the surface on which the wires are laid. In our method, it is possible to simply define the grid so that no collisions occur. Although the simple examples presented above used regular grids, we have not used symmetries to solve the problem. When many coils are designed and built, for instance, to produce homogeneous magnetic fields in each of the three dimensions, they can all share the same grid. The grid can, for example, be constructed out of cable channels into which the wires are laid.

A limitation associated with the finite size of the channels is the strength of the magnetic field that can be created, which, for given available power, is limited by the thickness of the wire. At the same time, the finite size of the cable channels can be neglected in the calculations only as long as it is small compared to the distance between the coils and the volume of interest. Using enameled wire, rather than standard, PVC-insulated cable, can reduce the overall thickness.

In our solution to produce the desired field, one still needs a system of several coils, even in the simplified solution. The more complicated the desired field and the more tiles used, the larger the set of different currents that are needed across the individual loops, which quickly becomes impractical. There are several ways to tackle the problem.

The first way is to use only one current and adjust the number of windings. In the example, when one decides for 60 as the maximum, then the current is round(597/60) = 10. The 597, 360, 25, 13, and 5 would be created with 60, 36, 3, 1, and 1 windings, respectively. A discretization error of 10/597 = 1.7% is of the same order as the accuracy of the solution in representing the field (see Fig. 5). For more precise designs the numbers of windings get larger, which is troublesome to construct and causes the coils to have larger inductances.

A second way is to use a current divider. Connect the different loops in parallel, each with an appropriately chosen resistance in series. This way the ratios between the currents in each loop can be tuned precisely. However, a practical realization will most likely involve routing all loops out of the system where the current divider is installed. For more complicated coil systems with tens of different currents this may be impractical.

Yet another way is to split the loops into decades of currents. In the coil, we use as an example the currents 597, 360, 13, 7, 5 (in arbitrary units) may be constructed from a set of wires with three relative currents of 100, 10, and 1, in the following way:

In this way, one can reach better than 1% accuracy in reproducing the solution in practice with only 3 different currents to control, even for complicated designs. Those can be either separately controlled or split with a current divider.

We do not consider any of the above ways superior. It depends on the particular application which one is the best suited.

## VI. AN EXAMPLE OF AN APPLICATION

The method presented may prove useful to precision physics experiments. Among them is the measurement of the electric dipole moment of the neutron, an observable providing direct insight into the fundamental $CP$ symmetry breaking in the strong interaction.^{10} In the experiment at the Paul Scherrer Institute, the electric dipole moment of the neutron is measured with a precise spectroscopy of polarized neutrons in a combination of electric and magnetic fields.^{11} The precision of the measurement sets very strict requirements for the magnetic field stability, yet the experiment is only one of the many in the experimental hall. Some of the other setups use strong magnets, often ramped up and down daily. The neutron electric dipole moment measurement would not be possible without counter measures against the external magnetic field changes.

The experiment uses an active magnetic field stabilization system. The apparatus is located in the middle of a set of large coils, connected in a feedback loop with magnetic field sensors. The present system is described in detail in Ref. 7. The geometry of the coils, six in total, is very simple: each coil is a rectangle and they are all arranged in three perpendicular Helmholtz-like pairs. The simple geometry causes the volume in which the field is stabilized (the volume of interest) to be small relative to the size of the coils. It is only because there was sufficient space around the less than 3 m large apparatus was it possible to make these simple coils large enough, 6–8 m side length, for the experiment to fit in the volume of interest. Also, the variety of field shapes that can be compensated is limited to homogeneous ones (parallel currents in a Helmholtz-like pair) plus fields created when each coil is controlled separately.

Using the method presented here to design the active compensation coils offers improvements in two areas. First, the size of the coils could be decreased, or the size of the experimental set-up increased, without loss of performance. Second, more coils could be built, extending the range of possible fields to be compensated. These could either be generic: e.g., pure independent gradients, or coils dedicated to counteract a particular known disturbance, designed based on a field map as the goal field in the described method.

## VII. CONCLUSION

Coil design is a complicated and very technical problem, especially when high accuracy is required. We have presented a method that is simple in terms of both the underlying math and the computational effort. We believe that the design method can find its niche in practical applications, where spatial constraints play a significant role and where a percent level is acceptable in the accuracy of the field produced.

For the sake of clarity, we explained the method for a simple example. We designed a coil for a homogeneous field with only few tiles. However, the method is much more powerful. We conclude by presenting, in Fig. 8, a showcase design with *N* = 6 × (5 × 5) = 150 tiles for compensating a nearby dipole source.

We publish the software implementation of the coil design, including examples, as open-source.^{6}

## ACKNOWLEDGMENTS

The authors would like to thank Georg Bison, Allard Schnabel and other members of the nEDM at PSI collaboration for fruitful discussions. This work was supported by the Swiss National Science Foundation under Grant Nos. 200020_162574 and 200020_172639.