Planar Coil Optimization in a Magnetically Shielded Cylinder

Hybrid magnetic shields with both active field generating components and high-permeability magnetic shielding are increasingly needed for a variety of technologies and experiments that require precision-controlled magnetic field environments. However, the fields generated by the active components interact with the passive magnetic shield, distorting the desired field profiles. Consequently, optimization of the active components needed to generate user-specified target fields must include coupling to the high-permeability passive components. Here, we consider the optimization of planar active systems, on which an arbitrary static current flows, coupled to a closed high-permeability cylindrical shield. We modify the Green's function for the magnetic vector potential to match boundary conditions on the shield's interior surface, enabling us to construct an inverse optimization problem to design planar coils that generate user-specified magnetic fields inside high-permeability shields. We validate our methodology by designing two bi-planar hybrid active--passive systems, which generate a constant transverse field, $\mathbf{B}=\mathbf{\hat{x}}$, and a linear field gradient, $\mathbf{B}=(-x~\mathbf{\hat{x}}-y~\mathbf{\hat{y}}+2z~\mathbf{\hat{z}})$, respectively. For both systems, the inverse-optimized magnetic field profiles agree well with forward numerical simulations. Our design methodology is accurate and flexible, facilitating the miniaturization of high-performance hybrid magnetic field generating technologies with strict design constraints and spatial limitations.

Hybrid magnetic shields with both active field generating components and high-permeability magnetic shielding are increasingly needed for various technologies and experiments that require precision-controlled magnetic field environments. However, the fields generated by the active components interact with the passive magnetic shield, distorting the desired field profiles. Consequently, optimization of the active components needed to generate user-specified target fields must include coupling to the high-permeability passive components. Here, we consider the optimization of planar active systems, on which an arbitrary static current flows, coupled to a closed high-permeability cylindrical shield. We modify the Green's function for the magnetic vector potential to match boundary conditions on the shield's interior surface, enabling us to construct an inverse optimization problem to design planar coils that generate user-specified magnetic fields inside high-permeability shields. We validate our methodology by designing two bi-planar hybrid active-passive systems, which generate a constant transverse field, B =x, and a linear field gradient, B = (−xx − yŷ + 2zẑ), respectively. For both systems, the inverse-optimized magnetic field profiles agree well with forward numerical simulations. Our design methodology is accurate and flexible, facilitating the miniaturization of high-performance hybrid magnetic field generating technologies with strict design constraints and spatial limitations.

I. INTRODUCTION
Tailored, high-precision, low magnetic field environments are required for many applications, devices, and experiments. Examples include magnetic field control in quantum sensing of gravity for underground surveying and mapping [1][2][3][4], magnetic field cancellation for atomic magnetometry [5] with applications in medical imaging such as magnetoencephalography [6][7][8][9] and neonatal/fetal magnetocardiography [10][11][12], and noise suppression in fundamental physics experiments [13][14][15]. Usually, these systems are enclosed by a shield formed from high-permeability material, which excludes stray external magnetic fields that can limit the accuracy and sensitivity of the measurements. In particular, cylindrical shield geometries are often used because the dimensions and spacings of multiple cylindrical shield layers can be optimized to generate a large interior shielded region [16,17]. The magnetic field interior to the shield can then be adjusted by active field generating components to either cancel background fields further or to define a specific field environment. However, the surrounding passive shielding material deforms the magnetic field profiles generated by the active components, making it hard to design wire patterns which accurately generate specified target magnetic field profiles [18]. Boundary element methods (BEMs) can be used to optimize magnetic fields generated by surface currents on a triangular mesh [19][20][21][22] to generate arbitrary target magnetic fields. BEMs are extremely powerful and flexible since they can be used to define active systems with complex geometries inside passive shields. They are, however, limited by computational power with results depending on mesh size and on the distance of the active components from the shielding material. Alternatively, analytical methods for optimizing hybrid active-passive systems are advantageous since they provide intuition and understanding of how magnetic fields are distorted by the presence of high-permeability material. Analytical formulations also provide a fast and efficient route to determine the best system for generating a bespoke user-specified magnetic field profile. Currently, however, analytical models for hybrid active-passive systems are restricted to a limited number of scenarios. Simple discrete coil geometries have been formulated in cylindrical high-permeability magnetic shields, where the magnetic field is decomposed into azimuthal Fourier modes [23][24][25] and matched at the shield boundary. Planar highpermeability materials have been incorporated into optimization procedures using the method of mirror images [26,27] to determine the total magnetic field generated by a static current source inside a magnetically shielded room [28][29][30][31]. Similar work in magnetic resonance imaging (MRI) has investigated the interaction of switched magnetic field gradients with high-conductivity materials used for passive shielding [32]. More recently, solutions for the interaction between a static current source and a finite closed cylindrical high-permeability shield have been formulated for the specific case of cylindrical co-axial surface current geometries [33]. These formulations incorporate the effect of the shielding material by explicitly solving Maxwell's equations by matching the required boundary conditions on the surface interface to find the total magnetic field generated by the system. Building on this previous work, here we determine the effect of a closed finite length high-permeability cylinder on the magnetic field generated by an arbitrary static current distribution on an interior circular plane, oriented perpendicular to the axis of the cylinder. We calculate the effect of the high-permeability cylinder on the magnetic field produced by current flow on the plane, which enables us to determine optimal current paths to generate user-specified target fields inside the cylinder. We first derive the vector potential generated by an arbitrary planar current source. Then, we calculate a pseudo-current density induced on the cylindrical surface of the high-permeability material in response to the planar current source, and, hence, derive a Green's function for our system. Finally, we implement a Fourier decomposition of the current paths to calculate the total magnetic field in terms of a set of weighted Fourier coefficients. This formulation allows the incorporation of a quadratic optimization procedure to determine globally optimal designs that generate user-specified target magnetic field profiles. Through the use of this optimization procedure, we design two example bi-planar coil systems optimized for operation inside a finite-length closed high-permeability cylindrical magnetic shield. In both cases, we confirm that our analytical model agrees well with the result of numerical finite element simulations. Our work extends the range of coil geometries which can be efficiently optimized to generate static user-specified target magnetic fields in the presence of high-permeability materials, expanding design flexibility for systems which require precision-controlled magnetic field environments.

