Identification of important error fields in stellarators using Hessian matrix method

TL;DR

This paper enhances a Hessian matrix method to quickly identify critical coil deviations affecting error fields in stellarators, aiming to optimize coil tolerances and reduce costs.

Contribution

The paper introduces an improved analytical Hessian matrix approach to pinpoint sensitive coil deviations impacting stellarator magnetic properties.

Findings

Successfully identified coil deviations affecting magnetic island size.

Determined coil perturbations that deteriorate quasi-symmetry.

Applied method to CFQS configuration with promising results.

Abstract

Error fields are predominantly attributed to inevitable coil imperfections. Controlling error fields during coil fabrication and assembly is crucial for stellarators. Excessively tight coil tolerance increases time and cost, and, in part, led to the cancellation of NCSX and delay of W7-X. In this paper, we improve the recently proposed Hessian matrix method to rapidly identify important coil deviations. Two of the most common figures of merit, magnetic island size and quasi-symmetry, are analytically differentiated over coil parameters. By extracting the eigenvectors of the Hessian matrix, we can directly identify sensitive coil deviations in the order of the eigenvalues. The new method is applied to the upcoming CFQS configuration. Important perturbations that enlarge n/m=4/11 islands and deteriorate quasi-axisymmetry of the magnetic field are successfully determined. The results suggest each modular coil should have separate tolerance and some certain perturbation combinations will produce significant error fields. By relaxing unnecessary coil tolerance, this method will hopefully lead to a substantial reduction in time and cost.

Equations

$$\delta F\approx \frac{1}{2}{\delta \mathbf{x}}^{\mathrm{T}}\cdot {\mathbf{H}}_0\cdot \delta \mathbf{x},$$

$$\delta F\approx \frac{1}{2}{\lambda }_i{\xi }^2.$$

$$w=4\sqrt{\frac{b_{mn}}{m{\iota }_{mn}^{\prime }}},$$

$$b_{mn}=\frac{1}{2{\pi }^2}{\int }_0^{2\pi }{\int }_0^{2\pi }\frac{\mathbf{B}\cdot \nabla {\psi }_t}{{\mathbf{B}}_0\cdot \nabla \zeta }{e}^{-i(m\theta -n\zeta )}\text{d}\theta \text{d}\zeta .$$

$$F_{RP}=\sum _{m,n}(b_{mn}-b_{mn}^o{)}^2,$$

$$\frac{\partial b_{mn}}{\partial X_i}$$

$$\frac{{\partial }^2{b}_{mn}}{\partial X_i\partial X_j}$$

$$F_{QS}=\sum _{n/m\ne N/M}{\left(\frac{B_{m,n}}{B_{0,0}}\right)}^2,$$

$$B_{m,n}=\frac{1}{2{\pi }^2}{\int }_0^{2\pi }{\int }_0^{2\pi }\sqrt{\mathbf{B}\cdot \mathbf{B}}{e}^{-i(m{\theta }_B-n{N}_p{\zeta }_B)}\text{d}{\theta }_B\text{d}{\zeta }_B.$$

$$\frac{\partial F_{QS}}{\partial X_i}=\sum _{n/m\ne N/M}2\frac{B_{m,n}}{{B_{0,0}}^2}\left[\frac{\partial B_{m,n}}{\partial X_i}-\frac{\partial B_{0,0}}{\partial X_i}\frac{B_{m,n}}{B_{0,0}}\right],$$

$$\displaystyle\displaystyle\frac{\partial^{2}F_{QS}}{\partial X_{i}\partial X_{j}}=\sum_{n/m\neq N/M}\frac{2}{{B_{0,0}}^{2}}\bigg{[}$$

$$+B_{m,n}\frac{{\partial }^2{B}_{m,n}}{\partial X_i\partial X_j}-\frac{{B_{m,n}}^2}{B_{0,0}}\frac{{\partial }^2{B}_{0,0}}{\partial X_i\partial X_j}$$

$$\displaystyle-\frac{2{B_{m,n}}}{{B_{0,0}}}\left(\frac{\partial B_{m,n}}{\partial X_{i}}\frac{\partial B_{0,0}}{\partial X_{j}}+\frac{\partial B_{0,0}}{\partial X_{i}}\frac{\partial B_{m,n}}{\partial X_{j}}\right)\bigg{]}\ .$$

$$\mathbf{X}\equiv \{x_{k,n}^c,x_{k,n}^s,y_{k,n}^c,y_{k,n}^s,z_{k,n}^c,z_{k,n}^s,I_k\},$$

$${\xi }_{tol}\approx \sqrt{2{\Delta F}_{tol}/{\lambda }_1}.$$

Full text

K****eywords Stellarator coils  $\cdot$ error field  $\cdot$ Hessian matrix  $\cdot$ sensitivity analysis

1 Introduction

