XGCa
|
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... | |
The module "ptb_3db" contains the framework for handling perturbed 3D magnetic fields in XGCa and XGC1.
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
[in,out] | psn | Field 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.
psn | (in) type(psn_type), field data |
n | (in) integer, the number of mesh vertices |
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.
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 |
subroutine ptb_3db_module::ptb_3db_da_bd_init | ( | type(grid_type), intent(in) | grid | ) |
Initializes the solver boundary for divergence cleaning.
grid | (in) type(grid_type), grid info |
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.
[in] | grid | grid data, type(grid_type) |
[out] | damp_fac | global damping factor, real(8) |
[in] | delta_dA | change of the vector potential, real(8) |
[in] | rhs | RHS vector used for the Ampere-solver, real(8) |
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 | ||
) |
subroutine ptb_3db_module::ptb_3db_div_clean_bd_init | ( | type(grid_type), intent(in) | grid | ) |
Initializes the solver boundary for divergence cleaning.
grid | (in) type(grid_type), grid info |
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.
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 |
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).
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 |
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.
r | (in) real(8), major radius |
z | (in) real(8), Z coordinate |
psi | (in) real(8), normalized poloidal flux corresponding to r and z |
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.
[in] | grid | type(grid_type), XGC grid information |
[in,out] | psn | XGC field data object, type(psn_type) |
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.
grid | (in) type(grid_type), XGC grid information |
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.
[in] | grid | type(grid_type), XGC grid information |
[in,out] | psn | XGC field data object, type(psn_type) |
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.
grid | (in) type(grid_type), XGC grid information |
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.
[in] | grid | XGC grid data object, type(grid_type) |
[in,out] | psn | XGC field data object, type(psn_type) |
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.
[in] | grid | type(grid_type), grid data |
[in] | psn | type(psn_type), field data (currents) |
[in] | ipc | RK2 stage index, integer |
[in] | adiabatic | whether the function is called after an adiabatic or full n=0 Poisson solve |
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.