II. THEORY
High-permeability magnetic shields, usually constructed from materials such as mumetal, are used to attenuate external background magnetic fields in order to shield magnetically sensitive equipment from spurious signals. These high-permeability materials generate an induced magnetization, M, in response to a incident field at their surface, S. This ensures continuity of the parallel and tangential components of the magnetic field, B, and magnetic field strength, H, respectively, at the interface between air and a material. Considering this interface while working in the magnetostatic regime with no surface currents, the boundary conditions take the form and wheren is the unit vector normal to the boundary and µ r is the relative permeability of the material. To design magnetic fields effectively in shielded environments, we need to determine the total field generated by an active coil structure and the high-permeability material. The total magnetic field in free space is related to the magnetic field strength and the magnetization by The induced magnetization can be formulated in terms of a pseudo-current densityj, confined to the material's surface with the magnetic field strength related to a current density j c , on an active structure, using Ampère's law, By using (3)- (5), and the relation between the vector potential and the magnetic field, B = ∇ ∧ A, the total vector potential in free space generated by the system can be cast as the Poisson equation, resulting in the integral solution in terms of an arbitrary current density where G(r, r ) is the associated Green's function for the system [26].
Here, we consider the specific example of a closed cylindrical magnetic shield surrounded by free space, as shown in Fig. 1. This cylinder has relative permeability, µ r 1, radius ρ s , and length L s , with planar end caps located at z = ±L s /2. Inside this cylinder, an arbitrary current flows on a circular planar surface of radius ρ c < ρ s , centered on the z-axis and lying in the z = z plane where |z | < L s /2. Previous work has shown that the magnetic field generated FIG. 1: A cylindrical hollow high-permeability magnetic shield of length L s , radius ρ s , and with planar end caps located at z = ±L s /2 encloses an interior circular current source (bounded by the upper dashed curve) of radius ρ c < ρ s , which lies in the z = z plane.
by a high-permeability material in response to an applied static field deviates from that of a perfect magnetic conductor on the scale of O µ −1 r [34,35]. More specifically, in the case of a closed finite length cylinder with permeability µ r = 20000, equal to industrial standard mumetal, a material thickness of d = 1 mm is sufficient to provide a response similar to that of bulk material, differing from the perfect magnetic conducting response only by approximately 0.005% [33]. Therefore, we can approximate sufficiently thick materials with high-permeability µ r > 20000, such as mumetal, to that of a perfect magnetic conductor without compromising the accuracy of a theoretical model. Assuming the surface of the shield shown in Fig. 1 is a perfect magnetic conductor, the boundary condition at its surface can be written as where the magnetic field B = B ρρ + B φφ + B zẑ , is expressed in cylindrical polar coordinates. As previously formulated for cylindrical current flow [33], the response generated by a finite closed cylinder with high-permeability, µ r 1, can be determined by matching the orthogonal modes in the magnetic field and generating a relation between the initial current flow and the resulting magnetization induced on the surface of the cylinder. Here we apply the same methodology to generate an analytical expression for the response of a finite closed high-permeability cylinder to a planar current source using a similar decomposition. The vector potential generated by a current source in cylindrical coordinates can be expressed as, A z (r) = µ 0 r d 3 r G(r, r )j z (r ).
Since there is no current flow in the z-direction for a planar current source perpendicular to the axis of the cylindrical shield, the continuity equation can be used to express the planar current flow in terms of a single scalar streamfunction defined on the planar surface where ϕ(r ) = ϕ(ρ , φ ). To exploit the radial symmetries of the system, we choose to decompose the Green's function in cylindrical coordinates in terms of Bessel functions of the first kind, allowing the vector potential to be expressed in terms of cylindrical harmonics defined on a circular plane. Using (9)-(11), (12), and (13) we cast the vector potential in a simplified form where J m (z) is the derivative of J m (z) with respect to z, and ϕ m (k) is defined as the m th order Hankel transform We now consider the vector potential generated by the magnetic shield. To do this, we first introduce a pseudo-current density induced on an infinite cylinder,j =j φ (φ , z )φ +j z (φ , z )ẑ. Next, we seek a Fourier representation of this pseudo-current density, which satisfies the boundary condition over the entire domain of the shield. In particular, we wish to equate the shared azimuthal Fourier modes at the radial boundary of the shield cylinder. This is achieved using a combination of methods that must be applied sequentially, since each method satisfies the condition at an orthogonal boundary. The radial condition is satisfied by equating the magnetic field generated by the cylindrical pseudo-current density and planar current flow, generating a relation between the response of an infinite cylindrical shield and the initial current source. Then, the boundary condition at the end caps can simultaneously be satisfied by applying the method of mirror images [26]. These methods can be combined because the infinite pseudo-current density and planar current flow are spatially orthogonal to the end caps, meaning that any reflections generated by the application of the method of mirror images continue to satisfy the radial condition. The components of the vector potential generated by a pseudo-current density induced on an infinite cylinder [32] in the region ρ < ρ s are given by where the Fourier transforms of the pseudo-currents are defined by The corresponding inverse transforms are given bỹ Therefore, adding the contributions from the planar current flow, (14)- (15), and the infinite pseudo-current density, (18)- (20), while using (21)- (24), and applying the method of mirror images for two infinitely large parallel planes, we can write the total magnetic field generated by the system in the region ρ < ρ s as where j mp φ (k) is the p th reflected Fourier transformed azimuthal pseudo-current density induced on the cylindrical surface of the magnetic shield with I m (z) and K m (z) defined as the derivatives with respect to z of the modified Bessel functions of the first and second kind I m (z) and K m (z), respectively. By applying the boundary condition at the radial surface, (8), we can match the shared m th azimuthal Fourier mode generated by each p th reflected pseudo-current and streamfunction, resulting in the relation Physically, due to the formulation of the response in terms of a pseudo-current density, there must be a unique solution that is independent of the axial position that satisfies the boundary condition over the infinite domain of the cylindrical shield. Therefore, we perform an inverse Fourier transform with respect to z to generate an integral representation of the p th reflected Fourier pseudo-current density, j mp φ (k), in terms of the m th order Hankel transform of the streamfunction defined on the planar surface, ϕ m (k), This expression for j mp φ (k) can now be substituted into (25)- (27) to determine the total magnetic field in terms of ϕ m (k). Next, we must choose an appropriate expansion of the streamfunction, ϕ(ρ , φ ). Although the choice of orthogonal basis for the expansion of the streamfunction is somewhat arbitrary, a choice of basis which considers the symmetries between the Hankel transform, coordinate system, and the integral representation of the pseudo-current yields a simpler solution. Here, we choose to decompose the radial component of the planar current flow into a Fourier-Bessel series while using a Fourier series representation of the azimuthal dependence, where (W nm , Q nm ) are Fourier coefficients and ρ nm is the n th zero of the m th Bessel function of the first kind, J m (ρ nm ) = 0. Therefore, using (25)- (27), (29), and (30), the total magnetic field generated by an arbitrary current flow on the planar surface inside the closed finite high-permeability cylinder can be written as where and withp = pπ/L s . A full derivation of these expressions is given in Appendix A. Solving for the unknown Fourier coefficients, (W nm , Q nm ), to generate a desired magnetic field using the system of governing equations (31)- (33) is an ill-conditioned problem due to the formulation of the vector potential through the integral representation in (9)- (11). As in previous work by Carlson et al. [36], this may be solved by a least squares minimization with the addition of a penalty term that acts as a regularization parameter. The choice of regularization parameter is arbitrary since it exists only to facilitate the inversion. Well-regularized solutions yield more simplistic coil designs at a cost to the field fidelity [37]. Here, we choose the regularization parameter to be the power, P , dissipated by a circular planar current source of thickness, t, and resistivity, , Substituting (30) into (12) and then (41), and integrating over the planar surface we find that, for m = 0, and, for m ∈ Z + , where iFj is the regularized hypergeometric function: see Appendix B for a full derivation. We can now formulate a cost function for the least squares minimization, where β is a weighting parameter chosen to adjust the physical constraints of the system. The cost function is minimized using a least squares fitting to calculate the optimal Fourier coefficients to generate the desired magnetic field at K target points. The minimization is achieved by taking the differential of the cost function with respect to the Fourier coefficients, which enables the optimal Fourier coefficients to be found for any given physical target magnetic field profile through matrix inversion. The inversion process yields the optimal continuous streamfunction defined on the planar surface for a finite number of Fourier coefficients. In the ideal case, the number of Fourier coefficients would be infinite. However, a finite number of terms can provide accurate solutions in well-regularized problems. This number can be approximated by ensuring that the distance between the planar current-carrying surface and the closest target field point is much larger than the smallest spatial frequency. This is due to the decreased contribution of higher-order terms at distances much larger than their spatial frequency. The final objective is to design a coil that generates the desired magnetic field to a specified accuracy. To do this, a discrete approximation of the field profile may be found by contouring the streamfunction at N ϕ discrete levels. The contours of the streamfunction generate streamlines where wires should be laid to replicate the desired target magnetic field. This is achieved by discretizing the streamfunction into N ϕ contours, where ϕ j = min ϕ + (j − 1/2)∆ϕ, j = 1, ..., N ϕ , separated by ∆ϕ = max ϕ−min ϕ

