1 Introduction
In plasma physics the magnetic hydro dynamic (MHD) equilibrium equation
$\unicode[STIX]{x1D735}P=J\times B$
is one of the most assessed equations. The analytical solution to the toroidal axisymmetric Grad Shafranov equation (Shafranov Reference Shafranov1960; Mukhovatov & Shafranov Reference Mukhovatov and Shafranov1971) was at the basis of the breakthrough realized by the tokamak configuration along the path of fusion energy achievements. In axisymmetry the Grad Shafranov equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn1.gif?pub-status=live)
where
$J_{\unicode[STIX]{x1D711}}$
is the toroidal plasma current,
$\unicode[STIX]{x1D713}$
is the flux function and (
$R,z,\unicode[STIX]{x1D711}$
) are the standard cylindrical coordinates; in the vacuum
$(J_{\unicode[STIX]{x1D711}}=0)$
coincides with the Laplace equation for the toroidal component of the vector potential (Bateman Reference Bateman1978). The analytical solution to the related Laplace equation in several different geometries was deeply studied during the second half of the 19th century (Neumann Reference Neumann1864; Böcher Reference Böcher1894; Wangerin Reference Wangerin1909), pointing to the possibility of finding an analytical solution also for oblate and prolate toroidal geometries. A good summary of all these works can be found in Moon and Spencer’s book (Reference Moon and Spencer1971). During the first half of the 20th century the problem of finding an analytical solution to the inductance of a torus was first tackled (Fock Reference Fock1932), and a first analytical solution to the Laplace equation, in elliptical oblate coordinates, was found by Lebedev (Reference Lebedev1937). More recently Chance (Reference Chance1997) has solved the Laplace equation for a general azimuthal geometry by using a Green function technique, even in presence of internal currents. Eventually, a team of Japanese scientists found the analytical solution to the vacuum Grad Shafranov equation for an oblate toroidal system (AikawaI & Takahara Reference AikawaI and Takahara1975; Honma et al.
Reference Honma, Kito, Kaji and Seki1979). In the meanwhile, a large amount of scientific papers have been written to solve the tokamak equilibrium problem, adopting several different approaches. Analytical solutions to the nonlinear Grad Shafranov equation have been studied and found for simple cases of the source term (Soloviev Reference Soloviev1968; Greene Reference Greene1988; Zheng, Wootton & Solano Reference Zheng, Wootton and Solano1996; Cerfon & Freidberg Reference Cerfon and Freidberg2010). The employment of a system of toroidal coordinates, which has the flux surfaces as the radial coordinate, has allowed us to develop robust codes with the mapping solution expanded on the momentum (Boozer Reference Boozer1980; Lao, Hirschmann & Wieland Reference Lao, Hirschmann and Wieland1981; Brambilla Reference Brambilla2003). Predictive equilibrium codes have been developed, through the numerical resolution (Feneberg & Lackner Reference Feneberg and Lackner1973; Cenacchi, Galvao & Taroni Reference Cenacchi, Galvao and Taroni1976; Blum, Le Foll & Thooris Reference Blum, Le Foll and Thooris1981; Albanese, Blum & De Barbieri Reference Albanese, Blum and De Barbieri1987; Albanese, Ambrosino & Mattei Reference Albanese, Ambrosino and Mattei2015) and/or the adoption of a semi-analytical approach to deal with (Alladio & Crisanti Reference Alladio and Crisanti1986; Alladio et al.
Reference Alladio, Crisanti, Marinucci, Micozzi and Tanga1991) the Grad Shafranov equation, and have been then employed to design the planned plasma configuration of practically all ongoing tokamak experiments. The linearization of the problem (Albanese & Villone Reference Albanese and Villone1998) has allowed us to easily optimize a given equilibrium configuration in terms of different aspects (currents in the poloidal coils, forces on the coils, interaction with the first wall
$\ldots$
). Moreover, the use of the linearized solver has allowed us to optimize the real time control of ongoing experiments (Sartori et al.
Reference Sartori, Ambrosino, Ariola, Cenedese, Crisanti and De Tommasi2005, Reference Sartori, Crisanti, Albanese, Ambrosino, Toigo and Hay2008; Albanese et al.
Reference Albanese, Ambrosino, Ariola, Artaserse, Bellizio and Coccorese2011) and to design an optimal control system for future devices (Albanese et al.
Reference Albanese, Ambrosino, Ariola, Cavinato, De Tommasi and Liu2009; Bachmann et al.
Reference Bachmann, Aiello, Albanese, Ambrosino, Arbeiter and Aubert2015). Although it is an ill-posed problem (different internal current distributions can produce the same magnetic signal), the magnetic measurements have been used to constrain the equilibrium codes, in order to reconstruct the experimental plasma boundary and to work out, as much as possible, the most important internal features of the confined plasma; i.e. the internal magnetic and kinetic inductances
$\unicode[STIX]{x1D6FD}_{\text{p}}$
and
$l_{\text{i}}$
(Shafranov Reference Shafranov1971), various plasma moments (Van Milligen Reference Van Milligen1990; Van Milligen & Lopez Fraguas Reference Van Milligen and Lopez Fraguas1994) (for instance the plasma column centroid and its elongation), the internal mapping of the magnetic surfaces and the identification of the plasma current density. As for the predictive equilibrium codes, both numerical and semi-analytical approaches have been used (Lao et al.
Reference Lao, John, Stambaugh, Kellman and Pfeiffer1985; Alladio & Crisanti Reference Alladio and Crisanti1986; Alladio et al.
Reference Alladio, Crisanti, Marinucci, Micozzi and Tanga1991). Eventually, the fact that the plasma equilibrium reconstruction is a mathematically ill-posed problem and that in some used codes it depends on the knowledge of external currents, and the fact that experimental noise must be accounted for, has always led to a poor reconstruction of the internal mapping of the magnetic topology. In several codes internal independent measurements have been used to constrain the solution. In particular, within the most used code, EFIT, the polarimetry (Li et al.
Reference Li, Lotte, Zwingmann, Gil and Imbeaux2011) and/or the motional Stark effect (Lao et al.
Reference Lao, St. John, Peng, Ferron, Strait and Taylor2005) have both been used to further constrain the problem and to improve the reconstructed poloidal flux mapping. This code has also been employed for real time reconstruction, which is to be used on the plasma shape control system (Ferron et al.
Reference Ferron, Walker, Lao, St John, Humphreys and Leuer1998). This latter aspect has been strongly improved thanks to the use of powerful modern computer parallelization (Yue et al.
Reference Yue, Xiao, Luo and Guo2013). Of course, the use of all these techniques and of modern computers has always led to result improvements. However, the details of the internal parameters, like the current density profile, are still far from what would be necessary for a deeper understanding of the physics and/or for an accurate real time control of the internal profiles. This represents a serious concern when considering that in a DEMO class machine the available diagnostics will be strongly reduced, with the consequent necessity of having very robust and flexible plasma models. Reconstructive equilibrium codes (Alladio & Crisanti Reference Alladio and Crisanti1986; Alladio et al.
Reference Alladio, Crisanti, Marinucci, Micozzi and Tanga1991), based on a semi-analytical approach, are in principle more flexible and robust than the completely numerical ones, as they use the multipolar analytical solution of the electromagnetic problem on the most suitable function basis. In particular this approach gives the possibility to fully separate the external world (passive structures,
$\text{coils},\ldots$
) from the internal plasma features. The most important advantage given by this approach is that any ‘unphysical’ solution, which could be worked out by a pure numerical code in the presence of measurement errors, is not intrinsically possible. This aspect was widely discussed in relation to the FTU experiment’s circular toroidal geometry in Alladio & Crisanti (Reference Alladio and Crisanti1986). In that paper the robustness of the approach was also verified on the old JET experiment, without the presence of internal divertor coils. The used reconstructive equilibrium code was the one that showed the possibility in JET of producing X point configurations and the consequent possibility of achieving an improved energy confinement regime (Alladio et al.
Reference Alladio, Crisanti, Lazzaro and Tanga1984). Unfortunately, so far, the exact analytical solution of the Grad Shafranov equation for the vacuum region has been found only in relation to the standard circular shaped and for the elliptical oblate toroidal geometries (AikawaI & Takahara Reference AikawaI and Takahara1975; Honma et al.
Reference Honma, Kito, Kaji and Seki1979). Both these geometries are not suitable to the present tokamak experiments, all of which are based on an elliptical prolate geometry. An interesting approach to analytically solving the Grad Shafranov equation through the separable variable technique and using not constant source functions has been used in Atanasiu et al. (Reference Atanasiu, Günter, Lackner and Miron2004) and Guazzotto & Freidberg (Reference Guazzotto and Freidberg2007). In both papers the solution has been worked out by using an
$(x,y)$
coordinate system. Unluckily this ‘not natural’ system of coordinates results in a solution with a large expansion that does not always converge. Guazzotto & Freidberg partially overcome the problem by linking some geometrical property of the plasma boundary to the solution expansion. Nonetheless, the non-robustness of the solution remains. In the present paper the mathematical problem is summarized and the analytical solution of the vacuum Grad Shafranov equation is worked out for the cap-cyclide coordinate system, which is the natural one for the modern elongated and D-shaped tokamak. Moreover, the paper addresses how the solution could be used to realize a robust reconstructive equilibrium code, capable of better fitting the needs of a future reactor relevant tokamak, like DEMO.
2 Cap-cyclide geometry
In cylindrical coordinates
$(R,Z)$
the vacuum Grad Shafranov equation (Shafranov Reference Shafranov1960) is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn2.gif?pub-status=live)
This former equation is formally similar to the Laplace equation, but for a sign (minus instead of plus) in the first derivative part.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn3.gif?pub-status=live)
Consequently, it is obvious that any analytical solution found for the Laplace equation will be correlated to the analytical solution of the Grad Shafranov equation. Since the end of the 19th century (Böcher Reference Böcher1894; Wangerin Reference Wangerin1909) analytical solutions to the Laplace equation have been worked out using a large set of system of coordinates, that allows writing of the solution in a ‘separable’ way, i.e.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn4.gif?pub-status=live)
where (
$\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D6FD},\unicode[STIX]{x1D6FE}$
) are independent orthogonal variables. For a subset of the system of coordinates a so-called ‘
$R$
-separable’ solution can be worked out. In particular this is the solution for the coordinate systems with a toroidal symmetry, such as the proper ones for the tokamak experiments
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn5.gif?pub-status=live)
${\mathcal{R}}(\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D6FD},\unicode[STIX]{x1D706})$
being a function determined by the coordinate system.
For the present toroidal elongated tokamak, the most appropriate coordinate system is the cap-cyclide one (Moon & Spencer Reference Moon and Spencer1971), (figure 1) with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn6.gif?pub-status=live)
where (
$x,y,z$
) are the standard Cartesian coordinates,
$\text{d}n,cn,sn$
are the Jacobi elliptical functions,
$a_{s}$
is a dimensional scale parameter and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn7.gif?pub-status=live)
The new variables (
$\unicode[STIX]{x1D707},\unicode[STIX]{x1D708},\unicode[STIX]{x1D6F1}$
) are defined between
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn8.gif?pub-status=live)
Here
$K$
and
$iK^{\prime }$
are respectively the real and the imaginary complete elliptical integrals
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn9.gif?pub-status=live)
and
$k$
and
$k_{1}$
are respectively the parameter and the complementary parameter of the elliptical integral.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_fig1g.gif?pub-status=live)
Figure 1. Cap-cyclide coordinate system. Surfaces at constant
$\unicode[STIX]{x1D707}$
and
$v$
.
The coordinate transformation of (2.5) can be also written in the complex plane as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn10.gif?pub-status=live)
Writing from here on out
$\text{d}n(\unicode[STIX]{x1D707})$
,
$cn(\unicode[STIX]{x1D707})$
,
$sn(\unicode[STIX]{x1D707})=\text{d}n\unicode[STIX]{x1D707}$
,
$cn\unicode[STIX]{x1D707}$
,
$sn\unicode[STIX]{x1D707}$
and
$\text{d}n(\unicode[STIX]{x1D708})$
,
$cn(\unicode[STIX]{x1D708})$
,
$sn(\unicode[STIX]{x1D708})=\text{d}n\unicode[STIX]{x1D708}$
,
$cn\unicode[STIX]{x1D708}$
,
$sn\unicode[STIX]{x1D708}$
, for
$\unicode[STIX]{x1D707}=\text{const.}$
the coordinate surfaces are ring-cyclides (see figure 1) with equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn11.gif?pub-status=live)
with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn12.gif?pub-status=live)
and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn13.gif?pub-status=live)
where (note the definitions reported in Moon & Spencer (Reference Moon and Spencer1971) present a small mistake)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn14.gif?pub-status=live)
For
$\unicode[STIX]{x1D708}=\text{const.}$
the coordinate surfaces are cap-cyclides (see figure 1) with equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn15.gif?pub-status=live)
where (note the definitions reported in Moon & Spencer (Reference Moon and Spencer1971) present a small mistake)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn16.gif?pub-status=live)
By varying the parameter
$k$
, the coordinates transformation of equation (2.5) describes a large set of quite different geometries (figure 2). For
$\unicode[STIX]{x1D707}\rightarrow 0$
, independently of
$k$
, the geometry resembles the standard toroidal geometry; for larger values of
$\unicode[STIX]{x1D707}$
the shape of the constant
$\unicode[STIX]{x1D707}$
surfaces depend on the value of the
$k$
parameter. For
$k\rightarrow 0$
(figure 2
a) the surfaces tends to assume a bean shape, like the old PBX-M (England et al.
Reference England, Bell, Hirshman, Kaita, Kugel and LeBlanc1997) experiment. For intermediate values of
$k$
(figure 2
b) the surfaces can be either
$D$
or purely elliptical prolate shaped, like the largest part of the present ongoing tokamak experiments. For
$k\rightarrow 1$
(figure 2
c) all the surfaces are similar to the standard toroidal ones, independently of the
$\unicode[STIX]{x1D707}$
value.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_fig2g.gif?pub-status=live)
Figure 2. (a) For
$k\rightarrow 0$
the constant
$\unicode[STIX]{x1D707}$
surfaces assume a bean shape; (b) for intermediate
$k$
the surfaces can be either D or purely elliptical prolate; (c) for
$k\rightarrow 1$
the surfaces tend to the standard toroidal ones.
For
$\unicode[STIX]{x1D707}\rightarrow K$
the surface at constant
$\unicode[STIX]{x1D707}$
degenerates in a circumference arc of radius
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn17.gif?pub-status=live)
For
$k\rightarrow 1$
,
$R_{0}=1/2a_{s}$
, as is clear from the previous discussion, is the centre of the standard toroidal coordinates. For the constant
$\unicode[STIX]{x1D707}$
contours it is possible to define the local major radius, the local minus radius and the local aspect ratio.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn18.gif?pub-status=live)
By using the expansion of the Jacobi elliptical functions in terms of the hyperbolic coordinates, when
$k\rightarrow 1$
, in the limit of large aspect ratio (
$\unicode[STIX]{x1D707}\rightarrow 0$
) we have that
$A(\unicode[STIX]{x1D707})\sim \cosh (\unicode[STIX]{x1D707})$
as in the standard toroidal geometry. We can define the local elongation as
$b(\unicode[STIX]{x1D707})/a(\unicode[STIX]{x1D707})$
, where
$b(\unicode[STIX]{x1D707})$
is the local maximum height of any contour; after some heavy algebra we have
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn19.gif?pub-status=live)
For
$\unicode[STIX]{x1D707}\rightarrow 0$
,
$b(\unicode[STIX]{x1D707})/a(\unicode[STIX]{x1D707})\rightarrow 1$
independently of the parameter
$k$
, i.e. the contours are circular as for the standard toroidal case. For
$\unicode[STIX]{x1D707}\rightarrow K$
,
$b(\unicode[STIX]{x1D707})/a(\unicode[STIX]{x1D707})\rightarrow \infty$
independently of the parameter
$k$
(but
$k=1$
, where
$b(\unicode[STIX]{x1D707})/a(\unicode[STIX]{x1D707})\rightarrow 1$
), this is because, in the centre of the ‘geometry’, the constant surface converges to a circumference arc and not to a point.
3 Geometry metric
For the Laplace equation (2.2) to admit a quasi-separable solution of type (2.4), in a rotational geometry like the cap-cyclides, it is necessary and sufficient that the three following conditions are satisfied (Moon & Spencer Reference Moon and Spencer1971)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn20.gif?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn21.gif?pub-status=live)
Here
${\mathcal{R}}$
is the funcion defined in (2.4) and
$g_{i}$
are the coordinate transformation metric factors
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn22.gif?pub-status=live)
where
$\unicode[STIX]{x1D6EC}$
and
$\unicode[STIX]{x1D6E4}$
are given in (2.6) and
$\unicode[STIX]{x1D6FA}^{2}=(1-sn^{2}\unicode[STIX]{x1D707}\,\text{d}n^{2}\unicode[STIX]{x1D708})(\text{d}n^{2}\unicode[STIX]{x1D708}-k^{2}sn^{2}\unicode[STIX]{x1D707})$
and
$\unicode[STIX]{x1D6F7}_{\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D6FD}}$
are elements of the Stäckel matrix (Morse & Feshbach Reference Morse and Feshbach1953)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn23.gif?pub-status=live)
Consequently, from the third condition of (3.2a,b )
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn24.gif?pub-status=live)
Eventually, by solving the equation (3.1) we find
$\unicode[STIX]{x1D6FC}_{1}=1$
.
4 Laplace equation
As a consequence of what was discussed in the previous paragraph, the Laplace equation (2.2) in the cap-cyclide coordinates admits a quasi-separable solution of the type of (2.4).
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn25.gif?pub-status=live)
where the equations for the radial
$M(\unicode[STIX]{x1D707})$
and the angular part
$N(\unicode[STIX]{x1D708})$
can be written as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn26.gif?pub-status=live)
By substituting
$z_{1}=sn^{2}\unicode[STIX]{x1D707}$
and
$z_{2}=\text{d}n^{2}\unicode[STIX]{x1D708}$
the two equations for
$M(\unicode[STIX]{x1D707})$
and
$N(\unicode[STIX]{x1D708})$
can be written as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn27.gif?pub-status=live)
Here, the variable
$z$
represents both the variables
$z_{\text{1}}$
and
$z_{2}$
and, consequently, the equations for the functions
$M(\unicode[STIX]{x1D707})$
and
$N(\unicode[STIX]{x1D708})$
are formally similar. In this equation
$a_{1}=(1;1)$
,
$a_{2}=(1/k^{2};k^{2})$
,
$a_{3}=(0;0)$
,
$A_{0}=(-q^{2}/k^{2};-q^{2}k^{2})$
,
$A_{1}=(-p^{\prime 2};p^{\prime 2})$
,
$A_{2}=(1-q^{2};1-q^{2})$
, where the first term in the brackets refers to the equation for
$M(\unicode[STIX]{x1D707})$
, the second to the equation for
$N(\unicode[STIX]{x1D708})$
,
$p$
and
$q$
are integer numbers and
$p^{\prime }=p/k$
. Equation (4.3) is a particular case of the general Böcher (Reference Böcher1894), Moon & Spencer (Reference Moon and Spencer1971) equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn28.gif?pub-status=live)
In the particular case of (4.3) the general Böcher equation is named a Wangerin equation. This equation has two singularities of the first order for
$z=1,a_{2}$
and two singularities of the second order for
$z=0$
and
$z\rightarrow \infty$
. Consequently the Wangerin equation is characterized by poles
$\{1,1,2,2\}$
and the Laplace equation has the solution (Moon & Spencer Reference Moon and Spencer1971).
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn29.gif?pub-status=live)
5 Grad Shafranov equation
The Grad Shafranov operator of equation (2.1) in axisymmetry (
$\unicode[STIX]{x2202}/\unicode[STIX]{x2202}\unicode[STIX]{x1D711}=0$
) is similar to the toroidal part of the Laplace operator applied to the potential vector
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn30.gif?pub-status=live)
This expression for an orthogonal system of coordinates
$(u_{1},u_{2},u_{3})$
assumes the explicit form
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn31.gif?pub-status=live)
That in our case of the cap-cyclide coordinates becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn32.gif?pub-status=live)
where here
${\mathcal{R}}$
is the separation function defined in (3.5). By remembering that in our coordinates the Laplace equation in axisymmetry can be written as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn33.gif?pub-status=live)
Equation (5.3) is identical to equation (5.4), plus a non-differential term in
$A_{\unicode[STIX]{x1D711}}$
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn34.gif?pub-status=live)
Where
$\unicode[STIX]{x1D6F7}_{\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D6FD}}$
are elements of the Stäckel matrix. Consequently the set of separable solutions of (5.4) for the toroidal component of the potential vector is the same as the one of set (4.2), when assuming that the axisymmetry is
$(\unicode[STIX]{x2202}/\unicode[STIX]{x2202}\unicode[STIX]{x1D711}=0\rightarrow q^{2}=\unicode[STIX]{x1D6FC}_{3}=0)$
, plus the non-differential term of (5.5).
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn35.gif?pub-status=live)
But this set of equations, when assuming
$q^{2}=\unicode[STIX]{x1D6FC}_{3}=1$
, is identical to the one of (4.2); consequently the two-dimensional (2-D) solution of the toroidal component of the vectorial Laplace equation for the potential vector is the same as the 3-D scalar Laplace equation when assuming
$q^{2}=\unicode[STIX]{x1D6FC}_{3}=1$
. If now we remember that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn36.gif?pub-status=live)
We get that the solution for the Grad Shafranov equation (2.1) can be expressed in terms of the Wangerin functions with
$q^{2}=1$
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn37.gif?pub-status=live)
It is important to note that in (5.8) all terms, but
$A_{\text{p,p'}}$
and
$B_{\text{p,p'}}$
, are geometrical terms. The physics information, regarding the topological structure of the magnetic field, depends only on the unknown moments
$A_{\text{p,p'}}$
and
$B_{\text{p.p'}}$
, constant in the vacuum region. The moments
$A_{\text{p,p'}}$
describe the field produced by all the currents flowing internally to the vacuum region, whilst the moments
$B_{\text{p,p'}}$
describe the field produced by all the currents flowing externally to the vacuum region. Consequently, the contribution of all the currents flowing externally from the plasma region is intrinsically separated by the internal properties of the confined plasma, which are fully described by the moments
$A_{\text{p,p'}}$
. In a region with a current distribution, the solution of the Grad Shafranov equation remains formally the same as (2.1), but now the moments
$A_{\text{p,p'}}$
and
$B_{\text{p.p'}}$
are not constant anymore, and can be worked out by using a Green function technique, once the current distribution is known.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn38.gif?pub-status=live)
On a given tokamak experiment we have a set of magnetic probes, measuring the local value of the magnetic field
$B$
and of the flux function
$\unicode[STIX]{x1D713}$
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20190404135617933-0537:S0022377819000175:S0022377819000175_eqn39.gif?pub-status=live)
here the subscript
$l$
indicates the generic position of a magnetic probe. Clearly this set of magnetic probes allows ua to exactly determine the values of the moments
$A_{\text{p,p'}}$
and
$B_{\text{p.p' }}$
and, consequently, to have the exact analytical expression of the magnetic topology in the region between the plasma column and the closer poloidal coil. Moreover, this methodology allows us to discriminate between the contribution provided by the external equilibrium currents (given by the moments
$B_{\text{p,p'}}$
) and the one provided by the internal plasma currents distribution (given by the moments
$A_{\text{p,p'}}$
). As described in Alladio & Crisanti (Reference Alladio and Crisanti1986) this approach, alongside the Green technique of (5.9) and the Grad Shafranov equation, allows us to build up a robust semi -analytical reconstructive equilibrium code. The robustness of using this approach to deal with the mathematically ill-posed problem of reconstructing the internal plasma equilibrium only relying on external magnetic measurements has been verified both in the context of the FTU circular tokamak and in the ‘old’ JET, where the divertor coils were missing. The robustness of the solution obtained, against the magnetic measurements error, was widely verified at JET and it allowed us to discover the possibility for JET to realize an X point configuration (Alladio et al.
Reference Alladio, Crisanti, Lazzaro and Tanga1984; Alladio & Crisanti Reference Alladio and Crisanti1986). The drawback of the equilibrium code described in Alladio & Crisanti (Reference Alladio and Crisanti1986) was that it was worked out by using the pure toroidal geometry, with a circular poloidal section. When applied to the modern elongated tokamak, this reconstructive equilibrium code was not able to exactly discriminate between the contribution provided by the external coils and the one provided by the plasma. The solution written down in (5.10) will allow us to solve this problem and to build up a very robust semi-analytical equilibrium code, even capable of dealing with the next generation of tokamaks, where the measurements will be far from the plasma border and the signal produced by high internal magnetic moments will be from the very low up to becoming comparable with the intrinsic system noise.
6 Conclusions
Equilibrium reconstructive codes (Lao et al. Reference Lao, John, Stambaugh, Kellman and Pfeiffer1985; Alladio & Crisanti Reference Alladio and Crisanti1986; Alladio et al. Reference Alladio, Crisanti, Marinucci, Micozzi and Tanga1991) are probably the most important diagnostic tool in fusion-based experiments, and their importance will be even stronger in future DEMO class experiments, where the scarce availability of diagnostics will require very robust codes, based on fundamental physics principles. The equilibrium codes based on an analytical solution to the Grad Shafranov equation are very robust and, probably, the most suitable tools to accomplish the necessary quality in terms of robustness and reliability. Unfortunately, so far, an analytical solution to the Grad Shafranov equation was available only for circular and elliptical oblate shaped toroidal geometries, not adaptable to the present elongated tokamak. The analytical solution to the vacuum Grad Shafranov, in the cap-cyclide coordinates, as worked out in this paper, opens the possibility of realizing an equilibrium reconstructive code based on an elongated geometry. The Wangerin functions have been evaluated for the first time in AikawaI & Takahara (Reference AikawaI and Takahara1975). However, in order to use them for an equilibrium code, an independent evaluation and cross-checking will be necessary. Moreover, it will necessary to evaluate their derivatives, to be able to get the poloidal magnetic field and to open the possibility of finding the Green function for this geometry. Eventually, this will allow us to write a robust reconstructive equilibrium code, to be used on the present and future elongated tokamak experiments.
Acknowledgements
The author wants to thank all the ENEA and European colleagues that have helped with fruitful discussions. A particular thank you is due to Dr E. Giovannozzi, for the strong support received in using Matlab and Mathematica tools. This work has been carried out within the framework of the EUROfusion Consortium and has received funding from the Euratom research and training programme 2014-2018 under grant agreement no. 633053. The views and opinions expressed herein do not necessarily reflect those of the European Commission.