Dynamic error fields derivation by inverting a validated interpretative perturbations model

The need for a model able to provide the dynamic evolution of the intrinsic magnetic error fields amplitude in tokamaks is of a significant importance in order to develop strategies to correct its destabilizing effect. The present paper specifically aims to deliver this kind of model and adjacent code. A previously built direct perturbations interpretative model (Miron (JET Contributors) 2021 Nucl. Fusion 61 106016) is used to calculate the plasma instabilities amplitude and frequency aiming to match the same quantities provided by the device diagnostics data analysis specific code for various discharges at JET. The mentioned good match ensures our model validity against the experimental results. Based on its proven structural validity, our model considering the error field quantities is inverted. This time the initial direct model output results, namely the plasma perturbations amplitude and frequency, become input data for our inverse model aiming to solve the corresponding error field modes system of equations searched as unknowns. The error fields basically satisfy the same outer plasma perturbed equations as the plasma perturbations do (vacuum/wall/plasma column external structures). Obviously they do not satisfy the perturbed plasma equations. A clear and explicit dynamic error field solution is finally provided. It has been demonstrated that whenever the error fields are responsible for the mode locking effect and plasma rotation damping, the calculated error fields shows its expected disruptive resonant behavior during various discharges at JET.


Introduction
Inherent misalignments or imperfections of the tokamaks poloidal field coils give rise to nonaxisymmetric, small (of 10 −5 -10 −4 order of the equilibrium magnetic field) magnetic error fields. The error fields affect the plasma stability the way any external perturbations do by breaking the toroidal symmetry of the equilibrium magnetic field through a resonant as well as a nonresonant coupling mechanism [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20] with the local plasma perturbations that finally drives to the occurrence of braking torques such as the electromagnetic [21][22][23] or neoclassical viscous torque [24][25][26][27][28][29] that are subsequently responsible for the deceleration or the drop of the plasma rotation that usually leads, if no countermeasures are taken, to the disruption of the plasma. The need to deliver a model that is able to provide the dynamics of the error field even under some used approximations is justified.
A key role for our goal is played by the perturbations model from [30] this paper relies upon. The mentioned model has been validated at JET regarding the quasi-analytic retrieval of the measured amplitude and frequency of the involved perturbations during various discharges at JET. This achievement has provided an explicit time dependent formula of the dynamic perturbations behavior. On the other hand this direct model has not explicitly considered the error field within our calculations. The eventually considered error field effect has been successfully embedded into the whole perturbation behavior especially by means of the perturbations initially conditions that have been chosen and the used plasma quantities diagnostics data. The error fields did not need to be explicitly defined at that time in order to preserve the validity of the model. Remember that the model from [30] is an interpretative, not a predictive model. The situation obviously changes when the error fields are the very quantities to be derived. This time the error fields have to be defined and inserted into the perturbed system of equations of a modified inverse model whose perturbed plasma solutions are different compared to the ones of the initial direct validated model we have started with. Fortunately, within the process of deriving of the error field solutions of the modified inverse model all we really need is the plasma perturbations amplitude evolution as input data for our new error fields system of equations. This quantity is provided by the JET data analysis tables. But the use of such a data table would drive our approach fully numerical instead of a quasi-analytic approach. The answer is to use an analytic expression that matches as good as possible the plasma perturbation amplitude data. Now it is the time our plasma perturbation solution of the initial direct validated model is used. It does not matter that it has a different formula compared to the similar solution of the inverse model where the error fields are explicitly taken into account that should also match the same experimental data. The only thing that matters is to match as good as possible the data analysis amplitude evolution. This is the way we rely upon the direct model from [30] and its plasma perturbation solution in order to invert and modify it to derive the explicitly introduced magnetic error fields. The agreement level between the experimental and the modelled results we expect to obtain is simply a reasonable one. A theoretical (mathematical) model is used, based on many justified physical approaches in order to be able to advance with the calculus starting from the beginning equations. Because of all these approaches it is expected to not deliver sharply accurate results, but good results that reasonably match the shape and the scale of the experimental ones.
We have organized the paper as follows: In section 2 the final form of the complete system of the perturbed equations within the whole plasma/vacuum/external structures space are recalled from [30]. The geometrical approach is described. In section 3 the magnetic error fields are introduced and the corresponding system of equations is derived. Section 4 is devoted to the quasi-analytic derivation of the error field solution as a function of the perturbed magnetic flux. The calculated results and the comparison against the data analysis ones are shown in section 5. Section 6 is dedicated to the conclusions.