Nϕ
, and the total current through each wire is I = ∆ϕ. The number of contours should be maximized, limited only by manufacturing since the accuracy of the theoretical model depends on the quality of the discrete approximation of the streamfunction. The approximation to the streamfunction can be determined by using the elemental Biot-Savart law to calculate the error as N ϕ is increased. In the case that multiple current-carrying planes are designed, the contours should be defined evenly between the global maximum and minimum of the streamfunction across all the planes so that the current through each wire is equal.

III. RESULTS
We now verify our analytical model by designing hybrid active-passive systems composed of bi-planar coils inside a closed high-permeability magnetic shield and compare the resulting magnetic field profiles with forward numerical simulations of each optimized system. We consider two distinct systems, each containing current confined to two disks of radius ρ c = 0.45 m and symmetrically placed at z = ±0.45 m. Both systems are interior to a perfect closed cylindrical magnetic shield of radius ρ s = 0.5 m and length L s = 1 m centered on the origin, as shown in Fig. 2. The first system is designed to generate a constant transverse field, B = B 0x , and the second system creates a linear field gradient, B = G(−xx − yŷ + 2zẑ), within an optimization region defined by −z /2 ≤ z ≤ z /2 and 0 ≤ ρ ≤ ρ c /4. The magnetic field profiles chosen, i.e. the uniform transverse and linear field gradient fields, are examples of tesseral (m = 0) and zonal (m = 0) harmonics, respectively, which exhibit m-fold and complete azimuthal symmetry [38], which facilitates analysis of the shield's particular response to tesseral and zonal harmonic fields generated by planar current sources. The two systems that we consider here are chosen to illustrate targeted magnetic field compensation in situations where compact systems are required, but space inside the central cylindrical cavity of the system is limited by the presence of experimental equipment (e.g. magnetic sensors). Coil designs generating other field profiles, for example using combinations of planar coils with varying sizes, can be found in the Appendix C. All designs can be replicated using our open-access Python code [42]. In Fig. 3a and Fig. 5a, we show the optimized contoured streamfunctions for the constant transverse and linear field gradient systems, respectively. Due to the symmetric placement of the bi-planar coil systems and optimized field regions within the shield (see Fig. 2), the magnitudes of the streamfunctions defined on both coils are identical. The current flow directions, however, are opposite, due to the form of the desired fields. Figures. 3b-d show the transverse magnetic field, B 0 = 1 µT, along the x-axis and z-axis, respectively, generated by the bi-planar coil design shown in Fig. 3a calculated in three different ways: analytically using (31)-(33) (solid red curves); numerically using COMSOL Multiphysics ® with the shield treated as a perfect magnetic conductor (blue dotted curves); numerically in free space, i.e. excluding the high-permeability material and evaluating the magnetic fields through the Biot-Savart law for discretized bi-planar coils with N ϕ = 250 (dashed green curve in b). Furthermore, color maps of the transverse and axial magnetic field components in the xz-plane are presented in Fig. 4a and Fig. 4b, respectively, calculated numerically using COMSOL Multiphysics ® with the shield treated as a perfect magnetic conductor. Similarly, Fig. 5b-c show the linear axial field gradient, G = 1 µT/m, along the x-axis and z-axis, respectively, generated by the bi-planar coil system shown in Fig. 5a in the same three cases with similar color maps of the magnetic field components in the xz-plane shown in Fig. 6a and Fig. 6b. For both designs, the analytical field profiles agree well with numerical simulations. The maximum errors between the analytical model and numerical simulation are 0.052% and 0.043%, for the constant transverse field and the gradient of the linear gradient field, respectively, along the z-axis of the optimization region. To quantify the performance of our optimization procedure, we can analyze the deviations between the magnetic fields generated by our theoretical model and the desired target fields. Examining the fidelity of the fields generated by both systems along x-axis of the optimization region, the maximum absolute deviations from the constant transverse and linear axial gradient fields are 6.78% and 0.380%, respectively. Along the z-axis of the optimization region, the maximum deviations are 7.50% and 0.306%, respectively. Clearly, the axial field gradient field is generated more accurately than the uniform transverse field. This can be seen from Fig. 3d, which reveals small oscillations in the transverse field over the z-axis of the optimization region. To understand the difference between the fidelity of the field profiles generated by the two systems, we must analyze how the passive magnetic shield affects the fields from the planar current distributions, decomposing its response into zonal (m = 0) and tesseral (m = 0) harmonic components. These harmonic responses relate to the variations in the induced pseudo-current density, corresponding to the planar streamfunction (29), required to satisfy the boundary condition at the shield wall. Schematic approximations of the surface currents for the m = 0 zonal and m = 1 tesseral harmonic responses can be seen in Fig. 7. For the boundary condition on the cylindrical surface to be satisfied, the induced azimuthal pseudo-currents must mirror azimuthal current paths on the planar coil surfaces, so that the associated fields cancel in the region ρ s > ρ > ρ c . Since zonal responses are composed of simple circular loops, which can be formed in either a cylindrical or planar basis, the response of the magnetic shield enhances the magnetic field in the optimization region. Consequently, the shield amplifies the axial magnetic field by a factor of 2.39 at the shield's center, as shown in Fig. 5b-c. The uniform field gradient therefore exhibits superior fidelity because the continuum response of the passive shield approximates a distributed cylindrical coil. The resulting system, composed of both the coil and shield, completely encloses the interior region, producing a high-fidelity magnetic field gradient. The tesseral responses are more complicated, with the cylindrical surface of the magnetic shield acting to oppose the magnetic field generated by the planar system in the region ρ < ρ s due to the formation of saddle-type currents [38]. These currents result from the required continuity of the induced azimuthal pseudo-current density in a cylindrical basis. Due to the restricted current flow on the planar surface, the only way to mitigate the opposing field generated by the cylindrical surface of the shield is to minimize the magnetic field at the boundary, resulting in field coil designs that are oscillatory, as seen in Fig. 3a. The response of the passive shield, for this configuration, reduces the magnetic field by a factor of 2.41 at the shield's center, as shown in Fig. 3b-c. The shield's response not only makes the optimization of tesseral harmonic fields more difficult in regions close to the cylindrical surface of the magnetic shield, but also has an effect on the level of fidelity that can be achieved over any region when the coils are in close proximity to the magnetic shield. To demonstrate this, we investigate the fidelity of the constant transverse field generated by optimized bi-planar coils, similar to those in Fig. 3, within a magnetic shield whose radius varies over the range ρ s = [ρ c , 3ρ c ]. In Fig. 8a we plot field profiles at selected shield radii to show how the uniformity of the transverse field along the x-axis improves as the radius of the shield increases, as expected, due to the reduced inhomogeneities introduced by the cylindrical surface of the magnetic shield at distances further from the field coils. To evaluate the deviation in the field uniformity, we can calculate the mean RMS error, which represents the averaged deviation of the transverse field from perfect uniformity. In Fig 8b, we evaluate the mean RMS error in the transverse field along the z-axis of the optimization region as the radius of the magnetic shield increases. The uniformity of the optimized fields in Fig. 8a and Fig. 8b initially decreases rapidly as the radius of the magnetic shield increases, but approaches the wide shield limit when ρ s ∼ 2ρ c , at which point the field generated by the cylindrical surface of the magnetic shield becomes negligible. In the wide shield limit the predominant contribution from the magnetic shield can then be approximated through the use of the method of mirror images from two infinite parallel perfect magnetic conductors. Although the pseudo-currents do not exist physically, analysis of the magnetic fields generated by the shield in response to the planar coils gives insight into how such coils should be designed in order to best realize a specified target field. For example, generating magnetic fields that require coil geometries with m-fold azimuthal symmetry causes the cylindrical surface of the magnetic shield to oppose the magnetic field generated by the planar coil. Consequently, if a tesseral field is required over a large radial region, the distance from this region to the planar coils and from the cylindrical shield surface should be minimized and maximized, respectively. The opposite is true for the fields generated by the planar end caps of the shield. When the bi-planar coils are located near the end caps, where −L s /2 < −z < z < z < L s /2, but radially distant from the cylindrical surface of the magnetic shield, the reflected pseudo-currents generated by the method of mirror images provide a field similar to that of the planar coils, resulting in a field that is magnified. Consequently, bi-planar coils are desirable in magnetic shields with large radii and small aspect ratios, L s /ρ s < 1. In comparison, the magnetic field generated by the shield in response to coils contained on a cylindrical surface [33] shows the opposite effect. Tesseral fields generated by cylindrical coils are enhanced by the cylindrical surface and are reduced by the planar end caps, favoring long magnetic shields with L s /ρ s > 1. Due to the conflicting conditions on the geometries (i.e. planar or cylindrical) of the coils and the magnetic shield, a system composed of coupled planar and cylindrical coils may have advantages for generating desired tesseral field profiles with the greatest accuracy. Intrinsic inaccuracies in calculating the optimum current density and approximating the current continuum by discrete wires are not the only sources of error that must be considered when designing fields using our method. Particularly for the tesseral harmonic field generating coils, wire patterns may be highly meandering, making them hard to manufacture and, due to their high resistance, power consuming. When manufacturing these systems, the achievable field fidelity is limited by the discretization of the continuum current flow pattern and by the creation of magnetic dipoles via the interaction between separate current streamlines and the wires that follow them. Fortunately, new ways to realize the current continuum are emerging due to advances in foldable PCBs [40] and 3D-printing technologies [41], which enable more complex coils to be made accurately. It should also be noted that approximating a perfect magnetic shield by a material of finite permeability µ r = 20000 and thickness d = 1 mm only introduces small deviations in the model of 0.005% [33], which is much less than the error introduced in the desired field by the coil designs. As a result, given an accurate representation of the current continuum, our methodology can be used reliably to generate target magnetic fields in high-permeability environments. The python code used to design arbitrary the planar coils in a magnetically shielded cylinder is openly available from GitHub and can be cited at http://doi.org/10.5281/zenodo.4442661 [42]. Verification using COMSOL Multiphysics ® requires a valid license.  Fig. 2. Blue and red shaded regions correspond to current counter flows and their intensity shows the streamfunction magnitude from low (white) to high (intense color). Solid and dashed black curves represent discrete wires with opposite senses of current flow, approximating the current continuum with N ϕ = 12 global contour levels. Streamfunction on the lower light gray plane in Fig. 2 is geometrically identical but the current direction is reversed. (b) B x , calculated versus transverse position, x, for y = z = 0, from the current continuum in (a) in three ways: analytically using (31)-(33) (solid red curve); numerically using COMSOL Multiphysics ® Version 5.3a, modeling the high-permeability cylinder as a perfect magnetic conductor (blue dotted curve); numerically without the high-permeability cylinder and using the Biot-Savart law with N ϕ = 100 contour levels (dashed green curve). Black lines enclose the regions where the calculated field deviates from the target field by 5% (dashed) and 1% (dot-dashed). (c) B x , calculated versus axial position, z, for x = y = 0, from the current continuum in the same three ways as for (b). (d) Enlarged section of (c) showing good agreement between the numerical and analytical results throughout the optimization region. FIG. 5: Wire layouts (a) and performance (b-c) of a hybrid active-passive system optimized to generate a linear variation of B z with z-position. Current flows on the surface of two disks of radius ρ c = 0.45 m, which are separated symmetrically from the origin and lie at z = ±0.45 m. The wire layouts are optimized to generate a linear axial field gradient, dB z /dz = 2 µT/m, along the z-axis of the cylinder through the harmonic B = G(−xx − yŷ + 2zẑ). The current-carrying planes are placed symmetrically inside a perfect closed magnetic shield of radius ρ s = 0.5 m and length L s = 1 m and the magnetic field is optimized between ρ = [0, ρ c /4] and z = ±z /2, as shown in Fig. 2. The least squares optimization was performed with parameters N = 50, M = 0, β = 1.77 × 10 −9 T 2 /W, t = 0.5 mm, and ρ = 1.68 × 10 −8 Ωm. (a) Color map of the optimal current streamfunction calculated for the upper current-carrying plane in Fig. 2. Intensity of red shaded regions shows the streamfunction magnitude from low (white) to high (intense color). Solid and dashed black curves represent discrete wires with opposite senses of current flow, approximating the current continuum with N ϕ = 12 global contour levels. Streamfunction on the lower light gray plane in Fig. 2 is geometrically identical but the current direction is reversed. (b) Transverse magnetic field, B x , calculated versus transverse position, x, for y = z = 0, from the current continuum in (a) in three ways: analytically using (31)-(33) (solid red curve); numerically using COMSOL Multiphysics ® Version 5.3a, modelling the high-permeability cylinder as a perfect magnetic conductor (blue dotted curve); numerically without the high-permeability cylinder and using the Biot-Savart law with N ϕ = 100 contour levels (dashed green curve). Black lines enclose the regions where the calculated field deviates from the target field by 5% (dashed) and 1% (dot-dashed). (c) Axial magnetic field, B z , calculated versus axial position, z, for x = y = 0, from the current continuum in (a) in the same three ways as for (b).  The discrete current loops (black) are located at z = z 1 and z = z 2 with their p th reflected pseudo-image-currents (red) generated by the surface of the cylindrical shield at ρ = ρ s and by the planar end caps at z = ±L s /2, in accordance with the method of images described by (29). Images of the two planar (black) coils resulting from the end caps are located at z = (−1) p z 1 + pL s and z = (−1) p z 2 + pL s , respectively, where p ∈ Z (two such image coils are shown blue). For all odd reflections the axial current direction is reversed.

