Analytic local equilibria via generalised forms of the Miller model
Justin Ball's MIT Masters thesis and/or "GS2 analytic geometry specification" provide a good overview of how this module works.
There are four equilibrium models in this module:
These can be chosen via the geoType
input parameter (0-4
,
respectively)
All of these models allow you to specify the amplitudes of two shaping modes, along with their phases and radial derivatives. The generalised elongation and Fourier models also allow you to specify the exact mode numbers.
Type | Visibility | Attributes | Name | Initial | |||
---|---|---|---|---|---|---|---|
character(len=:), | public, | allocatable | :: | type_name |
Name of the specific implementation |
||
logical, | public | :: | initialised | = | .false. |
Has this instance been initialised |
|
logical, | public | :: | has_full_theta_range | = | .false. |
Does this implementation have the full or half theta grid? If only half then we assume up-down symmetry and use this to set gradients etc. |
|
integer, | public | :: | nr |
Typically the number of radial and theta grid points. |
|||
integer, | public | :: | nt |
Typically the number of radial and theta grid points. |
|||
real, | public | :: | psi_0 | ||||
real, | public | :: | psi_a | ||||
real, | public | :: | B_T | ||||
real, | public | :: | beta_0 | ||||
real, | public | :: | R_mag | ||||
real, | public | :: | Z_mag | ||||
real, | public | :: | aminor | ||||
real, | public, | allocatable, dimension(:,:) | :: | R_psi |
2D map of co-ordinates to R, Z or B |
||
real, | public, | allocatable, dimension(:,:) | :: | Z_psi |
2D map of co-ordinates to R, Z or B |
||
real, | public, | allocatable, dimension(:,:) | :: | B_psi |
2D map of co-ordinates to R, Z or B |
||
real, | public, | allocatable, dimension(:,:) | :: | eqth |
2D maps of theta and eqpsi |
||
real, | public, | allocatable, dimension(:,:) | :: | eqpsi_2d |
2D maps of theta and eqpsi |
||
real, | public, | allocatable, dimension (:,:,:) | :: | dpcart |
Minor radius gradient in cylindrical coordinates |
||
real, | public, | allocatable, dimension (:,:,:) | :: | dtcart |
Theta gradient in cylindrical coordinates |
||
real, | public, | allocatable, dimension (:,:,:) | :: | dbcart |
B gradient in cylindrical coordinates |
||
real, | public, | allocatable, dimension (:,:,:) | :: | dpbish |
Minor radius gradient in Bishop coordinates |
||
real, | public, | allocatable, dimension (:,:,:) | :: | dtbish |
Theta gradient in Bishop coordinates |
||
real, | public, | allocatable, dimension (:,:,:) | :: | dbbish |
B gradient in Bishop coordinates |
||
real, | public, | allocatable, dimension(:) | :: | eqpsi |
Equilibrium psi (not normalised) and pressure grids |
||
real, | public, | allocatable, dimension(:) | :: | pressure |
Equilibrium psi (not normalised) and pressure grids |
||
real, | public, | allocatable, dimension(:) | :: | psi_bar |
Normalised psi, fp, q, beta, diameter and R_centre values on grid |
||
real, | public, | allocatable, dimension(:) | :: | fp |
Normalised psi, fp, q, beta, diameter and R_centre values on grid |
||
real, | public, | allocatable, dimension(:) | :: | qsf |
Normalised psi, fp, q, beta, diameter and R_centre values on grid |
||
real, | public, | allocatable, dimension(:) | :: | beta |
Normalised psi, fp, q, beta, diameter and R_centre values on grid |
||
real, | public, | allocatable, dimension(:) | :: | diam |
Normalised psi, fp, q, beta, diameter and R_centre values on grid |
||
real, | public, | allocatable, dimension(:) | :: | rc |
Normalised psi, fp, q, beta, diameter and R_centre values on grid |
||
real, | public | :: | thetaShift |
Used to translate definition of the poloidal angle such that the location of the maximium magnetic field is at pi (see section 3.2.3 of Ball MIT Masters thesis or section 2.4 of "GS2 analytic geometry specification") |
|||
type(flux_surface_type), | public | :: | surf |
procedure , public , :: initialise Subroutine | |
procedure , public , :: gradient Subroutine | Calculate the derivative of rp w.r.t. R and Z and return the modulus sqrt(drp/dR ^ 2 + drp/dZ^2). I.e. return |grad rp|. Parameters are: - rgrid: Submodule radial grid (r_circ for EFIT, rp otherwise) - theta: grid of theta - gradf(:,1): |grad psi| - gradf(:,2): 0.0 - char: if char = 'R', return |grad pressure| instead - char: if char = 'T', then return theta gradient in bishop form - gradf(:,1): (dtheta/dZ d rp/dR - dtheta/dR drp/dZ)/ |grad rp| - gradf(:,2): (dtheta/dR d rp/dR + dtheta/dZ drp/dZ)/ |grad rp| - rp: value of rp on the flux surface where we want the grad - nth: number of theta points - ntgrid: lower index of theta array is -ntgrid |
procedure , public , :: calculate_gradients Subroutine | |
procedure , public , :: rhofun => rhofun_null Function | |
procedure , public , :: derm Subroutine | |
procedure , public , :: invR Function | |
procedure , public , :: eqdbish Subroutine | |
procedure , public , :: eqdcart Subroutine | |
procedure , public , :: dealloc_common_arrays Subroutine | |
procedure , public , :: alloc_common_arrays Subroutine | |
procedure , public , :: hahm_burrell Subroutine | |
procedure , public , :: dealloc_arrays Subroutine | |
procedure , public , :: alloc_arrays Subroutine | |
procedure , public , :: alloc_special_arrays => alloc_special_arrays_null Subroutine | |
procedure , public , :: dealloc_special_arrays => dealloc_special_arrays_null Subroutine | |
procedure , public , :: read_and_set Subroutine | |
procedure , public , :: finalise Subroutine | |
procedure , public , :: finish_setup Subroutine | |
procedure , public , :: qfun Function | |
procedure , public , :: dbtori Function | |
procedure , public , :: btori Function | |
procedure , public , :: diameter Function | |
procedure , public , :: rcenter Function | |
procedure , public , :: zpos => Zpos Function | |
procedure , public , :: rpos => Rpos Function | |
procedure , public , :: eqitem Function | |
procedure , public , :: rfun Function | |
procedure , public , :: psi Function | |
procedure , public , :: betafun Function | |
procedure , public , :: dpfun Function | |
procedure , public , :: pfun Function | |
procedure , public , :: geo_FourierSeries Subroutine | |
procedure , public , :: geo_generalizedElongation Subroutine | |
procedure , public , :: geo_global Subroutine | |
procedure , public , :: geo_Miller Subroutine | |
procedure , public , :: geo_MillerExtendedHarmonic Subroutine |
Calculates fstar
which is f
interpolated at the location theta_in
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r_in | |||
real, | intent(in) | :: | theta_in |
Desired theta location |
||
real, | intent(in), | dimension(:, :) | :: | f |
Array to interpolate |
|
character(len=1), | intent(in) | :: | char |
FIXME: Unused |
Return the major radius at
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta |
Return the axial coordinate at
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta |
FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta |
FIXME : Add documentation | Should this be equivalent to psi as in others?
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta |
FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | rp |
Return the major radius of the centre of the flux surface
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | rp |
FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | pbar |
FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | pbar |
FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | pbar |
FIXME : Add documentation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | pbar |
Set module-level parameters and initialise the leq module
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(inout) | :: | self | |||
type(geo_input_type), | intent(in) | :: | inputs |
Finalise leq module
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(inout) | :: | self |
Finish the leq setup
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(inout) | :: | self |
Calculate the cylindrical (r, \theta) for the generalised Miller model
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta | |||
real, | intent(out) | :: | Rpos | |||
real, | intent(out) | :: | Zpos |
Calculate the cylindrical at for the global MHD model
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta | |||
real, | intent(out) | :: | Rpos | |||
real, | intent(out) | :: | Zpos |
Calculate the cylindrical at for the generalised elongation model
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta | |||
real, | intent(out) | :: | Rpos | |||
real, | intent(out) | :: | Zpos |
Calculate the cylindrical at for the Fourier model
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta | |||
real, | intent(out) | :: | Rpos | |||
real, | intent(out) | :: | Zpos |
Calculate the cylindrical at for the Miller Extended Harmonic model
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(leq_type), | intent(in) | :: | self | |||
real, | intent(in) | :: | r | |||
real, | intent(in) | :: | theta | |||
real, | intent(out) | :: | Rpos | |||
real, | intent(out) | :: | Zpos |