Perturbed equations
In order to derive the magnetic error field modes associated to the perturbations we are starting with the complete perturbed system of equations providing the dynamics of the perturbations from (49) and (79) of [30]: (1) encloses the ideal plasma perturbed system of equations whereas (2) provides the outer plasma (including various external resistive structures) perturbed equations, both of them localized at a magnetic surface of radial coordinate r s where the perturbation is developing. The system of equation (1) is the result of the linearization of the perturbed ideal plasma momentum equation and of the adiabatic equation of state that are subsequently expanded in Fourier harmonics with j, m and k, n the poloidal and the toroidal numbers, respectively. The linearized system (2) is the summation of all the perturbed vacuum equations crossing the plasma outer inhomogeneously resistive structures such as the external wall and coils system. P jk mni and Q jk mni gather the contributions of the perturbed magnetic force, inertial term, perturbed pressure gradient and perturbed parallel stress tensor on the perturbed plasma moving perpendicular and tangential to the corresponding magnetic surface, respectively. In a similar manner F jk mni and G jk mni collect the external plasma structures contributions to the local perturbed plasma moving in the above mentioned directions. All these coefficients are defined in (50), (80) and (81) of [30]. The local perturbation is searched for as Φ jk s where −∂Φ jk s /∂t is the Fourier term of the perturbed electric scalar potential at r s , from the displacement parametrization Φ ′jk s ≡ ∂Φ jk s /∂r is the radial derivative of the perturbed magnetic flux Φ jk s . ξ || is the displacement parallel to the magnetic field whose unit vector isn = B/B. The perturbed magnetic field is b = ∇ × (ξ × B). E mn (t) is found in (60) of [30] providing the plasma column external coils active contribution. Local flux coordinates (r, θ, φ) have been used related to the cylindrical polar coordinate system (R, φ, Z) by within the low inverse aspect ratio approximation. ! means factorial and ′ is the partial derivative ∂/∂r. The notation comprises the Shafranov shift aϵ∆(r), ellipticity aϵE(r) and triangularity aϵT(r) shaping parameters where a low plasma inverse aspect ratio ϵ = a/R 0 << 1 is considered, a and R 0 being the tokamak plasma minor and major radius, respectively. δ i,j is the Kronecker delta. The shaping parameters satisfy (7) of [30] The JET diagnostic data tables usually provides the shaping parameters quantities either at the boundary or at the magnetic axis. Therefore in order to derive it at the local magnetic surface where the perturbation occurs, we need to solve the corresponding differential boundary equation (8) listed above.

Introduction of the error field
The system of equation (2) basically describes the perturbation evolution in vacuum and plasma external structures. The perturbed quantities Φ jk s and Φ ′jk s are localized at the vicinity of the inertial layer or of the magnetic island lower boundary. It has been assumed that the equilibrium plasma current is mostly concentrated inside the mode rational surface and is negligible for r > r s [31]. Therefore the vacuum is basically prolonged until the perturbation inertial layer. The system (2) is the result of the outer plasma equations coupled with the perturbed boundary equation across the inertial layer. If we parametrize the error field in the absence of the plasma in the same manner as the perturbation by means of a vacuum flux Ψ EF , equation (72) of [30] showing the response of the plasma in the presence of a resonant error field becomes [21] s relates the both parametrizations and is valid within the low inverse aspect ratio approximation. Ψ mn s is the perturbed flux in the presence of both the error field and plasma. A similar but formal parametrization is also used for the error field magnetic flux. The inductive effect of the error field has been neglected in (9), the mode coupling effect being our primary interest. By introducing the error field using (57), (61), (67), (68), (71), (78) of [30] and (9) q s and s s are the safety factor and the shear at r s . β jkmn i is found in (82) of [30]. With the notation of an approximate fully jump across the layer (not only the real part) Equation (10) now is where (13) is our error field system of equations to be solved.

