WarpX
|
Functions | |
template<typename T_PC > | |
static void | deposit_charge (typename T_PC::ParIterType &pti, typename T_PC::RealVector const &wp, amrex::Real const charge, int const *const ion_lev, amrex::MultiFab *rho, amrex::FArrayBox &local_rho, int const particle_shape, std::array< amrex::Real, 3 > const &dx, std::array< amrex::Real, 3 > const &xyzmin, int const n_rz_azimuthal_modes=0, std::optional< amrex::IntVect > num_rho_deposition_guards=std::nullopt, std::optional< int > depos_lev=std::nullopt, std::optional< amrex::IntVect > rel_ref_ratio=std::nullopt, long const offset=0, std::optional< long > np_to_deposit=std::nullopt, int const icomp=0, int const nc=1, bool const do_device_synchronize=true) |
AMREX_GPU_HOST_DEVICE AMREX_INLINE void | compute_weights (const amrex::ParticleReal xp, const amrex::ParticleReal yp, const amrex::ParticleReal zp, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, int &i, int &j, int &k, amrex::Real W[AMREX_SPACEDIM][2], int nodal=1) noexcept |
Compute weight of each surrounding node (or cell-centered nodes) in interpolating a nodal ((or a cell-centered node) field to the given coordinates. If nodal=1, then the calculations will be done with respect to the nodes (default). If nodal=0, then the calculations will be done with respect to the cell-centered nodal) More... | |
AMREX_GPU_HOST_DEVICE AMREX_INLINE amrex::Real | interp_field_nodal (int i, int j, int k, const amrex::Real W[AMREX_SPACEDIM][2], amrex::Array4< const amrex::Real > const &scalar_field) noexcept |
Interpolate nodal field value based on surrounding indices and weights. More... | |
AMREX_GPU_HOST_DEVICE AMREX_INLINE amrex::Real | doGatherScalarFieldNodal (const amrex::ParticleReal xp, const amrex::ParticleReal yp, const amrex::ParticleReal zp, amrex::Array4< const amrex::Real > const &scalar_field, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &lo) noexcept |
Scalar field gather for a single particle. The field has to be defined at the cell nodes (see https://amrex-codes.github.io/amrex/docs_html/Basics.html#id2) More... | |
AMREX_GPU_HOST_DEVICE AMREX_INLINE amrex::GpuArray< amrex::Real, 3 > | doGatherVectorFieldNodal (const amrex::ParticleReal xp, const amrex::ParticleReal yp, const amrex::ParticleReal zp, amrex::Array4< const amrex::Real > const &vector_field_x, amrex::Array4< const amrex::Real > const &vector_field_y, amrex::Array4< const amrex::Real > const &vector_field_z, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &lo) noexcept |
Vector field gather for a single particle. The field has to be defined at the cell nodes (see https://amrex-codes.github.io/amrex/docs_html/Basics.html#id2) More... | |
template<typename T_PC > | |
static std::tuple< amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal > | MinAndMaxPositions (T_PC const &pc) |
template<typename T_PC , int T_RealSoAWeight> | |
static std::tuple< amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal, amrex::ParticleReal > | MeanAndStdPositions (T_PC const &pc) |
|
noexcept |
Compute weight of each surrounding node (or cell-centered nodes) in interpolating a nodal ((or a cell-centered node) field to the given coordinates. If nodal=1, then the calculations will be done with respect to the nodes (default). If nodal=0, then the calculations will be done with respect to the cell-centered nodal)
This currently only does linear order.
xp,yp,zp | Particle position coordinates |
plo | Index lower bounds of domain. |
dxi | inverse 3D cell spacing |
i,j,k | Variables to store indices of position on grid (nodal or cell-centered, depending of the value of nodal ) |
W | 2D array of weights to store each neighbouring node (or cell-centered node) |
nodal | Int that tells if the weights are calculated in respect to the nodes (nodal=1) of the cell-centered nodes (nodal=0) |
|
static |
Perform charge deposition for the particles on a tile.
T_PC | a type of amrex::ParticleContainer |
pti | an amrex::ParIter pointing to the tile to operate on |
wp | vector of the particle weights for those particles. |
charge | charge of the particle species |
ion_lev | pointer to array of particle ionization level. This is required to have the charge of each macroparticle since q is a scalar. For non-ionizable species, ion_lev is a null pointer. |
rho | MultiFab of the charge density |
local_rho | temporary FArrayBox for deposition with OpenMP |
particle_shape | shape factor in each direction |
dx | cell spacing at level lev |
xyzmin | lo corner of the current tile in physical coordinates. |
n_rz_azimuthal_modes | number of azimuthal modes in use, irrelevant outside RZ geometry (default: 0) |
num_rho_deposition_guards | number of ghost cells to use for rho (default: rho.nGrowVect()) |
depos_lev | the level to deposit the particles to (default: lev) |
rel_ref_ratio | mesh refinement ratio between lev and depos_lev (default: 1) |
offset | index to start at when looping over particles to deposit (default: 0) |
np_to_deposit | number of particles to deposit (default: pti.numParticles()) |
icomp | component in MultiFab to start depositing to |
nc | number of components to deposit |
do_device_synchronize | call amrex::Gpu::synchronize() for tiny profiler regions (default: true) |
|
noexcept |
Scalar field gather for a single particle. The field has to be defined at the cell nodes (see https://amrex-codes.github.io/amrex/docs_html/Basics.html#id2)
xp,yp,zp | Particle position coordinates |
scalar_field | Array4 of the nodal scalar field, either full array or tile. |
dxi | inverse 3D cell spacing |
lo | Index lower bounds of domain. |
|
noexcept |
Vector field gather for a single particle. The field has to be defined at the cell nodes (see https://amrex-codes.github.io/amrex/docs_html/Basics.html#id2)
xp,yp,zp | Particle position coordinates |
vector_field_x,vector_field_y,vector_field_z | Array4 of nodal scalar fields, either full array or tile. |
dxi | inverse 3D cell spacing |
lo | Index lower bounds of domain. |
|
noexcept |
Interpolate nodal field value based on surrounding indices and weights.
i,j,k | Indices of position on grid |
W | 2D array of weights for each neighbouring node |
scalar_field | Array4 of the nodal scalar field, either full array or tile. |
|
static |
Compute the mean and std of the particle position in each dimension
T_PC | a type of amrex::ParticleContainer |
T_RealSoAWeight |
pc | the particle container to operate on |
|
static |
Compute the min and max of the particle position in each dimension
T_PC | a type of amrex::ParticleContainer |
pc | the particle container to operate on |