XGCa
Public Member Functions | Public Attributes | Static Public Attributes | List of all members
Simulation< Device > Class Template Reference

#include <sml.hpp>

Public Member Functions

 Simulation (NLReader::NamelistReader &nlr, bool reduced_deltaf_in=false, double transit_time_in=1.0, double psi_norm=1.0)
 

Public Attributes

bool reduced_deltaf
 True if any species is reduced_deltaf. Will be removed. More...
 
bool separate_n0
 Use separate n0 solver even in reduced_deltaf. More...
 
int bounce
 Bounce routine switch 0 for off, 1 for inner boundary, 2 for both boundaries. More...
 
int bounce_zero_weight
 If ==1 and bounce>0, set particle weights to zero after bouncing from the outer boundary. More...
 
double bounce_buffer
 Buffer width between sml_outpsi and where the particle actually bounces (must be >=0) More...
 
double dt
 Time step (s) More...
 
double transit_time
 Torodial transit time of an ion with the characteristic energy. More...
 
int sheath_mode
 Type of sheath (0 is none) More...
 
bool sheath_adjust
 Whether to adjust the sheath. More...
 
bool sheath_global_balance
 Whether to use the global loss balance functionality in sheath_mode=2. More...
 
bool ptb_3db_on
 Toggle for ptb_3db. More...
 
bool drift_on
 Toggle for using drift. More...
 
bool neutrals
 Toggle for using neutrals. More...
 
bool diag_heat_on
 Use heat diagnostics in sheath. More...
 
int mstep
 Max number of time steps. More...
 
int nrk
 Order of Runge-Kutta time integration of particles + fields. More...
 
bool electron_on
 Use kinetic electrons. More...
 
bool limit_marker_den
 Whether to limit marker density. More...
 
int nthreads
 Number of OMP threads on the host. More...
 
bool field_solver
 If false, charge deposition and field solve is skipped. More...
 
bool thermal_bath_on
 Switch on thermal-bath and coarse-graining operation, see ccm_param parameters. More...
 
bool dwdt_fix_bg
 
bool dwdt_exb_only
 
bool em_mixed_variable
 Switch for use of mixed-variable formulation. More...
 
bool em_control_variate
 Switch for use of control variate method. More...
 
PullbackMethod em_pullback_method
 Electrostatic: mixed-variable pullback with dA_s/dt=0,. More...
 
double em_pullback_dampfac
 Damping term gamma on -b.grad(phi) in pullback mode 4. More...
 
bool em_control_variate_final_cv
 
bool em_control_variate_flag
 
int em_control_variate_niter
 Number of iterations for Ampere solve with control-variate method. More...
 
bool em_b_para_eff
 Effective dB_|| via a modified grad-B drift (Joiner et al., PoP, 2010) More...
 
bool em_es_step
 
bool em_dAsdt_filter_on
 Switch for applying Fourier filters on RHS of dAs/dt equation (pullback mode 4) More...
 
bool no_turb
 Set all non-axisymmetric field perturbations to zero (electromagnetic version only) More...
 
bool em_n0
 Include n=0 electromagnetic mode. More...
 
bool em_dAsdt_hypvis
 Use radial hyperviscosity in dA_s/dt (push_As) More...
 
bool exclude_private
 Whether to exclude private region in ion charge deposition. More...
 
bool grad_psitheta
 
bool iter_solver
 
bool update_poisson_solver
 Whether poisson solver is updated. More...
 
bool iter_solver_converge
 
bool positive_phi00_sol
 
bool heuristic_priv_pot
 
int add_pot0
 
bool poisson_use_bc
 (XGCa only) False: 00-bc is phi=0 everywhere; true: 00-boundary More...
 
bool dpot_bd_apply
 Damp (n=0,m>0) potential towards the magnetic axis. More...
 
double dpot_bd_width
 Decay length (normalized flux) for (n=0,m>0) potential towards the magnetic axis. More...
 
bool poisson_bias
 Whether to use a (0,0) bias potential on top of phi_00. More...
 
bool poisson_adia_wall
 EXPERIMENTAL! DON'T USE UNLESS ADVISED BY AN EXPERT. More...
 
double dpot_te_limit
 Max absolute value of dpot/temp in getf0. More...
 
double dpot_te_limit_n0
 Limits the magnitude of the normalized axisymmetric potential e*dphi_0/T_e. More...
 
bool em_use_dpot_te_limit
 In EM simulation: whether to force usage of the min-max limiter on the turbulent potential fluctuation. More...
 
bool no_fp_in_f
 
double f0_grid_alpha
 
bool f0_update_analytic
 
bool f0_update_analytic_local
 If .false. --> flux-surface average update, .true. --> local. More...
 
double f0_update_analytic_alpha
 Separate alpha from sml_f0_grid_alpha for Maxwellian contribution. More...
 
double f0_update_analytic_damp_width
 
int f0_grid_alpha_start
 
bool ignore_f0g
 
bool symmetric_f
 Enforce axisymmetry of the total distribution function on the grid. More...
 
double update_g_alpha
 
bool diff_on
 
bool resamp_on
 Do resampling. More...
 
bool resamp_for_final_restart_write
 Perform resampling before dumping the final restart file. More...
 
bool resamp_restart_read
 Whether to read a restart file written from a simulation with different grid. More...
 
double time
 Current simulation time. More...
 
int gstep
 Current time step. More...
 
bool multirate_timestepping
 Use multirate timestepping. More...
 
int mr_factor
 
double mr_psi_max
 maximum normalized psi value of each multirate region More...
 
bool use_em_bounds
 
bool em_exclude_private
 
double bd_ext_delta_Ai
 
double bd_ext_delta_Ao
 
bool cce_coupling_on
 Core-edge coupling flag. More...
 
bool use_unfused_electron_push_kernel
 
bool loop_voltage_on
 Inductive current drive: loop voltage (from Faraday's law curl(E)=-dB/dt) More...
 
double loop_voltage_psimax
 Inductive current drive: outer boundary (in pol. flux) of the loop voltage. More...
 
bool current_drive_on
 
bool zero_inner_bd_turb
 

Static Public Attributes

static constexpr bool is_XGCa = true
 Equivalent to the preprocessor flag for now. More...
 
static constexpr bool explicit_electromagnetic = false
 Equivalent to the preprocessor flag for now. More...
 
static constexpr bool old_f0_update_analytic = false
 Equivalent to the preprocessor flag for now. More...
 
static constexpr int nhybrid = 2
 Number of iterations in electrostatic electron weight evolution scheme. More...
 

Constructor & Destructor Documentation

◆ Simulation()

template<class Device >
Simulation< Device >::Simulation ( NLReader::NamelistReader nlr,
bool  reduced_deltaf_in = false,
double  transit_time_in = 1.0,
double  psi_norm = 1.0 
)

Constructor for simulation

< Whether to solve for fields. Scatter and solve are skipped if set to false

< Split Poisson equation into axisymmetric and

< If .true., use Dirichlet boundary condition \(\delta\phi=0\), ! if .false., set the magnetic axis to the flux-surface average of the ! first flux-surface (provides continuity of the m=0 mode). ! (requires the inner boundaries be set to <0)

Here is the call graph for this function:

Member Data Documentation

◆ add_pot0

template<class Device >
int Simulation< Device >::add_pot0

Additional electrostatic potential: (0) for off, (1) for reading data from file !! (2) for simple neoclassical axisymmetric field (not operational)

◆ bd_ext_delta_Ai

template<class Device >
double Simulation< Device >::bd_ext_delta_Ai

◆ bd_ext_delta_Ao

template<class Device >
double Simulation< Device >::bd_ext_delta_Ao

◆ bounce

template<class Device >
int Simulation< Device >::bounce

Bounce routine switch 0 for off, 1 for inner boundary, 2 for both boundaries.

◆ bounce_buffer

template<class Device >
double Simulation< Device >::bounce_buffer

Buffer width between sml_outpsi and where the particle actually bounces (must be >=0)

◆ bounce_zero_weight

template<class Device >
int Simulation< Device >::bounce_zero_weight

If ==1 and bounce>0, set particle weights to zero after bouncing from the outer boundary.

◆ cce_coupling_on

template<class Device >
bool Simulation< Device >::cce_coupling_on

Core-edge coupling flag.

◆ current_drive_on

template<class Device >
bool Simulation< Device >::current_drive_on

Whether to use dynamic current drive via adjustments of the loop voltage ! The controller output is implemented as ! \(V_{loop} = 2\pi R \eta_\parallel P \left[\Delta j_\parallel + I \int_0^t \Delta j_\parallel \mathrm{d}\tau + D \frac{\partial \Delta j_\parallel}{\partial t} \right] \) ! Time integral and derivative are with respect to the ion time step instead of time, so \(1/I\) is the integration time scale, ! and \(D\) is the derivative time scale

◆ diag_heat_on

template<class Device >
bool Simulation< Device >::diag_heat_on

Use heat diagnostics in sheath.

◆ diff_on

template<class Device >
bool Simulation< Device >::diff_on

◆ dpot_bd_apply

template<class Device >
bool Simulation< Device >::dpot_bd_apply

Damp (n=0,m>0) potential towards the magnetic axis.

◆ dpot_bd_width

template<class Device >
double Simulation< Device >::dpot_bd_width

Decay length (normalized flux) for (n=0,m>0) potential towards the magnetic axis.

◆ dpot_te_limit

template<class Device >
double Simulation< Device >::dpot_te_limit

Max absolute value of dpot/temp in getf0.

◆ dpot_te_limit_n0

template<class Device >
double Simulation< Device >::dpot_te_limit_n0

Limits the magnitude of the normalized axisymmetric potential e*dphi_0/T_e.

◆ drift_on

template<class Device >
bool Simulation< Device >::drift_on

Toggle for using drift.

◆ dt

template<class Device >
double Simulation< Device >::dt

Time step (s)

◆ dwdt_exb_only

template<class Device >
bool Simulation< Device >::dwdt_exb_only

◆ dwdt_fix_bg

template<class Device >
bool Simulation< Device >::dwdt_fix_bg

◆ electron_on

template<class Device >
bool Simulation< Device >::electron_on

Use kinetic electrons.

◆ em_b_para_eff

template<class Device >
bool Simulation< Device >::em_b_para_eff

Effective dB_|| via a modified grad-B drift (Joiner et al., PoP, 2010)

◆ em_control_variate

template<class Device >
bool Simulation< Device >::em_control_variate

Switch for use of control variate method.

◆ em_control_variate_final_cv

template<class Device >
bool Simulation< Device >::em_control_variate_final_cv

◆ em_control_variate_flag

template<class Device >
bool Simulation< Device >::em_control_variate_flag

◆ em_control_variate_niter

template<class Device >
int Simulation< Device >::em_control_variate_niter

Number of iterations for Ampere solve with control-variate method.

◆ em_dAsdt_filter_on

template<class Device >
bool Simulation< Device >::em_dAsdt_filter_on

Switch for applying Fourier filters on RHS of dAs/dt equation (pullback mode 4)

◆ em_dAsdt_hypvis

template<class Device >
bool Simulation< Device >::em_dAsdt_hypvis

Use radial hyperviscosity in dA_s/dt (push_As)

◆ em_es_step

template<class Device >
bool Simulation< Device >::em_es_step

If this is true, A_s is fixed, A_h=0, the Poisson equation is solved in the ES form (with adiabatic response on the LHS), and the electron weight update uses the ES method, i.e. the ES algorithm is used in a constant perturbed magnetic field. This parameter is not an input but set dynamically through a control function (intended for RMP penetration calculation).

◆ em_exclude_private

template<class Device >
bool Simulation< Device >::em_exclude_private

◆ em_mixed_variable

template<class Device >
bool Simulation< Device >::em_mixed_variable

Switch for use of mixed-variable formulation.

◆ em_n0

template<class Device >
bool Simulation< Device >::em_n0

Include n=0 electromagnetic mode.

◆ em_pullback_dampfac

template<class Device >
double Simulation< Device >::em_pullback_dampfac

Damping term gamma on -b.grad(phi) in pullback mode 4.

◆ em_pullback_method

template<class Device >
PullbackMethod Simulation< Device >::em_pullback_method

Electrostatic: mixed-variable pullback with dA_s/dt=0,.

◆ em_use_dpot_te_limit

template<class Device >
bool Simulation< Device >::em_use_dpot_te_limit

In EM simulation: whether to force usage of the min-max limiter on the turbulent potential fluctuation.

◆ exclude_private

template<class Device >
bool Simulation< Device >::exclude_private

Whether to exclude private region in ion charge deposition.

◆ explicit_electromagnetic

template<class Device >
constexpr bool Simulation< Device >::explicit_electromagnetic = false
staticconstexpr

Equivalent to the preprocessor flag for now.

◆ f0_grid_alpha

template<class Device >
double Simulation< Device >::f0_grid_alpha

Fraction of the particle weight that is transferred to the background distribution function (particle noise control) every sml_f_source_period time steps. As a rule of thumb, sml_f0_grid_alpha=sml_dt*sml_f_source_period can be used.

◆ f0_grid_alpha_start

template<class Device >
int Simulation< Device >::f0_grid_alpha_start

For delayed onset of total-f particle noise control. (must be >=1) Values >1 cause a linear ramp-up of sml_f0_grid_alpha from 0 to its input value from time step 1 to sml_f0_grid_alpha_start.

◆ f0_update_analytic

template<class Device >
bool Simulation< Device >::f0_update_analytic

Switch on/off the update of the analytic part of the distribution function

◆ f0_update_analytic_alpha

template<class Device >
double Simulation< Device >::f0_update_analytic_alpha

Separate alpha from sml_f0_grid_alpha for Maxwellian contribution.

◆ f0_update_analytic_damp_width

template<class Device >
double Simulation< Device >::f0_update_analytic_damp_width

For width of \(\exp^{-(x/w)^2}\) damping factor for updating density/temperature of the analytical f0

◆ f0_update_analytic_local

template<class Device >
bool Simulation< Device >::f0_update_analytic_local

If .false. --> flux-surface average update, .true. --> local.

◆ field_solver

template<class Device >
bool Simulation< Device >::field_solver

If false, charge deposition and field solve is skipped.

◆ grad_psitheta

template<class Device >
bool Simulation< Device >::grad_psitheta

If true, gradiant operator is calculated as \((\hat{\boldsymbol{\psi}}\cdot\nabla,\hat{\boldsymbol{\theta}}^\ast)\cdot\nabla\) instead of \((R,Z)\) coordinates.

◆ gstep

template<class Device >
int Simulation< Device >::gstep

Current time step.

◆ heuristic_priv_pot

template<class Device >
bool Simulation< Device >::heuristic_priv_pot

Override the Poisson solver in the private region and replace the solution by !! the separatrix potential scaled with the initial electron temperature

◆ ignore_f0g

template<class Device >
bool Simulation< Device >::ignore_f0g

Ignore f0g in f0 calculation. When true, it skips the update of den_f0_h, instead using the initialized value (which is zero). This is helpful to density conservation of heating when f0g is noisy due to insufficient number of particles.

◆ is_XGCa

template<class Device >
constexpr bool Simulation< Device >::is_XGCa = true
staticconstexpr

Equivalent to the preprocessor flag for now.

◆ iter_solver

template<class Device >
bool Simulation< Device >::iter_solver

Split Poisson equation into axisymmetric and !! non-axisymmetric parts. The axisymmetric part has the !! flux-surface averaged potential on the RHS. Solve !! iteratively by lagging the RHS potential by one !! iteration. Stop after sml_iter_solver_niter iterations. !! If .false., split Poisson equation into !! flux-surface averaged part and rest (assuming that the !! polarization term commutes with the flux-surface !! average) and the rest; solve in two steps.

◆ iter_solver_converge

template<class Device >
bool Simulation< Device >::iter_solver_converge

If .true., keep taking iterations in PETSc until termination criteria is met, !! i.e. until a residual tolerance is reached or it is determined that the iterations !! have diverged or exceeded sml_iter_solver_max_it. If .false., take a fixed number !! of iterations set from sml_iter_solver_niter.

◆ limit_marker_den

template<class Device >
bool Simulation< Device >::limit_marker_den

Whether to limit marker density.

◆ loop_voltage_on

template<class Device >
bool Simulation< Device >::loop_voltage_on

Inductive current drive: loop voltage (from Faraday's law curl(E)=-dB/dt)

◆ loop_voltage_psimax

template<class Device >
double Simulation< Device >::loop_voltage_psimax

Inductive current drive: outer boundary (in pol. flux) of the loop voltage.

◆ mr_factor

template<class Device >
int Simulation< Device >::mr_factor

multirate timestepping for core region ions; mr_ratio (> 1) is the factor used to accelerate the core region physics.

◆ mr_psi_max

template<class Device >
double Simulation< Device >::mr_psi_max

maximum normalized psi value of each multirate region

◆ mstep

template<class Device >
int Simulation< Device >::mstep

Max number of time steps.

◆ multirate_timestepping

template<class Device >
bool Simulation< Device >::multirate_timestepping

Use multirate timestepping.

◆ neutrals

template<class Device >
bool Simulation< Device >::neutrals

Toggle for using neutrals.

◆ nhybrid

template<class Device >
constexpr int Simulation< Device >::nhybrid = 2
staticconstexpr

Number of iterations in electrostatic electron weight evolution scheme.

◆ no_fp_in_f

template<class Device >
bool Simulation< Device >::no_fp_in_f

If .true. the distribution function used for the source routines will not include particle information (only for testing)

◆ no_turb

template<class Device >
bool Simulation< Device >::no_turb

Set all non-axisymmetric field perturbations to zero (electromagnetic version only)

◆ nrk

template<class Device >
int Simulation< Device >::nrk

Order of Runge-Kutta time integration of particles + fields.

◆ nthreads

template<class Device >
int Simulation< Device >::nthreads

Number of OMP threads on the host.

◆ old_f0_update_analytic

template<class Device >
constexpr bool Simulation< Device >::old_f0_update_analytic = false
staticconstexpr

Equivalent to the preprocessor flag for now.

◆ poisson_adia_wall

template<class Device >
bool Simulation< Device >::poisson_adia_wall

EXPERIMENTAL! DON'T USE UNLESS ADVISED BY AN EXPERT.

◆ poisson_bias

template<class Device >
bool Simulation< Device >::poisson_bias

Whether to use a (0,0) bias potential on top of phi_00.

◆ poisson_use_bc

template<class Device >
bool Simulation< Device >::poisson_use_bc

(XGCa only) False: 00-bc is phi=0 everywhere; true: 00-boundary

◆ positive_phi00_sol

template<class Device >
bool Simulation< Device >::positive_phi00_sol

EXPERIMENTAL! DON'T USE UNLESS ADVISED BY AN EXPERT !! If true, makes sure that \(\langle\phi\rangle > 0\) in SOL and !! private flux region by adding a constant to \(\phi\)

◆ ptb_3db_on

template<class Device >
bool Simulation< Device >::ptb_3db_on

Toggle for ptb_3db.

◆ reduced_deltaf

template<class Device >
bool Simulation< Device >::reduced_deltaf

True if any species is reduced_deltaf. Will be removed.

◆ resamp_for_final_restart_write

template<class Device >
bool Simulation< Device >::resamp_for_final_restart_write

Perform resampling before dumping the final restart file.

◆ resamp_on

template<class Device >
bool Simulation< Device >::resamp_on

Do resampling.

◆ resamp_restart_read

template<class Device >
bool Simulation< Device >::resamp_restart_read

Whether to read a restart file written from a simulation with different grid.

◆ separate_n0

template<class Device >
bool Simulation< Device >::separate_n0

Use separate n0 solver even in reduced_deltaf.

◆ sheath_adjust

template<class Device >
bool Simulation< Device >::sheath_adjust

Whether to adjust the sheath.

◆ sheath_global_balance

template<class Device >
bool Simulation< Device >::sheath_global_balance

Whether to use the global loss balance functionality in sheath_mode=2.

◆ sheath_mode

template<class Device >
int Simulation< Device >::sheath_mode

Type of sheath (0 is none)

◆ symmetric_f

template<class Device >
bool Simulation< Device >::symmetric_f

Enforce axisymmetry of the total distribution function on the grid.

◆ thermal_bath_on

template<class Device >
bool Simulation< Device >::thermal_bath_on

Switch on thermal-bath and coarse-graining operation, see ccm_param parameters.

◆ time

template<class Device >
double Simulation< Device >::time

Current simulation time.

◆ transit_time

template<class Device >
double Simulation< Device >::transit_time

Torodial transit time of an ion with the characteristic energy.

◆ update_g_alpha

template<class Device >
double Simulation< Device >::update_g_alpha

Fraction of the numerical marker particle density that is transferred to the g2

◆ update_poisson_solver

template<class Device >
bool Simulation< Device >::update_poisson_solver

Whether poisson solver is updated.

◆ use_em_bounds

template<class Device >
bool Simulation< Device >::use_em_bounds

◆ use_unfused_electron_push_kernel

template<class Device >
bool Simulation< Device >::use_unfused_electron_push_kernel

◆ zero_inner_bd_turb

template<class Device >
bool Simulation< Device >::zero_inner_bd_turb

The documentation for this class was generated from the following files: