8 #ifndef WARPX_CHARGEDEPOSITION_H_
9 #define WARPX_CHARGEDEPOSITION_H_
37 template <
int depos_order>
39 const amrex::ParticleReal *
const wp,
43 const std::array<amrex::Real, 3>&
dx,
44 const std::array<amrex::Real, 3> xyzmin,
47 int n_rz_azimuthal_modes)
49 using namespace amrex;
53 const bool do_ionization = ion_lev;
54 const amrex::Real dzi = 1.0_rt/
dx[2];
55 #if defined(WARPX_DIM_1D_Z)
56 const amrex::Real invvol = dzi;
58 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
59 const amrex::Real dxi = 1.0_rt/
dx[0];
60 const amrex::Real invvol = dxi*dzi;
61 #elif defined(WARPX_DIM_3D)
62 const amrex::Real dxi = 1.0_rt/
dx[0];
63 const amrex::Real dyi = 1.0_rt/
dx[1];
64 const amrex::Real invvol = dxi*dyi*dzi;
67 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) || defined(WARPX_DIM_3D)
68 const amrex::Real
xmin = xyzmin[0];
70 #if defined(WARPX_DIM_3D)
71 const amrex::Real ymin = xyzmin[1];
73 const amrex::Real zmin = xyzmin[2];
86 amrex::Real wq = q*wp[ip]*invvol;
91 amrex::ParticleReal xp, yp, zp;
92 GetPosition(ip, xp, yp, zp);
96 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) || defined(WARPX_DIM_3D)
99 #if defined(WARPX_DIM_RZ)
100 const amrex::Real rp = std::sqrt(xp*xp + yp*yp);
101 amrex::Real costheta;
102 amrex::Real sintheta;
111 const amrex::Real
x = (rp -
xmin)*dxi;
113 const amrex::Real
x = (xp -
xmin)*dxi;
118 amrex::Real sx[depos_order + 1] = {0._rt};
120 if (rho_type[0] ==
NODE) {
121 i = compute_shape_factor(sx,
x);
122 }
else if (rho_type[0] == CELL) {
123 i = compute_shape_factor(sx,
x - 0.5_rt);
126 #if defined(WARPX_DIM_3D)
128 const amrex::Real y = (yp - ymin)*dyi;
129 amrex::Real sy[depos_order + 1] = {0._rt};
131 if (rho_type[1] ==
NODE) {
132 j = compute_shape_factor(sy, y);
133 }
else if (rho_type[1] == CELL) {
134 j = compute_shape_factor(sy, y - 0.5_rt);
138 const amrex::Real z = (zp - zmin)*dzi;
139 amrex::Real sz[depos_order + 1] = {0._rt};
141 if (rho_type[WARPX_ZINDEX] ==
NODE) {
142 k = compute_shape_factor(sz, z);
143 }
else if (rho_type[WARPX_ZINDEX] == CELL) {
144 k = compute_shape_factor(sz, z - 0.5_rt);
148 #if defined(WARPX_DIM_1D_Z)
149 for (
int iz=0;
iz<=depos_order;
iz++){
151 &rho_arr(lo.
x+k+
iz, 0, 0, 0),
155 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
156 for (
int iz=0;
iz<=depos_order;
iz++){
157 for (
int ix=0;
ix<=depos_order;
ix++){
159 &rho_arr(lo.
x+
i+
ix, lo.
y+k+
iz, 0, 0),
161 #if defined(WARPX_DIM_RZ)
163 for (
int imode=1 ; imode < n_rz_azimuthal_modes ; imode++) {
172 #elif defined(WARPX_DIM_3D)
173 for (
int iz=0;
iz<=depos_order;
iz++){
174 for (
int iy=0;
iy<=depos_order;
iy++){
175 for (
int ix=0;
ix<=depos_order;
ix++){
212 template <
int depos_order>
214 const amrex::ParticleReal *
const wp,
218 const long np_to_deposit,
219 const std::array<amrex::Real, 3>&
dx,
220 const std::array<amrex::Real, 3> xyzmin,
223 const int n_rz_azimuthal_modes,
231 using namespace amrex;
235 #if !defined(AMREX_USE_GPU)
241 const bool do_ionization = ion_lev;
242 const amrex::Real dzi = 1.0_rt/
dx[2];
243 #if defined(WARPX_DIM_1D_Z)
244 const amrex::Real invvol = dzi;
246 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
247 const amrex::Real dxi = 1.0_rt/
dx[0];
248 const amrex::Real invvol = dxi*dzi;
249 #elif defined(WARPX_DIM_3D)
250 const amrex::Real dxi = 1.0_rt/
dx[0];
251 const amrex::Real dyi = 1.0_rt/
dx[1];
252 const amrex::Real invvol = dxi*dyi*dzi;
255 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) || defined(WARPX_DIM_3D)
256 const amrex::Real
xmin = xyzmin[0];
258 #if defined(WARPX_DIM_3D)
259 const amrex::Real ymin = xyzmin[1];
261 const amrex::Real zmin = xyzmin[2];
264 auto rho_box = rho_fab.
box();
271 #if defined(AMREX_USE_CUDA) || defined(AMREX_USE_HIP)
274 const auto domain = geom.
Domain();
279 sample_tbox_x.
grow(depos_order);
281 const auto npts = sample_tbox_x.
numPts();
283 const int nblocks = a_bins.
numBins();
285 const int threads_per_block = 256;
287 std::size_t shared_mem_bytes = npts*
sizeof(amrex::Real);
292 "Tile size too big for GPU shared memory charge deposition");
295 #if defined(AMREX_USE_CUDA) || defined(AMREX_USE_HIP)
307 #if defined(AMREX_USE_CUDA) || defined(AMREX_USE_HIP)
308 const int bin_id = blockIdx.x;
309 const unsigned int bin_start = offsets_ptr[bin_id];
310 const unsigned int bin_stop = offsets_ptr[bin_id+1];
312 if (bin_start == bin_stop) {
return; }
316 ParticleReal xp, yp, zp;
317 GetPosition(permutation[bin_start], xp, yp, zp);
318 #if defined(WARPX_DIM_3D)
319 IntVect iv =
IntVect(
int(amrex::Math::floor((xp-plo[0])*dxiarr[0])),
320 int(amrex::Math::floor((yp-plo[1])*dxiarr[1])),
321 int(amrex::Math::floor((zp-plo[2])*dxiarr[2])));
322 #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
324 int(amrex::Math::floor((xp-plo[0])*dxiarr[0])),
325 int(amrex::Math::floor((zp-plo[1])*dxiarr[1])));
326 #elif defined(WARPX_DIM_1D_Z)
327 IntVect iv =
IntVect(
int(amrex::Math::floor((zp-plo[0])*dxiarr[0])));
329 iv += domain.smallEnd();
334 tbx.
grow(depos_order);
337 amrex::Real*
const shared = gsm.
dataPtr();
342 volatile amrex::Real* vs = shared;
343 for (
int i = threadIdx.x;
i < tbx.
numPts();
i += blockDim.x) {
351 #if defined(AMREX_USE_CUDA) || defined(AMREX_USE_HIP)
353 for (
unsigned int ip_orig = bin_start + threadIdx.x; ip_orig < bin_stop; ip_orig += blockDim.x)
356 const unsigned int ip = permutation[ip_orig];
359 amrex::Real wq = q*wp[ip]*invvol;
364 amrex::ParticleReal xp, yp, zp;
365 GetPosition(ip, xp, yp, zp);
369 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) || defined(WARPX_DIM_3D)
372 #if defined(WARPX_DIM_RZ)
373 const amrex::Real rp = std::sqrt(xp*xp + yp*yp);
374 amrex::Real costheta;
375 amrex::Real sintheta;
384 const amrex::Real
x = (rp -
xmin)*dxi;
386 const amrex::Real
x = (xp -
xmin)*dxi;
391 amrex::Real sx[depos_order + 1] = {0._rt};
393 if (rho_type[0] ==
NODE) {
394 i = compute_shape_factor(sx,
x);
395 }
else if (rho_type[0] == CELL) {
396 i = compute_shape_factor(sx,
x - 0.5_rt);
399 #if defined(WARPX_DIM_3D)
401 const amrex::Real y = (yp - ymin)*dyi;
402 amrex::Real sy[depos_order + 1] = {0._rt};
404 if (rho_type[1] ==
NODE) {
405 j = compute_shape_factor(sy, y);
406 }
else if (rho_type[1] == CELL) {
407 j = compute_shape_factor(sy, y - 0.5_rt);
411 const amrex::Real z = (zp - zmin)*dzi;
412 amrex::Real sz[depos_order + 1] = {0._rt};
414 if (rho_type[WARPX_ZINDEX] ==
NODE) {
415 k = compute_shape_factor(sz, z);
416 }
else if (rho_type[WARPX_ZINDEX] == CELL) {
417 k = compute_shape_factor(sz, z - 0.5_rt);
421 #if defined(WARPX_DIM_1D_Z)
422 for (
int iz=0;
iz<=depos_order;
iz++){
424 &buf(lo.
x+k+
iz, 0, 0, 0),
428 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
429 for (
int iz=0;
iz<=depos_order;
iz++){
430 for (
int ix=0;
ix<=depos_order;
ix++){
432 &buf(lo.
x+
i+
ix, lo.
y+k+
iz, 0, 0),
434 #if defined(WARPX_DIM_RZ)
436 for (
int imode=1 ; imode < n_rz_azimuthal_modes ; imode++) {
445 #elif defined(WARPX_DIM_3D)
446 for (
int iz=0;
iz<=depos_order;
iz++){
447 for (
int iy=0;
iy<=depos_order;
iy++){
448 for (
int ix=0;
ix<=depos_order;
ix++){
458 #if defined(AMREX_USE_CUDA) || defined(AMREX_USE_HIP)
461 addLocalToGlobal(tbx, rho_arr, buf);
#define AMREX_D_DECL(a, b, c)
void doChargeDepositionShapeN(const GetParticlePosition< PIdx > &GetPosition, const amrex::ParticleReal *const wp, const int *ion_lev, amrex::FArrayBox &rho_fab, long np_to_deposit, const std::array< amrex::Real, 3 > &dx, const std::array< amrex::Real, 3 > xyzmin, amrex::Dim3 lo, amrex::Real q, int n_rz_azimuthal_modes)
Definition: ChargeDeposition.H:38
void doChargeDepositionSharedShapeN(const GetParticlePosition< PIdx > &GetPosition, const amrex::ParticleReal *const wp, const int *ion_lev, amrex::FArrayBox &rho_fab, const amrex::IntVect &ix_type, const long np_to_deposit, const std::array< amrex::Real, 3 > &dx, const std::array< amrex::Real, 3 > xyzmin, const amrex::Dim3 lo, const amrex::Real q, const int n_rz_azimuthal_modes, const amrex::DenseBins< WarpXParticleContainer::ParticleTileType::ParticleTileDataType > &a_bins, const amrex::Box &box, const amrex::Geometry &geom, const amrex::IntVect &a_tbox_max_size, const amrex::IntVect bin_size)
Definition: ChargeDeposition.H:213
#define WARPX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition: TextMsg.H:13
AMREX_FORCE_INLINE Array4< Real const > array() const noexcept
const Box & box() const noexcept
AMREX_GPU_HOST_DEVICE IntVect type() const noexcept
AMREX_GPU_HOST_DEVICE Long numPts() const noexcept
AMREX_GPU_HOST_DEVICE Box & grow(int i) noexcept
GpuArray< Real, AMREX_SPACEDIM > InvCellSizeArray() const noexcept
index_type * offsetsPtr() noexcept
index_type * permutationPtr() noexcept
Long numBins() const noexcept
GpuArray< Real, AMREX_SPACEDIM > ProbLoArray() const noexcept
const Box & Domain() const noexcept
static std::size_t sharedMemPerBlock() noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void AddNoRet(T *sum, T value) noexcept
gpuStream_t gpuStream() noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(Box const &box) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Box convert(const Box &b, const IntVect &typ) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(Box const &box) noexcept
void launch(T const &n, L &&f) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... >, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int getTileIndex(const IntVect &iv, const Box &box, const bool a_do_tiling, const IntVect &a_tile_size, Box &tbx)
i
Definition: check_interp_points_and_weights.py:174
tuple dx
lab frame
Definition: stencil.py:429
xmin
Definition: stencil.py:424
Definition: ShapeFactors.H:29
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE T real() const noexcept
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE T imag() const noexcept
AMREX_GPU_DEVICE T * dataPtr() noexcept