XGC1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Public Attributes | List of all members
ptb_3db_module Module Reference

The module "ptb_3db" contains the framework for handling perturbed 3D magnetic fields in XGCa and XGC1. More...

Public Member Functions

subroutine read_3db (grid, psn)
 read_3db allocates memory for the perturbed field and reads the perturbed field from a file. Supported file formats are (R,Z)-ASCII and M3D-C1. More...
 
subroutine read_3db_test (grid, psn)
 Sets up an analytical magnetic field perturbation with fixed toroidal mode number and (currently) 2 poloidal mode numbers. The number of poloidal mode numbers can be adjusted with few code changes. The most important parameters are the value of the poloidal mode numbers, the centers of the Gaussian envelopes of each mode number, the width of the Gaussian, and the amplitude. More...
 
subroutine read_3db_m3dc1_mesh (grid)
 Allocates memory for and reads perturbed field data from M3D-C1 file. For M3D-C1 input, the fusion-io library by N. Ferraro is used. The interpolation used in the fusion-io library keeps the field divergence free. More...
 
subroutine read_3db_m3dc1_rz (grid)
 Allocates memory for and reads perturbed field data from M3D-C1 file. This version of the routine reads the vector potential A on an (R,Z) grid. The perturbed field on the XGC mesh can be evaluated in a separate subroutine. Eventually, the perturbed magnetic field is to be calculated directy from A using spline interpolation. For M3D-C1 input, the fusion-io library by N. Ferraro is used. More...
 
subroutine extract_resonant_vector_potential_from_m3dc1 (grid, psn, itor)
 Extension of RZ reader for EM Use RZ reader to read the full vector potential from M3D-C1. For EM, proeceed to extract the (near) resonant components of \(A_\parallel\) (A_parallel) on the XGC mesh. Interpolate those back to the RZ grid and subtract from the raw vector potential from M3D-C1. Then re-initialize the spline interpolation for A. This enables pushing particles with the full M3D-C1 field while still calculating the plasma response to the near resonant modes. As a stop-gap measure, the such filtered vector potential is used to calculate the perturbed field on the XGC mesh because the spline interpolation of the vector potential on the RZ grid is not available in C++ yet. This currently works only for single toroidal mode number input. Extending it to multiple mode numbers seems straightforward. More...
 
subroutine read_3db_m3dc1_eem (grid, psn)
 Allocates memory for and reads perturbed field data from M3D-C1 file for use with the explicit EM method. This reader reads only the parallel component of the perturbed vacuum vector potential for a single toroidal mode number. More...
 
subroutine init_apar_3db (grid, psn)
 Sets the symplectic vector potential equal to the vacuum RMP vector potential. More...
 
subroutine psn_mem_alloc_3db (psn, n)
 Allocates memory for the perturbed current in a psn-object. More...
 
real(kind=8) function ptb_3db_screening_factor (r, z, psi)
 Evaluates the screening factor for artificial exponential damping of the vacuum perturbed field towards the magnetic axis. The screening shape function is simply exp((psi-psi0)/w). Since psi decreases towards the magnetic axis, psi-psi0 becomes negative at some point. More...
 
subroutine ptb_3db_calc_db (grid, dA_phi, dB_re, dB_im)
 Evaluates the perturbed field from the perturbed vector potential on the XGC grid. More...
 
subroutine update_3db ()
 Updates the ramp-up factor of the vacuum RMP field. The external (or vacuum) RMP field can be ramped up linearly instead of a hard switch from off to on. The ramp-up factor ptb_3db_rampup_fac determines the amplitude of the external RMP field throughout XGC. This routine is used only in electrostatic simulations. In electromagnetic simulations, push_As updates the ramp-up factor. More...
 
subroutine ptb_3db_damping_factor_em (grid, Apar, dApar)
 Calculates an adaptive damping factor to limit the time rate of change of the perturbed vector potential The routine calculates the poloidal Fourier spectrum of the current increment of the vector potential (dApar –> effectively the Hamiltonian vector potential Ah) and the vector potential itself (Apar). Then it compares the amplitude of the increment and the current value and computes a damping factor for each poloidal mode number such that the increment is limited to be a fraction of the current value. Using this damping factor, the increment is transformed back to real space and returned to the calling routine. More...
 

