aboutsummaryrefslogtreecommitdiff
path: root/Source/FieldSolver/WarpXPushFieldsEM.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'Source/FieldSolver/WarpXPushFieldsEM.cpp')
-rw-r--r--Source/FieldSolver/WarpXPushFieldsEM.cpp91
1 files changed, 77 insertions, 14 deletions
diff --git a/Source/FieldSolver/WarpXPushFieldsEM.cpp b/Source/FieldSolver/WarpXPushFieldsEM.cpp
index 11d598b98..b3d8e9d76 100644
--- a/Source/FieldSolver/WarpXPushFieldsEM.cpp
+++ b/Source/FieldSolver/WarpXPushFieldsEM.cpp
@@ -175,18 +175,19 @@ WarpX::EvolveB (int lev, PatchType patch_type, amrex::Real a_dt)
warpx_push_bz_nodal(j,k,l,Bzfab,Exfab,Eyfab,dtsdx,dtsdy);
});
} else if (WarpX::maxwell_fdtd_solver_id == 0) {
+ const long nmodes = n_rz_azimuthal_modes;
amrex::ParallelFor(tbx, tby, tbz,
[=] AMREX_GPU_DEVICE (int j, int k, int l)
{
- warpx_push_bx_yee(j,k,l,Bxfab,Eyfab,Ezfab,dtsdy,dtsdz);
+ warpx_push_bx_yee(j,k,l,Bxfab,Eyfab,Ezfab,dtsdx,dtsdy,dtsdz,dxinv,xmin,nmodes);
},
[=] AMREX_GPU_DEVICE (int j, int k, int l)
{
- warpx_push_by_yee(j,k,l,Byfab,Exfab,Ezfab,dtsdx,dtsdz);
+ warpx_push_by_yee(j,k,l,Byfab,Exfab,Ezfab,dtsdx,dtsdz,nmodes);
},
[=] AMREX_GPU_DEVICE (int j, int k, int l)
{
- warpx_push_bz_yee(j,k,l,Bzfab,Exfab,Eyfab,dtsdx,dtsdy,dxinv,xmin);
+ warpx_push_bz_yee(j,k,l,Bzfab,Exfab,Eyfab,dtsdx,dtsdy,dxinv,xmin,nmodes);
});
} else if (WarpX::maxwell_fdtd_solver_id == 1) {
Real betaxy, betaxz, betayx, betayz, betazx, betazy;
@@ -366,18 +367,19 @@ WarpX::EvolveE (int lev, PatchType patch_type, amrex::Real a_dt)
warpx_push_ez_nodal(j,k,l,Ezfab,Bxfab,Byfab,jzfab,mu_c2_dt,dtsdx_c2,dtsdy_c2);
});
} else {
+ const long nmodes = n_rz_azimuthal_modes;
amrex::ParallelFor(tex, tey, tez,
[=] AMREX_GPU_DEVICE (int j, int k, int l)
{
- warpx_push_ex_yee(j,k,l,Exfab,Byfab,Bzfab,jxfab,mu_c2_dt,dtsdy_c2,dtsdz_c2);
+ warpx_push_ex_yee(j,k,l,Exfab,Byfab,Bzfab,jxfab,mu_c2_dt,dtsdx_c2,dtsdy_c2,dtsdz_c2,dxinv,xmin,nmodes);
},
[=] AMREX_GPU_DEVICE (int j, int k, int l)
{
- warpx_push_ey_yee(j,k,l,Eyfab,Bxfab,Bzfab,jyfab,mu_c2_dt,dtsdx_c2,dtsdz_c2,xmin);
+ warpx_push_ey_yee(j,k,l,Eyfab,Bxfab,Bzfab,jyfab,Exfab,mu_c2_dt,dtsdx_c2,dtsdz_c2,xmin,nmodes);
},
[=] AMREX_GPU_DEVICE (int j, int k, int l)
{
- warpx_push_ez_yee(j,k,l,Ezfab,Bxfab,Byfab,jzfab,mu_c2_dt,dtsdx_c2,dtsdy_c2,dxinv,xmin);
+ warpx_push_ez_yee(j,k,l,Ezfab,Bxfab,Byfab,jzfab,mu_c2_dt,dtsdx_c2,dtsdy_c2,dxinv,xmin,nmodes);
});
}
@@ -647,6 +649,8 @@ WarpX::ApplyInverseVolumeScalingToCurrentDensity (MultiFab* Jx, MultiFab* Jy, Mu
const Real rmin = xyzmin[0];
const int irmin = lo.x;
+ const long nmodes = n_rz_azimuthal_modes;
+
// Rescale current in r-z mode since the inverse volume factor was not
// included in the current deposition.
amrex::ParallelFor(tbr, tbt, tbz,
@@ -656,13 +660,29 @@ WarpX::ApplyInverseVolumeScalingToCurrentDensity (MultiFab* Jx, MultiFab* Jy, Mu
// to the cells above the axis.
// Note that Jr(i==0) is at 1/2 dr.
if (rmin == 0. && 0 <= i && i < ngJ) {
- Jr_arr(i,j,0) -= Jr_arr(-1-i,j,0);
+ Jr_arr(i,j,0,0) -= Jr_arr(-1-i,j,0,0);
}
// Apply the inverse volume scaling
// Since Jr is not cell centered in r, no need for distinction
// between on axis and off-axis factors
const amrex::Real r = std::abs(rmin + (i - irmin + 0.5)*dr);
- Jr_arr(i,j,0) /= (2.*MathConst::pi*r);
+ Jr_arr(i,j,0,0) /= (2.*MathConst::pi*r);
+
+ for (int imode=1 ; imode < nmodes ; imode++) {
+ const Real ifact = ( (imode%2) == 0 ? +1. : -1.);
+ // Wrap the current density deposited in the guard cells around
+ // to the cells above the axis.
+ // Note that Jr(i==0) is at 1/2 dr.
+ if (rmin == 0. && 0 <= i && i < ngJ) {
+ Jr_arr(i,j,0,2*imode-1) -= ifact*Jr_arr(-1-i,j,0,2*imode-1);
+ Jr_arr(i,j,0,2*imode) -= ifact*Jr_arr(-1-i,j,0,2*imode);
+ }
+ // Apply the inverse volume scaling
+ // Since Jr is not cell centered in r, no need for distinction
+ // between on axis and off-axis factors
+ Jr_arr(i,j,0,2*imode-1) /= (2.*MathConst::pi*r);
+ Jr_arr(i,j,0,2*imode) /= (2.*MathConst::pi*r);
+ }
},
[=] AMREX_GPU_DEVICE (int i, int j, int k)
{
@@ -670,16 +690,37 @@ WarpX::ApplyInverseVolumeScalingToCurrentDensity (MultiFab* Jx, MultiFab* Jy, Mu
// to the cells above the axis.
// Jt is located on the boundary
if (rmin == 0. && 0 < i && i <= ngJ) {
- Jt_arr(i,j,0) += Jt_arr(-i,j,0);
+ Jt_arr(i,j,0,0) += Jt_arr(-i,j,0,0);
}
// Apply the inverse volume scaling
// Jt is forced to zero on axis.
const amrex::Real r = std::abs(rmin + (i - irmin)*dr);
if (r == 0.) {
- Jt_arr(i,j,0) = 0.;
+ Jt_arr(i,j,0,0) = 0.;
} else {
- Jt_arr(i,j,0) /= (2.*MathConst::pi*r);
+ Jt_arr(i,j,0,0) /= (2.*MathConst::pi*r);
+ }
+
+ for (int imode=1 ; imode < nmodes ; imode++) {
+ const Real ifact = ( (imode%2) == 0 ? +1. : -1.);
+ // Wrap the current density deposited in the guard cells around
+ // to the cells above the axis.
+ // Jt is located on the boundary
+ if (rmin == 0. && 0 < i && i <= ngJ) {
+ Jt_arr(i,j,0,2*imode-1) += ifact*Jt_arr(-i,j,0,2*imode-1);
+ Jt_arr(i,j,0,2*imode) += ifact*Jt_arr(-i,j,0,2*imode);
+ }
+
+ // Apply the inverse volume scaling
+ // Jt is forced to zero on axis.
+ if (r == 0.) {
+ Jt_arr(i,j,0,2*imode-1) = 0.;
+ Jt_arr(i,j,0,2*imode) = 0.;
+ } else {
+ Jt_arr(i,j,0,2*imode-1) /= (2.*MathConst::pi*r);
+ Jt_arr(i,j,0,2*imode) /= (2.*MathConst::pi*r);
+ }
}
},
[=] AMREX_GPU_DEVICE (int i, int j, int k)
@@ -688,17 +729,39 @@ WarpX::ApplyInverseVolumeScalingToCurrentDensity (MultiFab* Jx, MultiFab* Jy, Mu
// to the cells above the axis.
// Jz is located on the boundary
if (rmin == 0. && 0 < i && i <= ngJ) {
- Jz_arr(i,j,0) += Jz_arr(-i,j,0);
+ Jz_arr(i,j,0,0) += Jz_arr(-i,j,0,0);
}
// Apply the inverse volume scaling
const amrex::Real r = std::abs(rmin + (i - irmin)*dr);
if (r == 0.) {
// Verboncoeur JCP 164, 421-427 (2001) : corrected volume on axis
- Jz_arr(i,j,0) /= (MathConst::pi*dr/3.);
+ Jz_arr(i,j,0,0) /= (MathConst::pi*dr/3.);
} else {
- Jz_arr(i,j,0) /= (2.*MathConst::pi*r);
+ Jz_arr(i,j,0,0) /= (2.*MathConst::pi*r);
}
+
+ for (int imode=1 ; imode < nmodes ; imode++) {
+ const Real ifact = ( (imode%2) == 0 ? +1. : -1.);
+ // Wrap the current density deposited in the guard cells around
+ // to the cells above the axis.
+ // Jz is located on the boundary
+ if (rmin == 0. && 0 < i && i <= ngJ) {
+ Jz_arr(i,j,0,2*imode-1) += ifact*Jz_arr(-i,j,0,2*imode-1);
+ Jz_arr(i,j,0,2*imode) += ifact*Jz_arr(-i,j,0,2*imode);
+ }
+
+ // Apply the inverse volume scaling
+ if (r == 0.) {
+ // Verboncoeur JCP 164, 421-427 (2001) : corrected volume on axis
+ Jz_arr(i,j,0,2*imode-1) /= (MathConst::pi*dr/3.);
+ Jz_arr(i,j,0,2*imode) /= (MathConst::pi*dr/3.);
+ } else {
+ Jz_arr(i,j,0,2*imode-1) /= (2.*MathConst::pi*r);
+ Jz_arr(i,j,0,2*imode) /= (2.*MathConst::pi*r);
+ }
+ }
+
});
}
}