From 2c2f48ca24442c2efa62d3f64d3642e4d6335c58 Mon Sep 17 00:00:00 2001 From: Axel Huebl Date: Wed, 2 Oct 2024 08:53:11 -0700 Subject: [PATCH 1/6] MultiFabRegister: `throw` in get Throw a runtime exception instead of returning a `nullptr` if a field is requested via the getter. --- Source/ablastr/fields/MultiFabRegister.cpp | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/Source/ablastr/fields/MultiFabRegister.cpp b/Source/ablastr/fields/MultiFabRegister.cpp index 2c384a90089..722fa0ce106 100644 --- a/Source/ablastr/fields/MultiFabRegister.cpp +++ b/Source/ablastr/fields/MultiFabRegister.cpp @@ -350,8 +350,7 @@ namespace ablastr::fields ) { if (m_mf_register.count(internal_name) == 0) { - // FIXME: temporary, throw a std::runtime_error - // throw std::runtime_error("MultiFabRegister::get name does not exist in register: " + key); + throw std::runtime_error("MultiFabRegister::get name does not exist in register: " + internal_name); return nullptr; } amrex::MultiFab & mf = m_mf_register.at(internal_name).m_mf; @@ -365,8 +364,7 @@ namespace ablastr::fields ) const { if (m_mf_register.count(internal_name) == 0) { - // FIXME: temporary, throw a std::runtime_error - // throw std::runtime_error("MultiFabRegister::get name does not exist in register: " + internal_name); + throw std::runtime_error("MultiFabRegister::get name does not exist in register: " + internal_name); return nullptr; } amrex::MultiFab const & mf = m_mf_register.at(internal_name).m_mf; From 6a5ea04410bcb62676d754e43d4ed3072f05fc7d Mon Sep 17 00:00:00 2001 From: Edoardo Zoni Date: Fri, 18 Oct 2024 17:51:25 -0700 Subject: [PATCH 2/6] Add `skip_level_0` argument for coarse-patch fields --- Source/Evolve/WarpXEvolve.cpp | 32 +++++++------ .../LabFrameExplicitES.cpp | 4 +- Source/FieldSolver/WarpXPushFieldsEM.cpp | 10 ++-- .../FieldSolver/WarpXPushFieldsHybridPIC.cpp | 4 +- Source/Parallelization/WarpXComm.cpp | 21 +++++--- Source/Utils/WarpXMovingWindow.cpp | 6 ++- Source/WarpX.cpp | 1 + Source/ablastr/fields/MultiFabRegister.H | 36 +++++++++----- Source/ablastr/fields/MultiFabRegister.cpp | 48 +++++++++++++++---- 9 files changed, 114 insertions(+), 48 deletions(-) diff --git a/Source/Evolve/WarpXEvolve.cpp b/Source/Evolve/WarpXEvolve.cpp index a685afd28e7..01920169d67 100644 --- a/Source/Evolve/WarpXEvolve.cpp +++ b/Source/Evolve/WarpXEvolve.cpp @@ -663,6 +663,8 @@ WarpX::OneStep_multiJ (const amrex::Real cur_time) using warpx::fields::FieldType; + bool const skip_level_0 = true; + const int rho_mid = spectral_solver_fp[0]->m_spectral_index.rho_mid; const int rho_new = spectral_solver_fp[0]->m_spectral_index.rho_new; @@ -796,8 +798,8 @@ WarpX::OneStep_multiJ (const amrex::Real cur_time) PSATDBackwardTransformEBavg( m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_fp, finest_level), m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_fp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level) + m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level, skip_level_0), + m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level, skip_level_0) ); } @@ -868,11 +870,13 @@ WarpX::OneStep_sub1 (Real cur_time) using warpx::fields::FieldType; + bool const skip_level_0 = true; + // i) Push particles and fields on the fine patch (first fine step) PushParticlesandDeposit(fine_lev, cur_time, DtType::FirstHalf); RestrictCurrentFromFineToCoarsePatch( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level), fine_lev); + m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level, skip_level_0), fine_lev); RestrictRhoFromFineToCoarsePatch(fine_lev); if (use_filter) { ApplyFilterJ( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), fine_lev); @@ -883,7 +887,7 @@ WarpX::OneStep_sub1 (Real cur_time) ApplyFilterandSumBoundaryRho( m_fields.get_mr_levels(FieldType::rho_fp, finest_level), - m_fields.get_mr_levels(FieldType::rho_cp, finest_level), + m_fields.get_mr_levels(FieldType::rho_cp, finest_level, skip_level_0), fine_lev, PatchType::fine, 0, 2*ncomps); EvolveB(fine_lev, PatchType::fine, 0.5_rt*dt[fine_lev], DtType::FirstHalf); @@ -914,12 +918,12 @@ WarpX::OneStep_sub1 (Real cur_time) StoreCurrent(coarse_lev); AddCurrentFromFineLevelandSumBoundary( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::current_buf, finest_level), coarse_lev); + m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level, skip_level_0), + m_fields.get_mr_levels_alldirs(FieldType::current_buf, finest_level, skip_level_0), coarse_lev); AddRhoFromFineLevelandSumBoundary( m_fields.get_mr_levels(FieldType::rho_fp, finest_level), - m_fields.get_mr_levels(FieldType::rho_cp, finest_level), - m_fields.get_mr_levels(FieldType::rho_buf, finest_level), + m_fields.get_mr_levels(FieldType::rho_cp, finest_level, skip_level_0), + m_fields.get_mr_levels(FieldType::rho_buf, finest_level, skip_level_0), coarse_lev, 0, ncomps); EvolveB(fine_lev, PatchType::coarse, dt[fine_lev], DtType::FirstHalf); @@ -950,7 +954,7 @@ WarpX::OneStep_sub1 (Real cur_time) PushParticlesandDeposit(fine_lev, cur_time + dt[fine_lev], DtType::SecondHalf); RestrictCurrentFromFineToCoarsePatch( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level), fine_lev); + m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level, skip_level_0), fine_lev); RestrictRhoFromFineToCoarsePatch(fine_lev); if (use_filter) { ApplyFilterJ( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), fine_lev); @@ -958,7 +962,7 @@ WarpX::OneStep_sub1 (Real cur_time) SumBoundaryJ( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), fine_lev, Geom(fine_lev).periodicity()); ApplyFilterandSumBoundaryRho( m_fields.get_mr_levels(FieldType::rho_fp, finest_level), - m_fields.get_mr_levels(FieldType::rho_cp, finest_level), + m_fields.get_mr_levels(FieldType::rho_cp, finest_level, skip_level_0), fine_lev, PatchType::fine, 0, ncomps); EvolveB(fine_lev, PatchType::fine, 0.5_rt*dt[fine_lev], DtType::FirstHalf); @@ -988,13 +992,13 @@ WarpX::OneStep_sub1 (Real cur_time) RestoreCurrent(coarse_lev); AddCurrentFromFineLevelandSumBoundary( m_fields.get_mr_levels_alldirs(FieldType::current_fp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level), - m_fields.get_mr_levels_alldirs(FieldType::current_buf, finest_level), + m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level, skip_level_0), + m_fields.get_mr_levels_alldirs(FieldType::current_buf, finest_level, skip_level_0), coarse_lev); AddRhoFromFineLevelandSumBoundary( m_fields.get_mr_levels(FieldType::rho_fp, finest_level), - m_fields.get_mr_levels(FieldType::rho_cp, finest_level), - m_fields.get_mr_levels(FieldType::rho_buf, finest_level), + m_fields.get_mr_levels(FieldType::rho_cp, finest_level, skip_level_0), + m_fields.get_mr_levels(FieldType::rho_buf, finest_level, skip_level_0), coarse_lev, ncomps, ncomps); EvolveE(fine_lev, PatchType::coarse, dt[fine_lev]); diff --git a/Source/FieldSolver/ElectrostaticSolvers/LabFrameExplicitES.cpp b/Source/FieldSolver/ElectrostaticSolvers/LabFrameExplicitES.cpp index e973ae66975..c5f8353167c 100644 --- a/Source/FieldSolver/ElectrostaticSolvers/LabFrameExplicitES.cpp +++ b/Source/FieldSolver/ElectrostaticSolvers/LabFrameExplicitES.cpp @@ -31,8 +31,10 @@ void LabFrameExplicitES::ComputeSpaceChargeField ( using ablastr::fields::MultiLevelVectorField; using warpx::fields::FieldType; + bool const skip_level_0 = true; + const MultiLevelScalarField rho_fp = fields.get_mr_levels(FieldType::rho_fp, max_level); - const MultiLevelScalarField rho_cp = fields.get_mr_levels(FieldType::rho_cp, max_level); + const MultiLevelScalarField rho_cp = fields.get_mr_levels(FieldType::rho_cp, max_level, skip_level_0); const MultiLevelScalarField phi_fp = fields.get_mr_levels(FieldType::phi_fp, max_level); const MultiLevelVectorField Efield_fp = fields.get_mr_levels_alldirs(FieldType::Efield_fp, max_level); diff --git a/Source/FieldSolver/WarpXPushFieldsEM.cpp b/Source/FieldSolver/WarpXPushFieldsEM.cpp index fd786dc65ba..1ac456cd35a 100644 --- a/Source/FieldSolver/WarpXPushFieldsEM.cpp +++ b/Source/FieldSolver/WarpXPushFieldsEM.cpp @@ -721,6 +721,8 @@ WarpX::PushPSATD () "PushFieldsEM: PSATD solver selected but not built"); #else + bool const skip_level_0 = true; + const int rho_old = spectral_solver_fp[0]->m_spectral_index.rho_old; const int rho_new = spectral_solver_fp[0]->m_spectral_index.rho_new; @@ -852,8 +854,8 @@ WarpX::PushPSATD () if (WarpX::fft_do_time_averaging) { auto Efield_avg_fp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_fp, finest_level); auto Bfield_avg_fp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_fp, finest_level); - auto Efield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level); - auto Bfield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level); + auto Efield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level, skip_level_0); + auto Bfield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level, skip_level_0); PSATDBackwardTransformEBavg(Efield_avg_fp, Bfield_avg_fp, Efield_avg_cp, Bfield_avg_cp); } if (WarpX::do_dive_cleaning) { PSATDBackwardTransformF(); } @@ -1098,6 +1100,8 @@ WarpX::EvolveG (int lev, PatchType patch_type, amrex::Real a_dt, DtType /*a_dt_t WARPX_PROFILE("WarpX::EvolveG()"); + bool const skip_level_0 = true; + // Evolve G field in regular cells if (patch_type == PatchType::fine) { @@ -1108,7 +1112,7 @@ WarpX::EvolveG (int lev, PatchType patch_type, amrex::Real a_dt, DtType /*a_dt_t } else // coarse patch { - ablastr::fields::MultiLevelVectorField const& Bfield_cp_new = m_fields.get_mr_levels_alldirs(FieldType::Bfield_cp, finest_level); + ablastr::fields::MultiLevelVectorField const& Bfield_cp_new = m_fields.get_mr_levels_alldirs(FieldType::Bfield_cp, finest_level, skip_level_0); m_fdtd_solver_cp[lev]->EvolveG( m_fields.get(FieldType::G_cp, lev), Bfield_cp_new[lev], a_dt); diff --git a/Source/FieldSolver/WarpXPushFieldsHybridPIC.cpp b/Source/FieldSolver/WarpXPushFieldsHybridPIC.cpp index 5220419f822..6e1ac84d4b2 100644 --- a/Source/FieldSolver/WarpXPushFieldsHybridPIC.cpp +++ b/Source/FieldSolver/WarpXPushFieldsHybridPIC.cpp @@ -195,11 +195,13 @@ void WarpX::HybridPICDepositInitialRhoAndJ () { using warpx::fields::FieldType; + bool const skip_level_0 = true; + ablastr::fields::MultiLevelScalarField rho_fp_temp = m_fields.get_mr_levels(FieldType::hybrid_rho_fp_temp, finest_level); ablastr::fields::MultiLevelVectorField current_fp_temp = m_fields.get_mr_levels_alldirs(FieldType::hybrid_current_fp_temp, finest_level); mypc->DepositCharge(rho_fp_temp, 0._rt); mypc->DepositCurrent(current_fp_temp, dt[0], 0._rt); - SyncRho(rho_fp_temp, m_fields.get_mr_levels(FieldType::rho_cp, finest_level), m_fields.get_mr_levels(FieldType::rho_buf, finest_level)); + SyncRho(rho_fp_temp, m_fields.get_mr_levels(FieldType::rho_cp, finest_level, skip_level_0), m_fields.get_mr_levels(FieldType::rho_buf, finest_level, skip_level_0)); SyncCurrent("hybrid_current_fp_temp"); for (int lev=0; lev <= finest_level; ++lev) { // SyncCurrent does not include a call to FillBoundary, but it is needed diff --git a/Source/Parallelization/WarpXComm.cpp b/Source/Parallelization/WarpXComm.cpp index d64632d964a..3bf0b887f6f 100644 --- a/Source/Parallelization/WarpXComm.cpp +++ b/Source/Parallelization/WarpXComm.cpp @@ -836,6 +836,8 @@ WarpX::FillBoundaryE_avg(int lev, IntVect ng) void WarpX::FillBoundaryE_avg (int lev, PatchType patch_type, IntVect ng) { + bool const skip_level_0 = true; + if (patch_type == PatchType::fine) { if (do_pml && pml[lev]->ok()) @@ -865,7 +867,7 @@ WarpX::FillBoundaryE_avg (int lev, PatchType patch_type, IntVect ng) WARPX_ABORT_WITH_MESSAGE("Averaged Galilean PSATD with PML is not yet implemented"); } - ablastr::fields::MultiLevelVectorField Efield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level); + ablastr::fields::MultiLevelVectorField Efield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level, skip_level_0); const amrex::Periodicity& cperiod = Geom(lev-1).periodicity(); if ( safe_guard_cells ) { @@ -896,6 +898,8 @@ WarpX::FillBoundaryB_avg (int lev, PatchType patch_type, IntVect ng) { using ablastr::fields::Direction; + bool const skip_level_0 = true; + if (patch_type == PatchType::fine) { if (do_pml && pml[lev]->ok()) @@ -925,7 +929,7 @@ WarpX::FillBoundaryB_avg (int lev, PatchType patch_type, IntVect ng) WARPX_ABORT_WITH_MESSAGE("Averaged Galilean PSATD with PML is not yet implemented"); } - ablastr::fields::MultiLevelVectorField Bfield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level); + ablastr::fields::MultiLevelVectorField Bfield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level, skip_level_0); const amrex::Periodicity& cperiod = Geom(lev-1).periodicity(); if ( safe_guard_cells ){ @@ -1077,6 +1081,8 @@ WarpX::SyncCurrent (const std::string& current_fp_string) WARPX_PROFILE("WarpX::SyncCurrent()"); + bool const skip_level_0 = true; + ablastr::fields::MultiLevelVectorField const& J_fp = m_fields.get_mr_levels_alldirs(current_fp_string, finest_level); // If warpx.do_current_centering = 1, center currents from nodal grid to staggered grid @@ -1192,7 +1198,7 @@ WarpX::SyncCurrent (const std::string& current_fp_string) } }); // Now it's safe to apply filter and sumboundary on J_cp - ablastr::fields::MultiLevelVectorField const& J_cp = m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level); + ablastr::fields::MultiLevelVectorField const& J_cp = m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level, skip_level_0); if (use_filter) { ApplyFilterJ(J_cp, lev+1, idim); @@ -1207,14 +1213,14 @@ WarpX::SyncCurrent (const std::string& current_fp_string) // filtering depends on the level. This is also done before any // same-level communication because it's easier this way to // avoid double counting. - ablastr::fields::MultiLevelVectorField const& J_cp = m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level); + ablastr::fields::MultiLevelVectorField const& J_cp = m_fields.get_mr_levels_alldirs(FieldType::current_cp, finest_level, skip_level_0); J_cp[lev][Direction{idim}]->setVal(0.0); ablastr::coarsen::average::Coarsen(*J_cp[lev][Direction{idim}], *J_fp[lev][Direction{idim}], refRatio(lev-1)); if (m_fields.has(FieldType::current_buf, Direction{idim}, lev)) { - ablastr::fields::MultiLevelVectorField const& J_buffer = m_fields.get_mr_levels_alldirs(FieldType::current_buf, finest_level); + ablastr::fields::MultiLevelVectorField const& J_buffer = m_fields.get_mr_levels_alldirs(FieldType::current_buf, finest_level, skip_level_0); IntVect const& ng = J_cp[lev][Direction{idim}]->nGrowVect(); AMREX_ASSERT(ng.allLE(J_buffer[lev][Direction{idim}]->nGrowVect())); @@ -1241,14 +1247,15 @@ WarpX::SyncCurrent (const std::string& current_fp_string) void WarpX::SyncRho () { + bool const skip_level_0 = true; const ablastr::fields::MultiLevelScalarField rho_fp = m_fields.has(FieldType::rho_fp, 0) ? m_fields.get_mr_levels(FieldType::rho_fp, finest_level) : ablastr::fields::MultiLevelScalarField{static_cast(finest_level+1)}; const ablastr::fields::MultiLevelScalarField rho_cp = m_fields.has(FieldType::rho_cp, 1) ? - m_fields.get_mr_levels(FieldType::rho_cp, finest_level) : + m_fields.get_mr_levels(FieldType::rho_cp, finest_level, skip_level_0) : ablastr::fields::MultiLevelScalarField{static_cast(finest_level+1)}; const ablastr::fields::MultiLevelScalarField rho_buf = m_fields.has(FieldType::rho_buf, 1) ? - m_fields.get_mr_levels(FieldType::rho_buf, finest_level) : + m_fields.get_mr_levels(FieldType::rho_buf, finest_level, skip_level_0) : ablastr::fields::MultiLevelScalarField{static_cast(finest_level+1)}; SyncRho(rho_fp, rho_cp, rho_buf); diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index d5cebd69254..c119954ad13 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -143,6 +143,8 @@ WarpX::MoveWindow (const int step, bool move_j) using ablastr::fields::Direction; using warpx::fields::FieldType; + bool const skip_level_0 = true; + if (step == start_moving_window_step) { amrex::Print() << Utils::TextMsg::Info("Starting moving window"); } @@ -276,8 +278,8 @@ WarpX::MoveWindow (const int step, bool move_j) shiftMF(*m_fields.get(FieldType::Bfield_aux, Direction{dim}, lev), geom[lev], num_shift, dir, lev, do_update_cost); shiftMF(*m_fields.get(FieldType::Efield_aux, Direction{dim}, lev), geom[lev], num_shift, dir, lev, do_update_cost); if (fft_do_time_averaging) { - ablastr::fields::MultiLevelVectorField Efield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level); - ablastr::fields::MultiLevelVectorField Bfield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level); + ablastr::fields::MultiLevelVectorField Efield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Efield_avg_cp, finest_level, skip_level_0); + ablastr::fields::MultiLevelVectorField Bfield_avg_cp = m_fields.get_mr_levels_alldirs(FieldType::Bfield_avg_cp, finest_level, skip_level_0); shiftMF(*Bfield_avg_cp[lev][dim], geom[lev-1], num_shift_crse, dir, lev, do_update_cost, m_p_ext_field_params->B_external_grid[dim], use_Bparser, Bfield_parser); shiftMF(*Efield_avg_cp[lev][dim], geom[lev-1], num_shift_crse, dir, lev, do_update_cost, diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index d1e3108e32a..c996f847c57 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -2593,6 +2593,7 @@ WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm } // Create the MultiFabs for the current + amrex::Print() << "allocate current_cp" << std::endl; m_fields.alloc_init(FieldType::current_cp, Direction{0}, lev, amrex::convert(cba, jx_nodal_flag), dm, ncomps, ngJ, 0.0_rt); m_fields.alloc_init(FieldType::current_cp, Direction{1}, lev, amrex::convert(cba, jy_nodal_flag), dm, ncomps, ngJ, 0.0_rt); m_fields.alloc_init(FieldType::current_cp, Direction{2}, lev, amrex::convert(cba, jz_nodal_flag), dm, ncomps, ngJ, 0.0_rt); diff --git a/Source/ablastr/fields/MultiFabRegister.H b/Source/ablastr/fields/MultiFabRegister.H index 21df20c1678..3908c4eba0b 100644 --- a/Source/ablastr/fields/MultiFabRegister.H +++ b/Source/ablastr/fields/MultiFabRegister.H @@ -479,24 +479,28 @@ namespace ablastr::fields [[nodiscard]] MultiLevelScalarField get_mr_levels ( T name, - int finest_level + int finest_level, + bool skip_level_0=false ) { return internal_get_mr_levels( getExtractedName(name), - finest_level + finest_level, + skip_level_0 ); } template [[nodiscard]] ConstMultiLevelScalarField get_mr_levels ( T name, - int finest_level + int finest_level, + bool skip_level_0=false ) const { return internal_get_mr_levels( getExtractedName(name), - finest_level + finest_level, + skip_level_0 ); } //@} @@ -550,24 +554,28 @@ namespace ablastr::fields [[nodiscard]] MultiLevelVectorField get_mr_levels_alldirs ( T name, - int finest_level + int finest_level, + bool skip_level_0=false ) { return internal_get_mr_levels_alldirs( getExtractedName(name), - finest_level + finest_level, + skip_level_0 ); } template [[nodiscard]] ConstMultiLevelVectorField get_mr_levels_alldirs ( T name, - int finest_level + int finest_level, + bool skip_level_0=false ) const { return internal_get_mr_levels_alldirs( getExtractedName(name), - finest_level + finest_level, + skip_level_0 ); } //@} @@ -762,12 +770,14 @@ namespace ablastr::fields [[nodiscard]] MultiLevelScalarField internal_get_mr_levels ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ); [[nodiscard]] ConstMultiLevelScalarField internal_get_mr_levels ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ) const; [[nodiscard]] VectorField internal_get_alldirs ( @@ -782,12 +792,14 @@ namespace ablastr::fields [[nodiscard]] MultiLevelVectorField internal_get_mr_levels_alldirs ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ); [[nodiscard]] ConstMultiLevelVectorField internal_get_mr_levels_alldirs ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ) const; void diff --git a/Source/ablastr/fields/MultiFabRegister.cpp b/Source/ablastr/fields/MultiFabRegister.cpp index 722fa0ce106..709a0280bcf 100644 --- a/Source/ablastr/fields/MultiFabRegister.cpp +++ b/Source/ablastr/fields/MultiFabRegister.cpp @@ -417,14 +417,22 @@ namespace ablastr::fields MultiLevelScalarField MultiFabRegister::internal_get_mr_levels ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ) { MultiLevelScalarField field_on_level; field_on_level.reserve(finest_level+1); for (int lvl = 0; lvl <= finest_level; lvl++) { - field_on_level.push_back(internal_get(name, lvl)); + if (lvl == 0 && skip_level_0) + { + field_on_level.push_back(nullptr); + } + else + { + field_on_level.push_back(internal_get(name, lvl)); + } } return field_on_level; } @@ -432,14 +440,22 @@ namespace ablastr::fields ConstMultiLevelScalarField MultiFabRegister::internal_get_mr_levels ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ) const { ConstMultiLevelScalarField field_on_level; field_on_level.reserve(finest_level+1); for (int lvl = 0; lvl <= finest_level; lvl++) { - field_on_level.push_back(internal_get(name, lvl)); + if (lvl == 0 && skip_level_0) + { + field_on_level.push_back(nullptr); + } + else + { + field_on_level.push_back(internal_get(name, lvl)); + } } return field_on_level; } @@ -481,7 +497,8 @@ namespace ablastr::fields MultiLevelVectorField MultiFabRegister::internal_get_mr_levels_alldirs ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ) { MultiLevelVectorField field_on_level; @@ -495,7 +512,14 @@ namespace ablastr::fields // insert components for (Direction const & dir : m_all_dirs) { - field_on_level[lvl][dir] = internal_get(name, dir, lvl); + if (lvl == 0 && skip_level_0) + { + field_on_level[lvl][dir] = nullptr; + } + else + { + field_on_level[lvl][dir] = internal_get(name, dir, lvl); + } } } return field_on_level; @@ -504,7 +528,8 @@ namespace ablastr::fields ConstMultiLevelVectorField MultiFabRegister::internal_get_mr_levels_alldirs ( std::string const & name, - int finest_level + int finest_level, + bool skip_level_0 ) const { ConstMultiLevelVectorField field_on_level; @@ -518,7 +543,14 @@ namespace ablastr::fields // insert components for (Direction const & dir : m_all_dirs) { - field_on_level[lvl][dir] = internal_get(name, dir, lvl); + if (lvl == 0 && skip_level_0) + { + field_on_level[lvl][dir] = nullptr; + } + else + { + field_on_level[lvl][dir] = internal_get(name, dir, lvl); + } } } return field_on_level; From 29e08fe3919308c882e13aed676c39c844b22edd Mon Sep 17 00:00:00 2001 From: Edoardo Zoni Date: Mon, 21 Oct 2024 11:35:05 -0700 Subject: [PATCH 3/6] Fix bug in `WarpX::SyncCurrent` --- Source/Parallelization/WarpXComm.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Source/Parallelization/WarpXComm.cpp b/Source/Parallelization/WarpXComm.cpp index 3bf0b887f6f..318d363dfeb 100644 --- a/Source/Parallelization/WarpXComm.cpp +++ b/Source/Parallelization/WarpXComm.cpp @@ -1088,7 +1088,7 @@ WarpX::SyncCurrent (const std::string& current_fp_string) // If warpx.do_current_centering = 1, center currents from nodal grid to staggered grid if (do_current_centering) { - ablastr::fields::MultiLevelVectorField const& J_fp_nodal = m_fields.get_mr_levels_alldirs(FieldType::current_fp_nodal, finest_level+1); + ablastr::fields::MultiLevelVectorField const& J_fp_nodal = m_fields.get_mr_levels_alldirs(FieldType::current_fp_nodal, finest_level); AMREX_ALWAYS_ASSERT_WITH_MESSAGE(finest_level <= 1, "warpx.do_current_centering=1 not supported with more than one fine levels"); From a88403e11b48f102b86a2d1538236744636e634a Mon Sep 17 00:00:00 2001 From: Edoardo Zoni Date: Mon, 21 Oct 2024 11:39:01 -0700 Subject: [PATCH 4/6] Remove `return` statement after throwing runtime error --- Source/ablastr/fields/MultiFabRegister.cpp | 2 -- 1 file changed, 2 deletions(-) diff --git a/Source/ablastr/fields/MultiFabRegister.cpp b/Source/ablastr/fields/MultiFabRegister.cpp index 709a0280bcf..a1266deeab0 100644 --- a/Source/ablastr/fields/MultiFabRegister.cpp +++ b/Source/ablastr/fields/MultiFabRegister.cpp @@ -351,7 +351,6 @@ namespace ablastr::fields { if (m_mf_register.count(internal_name) == 0) { throw std::runtime_error("MultiFabRegister::get name does not exist in register: " + internal_name); - return nullptr; } amrex::MultiFab & mf = m_mf_register.at(internal_name).m_mf; @@ -365,7 +364,6 @@ namespace ablastr::fields { if (m_mf_register.count(internal_name) == 0) { throw std::runtime_error("MultiFabRegister::get name does not exist in register: " + internal_name); - return nullptr; } amrex::MultiFab const & mf = m_mf_register.at(internal_name).m_mf; From eddd8f7398352de84daddcb526a6198feae34541 Mon Sep 17 00:00:00 2001 From: Edoardo Zoni Date: Mon, 21 Oct 2024 14:07:09 -0700 Subject: [PATCH 5/6] Fix bug in `PML::PushPSATD` --- Source/BoundaryConditions/PML.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/Source/BoundaryConditions/PML.cpp b/Source/BoundaryConditions/PML.cpp index f45ca222e69..102f83ba478 100644 --- a/Source/BoundaryConditions/PML.cpp +++ b/Source/BoundaryConditions/PML.cpp @@ -1298,16 +1298,16 @@ PML::PushPSATD (ablastr::fields::MultiFabRegister& fields, const int lev) { ablastr::fields::VectorField pml_E_fp = fields.get_alldirs(FieldType::pml_E_fp, lev); ablastr::fields::VectorField pml_B_fp = fields.get_alldirs(FieldType::pml_B_fp, lev); - ablastr::fields::ScalarField pml_F_fp = fields.get(FieldType::pml_F_fp, lev); - ablastr::fields::ScalarField pml_G_fp = fields.get(FieldType::pml_G_fp, lev); + ablastr::fields::ScalarField pml_F_fp = (fields.has(FieldType::pml_F_fp, lev)) ? fields.get(FieldType::pml_F_fp, lev) : nullptr; + ablastr::fields::ScalarField pml_G_fp = (fields.has(FieldType::pml_G_fp, lev)) ? fields.get(FieldType::pml_G_fp, lev) : nullptr; // Update the fields on the fine and coarse patch PushPMLPSATDSinglePatch(lev, *spectral_solver_fp, pml_E_fp, pml_B_fp, pml_F_fp, pml_G_fp, m_fill_guards_fields); if (spectral_solver_cp) { ablastr::fields::VectorField pml_E_cp = fields.get_alldirs(FieldType::pml_E_cp, lev); ablastr::fields::VectorField pml_B_cp = fields.get_alldirs(FieldType::pml_B_cp, lev); - ablastr::fields::ScalarField pml_F_cp = fields.get(FieldType::pml_F_cp, lev); - ablastr::fields::ScalarField pml_G_cp = fields.get(FieldType::pml_G_cp, lev); + ablastr::fields::ScalarField pml_F_cp = (fields.has(FieldType::pml_F_cp, lev)) ? fields.get(FieldType::pml_F_cp, lev) : nullptr; + ablastr::fields::ScalarField pml_G_cp = (fields.has(FieldType::pml_G_cp, lev)) ? fields.get(FieldType::pml_G_cp, lev) : nullptr; PushPMLPSATDSinglePatch(lev, *spectral_solver_cp, pml_E_cp, pml_B_cp, pml_F_cp, pml_G_cp, m_fill_guards_fields); } } From 620d6eea510b2d1a8ac0651ab447a41ecf17689a Mon Sep 17 00:00:00 2001 From: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> Date: Mon, 21 Oct 2024 16:24:53 -0700 Subject: [PATCH 6/6] Remove debugging info left over --- Source/WarpX.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index c996f847c57..d1e3108e32a 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -2593,7 +2593,6 @@ WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm } // Create the MultiFabs for the current - amrex::Print() << "allocate current_cp" << std::endl; m_fields.alloc_init(FieldType::current_cp, Direction{0}, lev, amrex::convert(cba, jx_nodal_flag), dm, ncomps, ngJ, 0.0_rt); m_fields.alloc_init(FieldType::current_cp, Direction{1}, lev, amrex::convert(cba, jy_nodal_flag), dm, ncomps, ngJ, 0.0_rt); m_fields.alloc_init(FieldType::current_cp, Direction{2}, lev, amrex::convert(cba, jz_nodal_flag), dm, ncomps, ngJ, 0.0_rt);