Fix the order of applying excitation field and boundary conditions

Source/Evolve/WarpXEvolve.cpp

@@ -436,8 +436,6 @@ WarpX::OneStep_nosub (Real cur_time)
 #ifndef WARPX_MAG_LLG
         EvolveB(0.5_rt * dt[0], DtType::FirstHalf); // We now have B^{n+1/2}
         FillBoundaryB(guard_cells.ng_FieldSolver);
-        // ApplyExternalFieldExcitation
-        ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal); // apply B external excitation; soft source to be fixed
 #endif
 
 #ifdef WARPX_MAG_LLG
@@ -468,11 +466,9 @@ WarpX::OneStep_nosub (Real cur_time)
         }
 
         FillBoundaryE(guard_cells.ng_FieldSolver);
-        // ApplyExternalFieldExcitation
-        ApplyExternalFieldExcitationOnGrid(ExternalFieldType::EfieldExternal); // apply E external excitation; soft source to be fixed
 
-        EvolveF(0.5_rt * dt[0], DtType::SecondHalf);
-        EvolveG(0.5_rt * dt[0], DtType::SecondHalf);
+            EvolveF(0.5_rt * dt[0], DtType::SecondHalf);
+            EvolveG(0.5_rt * dt[0], DtType::SecondHalf);
 #ifndef WARPX_MAG_LLG
         EvolveB(0.5_rt * dt[0], DtType::SecondHalf); // We now have B^{n+1}
 
@@ -497,8 +493,6 @@ WarpX::OneStep_nosub (Real cur_time)
             // outdated.
             if (safe_guard_cells) {
                 FillBoundaryB(guard_cells.ng_alloc_EB);
-                // ApplyExternalFieldExcitation
-                ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal); // redundant for hs; need to fix the way to increment ss
             }
 #ifdef WARPX_MAG_LLG
             if (WarpX::em_solver_medium == MediumForEM::Macroscopic) {
@@ -518,8 +512,6 @@ WarpX::OneStep_nosub (Real cur_time)
             if ( safe_guard_cells ){
                 FillBoundaryH(guard_cells.ng_alloc_EB);
                 FillBoundaryM(guard_cells.ng_alloc_EB);
-               // ApplyExternalFieldExcitation
-               ApplyExternalFieldExcitationOnGrid(ExternalFieldType::HfieldExternal); // redundant for hs; need to fix the way to increment ss
             }
 #endif //
     } // !PSATD

Source/FieldSolver/WarpXPushFieldsEM.cpp

@@ -419,14 +419,15 @@ WarpX::PushPSATD ()
     // Evolve the fields in the PML boxes
     for (int lev = 0; lev <= finest_level; ++lev)
     {
-        if (do_pml && pml[lev]->ok())
-        {
-            pml[lev]->PushPSATD(lev);
-        }
         ApplyEfieldBoundary(lev, PatchType::fine);
         if (lev > 0) ApplyEfieldBoundary(lev, PatchType::coarse);
         ApplyBfieldBoundary(lev, PatchType::fine, DtType::FirstHalf);
         if (lev > 0) ApplyBfieldBoundary(lev, PatchType::coarse, DtType::FirstHalf);
+
+        if (do_pml && pml[lev]->ok())
+        {
+            pml[lev]->PushPSATD(lev);
+        }
     }
 #endif
 }
@@ -470,6 +471,11 @@ WarpX::EvolveB (int lev, PatchType patch_type, amrex::Real a_dt, DtType a_dt_typ
                                         m_face_areas[lev], lev, a_dt);
     }
 
+    // Apply external field excitation prior to applying boundary conditions such that excitation does not overwrite B.C.
+    ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal); // apply B external excitation
+
+    ApplyBfieldBoundary(lev, patch_type, a_dt_type);
+
     // Evolve B field in PML cells
     if (do_pml && pml[lev]->ok()) {
         if (patch_type == PatchType::fine) {
@@ -481,7 +487,6 @@ WarpX::EvolveB (int lev, PatchType patch_type, amrex::Real a_dt, DtType a_dt_typ
         }
     }
 
-    ApplyBfieldBoundary(lev, patch_type, a_dt_type);
 }
 
 
