Contacts: Guoqiang Li, Orso Meneghini

Short Description

PEST3 resistive MHD stability code


MHD, resistive, stability

Long Description

PEST3 is the nonideal (Delta’) version of the Princeton Equilibrium and STability (PEST) suite of codes, which has been developed since the late 1970s to determine the linear stability properties of toroidal, magnetically confined plasma equilibria.

Here is a short overview of the currently existing PEST versions:

  • PEST1: Determines the ideal MHD growth rate of a displacement perturbation, similarly to ERATO/GATO and other toroidal stability codes [1].

  • PEST2: Only solves for the normal displacement component, the two other components having been eliminated analytically by assuming the plasma to have zero inertia for displacements lying on the flux surfaces. This reduces the size of the problem and allows accurate determination of the marginal stability point. However, the obtained growth rate is ‘fake’ (not physically relevant). PEST2 is similar to DCON [2].

  • PEST3 (this version): Assumes zero-inertia for all components of the displacement. The result is a displacement field that is strongly singular (i.e. non-square integrable) about rational surfaces. The leading coefficients of the dominant and recessive components of the normal displacement are extracted to form a matching (Delta’-)matrix, which can be used in a postprocessing stage to determine the nonideal stability of linear or weakly nonlinear modes. Modes that fall in this category are the nonlinear resistive and neoclassical tearing modes for instance, provided the island width is small compared to the minor radius. In addition to computing the matching data, PEST3 also returns a yes/no answer to the ideal stability question.

PEST3 uses a variational scheme to calculate the Delta’-matrix described in PEST3 the reference. In the past few years, the code has changed significantly, although only cosmetically. PEST3 is has been rewritten in Fortran 90 to ensure full portability across UNIX workstations. Arrays are allocated dynamically, and more significantly perhaps, is the ability to embed or extend PEST3 by hooking to its API calls. Thus, it has now become effortless to integrate PEST3 into a foreign environment (e.g. Python, Matlab,…). A C++ driver is provided to illustrate the use of the API calls. Finally, PEST3 by interfacing to I2MEX, has overcome its long standing dependence on the JSOLVER equilibrium format; at the time this README is being written, there are 7 supported formats including direct access to TRANSP UFILES, G-EQDSK files stored locally or in a MDSPlus database, with the list likely to grow.

Typical workflows

This module is used to:

  • Loading of existing equilibrium file

  • Setup of toroidal/poloidal mode number, wall distance and polodial modes span

  • Plotting of 2D mode structures

Supported devices

  • Device independent

Relevant publications

  • PEST1: Grimm, R C, Green J M, and Johnson J L, Meth. Comput. Phys. 9, 253 (1976)

  • PEST2: Grimm, R C, Dewar R L, and Manickam J, J. Comput. Phys. 49, 94 (1983)

  • PEST3: Pletzer, A, Bondeson A and Dewar R L, J. Comput. Phys. 115, 530 (1994)

Technical info

Usage: pest3 [Options]
-i<if> Format (1). Possible values are:
        -1 read binary equilibrium input file (INP1) from CHEASE
         0 input equilibrium from TRANSP Ufiles (local or MDS+)
         1 read netCDF equilibrium input file (inp1.cdf) from CHEASE
         2 read netCDF equilibrium input file (eqdsk.cdf) from JSOLVER
         3 read G-EQDSK equilibrium input file from EFIT. This can also
           be an MDSPlus path, for instance:
         4 read G-EQDSK equilibrium input file from EFIT and rerun thru
           ESC-Q. Can also be an MDSPlus path (see above).
         5 read Menard netCDF input file in psipqgRZ MKS format
         6 read Belova netCDF input file format (freeqbe)
         8 read L. Don Pearlstein's inverse equilibrium file
-n<ns> Set the toroidal mode number (2)
-m<ms> Set the resonant poloidal mode (3). Alternatively,
         the user can also give a string of the form "..x.x" to
         specify the rational surfaces where the Delta' are to be
         computed: an 'x' for a selected surface and a '.' for surfaces
         that should be skipped. Each symbol ('.' or 'x') corresponds
         to a surface where ns*q is an integer, counting from the axis.
-k<Ns> Set number of radial finite elements (100)
         To perform a convergence study type -k"Ns1 Ns2 ..." where
         Ns1, Ns2, .. are a range of radial number of finite elements.
         The matching data so obtained are the extrapolation to inifinte
         resolution (Ns->oo) and an estimate of the finite element error
         is returned assuming a convergence in 1/Ns^2. Note that that the
         Frobebius coefficients of the large solution are only computed
         once, for Ns1. A typical choice is -k"100 70 140 200".
-l<Lmax1> Poloidal Fourier modes span -|Lmax1|,...+|Lmax1| (10)
-b<Bw> Relative wall distance from plasma surface normalized to
         minor radius (0). For instance -b0.2 for wall position at 1.2
         times minor radius. Bw=0 is wall on plasma surface whereas Bw>=10 is
         wall at infinity.
-S     Symmetric large solution support. Default is asymmetric.
         This option is typically used for higher m modes, which tend to
         be localized about their resonant surface
-d<rw> Relative extension of large solution to neighboring rational
         surfaces (0.9). Must have 0.05 < rw < 0.99. If -S is in addition
         selected, then rw is the smallest of the distance between the left
         and right rational surfaces.
-f<File> Read equilibrium data from <File>. <File> can be a local
         netCDF file (e.g. inp1.cdf) or, if MDSPlus access has been enabled
         a path such as: MDS+:TRANSPGRID.PPPL.GOV:TRANSP(TFTR.88,37065Y82)
-T<time> Set MDSPlus equilibrium time to <time> in [sec] (0).
         Will correct to the min/max value if selected <time> is outside
-q<newQ> Rescale equilibrium by adjusting the safety factor q=newQ
         at the plasma edge. If negative, then set q on axis to abs(newQ).
         Leave equilibrium unchanged if newQ=0 (default).
-s<is> Solver (0). Possible values are:
         0 use PEST3's built-in linear system solver
         1 compute lowest delta-W eigenvalue using LAPACK routine
-u     Uniform radial mesh. Default is mesh packing about rational surfaces.
         In general, automatic mesh packing is desirable. This option is useful
         at high finite element resolution by preventing the clutering of nodes
         about the resonant surface.
-a<xa> Extrapolate some profiles over distance xa to axis (0.01)
-x<xe> Extrapolate some profiles over distance xe to edge (0.00)
-P<p1> Use p1 poloidal rays for mapping (129)(direct equilibrium only).
-R<rs> Use rs radial surfaces for mapping (401)(direct equilibrium only).
-E<ed> Separatrix retreat factor (0.999) (direct equilibrium only).
-A     Akima spline Jacobian and friends (default is Jacobian derived
         from X and Z).
-M<mth> Set the packing method (0).
         0 Standard with axis, surfaces and edge packed.
         Standard zones are either left default or set with widmsh and alfmsh
         input.  See widmsh and alfmsh.
         1 Edge packing with power law dzx/dz=a*(c0-zx)^alfa. This can resolve
         the peeling mode very near the separatrix.  Power alfa is set with
-p<as> For packing at edge and boundary set alfmsh vector via
         -p"a0 a1" where a0 is the axis packing power and a1 is the
         edge packing power, with both >=1.0 (1.0).  Surfaces are packed
         with alpha(mu) at the surface if packing is turned on.
         A value of 1.0 produces no packing.  If one number is entered
         the axis and edge both get that alfmsh.
-K<pkpts> Set the packpts vector for the surfaces, pkpts
         via -K"p1 p2 ... ".  If one number is entered then all packing
         locations get that p1. The number of entries can also be the the
         number of resonant surfaces in mm (no axis and edge numbers).
         These entries should be in the same order as mm, from axis to
         edge, with the first surface packpts always coming first.
-W<ws> Set the widmsh vector for the surfaces, axis and edge with ws
         via -W"w1 w2 ... ".  If one number is entered then all packing
         locations get that widmsh. The number of entries can also be the
         number of resonant surfaces in mm plus 2 for the axis and edge.
         These entries should be in the same order as mm, from axis to
         edge, with the axis widmsh always coming first.


List of contributors sorted by number of lines authored:

290 Orso Meneghini
 22 Fusion Bot


List of usernames sorted by number of module imports: kimkyungjin, meneghini, halfmoonm, kyungjin, marinoni, MHalfmoon, chenji, knolkerm, eldond, fila, halfmoon, lig, lmorton, mcclenaghanj, nelsonand, smithsp