Plasma performance in magnetically confined fusion devices depends on the quality of the magnetic field. Magnetic field irregularities, namely ‘error fields’, can lead to the destruction of magnetic surfaces [1] and locked modes [2]. One of the main sources of error field is inevitable coil deviations. The stellarator generally has more complicated coils than axisymmetric devices (like tokamaks) and the confining magnetic field predominantly arises from carefully shaped external coils. Therefore, controlling error fields during coil fabrication and assembly is crucial for stellarator construction. This is even more challenging than expected, as high accuracy requirement on modular coils was the largest driver of cost growth for the National Compact Stellarator Experiment (NCSX) and partly led to the project cancellation [3, 4]. The high demand for coil tolerance was also identified during the assembly of W7-X [5].

As-built coil geometries might differ from designed models in location, orientation or even shape. Coil deviations, as measured between designed models and final built shapes, have different effects on plasma performance. The cost and complexity of a device can increase dramatically with tight construction tolerance. Thus, it is necessary to carry out an error field sensitivity analysis for identifying important coil deviations and determining acceptable tolerance prior to machine construction.

The Large Helical Devices (LHD) evaluated the irregular magnetic fields produced by global deformations and local irregularities on coils. Poloidal field (PF) and helical coils were represented with Fourier coefficients and the vacuum flux surfaces were plotted out when perturbing different Fourier modes. It was found that horizontal shifts ($n=1$) had the most significant influence on the destruction of magnetic surfaces [6]. On the Columbia Non-neutral Torus (CNT), Kremer [7] analyzed the volume of magnetic surfaces under hundreds of random coil displacements, including shifts and rotations of coils and he suggested the PF coils could be misplaced by a maximum distance of 1 cm while inter-linked (IL) coils could only be perturbed up to 2 mm. This was also observed later by Hammond et al. [8] when the gradient of the rotational transform on the magnetic axis with respect to several defined rigid displacements were calculated by the finite difference. NCSX investigated the impacts of systematic coil geometric perturbations and tolerance schemes on magnetic island size [9]. Local errors including coil-plasma spacing, short wavelet orthogonal displacements, and broad coil deformations causing coil length errors were also examined [10]. The results from NCSX suggested that modular coils required more tight tolerance (about 1.5 mm) than PF & TF coils, and particularly the inboard regions of modular coils had more significant effects on flux quality, while errors in other regions might approach 3 mm or even larger coil tolerance. As the largest stellarator, the Wendelstein 7-X (W7-X) performed extensive studies on coil tolerance [11]. The primary criteria were the resonant magnetic perturbations, $B_{11},B_{22},B_{33}\ \&\ B_{44}$ since W7-X has a $m/n=5/5$ island chain outside the last closed flux surface (LCFS). Numerical investigations [12, 13] showed that the resonant magnetic field perturbations were the most sensitive to rotations of coils and modules, after the effects of manufacturing errors, shifts and rotations of individual coils and modules were compared. With numerous efforts and advanced manufacturing techniques, W7-X superconducting coil system was built and assembled with impressively high accuracy. The maximum deviation for non-planar coils from the average shape is of the order of 2 mm [14] and the average alignment deviation for all 70 main field coils after the assembly is 1.2 mm [5].

More recently, there are several new approaches proposed to address the challenge of tight coil tolerance for stellarators. A shape gradient method was developed by Landreman & Paul [15] to compute local coil tolerance after calculating all the derivatives of a figure of merit with respect to coil parameters. Coil shapes with respect to plasma boundary are analytically differentiated by Hudson et al. [16] and it can be used to find simpler coils. Lobsien et al. [17] adopted a stochastic optimization method to find more robust coils with higher tolerance could by evaluating numerous randomized coil perturbations using a Monte Carlo sampling approach in optimization runs.

Most of the error field sensitivity studies to date are carried out by perturbing (shifting or rotating) coils in certain potential directions and then evaluate the changes in the figure of merit, which is normally the quality of magnetic surfaces, i.e., magnetic islands, as the resonant magnetic perturbation will be amplified, destruct the flux surfaces and degrade plasma confinement. This approach has been successfully applied to several devices. However, it requires massive computation resources, perhaps as well as man-hours, to scan all possible individual coil deviations and compare different combination scenarios. In this paper, we propose a new method to rapidly identify important coil deviations that could possibly appear during coil fabrication and assembly.

The Hessian matrix can be used for sensitivity analysis [18]. In [19], a Hessian matrix method was described to determine error field sensitivity to coil deviations. The figure of merit used was the root-mean-squared (RMS) normal field error on the target plasma boundary, which comes from the coil design code FOCUS [20]. A quadratic approximation indicates coil perturbation in the direction of the eigenvector corresponding to the largest eigenvalue has the most significant effect on the figure of merit. The Hessian matrix method was then applied to a CNT-like configuration as a proof of principle and the results were consistent with previous observations. The RMS normal field error represents the discrepancy between the desired magnetic field and the one produced by coils, but it does not have any particular physics meaning. In practical situations, what should be considered when evaluating the error field is plasma confinement performance. In this paper, we will implement the Hessian matrix method over the magnetic island size, which is used in most error field studies. Additionally, we would also evaluate the so-called ‘quasi-symmetry’ of the magnetic field [21], which is the symmetry of magnetic field strength $|\mathbf{B}|$. Quasi-symmetry has been both theoretically predicted and experimentally confirmed to reduce neoclassical transport and is one of the key qualities used in today’s stellarator optimizations.

