XGCa
 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_rz (grid)
 Allocates memory for and reads perturbed field data from (R,Z)-ASCII file. Does not need any input variables. (R,Z)-ASCII files are expected in the following format Line 1: Header (description of file contents) Line 2: '(3i8)' —> ptb_3db_ntor – ptb_3db_nr – ptb_3db_nz Line 3: '(4(f,1x))' —> ptb_3db_r_min – ptb_3db_r_max – ptb_3db_z_min – ptb_3db_z_max Line 4: empty Then for each toroidal mode number in the file: '(i8)' —> tor. mode number n Linewise (all R at with the same Z, then next Z from small to large R and Z) '(6f)' —> Re[dB_R_n] – Re[dB_Z_n] – Re[dB_phi_n] – Im[dB_R_n] – Im[dB_Z_n] – Im[dB_phi_n] Empty line: '#...' Next toroidal mode number. More...
 
subroutine read_3db_m3dc1 (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_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 ptb_3db_damping_factor (grid, damp_fac, delta_dA, rhs)
 Calculates Fourier decomposition of the change of dA_phi and determines the damping factor for a given iteration of the 3D-field Picard iteration. More...
 
subroutine update_3db (grid, psn, ipc, adiabatic)
 Updates the perturbed electric and magnetic field. In self-consistent EM RMP simulations (ptb_3d_mode==2), the Poisson equation and Ampere's law are solved only periodically, and a damping function is applied to remove sudden changes of the solution over time. More...
 
subroutine ptb_3db_damping_factor_em (grid, Apar, dApar, damp_fac)
 
subroutine ptb_3db_div_clean_bd_init (grid)
 Initializes the solver boundary for divergence cleaning. More...
 
subroutine ptb_3db_div_cleaning (grid, db_re, db_im)
 Performs divergence cleaning of the perturbed field that may have been introduced by interpolation to the XGCa grid. Solves d_R(R*d_R(phi_b)) + d_Z(R*d_Z(phi_b)) - n^2/R phi_b = r*div(dB) for all toroidal mode numbers in the system. More...
 
subroutine ptb_3db_r_div_b (grid, db_re, db_im, itor, r_div_db)
 Evaluates the divergence of a single Fourier component of the perturbed magnetic field db (sine+cosine parts for a single n). More...
 
subroutine ptb_3db_da_bd_init (grid)
 Initializes the solver boundary for divergence cleaning. 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, (1) self-consistent RMP field (ES version), (2) self-consistent RMP field (EM version) More...
 
character(len=250) ptb_3db_filename ='3db.dat'
 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. 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_force_symmetric_f =.false.
 Whether to force f0_f and f0_f0g to be axisymmetric. More...
 
integer ptb_3db_field_update_period =100
 Period in time steps between two subsequent updates of the perturbed field (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...
 
integer ptb_3db_rampup_interval =100
 Number of time steps over which the perturbed field is constant during ramp-up The rampup, thus, proceeds in ptb_3db_rampup_time/ptb_3db_rampup_interval steps. 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...
 
logical ptb_3db_use_gyro_avg =.false.
 If .true., use gyro-averaged perturbed magnetic field for ions. 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 =1
 Number of poloidal mode numbers. More...
 
integer, dimension(:), allocatable ptb_3db_mpol
 Array to store the poloidal mode numbers. More...
 
integer ptb_3db_mpol_min =1
 smallest poloidal mode number 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_re
 perturbed magnetic field on XGC mesh, real part 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, target 
ptb_3db_bfield_im
 perturbed magnetic field on XGC mesh, imaginary part More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_bfield_re_vac_rz
 perturbed vacuum field on (R,Z) mesh, real part, dimensions: (R-index,Z-index, R-Z-phi components, tor. mode number) More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_bfield_im_vac_rz
 perturbed vacuum field on (R,Z) mesh, imaginary part More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_bfield_re_rz
 perturbed field on (R,Z) mesh, real part More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_bfield_im_rz
 perturbed field on (R,Z) 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_vecpot_re_rz
 perturbed vector potential on (R,Z) mesh, real part More...
 
real(kind=8), dimension(:,:,:,:),
allocatable 
ptb_3db_vecpot_im_rz
 perturbed vector potential on (R,Z) mesh, imaginary part More...
 
real(kind=8), dimension(:,:,:),
allocatable 
ptb_3db_da_phi
 covariant phi-component of the perturbed, induced vector potential; for basis grad(R), grad(phi), grad(Z) –> dpsi = R*dA More...
 
real(kind=8), dimension(:,:,:),
allocatable 
ptb_3db_da_phi_vac
 covariant phi-component of the perturbed, induced vector potential: for basis grad(R), grad(phi), grad(Z) –> dpsi = R*dA More...
 
real(kind=8), dimension(:,:,:),
allocatable 
ptb_3db_da_par_vac
 Parallel component of the perturbed vacuum vector potential (B.dA)/|B| for use with explicit EM method. 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_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_spl_rz_im
 Spline for interpolation of the imaginary part of (R,Z) data to XGC mesh. 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...
 
integer ptb_3db_jpar_avg_period =100
 Averaging period for the perturbed curreint in time steps. More...
 
real(kind=8) ptb_3db_diag_mn_inpsi =0.1D0
 inner boundary for Fourier diagnostic More...
 
real(kind=8) ptb_3db_diag_mn_outpsi =0.995D0
 outer boundary for Fourier diagnostic (must be smaller than sml_outpsi) More...
 
real(kind=8) ptb_3db_solve_mn_inpsi =0.1D0
 inner boundary for Fourier filter used by the Ampere solver More...
 
real(kind=8) ptb_3db_solve_mn_outpsi =0.995D0
 outer boundary for Fourier filter used by the Ampere solver (must be smaller than sml_outpsi) More...
 
integer ptb_3db_solve_mpol_min =1
 Minimal poloidal mode number to be retained in the Amper solver. More...
 
integer ptb_3db_solve_mpol_max =25
 Maximal poloidal mode number to be retained in the Amper solver. More...
 
logical ptb_3db_div_clean_on =.false.
 .true.: Use divergence cleaning, .false.: no divergence cleaning More...
 
type(xgc_solver) ptb_3db_div_clean_solver
 Poisson solver for divergence cleaning. More...
 
type(xgc_solver) ptb_3db_da_solver
 Poisson solver for solving Ampere's law. More...
 
type(boundary2_type) ptb_3db_div_clean_bd
 Boundary for divergence cleaning solver (RHS boundary) More...
 
type(boundary2_type) ptb_3db_div_clean_bd2
 Boundary for divergence cleaning solver (Solver boundary) More...
 
integer ptb_3db_div_clean_bc_inode1
 index of first local boundary condition More...
 
integer ptb_3db_div_clean_bc_inode2
 index of last local boundary condition More...
 
real(kind=8), dimension(:,:),
allocatable 
ptb_3db_div_clean_bc_weights
 Weights for the calculation of boundary conditions. More...
 
type(boundary2_type) ptb_3db_da_bd
 Boundary for Ampere's law solver. More...
 
type(boundary2_type) ptb_3db_da_bd2
 Boundary for Ampere's law solver. More...
 
integer ptb_3db_da_bc_inode1
 index of first local boundary condition More...
 
integer ptb_3db_da_bc_inode2
 index of last local boundary condition More...
 
real(kind=8), dimension(:,:),
allocatable 
ptb_3db_da_bc_weights
 Weights for the calculation of boundary conditions. More...
 
real(kind=8) ptb_3db_da_inpsi =0D0
 Inner boundary for perturbed Ampere solver. More...
 
real(kind=8) ptb_3db_da_outpsi =1D0
 Outer boundary for perturbed Ampere solver. 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::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)
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,ptb_3db_num_ntor,grid%nnode), intent(out)  dB_re,
real (kind=8), dimension(3,ptb_3db_num_ntor,grid%nnode), 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_da_bd_init ( type(grid_type), intent(in)  grid)

