ceq Module

Reads in the geometry using a CHEASE output file. The CHEASE output file is read using the helper module read_chease.


Contents


Variables

TypeVisibility AttributesNameInitial
integer, private :: nr
integer, private :: nt
real, private, allocatable, dimension (:):: rho_d
real, private, allocatable, dimension (:):: eqpsi
real, private, allocatable, dimension (:):: psi_bar
real, private, allocatable, dimension (:):: fp
real, private, allocatable, dimension (:):: beta
real, private, allocatable, dimension (:):: pressure
real, private, allocatable, dimension (:):: diam
real, private, allocatable, dimension (:):: rc
real, private, allocatable, dimension (:):: qsf
real, private, allocatable, dimension (:):: rho_b
real, private, allocatable, dimension (:,:):: R_psi
real, private, allocatable, dimension (:,:):: Z_psi
real, private, allocatable, dimension (:,:,:):: drm
real, private, allocatable, dimension (:,:,:):: dzm
real, private, allocatable, dimension (:,:,:):: dbtm
real, private, allocatable, dimension (:,:,:):: dpm
real, private, allocatable, dimension (:,:,:):: dtm
real, private, allocatable, dimension (:,:,:):: dpcart
real, private, allocatable, dimension (:,:,:):: dbcart
real, private, allocatable, dimension (:,:,:):: dtcart
real, private, allocatable, dimension (:,:,:):: dbtcart
real, private, allocatable, dimension (:,:,:):: dpbish
real, private, allocatable, dimension (:,:,:):: dtbish
real, private, allocatable, dimension (:,:,:):: dbtbish
real, private :: psi_0
real, private :: psi_a
real, private :: B_T
real, private :: beta_0
real, private :: R_mag
real, private :: Z_mag
real, private :: aminor
logical, private :: init_rcenter =.true.
logical, private :: init_diameter =.true.
logical, private :: init_btori =.true.
logical, private :: init_dbtori =.true.
logical, private :: init_q =.true.
logical, private :: init_pressure =.true.
logical, private :: init_dpressure =.true.
logical, private :: init_beta =.true.
logical, private :: init_psi =.true.
logical, private :: init_invR =.true.
logical, private, parameter:: debug =.true.

Derived Types

type, public :: ceq_testing_type

Expose private routines for testing only!

Type-Bound Procedures

procedure, public, nopass :: alloc_arrays
procedure, public, nopass :: chease_chi_index
procedure, public, nopass :: dealloc_arrays
procedure, public, nopass :: derm
procedure, public, nopass :: eqdbish
procedure, public, nopass :: eqdcart
procedure, public, nopass :: mod2pi
procedure, public, nopass :: set_ceq_from_chease
procedure, public, nopass :: set_dpcart
procedure, public, nopass :: set_drm
procedure, public, nopass :: set_dzm
procedure, public, nopass :: set_eqpsi
procedure, public, nopass :: set_nrnt
procedure, public, nopass :: set_R_psi
procedure, public, nopass :: set_Z_psi
procedure, public, nopass :: set_cart_arrays
procedure, public, nopass :: set_bish_arrays
procedure, public, nopass :: get_cart_arrays
procedure, public, nopass :: get_bish_arrays

Functions

private function chease_chi_index(nchi, itheta)

Convert from theta-index which is -pi to pi, to chi-index, which is 0 to 2*pi. Assumes nchi is even?

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: nchi
integer, intent(in) :: itheta

Return Value integer

public function initialize_invR(init)

If init == 1, set flag to initialize invR spline when invR is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer :: init

Return Value integer

public function invR(r, theta)

Return the inverse of the major radius at

Arguments

Type IntentOptional AttributesName
real, intent(in) :: r
real, intent(in) :: theta

Return Value real

public function Rpos(r, theta)

Return the major radius at

Arguments

Type IntentOptional AttributesName
real, intent(in) :: r
real, intent(in) :: theta

Return Value real

public function Zpos(r, theta)

Return the vertical coordinate at

Arguments

Type IntentOptional AttributesName
real, intent(in) :: r
real, intent(in) :: theta

Return Value real

public function initialize_psi(init)

If init == 1, set flag to initialize psi spline when psi is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function psi(r)

Return the flux surface label psi

Arguments

Type IntentOptional AttributesName
real, intent(in) :: r

Return Value real

private function mod2pi(theta)

Put theta into the range

Arguments

Type IntentOptional AttributesName
real, intent(in) :: theta

Return Value real

public function initialize_diameter(init)

If init == 1, set flag to initialize diameter spline when diameter is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function diameter(rp)

Return the diameter of the flux surface at a given major radius

Arguments

Type IntentOptional AttributesName
real, intent(in) :: rp

Return Value real

public function initialize_rcenter(init)

If init == 1, set flag to initialize rcenter spline when rcenter is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function rcenter(rp)

Return the major radius of the centre of a given flux surface

Arguments

Type IntentOptional AttributesName
real, intent(in) :: rp

Return Value real

public function initialize_dbtori(init)

If init == 1, set flag to initialize dbtori spline when dbtori is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function dbtori(pbar)

at the given normalised

Arguments

Type IntentOptional AttributesName
real, intent(in) :: pbar

Return Value real

public function initialize_btori(init)

If init == 1, set flag to initialize btori spline when btori is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function btori(pbar)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
real, intent(in) :: pbar

Return Value real

public function initialize_q(init)

If init == 1, set flag to initialize q spline when qfun is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function qfun(pbar)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
real, intent(in) :: pbar

Return Value real

public function initialize_pressure(init)

If init == 1, set flag to initialize pressure spline when pfun is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function pfun(pbar)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
real, intent(in) :: pbar

Return Value real

public function initialize_dpressure(init)

If init == 1, set flag to initialize dpressure spline when dpfun is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function dpfun(pbar)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
real, intent(in) :: pbar

Return Value real

public function initialize_beta(init)

If init == 1, set flag to initialize beta spline when betafun is called FIXME: this always returns 1, change to subroutine

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: init

Return Value integer

public function betafun(pbar)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
real, intent(in) :: pbar

Return Value real


Subroutines

private subroutine set_nrnt(nradius, ntheta)

Set module-level nt, nr PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: nradius
integer, intent(in) :: ntheta

private subroutine set_drm(drm_in)

Set module-level drm PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr,nt,2):: drm_in

private subroutine set_dzm(dzm_in)

Set module-level dzm PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr,nt,2):: dzm_in

private subroutine set_dpcart(dpcart_in)

Set module-level dpcart PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr,nt,2):: dpcart_in

private subroutine set_eqpsi(eqpsi_in)

Set module-level eqpsi PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr):: eqpsi_in

private subroutine set_R_psi(R_psi_in)

Set module-level R_psi PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr, nt):: R_psi_in

private subroutine set_Z_psi(Z_psi_in)

Set module-level Z_psi PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr, nt):: Z_psi_in

private subroutine set_cart_arrays(dpcart_in, dbtcart_in, dtcart_in)

Set module-level dpcart, dbtcart and dtcart PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr, nt, 2):: dpcart_in
real, intent(in), dimension(nr, nt, 2):: dbtcart_in
real, intent(in), dimension(nr, nt, 2):: dtcart_in

private subroutine set_bish_arrays(dpbish_in, dbtbish_in, dtbish_in)

Set module-level dpbish, dbtbish and dtbish PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(nr, nt, 2):: dpbish_in
real, intent(in), dimension(nr, nt, 2):: dbtbish_in
real, intent(in), dimension(nr, nt, 2):: dtbish_in

private subroutine get_cart_arrays(dpcart_out, dbtcart_out, dtcart_out)

Get module-level dpcart, dbtcart and dtcart PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(out), dimension(nr, nt, 2):: dpcart_out
real, intent(out), dimension(nr, nt, 2):: dbtcart_out
real, intent(out), dimension(nr, nt, 2):: dtcart_out

private subroutine get_bish_arrays(dpbish_out, dbtbish_out, dtbish_out)

Get module-level dpbish, dbtbish and dtbish PURELY FOR TESTING FIXME: Remove when we have more control over module-level variables for testing

Arguments

Type IntentOptional AttributesName
real, intent(out), dimension(nr, nt, 2):: dpbish_out
real, intent(out), dimension(nr, nt, 2):: dbtbish_out
real, intent(out), dimension(nr, nt, 2):: dtbish_out

public subroutine ceqin(eqfile, psi_0_out, psi_a_out, rmaj, B_T0, avgrmid, initeq, in_nt, nthg)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
character(len=*), intent(in) :: eqfile
real, intent(out) :: psi_0_out
real, intent(out) :: psi_a_out
real, intent(out) :: rmaj
real, intent(out) :: B_T0
real, intent(out) :: avgrmid
integer, intent(in) :: initeq
logical, intent(in) :: in_nt
integer, intent(out) :: nthg

private subroutine set_ceq_from_chease()

Set the ceq module-level variables from the read_chease module variables.

Read more…

Arguments

None

private subroutine alloc_arrays(nr, nt)

Allocate module-level arrays

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: nr

Number of points in radial direction

integer, intent(in) :: nt

Number of points in theta direction

private subroutine dealloc_arrays()

Deallocate module-level arrays

Arguments

None

public subroutine ceq_finish()

Finish this module; deallocate arrays

Arguments

None

public subroutine ceq_init()

Initialise CHEASE equilibrium module

Read more…

Arguments

None

private subroutine derm(f, dfm, char)

Calculate the derivative of f w.r.t to the radial and poloidal indexes (i.e. calculate the finite differences without dividing by delta psi and delta theta). - dfm(:,:,1) is the psi derivative - dfm(:,:,2) is the theta derivative - char specifies the periodicity in theta - 'E', 'O' mean continuous at theta = +/- pi - 'T' means a jump of 2 pi at theta = +/- pi

Arguments

Type IntentOptional AttributesName
real, intent(in) :: f(:,:)
real, intent(out) :: dfm(:,:,:)
character(len=1), intent(in) :: char

public subroutine gradient(rgrid, theta, grad, char, rp, nth_used, ntm)

Calculate the derivative of a selected field w.r.t. R and Z

Arguments

Type IntentOptional AttributesName
real, intent(in) :: rgrid(-ntm:)

Minor radius grid of flux surface (?)

real, intent(in) :: theta(-ntm:)

Theta grid

real, intent(out) :: grad(-ntm:,:)

grad(:,1): derivative w.r.t. R as function of theta grad(:,2): derivative w.r.t Z as function of theta

character(len=1), intent(in) :: char

Select field to take gradient of: - 'D': Bt - 'P': psi - 'R': pressure - 'T': theta - 'B': for debugging with bishop==0 only

real, intent(in) :: rp

Radius for taking pressure derivative (FIXME: double-check)

integer, intent(in) :: nth_used

Number of theta points

integer, intent(in) :: ntm

Lower index of theta array

public subroutine bgradient(rgrid, theta, grad, char, rp, nth_used, ntm)

Calculate the derivative of a selected field w.r.t. the Bishop coordinate system

Arguments

Type IntentOptional AttributesName
real, intent(in) :: rgrid(-ntm:)

Minor radius grid of flux surface (?)

real, intent(in) :: theta(-ntm:)

Theta grid

real, intent(out) :: grad(-ntm:,:)

Computed derivative, second dimension is the Bishop coordinate

character(len=1), intent(in) :: char

Select field to take gradient of: - 'D': Bt - 'P': psi - 'R': pressure - 'T': theta - 'B': for debugging with bishop==0 only

real, intent(in) :: rp

Radius for taking pressure derivative (?)

integer, intent(in) :: nth_used

Number of theta points

integer, intent(in) :: ntm

Lower index of theta array

public subroutine eqitem(r, theta_in, f, fstar)

Calculates fstar which is f interpolated at the location (r,theta). Here r is the normalised poloidal flux coordinate rp (= psi_pN + psi_0N) and theta_in is theta. f is a grid of values of f as a function of psi_p,theta

Arguments

Type IntentOptional AttributesName
real, intent(in) :: r
real, intent(in) :: theta_in
real, intent(inout), dimension(:,:):: f
real, intent(out) :: fstar

private subroutine eqdcart(dfm, dfcart)

Converts derivatives w.r.t. (psi_index,theta_index) to derivatives w.r.t. (R,Z). Note that dfcart(1, :, :) is not valid

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension (:,:,:):: dfm

dfm(:,:,1) is deriv w.r.t. psi_index (i.e. (df/dpsi)_theta * delta psi); dfm(:,:,2) is deriv w.r.t. theta_index

real, intent(out), dimension (:,:,:):: dfcart

dfcart(:,:,1) is deriv w.r.t. R; dfcart(:,:,2) is deriv w.r.t. Z

public subroutine eqdbish(dcart, dbish)

Convert gradients of a function f w.r.t. R,Z into bishop form Note that dbish(1, :, :) is not valid

Arguments

Type IntentOptional AttributesName
real, intent(in), dimension(:, :, :):: dcart

dcart(:,:,1) is gradient of f w.r.t. R; dcart(:,:,2) is gradient of f w.r.t. Z

real, intent(out), dimension(:, :, :):: dbish

dbish(:,:,1) is set to (df/dR dpsi/dR + df/dZ dpsi/dZ)/|grad psi|; dbish(:,:,2) is set to (-df/dR dpsi/dZ + df/dZ dpsi/dR)/|grad psi|. Note that in the special case where f=psi: dbish(:,:,1) is |grad psi|; dbish(:,:,2) is 0

public subroutine Hahm_Burrell(irho, a)

FIXME : Add documentation

Arguments

Type IntentOptional AttributesName
integer, intent(in) :: irho
real, intent(in) :: a