This paper is organized as follows. In Sec. 2, we have a brief view of the Hessian matrix method. The two new figures of merit, island size and quasi-symmetry, are also described in Sec. 2. In Sec. 3, we apply the method to an actual device, the Chinese First Quasi-axisymmetric Stellarator (CFQS) [22], which is currently under construction. The results of analyzing the most important islands and evaluations on quasi-axisymmetry (QA) are presented. A method for computing the maximum allowable coil deviation under the worst scheme is discussed in Sec. 4 and a demonstration using the information from the analysis of magnetic islands and quasi-symmetry to improve coil designs is shown in Sec. 5. We will conclude in Sec. 6.

2 Hessian matrix method, magnetic island size and quasi-symmetry

2.1 Hessian matrix method

Coil optimization consists of varying coil parameters, in whatever representation, and minimizing an objective function which is the weighted summation of multiple penalty functions. Once optimal coils are found, a small change in coils shapes (and currents), which can be described as $\delta\mathbf{x}$ in parameter space, will cause a departure in the figure of merit away from the optimum. This change can be approximated by

[TABLE]

where the matrix $\mathbf{H}_{0}$ is known as the Hessian matrix (second-order derivatives). In equation (1), only the quadratic term is left since the linear term is zero at the optimum (stationary point) and we are only considering a small perturbation (neglecting higher-order terms). The Hessian matrix is symmetric and its eigenvalues are positive. By perturbing coils in the direction of eigenvector, $\delta\mathbf{x}=\xi\mathbf{v}_{i}$, the change of the figure of merit is

[TABLE]

Its eigenvectors can be ordered by the magnitude of associated eigenvalues $\lambda_{i}$. The first principal eigenvector, which is the one corresponding to the largest eigenvalue, describes the most sensitive perturbation. With the information in eigenvectors and eigenvalues of the Hessian matrix, we can easily identify the sensitivity of the error field to coil deviations. More details about the Hessian matrix method can be found in [19].

2.2 Magnetic island size

For magnetic fields with a continuously nested family of flux surfaces, small perturbations may break the surfaces. The destroyed magnetic surface has a rational rotational transform ($\iota=n/m$) which is in resonance with perturbations. Magnetic islands are then formed. Plasma confinement across islands is different from within magnetic surfaces and magnetic islands can give rise to plasma instabilities [23]. For this reason, an important figure of merit to measure the quality of the magnetic field is the size of islands, especially in the core.

The width of a magnetic island in toroidal flux space is calculated as [24]

[TABLE]

where $m$ the poloidal mode number, and $\iota^{\prime}_{mn}=\textnormal{d}\iota/\textnormal{d}\psi_{t}|_{\iota=n/m}$. The resonant Fourier component, $b_{mn}$, is decomposed at the associated rational surface,

[TABLE]

In this calculation, we write the magnetic field as $\mathbf{B}=\mathbf{B_{0}}+\delta\mathbf{B}$, where $\mathbf{B_{0}}$ is the nearby magnetic field with perfect magnetic surfaces and $\delta\mathbf{B}$ is a small perturbation. Island width is proportional to the square root of the magnitude of $b_{mn}$, therefore even a small perturbation can generate sizable islands. This is particularly true for low poloidal number islands in low shear configurations.

For a given configuration, one could find the most important island chains that the rotational transform profile will cross. For instance, NCSX cares $\iota=3/7$, $3/6$ and $3/5$ [25]. Furthermore, some configurations might have special requirements on finite-size islands, e.g. W7-X has $\iota=5/5$ island divertors [26]. To manipulate islands, pre-selected resonant components should be maintained to target values, which is to minimize the following functional,

[TABLE]

where $b^{o}_{mn}$ is the target value of $b_{mn}$ and the summation is over the chosen resonant components. To suppress island, $b^{o}_{mn}$ should be zero. Here, we have considered only the effect of the resonant perturbation in equation (3) since the rotational transform profile will not change substantially under small perturbations.

Accurately computing the derivatives of $b_{mn}$ with respect to coil parameters is not easy. As the magnetic field is changed under coil perturbations, the resonant harmonics, $b_{mn}$ in equation (4), varies in several ways: the perturbation field is changed, the flux surface is moved, the straight field line coordinates are reconstructed and the nearby perfect magnetic field is altered. In this paper, we shall only consider small perturbations, by which we assume $\mathbf{B_{0}}$ is not altered and omit the motion of flux surfaces. Hence, only the change in the perturbation field $\delta\mathbf{B}$ will be considered under small coil deformations. This linear approximation is used in W7-X error field studies [11, 12, 13, 14, 27] and on NCSX [9]. Rational surface information, the straight field line coordinates $(\psi_{t},\theta,\zeta)$ and the normal vector $\nabla{\psi_{t}}$, could be obtained by the so-called “quadratic-flux-minimizing surface” [28] or approximated by a closest flux surface in free-boundary MHD equilibrium calculations (like VMEC [29]). Derivatives of $b_{mn}$ are now calculated as

[TABLE]

where $X_{i}$ and $X_{j}$ are arbitrary terms in the coil parameter vector $\mathbf{X}$. In equation (6) and equation (7), we are using the linear approximation and have considered only the changes in the perturbation field. Derivatives of $F_{RP}$ could be then obtained. The first- and second-order functional derivatives of the magnetic field produced by external coils are calculated in FOCUS (Eq. (A.4) & (A.5) in [19]). We can then use this information to rapidly compute the Hessian matrix of the magnetic island metric.

2.3 Quasi-symmetry

The guiding center motion of collisionless particles was found to be governed by the magnetic field strength $B$ alone [30] in Boozer coordinates. It implies that configurations with symmetry in $B$ will have good neoclassical transport, although the magnetic field itself is not symmetric. The condition that $B$ only depends on one helicity is called ‘quasi-symmetry’. This leads to an important category of optimized stellarator designs, like quasi-axisymmetric stellarators [31, 32, 33], and quasi-helical stellarators [21, 34, 35].

The quality of quasi-symmetry (on one flux surface) can be evaluated using

[TABLE]

where $N,M$ are the target Fourier modes to be conserved, e.g. $N=0$ for quasi-axisymmetry (QA). Note that $F_{QS}=0$ indicates perfect quasi-symmetry of the magnetic field on the reference flux surface. The Fourier component of magnetic field strength, $B_{m,n}$, is then decomposed in Boozer angles $(\theta_{B},\zeta_{B})$,

[TABLE]

Here, $N_{p}$ is the number of periodicity. $B_{m,n}$ depends on the magnetic field $\mathbf{B}$, the flux surface and the Boozer coordinates on the surface. When calculating the derivatives of $B_{m,n}$, a finite difference method can be used. Each evaluation involves preparing the perturbed magnetic field, solving the free-boundary MHD equilibrium, constructing Boozer coordinates and calculating $B_{m,n}$. To calculate the Hessian matrix, the finite difference method might be computationally expensive especially when the dimension of the Hessian matrix is large. Instead of using the finite difference, in this paper we will apply the same linear approximation as the magnetic island to quickly construct the Hessian matrix.

If we apply the above linear approximation and consider only the change in the magnetic field $\mathbf{B}$, $B_{m,n}$ derivatives can be calculated as in equation (6) and equation (7) (they are not written down here for brevity). Consequently, the derivatives of $F_{QS}$ are computed as,

[TABLE]

[TABLE]

These equations are evaluated in FOCUS and the Hessian matrix of the quasi-symmetry metric could be computed.

3 Numerical application on CFQS

3.1 Brief introduction of CFQS and its coils

CFQS is a quasi-axisymmetric stellarator being built in China under the collaboration of Southwest Jiaotong University in China and the National Institute for Fusion Science in Japan. The main parameters of CFQS are as follows: the toroidal period number $N_{p}=2$, major radius $R_{0}=1.0$ m, minor radius $a=0.25$ m and magnetic field strength $B_{t}=1.0$ T [22]. Sixteen modular coils are designed to provide a relatively accurate magnetic field along with adequate coil intervals and acceptable coil curvature [36]. The target plasma boundary of CFQS and modular coils are shown in figure 1. Because of stellarator symmetry and periodicity, there are only four unique coil shapes, denoted as M1, M2, M3 and M4. Among the four coil types, the M4 is considered to be the most complex one, as it has the largest toroidal excursions.

Although multiple scenarios have been investigated with different beta values [37], CFQS and its modular coil system were optimized by targeting a zero-beta (vacuum) configuration. For simplicity, we only consider the vacuum configuration for the following calculations. The target plasma boundary at the bullet-shaped cross-section and Poincar$\acute{\rm e}$ plots by tracing field lines in the magnetic field produced by the designed coils are shown in figure 2. It should be noted that the magnetic field is calculated from single filaments located at the center of each coil. In this paper, we will ignore the effect of finite-build coils and only use coil filaments for the following calculation (although coils are plotted with finite size). Note that the target boundary is not exactly matched in figure 2. There are two reasons. First, the desired magnetic field is not the only target when designing CFQS coils. The average residual $\mathbf{B}\cdot\mathbf{n}/|\mathbf{B}|$ on the target boundary is about 1%. Secondly, the coils used in this paper are slightly different from the original ones, as FOCUS does not directly use the raw points in space. FOCUS adopts Fourier representation to describe coils and truncation errors exist when fitting the actual coil data. Hereafter, CFQS (modular) coils are referred to single filaments described with Fourier representation in FOCUS.

3.2 n/m=4/11 magnetic islands

CFQS is a low shear configuration. The rotational transform profile was carefully chosen to avoid low-poloidal-number islands. Rotational transform profiles from VMEC free-boundary calculation with the designed coils and from field line tracing are shown in figure 3. The $n/m=4/11$ island chain is found to have the lowest poloidal number inside the boundary. From figure 2, one can clearly observe the 4/11 islands, although their sizes are small. We shall select $b_{11,4}$ as our reference resonant perturbation and the present value of $b_{11,4}$ is served as the target value $b^{o}_{11,4}$. The metric described in equation (5) is now evaluating how far the 4/11 islands are away from the present size. With the designed coils, $F_{RP}=0$ is at a global minimum. The 4/11 island chain might not be serious (because of relatively large $m$), but it would be sufficient for demonstrating the new Hessian matrix method.

To approximate the $\iota=4/11$ rational surface, free-boundary VMEC was executed with relatively high radial resolution (256 surfaces with maximum Fourier modes Mpol=10, Ntor=20). Afterwards, the flux surface that is closest to $\iota=4/11$ surface was selected and the code BOOZ_XFORM [38] was used to convert VMEC coordinates into Boozer coordinates. FOCUS then read the reference flux surface parameterized in Boozer angles ($\theta_{B},\zeta_{B}$). The magnetic field $\mathbf{B}_{0}$ was calculated from the actual coils, as the reference configuration is zero-beta. The resonant harmonic $b_{11,4}$ and its derivatives (both first- and second-order ones) were then calculated by following equation (4), equation (6) and equation (7), where we fixed the flux surface and the Boozer coordinates on it.

When using Fourier representation in FOCUS, the vector of coil parameters consists of Fourier coefficients and coil currents,

[TABLE]

where $x^{c}_{k,n}$ ($x^{s}_{k,n}$) is the $n$-th cosine (sine) Fourier coefficient of the $k$-th coil and $I_{k}$ is the current. To accurately describe CFQS modular coils, we have performed a scan of the number of Fourier coefficients used in each spatial coordinate. According to appendix A, the optimal value of the maximum Fourier mode $N_{F}$ is 10. The total number of degrees of freedom is 1024 since each modular coil is considered to be independent (as they will be fabricated and assembled separately). Coil currents are normalized to the average current, $3.125\times 10^{5}$ A, while Fourier coefficients are normalized to a geometric quantity, $1.23$ m.

The Hessian matrix of the 4/11 islands is then computed by FOCUS. (The computation time is about 100 seconds with 64 cores.) Eigenvalue decomposition is carried out and figure 4 shows the spectrum of eigenvalues in descending order. Only the first (largest) eigenvalue has a significant value and all the others are almost zero. This implies that the Hessian matrix is positive semi-definite, which is no surprise as $F_{RP}$ is at an exact global minimum. The first principal eigenvector, which is associated with the largest eigenvalue, will have the highest impact on the size of 4/11 islands. It represents the most important perturbation scheme in all the possible combinations with the same magnitude in parameter space. figure 5 illustrates how the coils will deform if perturbed in the direction of the first principal eigenvector with a finite magnitude ($\xi=0.01$). It is interesting that the first principal eigenvector is preserving periodicity and stellarator symmetry, although all the 16 modular coils are considered independently. Coil currents are also free variables and the currents will be perturbed by amounts of $\xi\times(-2160.67,-1949.13,1877.62,2529.56)$ A in each type of coil, which are negligible compared to the actual coil current ($3.125\times 10^{5}$A). Results in figure 5 suggest that the M1 & M2 coils are much more sensitive than M3 & M4. Particularly, toroidal deformations at three regions, inboard upper corner, inboard lower corner of M1 and inboard part under the middle plane of M2, have the most significant effects on preserving 4/11 islands.

One can verify the effects of coil deviations by checking the size of 4/11 islands. Field line tracings are performed with the perturbed coils. Poincar$\acute{\rm e}$ plot of the negatively perturbed coils ($\mathbf{X}_{0}-\xi\mathbf{v}_{1}$) is shown in figure 6(a) and the positively perturbed one ($\mathbf{X}_{0}+\xi\mathbf{v}_{1}$) is in figure 6(b), while $\xi=0.01$ in both cases. As shown in figure 6, the 4/11 islands are either suppressed or enlarged when perturbing the coils in the two directions. Both eliminations and enlargements of the islands are deteriorating our figure of merit since we provide a non-zero $b^{o}_{11,4}$ in equation (5) and it will preserve a finite island chain.

3.3 Quasi-axisymmetry

Similarly, quasi-axisymmetry of the magnetic field could be our figure of merit. By assessing the non-axisymmetric terms ($n\neq 0$) in equation (8), we could evaluate the QA quality of the produced magnetic field. figure 7 shows the spectra of main non-zero non-axisymmetric terms from BOOZ_XFORM converting free-boundary VMEC results with the original CFQS coils. Because the flux surface is stellarator symmetric, odd (sin) components in $B_{m,n}$ are zero and hereafter we only consider even (cosine) components. The leading non-axisymmetric terms, $B_{0,1}$, $B_{1,1}$, $B_{1,2}$ and $B_{2,2}$, are significantly smaller than $B_{0,0}$ ($\sim 0.95$T). QA quality, as evaluated in equation (8), monotonically degrades from the axis to edge. We can select an arbitrary flux surface, e.g. the flux surface at $s=\Psi/\Psi_{edge}=0.5$, as the reference surface. After reading in flux surface information and Boozer coordinates, FOCUS calculates the magnetic field produced by coils on the reference surface and then Fourier coefficients $B_{m,n}$ are gauged by using equation (9). The maximum poloidal and toroidal mode numbers used here are Mpol=4, Ntor=8 ( $m\in[0,\mathrm{Mpol}]$ and $n\in[-\mathrm{Ntor},\mathrm{Ntor}]$).

Eigenvalues of the Hessian matrix regarding QA quality are shown in figure 8. Numerical details are the same as in Sec. 3.2. Unlike the 4/11 island size, the eigenvalues for QA quality on the half-toroidal-flux surface are more divergent. The four largest eigenvalues ($\lambda_{1}$, $\lambda_{2}$, $\lambda_{3}$, $\lambda_{4}$) are of values ($1.07$, $0.92$, $0.56$, $0.28$), while the smallest one $\lambda_{1024}=-6.42\times 10^{-2}$. Actually, there are 515 negative eigenvalues. This implies that the Hessian matrix is not positive semi-definite and the figure of merit (QA quality) used here is not at a local minimum. It is expected since the coils are not designed by targeting quasi-axisymmetry (on this particular flux surface) and there is a trade-off for finding easy-to-build coils. However, the quadratic approximation in equation (1) is still valid, as checked in appendix B. QA quality is marginally close to a minimum. Therefore, we can keep using the Hessian matrix method. Here, the target value for QA quality was set to be zero, i.e., CFQS would like to have as high QA quality as possible. This is not essential and we could modify it to include a non-zero target value like in equation (5).

Likewise, the effects of eigenvectors on coils can be studied by mapping the perturbed coils back to real space. This time, we need to consider more than one eigenvector. In figure 9, how coils deform under the first four principal eigenvectors with the same size perturbation $\xi=0.01$ (in positive directions) are displayed. Currents flowing in the coils are also varied but within a small range (less than 0.09% of the original current). Again, perturbed coils are still following stellarator symmetry and periodicity, so only one fourth of the coils are shown. The maximum departures in the four schemes are not identical, but close. Distributions of highly deformed regions are also different.

4 Coil tolerance calculation

Once the eigenvalues and eigenvectors of the Hessian matrix have been obtained, coil tolerance can be computed. The change in the target figure of merit is quadratically approximated, as shown in equation (1). Suppose we have an arbitrary perturbation $\displaystyle\Delta\mathbf{X}=\sum_{i}^{N}a_{i}\mathbf{v}_{i}$ where $a_{i}$ is the $i$-th component in the direction of $\mathbf{v}_{i}$. The resulting error field can be approximated by $\displaystyle\Delta F\approx\sum_{i}^{N}\frac{1}{2}\lambda_{i}{a_{i}}^{2}$. Since the eigenvalues are in descending order ($\lambda_{1}\geq\lambda_{2}\geq\cdots\geq\lambda_{N}$), for all the perturbations that have the same magnitude ($\displaystyle|\Delta\mathbf{X}|=\sqrt{\sum_{i}^{N}{a_{i}}^{2}}=\xi$), the one that only consists of the first principal eigenvector ($\Delta\mathbf{X}=\xi\mathbf{v}_{1}$) will generate the largest error field ($\Delta F=\frac{1}{2}\lambda_{1}\xi^{2}$).

If we have a maximum allowable deterioration on the figure of merit, i.e., ${\Delta F}_{tol}$, the associated tolerance of the perturbation magnitude in the direction of the first principal eigenvector is calculated as,

[TABLE]

Here, the perturbation magnitude is calculated in parameter space (e.g. Fourier space). More useful information for engineers would be the actual coil tolerance in real space (values in $mm$). We can obtain the worst coil deviation scheme and the allowable perturbation magnitude in real space by mapping the perturbed coils back to the real space. Discrepancies on values after the transformation are expected, especially when we use Fourier representation for coils. Other local representations, like piece-wise linear representation and splines, might be able to resolve this discrepancy. In this paper, we are going to stay with Fourier representation and conduct some initial attempts to calculate coil tolerance in real space.