Initializes the solver boundary for divergence cleaning.

Parameters
grid(in) type(grid_type), grid info

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::ptb_3db_damping_factor ( type(grid_type), intent(in)  grid,
real (kind=8), intent(out)  damp_fac,
real (kind=8), dimension(grid%nnode,2,ptb_3db_num_ntor), intent(in)  delta_dA,
real (kind=8), dimension(grid%nnode,2,ptb_3db_num_ntor), intent(in)  rhs 
)

Calculates Fourier decomposition of the change of dA_phi and determines the damping factor for a given iteration of the 3D-field Picard iteration.

Parameters
[in]gridgrid data, type(grid_type)
[out]damp_facglobal damping factor, real(8)
[in]delta_dAchange of the vector potential, real(8)
[in]rhsRHS vector used for the Ampere-solver, real(8)

Here is the call 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(in)  dApar,
real(kind=8), intent(out)  damp_fac 
)

Here is the call graph for this function:

Here is the caller graph for this function:

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

Initializes the solver boundary for divergence cleaning.

Parameters
grid(in) type(grid_type), grid info

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::ptb_3db_div_cleaning ( type(grid_type), intent(in)  grid,
real (kind=8), dimension(3,ptb_3db_num_ntor,grid%nnode), intent(inout)  db_re,
real (kind=8), dimension(3,ptb_3db_num_ntor,grid%nnode), intent(inout)  db_im 
)

