aboutsummaryrefslogtreecommitdiff
path: root/Source/Particles/Deposition/CurrentDeposition.H
diff options
context:
space:
mode:
Diffstat (limited to 'Source/Particles/Deposition/CurrentDeposition.H')
-rw-r--r--Source/Particles/Deposition/CurrentDeposition.H161
1 files changed, 57 insertions, 104 deletions
diff --git a/Source/Particles/Deposition/CurrentDeposition.H b/Source/Particles/Deposition/CurrentDeposition.H
index b4531d3b3..58621c76b 100644
--- a/Source/Particles/Deposition/CurrentDeposition.H
+++ b/Source/Particles/Deposition/CurrentDeposition.H
@@ -208,27 +208,34 @@ void doEsirkepovDepositionShapeN (const Real * const xp, const Real * const yp,
const Dim3 lo,
const amrex::Real q)
{
-
-#if (AMREX_SPACEDIM == 3)
-
const Real dxi = 1.0/dx[0];
- const Real dyi = 1.0/dx[1];
- const Real dzi = 1.0/dx[2];
const Real dtsdx0 = dt*dxi;
+ const Real xmin = xyzmin[0];
+#if (AMREX_SPACEDIM == 3)
+ const Real dyi = 1.0/dx[1];
const Real dtsdy0 = dt*dyi;
+ const Real ymin = xyzmin[1];
+#endif
+ const Real dzi = 1.0/dx[2];
const Real dtsdz0 = dt*dzi;
+ const Real zmin = xyzmin[2];
+
+#if (AMREX_SPACEDIM == 3)
const Real invdtdx = 1.0/(dt*dx[1]*dx[2]);
const Real invdtdy = 1.0/(dt*dx[0]*dx[2]);
const Real invdtdz = 1.0/(dt*dx[0]*dx[1]);
+#elif (AMREX_SPACEDIM == 2)
+ const Real invdtdx = 1.0/(dt*dx[2]);
+ const Real invdtdz = 1.0/(dt*dx[0]);
+ const Real invvol = 1.0/(dx[0]*dx[2]);
+#endif
- const Real xmin = xyzmin[0];
- const Real ymin = xyzmin[1];
- const Real zmin = xyzmin[2];
const Real clightsq = 1.0/PhysConst::c/PhysConst::c;
// Loop over particles and deposit into Jx_arr, Jy_arr and Jz_arr
ParallelFor( np_to_depose,
[=] AMREX_GPU_DEVICE (long ip) {
+
// --- Get particle quantities
const Real gaminv = 1.0/std::sqrt(1.0 + uxp[ip]*uxp[ip]*clightsq
+ uyp[ip]*uyp[ip]*clightsq
@@ -237,17 +244,19 @@ void doEsirkepovDepositionShapeN (const Real * const xp, const Real * const yp,
// wqx, wqy wqz are particle current in each direction
const Real wq = q*wp[ip];
const Real wqx = wq*invdtdx;
+#if (AMREX_SPACEDIM == 3)
const Real wqy = wq*invdtdy;
+#endif
const Real wqz = wq*invdtdz;
- // computes current position in grid units
+ // computes current and old position in grid units
const Real x_new = (xp[ip] - xmin)*dxi;
- const Real y_new = (yp[ip] - ymin)*dyi;
- const Real z_new = (zp[ip] - zmin)*dzi;
-
- // computes old position in grid units
const Real x_old = x_new - dtsdx0*uxp[ip]*gaminv;
+#if (AMREX_SPACEDIM == 3)
+ const Real y_new = (yp[ip] - ymin)*dyi;
const Real y_old = y_new - dtsdy0*uyp[ip]*gaminv;
+#endif
+ const Real z_new = (zp[ip] - zmin)*dzi;
const Real z_old = z_new - dtsdz0*uzp[ip]*gaminv;
// Shape factor arrays
@@ -255,34 +264,41 @@ void doEsirkepovDepositionShapeN (const Real * const xp, const Real * const yp,
// to possibly hold the factor for the old particle
// which can be at a different grid location.
Real AMREX_RESTRICT sx_new[depos_order + 3] = {0.};
- Real AMREX_RESTRICT sy_new[depos_order + 3] = {0.};
- Real AMREX_RESTRICT sz_new[depos_order + 3] = {0.};
Real AMREX_RESTRICT sx_old[depos_order + 3] = {0.};
+#if (AMREX_SPACEDIM == 3)
+ Real AMREX_RESTRICT sy_new[depos_order + 3] = {0.};
Real AMREX_RESTRICT sy_old[depos_order + 3] = {0.};
+#endif
+ Real AMREX_RESTRICT sz_new[depos_order + 3] = {0.};
Real AMREX_RESTRICT sz_old[depos_order + 3] = {0.};
// --- Compute shape factors
// Compute shape factors for position as they are now and at old positions
// [ijk]_new: leftmost grid point that the particle touches
const int i_new = compute_shape_factor<depos_order>(sx_new+1, x_new);
- const int j_new = compute_shape_factor<depos_order>(sy_new+1, y_new);
- const int k_new = compute_shape_factor<depos_order>(sz_new+1, z_new);
-
const int i_old = compute_shifted_shape_factor<depos_order>(sx_old, x_old, i_new);
+#if (AMREX_SPACEDIM == 3)
+ const int j_new = compute_shape_factor<depos_order>(sy_new+1, y_new);
const int j_old = compute_shifted_shape_factor<depos_order>(sy_old, y_old, j_new);
+#endif
+ const int k_new = compute_shape_factor<depos_order>(sz_new+1, z_new);
const int k_old = compute_shifted_shape_factor<depos_order>(sz_old, z_old, k_new);
// computes min/max positions of current contributions
int dil = 1, diu = 1;
if (i_old < i_new) dil = 0;
if (i_old > i_new) diu = 0;
+#if (AMREX_SPACEDIM == 3)
int djl = 1, dju = 1;
if (j_old < j_new) djl = 0;
if (j_old > j_new) dju = 0;
+#endif
int dkl = 1, dku = 1;
if (k_old < k_new) dkl = 0;
if (k_old > k_new) dku = 0;
+#if (AMREX_SPACEDIM == 3)
+
for (int k=dkl; k<=depos_order+2-dku; k++) {
for (int j=djl; j<=depos_order+2-dju; j++) {
Real sdxi = 0.;
@@ -293,7 +309,6 @@ void doEsirkepovDepositionShapeN (const Real * const xp, const Real * const yp,
}
}
}
-
for (int k=dkl; k<=depos_order+2-dku; k++) {
for (int i=dil; i<=depos_order+2-diu; i++) {
Real sdyj = 0.;
@@ -304,7 +319,6 @@ void doEsirkepovDepositionShapeN (const Real * const xp, const Real * const yp,
}
}
}
-
for (int j=djl; j<=depos_order+2-dju; j++) {
for (int i=dil; i<=depos_order+2-diu; i++) {
Real sdzk = 0.;
@@ -315,97 +329,36 @@ void doEsirkepovDepositionShapeN (const Real * const xp, const Real * const yp,
}
}
}
- }
- );
#elif (AMREX_SPACEDIM == 2)
- const Real dxi = 1.0/dx[0];
- const Real dzi = 1.0/dx[2];
- const Real dtsdx0 = dt*dxi;
- const Real dtsdz0 = dt*dzi;
- const Real invdtdx = 1.0/(dt*dx[2]);
- const Real invdtdz = 1.0/(dt*dx[0]);
- const Real invvol = 1.0/(dx[0]*dx[2]);
+ for (int k=dkl; k<=depos_order+2-dku; k++) {
+ Real sdxi = 0.;
+ for (int i=dil; i<=depos_order+1-diu; i++) {
+ sdxi += wqx*(sx_old[i] - sx_new[i])*(sz_new[k] + 0.5*(sz_old[k] - sz_new[k]));
+ amrex::Gpu::Atomic::Add( &Jx_arr(lo.x+i_new-1+i, lo.y+k_new-1+k, 0), sdxi);
+ }
+ }
+ for (int k=dkl; k<=depos_order+2-dku; k++) {
+ for (int i=dil; i<=depos_order+2-diu; i++) {
+ const Real sdyj = wq*uyp[ip]*gaminv*invvol*((sz_new[k] + 0.5*(sz_old[k] - sz_new[k]))*sx_new[i] +
+ (0.5*sz_new[k] + 1./3.*(sz_old[k] - sz_new[k]))*(sx_old[i] - sx_new[i]));
+ amrex::Gpu::Atomic::Add( &Jy_arr(lo.x+i_new-1+i, lo.y+k_new-1+k, 0), sdyj);
+ }
+ }
+ for (int i=dil; i<=depos_order+2-diu; i++) {
+ Real sdzk = 0.;
+ for (int k=dkl; k<=depos_order+1-dku; k++) {
+ sdzk += wqz*(sz_old[k] - sz_new[k])*(sx_new[i] + 0.5*(sx_old[i] - sx_new[i]));
+ amrex::Gpu::Atomic::Add( &Jz_arr(lo.x+i_new-1+i, lo.y+k_new-1+k, 0), sdzk);
+ }
+ }
- const Real xmin = xyzmin[0];
- const Real zmin = xyzmin[2];
- const Real clightsq = 1.0/PhysConst::c/PhysConst::c;
-
- // Loop over particles and deposit into Jx_arr, Jy_arr and Jz_arr
- ParallelFor( np_to_depose,
- [=] AMREX_GPU_DEVICE (long ip) {
- // --- Get particle quantities
- const Real gaminv = 1.0/std::sqrt(1.0 + uxp[ip]*uxp[ip]*clightsq
- + uyp[ip]*uyp[ip]*clightsq
- + uzp[ip]*uzp[ip]*clightsq);
-
- // wqx, wqy wqz are particle current in each direction
- const Real wq = q*wp[ip];
- const Real wqx = wq*invdtdx;
- const Real wqz = wq*invdtdz;
-
- // computes current position in grid units
- const Real x_new = (xp[ip] - xmin)*dxi;
- const Real z_new = (zp[ip] - zmin)*dzi;
-
- // computes old position in grid units
- const Real x_old = x_new - dtsdx0*uxp[ip]*gaminv;
- const Real z_old = z_new - dtsdz0*uzp[ip]*gaminv;
-
- // Shape factor arrays
- // Note that there are extra values above and below
- // to possibly hold the factor for the old particle
- // which can be at a different grid location.
- Real AMREX_RESTRICT sx_new[depos_order + 3] = {0.};
- Real AMREX_RESTRICT sz_new[depos_order + 3] = {0.};
- Real AMREX_RESTRICT sx_old[depos_order + 3] = {0.};
- Real AMREX_RESTRICT sz_old[depos_order + 3] = {0.};
-
- // --- Compute shape factors
- // Compute shape factors for position as they are now and at old positions
- // [ik]_new: leftmost grid point that the particle touches
- const int i_new = compute_shape_factor<depos_order>(sx_new+1, x_new);
- const int k_new = compute_shape_factor<depos_order>(sz_new+1, z_new);
-
- const int i_old = compute_shifted_shape_factor<depos_order>(sx_old, x_old, i_new);
- const int k_old = compute_shifted_shape_factor<depos_order>(sz_old, z_old, k_new);
-
- // computes min/max positions of current contributions
- int dil = 1, diu = 1;
- if (i_old < i_new) dil = 0;
- if (i_old > i_new) diu = 0;
- int dkl = 1, dku = 1;
- if (k_old < k_new) dkl = 0;
- if (k_old > k_new) dku = 0;
-
- for (int k=dkl; k<=depos_order+2-dku; k++) {
- Real sdxi = 0.;
- for (int i=dil; i<=depos_order+1-diu; i++) {
- sdxi += wqx*(sx_old[i] - sx_new[i])*(sz_new[k] + 0.5*(sz_old[k] - sz_new[k]));
- amrex::Gpu::Atomic::Add( &Jx_arr(lo.x+i_new-1+i, lo.y+k_new-1+k, 0), sdxi);
- }
- }
-
- for (int k=dkl; k<=depos_order+2-dku; k++) {
- for (int i=dil; i<=depos_order+2-diu; i++) {
- const Real sdyj = wq*uyp[ip]*gaminv*invvol*((sz_new[k] + 0.5*(sz_old[k] - sz_new[k]))*sx_new[i] +
- (0.5*sz_new[k] + 1./3.*(sz_old[k] - sz_new[k]))*(sx_old[i] - sx_new[i]));
- amrex::Gpu::Atomic::Add( &Jy_arr(lo.x+i_new-1+i, lo.y+k_new-1+k, 0), sdyj);
- }
- }
-
- for (int i=dil; i<=depos_order+2-diu; i++) {
- Real sdzk = 0.;
- for (int k=dkl; k<=depos_order+1-dku; k++) {
- sdzk += wqz*(sz_old[k] - sz_new[k])*(sx_new[i] + 0.5*(sx_old[i] - sx_new[i]));
- amrex::Gpu::Atomic::Add( &Jz_arr(lo.x+i_new-1+i, lo.y+k_new-1+k, 0), sdzk);
- }
- }
+#endif
}
);
-#endif
+
}