IV. CONCLUSION
Here we have formulated an analytical model to calculate the magnetic field generated by an arbitrary static current distribution confined to a plane whose normal is parallel to the axis of a finite closed high-permeability cylinder. Our formalism is based on a Green's function expansion of the vector potential generated by a planar coil system. To satisfy the boundary conditions at the surface of the cylinder, we assume that it is a perfect magnetic conductor and express the induced magnetization in terms of a pseudo-current density. Due to the shared azimuthal symmetries of the planar current flows and the induced pseudo-current on the cylindrical surface of the magnetic shield, the response of the magnetic shield can be expressed in terms of the planar current distribution. We used this formalism to enable a priori optimization of the planar current distribution required to generate specified target magnetic field profiles by combining a Fourier-Bessel decomposition of the current flow with a quadratic minimization procedure.
We demonstrated this optimization methodology by using it to design bi-planar coils that accurately generate a constant transverse field, B x , across the cylinder, and a linear gradient field, dB z /dz, along the axis of the cylinder. Predictions from our analytical model agree well with with subsequent forward finite element simulations, validating our methodology. We quantified the interaction of planar systems with a closed finite high-permeability cylinder and found that all coil designs without complete azimuthal symmetry induced magnetic fields in the cylindrical surface of the shield that opposed the field generated by the planar current source. Conversely, planar end caps have the opposite effect and amplify the field from the planar current source.
Our analytical model and optimization procedure extends the range of current geometries that can be tailored within finite closed cylindrical magnetic shields in order to accurately generate specified target field profiles. This enables the development of systems that require the best magnetic field control that can be achieved subject to size, weight, power and cost constraints. It may also enable hybrid active/passive planar coils to be retrofitted to existing cylindrical magnetic shields in order to improve their performance for a low cost and without disrupting existing experimental setups. In the future, combining planar and cylindrical coils could enable ultra-precise magnetic fields to be generated through even larger interior fractions of magnetic shields. The PYTHON code used to design arbitrary planar coils in a magnetically shielded cylinder is openly available from GitHub [42]. Verification using COMSOL Multiphysics requires a valid license. generality, applying a Bessel function product identity (5.54.1 in [39]), we may now write Similarly, we may also calculate the induced Fourier pseudo-current density. Substituting (30) into (29) and integrating over the azimuthal coordinate yields To simplify these expressions, we evaluate the integral overk first and then separate the integrand using partial fractions to yield Inspecting (49) we see that the first component is simply the orthogonality relation between Bessel functions, and the second has a known solution (6.541.1 in [39]). Hence, we can state that However, since ρ < ρ s , the first term in (50) can be neglected. Inserting (50) into (48) and integrating over ρ results in the expression Substituting (47) and (51) into (25)- (27), we can express the summation over the exponentials associated with the infinite reflections of the planar streamfuction as ∞ p=−∞ e −k|z−(−1) p z +pLs| = e −k|z−z | + 2 e 2kL − 1 e kL cosh (k (z + z )) + cosh (k (z − z )) , and Similarly, the summation over the infinite pseudo-current reflections can be written as a Fourier series expansion, The power dissipated in the surface of a circular planar current source of thickness t and resistivity is given by We can substitute the streamfunction (30) into the continuity relations (12) to obtain Inserting (66)-(67) into (65) and integrating over the azimuthal component results in which, when integrated over the radial component, can be expressed as for m = 0, and for m ∈ Z + , where iFj is the regularized hypergeometric function. Some specific evaluations of (70) are:  Fig. 2 [blue and red shaded regions correspond to current counter flows and their intensity shows the streamfunction magnitude from low (white) to high (intense color)]. Solid and dashed black curves represent discrete wires with opposite senses of current flow, approximating the current continuum with N ϕ = 12 global contour levels. Streamfunction on the lower light brown plane in Fig. 2 is geometrically identical but the current direction is reversed. (b) Transverse magnetic field, B x , calculated versus axial position, z, from the current continuum in (a) in three ways: analytically using (31)-(33) (solid red curve); numerically using COMSOL Multiphysics ® Version 5.5a, modeling the high-permeability cylinder as a perfect magnetic conductor (blue dotted curve); numerically without the high-permeability cylinder and using the Biot-Savart law with N ϕ = 100 contour levels (dashed green curve). Black lines enclose the regions where the axial field gradient deviates by 5% (dashed) and 1% (dot-dashed).
To further validate our design process, in Figs. 9-10 we show additional bi-planar hybrid active-passive systems and analyze their behavior. These systems are designed using our open-access Python code [42]. The coordinate axes and magnetic field plots are labeled in the same way as the systems presented in the main text. For both systems, the coils lie in the z = ±0.45 m planes, as in the main text and shown in Fig. 2, and the fields are optimized along the z-axis between z = ±z /2. The design in Fig. 9 generates linear transverse field gradient, dB x /dz, whose spatial variation matches the harmonic B = (zx + xẑ) between two identical plates of radius ρ c = 0.45 m, inside a magnetic shield of radius ρ s = 0.5 m and length L s = 1 m. The field gradient is amplified by a factor of 1.06 at the shield's center. If the shield is lengthened, the field gradient is diminished to 0.81 at the shield's center in the long shield limit, as expected from our analysis in the main body of the text. The design in Fig. 10 generates a uniform field, B z , in the FIG. 10: Wire layouts (a-b) and performance (c) of a hybrid active-passive linear system optimized to generate a constant axial field, B z . Current flows on the surface of two disks of radii ρ c1 = 0.95 m (on the upper current-carrying plane in Fig. 2) and ρ c2 = 0.35 m (on the lower current-carrying plane in Fig. 2), respectively, which are separated symmetrically from the origin and lie at z = ±0.45 m. The wire layouts are optimized to generate a constant axial field, B 0 = 1 µT, along the cylinder and parallel to its axis of symmetry. The current-carrying planes are placed symmetrically inside a perfect closed magnetic shield of radius ρ s = 1 m and length L s = 1 m and the magnetic field is optimized along the z-axis between z = ±z /2. The least squares optimization was performed with parameters N = 100, M = 0, β = 1.77 × 10 −8 T 2 /W, t = 0.5 mm, and ρ = 1.68 × 10 −8 Ωm. (a-b) Color maps of the optimal current streamfunction on the upper and lower planar disks, respectively [blue and red shaded regions correspond to current counter-flows and their intensity shows the streamfunction magnitude from low (white) to high (intense color)]. Solid and dashed black curves represent discrete wires with opposite senses of current flow, approximating the current continuum with N ϕ = 6 global contour levels. (c) Axial magnetic field, B z , versus axial position, z, calculated from the current continuum in (a-b) in three ways: analytically using (31)-(33) (solid red curve); numerically using COMSOL Multiphysics ® Version 5.5a, modeling the high-permeability cylinder as a perfect magnetic conductor (blue dotted curve); numerically without the high-permeability cylinder and using the Biot-Savart law with N ϕ = 100 contour levels (dashed green curve). Black lines enclose the regions where the axial field gradient deviates by 5% (dashed) and 1% (dot-dashed).
z-direction and extending between two planes of different sizes, with upper and lower coil sizes ρ c1 = 0.95 m and ρ c2 = 0.35 m, respectively, inside a magnetic shield of radius ρ s = 1 m and length L s = 1 m. Since the shield and coils have a relatively wide form factor, the field is highly uniform, with the region of axial field deviation below 1% extending over more than 70% of the distance between the planes.