Performs divergence cleaning of the perturbed field that may have been introduced by interpolation to the XGCa grid. Solves d_R(R*d_R(phi_b)) + d_Z(R*d_Z(phi_b)) - n^2/R phi_b = r*div(dB) for all toroidal mode numbers in the system.

Parameters
grid(in) type(grid_type), grid info
db_re(inout) real(8), cosine component of the perturbed field
db_im(inout) real(8), sine component of the perturbed field

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ptb_3db_module::ptb_3db_r_div_b ( type(grid_type), intent(in)  grid,
real (kind=8), dimension(3,ptb_3db_num_ntor,grid%nnode), intent(in)  db_re,
real (kind=8), dimension(3,ptb_3db_num_ntor,grid%nnode), intent(in)  db_im,
integer, intent(in)  itor,
real (kind=8), dimension(grid%nnode,2), intent(out)  r_div_db 
)

Evaluates the divergence of a single Fourier component of the perturbed magnetic field db (sine+cosine parts for a single n).

Parameters
grid(in) type(grid_type), grid info
db_re(in) real(8), cosine component of the perturbed field
db_im(in) real(8), sine component of the perturbed field
r_div_db(out) real(8), divergence of the perturbed field

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 ( 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_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_rz ( type(grid_type), intent(in)  grid)

Allocates memory for and reads perturbed field data from (R,Z)-ASCII file. Does not need any input variables. (R,Z)-ASCII files are expected in the following format Line 1: Header (description of file contents) Line 2: '(3i8)' —> ptb_3db_ntor – ptb_3db_nr – ptb_3db_nz Line 3: '(4(f,1x))' —> ptb_3db_r_min – ptb_3db_r_max – ptb_3db_z_min – ptb_3db_z_max Line 4: empty Then for each toroidal mode number in the file: '(i8)' —> tor. mode number n Linewise (all R at with the same Z, then next Z from small to large R and Z) '(6f)' —> Re[dB_R_n] – Re[dB_Z_n] – Re[dB_phi_n] – Im[dB_R_n] – Im[dB_Z_n] – Im[dB_phi_n] Empty line: '#...' Next toroidal mode number.

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 ( type(grid_type), intent(in)  grid,
type(psn_type), intent(inout)  psn,
integer, intent(in)  ipc,
logical, intent(in)  adiabatic 
)

Updates the perturbed electric and magnetic field. In self-consistent EM RMP simulations (ptb_3d_mode==2), the Poisson equation and Ampere's law are solved only periodically, and a damping function is applied to remove sudden changes of the solution over time.

Parameters
[in]gridtype(grid_type), grid data
[in]psntype(psn_type), field data (currents)
[in]ipcRK2 stage index, integer
[in]adiabaticwhether the function is called after an adiabatic or full n=0 Poisson solve

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

perturbed magnetic field on XGC mesh, imaginary part

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

perturbed field on (R,Z) mesh, imaginary part

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 ptb_3db_module::ptb_3db_bfield_im_vac_rz

perturbed vacuum field on (R,Z) mesh, imaginary part

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

perturbed magnetic field on XGC mesh, real part

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

perturbed field on (R,Z) mesh, real 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)

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

