7 #ifndef ABLASTR_NODALFIELDGATHER_H_
8 #define ABLASTR_NODALFIELDGATHER_H_
35 const amrex::ParticleReal yp,
36 const amrex::ParticleReal zp,
39 int&
i,
int& j,
int& k, amrex::Real W[AMREX_SPACEDIM][2],
int nodal=1) noexcept
41 using namespace amrex::literals;
43 #if !((nodal==0)||(nodal==1))
47 #if (defined WARPX_DIM_3D)
48 const amrex::Real
x = (xp - plo[0]) * dxi[0] +
static_cast<amrex::Real
>(nodal-1)*0.5_rt;
49 const amrex::Real y = (yp - plo[1]) * dxi[1] +
static_cast<amrex::Real
>(nodal-1)*0.5_rt;
50 const amrex::Real z = (zp - plo[2]) * dxi[2] +
static_cast<amrex::Real
>(nodal-1)*0.5_rt;
52 i =
static_cast<int>(amrex::Math::floor(
x));
53 j =
static_cast<int>(amrex::Math::floor(y));
54 k =
static_cast<int>(amrex::Math::floor(z));
60 W[0][0] = 1.0_rt - W[0][1];
61 W[1][0] = 1.0_rt - W[1][1];
62 W[2][0] = 1.0_rt - W[2][1];
64 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
66 # if (defined WARPX_DIM_XZ)
67 const amrex::Real
x = (xp - plo[0]) * dxi[0] +
static_cast<amrex::Real
>(nodal-1)*0.5_rt;
69 i =
static_cast<int>(amrex::Math::floor(
x));
71 # elif (defined WARPX_DIM_RZ)
72 const amrex::Real
r = (std::sqrt(xp*xp+yp*yp) - plo[0]) * dxi[0] +
static_cast<amrex::Real
>(nodal-1)*0.5_rt;
73 i =
static_cast<int>(amrex::Math::floor(
r));
77 const amrex::Real z = (zp - plo[1]) * dxi[1] +
static_cast<amrex::Real
>(nodal-1)*0.5_rt;
78 j =
static_cast<int>(amrex::Math::floor(z));
81 W[0][0] = 1.0_rt - W[0][1];
82 W[1][0] = 1.0_rt - W[1][1];
100 const amrex::Real W[AMREX_SPACEDIM][2],
103 #if (defined WARPX_DIM_3D)
104 amrex::Real
value = 0;
105 value += scalar_field(
i, j , k ) * W[0][0] * W[1][0] * W[2][0];
106 value += scalar_field(
i+1, j , k ) * W[0][1] * W[1][0] * W[2][0];
107 value += scalar_field(
i, j+1, k ) * W[0][0] * W[1][1] * W[2][0];
108 value += scalar_field(
i+1, j+1, k ) * W[0][1] * W[1][1] * W[2][0];
109 value += scalar_field(
i, j , k+1) * W[0][0] * W[1][0] * W[2][1];
110 value += scalar_field(
i+1, j , k+1) * W[0][1] * W[1][0] * W[2][1];
111 value += scalar_field(
i , j+1, k+1) * W[0][0] * W[1][1] * W[2][1];
112 value += scalar_field(
i+1, j+1, k+1) * W[0][1] * W[1][1] * W[2][1];
113 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
114 amrex::Real
value = 0;
115 value += scalar_field(
i, j , k) * W[0][0] * W[1][0];
116 value += scalar_field(
i+1, j , k) * W[0][1] * W[1][0];
117 value += scalar_field(
i, j+1, k) * W[0][0] * W[1][1];
118 value += scalar_field(
i+1, j+1, k) * W[0][1] * W[1][1];
120 const amrex::Real
value = 0;
138 const amrex::ParticleReal yp,
139 const amrex::ParticleReal zp,
146 amrex::Real W[AMREX_SPACEDIM][2];
164 const amrex::ParticleReal yp,
165 const amrex::ParticleReal zp,
174 amrex::Real W[AMREX_SPACEDIM][2];
#define AMREX_GPU_HOST_DEVICE
#define ABLASTR_ABORT_WITH_MESSAGE(MSG)
Definition: TextMsg.H:77
Definition: DepositCharge.H:21
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.
Definition: NodalFieldGather.H:99
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...
Definition: NodalFieldGather.H:34
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:/...
Definition: NodalFieldGather.H:137
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:/...
Definition: NodalFieldGather.H:163
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
i
Definition: check_interp_points_and_weights.py:174
ii
Definition: check_interp_points_and_weights.py:148
jj
Definition: check_interp_points_and_weights.py:160
value
Definition: updateAMReX.py:141