XGC1
|
#include <perturbed_B_field.hpp>
Public Member Functions | |
PerturbedBField (NLReader::NamelistReader &nlr, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const Smoothing &smoothing, bool explicit_em) | |
PerturbedBField () | |
void | initialize (NLReader::NamelistReader &nlr, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const Smoothing &smoothing, bool explicit_em) |
void | read_3db_xgc_mesh (const std::string &input_dir, const std::string &filename, const Grid< DeviceType > &grid) |
void | copy_bfield_from_fortran () const |
void | copy_to_device_if_on (int gstep) |
void | save_Ah (const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &Ah) const |
void | replace_B_grad_phi (const View< double **, CLayout, HostType > &view, int gstep, double dt) const |
void | ramp_up_As (int ipc, int gstep, double dt, double dt_full, const View< double *, CLayout, HostType > &As_local) |
void | update_rampup_and_write (bool explicit_electromagnetic, int gstep, double time) |
void | set_to_As_vac (const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &As) const |
void | subtract_As_vac_fac (const View< double *, CLayout, HostType > &rhs2) const |
void | save_rhs_poisson (int iflag, int ipc, bool em_es_step, const View< double *, CLayout, HostType > &rhs) const |
void | save_rhs_ampere (const Simulation< DeviceType > &sml, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, Plasma &plasma, const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &dpot_h, const View< double *, CLayout, HostType > &rhs, int ipc) const |
void | write_checkpoint_files (const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const XGC_IO_Stream &stream) const |
void | read_checkpoint_files (const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const XGC_IO_Stream &stream, const XGC_IO_Stream &em_stream) |
KOKKOS_INLINE_FUNCTION bool | is_triggered (int gstep) const |
KOKKOS_INLINE_FUNCTION bool | rampup_is_triggered (int gstep) const |
template<Order OT> | |
KOKKOS_INLINE_FUNCTION void | get_delta_b (const Grid< Device > &grid, const Simd< double > &phi, const SimdGridWeights< OT, PhiInterpType::None > &grid_wts, SimdVector &tdb) const |
void | read_3db_xgc_mesh (const std::string &input_dir, const std::string &filename, const Grid< DeviceType > &grid) |
void | initialize (NLReader::NamelistReader &nlr, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const Smoothing &smoothing, bool explicit_em) |
void | save_Ah (const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &Ah) const |
void | replace_B_grad_phi (const View< double **, CLayout, HostType > &view, int gstep, double dt) const |
void | update_rampup_and_write (bool explicit_electromagnetic, int gstep, double time) |
void | ramp_up_As (int ipc, int gstep, double dt, double dt_full, const View< double *, CLayout, HostType > &As_local) |
void | set_to_As_vac (const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &As) const |
void | subtract_As_vac_fac (const View< double *, CLayout, HostType > &rhs2) const |
void | save_rhs_poisson (int iflag, int ipc, bool em_es_step, const View< double *, CLayout, HostType > &rhs) const |
void | save_rhs_ampere (const Simulation< DeviceType > &sml, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, Plasma &plasma, const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &dpot_h, const View< double *, CLayout, HostType > &rhs, int ipc) const |
void | write_checkpoint_files (const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const XGC_IO_Stream &stream) const |
void | read_checkpoint_files (const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const XGC_IO_Stream &stream, const XGC_IO_Stream &em_stream) |
PerturbedBField (NLReader::NamelistReader &nlr, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const Smoothing &smoothing, bool explicit_em) | |
Public Attributes | |
bool | rampup_vac |
(.true.) Ramp up perturbed field slowly, (.false.) turn on perturbed field abruptly More... | |
int | num_ntor |
Number of toroidal mode numbers (<= sml_nphi_total) More... | |
int | rampup_time |
Number of time steps over which the perturbed field is ramped up. More... | |
int | start_time |
Time step in which perturbed field is switched on. More... | |
int | mode |
Mode of RMP operation (details in XGC_core/module_ptb_3db.F90 --> ptb_3db_mode) More... | |
bool | full_spec_on |
whether to retain the full RMP spectrum from M3D-C1 or only |m/q-n|<=sml_mode_select_mres_q ; Always on for ES More... | |
double | rampup_fac |
The current relative amplitude of the vacuum vector potential (if ptb_3db_mode==2) - this is updated during the simulation. More... | |
double | rampup_fac0 |
The previous (backup) relative amplitude of the vacuum vector potential (if ptb_3db_mode==2) More... | |
int | dasdt_opt |
Selector for how to compute dAs/dt for RMP penetration. More... | |
int | mstep_es |
int | mstep_em |
int | es_to_em_dt_ratio |
Ratio of EM to ES time step size in RMP penetration calculation. More... | |
bool | active |
Whether the field is in device memory and should be used in the push etc. Could be renamed for clarity. More... | |
View< int *, CLayout, Device > | ntor |
Array to store the toroidal mode numbers. More... | |
View< Field< VarType::Vector, PhiInterpType::None > **, CLayout, Device > | bfield_im_vac |
perturbed vacuum field on XGC mesh, imaginary part More... | |
View< Field< VarType::Vector, PhiInterpType::None > **, CLayout, Device > | bfield_re_vac |
perturbed vacuum field on XGC mesh, real part, dimensions: (grid vertex, R-Z-phi components, tor. mode number) More... | |
Private Attributes | |
View< double **, CLayout, HostType > | jpar0 |
View< double ***, CLayout, HostType > | cden0 |
View< double **, CLayout, HostType > | As_vac |
View< double ***, CLayout, HostType > | Ah_rmp |
PerturbedBField< Device >::PerturbedBField | ( | NLReader::NamelistReader & | nlr, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const MagneticField< DeviceType > & | magnetic_field, | ||
const Grid< DeviceType > & | grid, | ||
const Smoothing & | smoothing, | ||
bool | explicit_em | ||
) |
|
inline |
PerturbedBField< DeviceType >::PerturbedBField | ( | NLReader::NamelistReader & | nlr, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const MagneticField< DeviceType > & | magnetic_field, | ||
const Grid< DeviceType > & | grid, | ||
const Smoothing & | smoothing, | ||
bool | explicit_em | ||
) |
Constructor for perturbed magnetic field class
< Selector for how to compute dAs/dt for RMP penetration
|
inline |
|
inline |
KOKKOS_INLINE_FUNCTION void PerturbedBField< Device >::get_delta_b | ( | const Grid< Device > & | grid, |
const Simd< double > & | phi, | ||
const SimdGridWeights< OT, PhiInterpType::None > & | grid_wts, | ||
SimdVector & | tdb | ||
) | const |
Get magnetic field perturbation vector for a vector of particle locations
[in] | grid | The grid object |
[in] | phi | Phi coordinate |
[in] | grid_wts | Triangle the particles are in |
[out] | tdb | Magnetic field perturbation vector |
void PerturbedBField< DeviceType >::initialize | ( | NLReader::NamelistReader & | nlr, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const MagneticField< DeviceType > & | magnetic_field, | ||
const Grid< DeviceType > & | grid, | ||
const Smoothing & | smoothing, | ||
bool | explicit_em | ||
) |
< Inner boundary of RMP representation in A_s
< Outer boundary of RMP representation in A_s
< Blending width between A_s and full A representation
void PerturbedBField< Device >::initialize | ( | NLReader::NamelistReader & | nlr, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const MagneticField< DeviceType > & | magnetic_field, | ||
const Grid< DeviceType > & | grid, | ||
const Smoothing & | smoothing, | ||
bool | explicit_em | ||
) |
KOKKOS_INLINE_FUNCTION bool PerturbedBField< Device >::is_triggered | ( | int | gstep | ) | const |
void PerturbedBField< DeviceType >::ramp_up_As | ( | int | ipc, |
int | gstep, | ||
double | dt, | ||
double | dt_full, | ||
const View< double *, CLayout, HostType > & | As_local | ||
) |
void PerturbedBField< Device >::ramp_up_As | ( | int | ipc, |
int | gstep, | ||
double | dt, | ||
double | dt_full, | ||
const View< double *, CLayout, HostType > & | As_local | ||
) |
KOKKOS_INLINE_FUNCTION bool PerturbedBField< Device >::rampup_is_triggered | ( | int | gstep | ) | const |
void PerturbedBField< DeviceType >::read_3db_xgc_mesh | ( | const std::string & | input_dir, |
const std::string & | filename, | ||
const Grid< DeviceType > & | grid | ||
) |
void PerturbedBField< Device >::read_3db_xgc_mesh | ( | const std::string & | input_dir, |
const std::string & | filename, | ||
const Grid< DeviceType > & | grid | ||
) |
void PerturbedBField< DeviceType >::read_checkpoint_files | ( | const Grid< DeviceType > & | grid, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const XGC_IO_Stream & | stream, | ||
const XGC_IO_Stream & | em_stream | ||
) |
void PerturbedBField< Device >::read_checkpoint_files | ( | const Grid< DeviceType > & | grid, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const XGC_IO_Stream & | stream, | ||
const XGC_IO_Stream & | em_stream | ||
) |
void PerturbedBField< DeviceType >::replace_B_grad_phi | ( | const View< double **, CLayout, HostType > & | view, |
int | gstep, | ||
double | dt | ||
) | const |
void PerturbedBField< Device >::replace_B_grad_phi | ( | const View< double **, CLayout, HostType > & | view, |
int | gstep, | ||
double | dt | ||
) | const |
void PerturbedBField< DeviceType >::save_Ah | ( | const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > & | Ah | ) | const |
void PerturbedBField< Device >::save_Ah | ( | const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > & | Ah | ) | const |
void PerturbedBField< DeviceType >::save_rhs_ampere | ( | const Simulation< DeviceType > & | sml, |
const MagneticField< DeviceType > & | magnetic_field, | ||
const Grid< DeviceType > & | grid, | ||
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
Plasma & | plasma, | ||
const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > & | dpot_h, | ||
const View< double *, CLayout, HostType > & | rhs, | ||
int | ipc | ||
) | const |
void PerturbedBField< Device >::save_rhs_ampere | ( | const Simulation< DeviceType > & | sml, |
const MagneticField< DeviceType > & | magnetic_field, | ||
const Grid< DeviceType > & | grid, | ||
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
Plasma & | plasma, | ||
const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > & | dpot_h, | ||
const View< double *, CLayout, HostType > & | rhs, | ||
int | ipc | ||
) | const |
void PerturbedBField< DeviceType >::save_rhs_poisson | ( | int | iflag, |
int | ipc, | ||
bool | em_es_step, | ||
const View< double *, CLayout, HostType > & | rhs | ||
) | const |
void PerturbedBField< Device >::save_rhs_poisson | ( | int | iflag, |
int | ipc, | ||
bool | em_es_step, | ||
const View< double *, CLayout, HostType > & | rhs | ||
) | const |
void PerturbedBField< DeviceType >::set_to_As_vac | ( | const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > & | As | ) | const |
void PerturbedBField< Device >::set_to_As_vac | ( | const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > & | As | ) | const |
void PerturbedBField< DeviceType >::subtract_As_vac_fac | ( | const View< double *, CLayout, HostType > & | rhs2 | ) | const |
void PerturbedBField< Device >::subtract_As_vac_fac | ( | const View< double *, CLayout, HostType > & | rhs2 | ) | const |
void PerturbedBField< DeviceType >::update_rampup_and_write | ( | bool | explicit_electromagnetic, |
int | gstep, | ||
double | time | ||
) |
void PerturbedBField< Device >::update_rampup_and_write | ( | bool | explicit_electromagnetic, |
int | gstep, | ||
double | time | ||
) |
void PerturbedBField< DeviceType >::write_checkpoint_files | ( | const Grid< DeviceType > & | grid, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const XGC_IO_Stream & | stream | ||
) | const |
void PerturbedBField< Device >::write_checkpoint_files | ( | const Grid< DeviceType > & | grid, |
const DomainDecomposition< DeviceType > & | pol_decomp, | ||
const XGC_IO_Stream & | stream | ||
) | const |
bool PerturbedBField< Device >::active |
Whether the field is in device memory and should be used in the push etc. Could be renamed for clarity.
|
private |
|
private |
View<Field<VarType::Vector,PhiInterpType::None>**,CLayout,Device> PerturbedBField< Device >::bfield_im_vac |
perturbed vacuum field on XGC mesh, imaginary part
View<Field<VarType::Vector,PhiInterpType::None>**,CLayout,Device> PerturbedBField< Device >::bfield_re_vac |
perturbed vacuum field on XGC mesh, real part, dimensions: (grid vertex, R-Z-phi components, tor. mode number)
|
private |
int PerturbedBField< Device >::dasdt_opt |
Selector for how to compute dAs/dt for RMP penetration.
int PerturbedBField< Device >::es_to_em_dt_ratio |
Ratio of EM to ES time step size in RMP penetration calculation.
bool PerturbedBField< Device >::full_spec_on |
whether to retain the full RMP spectrum from M3D-C1 or only |m/q-n|<=sml_mode_select_mres_q ; Always on for ES
|
private |
int PerturbedBField< Device >::mode |
Mode of RMP operation (details in XGC_core/module_ptb_3db.F90 --> ptb_3db_mode)
int PerturbedBField< Device >::mstep_em |
Number of electromagnetic time steps for RMP penetration calculation with damped Newton iteration
int PerturbedBField< Device >::mstep_es |
Number of electrostatic time steps for RMP penetration calculation with damped Newton iteration
View<int*,CLayout,Device> PerturbedBField< Device >::ntor |
Array to store the toroidal mode numbers.
int PerturbedBField< Device >::num_ntor |
Number of toroidal mode numbers (<= sml_nphi_total)
double PerturbedBField< Device >::rampup_fac |
The current relative amplitude of the vacuum vector potential (if ptb_3db_mode==2) - this is updated during the simulation.
double PerturbedBField< Device >::rampup_fac0 |
The previous (backup) relative amplitude of the vacuum vector potential (if ptb_3db_mode==2)
int PerturbedBField< Device >::rampup_time |
Number of time steps over which the perturbed field is ramped up.
bool PerturbedBField< Device >::rampup_vac |
(.true.) Ramp up perturbed field slowly, (.false.) turn on perturbed field abruptly
int PerturbedBField< Device >::start_time |
Time step in which perturbed field is switched on.