perturbed vacuum field on (R,Z) mesh, real part, dimensions: (R-index,Z-index, 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_da_bc_inode1

index of first local boundary condition

integer ptb_3db_module::ptb_3db_da_bc_inode2

index of last local boundary condition

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

Weights for the calculation of boundary conditions.

type(boundary2_type) ptb_3db_module::ptb_3db_da_bd

Boundary for Ampere's law solver.

type(boundary2_type) ptb_3db_module::ptb_3db_da_bd2

Boundary for Ampere's law solver.

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

Inner boundary for perturbed Ampere solver.

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

Outer boundary for perturbed Ampere solver.

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

Parallel component of the perturbed vacuum vector potential (B.dA)/|B| for use with explicit EM method.

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

covariant phi-component of the perturbed, induced vector potential; for basis grad(R), grad(phi), grad(Z) –> dpsi = R*dA

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

covariant phi-component of the perturbed, induced vector potential: for basis grad(R), grad(phi), grad(Z) –> dpsi = R*dA

type(xgc_solver) ptb_3db_module::ptb_3db_da_solver

Poisson solver for solving Ampere's law.

real (kind=8) ptb_3db_module::ptb_3db_diag_mn_inpsi =0.1D0

inner boundary for Fourier diagnostic

real (kind=8) ptb_3db_module::ptb_3db_diag_mn_outpsi =0.995D0

outer boundary for Fourier diagnostic (must be smaller than sml_outpsi)

integer ptb_3db_module::ptb_3db_div_clean_bc_inode1

index of first local boundary condition

integer ptb_3db_module::ptb_3db_div_clean_bc_inode2

index of last local boundary condition

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

Weights for the calculation of boundary conditions.

type(boundary2_type) ptb_3db_module::ptb_3db_div_clean_bd

Boundary for divergence cleaning solver (RHS boundary)

type(boundary2_type) ptb_3db_module::ptb_3db_div_clean_bd2

Boundary for divergence cleaning solver (Solver boundary)

logical ptb_3db_module::ptb_3db_div_clean_on =.false.

.true.: Use divergence cleaning, .false.: no divergence cleaning

type(xgc_solver) ptb_3db_module::ptb_3db_div_clean_solver

Poisson solver for divergence cleaning.

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

grid spacing of (R,Z) grid

integer ptb_3db_module::ptb_3db_field_update_period =100

Period in time steps between two subsequent updates of the perturbed field (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.

character(len=250) ptb_3db_module::ptb_3db_filename ='3db.dat'

File containing the perturbed field data.

logical ptb_3db_module::ptb_3db_force_symmetric_f =.false.

Whether to force f0_f and f0_f0g to be axisymmetric.

integer ptb_3db_module::ptb_3db_jpar_avg_period =100

Averaging period for the perturbed curreint in time steps.

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, (1) self-consistent RMP field (ES version), (2) self-consistent RMP field (EM version)

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

Array to store the poloidal mode numbers.

integer ptb_3db_module::ptb_3db_mpol_min =1

smallest poloidal mode number

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 =1

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
integer ptb_3db_module::ptb_3db_rampup_interval =100

Number of time steps over which the perturbed field is constant during ramp-up The rampup, thus, proceeds in ptb_3db_rampup_time/ptb_3db_rampup_interval steps.

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.

real (kind=8) ptb_3db_module::ptb_3db_solve_mn_inpsi =0.1D0

inner boundary for Fourier filter used by the Ampere solver

real (kind=8) ptb_3db_module::ptb_3db_solve_mn_outpsi =0.995D0

outer boundary for Fourier filter used by the Ampere solver (must be smaller than sml_outpsi)

integer ptb_3db_module::ptb_3db_solve_mpol_max =25

Maximal poloidal mode number to be retained in the Amper solver.

integer ptb_3db_module::ptb_3db_solve_mpol_min =1

Minimal poloidal mode number to be retained in the Amper solver.

type (ezspline2_r8), dimension(:,:), allocatable ptb_3db_module::ptb_3db_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_spl_rz_re

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

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.

logical ptb_3db_module::ptb_3db_use_gyro_avg =.false.

If .true., use gyro-averaged perturbed magnetic field for ions.

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

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

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_rz

perturbed vector potential on (R,Z) mesh, real 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: