XGC1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator 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)
 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...
 

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::psn_mem_alloc_3db ( type(psn_type)  psn,
integer  n,
integer  nrho 
)

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
nrho(in) integer, size of rho_i grid

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_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

Parameters
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

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_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.

Parameters
grid[in]type(grid_type), grid info
psn[in]type(psn_type), field information (currents)
[out]dA_phireal(8), new perturbed vector potential dA_phi
[out]rhsRHS vector used in the solve

Here is the call graph for this function:

Here is the caller graph for this function:

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 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

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

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:

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_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.

Parameters
grid(in) type(grid_type), grid information
ntor(in) integer, toroidal mode number

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)

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
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 ( 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_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)

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(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.

Parameters
grid(in) type(grid_type), grid data
psn(in) type(psn_type), field data (currents)

Here is the call 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_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.


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