For CFQS, as stated above, the rotational transform profile is carefully chosen to avoid low-rational-number islands. Therefore, resonant magnetic islands are not the first priority for the error field study. Meanwhile, the primary goal of the device is to demonstrate the properties of quasi-axisymmetry. We shall only use the QA quality on the half-toroidal-flux surface to calculate coil tolerances.

If we assume the maximum acceptable QA deterioration is 10% of the present value ($F_{QA}=2.95\times 10^{-4}$), the magnitude of the allowable perturbation calculated using equation (13) in the first principal eigenvector is $7.42\times 10^{-3}$. Note that for perturbations having the same magnitudes in Fourier space , coil deviations in real space are not of the same amount (albeit they are of the same order), as shown in figure 9. In this case, the magnitudes of perturbation in the other three principal eigenvectors are $8.02\times 10^{-3}$, $1.03\times 10^{-2}$ and $1.45\times 10^{-2}$. Since the differences of the four largest eigenvalues are not substantial, we will try to compare all the four associated eigenvectors. If we only consider perturbations in the direction of eigenvectors, the allowable coil deviations for the first four eigenvectors are 5.74, 4.57, 5.54 and 10.02 mm. Detailed coil deformations can be visualized in figure 9 (using allowable perturbation magnitudes). A maximum allowable deviation of about 5 mm is much more generous than coil tolerance in NCSX and W7-X, especially considering that CFQS has a smaller size.

Here, the value of allowable 10% of QA quality deterioration is chosen just for demonstration. Serious numbers can be assessed by comprehensive calculations, like evaluating the $\epsilon_{eff}$ or do particle transport simulations, and it is beyond the scope of this paper. In the above calculations, we were evaluating the first four principal eigenvectors and the contributions from coil currents were all small. For obtaining coil tolerance in real space, coil currents should probably be excluded from free coil parameters.

5 Improve coil designs towards better physics properties

The two figures of merit mentioned in this paper, magnetic islands and quasi-symmetry, are two important physics properties. In addition to quantifying error fields, they can be used to optimize coils for better physics properties. For instance, the 4/11 islands in CFQS can be eliminated if the coil shapes are varied following a certain direction, as shown in figure 6(a). More generally, specific islands could be reduced if we minimize equation (5) by targeting zero island width. This idea was successfully applied to NCSX by Hudson et al. [25].

Here, we will show that we can also improve the quasi-axisymmetry quality of CFQS. As shown in figure 8, there are several negative eigenvalues of the Hessian matrix. If the coils are varied following the eigenvectors associated with negative eigenvalues, the change of the functional in equation (8) will be negative as well. This indicates that the QA quality on the half-toroidal-flux surface will be improved. Instead of varying coils manually, we can directly optimize the QA quality as one of the cost functions in FOCUS. A simple optimization was carried out by targeting only the cost function of QA quality on the half-toroidal-flux surface. The actual CFQS coils were used as the initial guess. After 10 iterations using a conjugate gradient method, a new coil set was found. $F_{QA}$ was reduced from $2.95\times 10^{-4}$ to $4.42\times 10^{-5}$ on the half-toroidal-flux surface. Free_boundary VMEC with the new coils and BOOZ_XFORM calculations were then conducted. In figure 10, we show the leading non-axisymmetric terms from the optimized coils compared with the original coils. All the main leading non-axisymmetric terms, except $B_{0,1}$, have been reduced. It indicates that the optimized coils have a better quasi-axisymmetric magnetic field across the whole plasma, although we are trying to improve QA quality on only one flux surface. The optimized coils are shown in figure 11. Compared to the original coils, the optimized ones do not have unrealistic increase in geometry complexity (just naively comparing by eye). Note that this is not a comprehensive optimization and the only objective function, $F_{QS}$ is not at a local minimum after 10 iterations. We are also not attempting to compare which coil set is better. The results demonstrate that the QA quality could be directly optimized during coil designs. In future studies, we could actually employ other powerful optimization algorithms, like the modified Newton method (MN) [39] since the gradient and Hessian have been analytically calculated.

6 Conclusions

The recently developed Hessian matrix method for analyzing error field sensitivity has been extended to identify important coil deviations that will open magnetic islands or degrade quasi-symmetry. To manipulate islands, a new functional that evaluates resonant perturbations is incorporated. The quasi-symmetry quality on one flux surface is assessed by computing Fourier components in Boozer coordinates. With the two new figures of merit for quantifying error fields, we have applied the Hessian matrix method to the CFQS device which is under construction. We have successfully found the most effective coil perturbation scheme that will enlarge the 4/11 islands, as shown in figure 6. In addition, the important coil perturbations deteriorating the quasi-axisymmetry on the half-toroidal-flux surface are also identified (figure 9). These results will provide insights for the upcoming coil manufacture and assembly: avoiding critical coil deformations and relaxing coil tolerance at insensitive parts. A simple estimation indicates that the maximum allowable coil deviation is about 5 mm if 10% of QA quality deterioration is acceptable. We have also demonstrated the ability to improve coil designs towards better QA quality. This could be used in future coil optimization studies.

