7 #ifndef WARPX_SHAREDDEPOSITIONUTILS_H_
8 #define WARPX_SHAREDDEPOSITIONUTILS_H_
28 const int nTiles = nCells / tilesize;
29 const int remainder = nCells % tilesize;
30 maxTilesize = tilesize +
int(std::ceil((amrex::Real) remainder / nTiles));
40 #if defined(AMREX_USE_HIP) || defined(AMREX_USE_CUDA)
46 using namespace amrex::literals;
50 for (
int icell = threadIdx.x; icell < bx.numPts(); icell += blockDim.x)
52 int k = icell / (len.x*len.y);
53 int j = (icell - k*(len.x*len.y)) / len.x;
54 int i = (icell - k*(len.x*len.y)) - j*len.x;
58 if (amrex::Math::abs(local(
i, j, k)) > 0.0_rt) {
65 #if defined(AMREX_USE_HIP) || defined(AMREX_USE_CUDA)
66 template <
int depos_order>
69 const amrex::ParticleReal *
const wp,
70 const amrex::ParticleReal *
const uxp,
71 const amrex::ParticleReal *
const uyp,
72 const amrex::ParticleReal *
const uzp,
76 const amrex::Real relative_time,
77 const std::array<amrex::Real,3>&
dx,
78 const std::array<amrex::Real,3>& xyzmin,
81 const int n_rz_azimuthal_modes,
82 const unsigned int ip,
83 const int zdir,
const int NODE,
const int CELL,
const int dir)
85 using namespace amrex::literals;
87 #if !defined(WARPX_DIM_RZ)
93 const bool do_ionization = ion_lev;
94 const amrex::Real dzi = 1.0_rt/
dx[2];
95 #if defined(WARPX_DIM_1D_Z)
96 const amrex::Real invvol = dzi;
98 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
99 const amrex::Real dxi = 1.0_rt/
dx[0];
100 const amrex::Real invvol = dxi*dzi;
101 #elif defined(WARPX_DIM_3D)
102 const amrex::Real dxi = 1.0_rt/
dx[0];
103 const amrex::Real dyi = 1.0_rt/
dx[1];
104 const amrex::Real invvol = dxi*dyi*dzi;
107 #if (AMREX_SPACEDIM >= 2)
108 const amrex::Real
xmin = xyzmin[0];
110 #if defined(WARPX_DIM_3D)
111 const amrex::Real ymin = xyzmin[1];
113 const amrex::Real zmin = xyzmin[2];
118 const amrex::Real gaminv = 1.0_rt/std::sqrt(1.0_rt + uxp[ip]*uxp[ip]*clightsq
119 + uyp[ip]*uyp[ip]*clightsq
120 + uzp[ip]*uzp[ip]*clightsq);
121 amrex::Real wq = q*wp[ip];
126 amrex::ParticleReal xp, yp, zp;
127 GetPosition(ip, xp, yp, zp);
129 const amrex::Real vx = uxp[ip]*gaminv;
130 const amrex::Real vy = uyp[ip]*gaminv;
131 const amrex::Real vz = uzp[ip]*gaminv;
133 #if defined(WARPX_DIM_RZ)
136 const amrex::Real xpmid = xp + relative_time*vx;
137 const amrex::Real ypmid = yp + relative_time*vy;
138 const amrex::Real rpmid = std::sqrt(xpmid*xpmid + ypmid*ypmid);
139 amrex::Real costheta;
140 amrex::Real sintheta;
142 costheta = xpmid/rpmid;
143 sintheta = ypmid/rpmid;
149 const amrex::Real wqx = wq*invvol*(+vx*costheta + vy*sintheta);
150 const amrex::Real wqy = wq*invvol*(-vx*sintheta + vy*costheta);
152 const amrex::Real wqx = wq*invvol*vx;
153 const amrex::Real wqy = wq*invvol*vy;
155 const amrex::Real wqz = wq*invvol*vz;
157 amrex::Real pcurrent = 0.0;
160 }
else if (dir == 1) {
162 }
else if (dir == 2) {
168 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) || defined(WARPX_DIM_3D)
172 #if defined(WARPX_DIM_RZ)
174 const double xmid = (rpmid -
xmin)*dxi;
176 const double xmid = ((xp -
xmin) + relative_time*vx)*dxi;
183 double sx_node[depos_order + 1] = {0.};
184 double sx_cell[depos_order + 1] = {0.};
187 if (j_type[0] == NODE) {
188 j_node = compute_shape_factor(sx_node, xmid);
190 if (j_type[0] == CELL) {
191 j_cell = compute_shape_factor(sx_cell, xmid - 0.5);
194 amrex::Real sx_j[depos_order + 1] = {0._rt};
195 for (
int ix=0;
ix<=depos_order;
ix++)
197 sx_j[
ix] = ((j_type[0] ==
NODE) ? amrex::Real(sx_node[ix]) : amrex::Real(sx_cell[ix]));
200 int const j_j = ((j_type[0] ==
NODE) ? j_node : j_cell);
203 #if defined(WARPX_DIM_3D)
206 const double ymid = ((yp - ymin) + relative_time*vy)*dyi;
207 double sy_node[depos_order + 1] = {0.};
208 double sy_cell[depos_order + 1] = {0.};
211 if (j_type[1] == NODE) {
212 k_node = compute_shape_factor(sy_node, ymid);
214 if (j_type[1] == CELL) {
215 k_cell = compute_shape_factor(sy_cell, ymid - 0.5);
217 amrex::Real sy_j[depos_order + 1] = {0._rt};
218 for (
int iy=0;
iy<=depos_order;
iy++)
220 sy_j[
iy] = ((j_type[1] ==
NODE) ? amrex::Real(sy_node[iy]) : amrex::Real(sy_cell[iy]));
222 int const k_j = ((j_type[1] ==
NODE) ? k_node : k_cell);
227 const double zmid = ((zp - zmin) + relative_time*vz)*dzi;
228 double sz_node[depos_order + 1] = {0.};
229 double sz_cell[depos_order + 1] = {0.};
232 if (j_type[zdir] == NODE) {
233 l_node = compute_shape_factor(sz_node, zmid);
235 if (j_type[zdir] == CELL) {
236 l_cell = compute_shape_factor(sz_cell, zmid - 0.5);
238 amrex::Real sz_j[depos_order + 1] = {0._rt};
239 for (
int iz=0;
iz<=depos_order;
iz++)
241 sz_j[
iz] = ((j_type[zdir] ==
NODE) ? amrex::Real(sz_node[iz]) : amrex::Real(sz_cell[iz]));
243 int const l_j = ((j_type[zdir] ==
NODE) ? l_node : l_cell);
246 #if defined(WARPX_DIM_1D_Z)
247 for (
int iz=0;
iz<=depos_order;
iz++){
249 &j_buff(lo.
x+l_j+iz, 0, 0, 0),
253 #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
254 for (
int iz=0;
iz<=depos_order;
iz++){
255 for (
int ix=0;
ix<=depos_order;
ix++){
257 &j_buff(lo.
x+j_j+ix, lo.
y+l_j+iz, 0, 0),
258 sx_j[ix]*sz_j[iz]*pcurrent);
259 #if defined(WARPX_DIM_RZ)
261 for (
int imode=1 ; imode < n_rz_azimuthal_modes ; imode++) {
270 #elif defined(WARPX_DIM_3D)
271 for (
int iz=0;
iz<=depos_order;
iz++){
272 for (
int iy=0;
iy<=depos_order;
iy++){
273 for (
int ix=0;
ix<=depos_order;
ix++){
275 &j_buff(lo.
x+j_j+ix, lo.
y+k_j+iy, lo.
z+l_j+iz),
276 sx_j[ix]*sy_j[iy]*sz_j[iz]*pcurrent);
#define AMREX_FORCE_INLINE
#define AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE int getMaxTboxAlongDim(int nCells, int tilesize)
Definition: SharedDepositionUtils.H:26
static constexpr auto c
vacuum speed of light [m/s]
Definition: constant.H:44
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void AddNoRet(T *sum, T value) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 length(Array4< T > const &a) noexcept
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