Public Attributes

logical ptb_3db_on =.false.
 Switch 3D perturbed field on (.true.) and off (.false.) More...
 
integer ptb_3db_mode =0
 Mode of operation: (0) static RMP field (ES), or fully self-consistent RMP field (EM version) (2) self-consistent RMP field with time-averaged field-equations (EM version) More...
 
character(len=250) ptb_3db_filename ='C1.h5'
 File containing the perturbed field data. More...
 
integer ptb_3db_file_format =1
 Format of input file with perturbed field data: (0) toroidal Fourier coefficients in R,Z, (1) read from M3D-C1 file using fusion-io (2) generic test data (no file) More...
 
logical ptb_3db_single_mode_input =.true.
 For M3D-C1 input; If .true., M3D-C1 file contains only one toroidal mode number. More...
 
integer ptb_3db_single_mode_ntor =3
 The mode number in the M3D-C1 file in case of single mode input. More...
 
real(kind=8) ptb_3db_mult_factor =1D0
 Factor for scaling the perturbed field. More...
 
integer ptb_3db_m3dc1_timeslice =1
 Time slice that is to be read from M3D-C1 file. More...
 
logical ptb_3db_em_full_spec =.false.
 If false, only modes |m/q-n|<=sml_mode_select_mres_q are retained, if true, the remaining modes are retained in toroidal Fourier representation but without RMP response (EM only). More...
 
real(kind=8) ptb_3db_em_inpsi =-1D0
 Inner boundary of RMP representation in A_s. More...
 
real(kind=8) ptb_3db_em_outpsi =1.1D0
 Outer boundary of RMP representation in A_s. More...
 
real(kind=8) ptb_3db_em_bd_width =0.05D0
 Blending width between A_s and full A representation. More...
 
integer ptb_3db_dasdt_opt =1
 Selector for how to compute dAs/dt for RMP penetration. More...
 
integer ptb_3db_mstep_es =0
 Number of ES time steps in RMP penetration calculation with ptb_3db_mode==2. More...
 
integer ptb_3db_mstep_em =1
 Number of EM time steps in RMP penetration calculation with ptb_3db_mode==2. More...
 
integer ptb_3db_es_to_em_dt_ratio =1
 Ratio of ES to EM time step size in RMP penetration calculation with ptb_3db_mode==2 (also: the number of time steps over which the perturbed current is averaged) More...
 
real(kind=8) ptb_3db_update_alpha =1.0D0
 Damping factor for update of the perturbed field. More...
 
integer ptb_3db_start_time =1
 Time step in which perturbed field is switched on. More...
 
logical ptb_3db_rampup_vac =.true.
 (.true.) Ramp up perturbed field slowly, (.false.) turn on perturbed field abruptly More...
 
integer ptb_3db_rampup_time =1000
 Number of time steps over which the perturbed field is ramped up. More...
 
real(kind=8) ptb_3db_rampup_fac =0D0
 The current relative amplitude of the vacuum vector potential (if ptb_3db_mode==2) More...
 
real(kind=8) ptb_3db_rampup_fac0 =0D0
 The previous (backup) relative amplitude of the vacuum vector potential (if ptb_3db_mode==2) More...
 
logical ptb_3db_screening_on =.false.
 Artificial screening for static RMP field: (.false.) no screening, (.true.) exponential damping towards the axis. More...
 
real(kind=8) ptb_3db_screening_width1 =0.02D0
 Decay length of the exponential damping function in units of normalized poloidal flux. More...
 
real(kind=8) ptb_3db_screening_width2 =0.02D0
 Width of tanh for damping in the core in units of normalized poloidal flux. More...
 
real(kind=8) ptb_3db_screening_psi1 =0.975D0
 Center of the tanh function for screening in the core (<1) More...
 
real(kind=8) ptb_3db_screening_fac1 =0.9D0
 Reduction of deltaB at the separatrix (<=1) More...
 
real(kind=8) ptb_3db_screening_fac2 =0.8D0
 Fraction of vacuum RMP in the core (<=1) More...
 