Error field solution
We intend to solve the system of the error field equation (13) by using our validated model results at JET [30]. This time the perturbed magnetic fluxes associated to the perturbations, Φ jk s and Φ ′jk s , play a parametric role for (13). Ideally the latter quantities should be taken from the JET mode analysis MHD python code [32] perturbations amplitude results but the solving of the system (13) would become a pure numerical matter. By relying on our direct model validity, our model matching perturbed solutions will be used instead. This approach allows us to obtain a quasi-analytic solution for the involved magnetic error fields.
We simply search the solution as a sum of a general solution of the homogenous system of equations and a particular solution of the non-homogeneous one. The solution of the homogenous system associated to (13) simply relies on solving the eigenvalues problem for the 2L × 2L matrix O L and I L are the L × L null and identity matrices, respectively. The homogenous solution is where m and n span the m 1 < m < m 2 and n 1 < n < n 2 poloidal and toroidal number intervals, respectively. λ i and v i are the eigenvalues and the eigenfunctions of C, respectively.ĉ i are constants to be found from the initial conditions. We are looking for a particular solution of (13) so that having the matrices Ω MP is the external coils signal rotation frequency and E DC and E AC are the static and rotational contributions of the coils current intensities defined in (60) of [30]. This choice allows us to turn (13) into with the particular solution Λ jk p to be found. As in the [30] case the assumption of considering the diagnostic quantities almost static on the perturbations time scale has been used as a fair approximation, the growth rates of the former being significantly lower compared to the ones of the perturbations themselves. By using this approach the coefficients of the system of the perturbed equation (20) bypass the time derivation or the Laplace integration. This assumption keeps the linear approach and allows us to provide a quasi-analytic solution for the magnetic error field mode. After Laplace transforming (20) we are choosing Λ jk p as to satisfy where the contribution of all the terms containing the initial value of the involved perturbed fluxes and its subsequent time derivatives following the Laplace transform operation is not considered. This particular choice also provides a solution for (20). The circumflex sign indicates the Laplace transformf(τ ) ≡ L(f(t)). At this point we make use of our validated model [30] in order to replace the perturbed flux Φ jk s from (21) with our obtained solution (93) of [30] whose amplitude and frequency are expected to match the ones provided by the JET mode analysis MHD python code and the mode spectrogram, respectively. The solution of the Laplace transformed equation (21) isΛ and where δ ab is the Kronecker symbol, sgn is the sign of the L index permutations and {τ p } p=1,...,7L are the roots of the ∆ = 0 equation from (96) of [30]. If we carefully look at (22) we observe that the order of the numerator polynomial in τ is equal to the one of the denominator. In order to use the partial fraction decomposition method that allows us to easily apply the inverse Laplace transform and find the solution, we need a higher order of the denominator polynomial compared to numerator. An extra step is needed, the use of long division of fractions. A straightforward calculus leads to By defining the operator we get and The used coefficients are and c l 2L and c 2L are the polynomial coefficient of the highest degree in τ of D l and D respectively. We get rid of the the first term of (26) (not dependent of τ ) and choose our particular solution accordingly. This choice also provides a solution for (20). The rest of the solution contains terms whose order of the numerator polynomial in τ is at least one degree lower compared to the denominator one. Now (26) can be inverse Laplace transformed. With the notations {s i } i =1,...,2L for the 2L roots of the D(τ ) = 0 equation, our chosen particular solution becomes By using (16), (17)- (19) and (35) the general error field solution is finally with the {ĉ i } i =1,...,2L constants only left to be derived from the initial conditions. A suitable assumption is to consider the only presence of the error field at t = 0. Recall the one-toone inverse relations between l and the mode numbers (m, n) [33], namely m = l + m 1 − 1 + (m 2 − m 1 + 1)[(l − 1)/(m 2 − m 1 + 1)] and n = n 1 + [(l − 1)/(m 2 − m 1 + 1)]. The square bracket here means the integer part.