@@ -519,6 +524,11 @@ WarpX::EvolveE (int lev, PatchType patch_type, amrex::Real a_dt)
                                        F_cp[lev], lev, a_dt );
     }
 
+    // Apply external field excitation prior to applying boundary conditions such that excitation does not overwrite B.C.
+    ApplyExternalFieldExcitationOnGrid(ExternalFieldType::EfieldExternal); // apply E external excitation
+
+    ApplyEfieldBoundary(lev, patch_type);
+
     // Evolve E field in PML cells
     if (do_pml && pml[lev]->ok()) {
         if (patch_type == PatchType::fine) {
@@ -546,8 +556,6 @@ WarpX::EvolveE (int lev, PatchType patch_type, amrex::Real a_dt)
         }
     }
 
-    ApplyEfieldBoundary(lev, patch_type);
-
 }
 
 
@@ -681,6 +689,11 @@ WarpX::MacroscopicEvolveE (int lev, PatchType patch_type, amrex::Real a_dt) {
         amrex::Abort("Macroscopic EvolveE is not implemented for lev > 0, yet.");
     }
 
+    // Apply external field excitation prior to the boundary condition such that excitation does not overwrite B.C.
+    ApplyExternalFieldExcitationOnGrid(ExternalFieldType::EfieldExternal); // apply E external excitation
+
+    ApplyEfieldBoundary(lev, patch_type);
+
     // Evolve E field in PML cells
     if (do_pml && pml[lev]->ok()) {
         if (patch_type == PatchType::fine) {
@@ -716,7 +729,6 @@ WarpX::MacroscopicEvolveE (int lev, PatchType patch_type, amrex::Real a_dt) {
         }
     }
 
-    ApplyEfieldBoundary(lev, patch_type);
 }
 
 #ifdef WARPX_MAG_LLG
@@ -751,6 +763,9 @@ WarpX::MacroscopicEvolveHM (int lev, PatchType patch_type, amrex::Real a_dt) {
         amrex::Abort("Macroscopic EvolveHM is not implemented for lev > 0 yet");
     }
 
+    // Apply external field excitation prior to applying boundary conditions such that excitation does not overwrite B.C.
+    ApplyExternalFieldExcitationOnGrid(ExternalFieldType::HfieldExternal); // apply H external excitation
+
     // Evolve H field in PML cells
     if (do_pml && pml[lev]->ok()) {
         if (patch_type == PatchType::fine) {
@@ -794,6 +809,9 @@ WarpX::MacroscopicEvolveHM_2nd (int lev, PatchType patch_type, amrex::Real a_dt)
         amrex::Abort("Macroscopic EvolveHM_2nd is not implemented for lev > 0 yet");
     }
 
+    // Apply external field excitation prior to applying boundary conditions such that excitation does not overwrite B.C.
+    ApplyExternalFieldExcitationOnGrid(ExternalFieldType::HfieldExternal); // apply H external excitation
+
     // Evolve H field in PML cells
     if (do_pml && pml[lev]->ok()) {
         if (patch_type == PatchType::fine) {

Source/WarpX.cpp

@@ -1512,12 +1512,11 @@ WarpX::AllocLevelMFs (int lev, const BoxArray& ba, const DistributionMapping& dm
 {
     // Declare nodal flags
     IntVect Ex_nodal_flag, Ey_nodal_flag, Ez_nodal_flag;
+    IntVect Bx_nodal_flag, By_nodal_flag, Bz_nodal_flag;
 #ifdef WARPX_MAG_LLG
     IntVect Mx_nodal_flag, My_nodal_flag, Mz_nodal_flag;
     IntVect Hx_nodal_flag, Hy_nodal_flag, Hz_nodal_flag;
     IntVect Hx_bias_nodal_flag, Hy_bias_nodal_flag, Hz_bias_nodal_flag;
-#else
-    IntVect Bx_nodal_flag, By_nodal_flag, Bz_nodal_flag;
 #endif
     IntVect jx_nodal_flag, jy_nodal_flag, jz_nodal_flag;
     IntVect rho_nodal_flag;