From 641ef5275815401abc8cffa17e28d313403cfa34 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Tue, 4 Jun 2024 17:10:42 -0700 Subject: [PATCH] Clean up of ParticleExtrema --- .../ReducedDiags/ParticleExtrema.cpp | 237 +++++------------- 1 file changed, 64 insertions(+), 173 deletions(-) diff --git a/Source/Diagnostics/ReducedDiags/ParticleExtrema.cpp b/Source/Diagnostics/ReducedDiags/ParticleExtrema.cpp index 0cc5429be7a..d93e570a2e6 100644 --- a/Source/Diagnostics/ReducedDiags/ParticleExtrema.cpp +++ b/Source/Diagnostics/ReducedDiags/ParticleExtrema.cpp @@ -51,7 +51,7 @@ #include #include -using namespace amrex; +using namespace amrex::literals; using namespace warpx::fields; // constructor @@ -59,7 +59,7 @@ ParticleExtrema::ParticleExtrema (const std::string& rd_name) : ReducedDiags{rd_name} { // read species name - const ParmParse pp_rd_name(rd_name); + const amrex::ParmParse pp_rd_name(rd_name); pp_rd_name.get("species",m_species_name); // get WarpX class object @@ -122,7 +122,7 @@ ParticleExtrema::ParticleExtrema (const std::string& rd_name) m_data.resize(all_diag_names.size()); - if (ParallelDescriptor::IOProcessor()) + if (amrex::ParallelDescriptor::IOProcessor()) { if ( m_write_header ) { @@ -165,16 +165,7 @@ void ParticleExtrema::ComputeDiags (int step) const auto species_names = mypc.GetSpeciesNames(); // inverse of speed of light squared - Real constexpr inv_c2 = 1.0_rt / (PhysConst::c * PhysConst::c); - - // If 2D-XZ, p.pos(1) is z, rather than p.pos(2). -#if (defined WARPX_DIM_3D) - int const index_z = 2; -#elif (defined WARPX_DIM_XZ || defined WARPX_DIM_RZ) - int const index_z = 1; -#elif (defined WARPX_DIM_1D_Z) - int const index_z = 0; -#endif + amrex::Real constexpr inv_c2 = 1.0_rt / (PhysConst::c * PhysConst::c); // loop over species for (int i_s = 0; i_s < nSpecies; ++i_s) @@ -193,172 +184,72 @@ void ParticleExtrema::ComputeDiags (int step) } using PType = typename WarpXParticleContainer::SuperParticleType; - - // xmin -#if (defined WARPX_DIM_RZ) - Real xmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(0)*std::cos(p.rdata(PIdx::theta)); }); -#elif (defined WARPX_DIM_1D_Z) - Real xmin = 0.0_rt; -#else - Real xmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(0); }); -#endif - - // xmax -#if (defined WARPX_DIM_RZ) - Real xmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(0)*std::cos(p.rdata(PIdx::theta)); }); -#elif (defined WARPX_DIM_1D_Z) - Real xmax = 0.0_rt; -#else - Real xmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(0); }); -#endif - - // ymin -#if (defined WARPX_DIM_RZ) - Real ymin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(0)*std::sin(p.rdata(PIdx::theta)); }); -#elif (defined WARPX_DIM_XZ || WARPX_DIM_1D_Z) - Real ymin = 0.0_rt; -#else - Real ymin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(1); }); -#endif - - // ymax -#if (defined WARPX_DIM_RZ) - Real ymax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(0)*std::sin(p.rdata(PIdx::theta)); }); -#elif (defined WARPX_DIM_XZ || WARPX_DIM_1D_Z) - Real ymax = 0.0_rt; -#else - Real ymax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(1); }); -#endif - - // zmin - Real zmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(index_z); }); - - // zmax - Real zmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.pos(index_z); }); - - // uxmin - Real uxmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::ux); }); - - // uxmax - Real uxmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::ux); }); - - // uymin - Real uymin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::uy); }); - - // uymax - Real uymax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::uy); }); - - // uzmin - Real uzmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::uz); }); - - // uzmax - Real uzmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::uz); }); - - // gmin - Real gmin = 0.0_rt; - if ( is_photon ) { - gmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) + using OpMin = amrex::ReduceOpMin; + using OpMax = amrex::ReduceOpMax; + + amrex::ReduceOps reduce_ops; + auto posminmax = amrex::ParticleReduce>( + myspc, + [=] AMREX_GPU_DEVICE(const PType& p) noexcept -> amrex::GpuTuple { - const Real ux = p.rdata(PIdx::ux); - const Real uy = p.rdata(PIdx::uy); - const Real uz = p.rdata(PIdx::uz); - const Real us = ux*ux + uy*uy + uz*uz; - return std::sqrt(us*inv_c2); - }); - } else { - gmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) + amrex::ParticleReal x, y, z; + get_particle_position(p, x, y, z); + amrex::Real const w = p.rdata(PIdx::w); + return {w, x, y, z, w, x, y, z}; + }, + reduce_ops); + + amrex::Real wmin = amrex::get<0>(posminmax); + amrex::Real xmin = amrex::get<1>(posminmax); + amrex::Real ymin = amrex::get<2>(posminmax); + amrex::Real zmin = amrex::get<3>(posminmax); + amrex::Real wmax = amrex::get<4>(posminmax); + amrex::Real xmax = amrex::get<5>(posminmax); + amrex::Real ymax = amrex::get<6>(posminmax); + amrex::Real zmax = amrex::get<7>(posminmax); + + amrex::Real gfactor = (is_photon ? 0._rt : 1._rt); + auto uminmax = amrex::ParticleReduce>( + myspc, + [=] AMREX_GPU_DEVICE(const PType& p) noexcept -> amrex::GpuTuple { - const Real ux = p.rdata(PIdx::ux); - const Real uy = p.rdata(PIdx::uy); - const Real uz = p.rdata(PIdx::uz); - const Real us = ux*ux + uy*uy + uz*uz; - return std::sqrt(1.0_rt + us*inv_c2); - }); - } - - // gmax - Real gmax = 0.0_rt; - if ( is_photon ) { - gmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { - const Real ux = p.rdata(PIdx::ux); - const Real uy = p.rdata(PIdx::uy); - const Real uz = p.rdata(PIdx::uz); - const Real us = ux*ux + uy*uy + uz*uz; - return std::sqrt(us*inv_c2); - }); - } else { - gmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { - const Real ux = p.rdata(PIdx::ux); - const Real uy = p.rdata(PIdx::uy); - const Real uz = p.rdata(PIdx::uz); - const Real us = ux*ux + uy*uy + uz*uz; - return std::sqrt(1.0_rt + us*inv_c2); - }); - } - - // wmin - Real wmin = ReduceMin( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::w); }); - - // wmax - Real wmax = ReduceMax( myspc, - [=] AMREX_GPU_HOST_DEVICE (const PType& p) - { return p.rdata(PIdx::w); }); - - ParallelDescriptor::ReduceRealMin({xmin,ymin,zmin,uxmin,uymin,uzmin,gmin,wmin}); - ParallelDescriptor::ReduceRealMax({xmax,ymax,zmax,uxmax,uymax,uzmax,gmax,wmax}); + amrex::Real const ux = p.rdata(PIdx::ux); + amrex::Real const uy = p.rdata(PIdx::uy); + amrex::Real const uz = p.rdata(PIdx::uz); + amrex::Real const g = std::sqrt(gfactor + (ux*ux + uy*uy + uz*uz)*inv_c2); + return {g, ux, uy, uz, g, ux, uy, uz}; + }, + reduce_ops); + + amrex::Real gmin = amrex::get<0>(uminmax); + amrex::Real uxmin = amrex::get<1>(uminmax); + amrex::Real uymin = amrex::get<2>(uminmax); + amrex::Real uzmin = amrex::get<3>(uminmax); + amrex::Real gmax = amrex::get<4>(uminmax); + amrex::Real uxmax = amrex::get<5>(uminmax); + amrex::Real uymax = amrex::get<6>(uminmax); + amrex::Real uzmax = amrex::get<7>(uminmax); + + amrex::ParallelDescriptor::ReduceRealMin({xmin,ymin,zmin,uxmin,uymin,uzmin,gmin,wmin}); + amrex::ParallelDescriptor::ReduceRealMax({xmax,ymax,zmax,uxmax,uymax,uzmax,gmax,wmax}); #if (defined WARPX_QED) // get number of level (int) const auto level_number = WarpX::GetInstance().finestLevel(); // compute chimin and chimax - Real chimin_f = 0.0_rt; - Real chimax_f = 0.0_rt; + amrex::Real chimin_f = 0.0_rt; + amrex::Real chimax_f = 0.0_rt; if (myspc.DoQED()) { // declare chi arrays - std::vector chimin, chimax; + std::vector chimin, chimax; chimin.resize(level_number+1,0.0_rt); chimax.resize(level_number+1,0.0_rt); @@ -383,8 +274,8 @@ void ParticleExtrema::ComputeDiags (int step) const MultiFab & Bz = warpx.getField(FieldType::Bfield_aux, lev,2); // declare reduce_op - ReduceOps reduce_op; - ReduceData reduce_data(reduce_op); + amrex::ReduceOps reduce_op; + amrex::ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; // Loop over boxes @@ -429,7 +320,7 @@ void ParticleExtrema::ComputeDiags (int step) [=] AMREX_GPU_DEVICE (int i) -> ReduceTuple { // get external fields - ParticleReal xp, yp, zp; + amrex::ParticleReal xp, yp, zp; GetPosition(i, xp, yp, zp); amrex::ParticleReal ex = Ex_external_particle; amrex::ParticleReal ey = Ey_external_particle; @@ -449,7 +340,7 @@ void ParticleExtrema::ComputeDiags (int step) dx_arr, xyzmin_arr, lo, n_rz_azimuthal_modes, nox, galerkin_interpolation); // compute chi - Real chi = 0.0_rt; + amrex::Real chi = 0.0_rt; if ( is_photon ) { chi = QedUtils::chi_photon(ux[i]*m, uy[i]*m, uz[i]*m, ex, ey, ez, bx, by, bz); @@ -466,8 +357,8 @@ void ParticleExtrema::ComputeDiags (int step) } chimin_f = *std::min_element(chimin.begin(), chimin.end()); chimax_f = *std::max_element(chimax.begin(), chimax.end()); - ParallelDescriptor::ReduceRealMin(chimin_f, ParallelDescriptor::IOProcessorNumber()); - ParallelDescriptor::ReduceRealMax(chimax_f, ParallelDescriptor::IOProcessorNumber()); + amrex::ParallelDescriptor::ReduceRealMin(chimin_f, amrex::ParallelDescriptor::IOProcessorNumber()); + amrex::ParallelDescriptor::ReduceRealMax(chimax_f, amrex::ParallelDescriptor::IOProcessorNumber()); } #endif