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) |
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) |
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 | psn_mem_alloc_3db (psn, n, nrho) |
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 | update_3db (grid, psn) |
Updates the perturbed magnetic field. In self-consistent RMP simulations (ptb_3d_mode==1), the current calculated with XGC is used to update the perturbed magnetic field. 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 | ptb_3db_div_clean_init (grid) |
Initializes a solver for divergence cleaning of the perturbed magnetic field This equation is being solved: div.(R*phi_b) - (ntor^2/R)*phi_b = R*div.delta_B. More... | |
subroutine | ptb_3db_create_1field_solver (grid, solver, prefix) |
create 1 field solver with A00 built (Amat), just setoperator & setup like 2 field More... | |
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_solver_init (grid) |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!! We solve perturbed Ampere's law in the following form: !!!! dB = - grad(phi).grad(dA_phi) = grad(phi).grad(dpsi) !!!! –> grad(phi).curl(dB) = - grad(phi).curl(grad(phi).grad(dA_phi)) !!!! = div.(grad(phi) x (grad(phi) x grad(dA_phi)) !!!! = -div.[g^{phi,phi} * (grad(dA_phi) - d(dA_phi)/dphi grad(phi))] !!!! = -(1/R) * [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] !!!! = mu_0 j.grad(phi) !!!! –> [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] = -mu_0 R j^phi !!!! = -mu_0 (B_phi/B) j_parallel !!!! –> We need to set up one Helmholtz solver with alpha=1/R and beta=0 !!!! –> The RHS is -mu_0*(B_phi/B)*j_parallel !!!! –> j_parallel = psnijpar + psnejpar !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! More... | |
subroutine | ptb_3db_da_bd_init (grid) |
Initializes the solver boundary for divergence cleaning. More... | |
subroutine | ptb_3db_da_solve (grid, psn, dA_phi, rhs) |
Solves perturbed Ampere's law using psnijpar and psnejpar. 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. 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 =0.1D0 |
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_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::psn_mem_alloc_3db | ( | type(psn_type) | psn, |
integer | n, | ||
integer | nrho | ||
) |
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 |
nrho | (in) integer, size of rho_i grid |
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_create_1field_solver | ( | type(grid_type), intent(in) | grid, |
type(xgc_solver), intent(inout) | solver, | ||
character(*), intent(in) | prefix | ||
) |
create 1 field solver with A00 built (Amat), just setoperator & setup like 2 field
grid | (in) type(grid_type), grid info |
solver | (inout) type(xgc_solver_type), XGC data structure for the solver |
prefix | (in) char(*), prefix for the solver |
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_da_solve | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(in) | psn, | ||
real (kind=8), dimension(grid%nnode,2,ptb_3db_num_ntor), intent(out) | dA_phi, | ||
real (kind=8), dimension(grid%nnode,2,ptb_3db_num_ntor), intent(inout) | rhs | ||
) |
Solves perturbed Ampere's law using psnijpar and psnejpar.
grid[in] | type(grid_type), grid info | |
psn[in] | type(psn_type), field information (currents) | |
[out] | dA_phi | real(8), new perturbed vector potential dA_phi |
[out] | rhs | RHS vector used in the solve |
subroutine ptb_3db_module::ptb_3db_da_solver_init | ( | type(grid_type), intent(in) | grid | ) |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!! We solve perturbed Ampere's law in the following form: !!!! dB = - grad(phi).grad(dA_phi) = grad(phi).grad(dpsi) !!!! –> grad(phi).curl(dB) = - grad(phi).curl(grad(phi).grad(dA_phi)) !!!! = div.(grad(phi) x (grad(phi) x grad(dA_phi)) !!!! = -div.[g^{phi,phi} * (grad(dA_phi) - d(dA_phi)/dphi grad(phi))] !!!! = -(1/R) * [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] !!!! = mu_0 j.grad(phi) !!!! –> [d_R((1/R)*d_R(dA_phi)) + d_Z((1/R)*d_Z(dA_phi))] = -mu_0 R j^phi !!!! = -mu_0 (B_phi/B) j_parallel !!!! –> We need to set up one Helmholtz solver with alpha=1/R and beta=0 !!!! –> The RHS is -mu_0*(B_phi/B)*j_parallel !!!! –> j_parallel = psnijpar + psnejpar !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
grid | (in) type(grid_type), grid information |
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_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_clean_init | ( | type(grid_type), intent(in) | grid | ) |
Initializes a solver for divergence cleaning of the perturbed magnetic field This equation is being solved: div.(R*phi_b) - (ntor^2/R)*phi_b = R*div.delta_B.
grid | (in) type(grid_type), grid information |
ntor | (in) integer, toroidal mode number |
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 | ) |
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.
grid | (in) type(grid_type), XGC grid information |
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_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 | ) |
subroutine ptb_3db_module::update_3db | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(in) | psn | ||
) |
Updates the perturbed magnetic field. In self-consistent RMP simulations (ptb_3d_mode==1), the current calculated with XGC is used to update the perturbed magnetic field.
grid | (in) type(grid_type), grid data |
psn | (in) type(psn_type), field data (currents) |
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_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.
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 =0.1D0 |
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.