integer ptb_3db_num_ntor =1
 Number of toroidal mode numbers (<= sml_nphi_total) More...
 
integer, dimension(:),
allocatable, target 
ptb_3db_ntor
 Array to store the toroidal mode numbers. More...
 
integer ptb_3db_ntor_min =3
 smallest toroidal mode number More...
 
integer ptb_3db_num_mpol =30
 Number of poloidal mode numbers. More...
 
integer, dimension(:), allocatable ptb_3db_mpol
 Array to store the poloidal mode numbers. More...
 
integer ptb_3db_nr =1000
 (R,Z) grid size in R-direction More...
 
integer ptb_3db_nz =1000
 (R,Z) grid size in Z-direction More...
 
real(kind=8) ptb_3db_r_min
 
real(kind=8) ptb_3db_r_max
 Range of R in (R,Z) grid. More...
 
real(kind=8) ptb_3db_z_min
 
real(kind=8) ptb_3db_z_max
 Range of Z in (R,Z) grid. More...
 
real(kind=8) ptb_3db_dr
 
real(kind=8) ptb_3db_dz
 grid spacing of (R,Z) grid More...
 
real(kind=8), dimension(:,:,:),
allocatable, target 
ptb_3db_bfield_re_vac
 perturbed vacuum field on XGC mesh, real part, dimensions: (grid vertex, R-Z-phi components, tor. mode number) More...
 
real(kind=8), dimension(:,:,:),
allocatable, target 
ptb_3db_bfield_im_vac
 perturbed vacuum field on XGC mesh, imaginary part More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_vecpot_re_vac_rz
 perturbed vacuum vector potential on (R,Z) mesh, real part More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_vecpot_im_vac_rz
 perturbed vacuum vector potential on (R,Z) mesh, imaginary part More...
 
real(kind=8), dimension(:),
allocatable 
ptb_3db_rgrid
 
real(kind=8), dimension(:),
allocatable 
ptb_3db_zgrid
 R and Z grid in case of (R,Z) input data. More...
 
type(ezspline2_r8), dimension(:,:),
allocatable 
ptb_3db_vecpot_spl_rz_re
 Spline for interpolation of the real part of (R,Z) data to XGC mesh. More...
 
type(ezspline2_r8), dimension(:,:),
allocatable 
ptb_3db_vecpot_spl_rz_im
 Spline for interpolation of the imaginary part of (R,Z) data to XGC mesh. More...
 
type(ezspline2_r8) ptb_3db_spl_tmp
 Dummy 2D spline used for setting up the 2D interpolations with bicub_mod. More...
 
type(bicub_type), dimension(:,:,:),
allocatable 
ptb_3db_bicub_da
 Spline for interpolation of the perturbed vector potential with bicub_mod. More...
 

Detailed Description

The module "ptb_3db" contains the framework for handling perturbed 3D magnetic fields in XGCa and XGC1.

Member Function/Subroutine Documentation

subroutine ptb_3db_module::extract_resonant_vector_potential_from_m3dc1 ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn,
integer, intent(in)  itor 
)

Extension of RZ reader for EM Use RZ reader to read the full vector potential from M3D-C1. For EM, proeceed to extract the (near) resonant components of \(A_\parallel\) (A_parallel) on the XGC mesh. Interpolate those back to the RZ grid and subtract from the raw vector potential from M3D-C1. Then re-initialize the spline interpolation for A. This enables pushing particles with the full M3D-C1 field while still calculating the plasma response to the near resonant modes. As a stop-gap measure, the such filtered vector potential is used to calculate the perturbed field on the XGC mesh because the spline interpolation of the vector potential on the RZ grid is not available in C++ yet. This currently works only for single toroidal mode number input. Extending it to multiple mode numbers seems straightforward.

Parameters
[in]gridXGC mesh object, type(grid_type)
[in,out]psnXGC field object, type(psn_type)
[in]itorIndex of targeted mode number, integer

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::init_apar_3db ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn 
)

Sets the symplectic vector potential equal to the vacuum RMP vector potential.

[in] grid type(grid_type), XGC grid information

Parameters
[in,out]psnField data object, type(psn_type)

Here is the caller graph for this function:

