SPEC 3.20
Stepped Pressure Equilibrium Code

The namelist numericlist controls internal resolution parameters that the user rarely needs to consider. More...

Collaboration diagram for numericlist:

Variables

integer inputlist::linitialize = 0
 Used to initialize geometry using a regularization / extrapolation method. More...
 
integer inputlist::lautoinitbn = 1
 Used to initialize \(B_{ns}\) using an initial fixed-boundary calculation. More...
 
integer inputlist::lzerovac = 0
 Used to adjust vacuum field to cancel plasma field on computational boundary. More...
 
integer inputlist::ndiscrete = 2
 resolution of the real space grid on which fast Fourier transforms are performed is given by Ndiscrete*Mpol*4 More...
 
integer inputlist::nquad = -1
 Resolution of the Gaussian quadrature. More...
 
integer inputlist::impol = -4
 Fourier resolution of straight-fieldline angle on interfaces. More...
 
integer inputlist::intor = -4
 Fourier resolution of straight-fieldline angle on interfaces;. More...
 
integer inputlist::lsparse = 0
 controls method used to solve for rotational-transform on interfaces More...
 
integer inputlist::lsvdiota = 0
 controls method used to solve for rotational-transform on interfaces; only relevant if Lsparse = 0 More...
 
integer inputlist::imethod = 3
 controls iterative solution to sparse matrix arising in real-space transformation to the straight-fieldline angle; only relevant if Lsparse.eq.2; More...
 
integer inputlist::iorder = 2
 controls real-space grid resolution for constructing the straight-fieldline angle; only relevant if Lsparse>0 More...
 
integer inputlist::iprecon = 0
 controls iterative solution to sparse matrix arising in real-space transformation to the straight-fieldline angle; only relevant if Lsparse.eq.2; More...
 
real inputlist::iotatol = -1.0
 tolerance required for iterative construction of straight-fieldline angle; only relevant if Lsparse.ge.2
 
integer inputlist::lextrap = 0
 geometry of innermost interface is defined by extrapolation
 
integer inputlist::mregular = -1
 maximum regularization factor More...
 
integer inputlist::lrzaxis = 1
 controls the guess of geometry axis in the innermost volume or initialization of interfaces More...
 
integer inputlist::ntoraxis = 3
 the number of \(n\) harmonics used in the Jacobian \(m=1\) harmonic elimination method; only relevant if Lrzaxis.ge.1 .
 

Detailed Description

The namelist numericlist controls internal resolution parameters that the user rarely needs to consider.

Variable Documentation

◆ linitialize

integer inputlist::linitialize = 0

Used to initialize geometry using a regularization / extrapolation method.

  • if Linitialize = \(-I\) , where \(I\) is a positive integer, the geometry of the \(i=1,N_V-I\) surfaces constructed by an extrapolation
  • if Linitialize = 0, the geometry of the interior surfaces is provided after the namelists in the input file
  • if Linitialize = 1, the interior surfaces will be intialized as \(R_{l,m,n} = R_{N,m,n} \psi_{t,l}^{m/2}\), where \(R_{N,m,n}\) is the plasma boundary and \(\psi_{t,l}\) is the given toroidal flux enclosed by the \(l\)-th interface, normalized to the total enclosed toroidal flux; a similar extrapolation is used for \(Z_{l,m,n}\)
  • Note that the Fourier harmonics of the boundary is always given by the Rbc and Zbs given in physicslist.
  • if Linitialize = 2, the interior surfaces and the plasma boundary will be intialized as \(R_{l,m,n} = R_{W,m,n} \psi_{t,l}^{m/2}\), where \(R_{W,m,n}\) is the computational boundary and \(\psi_{t,l}\) is the given toroidal flux enclosed by the \(l\)-th interface, normalized to the total enclosed toroidal flux; a similar extrapolation is used for \(Z_{l,m,n}\)
  • Note that, for free-boundary calculations, the Fourier harmonics of the computational boundary are always given by the Rwc and Zws given in physicslist.
  • if Linitialize = 1, 2, it is not required to provide the geometry of the interfaces after the namelists

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), preset(), rzaxis(), and allglobal::wrtend().

◆ lautoinitbn

integer inputlist::lautoinitbn = 1

Used to initialize \(B_{ns}\) using an initial fixed-boundary calculation.

  • only relevant if Lfreebound = 1
  • user-supplied Bns will only be considered if LautoinitBn = 0

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), spec(), and allglobal::wrtend().

◆ lzerovac

integer inputlist::lzerovac = 0

Used to adjust vacuum field to cancel plasma field on computational boundary.

  • only relevant if Lfreebound = 1

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), spec(), and allglobal::wrtend().

◆ ndiscrete

integer inputlist::ndiscrete = 2

resolution of the real space grid on which fast Fourier transforms are performed is given by Ndiscrete*Mpol*4

  • constraint Ndiscrete>0

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), preset(), and allglobal::wrtend().

◆ nquad

integer inputlist::nquad = -1

Resolution of the Gaussian quadrature.

  • The resolution of the Gaussian quadrature, \(\displaystyle \int \!\! f(s) ds = \sum_k \omega_k f(s_k)\), in each volume is given by Iquad \(_v\),
  • Iquad \(_v\) is set in preset()

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), preset(), and allglobal::wrtend().

◆ impol

integer inputlist::impol = -4

Fourier resolution of straight-fieldline angle on interfaces.

  • the rotational-transform on the interfaces is determined by a transformation to the straight-fieldline angle, with poloidal resolution given by iMpol
  • if iMpol<=0, then iMpol = Mpol - iMpol

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), preset(), and allglobal::wrtend().

◆ intor

integer inputlist::intor = -4

Fourier resolution of straight-fieldline angle on interfaces;.

  • the rotational-transform on the interfaces is determined by a transformation to the straight-fieldline angle, with toroidal resolution given by iNtor
  • if iNtor<=0 then iNtor = Ntor - iNtor
  • if Ntor==0, then the toroidal resolution of the angle transformation is set lNtor = 0

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), preset(), and allglobal::wrtend().

◆ lsparse

integer inputlist::lsparse = 0

controls method used to solve for rotational-transform on interfaces

  • if Lsparse = 0, the transformation to the straight-fieldline angle is computed in Fourier space using a dense matrix solver, F04AAF
  • if Lsparse = 1, the transformation to the straight-fieldline angle is computed in real space using a dense matrix solver, F04ATF
  • if Lsparse = 2, the transformation to the straight-fieldline angle is computed in real space using a sparse matrix solver, F11DEF
  • if Lsparse = 3, the different methods for constructing the straight-fieldline angle are compared

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), tr00ab(), and allglobal::wrtend().

◆ lsvdiota

integer inputlist::lsvdiota = 0

controls method used to solve for rotational-transform on interfaces; only relevant if Lsparse = 0

  • if Lsvdiota = 0, use standard linear solver to construct straight fieldline angle transformation
  • if Lsvdiota = 1, use SVD method to compute rotational-transform

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), final_diagnostics(), sphdf5::mirror_input_to_outfile(), tr00ab(), and allglobal::wrtend().

◆ imethod

integer inputlist::imethod = 3

controls iterative solution to sparse matrix arising in real-space transformation to the straight-fieldline angle; only relevant if Lsparse.eq.2;

See also
tr00ab() for details
  • if imethod = 1, the method is RGMRES
  • if imethod = 2, the method is CGS
  • if imethod = 3, the method is BICGSTAB

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), tr00ab(), and allglobal::wrtend().

◆ iorder

integer inputlist::iorder = 2

controls real-space grid resolution for constructing the straight-fieldline angle; only relevant if Lsparse>0

determines order of finite-difference approximation to the derivatives

  • if iorder = 2,
  • if iorder = 4,
  • if iorder = 6,

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), tr00ab(), and allglobal::wrtend().

◆ iprecon

integer inputlist::iprecon = 0

controls iterative solution to sparse matrix arising in real-space transformation to the straight-fieldline angle; only relevant if Lsparse.eq.2;

See also
tr00ab() for details
  • if iprecon = 0, the preconditioner is ‘N’
  • if iprecon = 1, the preconditioner is ‘J’
  • if iprecon = 2, the preconditioner is ‘S’

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), tr00ab(), and allglobal::wrtend().

◆ mregular

integer inputlist::mregular = -1

maximum regularization factor

  • if Mregular.ge.2, then regumm \(_i\) = Mregular \(/ 2 \) where m \(_i > \) Mregular

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), preset(), and allglobal::wrtend().

◆ lrzaxis

integer inputlist::lrzaxis = 1

controls the guess of geometry axis in the innermost volume or initialization of interfaces

  • if iprecon = 1, the centroid is used
  • if iprecon = 2, the Jacobian \(m=1\) harmonic elimination method is used

Referenced by allglobal::broadcast_inputs(), allglobal::check_inputs(), sphdf5::mirror_input_to_outfile(), rzaxis(), and allglobal::wrtend().