From 0fe1905133033f128b235a14b04135c064800521 Mon Sep 17 00:00:00 2001 From: MaxThevenet Date: Tue, 29 Oct 2019 20:46:45 -0700 Subject: Reduce number of guard cells exchanged in Moving window and EvolveEM --- Source/Utils/WarpXMovingWindow.cpp | 54 +++++++++++++++++++++++--------------- 1 file changed, 33 insertions(+), 21 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index c577da7f3..59810d817 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -99,32 +99,34 @@ WarpX::MoveWindow (bool move_j) for (int dim = 0; dim < 3; ++dim) { // Fine grid - shiftMF(*Bfield_fp[lev][dim], geom[lev], num_shift, dir, B_external_grid[dim]); - shiftMF(*Efield_fp[lev][dim], geom[lev], num_shift, dir, E_external_grid[dim]); + shiftMF(*Bfield_fp[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE, B_external_grid[dim]); + shiftMF(*Efield_fp[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE, E_external_grid[dim]); if (move_j) { - shiftMF(*current_fp[lev][dim], geom[lev], num_shift, dir); + shiftMF(*current_fp[lev][dim], geom[lev], num_shift, dir, guard_cells.ngJ); } if (do_pml && pml[lev]->ok()) { const std::array& pml_B = pml[lev]->GetB_fp(); const std::array& pml_E = pml[lev]->GetE_fp(); - shiftMF(*pml_B[dim], geom[lev], num_shift, dir); - shiftMF(*pml_E[dim], geom[lev], num_shift, dir); + IntVect ng_exchange = pml_B[dim]->nGrowVect(); + shiftMF(*pml_B[dim], geom[lev], num_shift, dir, ng_exchange); + shiftMF(*pml_E[dim], geom[lev], num_shift, dir, ng_exchange); } if (lev > 0) { // Coarse grid - shiftMF(*Bfield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, B_external_grid[dim]); - shiftMF(*Efield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, E_external_grid[dim]); - shiftMF(*Bfield_aux[lev][dim], geom[lev], num_shift, dir); - shiftMF(*Efield_aux[lev][dim], geom[lev], num_shift, dir); + shiftMF(*Bfield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, guard_cells.ngE, B_external_grid[dim]); + shiftMF(*Efield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, guard_cells.ngE, E_external_grid[dim]); + shiftMF(*Bfield_aux[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE); + shiftMF(*Efield_aux[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE); if (move_j) { - shiftMF(*current_cp[lev][dim], geom[lev-1], num_shift_crse, dir); + shiftMF(*current_cp[lev][dim], geom[lev-1], num_shift_crse, dir, guard_cells.ngJ); } if (do_pml && pml[lev]->ok()) { const std::array& pml_B = pml[lev]->GetB_cp(); const std::array& pml_E = pml[lev]->GetE_cp(); - shiftMF(*pml_B[dim], geom[lev-1], num_shift_crse, dir); - shiftMF(*pml_E[dim], geom[lev-1], num_shift_crse, dir); + IntVect ng_exchange = pml_B[dim]->nGrowVect(); + shiftMF(*pml_B[dim], geom[lev-1], num_shift_crse, dir, ng_exchange); + shiftMF(*pml_E[dim], geom[lev-1], num_shift_crse, dir, ng_exchange); } } } @@ -132,19 +134,21 @@ WarpX::MoveWindow (bool move_j) // Shift scalar component F for dive cleaning if (do_dive_cleaning) { // Fine grid - shiftMF(*F_fp[lev], geom[lev], num_shift, dir); + shiftMF(*F_fp[lev], geom[lev], num_shift, dir, guard_cells.ngF); if (do_pml && pml[lev]->ok()) { MultiFab* pml_F = pml[lev]->GetF_fp(); - shiftMF(*pml_F, geom[lev], num_shift, dir); + IntVect ng_exchange = pml_F->nGrowVect(); + shiftMF(*pml_F, geom[lev], num_shift, dir, ng_exchange); } if (lev > 0) { // Coarse grid - shiftMF(*F_cp[lev], geom[lev-1], num_shift_crse, dir); + shiftMF(*F_cp[lev], geom[lev-1], num_shift_crse, dir, guard_cells.ngF); if (do_pml && pml[lev]->ok()) { MultiFab* pml_F = pml[lev]->GetF_cp(); - shiftMF(*pml_F, geom[lev-1], num_shift_crse, dir); + IntVect ng_exchange = pml_F->nGrowVect(); + shiftMF(*pml_F, geom[lev-1], num_shift_crse, dir, ng_exchange); } - shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir); + shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir, guard_cells.ngRho); } } @@ -152,10 +156,10 @@ WarpX::MoveWindow (bool move_j) if (move_j) { if (rho_fp[lev]){ // Fine grid - shiftMF(*rho_fp[lev], geom[lev], num_shift, dir); + shiftMF(*rho_fp[lev], geom[lev], num_shift, dir, guard_cells.ngRho); if (lev > 0){ // Coarse grid - shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir); + shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir, guard_cells.ngRho); } } } @@ -204,7 +208,7 @@ WarpX::MoveWindow (bool move_j) void WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, - amrex::Real external_field) + IntVect ng_exchange, amrex::Real external_field) { BL_PROFILE("WarpX::shiftMF()"); const BoxArray& ba = mf.boxArray(); @@ -216,7 +220,15 @@ WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, MultiFab tmpmf(ba, dm, nc, ng); MultiFab::Copy(tmpmf, mf, 0, 0, nc, ng); - tmpmf.FillBoundary(geom.periodicity()); + + // Not sure why this is needed, but it is... + ng_exchange[0] = 1; + ng_exchange[1] = num_shift; // 2 + Print()<<"ng_exchange "< Date: Wed, 30 Oct 2019 11:40:55 -0700 Subject: no need to pass ng in ShiftMF --- Source/Evolve/WarpXEvolveEM.cpp | 13 ++++---- Source/Parallelization/GuardCellManager.H | 3 +- Source/Parallelization/GuardCellManager.cpp | 9 +++-- Source/Utils/WarpXMovingWindow.cpp | 52 +++++++++++++---------------- Source/WarpX.H | 2 +- 5 files changed, 39 insertions(+), 40 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index c98a4688e..8707b9fde 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -93,10 +93,9 @@ WarpX::EvolveEM (int numsteps) // FillBoundaryE(guard_cells.ngE); // This is probably overkill // FillBoundaryB(guard_cells.ngE); - IntVect my_nc; - my_nc = guard_cells.ng_FieldGather+guard_cells.ng_NCIFilter; - FillBoundaryE(my_nc); - FillBoundaryB(my_nc); + IntVect ng_gather = guard_cells.ng_FieldGather+guard_cells.ng_NCIFilter; + FillBoundaryE(ng_gather); + FillBoundaryB(ng_gather); FillBoundaryAux(guard_cells.ng_Aux); UpdateAuxilaryData(); } @@ -328,12 +327,12 @@ WarpX::OneStep_nosub (Real cur_time) FillBoundaryB(guard_cells.ngE); #else EvolveF(0.5*dt[0], DtType::FirstHalf); - FillBoundaryF(guard_cells.ngF); + FillBoundaryF(guard_cells.ng_FieldSolver); EvolveB(0.5*dt[0]); // We now have B^{n+1/2} - FillBoundaryB(guard_cells.ngE_FieldSolver); + FillBoundaryB(guard_cells.ng_FieldSolver); EvolveE(dt[0]); // We now have E^{n+1} - FillBoundaryE(guard_cells.ngE_FieldSolver); + FillBoundaryE(guard_cells.ng_FieldSolver); EvolveF(0.5*dt[0], DtType::SecondHalf); FillBoundaryF(guard_cells.ngF); EvolveB(0.5*dt[0]); // We now have B^{n+1} diff --git a/Source/Parallelization/GuardCellManager.H b/Source/Parallelization/GuardCellManager.H index 706b5df79..e241eed75 100644 --- a/Source/Parallelization/GuardCellManager.H +++ b/Source/Parallelization/GuardCellManager.H @@ -30,8 +30,7 @@ public: int ngF_int = 0; // Guard cells to exchange data - amrex::IntVect ngB_FieldSolver = amrex::IntVect::TheZeroVector(); - amrex::IntVect ngE_FieldSolver = amrex::IntVect::TheZeroVector(); + amrex::IntVect ng_FieldSolver = amrex::IntVect::TheZeroVector(); amrex::IntVect ng_FieldGather = amrex::IntVect::TheZeroVector(); amrex::IntVect ngJ_CurrentDepo = amrex::IntVect::TheZeroVector(); amrex::IntVect ng_MovingWindow = amrex::IntVect::TheZeroVector(); diff --git a/Source/Parallelization/GuardCellManager.cpp b/Source/Parallelization/GuardCellManager.cpp index 166f0d58d..a67e15eb7 100644 --- a/Source/Parallelization/GuardCellManager.cpp +++ b/Source/Parallelization/GuardCellManager.cpp @@ -1,4 +1,5 @@ #include "GuardCellManager.H" +#include using namespace amrex; @@ -121,11 +122,15 @@ guardCellManager::Init( ngF = IntVect(AMREX_D_DECL(ngF_int, ngF_int, ngF_int)); #endif + Print()<<"rrr ngE : "<(aux_is_nodal and !do_nodal)); // Guard cells for field solver - ngB_FieldSolver = IntVect(AMREX_D_DECL(1,1,1)); - ngE_FieldSolver = IntVect(AMREX_D_DECL(1,1,1)); + ng_FieldSolver = IntVect(AMREX_D_DECL(1,1,1)); ng_MovingWindow = IntVect(AMREX_D_DECL(0,0,0)); // Multiplied by number of cells moved at each timestep ng_MovingWindow[moving_window_dir] = 1; int FGcell[4] = {0,1,1,2}; // Index is nox diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index 59810d817..db58ad53c 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -99,34 +99,32 @@ WarpX::MoveWindow (bool move_j) for (int dim = 0; dim < 3; ++dim) { // Fine grid - shiftMF(*Bfield_fp[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE, B_external_grid[dim]); - shiftMF(*Efield_fp[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE, E_external_grid[dim]); + shiftMF(*Bfield_fp[lev][dim], geom[lev], num_shift, dir, B_external_grid[dim]); + shiftMF(*Efield_fp[lev][dim], geom[lev], num_shift, dir, E_external_grid[dim]); if (move_j) { - shiftMF(*current_fp[lev][dim], geom[lev], num_shift, dir, guard_cells.ngJ); + shiftMF(*current_fp[lev][dim], geom[lev], num_shift, dir); } if (do_pml && pml[lev]->ok()) { const std::array& pml_B = pml[lev]->GetB_fp(); const std::array& pml_E = pml[lev]->GetE_fp(); - IntVect ng_exchange = pml_B[dim]->nGrowVect(); - shiftMF(*pml_B[dim], geom[lev], num_shift, dir, ng_exchange); - shiftMF(*pml_E[dim], geom[lev], num_shift, dir, ng_exchange); + shiftMF(*pml_B[dim], geom[lev], num_shift, dir); + shiftMF(*pml_E[dim], geom[lev], num_shift, dir); } if (lev > 0) { // Coarse grid - shiftMF(*Bfield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, guard_cells.ngE, B_external_grid[dim]); - shiftMF(*Efield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, guard_cells.ngE, E_external_grid[dim]); - shiftMF(*Bfield_aux[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE); - shiftMF(*Efield_aux[lev][dim], geom[lev], num_shift, dir, guard_cells.ngE); + shiftMF(*Bfield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, B_external_grid[dim]); + shiftMF(*Efield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, E_external_grid[dim]); + shiftMF(*Bfield_aux[lev][dim], geom[lev], num_shift, dir); + shiftMF(*Efield_aux[lev][dim], geom[lev], num_shift, dir); if (move_j) { - shiftMF(*current_cp[lev][dim], geom[lev-1], num_shift_crse, dir, guard_cells.ngJ); + shiftMF(*current_cp[lev][dim], geom[lev-1], num_shift_crse, dir); } if (do_pml && pml[lev]->ok()) { const std::array& pml_B = pml[lev]->GetB_cp(); const std::array& pml_E = pml[lev]->GetE_cp(); - IntVect ng_exchange = pml_B[dim]->nGrowVect(); - shiftMF(*pml_B[dim], geom[lev-1], num_shift_crse, dir, ng_exchange); - shiftMF(*pml_E[dim], geom[lev-1], num_shift_crse, dir, ng_exchange); + shiftMF(*pml_B[dim], geom[lev-1], num_shift_crse, dir); + shiftMF(*pml_E[dim], geom[lev-1], num_shift_crse, dir); } } } @@ -134,21 +132,19 @@ WarpX::MoveWindow (bool move_j) // Shift scalar component F for dive cleaning if (do_dive_cleaning) { // Fine grid - shiftMF(*F_fp[lev], geom[lev], num_shift, dir, guard_cells.ngF); + shiftMF(*F_fp[lev], geom[lev], num_shift, dir); if (do_pml && pml[lev]->ok()) { MultiFab* pml_F = pml[lev]->GetF_fp(); - IntVect ng_exchange = pml_F->nGrowVect(); - shiftMF(*pml_F, geom[lev], num_shift, dir, ng_exchange); + shiftMF(*pml_F, geom[lev], num_shift, dir); } if (lev > 0) { // Coarse grid - shiftMF(*F_cp[lev], geom[lev-1], num_shift_crse, dir, guard_cells.ngF); + shiftMF(*F_cp[lev], geom[lev-1], num_shift_crse, dir); if (do_pml && pml[lev]->ok()) { MultiFab* pml_F = pml[lev]->GetF_cp(); - IntVect ng_exchange = pml_F->nGrowVect(); - shiftMF(*pml_F, geom[lev-1], num_shift_crse, dir, ng_exchange); + shiftMF(*pml_F, geom[lev-1], num_shift_crse, dir); } - shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir, guard_cells.ngRho); + shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir); } } @@ -156,10 +152,10 @@ WarpX::MoveWindow (bool move_j) if (move_j) { if (rho_fp[lev]){ // Fine grid - shiftMF(*rho_fp[lev], geom[lev], num_shift, dir, guard_cells.ngRho); + shiftMF(*rho_fp[lev], geom[lev], num_shift, dir); if (lev > 0){ // Coarse grid - shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir, guard_cells.ngRho); + shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir); } } } @@ -208,7 +204,7 @@ WarpX::MoveWindow (bool move_j) void WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, - IntVect ng_exchange, amrex::Real external_field) + amrex::Real external_field) { BL_PROFILE("WarpX::shiftMF()"); const BoxArray& ba = mf.boxArray(); @@ -222,11 +218,11 @@ WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, MultiFab::Copy(tmpmf, mf, 0, 0, nc, ng); // Not sure why this is needed, but it is... - ng_exchange[0] = 1; - ng_exchange[1] = num_shift; // 2 - Print()<<"ng_exchange "< Date: Thu, 31 Oct 2019 16:40:10 -0700 Subject: Each call to FillBoundary in regular PIC loop (FDTD) uses fewer guard cells --- Source/Evolve/WarpXEvolveEM.cpp | 39 +++++++++++++---------- Source/Parallelization/GuardCellManager.H | 1 + Source/Parallelization/GuardCellManager.cpp | 27 +++++++++++----- Source/Parallelization/WarpXComm.cpp | 49 ++++++++++++++++++----------- Source/Utils/WarpXMovingWindow.cpp | 48 +++++++++++++++------------- Source/WarpX.H | 10 +++--- 6 files changed, 104 insertions(+), 70 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index 8707b9fde..3247ca64f 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -75,8 +75,10 @@ WarpX::EvolveEM (int numsteps) // Particles have p^{n} and x^{n}. // is_synchronized is true. if (is_synchronized) { - FillBoundaryE(guard_cells.ng_FieldGather); - FillBoundaryB(guard_cells.ng_FieldGather); + // FillBoundaryE(guard_cells.ng_FieldGather, guard_cells.ngExtra); + // FillBoundaryB(guard_cells.ng_FieldGather, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); + FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); UpdateAuxilaryData(); // on first step, push p by -0.5*dt for (int lev = 0; lev <= finest_level; ++lev) { @@ -90,12 +92,14 @@ WarpX::EvolveEM (int numsteps) // Beyond one step, we have E^{n} and B^{n}. // Particles have p^{n-1/2} and x^{n}. // This is probably overkill - // FillBoundaryE(guard_cells.ngE); // This is probably overkill - // FillBoundaryB(guard_cells.ngE); - IntVect ng_gather = guard_cells.ng_FieldGather+guard_cells.ng_NCIFilter; - FillBoundaryE(ng_gather); - FillBoundaryB(ng_gather); + // IntVect ng_gather = guard_cells.ng_FieldGather+guard_cells.ng_NCIFilter; + // FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); + // FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); + // FillBoundaryE(ng_gather, guard_cells.ngExtra); + // FillBoundaryB(ng_gather, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_FieldGatherAndNCIFilter, guard_cells.ngExtra); + FillBoundaryB(guard_cells.ng_FieldGatherAndNCIFilter, guard_cells.ngExtra); FillBoundaryAux(guard_cells.ng_Aux); UpdateAuxilaryData(); } @@ -195,9 +199,9 @@ WarpX::EvolveEM (int numsteps) if (to_make_plot || do_insitu || to_make_slice_plot) { // This is probably overkill - FillBoundaryE(guard_cells.ngE); + FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); // This is probably overkill - FillBoundaryB(guard_cells.ngE); + FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); FillBoundaryAux(guard_cells.ng_Aux); UpdateAuxilaryData(); @@ -250,9 +254,9 @@ WarpX::EvolveEM (int numsteps) if (write_plot_file || do_insitu) { // This is probably overkill - FillBoundaryE(guard_cells.ngE); + FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); // This is probably overkill - FillBoundaryB(guard_cells.ngE); + FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); FillBoundaryAux(guard_cells.ng_Aux); UpdateAuxilaryData(); @@ -323,23 +327,24 @@ WarpX::OneStep_nosub (Real cur_time) #ifdef WARPX_USE_PSATD PushPSATD(dt[0]); if (do_pml) DampPML(); - FillBoundaryE(guard_cells.ngE); - FillBoundaryB(guard_cells.ngE); + FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); + FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); #else + Print()<<"Push fields\n"; EvolveF(0.5*dt[0], DtType::FirstHalf); FillBoundaryF(guard_cells.ng_FieldSolver); EvolveB(0.5*dt[0]); // We now have B^{n+1/2} - FillBoundaryB(guard_cells.ng_FieldSolver); - + // FillBoundaryB(guard_cells.ng_FieldSolver, guard_cells.ngExtra); + FillBoundaryB(guard_cells.ng_FieldSolver, IntVect::TheZeroVector()); EvolveE(dt[0]); // We now have E^{n+1} - FillBoundaryE(guard_cells.ng_FieldSolver); + // FillBoundaryE(guard_cells.ng_FieldSolver, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_FieldSolver, IntVect::TheZeroVector()); EvolveF(0.5*dt[0], DtType::SecondHalf); FillBoundaryF(guard_cells.ngF); EvolveB(0.5*dt[0]); // We now have B^{n+1} if (do_pml) { DampPML(); } - #endif } diff --git a/Source/Parallelization/GuardCellManager.H b/Source/Parallelization/GuardCellManager.H index e241eed75..93f1e9fd1 100644 --- a/Source/Parallelization/GuardCellManager.H +++ b/Source/Parallelization/GuardCellManager.H @@ -32,6 +32,7 @@ public: // Guard cells to exchange data amrex::IntVect ng_FieldSolver = amrex::IntVect::TheZeroVector(); amrex::IntVect ng_FieldGather = amrex::IntVect::TheZeroVector(); + amrex::IntVect ng_FieldGatherAndNCIFilter = amrex::IntVect::TheZeroVector(); amrex::IntVect ngJ_CurrentDepo = amrex::IntVect::TheZeroVector(); amrex::IntVect ng_MovingWindow = amrex::IntVect::TheZeroVector(); amrex::IntVect ng_NCIFilter = amrex::IntVect::TheZeroVector(); diff --git a/Source/Parallelization/GuardCellManager.cpp b/Source/Parallelization/GuardCellManager.cpp index a67e15eb7..9f3e5b18d 100644 --- a/Source/Parallelization/GuardCellManager.cpp +++ b/Source/Parallelization/GuardCellManager.cpp @@ -110,7 +110,7 @@ guardCellManager::Init( ng_required = std::max( ng_required, ngE[i_dim] ); ng_required = std::max( ng_required, ngJ[i_dim] ); ng_required = std::max( ng_required, ngRho[i_dim] ); - ng_required = std::max( ng_required, ngF[i_dim] ); +v ng_required = std::max( ng_required, ngF[i_dim] ); // Set the guard cells to this max ngE[i_dim] = ng_required; ngJ[i_dim] = ng_required; @@ -122,11 +122,6 @@ guardCellManager::Init( ngF = IntVect(AMREX_D_DECL(ngF_int, ngF_int, ngF_int)); #endif - Print()<<"rrr ngE : "<(aux_is_nodal and !do_nodal)); // Guard cells for field solver @@ -134,12 +129,28 @@ guardCellManager::Init( ng_MovingWindow = IntVect(AMREX_D_DECL(0,0,0)); // Multiplied by number of cells moved at each timestep ng_MovingWindow[moving_window_dir] = 1; int FGcell[4] = {0,1,1,2}; // Index is nox - ng_FieldGather = IntVect(AMREX_D_DECL(FGcell[nox],FGcell[nox],FGcell[nox])); - ngJ_CurrentDepo = ng_FieldGather; + ng_FieldGather = IntVect(AMREX_D_DECL(FGcell[nox],FGcell[nox],FGcell[nox])) + 2*ngExtra; + + ng_FieldGather = ng_FieldGather.min(ngE); if (do_fdtd_nci_corr){ ng_NCIFilter = IntVect::TheZeroVector(); ng_NCIFilter[AMREX_SPACEDIM-1] = 4; } + + ng_FieldGatherAndNCIFilter = ng_FieldGather + ng_NCIFilter; + ng_FieldGatherAndNCIFilter = ng_FieldGatherAndNCIFilter.min(ngE); + ng_Aux = 2*ng_FieldGather+ng_NCIFilter; ng_Aux = ng_Aux.min(ngE); + + Print()<<"rrr ngE : "< #include #include +#include "WarpXAlgorithmSelection.H" #include #include @@ -321,20 +322,20 @@ WarpX::UpdateAuxilaryDataSameType () } void -WarpX::FillBoundaryB (IntVect ng) +WarpX::FillBoundaryB (IntVect ng, IntVect ng_extra_fine) { for (int lev = 0; lev <= finest_level; ++lev) { - FillBoundaryB(lev, ng); + FillBoundaryB(lev, ng, ng_extra_fine); } } void -WarpX::FillBoundaryE (IntVect ng) +WarpX::FillBoundaryE (IntVect ng, IntVect ng_extra_fine) { for (int lev = 0; lev <= finest_level; ++lev) { - FillBoundaryE(lev, ng); + FillBoundaryE(lev, ng, ng_extra_fine); } } @@ -348,9 +349,11 @@ WarpX::FillBoundaryF (IntVect ng) } void -WarpX::FillBoundaryE(int lev, IntVect ng) +WarpX::FillBoundaryE(int lev, IntVect ng, IntVect ng_extra_fine) { - FillBoundaryE(lev, PatchType::fine, ng); + Print()<<"FillBoundaryE ng_extra_fine "<< ng_extra_fine <<'\n'; + Print()<<"FillBoundaryE exchanges "<< ng+ng_extra_fine <<'\n'; + FillBoundaryE(lev, PatchType::fine, ng+ng_extra_fine); if (lev > 0) FillBoundaryE(lev, PatchType::coarse, ng); } @@ -370,11 +373,12 @@ WarpX::FillBoundaryE (int lev, PatchType patch_type, IntVect ng) } const auto& period = Geom(lev).periodicity(); + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + ng <= Efield_fp[lev][0]->nGrowVect(), + "Error: in FillBoundaryE, requested more guard cells than allocated"); Efield_fp[lev][0]->FillBoundary(0, Efield_fp[lev][0]->nComp(), ng, period); Efield_fp[lev][1]->FillBoundary(0, Efield_fp[lev][1]->nComp(), ng, period); Efield_fp[lev][2]->FillBoundary(0, Efield_fp[lev][2]->nComp(), ng, period); - // Vector mf{Efield_fp[lev][0].get(),Efield_fp[lev][1].get(),Efield_fp[lev][2].get()}; - // amrex::FillBoundary(mf, period); } else if (patch_type == PatchType::coarse) { @@ -389,18 +393,21 @@ WarpX::FillBoundaryE (int lev, PatchType patch_type, IntVect ng) } const auto& cperiod = Geom(lev-1).periodicity(); + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + ng <= Efield_cp[lev][0]->nGrowVect(), + "Error: in FillBoundaryE, requested more guard cells than allocated"); Efield_cp[lev][0]->FillBoundary(0, Efield_cp[lev][0]->nComp(), ng, cperiod); Efield_cp[lev][1]->FillBoundary(0, Efield_cp[lev][1]->nComp(), ng, cperiod); Efield_cp[lev][2]->FillBoundary(0, Efield_cp[lev][2]->nComp(), ng, cperiod); - // Vector mf{Efield_cp[lev][0].get(),Efield_cp[lev][1].get(),Efield_cp[lev][2].get()}; - // amrex::FillBoundary(mf, cperiod); } } void -WarpX::FillBoundaryB (int lev, IntVect ng) +WarpX::FillBoundaryB (int lev, IntVect ng, IntVect ng_extra_fine) { - FillBoundaryB(lev, PatchType::fine, ng); + Print()<<"FillBoundaryB ng_extra_fine "<< ng_extra_fine <<'\n'; + Print()<<"FillBoundaryB exchanges "<< ng+ng_extra_fine <<'\n'; + FillBoundaryB(lev, PatchType::fine, ng + ng_extra_fine); if (lev > 0) FillBoundaryB(lev, PatchType::coarse, ng); } @@ -419,11 +426,12 @@ WarpX::FillBoundaryB (int lev, PatchType patch_type, IntVect ng) pml[lev]->FillBoundaryB(patch_type); } const auto& period = Geom(lev).periodicity(); + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + ng <= Bfield_fp[lev][0]->nGrowVect(), + "Error: in FillBoundaryB, requested more guard cells than allocated"); Bfield_fp[lev][0]->FillBoundary(0, Bfield_fp[lev][0]->nComp(), ng, period); Bfield_fp[lev][1]->FillBoundary(0, Bfield_fp[lev][1]->nComp(), ng, period); Bfield_fp[lev][2]->FillBoundary(0, Bfield_fp[lev][2]->nComp(), ng, period); - // Vector mf{Bfield_fp[lev][0].get(),Bfield_fp[lev][1].get(),Bfield_fp[lev][2].get()}; - // amrex::FillBoundary(mf, period); } else if (patch_type == PatchType::coarse) { @@ -437,11 +445,12 @@ WarpX::FillBoundaryB (int lev, PatchType patch_type, IntVect ng) pml[lev]->FillBoundaryB(patch_type); } const auto& cperiod = Geom(lev-1).periodicity(); + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + ng <= Bfield_cp[lev][0]->nGrowVect(), + "Error: in FillBoundaryB, requested more guard cells than allocated"); Bfield_cp[lev][0]->FillBoundary(0, Bfield_cp[lev][0]->nComp(), ng, cperiod); Bfield_cp[lev][1]->FillBoundary(0, Bfield_cp[lev][1]->nComp(), ng, cperiod); Bfield_cp[lev][2]->FillBoundary(0, Bfield_cp[lev][2]->nComp(), ng, cperiod); - // Vector mf{Bfield_cp[lev][0].get(),Bfield_cp[lev][1].get(),Bfield_cp[lev][2].get()}; - // amrex::FillBoundary(mf, cperiod); } } @@ -465,8 +474,10 @@ WarpX::FillBoundaryF (int lev, PatchType patch_type, IntVect ng) } const auto& period = Geom(lev).periodicity(); + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + ng <= F_fp[lev]->nGrowVect(), + "Error: in FillBoundaryF, requested more guard cells than allocated"); F_fp[lev]->FillBoundary(0, F_fp[lev]->nComp(), ng, period); - // F_fp[lev]->FillBoundary(period); } else if (patch_type == PatchType::coarse && F_cp[lev]) { @@ -478,8 +489,10 @@ WarpX::FillBoundaryF (int lev, PatchType patch_type, IntVect ng) } const auto& cperiod = Geom(lev-1).periodicity(); + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + ng <= F_cp[lev]->nGrowVect(), + "Error: in FillBoundaryF, requested more guard cells than allocated"); F_cp[lev]->FillBoundary(0, F_cp[lev]->nComp(), ng, cperiod); - // F_cp[lev]->FillBoundary(cperiod); } } diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index db58ad53c..8770d2059 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -1,3 +1,4 @@ +#include "GuardCellManager.H" #include #include @@ -28,6 +29,9 @@ WarpX::MoveWindow (bool move_j) { if (do_moving_window == 0) return 0; + IntVect ng_extra = guard_cells.ngExtra; + IntVect ng_zero = IntVect::TheZeroVector(); + // Update the continuous position of the moving window, // and of the plasma injection moving_window_x += moving_window_v * dt[0]; @@ -99,32 +103,32 @@ WarpX::MoveWindow (bool move_j) for (int dim = 0; dim < 3; ++dim) { // Fine grid - shiftMF(*Bfield_fp[lev][dim], geom[lev], num_shift, dir, B_external_grid[dim]); - shiftMF(*Efield_fp[lev][dim], geom[lev], num_shift, dir, E_external_grid[dim]); + shiftMF(*Bfield_fp[lev][dim], geom[lev], num_shift, dir, ng_extra, B_external_grid[dim]); + shiftMF(*Efield_fp[lev][dim], geom[lev], num_shift, dir, ng_extra, E_external_grid[dim]); if (move_j) { - shiftMF(*current_fp[lev][dim], geom[lev], num_shift, dir); + shiftMF(*current_fp[lev][dim], geom[lev], num_shift, dir, ng_zero); } if (do_pml && pml[lev]->ok()) { const std::array& pml_B = pml[lev]->GetB_fp(); const std::array& pml_E = pml[lev]->GetE_fp(); - shiftMF(*pml_B[dim], geom[lev], num_shift, dir); - shiftMF(*pml_E[dim], geom[lev], num_shift, dir); + shiftMF(*pml_B[dim], geom[lev], num_shift, dir, ng_extra); + shiftMF(*pml_E[dim], geom[lev], num_shift, dir, ng_extra); } if (lev > 0) { // Coarse grid - shiftMF(*Bfield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, B_external_grid[dim]); - shiftMF(*Efield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, E_external_grid[dim]); - shiftMF(*Bfield_aux[lev][dim], geom[lev], num_shift, dir); - shiftMF(*Efield_aux[lev][dim], geom[lev], num_shift, dir); + shiftMF(*Bfield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, ng_zero, B_external_grid[dim]); + shiftMF(*Efield_cp[lev][dim], geom[lev-1], num_shift_crse, dir, ng_zero, E_external_grid[dim]); + shiftMF(*Bfield_aux[lev][dim], geom[lev], num_shift, dir, ng_zero); + shiftMF(*Efield_aux[lev][dim], geom[lev], num_shift, dir, ng_zero); if (move_j) { - shiftMF(*current_cp[lev][dim], geom[lev-1], num_shift_crse, dir); + shiftMF(*current_cp[lev][dim], geom[lev-1], num_shift_crse, dir, ng_zero); } if (do_pml && pml[lev]->ok()) { const std::array& pml_B = pml[lev]->GetB_cp(); const std::array& pml_E = pml[lev]->GetE_cp(); - shiftMF(*pml_B[dim], geom[lev-1], num_shift_crse, dir); - shiftMF(*pml_E[dim], geom[lev-1], num_shift_crse, dir); + shiftMF(*pml_B[dim], geom[lev-1], num_shift_crse, dir, ng_extra); + shiftMF(*pml_E[dim], geom[lev-1], num_shift_crse, dir, ng_extra); } } } @@ -132,19 +136,19 @@ WarpX::MoveWindow (bool move_j) // Shift scalar component F for dive cleaning if (do_dive_cleaning) { // Fine grid - shiftMF(*F_fp[lev], geom[lev], num_shift, dir); + shiftMF(*F_fp[lev], geom[lev], num_shift, dir, ng_zero); if (do_pml && pml[lev]->ok()) { MultiFab* pml_F = pml[lev]->GetF_fp(); - shiftMF(*pml_F, geom[lev], num_shift, dir); + shiftMF(*pml_F, geom[lev], num_shift, dir, ng_extra); } if (lev > 0) { // Coarse grid - shiftMF(*F_cp[lev], geom[lev-1], num_shift_crse, dir); + shiftMF(*F_cp[lev], geom[lev-1], num_shift_crse, dir, ng_zero); if (do_pml && pml[lev]->ok()) { MultiFab* pml_F = pml[lev]->GetF_cp(); - shiftMF(*pml_F, geom[lev-1], num_shift_crse, dir); + shiftMF(*pml_F, geom[lev-1], num_shift_crse, dir, ng_zero); } - shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir); + shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir, ng_zero); } } @@ -152,10 +156,10 @@ WarpX::MoveWindow (bool move_j) if (move_j) { if (rho_fp[lev]){ // Fine grid - shiftMF(*rho_fp[lev], geom[lev], num_shift, dir); + shiftMF(*rho_fp[lev], geom[lev], num_shift, dir, ng_zero); if (lev > 0){ // Coarse grid - shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir); + shiftMF(*rho_cp[lev], geom[lev-1], num_shift_crse, dir, ng_zero); } } } @@ -204,7 +208,7 @@ WarpX::MoveWindow (bool move_j) void WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, - amrex::Real external_field) + IntVect ng_extra, amrex::Real external_field) { BL_PROFILE("WarpX::shiftMF()"); const BoxArray& ba = mf.boxArray(); @@ -220,12 +224,12 @@ WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, // Not sure why this is needed, but it is... IntVect ng_mw = IntVect::TheUnitVector(); ng_mw[dir] = num_shift; + ng_mw += ng_extra; + ng_mw = ng_mw.min(ng); Print()<<"ng_mw "< Date: Thu, 31 Oct 2019 18:00:23 -0700 Subject: more systematic names for guard cell variables --- Source/Evolve/WarpXEvolveEM.cpp | 79 +++++++++++++++-------------- Source/Parallelization/GuardCellManager.H | 20 +++----- Source/Parallelization/GuardCellManager.cpp | 61 +++++++++------------- Source/Utils/WarpXMovingWindow.cpp | 2 +- Source/WarpX.H | 4 +- Source/WarpX.cpp | 24 ++++----- 6 files changed, 84 insertions(+), 106 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index 4e44aea48..2d28641d2 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -76,8 +76,8 @@ WarpX::EvolveEM (int numsteps) // is_synchronized is true. if (is_synchronized) { // Not called at each iteration, so exchange all guard cells - FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); - FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); + FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); UpdateAuxilaryData(); // on first step, push p by -0.5*dt for (int lev = 0; lev <= finest_level; ++lev) { @@ -92,10 +92,11 @@ WarpX::EvolveEM (int numsteps) // Particles have p^{n-1/2} and x^{n}. // E and B are up-to-date inside the domain only - FillBoundaryE(guard_cells.ng_FieldGather, guard_cells.ngExtra); - FillBoundaryB(guard_cells.ng_FieldGather, guard_cells.ngExtra); - // E and B are up-to-date inside the domain only - FillBoundaryAux(guard_cells.ng_Aux); + FillBoundaryE(guard_cells.ng_FieldGather, guard_cells.ng_Extra); + FillBoundaryB(guard_cells.ng_FieldGather, guard_cells.ng_Extra); + // E and B: enough guard cells to update Aux or call Field Gather in fp and cp + // Need to update Aux on lower levels, to interpolate to higher levels. + FillBoundaryAux(guard_cells.ng_UpdateAux); UpdateAuxilaryData(); } @@ -194,10 +195,10 @@ WarpX::EvolveEM (int numsteps) if (to_make_plot || do_insitu || to_make_slice_plot) { // This is probably overkill - FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); // This is probably overkill - FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); - FillBoundaryAux(guard_cells.ng_Aux); + FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); + FillBoundaryAux(guard_cells.ng_UpdateAux); UpdateAuxilaryData(); for (int lev = 0; lev <= finest_level; ++lev) { @@ -249,10 +250,10 @@ WarpX::EvolveEM (int numsteps) if (write_plot_file || do_insitu) { // This is probably overkill - FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); // This is probably overkill - FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); - FillBoundaryAux(guard_cells.ng_Aux); + FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); + FillBoundaryAux(guard_cells.ng_UpdateAux); UpdateAuxilaryData(); for (int lev = 0; lev <= finest_level; ++lev) { @@ -322,8 +323,8 @@ WarpX::OneStep_nosub (Real cur_time) #ifdef WARPX_USE_PSATD PushPSATD(dt[0]); if (do_pml) DampPML(); - FillBoundaryE(guard_cells.ngE, guard_cells.ngExtra); - FillBoundaryB(guard_cells.ngE, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); + FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); #else Print()<<"Push fields\n"; EvolveF(0.5*dt[0], DtType::FirstHalf); @@ -335,7 +336,7 @@ WarpX::OneStep_nosub (Real cur_time) // FillBoundaryE(guard_cells.ng_FieldSolver, guard_cells.ngExtra); FillBoundaryE(guard_cells.ng_FieldSolver, IntVect::TheZeroVector()); EvolveF(0.5*dt[0], DtType::SecondHalf); - FillBoundaryF(guard_cells.ngF); + FillBoundaryF(guard_cells.ng_alloc_F); EvolveB(0.5*dt[0]); // We now have B^{n+1} if (do_pml) { DampPML(); @@ -382,21 +383,21 @@ WarpX::OneStep_sub1 (Real curtime) EvolveB(fine_lev, PatchType::fine, 0.5*dt[fine_lev]); EvolveF(fine_lev, PatchType::fine, 0.5*dt[fine_lev], DtType::FirstHalf); - FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ngE); - FillBoundaryF(fine_lev, PatchType::fine, guard_cells.ngF); + FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); + FillBoundaryF(fine_lev, PatchType::fine, guard_cells.ng_alloc_F); EvolveE(fine_lev, PatchType::fine, dt[fine_lev]); - FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); EvolveB(fine_lev, PatchType::fine, 0.5*dt[fine_lev]); EvolveF(fine_lev, PatchType::fine, 0.5*dt[fine_lev], DtType::SecondHalf); if (do_pml) { DampPML(fine_lev, PatchType::fine); - FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); } - FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); // ii) Push particles on the coarse patch and mother grid. // Push the fields on the coarse patch and mother grid @@ -408,19 +409,19 @@ WarpX::OneStep_sub1 (Real curtime) EvolveB(fine_lev, PatchType::coarse, dt[fine_lev]); EvolveF(fine_lev, PatchType::coarse, dt[fine_lev], DtType::FirstHalf); - FillBoundaryB(fine_lev, PatchType::coarse, guard_cells.ngE); - FillBoundaryF(fine_lev, PatchType::coarse, guard_cells.ngF); + FillBoundaryB(fine_lev, PatchType::coarse, guard_cells.ng_alloc_EB); + FillBoundaryF(fine_lev, PatchType::coarse, guard_cells.ng_alloc_F); EvolveE(fine_lev, PatchType::coarse, dt[fine_lev]); - FillBoundaryE(fine_lev, PatchType::coarse, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::coarse, guard_cells.ng_alloc_EB); EvolveB(coarse_lev, PatchType::fine, 0.5*dt[coarse_lev]); EvolveF(coarse_lev, PatchType::fine, 0.5*dt[coarse_lev], DtType::FirstHalf); - FillBoundaryB(coarse_lev, PatchType::fine, guard_cells.ngE); - FillBoundaryF(coarse_lev, PatchType::fine, guard_cells.ngF); + FillBoundaryB(coarse_lev, PatchType::fine, guard_cells.ng_alloc_EB); + FillBoundaryF(coarse_lev, PatchType::fine, guard_cells.ng_alloc_F); EvolveE(coarse_lev, PatchType::fine, 0.5*dt[coarse_lev]); - FillBoundaryE(coarse_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(coarse_lev, PatchType::fine, guard_cells.ng_alloc_EB); // iii) Get auxiliary fields on the fine grid, at dt[fine_lev] UpdateAuxilaryData(); @@ -436,22 +437,22 @@ WarpX::OneStep_sub1 (Real curtime) EvolveB(fine_lev, PatchType::fine, 0.5*dt[fine_lev]); EvolveF(fine_lev, PatchType::fine, 0.5*dt[fine_lev], DtType::FirstHalf); - FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ngE); - FillBoundaryF(fine_lev, PatchType::fine, guard_cells.ngF); + FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); + FillBoundaryF(fine_lev, PatchType::fine, guard_cells.ng_alloc_F); EvolveE(fine_lev, PatchType::fine, dt[fine_lev]); - FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); EvolveB(fine_lev, PatchType::fine, 0.5*dt[fine_lev]); EvolveF(fine_lev, PatchType::fine, 0.5*dt[fine_lev], DtType::SecondHalf); if (do_pml) { DampPML(fine_lev, PatchType::fine); - FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); } - FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ngE); - FillBoundaryF(fine_lev, PatchType::fine, guard_cells.ngF); + FillBoundaryB(fine_lev, PatchType::fine, guard_cells.ng_alloc_EB); + FillBoundaryF(fine_lev, PatchType::fine, guard_cells.ng_alloc_F); // v) Push the fields on the coarse patch and mother grid // by only half a coarse step (second half) @@ -460,7 +461,7 @@ WarpX::OneStep_sub1 (Real curtime) AddRhoFromFineLevelandSumBoundary(coarse_lev, ncomps, ncomps); EvolveE(fine_lev, PatchType::coarse, dt[fine_lev]); - FillBoundaryE(fine_lev, PatchType::coarse, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::coarse, guard_cells.ng_alloc_EB); EvolveB(fine_lev, PatchType::coarse, dt[fine_lev]); EvolveF(fine_lev, PatchType::coarse, dt[fine_lev], DtType::SecondHalf); @@ -468,24 +469,24 @@ WarpX::OneStep_sub1 (Real curtime) if (do_pml) { DampPML(fine_lev, PatchType::coarse); // do it twice DampPML(fine_lev, PatchType::coarse); - FillBoundaryE(fine_lev, PatchType::coarse, guard_cells.ngE); + FillBoundaryE(fine_lev, PatchType::coarse, guard_cells.ng_alloc_EB); } - FillBoundaryB(fine_lev, PatchType::coarse, guard_cells.ngE); - FillBoundaryF(fine_lev, PatchType::coarse, guard_cells.ngF); + FillBoundaryB(fine_lev, PatchType::coarse, guard_cells.ng_alloc_EB); + FillBoundaryF(fine_lev, PatchType::coarse, guard_cells.ng_alloc_F); EvolveE(coarse_lev, PatchType::fine, 0.5*dt[coarse_lev]); - FillBoundaryE(coarse_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(coarse_lev, PatchType::fine, guard_cells.ng_alloc_EB); EvolveB(coarse_lev, PatchType::fine, 0.5*dt[coarse_lev]); EvolveF(coarse_lev, PatchType::fine, 0.5*dt[coarse_lev], DtType::SecondHalf); if (do_pml) { DampPML(coarse_lev, PatchType::fine); - FillBoundaryE(coarse_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryE(coarse_lev, PatchType::fine, guard_cells.ng_alloc_EB); } - FillBoundaryB(coarse_lev, PatchType::fine, guard_cells.ngE); + FillBoundaryB(coarse_lev, PatchType::fine, guard_cells.ng_alloc_EB); } void diff --git a/Source/Parallelization/GuardCellManager.H b/Source/Parallelization/GuardCellManager.H index 19b5ac14a..7e6a0d9f5 100644 --- a/Source/Parallelization/GuardCellManager.H +++ b/Source/Parallelization/GuardCellManager.H @@ -21,24 +21,16 @@ public: const int maxwell_fdtd_solver_id, const int max_level); - // Guard cells to initialize multifabs - amrex::IntVect ngE = amrex::IntVect::TheZeroVector(); - amrex::IntVect ngJ = amrex::IntVect::TheZeroVector(); - amrex::IntVect ngRho = amrex::IntVect::TheZeroVector(); - amrex::IntVect ngExtra = amrex::IntVect::TheZeroVector(); - amrex::IntVect ngF = amrex::IntVect::TheZeroVector(); - int ngF_int = 0; - - //amrex::IntVect ng_alloc_EB = amrex::IntVect::TheZeroVector(); - //amrex::IntVect ng_alloc_J = amrex::IntVect::TheZeroVector(); - //amrex::IntVect ng_alloc_Rho = amrex::IntVect::TheZeroVector(); - //amrex::IntVect ng_alloc_F = amrex::IntVect::TheZeroVector(); - //int ng_alloc_F_int = 0; + amrex::IntVect ng_alloc_EB = amrex::IntVect::TheZeroVector(); + amrex::IntVect ng_alloc_J = amrex::IntVect::TheZeroVector(); + amrex::IntVect ng_alloc_Rho = amrex::IntVect::TheZeroVector(); + amrex::IntVect ng_alloc_F = amrex::IntVect::TheZeroVector(); + int ng_alloc_F_int = 0; // Guard cells to exchange data amrex::IntVect ng_FieldSolver = amrex::IntVect::TheZeroVector(); amrex::IntVect ng_FieldGather = amrex::IntVect::TheZeroVector(); - amrex::IntVect ng_Aux = amrex::IntVect::TheZeroVector(); + amrex::IntVect ng_UpdateAux = amrex::IntVect::TheZeroVector(); // Extra guard cells for fine level of E and B amrex::IntVect ng_Extra = amrex::IntVect::TheZeroVector(); diff --git a/Source/Parallelization/GuardCellManager.cpp b/Source/Parallelization/GuardCellManager.cpp index 3777c7fcc..d970cd464 100644 --- a/Source/Parallelization/GuardCellManager.cpp +++ b/Source/Parallelization/GuardCellManager.cpp @@ -62,19 +62,19 @@ guardCellManager::Init( } #if (AMREX_SPACEDIM == 3) - ngE = IntVect(ngx,ngy,ngz); - ngJ = IntVect(ngJx,ngJy,ngJz); + ng_alloc_EB = IntVect(ngx,ngy,ngz); + ng_alloc_J = IntVect(ngJx,ngJy,ngJz); #elif (AMREX_SPACEDIM == 2) - ngE = IntVect(ngx,ngz); - ngJ = IntVect(ngJx,ngJz); + ng_alloc_EB = IntVect(ngx,ngz); + ng_alloc_J = IntVect(ngJx,ngJz); #endif - ngRho = ngJ+1; //One extra ghost cell, so that it's safe to deposit charge density + ng_alloc_Rho = ng_alloc_J+1; //One extra ghost cell, so that it's safe to deposit charge density // after pushing particle. - ngF_int = (do_moving_window) ? 2 : 0; + ng_alloc_F_int = (do_moving_window) ? 2 : 0; // CKC solver requires one additional guard cell - if (maxwell_fdtd_solver_id == 1) ngF_int = std::max( ngF_int, 1 ); - ngF = IntVect(AMREX_D_DECL(ngF_int, ngF_int, ngF_int)); + if (maxwell_fdtd_solver_id == 1) ng_alloc_F_int = std::max( ng_alloc_F_int, 1 ); + ng_alloc_F = IntVect(AMREX_D_DECL(ng_alloc_F_int, ng_alloc_F_int, ng_alloc_F_int)); #ifdef WARPX_USE_PSATD if (do_fft_mpi_dec == false){ @@ -94,22 +94,22 @@ guardCellManager::Init( for (int i_dim=0; i_dim(aux_is_nodal and !do_nodal)); + ng_Extra = IntVect(static_cast(aux_is_nodal and !do_nodal)); // Compute number of cells required for Field Solver ng_FieldSolver = IntVect(AMREX_D_DECL(1,1,1)); @@ -118,10 +118,10 @@ guardCellManager::Init( int FGcell[4] = {0,1,1,2}; // Index is nox IntVect ng_FieldGather_noNCI = IntVect(AMREX_D_DECL(FGcell[nox],FGcell[nox],FGcell[nox])); // Add one cell if momentum_conserving gather in a staggered-field simulation - ng_FieldGather_noNCI += ngExtra; + ng_FieldGather_noNCI += ng_Extra; // Not sure why, but need one extra guard cell when using MR - if (max_level >= 1) ng_FieldGather_noNCI += ngExtra; - ng_FieldGather_noNCI = ng_FieldGather_noNCI.min(ngE); + if (max_level >= 1) ng_FieldGather_noNCI += ng_Extra; + ng_FieldGather_noNCI = ng_FieldGather_noNCI.min(ng_alloc_EB); // If NCI filter, add guard cells in the z direction IntVect ng_NCIFilter = IntVect::TheZeroVector(); if (do_fdtd_nci_corr) @@ -131,20 +131,9 @@ guardCellManager::Init( ng_FieldGather = ng_FieldGather_noNCI + ng_NCIFilter; // Guard cells for auxiliary grid - ng_Aux = 2*ng_FieldGather_noNCI + ng_NCIFilter; + ng_UpdateAux = 2*ng_FieldGather_noNCI + ng_NCIFilter; // Make sure we do not exchange more guard cells than allocated. - ng_FieldGather = ng_FieldGather.min(ngE); - ng_Aux = ng_Aux.min(ngE); - - Print()<<"rrr ngE : "<& B, const std::array& dx, int ngrow); - const amrex::IntVect getngE() const { return guard_cells.ngE; }; - const amrex::IntVect getngF() const { return guard_cells.ngF; }; + const amrex::IntVect getngE() const { return guard_cells.ng_alloc_EB; }; + const amrex::IntVect getngF() const { return guard_cells.ng_alloc_F; }; protected: diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index 1915434b6..489692b35 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -717,18 +717,6 @@ WarpX::AllocLevelData (int lev, const BoxArray& ba, const DistributionMapping& d maxwell_fdtd_solver_id, maxLevel()); - IntVect ngE = guard_cells.ngE; - IntVect ngJ = guard_cells.ngJ; - IntVect ngRho = guard_cells.ngRho; - int ngF_int = guard_cells.ngF_int; - IntVect ngextra = guard_cells.ngExtra; - - Print()<<"ngE "<nSpeciesDepositOnMainGrid() && n_current_deposition_buffer == 0) { n_current_deposition_buffer = 1; // This forces the allocation of buffers and allows the code associated @@ -738,14 +726,22 @@ WarpX::AllocLevelData (int lev, const BoxArray& ba, const DistributionMapping& d } if (n_current_deposition_buffer < 0) { - n_current_deposition_buffer = ngJ.max(); + n_current_deposition_buffer = guard_cells.ng_alloc_J.max(); } if (n_field_gather_buffer < 0) { // Field gather buffer should be larger than current deposition buffers n_field_gather_buffer = n_current_deposition_buffer + 1; } - AllocLevelMFs(lev, ba, dm, ngE, ngJ, ngRho, ngF_int, ngextra, aux_is_nodal); + //IntVect ngE = guard_cells.ng_alloc_EB; + //IntVect ngJ = guard_cells.ng_alloc_J; + //IntVect ngRho = guard_cells.ng_alloc_Rho; + //int ngF_int = guard_cells.ng_alloc_F_int; + //IntVect ngextra = guard_cells.ngExtra; + + AllocLevelMFs(lev, ba, dm, guard_cells.ng_alloc_EB, guard_cells.ng_alloc_J, + guard_cells.ng_alloc_Rho, guard_cells.ng_alloc_F_int, + guard_cells.ng_Extra, aux_is_nodal); } void -- cgit v1.2.3 From 82c791de01622ff1556203902c0782e396719629 Mon Sep 17 00:00:00 2001 From: MaxThevenet Date: Thu, 31 Oct 2019 18:25:55 -0700 Subject: cleaning and use better naming --- Source/Evolve/WarpXEvolveEM.cpp | 21 ++++++++++++++------- Source/Parallelization/GuardCellManager.cpp | 10 +++++++++- Source/Utils/WarpXMovingWindow.cpp | 6 ++++-- 3 files changed, 27 insertions(+), 10 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index 2d28641d2..550f7a708 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -194,10 +194,11 @@ WarpX::EvolveEM (int numsteps) // slice gen // if (to_make_plot || do_insitu || to_make_slice_plot) { - // This is probably overkill + // This is probably overkill, but it's not called often FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); - // This is probably overkill + // This is probably overkill, but it's not called often FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); + // This is probably overkill, but it's not called often FillBoundaryAux(guard_cells.ng_UpdateAux); UpdateAuxilaryData(); @@ -249,10 +250,11 @@ WarpX::EvolveEM (int numsteps) if (write_plot_file || do_insitu) { - // This is probably overkill + // This is probably overkill, but it's not called often FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); - // This is probably overkill + // This is probably overkill, but it's not called often FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); + // This is probably overkill FillBoundaryAux(guard_cells.ng_UpdateAux); UpdateAuxilaryData(); @@ -314,6 +316,10 @@ WarpX::OneStep_nosub (Real cur_time) SyncRho(); + // At this point, J is up-to-date inside the domain, and E and B are + // up-to-date including enough guard cells for first step of the field + // solve. + // For extended PML: copy J from regular grid to PML, and damp J in PML if (do_pml && pml_has_particles) CopyJPML(); if (do_pml && do_pml_j_damping) DampJPML(); @@ -326,14 +332,13 @@ WarpX::OneStep_nosub (Real cur_time) FillBoundaryE(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); FillBoundaryB(guard_cells.ng_alloc_EB, guard_cells.ng_Extra); #else - Print()<<"Push fields\n"; EvolveF(0.5*dt[0], DtType::FirstHalf); FillBoundaryF(guard_cells.ng_FieldSolver); EvolveB(0.5*dt[0]); // We now have B^{n+1/2} - // FillBoundaryB(guard_cells.ng_FieldSolver, guard_cells.ngExtra); + FillBoundaryB(guard_cells.ng_FieldSolver, IntVect::TheZeroVector()); EvolveE(dt[0]); // We now have E^{n+1} - // FillBoundaryE(guard_cells.ng_FieldSolver, guard_cells.ngExtra); + FillBoundaryE(guard_cells.ng_FieldSolver, IntVect::TheZeroVector()); EvolveF(0.5*dt[0], DtType::SecondHalf); FillBoundaryF(guard_cells.ng_alloc_F); @@ -341,6 +346,8 @@ WarpX::OneStep_nosub (Real cur_time) if (do_pml) { DampPML(); } + // E and B are up-to-date in the domain, but all guard cells are + // outdated. #endif } diff --git a/Source/Parallelization/GuardCellManager.cpp b/Source/Parallelization/GuardCellManager.cpp index d970cd464..415a13c00 100644 --- a/Source/Parallelization/GuardCellManager.cpp +++ b/Source/Parallelization/GuardCellManager.cpp @@ -107,12 +107,16 @@ guardCellManager::Init( } } ng_alloc_F = IntVect(AMREX_D_DECL(ng_alloc_F_int, ng_alloc_F_int, ng_alloc_F_int)); -#endif +#endif ng_Extra = IntVect(static_cast(aux_is_nodal and !do_nodal)); // Compute number of cells required for Field Solver +#ifdef WARPX_USE_PSATD + ng_FieldSolver = ng_alloc_EB; +#else ng_FieldSolver = IntVect(AMREX_D_DECL(1,1,1)); +#endif // Compute number of cells required for Field Gather int FGcell[4] = {0,1,1,2}; // Index is nox @@ -136,4 +140,8 @@ guardCellManager::Init( // Make sure we do not exchange more guard cells than allocated. ng_FieldGather = ng_FieldGather.min(ng_alloc_EB); ng_UpdateAux = ng_UpdateAux.min(ng_alloc_EB); + // Only FillBoundary(ng_FieldGather) is called between consecutive + // field solves. So ng_FieldGather must have enough cells + // for the field solve too. + ng_FieldGather = ng_FieldGather.max(ng_FieldSolver); } diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index d59d0d8b0..e3819966d 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -223,11 +223,13 @@ WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, // Not sure why this is needed, but it is... IntVect ng_mw = IntVect::TheUnitVector(); + // Enough guard cells in the MW direction ng_mw[dir] = num_shift; + // Add the extra cell (if momentum-conserving gather with staggered field solve) ng_mw += ng_extra; + // Make sure we don't exceed number of guard cells allocated ng_mw = ng_mw.min(ng); - Print()<<"ng_mw "< Date: Tue, 5 Nov 2019 09:28:35 -0800 Subject: use new version of FillBoundary from amrex --- Source/Parallelization/WarpXComm.cpp | 40 ++++++++++++++++++------------------ Source/Utils/WarpXMovingWindow.cpp | 2 +- 2 files changed, 21 insertions(+), 21 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Parallelization/WarpXComm.cpp b/Source/Parallelization/WarpXComm.cpp index ea9408953..c23447e8f 100644 --- a/Source/Parallelization/WarpXComm.cpp +++ b/Source/Parallelization/WarpXComm.cpp @@ -375,9 +375,9 @@ WarpX::FillBoundaryE (int lev, PatchType patch_type, IntVect ng) ng <= Efield_fp[lev][0]->nGrowVect(), "Error: in FillBoundaryE, requested more guard cells than allocated"); Print()<<"FillBoundaryE exchanges "<< Efield_fp[lev][0]->nGrowVect() <<'\n'; - Efield_fp[lev][0]->FillBoundary(0, Efield_fp[lev][0]->nComp(), ng, period); - Efield_fp[lev][1]->FillBoundary(0, Efield_fp[lev][1]->nComp(), ng, period); - Efield_fp[lev][2]->FillBoundary(0, Efield_fp[lev][2]->nComp(), ng, period); + Efield_fp[lev][0]->FillBoundary(ng, period); + Efield_fp[lev][1]->FillBoundary(ng, period); + Efield_fp[lev][2]->FillBoundary(ng, period); } else if (patch_type == PatchType::coarse) { @@ -396,9 +396,9 @@ WarpX::FillBoundaryE (int lev, PatchType patch_type, IntVect ng) ng <= Efield_cp[lev][0]->nGrowVect(), "Error: in FillBoundaryE, requested more guard cells than allocated"); Print()<<"FillBoundaryE exchanges "<< Efield_cp[lev][0]->nGrowVect() <<'\n'; - Efield_cp[lev][0]->FillBoundary(0, Efield_cp[lev][0]->nComp(), ng, cperiod); - Efield_cp[lev][1]->FillBoundary(0, Efield_cp[lev][1]->nComp(), ng, cperiod); - Efield_cp[lev][2]->FillBoundary(0, Efield_cp[lev][2]->nComp(), ng, cperiod); + Efield_cp[lev][0]->FillBoundary(ng, cperiod); + Efield_cp[lev][1]->FillBoundary(ng, cperiod); + Efield_cp[lev][2]->FillBoundary(ng, cperiod); } } @@ -428,9 +428,9 @@ WarpX::FillBoundaryB (int lev, PatchType patch_type, IntVect ng) ng <= Bfield_fp[lev][0]->nGrowVect(), "Error: in FillBoundaryB, requested more guard cells than allocated"); Print()<<"FillBoundaryB exchanges "<< Bfield_fp[lev][0]->nGrowVect() <<'\n'; - Bfield_fp[lev][0]->FillBoundary(0, Bfield_fp[lev][0]->nComp(), ng, period); - Bfield_fp[lev][1]->FillBoundary(0, Bfield_fp[lev][1]->nComp(), ng, period); - Bfield_fp[lev][2]->FillBoundary(0, Bfield_fp[lev][2]->nComp(), ng, period); + Bfield_fp[lev][0]->FillBoundary(ng, period); + Bfield_fp[lev][1]->FillBoundary(ng, period); + Bfield_fp[lev][2]->FillBoundary(ng, period); } else if (patch_type == PatchType::coarse) { @@ -448,9 +448,9 @@ WarpX::FillBoundaryB (int lev, PatchType patch_type, IntVect ng) ng <= Bfield_cp[lev][0]->nGrowVect(), "Error: in FillBoundaryB, requested more guard cells than allocated"); Print()<<"FillBoundaryB exchanges "<< Bfield_cp[lev][0]->nGrowVect() <<'\n'; - Bfield_cp[lev][0]->FillBoundary(0, Bfield_cp[lev][0]->nComp(), ng, cperiod); - Bfield_cp[lev][1]->FillBoundary(0, Bfield_cp[lev][1]->nComp(), ng, cperiod); - Bfield_cp[lev][2]->FillBoundary(0, Bfield_cp[lev][2]->nComp(), ng, cperiod); + Bfield_cp[lev][0]->FillBoundary(ng, cperiod); + Bfield_cp[lev][1]->FillBoundary(ng, cperiod); + Bfield_cp[lev][2]->FillBoundary(ng, cperiod); } } @@ -477,7 +477,7 @@ WarpX::FillBoundaryF (int lev, PatchType patch_type, IntVect ng) AMREX_ALWAYS_ASSERT_WITH_MESSAGE( ng <= F_fp[lev]->nGrowVect(), "Error: in FillBoundaryF, requested more guard cells than allocated"); - F_fp[lev]->FillBoundary(0, F_fp[lev]->nComp(), ng, period); + F_fp[lev]->FillBoundary(ng, period); } else if (patch_type == PatchType::coarse && F_cp[lev]) { @@ -492,7 +492,7 @@ WarpX::FillBoundaryF (int lev, PatchType patch_type, IntVect ng) AMREX_ALWAYS_ASSERT_WITH_MESSAGE( ng <= F_cp[lev]->nGrowVect(), "Error: in FillBoundaryF, requested more guard cells than allocated"); - F_cp[lev]->FillBoundary(0, F_cp[lev]->nComp(), ng, cperiod); + F_cp[lev]->FillBoundary(ng, cperiod); } } @@ -509,12 +509,12 @@ void WarpX::FillBoundaryAux (int lev, IntVect ng) { const auto& period = Geom(lev).periodicity(); - Efield_aux[lev][0]->FillBoundary(0, Efield_aux[lev][0]->nComp(), ng, period); - Efield_aux[lev][1]->FillBoundary(0, Efield_aux[lev][1]->nComp(), ng, period); - Efield_aux[lev][2]->FillBoundary(0, Efield_aux[lev][2]->nComp(), ng, period); - Bfield_aux[lev][0]->FillBoundary(0, Bfield_aux[lev][0]->nComp(), ng, period); - Bfield_aux[lev][1]->FillBoundary(0, Bfield_aux[lev][1]->nComp(), ng, period); - Bfield_aux[lev][2]->FillBoundary(0, Bfield_aux[lev][2]->nComp(), ng, period); + Efield_aux[lev][0]->FillBoundary(ng, period); + Efield_aux[lev][1]->FillBoundary(ng, period); + Efield_aux[lev][2]->FillBoundary(ng, period); + Bfield_aux[lev][0]->FillBoundary(ng, period); + Bfield_aux[lev][1]->FillBoundary(ng, period); + Bfield_aux[lev][2]->FillBoundary(ng, period); } void diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index e3819966d..2c0c6eac9 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -230,7 +230,7 @@ WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, // Make sure we don't exceed number of guard cells allocated ng_mw = ng_mw.min(ng); // Fill guard cells. - tmpmf.FillBoundary(0, tmpmf.nComp(), ng_mw, geom.periodicity()); + tmpmf.FillBoundary(ng_mw, geom.periodicity()); // Make a box that covers the region that the window moved into const IndexType& typ = ba.ixType(); -- cgit v1.2.3 From 72c751d492c44c0bc4d318aeb195e05636c8c38a Mon Sep 17 00:00:00 2001 From: MaxThevenet Date: Fri, 3 Jan 2020 19:30:34 -0800 Subject: include Remi's suggestions --- Source/Evolve/WarpXEvolveEM.cpp | 1 + Source/Parallelization/GuardCellManager.H | 43 +++++++++++++++++++++++++---- Source/Parallelization/GuardCellManager.cpp | 8 ++---- Source/Utils/WarpXMovingWindow.cpp | 3 +- Source/WarpX.H | 6 ++-- Source/WarpX.cpp | 14 +++++----- 6 files changed, 52 insertions(+), 23 deletions(-) (limited to 'Source/Utils/WarpXMovingWindow.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index 366be28d2..e313aab01 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -5,6 +5,7 @@ #include #include #include +#include #ifdef WARPX_USE_PY #include #endif diff --git a/Source/Parallelization/GuardCellManager.H b/Source/Parallelization/GuardCellManager.H index 2e1cebff8..34cf549cf 100644 --- a/Source/Parallelization/GuardCellManager.H +++ b/Source/Parallelization/GuardCellManager.H @@ -3,11 +3,33 @@ #include +/** + * \brief This class computes and stores the number of guard cells needed for + * the allocation of the MultiFabs and required for each part of the PIC loop. + */ class guardCellManager{ public: - int Init( + /** + * \brief Initialize number of guard cells depending on the options used. + * + * \param do_subcycling bool, whether to use subcycling + * \param do_fdtd_nci_corr bool, whether to use Godfrey NCI corrector + * \param do_nodal bool, whether the field solver is nodal + * \param do_moving_window bool, whether to use moving window + * \param do_fft_mpi_dec bool, whether to do parallel FFTs for PSATD + * \param aux_is_nodal bool, true if the aux grid is nodal + * \param moving_window_dir direction of moving window + * \param nox order of current deposition + * \param nox_fft order of PSATD in x direction + * \param noy_fft order of PSATD in y direction + * \param noz_fft order of PSATD in z direction + * \param nci_corr_stencil stencil of NCI corrector + * \param maxwell_fdtd_solver_id if of Maxwell solver + * \param max_level max level of the simulation + */ + void Init( const bool do_subcycling, const bool do_fdtd_nci_corr, const bool do_nodal, @@ -21,21 +43,32 @@ public: const int maxwell_fdtd_solver_id, const int max_level); - // Guard cells allocated for each multifab + // Guard cells allocated for MultiFabs E and B amrex::IntVect ng_alloc_EB = amrex::IntVect::TheZeroVector(); + // Guard cells allocated for MultiFab J amrex::IntVect ng_alloc_J = amrex::IntVect::TheZeroVector(); + // Guard cells allocated for MultiFab Rho amrex::IntVect ng_alloc_Rho = amrex::IntVect::TheZeroVector(); + // Guard cells allocated for MultiFab F amrex::IntVect ng_alloc_F = amrex::IntVect::TheZeroVector(); - int ng_alloc_F_int = 0; - // Guard cells exchanged for specific in the PIC loop + // Guard cells exchanged for specific parts of the PIC loop + + // Number of guard cells of E and B that must exchanged before Field Solver amrex::IntVect ng_FieldSolver = amrex::IntVect::TheZeroVector(); + // Number of guard cells of F that must exchanged before Field Solver amrex::IntVect ng_FieldSolverF = amrex::IntVect::TheZeroVector(); + // Number of guard cells of E and B that must exchanged before Field Gather amrex::IntVect ng_FieldGather = amrex::IntVect::TheZeroVector(); + // Number of guard cells of E and B that must exchanged before updating the Aux grid amrex::IntVect ng_UpdateAux = amrex::IntVect::TheZeroVector(); + // Number of guard cells of all MultiFabs that must exchanged before moving window amrex::IntVect ng_MovingWindow = amrex::IntVect::TheZeroVector(); - // Extra guard cells for fine level of E and B + // When the auxiliary grid is nodal but the field solver is staggered + // (typically with momentum-conserving gather with FDTD Yee solver), + // An extra guard cell is needed on the fine grid to do the interpolation + // for E and B. amrex::IntVect ng_Extra = amrex::IntVect::TheZeroVector(); }; diff --git a/Source/Parallelization/GuardCellManager.cpp b/Source/Parallelization/GuardCellManager.cpp index 34454bd7e..99feca516 100644 --- a/Source/Parallelization/GuardCellManager.cpp +++ b/Source/Parallelization/GuardCellManager.cpp @@ -4,7 +4,7 @@ using namespace amrex; -int +void guardCellManager::Init( const bool do_subcycling, const bool do_fdtd_nci_corr, @@ -69,11 +69,9 @@ guardCellManager::Init( ng_alloc_J = IntVect(ngJx,ngJz); #endif - int nJ_buffer = ng_alloc_J.max(); // guard cells for J required for deposition only. - ng_alloc_Rho = ng_alloc_J+1; //One extra ghost cell, so that it's safe to deposit charge density // after pushing particle. - ng_alloc_F_int = (do_moving_window) ? 2 : 0; + int ng_alloc_F_int = (do_moving_window) ? 2 : 0; // CKC solver requires one additional guard cell if (maxwell_fdtd_solver_id == 1) ng_alloc_F_int = std::max( ng_alloc_F_int, 1 ); ng_alloc_F = IntVect(AMREX_D_DECL(ng_alloc_F_int, ng_alloc_F_int, ng_alloc_F_int)); @@ -155,6 +153,4 @@ guardCellManager::Init( if (do_moving_window){ ng_MovingWindow[moving_window_dir] = 1; } - - return nJ_buffer; } diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index 2c0c6eac9..7044d75d8 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -221,8 +221,7 @@ WarpX::shiftMF (MultiFab& mf, const Geometry& geom, int num_shift, int dir, MultiFab tmpmf(ba, dm, nc, ng); MultiFab::Copy(tmpmf, mf, 0, 0, nc, ng); - // Not sure why this is needed, but it is... - IntVect ng_mw = IntVect::TheUnitVector(); + IntVect ng_mw = IntVect::TheZeroVector(); // Enough guard cells in the MW direction ng_mw[dir] = num_shift; // Add the extra cell (if momentum-conserving gather with staggered field solve) diff --git a/Source/WarpX.H b/Source/WarpX.H index 4b2404198..80d4e3367 100644 --- a/Source/WarpX.H +++ b/Source/WarpX.H @@ -138,7 +138,7 @@ public: // do nodal static int do_nodal; - + const amrex::MultiFab& getcurrent (int lev, int direction) {return *current_fp[lev][direction];} const amrex::MultiFab& getEfield (int lev, int direction) {return *Efield_aux[lev][direction];} const amrex::MultiFab& getBfield (int lev, int direction) {return *Bfield_aux[lev][direction];} @@ -455,8 +455,8 @@ private: void AllocLevelMFs (int lev, const amrex::BoxArray& ba, const amrex::DistributionMapping& dm, const amrex::IntVect& ngE, const amrex::IntVect& ngJ, - const amrex::IntVect& ngRho, int ngF, const amrex::IntVect& ngextra, - const bool aux_is_nodal); + const amrex::IntVect& ngRho, const amrex::IntVect& ngF, + const amrex::IntVect& ngextra, const bool aux_is_nodal); amrex::Vector istep; // which step? amrex::Vector nsubsteps; // how many substeps on each level? diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index 470c8ea7e..c9bd0c5d7 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -720,7 +720,7 @@ WarpX::AllocLevelData (int lev, const BoxArray& ba, const DistributionMapping& d bool aux_is_nodal = (field_gathering_algo == GatheringAlgo::MomentumConserving); - int nJ_buffer = guard_cells.Init( + guard_cells.Init( do_subcycling, WarpX::use_fdtd_nci_corr, do_nodal, @@ -743,7 +743,7 @@ WarpX::AllocLevelData (int lev, const BoxArray& ba, const DistributionMapping& d } if (n_current_deposition_buffer < 0) { - n_current_deposition_buffer = nJ_buffer; + n_current_deposition_buffer = guard_cells.ng_alloc_J.max(); } if (n_field_gather_buffer < 0) { // Field gather buffer should be larger than current deposition buffers @@ -751,14 +751,14 @@ WarpX::AllocLevelData (int lev, const BoxArray& ba, const DistributionMapping& d } AllocLevelMFs(lev, ba, dm, guard_cells.ng_alloc_EB, guard_cells.ng_alloc_J, - guard_cells.ng_alloc_Rho, guard_cells.ng_alloc_F_int, + guard_cells.ng_alloc_Rho, guard_cells.ng_alloc_F, guard_cells.ng_Extra, aux_is_nodal); } void WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm, - const IntVect& ngE, const IntVect& ngJ, const IntVect& ngRho, int ngF, - const IntVect& ngextra, const bool aux_is_nodal) + const IntVect& ngE, const IntVect& ngJ, const IntVect& ngRho, + const IntVect& ngF, const IntVect& ngextra, const bool aux_is_nodal) { #if defined WARPX_DIM_RZ @@ -799,7 +799,7 @@ WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm if (do_dive_cleaning) { - F_fp[lev].reset (new MultiFab(amrex::convert(ba,IntVect::TheUnitVector()),dm,ncomps, ngF)); + F_fp[lev].reset (new MultiFab(amrex::convert(ba,IntVect::TheUnitVector()),dm,ncomps, ngF.max())); } #ifdef WARPX_USE_PSATD else @@ -884,7 +884,7 @@ WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm } if (do_dive_cleaning) { - F_cp[lev].reset (new MultiFab(amrex::convert(cba,IntVect::TheUnitVector()),dm,ncomps, ngF)); + F_cp[lev].reset (new MultiFab(amrex::convert(cba,IntVect::TheUnitVector()),dm,ncomps, ngF.max())); } #ifdef WARPX_USE_PSATD else -- cgit v1.2.3