In this paper, each coil is considered independently, but the leading eigenvectors shown above are surprisingly preserving stellarator symmetry and periodicity. This is not true for all the eigenvectors and some eigenvectors with smaller eigenvalues do not retain stellarator symmetry or periodicity. There is a possibility that the eigenvectors preserving stellarator symmetry (and/or periodicity) will be more effective than others if the figure of merit and present coils are both symmetric (and/or periodic). Rigorous mathematical explorations are remaining for future work. Here, we only employ Fourier representation to describe coil geometries. In general, any other representation could be used. Some representation might be able to resolve the discrepancy between the magnitude in parameter space and in real space as displayed in Sec. 4. So far we have only considered coil filaments. When coils are not too close to the plasma, this approximation is sufficient, as the magnetic field produced by coils is proportional to $r^{-3}$. However, coils are eventually built with finite thickness. A finite-build model would be still essential. We can use the Hessian matrix method as a pre-processing step since it is fast and the sensitivity information could provide valuable guidance for detailed engineering analyses. To analytically calculate the derivatives of resonant perturbations and quasi-symmetry, we assume the flux surface and the straight field-line coordinates system do not move under small perturbations. The fact that we successfully manipulated the island size (figure 6) proves this simplification is valid here. Additionally, as shown in Sec. 5, the global QA quality was improved by minimizing $F_{QS}$ in equation (8) using the same assumption. This linear approximation is not essential and might not be as accurate as nonlinear analysis since in reality flux surfaces and coordinates systems are likely changing under coil perturbations. By using the finite difference to calculate the Hessian matrix, the assumption can be relaxed and all the calculations would be self-consistent, although the finite difference suffers numerical noises and takes much more time. Careful validations on the linear approximation should be carried out in the future work. While we have discussed using the Hessian matrix method for only studying error fields from stellarator coils in vacuum configurations, it could also be applied to equilibria with finite plasma currents, or even in tokamaks and other general fields to obtain sensitivity information.

Acknowledgments

The authors gratefully appreciate fruitful discussions with A. Brooks, Y. Suzuki, M. Landreman S. Lazerson and N. Pomphrey. This work was supported by the U.S. Department of Energy under Contract No. DE-AC02-09CH11466 through the Princeton Plasma Physics Laboratory and the Max-Planck-Princeton Center for Fusion and Astro. Plasma Physics.HL and YX would like to acknowledge the National Science Foundation of China under Grant No. 11820101004 and the National Key R&D Program of China under Grant No. 2017YFE0301705.

Appendix A Optimal Fourier resolution for CFQS modular coils

Each CFQS modular coil is piece-wisely described by 48 points in space. However, FOCUS employs Fourier representation. We need to fit CFQS coil data with Fourier coefficients. First of all, we used cubic spline interpolation to smooth coils and each coil now has 256 points. Afterwards, FOCUS read coil data, performed Fourier transformations and calculated the residual $\mathbf{B}\cdot\mathbf{n}$ error on the last closed flux surface (obtained from free-boundary VMEC runs). figure 12 shows how the residual normal field error varies when we increase the resolution of Fourier coefficients. $N_{F}$ is the maximum Fourier mode truncated in each coil. The residual error is exponentially reduced with increasing $N_{F}$.Since $N_{F}=10$ has the minimum error, we will use $N_{F}=10$ for future calculations. Then the total number of degrees of freedom, as listed in equation (12), is $1024=16\times(6\times 10+3+1)$.

Actually, the number of required Fourier coefficients to adequately describe CFQS coils should be smaller than 1024. The spectrum of Fourier modes can be condensed by carefully choosing the parameterization angle. Besides, some type of coil with simpler geometry might need less Fourier coefficients. But we are not going to explore these issues in this paper.

Appendix B Quadratic approximation of QA quality when the gradient is non-zero

In Sec. 3.3, the Hessian matrix for QA quality on the half-toroidal-flux surface is not positive definite. It implies the figure of merit is not at a strict local minimum with CFQS coils. However, the quadratic approximation assumes the gradient is negligible. Therefore, we need to test if equation (1) is still valid.

Choose a perturbation direction, e.g. the first principal eigenvector of the Hessian matrix, apply different magnitudes of perturbation and then use FOCUS to compute the QA quality with each perturbed coil set. figure 13 shows the relationship between changes in the QA quality $\Delta F_{QS}$ and the perturbation magnitudes $\xi$ in both the first principal eigenvector $\mathbf{v}_{1}$ and the last principal eigenvector $\mathbf{v}_{1024}$. The two lines are consistent with the predicted quadratic lines, even in a relatively sizable range. The results imply that the QA quality is marginally close to a local minimum and the quadratic approximation is still valid.