subroutine ptb_3db_module::psn_mem_alloc_3db ( type(psn_type)  psn,
integer  n 
)

Allocates memory for the perturbed current in a psn-object.

Parameters
psn(in) type(psn_type), field data
n(in) integer, the number of mesh vertices

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::ptb_3db_calc_db ( type(grid_type), intent(in)  grid,
real (kind=8), dimension(grid%nnode,2,ptb_3db_num_ntor), intent(in)  dA_phi,
real (kind=8), dimension(3,grid%nnode,ptb_3db_num_ntor), intent(out)  dB_re,
real (kind=8), dimension(3,grid%nnode,ptb_3db_num_ntor), intent(out)  dB_im 
)

Evaluates the perturbed field from the perturbed vector potential on the XGC grid.

Parameters
grid(in) type(grid_type), grid data
dA_phi(in) real(8), perturbed vector potential (dA_phi, covariant phi-component)
dB_re(out) real(8), real (cos) part of the perturbed mag. field
dB_im(out) real(8), imaginary (sin) part of the perturbed mag. field

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::ptb_3db_damping_factor_em ( type(grid_type), intent(in)  grid,
real(kind=8), dimension(grid%nnode), intent(in)  Apar,
real(kind=8), dimension(grid%nnode), intent(inout)  dApar 
)

Calculates an adaptive damping factor to limit the time rate of change of the perturbed vector potential The routine calculates the poloidal Fourier spectrum of the current increment of the vector potential (dApar –> effectively the Hamiltonian vector potential Ah) and the vector potential itself (Apar). Then it compares the amplitude of the increment and the current value and computes a damping factor for each poloidal mode number such that the increment is limited to be a fraction of the current value. Using this damping factor, the increment is transformed back to real space and returned to the calling routine.

Parameters
[in]gridXGC mesh object, type(grid_type)
[in]AparCurrent symplectic vector potential As, double
[in,out]dAparIncrement of the vector potential (=Ah), double

Here is the call graph for this function:

Here is the caller graph for this function:

real (kind=8) function ptb_3db_module::ptb_3db_screening_factor ( real (kind=8), intent(in)  r,
real (kind=8), intent(in)  z,
real (kind=8), intent(in)  psi 
)

Evaluates the screening factor for artificial exponential damping of the vacuum perturbed field towards the magnetic axis. The screening shape function is simply exp((psi-psi0)/w). Since psi decreases towards the magnetic axis, psi-psi0 becomes negative at some point.

Parameters
r(in) real(8), major radius
z(in) real(8), Z coordinate
psi(in) real(8), normalized poloidal flux corresponding to r and z

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::read_3db ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn 
)

read_3db allocates memory for the perturbed field and reads the perturbed field from a file. Supported file formats are (R,Z)-ASCII and M3D-C1.

Parameters
[in]gridtype(grid_type), XGC grid information
[in,out]psnXGC field data object, type(psn_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::read_3db_m3dc1_eem ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn 
)

Allocates memory for and reads perturbed field data from M3D-C1 file for use with the explicit EM method. This reader reads only the parallel component of the perturbed vacuum vector potential for a single toroidal mode number.

Parameters
[in]gridtype(grid_type), XGC grid information
[in,out]psnXGC field data object, type(psn_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::read_3db_m3dc1_mesh ( type(grid_type), intent(in)  grid)

Allocates memory for and reads perturbed field data from M3D-C1 file. For M3D-C1 input, the fusion-io library by N. Ferraro is used. The interpolation used in the fusion-io library keeps the field divergence free.

Parameters
grid(in) type(grid_type), XGC grid information

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::read_3db_m3dc1_rz ( type(grid_type), intent(in)  grid)

Allocates memory for and reads perturbed field data from M3D-C1 file. This version of the routine reads the vector potential A on an (R,Z) grid. The perturbed field on the XGC mesh can be evaluated in a separate subroutine. Eventually, the perturbed magnetic field is to be calculated directy from A using spline interpolation. For M3D-C1 input, the fusion-io library by N. Ferraro is used.

Parameters
grid(in) type(grid_type), XGC grid information

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::read_3db_test ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn 
)

