diff --git a/Source/ablastr/fields/PoissonSolver.H b/Source/ablastr/fields/PoissonSolver.H
index 727280d630b..cc2f7b456fb 100644
--- a/Source/ablastr/fields/PoissonSolver.H
+++ b/Source/ablastr/fields/PoissonSolver.H
@@ -275,7 +275,7 @@ computePhi (
 #endif
         // Use the Multigrid (MLMG) solver if selected or on refined patches
         // but first scale rho appropriately
-        rho[lev]->mult(-1._rt / ablastr::constant::SI::ep0);  // TODO: when do we "un-multiply" this? We need to document this side-effect!
+        rho[lev]->mult(-1._rt / ablastr::constant::SI::ep0);
 
 #ifdef WARPX_DIM_RZ
         constexpr bool is_rz = true;
@@ -409,6 +409,8 @@ computePhi (
                 post_phi_calculation.value()(mlmg, lev);
             }
         }
+        rho[lev]->mult(-ablastr::constant::SI::ep0);  // Divide rho by epsilon again
+
     } // loop over lev(els)
 } // computePhi
 } // namespace ablastr::fields