Results
The above formula is ready to be tested for various discharges involving locked modes and magnetic error fields at JET. Although the error field penetration effect could be significantly due to the onset of an neoclassical torque that acts globally following the nonresonant mode coupling mechanism, we limit ourselves in this paper to the resonant approximation for an obvious desired easiness of our calculus.
The experimental data at JET has been used in the same manner as in [30], namely using the EFIT equilibrium reconstruction constrained by polarimetry measurements for the plasma safety factor, shear, pressure, pressure radial derivative, toroidal magnetic flux, plasma minor and major radius, poloidal beta, plasma boundary ellipticity and triangularity and the poloidal magnetic flux on the magnetic axis and at the plasma boundary, respectively. The local values of the above first five spatially distributed dynamic quantities are derived at the radial position r s where the perturbation occurs, where r s is found either using the spatial distribution of the safety factor values in the vicinity of the local magnetic surface or from the perturbation magnetic surface major radius provided by the JET mode analysis code, when available. The charge exchange recombination spectroscopy (CXRS) diagnostics provides the toroidal plasma rotation and the ions temperature. The density and temperature of the electrons are obtained via Thomson scattered laser light (HRTS) and the spectroscopy diagnostics using Bremstrahlung signal along the vertical channel (ZEFV) through the plasma delivers the effective ionic charge quantity. A first example of how the calculated error field looks like corresponds to the JET shot no. 99 902 (toroidal magnetic field B T = 3.28 T) where a locked mode occurence has been reported. The overview data for this shot is depicted in figure 1. The neutral beam injection (NBI) has a single step evolution whereas the ion cyclotron resonance heating (ICRH) is kept constant during the discharge ( figure 1(c)). The mode frequency spectrogram is shown in figure 1(d) [32]. The shaded areas in the figure are accompanied by a drop of the modes frequencies, normalized beta and of the absolute value of the plasma current indicating an unstable plasma behavior. The mode of the n = 1 toroidal number is analysed in figure 2. As it has been pointed out earlier, as a first step, we rely on the best obtained match between our calculated perturbed poloidal magnetic field amplitude at the level of the JET fast magnetic aquisition system diagnostic coils (of radial coordinate r f ) and the similar quantity provided by the JET mode analysis tool [32] in order to test the validity of our direct model used in [30]. This good match between the experimental (data analysis) and the calculated mode amplitude is clearly visible in figure 2(a). Based on this match we are confident to deliver a reliable corresponding error field amplitude introduced by (9) into the inverse model.
The mode locks and then unlocks within areas centred around 10.2 s and 12.8 s, respectively. The mode amplitude b mn θf = m(r m s /r m+1 f )|Ψ mn s | shows an increase during the locking intervals. No CXRS data is available for the 99 902 pulse thus no plasma toroidal rotation data is provided.
The local toroidal rotation velocity has been calculated using the mode frequency data instead, shown in figure 2(b), as an acceptable approximation. The n = 1 error field mode in figure 2(a), calculated using (36), increases during locking and seems to be enhanced by its corresponding plasma mode, most probably due to the resonant field amplification effect. It is thus probably safe to say that the error field is the source of the mode locking effect.
A similar reported n = 1 locked mode is also tested for the JET shot no. 77 635 (B T = 1.6 T) summarized in figure 3. According to the measurements the mode locks at 3.65 s and unlocks at 3.8 s within the shaded area. There is again a good match between the calculated and the measured perturbed poloidal magnetic field amplitude depicted in figure 4(a). This time the CXRS data is available and we are able to also check the reasonably good match between the calculated mode frequency Im[(∂Ψ mn s /∂t)/Ψ mn s ] and the frequency provided by the JET mode analysis data in figure 4(b) that is also depicted in the spectrogram from figure 3(d) (the red plot).  The locking effect is illustrated by the behaviour of the perturbed poloidal magnetic field whose abrupt increase during locking seems to be accompanied by a a similar behavior of the corresponding derived error field ( figure 4(a)). Although no pulse comments have been advanced related to the locking cause at the time, according to our calculations a clear disruptive mode resonance seems to occur indicating the error field as the cause for the mode locking effect despite the existent Error Field Correction Coil (EFCC) system on, aiming to compensate the intrinsic error fields. However after 3.8 s the perturbation seems to not completely unlock from the error field, the latter keeping a lower but comparable amplitude that follows the shape of the perturbation. Finally two examples of reported error field modes lock followed by an immediate disruption is shown for the case of two JET discharges. First the JET shot no. 84 308 (B T = 2.5 T, I p =2 A, P ICRH =3.1 MW), depicted in figure 5 is studied. No CXRS data is available therefore the modes spectrogram data from figure 5(b) is used. No particular mode number has been specified for the reported locked mode. The n = 2 mode evolution is checked (blue lines in figure 5(b). Due to a high rate of bad diagnostics data for the electrons density and temperature distribution in the vicinity of the n = 1 magnetic surface, no testing is possible for this mode of perturbation. The disrupting effect of the mode locking is clearly visible from the evolution of the poloidal component of the perturbation that matches the data analysis result (figure 5(a)) being accompanied by the corresponding resonant error field prior to the pulse disruption.
The locking area is clearly linked to the drop of the mode frequency within the 13-14 s interval (blue frequency). As it has been experimentally reported the locking effect followed by disruption is due to the error field. A similar n = 2 reported locked mode case leading to disruption is presented in figure 6 for the JET pulse no. 97 753 (B T = 2.8 T, I p =2.5 A, P NBI /P RF =15/36 MW). The good calculated vs. measured n = 2 mode amplitude match from figure 6(a) indulge us again to confidently derive the corresponding error field amplitude. The final disruptive stage is clearly felt by the mode amplitude linked to a sharp drop of the mode frequency. Contrary to the previous case, the error field is very low and has no influence on the pulse disruption. The JET data analysis software also reported no neighboring n = 1 locked mode therefore it seems that neither the resonant nor the nonresonant error field spectrum is to be blamed for the n = 2 mode locking effect. Indeed, according to the JET postpulse comments, a large sawtooth perturbation crash leads to a massive n = 2 locked mode finally contributing to the plasma collapse in this case.

Conclusions
The main goal of this paper was to provide a clear and explicit time dependent formula for the intrinsic resonant error field amplitude. A model has been proposed that not only provides the dynamic evolution of the intrinsic magnetic error field but also identifies whether the resonant error field are causing or not the mode locking effect eventually leading to disruption, for pulses whose locking mechanism cause may be difficult to be concluded. Beyond the provision of a straight calculus of the error field dynamic amplitude several JET discharges have been chosen in order to test the peculiar behavior of the intrinsic magnetic error fields during error field mode locking to check its resonant evolution. Due to the usual EFCC system on not many shots showing an significant error field effect are available to chose from. On the other hand, no dynamic measured error fields amplitudes are usually provided in order to be compared against our calculated results. This is the reason we are primarily checking the associated perturbations amplitude and frequency against the experimental results provided by data analysis. Our error field model is basically the same validated, inversely approached perturbations model we have started with plus the introduced error fields. Therefore we need to rely on the precision and correctness of the used mathematical calculus, on the structural consistency of our initial proposed model validated by the good match between the calculated and the experimental perturbations amplitude and frequency and finally on the specific and expected behavior of the error fields to be checked in various situations. To summarize, we believe that this dynamic model could be a handy tool to rapidly check whether the magnetic error field is the cause for disruptions during mode locking situations or not at least within the resonant approach. However this is a multimode model therefore the nonresonant approach is also possible to be tested. But this approach is the subject of a future work.