Sets up an analytical magnetic field perturbation with fixed toroidal mode number and (currently) 2 poloidal mode numbers. The number of poloidal mode numbers can be adjusted with few code changes. The most important parameters are the value of the poloidal mode numbers, the centers of the Gaussian envelopes of each mode number, the width of the Gaussian, and the amplitude.

Parameters
[in]gridXGC grid data object, type(grid_type)
[in,out]psnXGC field data object, type(psn_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::update_3db ( )

Updates the ramp-up factor of the vacuum RMP field. The external (or vacuum) RMP field can be ramped up linearly instead of a hard switch from off to on. The ramp-up factor ptb_3db_rampup_fac determines the amplitude of the external RMP field throughout XGC. This routine is used only in electrostatic simulations. In electromagnetic simulations, push_As updates the ramp-up factor.

Here is the call graph for this function:

Here is the caller graph for this function:

Member Data Documentation

real (kind=8), dimension(:,:,:), allocatable, target ptb_3db_module::ptb_3db_bfield_im_vac

perturbed vacuum field on XGC mesh, imaginary part

real (kind=8), dimension(:,:,:), allocatable, target ptb_3db_module::ptb_3db_bfield_re_vac

perturbed vacuum field on XGC mesh, real part, dimensions: (grid vertex, R-Z-phi components, tor. mode number)

type (bicub_type), dimension(:,:,:), allocatable ptb_3db_module::ptb_3db_bicub_da

Spline for interpolation of the perturbed vector potential with bicub_mod.

integer ptb_3db_module::ptb_3db_dasdt_opt =1

Selector for how to compute dAs/dt for RMP penetration.

real (kind=8) ptb_3db_module::ptb_3db_dr
real (kind=8) ptb_3db_module::ptb_3db_dz

grid spacing of (R,Z) grid

real (kind=8) ptb_3db_module::ptb_3db_em_bd_width =0.05D0

Blending width between A_s and full A representation.

logical ptb_3db_module::ptb_3db_em_full_spec =.false.

If false, only modes |m/q-n|<=sml_mode_select_mres_q are retained, if true, the remaining modes are retained in toroidal Fourier representation but without RMP response (EM only).

real (kind=8) ptb_3db_module::ptb_3db_em_inpsi =-1D0

Inner boundary of RMP representation in A_s.

real (kind=8) ptb_3db_module::ptb_3db_em_outpsi =1.1D0

Outer boundary of RMP representation in A_s.

integer ptb_3db_module::ptb_3db_es_to_em_dt_ratio =1

Ratio of ES to EM time step size in RMP penetration calculation with ptb_3db_mode==2 (also: the number of time steps over which the perturbed current is averaged)

integer ptb_3db_module::ptb_3db_file_format =1

Format of input file with perturbed field data: (0) toroidal Fourier coefficients in R,Z, (1) read from M3D-C1 file using fusion-io (2) generic test data (no file)

character(len=250) ptb_3db_module::ptb_3db_filename ='C1.h5'

File containing the perturbed field data.

integer ptb_3db_module::ptb_3db_m3dc1_timeslice =1

Time slice that is to be read from M3D-C1 file.

integer ptb_3db_module::ptb_3db_mode =0

Mode of operation: (0) static RMP field (ES), or fully self-consistent RMP field (EM version) (2) self-consistent RMP field with time-averaged field-equations (EM version)

integer, dimension(:), allocatable ptb_3db_module::ptb_3db_mpol

Array to store the poloidal mode numbers.

integer ptb_3db_module::ptb_3db_mstep_em =1

Number of EM time steps in RMP penetration calculation with ptb_3db_mode==2.

integer ptb_3db_module::ptb_3db_mstep_es =0

Number of ES time steps in RMP penetration calculation with ptb_3db_mode==2.

real (kind=8) ptb_3db_module::ptb_3db_mult_factor =1D0

Factor for scaling the perturbed field.

integer ptb_3db_module::ptb_3db_nr =1000

(R,Z) grid size in R-direction

integer, dimension(:), allocatable, target ptb_3db_module::ptb_3db_ntor

Array to store the toroidal mode numbers.

integer ptb_3db_module::ptb_3db_ntor_min =3

smallest toroidal mode number

integer ptb_3db_module::ptb_3db_num_mpol =30

Number of poloidal mode numbers.

integer ptb_3db_module::ptb_3db_num_ntor =1

Number of toroidal mode numbers (<= sml_nphi_total)

integer ptb_3db_module::ptb_3db_nz =1000

(R,Z) grid size in Z-direction

logical ptb_3db_module::ptb_3db_on =.false.

Switch 3D perturbed field on (.true.) and off (.false.)

real (kind=8) ptb_3db_module::ptb_3db_r_max

Range of R in (R,Z) grid.

real (kind=8) ptb_3db_module::ptb_3db_r_min
real (kind=8) ptb_3db_module::ptb_3db_rampup_fac =0D0

The current relative amplitude of the vacuum vector potential (if ptb_3db_mode==2)

real (kind=8) ptb_3db_module::ptb_3db_rampup_fac0 =0D0

The previous (backup) relative amplitude of the vacuum vector potential (if ptb_3db_mode==2)

integer ptb_3db_module::ptb_3db_rampup_time =1000

Number of time steps over which the perturbed field is ramped up.

logical ptb_3db_module::ptb_3db_rampup_vac =.true.

(.true.) Ramp up perturbed field slowly, (.false.) turn on perturbed field abruptly

real (kind=8), dimension(:), allocatable ptb_3db_module::ptb_3db_rgrid
real (kind=8) ptb_3db_module::ptb_3db_screening_fac1 =0.9D0

Reduction of deltaB at the separatrix (<=1)

real (kind=8) ptb_3db_module::ptb_3db_screening_fac2 =0.8D0

Fraction of vacuum RMP in the core (<=1)

logical ptb_3db_module::ptb_3db_screening_on =.false.

Artificial screening for static RMP field: (.false.) no screening, (.true.) exponential damping towards the axis.

real (kind=8) ptb_3db_module::ptb_3db_screening_psi1 =0.975D0

Center of the tanh function for screening in the core (<1)

real (kind=8) ptb_3db_module::ptb_3db_screening_width1 =0.02D0

Decay length of the exponential damping function in units of normalized poloidal flux.

real (kind=8) ptb_3db_module::ptb_3db_screening_width2 =0.02D0

Width of tanh for damping in the core in units of normalized poloidal flux.

logical ptb_3db_module::ptb_3db_single_mode_input =.true.

For M3D-C1 input; If .true., M3D-C1 file contains only one toroidal mode number.

integer ptb_3db_module::ptb_3db_single_mode_ntor =3

The mode number in the M3D-C1 file in case of single mode input.

type (ezspline2_r8) ptb_3db_module::ptb_3db_spl_tmp

Dummy 2D spline used for setting up the 2D interpolations with bicub_mod.

integer ptb_3db_module::ptb_3db_start_time =1

Time step in which perturbed field is switched on.

real (kind=8) ptb_3db_module::ptb_3db_update_alpha =1.0D0

Damping factor for update of the perturbed field.

real (kind=8), dimension(:,:,:,:), allocatable ptb_3db_module::ptb_3db_vecpot_im_vac_rz

perturbed vacuum vector potential on (R,Z) mesh, imaginary part

real (kind=8), dimension(:,:,:,:), allocatable ptb_3db_module::ptb_3db_vecpot_re_vac_rz

perturbed vacuum vector potential on (R,Z) mesh, real part

type (ezspline2_r8), dimension(:,:), allocatable ptb_3db_module::ptb_3db_vecpot_spl_rz_im

Spline for interpolation of the imaginary part of (R,Z) data to XGC mesh.

type (ezspline2_r8), dimension(:,:), allocatable ptb_3db_module::ptb_3db_vecpot_spl_rz_re

Spline for interpolation of the real part of (R,Z) data to XGC mesh.

real (kind=8) ptb_3db_module::ptb_3db_z_max

Range of Z in (R,Z) grid.

real (kind=8) ptb_3db_module::ptb_3db_z_min
real (kind=8), dimension(:), allocatable ptb_3db_module::ptb_3db_zgrid

R and Z grid in case of (R,Z) input data.


The documentation for this module was generated from the following file: