XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Functions
update_f0.cpp File Reference
#include "timer_macro.hpp"
#include "globals.hpp"
#include "sml.hpp"
#include "magnetic_field.hpp"
#include "grid.hpp"
#include "particles.hpp"
#include "plasma.hpp"
#include "domain_decomposition.hpp"
#include "quadratic_weights.hpp"
#include <math.h>
#include "send_particles.hpp"
#include "streamed_parallel_for.hpp"
#include "update_f0.hpp"
Include dependency graph for update_f0.cpp:

Functions

double * get_inv_grid_vol_loc ()
 
double * get_f0g_loc ()
 
double * get_f_loc ()
 
double * get_n_loc ()
 
double * get_temp_ev_loc ()
 
template<class Device >
KOKKOS_INLINE_FUNCTION void f0_update_f0g (const Grid< Device > &grid, const Species< Device > &species, const VelocityGrid &vgrid, int node, double phi, double mu_n, double vp_n, double df0g, double df2, bool update_f_and_n)
 
template<class Device >
KOKKOS_INLINE_FUNCTION void f0_update_f0g_quadratic (const Grid< Device > &grid, const Species< Device > &species, const VelocityGrid &vgrid, int node, double phi, double mu_n, double vp_n, double df0g, double df2, bool update_f_and_n)
 
template<class Device >
KOKKOS_INLINE_FUNCTION void f0_update_f0_n (const Grid< Device > &grid, const Species< Device > &species, const VelocityGrid &vgrid, int node, double mu_n, double vp_n)
 
template<class Device >
KOKKOS_INLINE_FUNCTION void update_f0_sp_c (const TmpSpecies< Device > &tmp_species, const Simulation< Device > &sml, const Grid< Device > &grid, const MagneticField< Device > &magnetic_field, const Species< Device > &species, const VelocityGrid &vgrid, const DomainDecomposition< Device > &pol_decomp, double alpha_in, bool update_f_and_n, bool f_source_on, int i_item)
 
template<class Device >
KOKKOS_INLINE_FUNCTION void update_f0_sp_c_pseudo_inv (const TmpSpecies< Device > &tmp_species, const Simulation< Device > &sml, const Grid< Device > &grid, const MagneticField< Device > &magnetic_field, const Species< Device > &species, const VelocityGrid &vgrid, const DomainDecomposition< Device > &pol_decomp, int isp_non_ad, double alpha_in, bool update_f_and_n, bool f_source_on, const Pseudo_inverse< Device > &pseudo_inv, int i_item)
 
void f0_update (const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, const Species< DeviceType > &species, const TmpSpecies< DeviceType > &tmp_species, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, int isp_non_ad, double alpha_in, bool f_source_on, Pseudo_inverse< DeviceType > &pseudo_inv)
 
void update_f0_sp (const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, Species< DeviceType > &species, const TmpSpecies< DeviceType > &tmp_species, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, int isp_non_ad, double alpha_in, bool f_source_on, DMWrapper &pseudo_inv_dm, Pseudo_inverse< DeviceType > &pseudo_inv)
 
void all_species_update_f0 (const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, Plasma &plasma, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, double alpha_in, bool f_source_on, DMWrapper &pseudo_inv_dm, Pseudo_inverse< DeviceType > &pseudo_inv)
 
void f0_update_f0g_pseudo_inv (const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, const Species< DeviceType > &species, const VelocityGrid &vgrid, const DomainDecomposition< DeviceType > &pol_decomp, int isp_non_ad, bool update_f_and_n, DMWrapper &pseudo_inv_dm, Pseudo_inverse< DeviceType > &pseudo_inv)
 

Function Documentation

void all_species_update_f0 ( const Simulation< DeviceType > &  sml,
const Grid< DeviceType > &  grid,
const MagneticField< DeviceType > &  magnetic_field,
Plasma plasma,
const VelocityGrid vgrid,
const DomainDecomposition< DeviceType > &  pol_decomp,
double  alpha_in,
bool  f_source_on,
DMWrapper pseudo_inv_dm,
Pseudo_inverse< DeviceType > &  pseudo_inv 
)

Here is the call graph for this function:

Here is the caller graph for this function:

void f0_update ( const Simulation< DeviceType > &  sml,
const Grid< DeviceType > &  grid,
const MagneticField< DeviceType > &  magnetic_field,
const Species< DeviceType > &  species,
const TmpSpecies< DeviceType > &  tmp_species,
const VelocityGrid vgrid,
const DomainDecomposition< DeviceType > &  pol_decomp,
int  isp_non_ad,
double  alpha_in,
bool  f_source_on,
Pseudo_inverse< DeviceType > &  pseudo_inv 
)

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device >
KOKKOS_INLINE_FUNCTION void f0_update_f0_n ( const Grid< Device > &  grid,
const Species< Device > &  species,
const VelocityGrid vgrid,
int  node,
double  mu_n,
double  vp_n 
)

Update velocity interpolation normalization factors bilinearly or quadratically depending on element order. Uses quadratic Lagrange polynomials if quadratic interpolation.

Parameters
[in]gridspatial grid object
[in]speciescontains species parameters
[in]vgridcontains the velocity grid dimensions
[in]nodegrid node of particle
[in]mu_nmu [(v_{perpendicular}/v_{thermal})^2] of particle
[in]vp_nparallel velocity of particle
Returns
void

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device >
KOKKOS_INLINE_FUNCTION void f0_update_f0g ( const Grid< Device > &  grid,
const Species< Device > &  species,
const VelocityGrid vgrid,
int  node,
double  phi,
double  mu_n,
double  vp_n,
double  df0g,
double  df2,
bool  update_f_and_n 
)

Interpolate weight of one particle to grid distribution bilinearly in velocity space.

Parameters
[in]gridspatial grid object
[in]speciescontains species parameters
[in]vgridcontains the velocity grid dimensions
[in]nodegrid node of particle
[in]phiphi of particle
[in]mu_nmu [(v_{perpendicular}/v_{thermal})^2] of particle
[in]vp_nparallel velocity of particle
[in]df0gparticle weight multiplied by fraction transferred to distribution function every source time step: alpha*w1*w0
[in]df2particle weight: w1*w0
[in]update_f_and_nwhether to add particle weights in the distribution function used for the source routines
Returns
void

Here is the caller graph for this function:

void f0_update_f0g_pseudo_inv ( const Grid< DeviceType > &  grid,
const MagneticField< DeviceType > &  magnetic_field,
const Species< DeviceType > &  species,
const VelocityGrid vgrid,
const DomainDecomposition< DeviceType > &  pol_decomp,
int  isp_non_ad,
bool  update_f_and_n,
DMWrapper pseudo_inv_dm,
Pseudo_inverse< DeviceType > &  pseudo_inv 
)

Loop over nodes and interpolate weights to grid distribution. This is pseudo-inverse version with particles -> velocity grid interpolation done with PETSc.

NOTE: This method is slower than the other but may be necessary for a PETSc collision operator.

Parameters
[in]gridspatial grid object
[in]magnetic_fieldmagnetic field object
[in]speciescontains species parameters
[in]vgridcontains the velocity grid dimensions
[in]pol_decompcontains poloidal decomposition info
[in]isp_non_adnon-adiabatic species index
[in]update_f_and_nwhether to add particle weights in the distribution function used for the source routines
[in]pseudo_inv_dmpseudo-inverse mesh object
[in]pseudo_invpseudo-inverse object (contains pseudo-inverse arrays)
Returns
void

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device >
KOKKOS_INLINE_FUNCTION void f0_update_f0g_quadratic ( const Grid< Device > &  grid,
const Species< Device > &  species,
const VelocityGrid vgrid,
int  node,
double  phi,
double  mu_n,
double  vp_n,
double  df0g,
double  df2,
bool  update_f_and_n 
)

Interpolate weight of one particle to grid distribution quadratically in velocity space. Uses quadratic Lagrange polynomials.

Parameters
[in]gridspatial grid object
[in]speciescontains species parameters
[in]vgridcontains the velocity grid dimensions
[in]nodegrid node of particle
[in]phiphi of particle
[in]mu_nmu [(v_{perpendicular}/v_{thermal})^2] of particle
[in]vp_nparallel velocity of particle
[in]df0gparticle weight multiplied by fraction transferred to distribution function every source time step: alpha*w1*w0
[in]df2particle weight: w1*w0
[in]update_f_and_nwhether to add particle weights in the distribution function used for the source routines
Returns
void

Here is the call graph for this function:

Here is the caller graph for this function:

double* get_f0g_loc ( )
double* get_f_loc ( )
double* get_inv_grid_vol_loc ( )

Here is the caller graph for this function:

double* get_n_loc ( )

Here is the caller graph for this function:

double* get_temp_ev_loc ( )
void update_f0_sp ( const Simulation< DeviceType > &  sml,
const Grid< DeviceType > &  grid,
const MagneticField< DeviceType > &  magnetic_field,
Species< DeviceType > &  species,
const TmpSpecies< DeviceType > &  tmp_species,
const VelocityGrid vgrid,
const DomainDecomposition< DeviceType > &  pol_decomp,
int  isp_non_ad,
double  alpha_in,
bool  f_source_on,
DMWrapper pseudo_inv_dm,
Pseudo_inverse< DeviceType > &  pseudo_inv 
)

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device >
KOKKOS_INLINE_FUNCTION void update_f0_sp_c ( const TmpSpecies< Device > &  tmp_species,
const Simulation< Device > &  sml,
const Grid< Device > &  grid,
const MagneticField< Device > &  magnetic_field,
const Species< Device > &  species,
const VelocityGrid vgrid,
const DomainDecomposition< Device > &  pol_decomp,
double  alpha_in,
bool  update_f_and_n,
bool  f_source_on,
int  i_item 
)

Loop over particles and interpolate weights to grid distribution.

Parameters
[in]tmp_speciescontains the particle data
[in]smlcontains simulation control parameters
[in]gridspatial grid object
[in]magnetic_fieldmagnetic field object
[in]speciescontains species parameters
[in]vgridcontains the velocity grid dimensions
[in]pol_decompcontains poloidal decomposition info
[in]alpha_infraction of particle weight transferred to distribution function every source time step
[in]update_f_and_nwhether to add particle weights in the distribution function used for the source routines
[in]f_source_onis this function ever called with this set to false?
[in]i_itemthe particle or particle vector index
Returns
void

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device >
KOKKOS_INLINE_FUNCTION void update_f0_sp_c_pseudo_inv ( const TmpSpecies< Device > &  tmp_species,
const Simulation< Device > &  sml,
const Grid< Device > &  grid,
const MagneticField< Device > &  magnetic_field,
const Species< Device > &  species,
const VelocityGrid vgrid,
const DomainDecomposition< Device > &  pol_decomp,
int  isp_non_ad,
double  alpha_in,
bool  update_f_and_n,
bool  f_source_on,
const Pseudo_inverse< Device > &  pseudo_inv,
int  i_item 
)

Loop over particles and interpolate weights to grid distribution. Pseudo-inverse version.

Parameters
[in]tmp_speciescontains the particle data
[in]smlcontains simulation control parameters
[in]gridspatial grid object
[in]magnetic_fieldmagnetic field object
[in]speciescontains species parameters
[in]vgridcontains the velocity grid dimensions
[in]pol_decompcontains poloidal decomposition info
[in]isp_non_adnon-adiabatic species index
[in]alpha_infraction of particle weight transferred to distribution function every source time step
[in]update_f_and_nwhether to add particle weights in the distribution function used for the source routines
[in]f_source_onis this function ever called with this set to false?
[in]pseudo_invpseudo-inverse object (contains pseudo-inverse arrays)
[in]i_itemthe particle or particle vector index
Returns
void

Here is the call graph for this function:

Here is the caller graph for this function: