diff options
Diffstat (limited to 'Source/Particles/RigidInjectedParticleContainer.cpp')
| -rw-r--r-- | Source/Particles/RigidInjectedParticleContainer.cpp | 135 |
1 files changed, 58 insertions, 77 deletions
diff --git a/Source/Particles/RigidInjectedParticleContainer.cpp b/Source/Particles/RigidInjectedParticleContainer.cpp index f049fdb7c..891ade76d 100644 --- a/Source/Particles/RigidInjectedParticleContainer.cpp +++ b/Source/Particles/RigidInjectedParticleContainer.cpp @@ -76,7 +76,7 @@ RigidInjectedParticleContainer::RemapParticles() // Note that the particles are already in the boosted frame. // This value is saved to advance the particles not injected yet - Cuda::ManagedDeviceVector<Real> xp, yp, zp; + Cuda::ManagedDeviceVector<ParticleReal> xp, yp, zp; for (WarpXParIter pti(*this, lev); pti.isValid(); ++pti) { @@ -136,7 +136,7 @@ RigidInjectedParticleContainer::BoostandRemapParticles() #pragma omp parallel #endif { - Cuda::ManagedDeviceVector<Real> xp, yp, zp; + Cuda::ManagedDeviceVector<ParticleReal> xp, yp, zp; for (WarpXParIter pti(*this, 0); pti.isValid(); ++pti) { @@ -207,11 +207,10 @@ RigidInjectedParticleContainer::BoostandRemapParticles() void RigidInjectedParticleContainer::PushPX(WarpXParIter& pti, - Cuda::ManagedDeviceVector<Real>& xp, - Cuda::ManagedDeviceVector<Real>& yp, - Cuda::ManagedDeviceVector<Real>& zp, - Cuda::ManagedDeviceVector<Real>& giv, - Real dt) + Cuda::ManagedDeviceVector<ParticleReal>& xp, + Cuda::ManagedDeviceVector<ParticleReal>& yp, + Cuda::ManagedDeviceVector<ParticleReal>& zp, + Real dt, DtType a_dt_type) { // This wraps the momentum and position advance so that inheritors can modify the call. @@ -221,33 +220,30 @@ RigidInjectedParticleContainer::PushPX(WarpXParIter& pti, auto& uzp = attribs[PIdx::uz]; // Save the position and momenta, making copies - Cuda::ManagedDeviceVector<Real> xp_save, yp_save, zp_save; + Cuda::ManagedDeviceVector<ParticleReal> xp_save, yp_save, zp_save; RealVector uxp_save, uyp_save, uzp_save; - Real* const AMREX_RESTRICT x = xp.dataPtr(); - Real* const AMREX_RESTRICT y = yp.dataPtr(); - Real* const AMREX_RESTRICT z = zp.dataPtr(); - Real* const AMREX_RESTRICT gi = giv.dataPtr(); - Real* const AMREX_RESTRICT ux = uxp.dataPtr(); - Real* const AMREX_RESTRICT uy = uyp.dataPtr(); - Real* const AMREX_RESTRICT uz = uzp.dataPtr(); - Real* const AMREX_RESTRICT Exp = attribs[PIdx::Ex].dataPtr(); - Real* const AMREX_RESTRICT Eyp = attribs[PIdx::Ey].dataPtr(); - Real* const AMREX_RESTRICT Ezp = attribs[PIdx::Ez].dataPtr(); - Real* const AMREX_RESTRICT Bxp = attribs[PIdx::Bx].dataPtr(); - Real* const AMREX_RESTRICT Byp = attribs[PIdx::By].dataPtr(); - Real* const AMREX_RESTRICT Bzp = attribs[PIdx::Bz].dataPtr(); + ParticleReal* const AMREX_RESTRICT x = xp.dataPtr(); + ParticleReal* const AMREX_RESTRICT y = yp.dataPtr(); + ParticleReal* const AMREX_RESTRICT z = zp.dataPtr(); + ParticleReal* const AMREX_RESTRICT ux = uxp.dataPtr(); + ParticleReal* const AMREX_RESTRICT uy = uyp.dataPtr(); + ParticleReal* const AMREX_RESTRICT uz = uzp.dataPtr(); + ParticleReal* const AMREX_RESTRICT Exp = attribs[PIdx::Ex].dataPtr(); + ParticleReal* const AMREX_RESTRICT Eyp = attribs[PIdx::Ey].dataPtr(); + ParticleReal* const AMREX_RESTRICT Ezp = attribs[PIdx::Ez].dataPtr(); + ParticleReal* const AMREX_RESTRICT Bxp = attribs[PIdx::Bx].dataPtr(); + ParticleReal* const AMREX_RESTRICT Byp = attribs[PIdx::By].dataPtr(); + ParticleReal* const AMREX_RESTRICT Bzp = attribs[PIdx::Bz].dataPtr(); if (!done_injecting_lev) { - if (!(WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags)) { - // If the old values are not already saved, create copies here. - xp_save = xp; - yp_save = yp; - zp_save = zp; - uxp_save = uxp; - uyp_save = uyp; - uzp_save = uzp; - } + // If the old values are not already saved, create copies here. + xp_save = xp; + yp_save = yp; + zp_save = zp; + uxp_save = uxp; + uyp_save = uyp; + uzp_save = uzp; // Scale the fields of particles about to cross the injection plane. // This only approximates what should be happening. The particles @@ -271,31 +267,16 @@ RigidInjectedParticleContainer::PushPX(WarpXParIter& pti, ); } - PhysicalParticleContainer::PushPX(pti, xp, yp, zp, giv, dt); + PhysicalParticleContainer::PushPX(pti, xp, yp, zp, dt, a_dt_type); if (!done_injecting_lev) { - Real* AMREX_RESTRICT x_save; - Real* AMREX_RESTRICT y_save; - Real* AMREX_RESTRICT z_save; - Real* AMREX_RESTRICT ux_save; - Real* AMREX_RESTRICT uy_save; - Real* AMREX_RESTRICT uz_save; - if (!(WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags)) { - x_save = xp_save.dataPtr(); - y_save = yp_save.dataPtr(); - z_save = zp_save.dataPtr(); - ux_save = uxp_save.dataPtr(); - uy_save = uyp_save.dataPtr(); - uz_save = uzp_save.dataPtr(); - } else { - x_save = pti.GetAttribs(particle_comps["xold"]).dataPtr(); - y_save = pti.GetAttribs(particle_comps["yold"]).dataPtr(); - z_save = pti.GetAttribs(particle_comps["zold"]).dataPtr(); - ux_save = pti.GetAttribs(particle_comps["uxold"]).dataPtr(); - uy_save = pti.GetAttribs(particle_comps["uyold"]).dataPtr(); - uz_save = pti.GetAttribs(particle_comps["uzold"]).dataPtr(); - } + ParticleReal* AMREX_RESTRICT x_save = xp_save.dataPtr(); + ParticleReal* AMREX_RESTRICT y_save = yp_save.dataPtr(); + ParticleReal* AMREX_RESTRICT z_save = zp_save.dataPtr(); + ParticleReal* AMREX_RESTRICT ux_save = uxp_save.dataPtr(); + ParticleReal* AMREX_RESTRICT uy_save = uyp_save.dataPtr(); + ParticleReal* AMREX_RESTRICT uz_save = uzp_save.dataPtr(); // Undo the push for particles not injected yet. // The zp are advanced a fixed amount. @@ -309,14 +290,14 @@ RigidInjectedParticleContainer::PushPX(WarpXParIter& pti, ux[i] = ux_save[i]; uy[i] = uy_save[i]; uz[i] = uz_save[i]; - gi[i] = 1./std::sqrt(1. + (ux[i]*ux[i] + uy[i]*uy[i] + uz[i]*uz[i])*inv_csq); x[i] = x_save[i]; y[i] = y_save[i]; if (rigid) { z[i] = z_save[i] + dt*vz_ave_boosted; } else { - z[i] = z_save[i] + dt*uz[i]*gi[i]; + const Real gi = 1./std::sqrt(1. + (ux[i]*ux[i] + uy[i]*uy[i] + uz[i]*uz[i])*inv_csq); + z[i] = z_save[i] + dt*uz[i]*gi; } } } @@ -333,7 +314,7 @@ RigidInjectedParticleContainer::Evolve (int lev, MultiFab* rho, MultiFab* crho, const MultiFab* cEx, const MultiFab* cEy, const MultiFab* cEz, const MultiFab* cBx, const MultiFab* cBy, const MultiFab* cBz, - Real t, Real dt) + Real t, Real dt, DtType a_dt_type) { // Update location of injection plane in the boosted frame @@ -347,7 +328,9 @@ RigidInjectedParticleContainer::Evolve (int lev, // particles have crossed the inject plane. const Real* plo = Geom(lev).ProbLo(); const Real* phi = Geom(lev).ProbHi(); - done_injecting[lev] = (zinject_plane_levels[lev] < plo[2] || zinject_plane_levels[lev] > phi[2]); + const int zdir = AMREX_SPACEDIM-1; + done_injecting[lev] = ((zinject_plane_levels[lev] < plo[zdir] && WarpX::moving_window_v + WarpX::beta_boost*PhysConst::c >= 0.) || + (zinject_plane_levels[lev] > phi[zdir] && WarpX::moving_window_v + WarpX::beta_boost*PhysConst::c <= 0.)); done_injecting_lev = done_injecting[lev]; PhysicalParticleContainer::Evolve (lev, @@ -358,7 +341,7 @@ RigidInjectedParticleContainer::Evolve (int lev, rho, crho, cEx, cEy, cEz, cBx, cBy, cBz, - t, dt); + t, dt, a_dt_type); } void @@ -414,8 +397,6 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt, Byp.assign(np,WarpX::B_external[1]); Bzp.assign(np,WarpX::B_external[2]); - m_giv[thread_num].resize(np); - // // copy data from particle container to temp arrays // @@ -424,7 +405,8 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt, int e_is_nodal = Ex.is_nodal() and Ey.is_nodal() and Ez.is_nodal(); FieldGather(pti, Exp, Eyp, Ezp, Bxp, Byp, Bzp, &exfab, &eyfab, &ezfab, &bxfab, &byfab, &bzfab, - Ex.nGrow(), e_is_nodal, 0, np, thread_num, lev, lev); + Ex.nGrow(), e_is_nodal, + 0, np, thread_num, lev, lev); // Save the position and momenta, making copies auto uxp_save = uxp; @@ -433,17 +415,16 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt, // This wraps the momentum advance so that inheritors can modify the call. // Extract pointers to the different particle quantities - const Real* const AMREX_RESTRICT zp = m_zp[thread_num].dataPtr(); - Real* const AMREX_RESTRICT gi = m_giv[thread_num].dataPtr(); - Real* const AMREX_RESTRICT uxpp = uxp.dataPtr(); - Real* const AMREX_RESTRICT uypp = uyp.dataPtr(); - Real* const AMREX_RESTRICT uzpp = uzp.dataPtr(); - const Real* const AMREX_RESTRICT Expp = Exp.dataPtr(); - const Real* const AMREX_RESTRICT Eypp = Eyp.dataPtr(); - const Real* const AMREX_RESTRICT Ezpp = Ezp.dataPtr(); - const Real* const AMREX_RESTRICT Bxpp = Bxp.dataPtr(); - const Real* const AMREX_RESTRICT Bypp = Byp.dataPtr(); - const Real* const AMREX_RESTRICT Bzpp = Bzp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT zp = m_zp[thread_num].dataPtr(); + ParticleReal* const AMREX_RESTRICT uxpp = uxp.dataPtr(); + ParticleReal* const AMREX_RESTRICT uypp = uyp.dataPtr(); + ParticleReal* const AMREX_RESTRICT uzpp = uzp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT Expp = Exp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT Eypp = Eyp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT Ezpp = Ezp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT Bxpp = Bxp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT Bypp = Byp.dataPtr(); + const ParticleReal* const AMREX_RESTRICT Bzpp = Bzp.dataPtr(); // Loop over the particles and update their momentum const Real q = this->charge; @@ -451,14 +432,14 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt, if (WarpX::particle_pusher_algo == ParticlePusherAlgo::Boris){ amrex::ParallelFor( pti.numParticles(), [=] AMREX_GPU_DEVICE (long i) { - UpdateMomentumBoris( uxpp[i], uypp[i], uzpp[i], gi[i], + UpdateMomentumBoris( uxpp[i], uypp[i], uzpp[i], Expp[i], Eypp[i], Ezpp[i], Bxpp[i], Bypp[i], Bzpp[i], q, m, dt); } ); } else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::Vay) { amrex::ParallelFor( pti.numParticles(), [=] AMREX_GPU_DEVICE (long i) { - UpdateMomentumVay( uxpp[i], uypp[i], uzpp[i], gi[i], + UpdateMomentumVay( uxpp[i], uypp[i], uzpp[i], Expp[i], Eypp[i], Ezpp[i], Bxpp[i], Bypp[i], Bzpp[i], q, m, dt); } ); @@ -469,10 +450,10 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt, // Undo the push for particles not injected yet. // It is assumed that PushP will only be called on the first and last steps // and that no particles will cross zinject_plane. - const Real* const AMREX_RESTRICT ux_save = uxp_save.dataPtr(); - const Real* const AMREX_RESTRICT uy_save = uyp_save.dataPtr(); - const Real* const AMREX_RESTRICT uz_save = uzp_save.dataPtr(); - const Real zz = zinject_plane_levels[lev]; + const ParticleReal* const AMREX_RESTRICT ux_save = uxp_save.dataPtr(); + const ParticleReal* const AMREX_RESTRICT uy_save = uyp_save.dataPtr(); + const ParticleReal* const AMREX_RESTRICT uz_save = uzp_save.dataPtr(); + const ParticleReal zz = zinject_plane_levels[lev]; amrex::ParallelFor( pti.numParticles(), [=] AMREX_GPU_DEVICE (long i) { if (zp[i] <= zz) { |