Multi-J scheme (#1828)

* Introduce new option skip_deposition

* Properly implement the option to skip deposition

* Skip deposition for electrostatic solver

* Correct typo

* Add Index Enumerator and Equations for F/G Without Averaging

* Define new functions for current deposition and charge deposition

* Disable interpolation between levels

* Correct compilation error in RZ mode

* Add argument for relative time

* Add Index Enumerator and Equations for F/G With Averaging

* [skip ci] Add new OneStep function

* Fix compilation errors

* Correct more compilation errors

* [skip ci] Fix compilation

* Split the PSATD push into separate functions

* Add guards for rho field

* [skip ci] Use new functions in OneStep

* [skip ci] Separate the inverse transform of E_avg, B_avg

* Add deposition of rho

* [skip ci] Prevent deposition of rho if unallocated

* Fix error in deposition function

* Add subcycling of current deposition

* [skip ci] Add inverse transform of averaged fields

* [skip ci] Move component of rho

* Add function to copy J

* Temporarily deactivate contribution from F

* [skip ci] Implement call to linear in J

* Add psatd time averaging for multiJ

* [skip ci] Fix implementation of averaging

* [skip ci] Implement divE cleaning

* Fix Bug for RZ Builds

* Fix Bug for RZ Builds

* Fix Bug in Init of PML Spectral Solvers

* Cleaning

* Cleaning

* Add div(B) Cleaning (G Equation)

* Multi-J Not Implemented with Galilean PSATD or PMLs

* Add 2D CI Test Using Multi-J Scheme

* Add More Inline Comments

* Add More Inline Comments & Doxygen

* Add Doxygen for Constructor of SpectralSolver

* More Doxygen in Class SpectralSolver

* Add Doxygen for New Functions in WarpXPushFieldsEM.cpp

* Add Doxygen for New Functions in WarpX/MultiParticleContainer

* do_dive/b_cleaning Must Be True With linear_in_J Option

* Cast multij_n_depose to Real in Divisions

* New Input Syntax

* Add const where Possible, Fix Warnings

* Docs for New Input Syntax, Fix Abort Messages

* Consistent Use of Idx, IdxAvg, IdxLin

* Improve Documentation of psatd.J_linear_in_time

* Use const Type Qualifier whenever Possible

* Simplify Initialization of Pointer ion_lev

* Improve Doxygen

* Update documentation

* Add Note on NCI to Docs

* Make warpx.do_multi_J_n_depositions Not Optional

* Simplify Logic in getRequiredNumberOfFields

* Use More const Type Qualifiers

Co-authored-by: Edoardo Zoni <ezoni@lbl.gov>

1 files changed, 66 insertions, 0 deletions

diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp
index 4990f7c6d..b7d34ec75 100644
--- a/Source/Particles/MultiParticleContainer.cpp
+++ b/Source/Particles/MultiParticleContainer.cpp
@@ -421,6 +421,72 @@ MultiParticleContainer::GetZeroChargeDensity (const int lev)
     return zero_rho;
 }
 
+void
+MultiParticleContainer::DepositCurrent (
+    amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3 > >& J,
+    const amrex::Real dt, const amrex::Real relative_t)
+{
+    // Reset the J arrays
+    for (int lev = 0; lev < J.size(); ++lev)
+    {
+        J[lev][0]->setVal(0.0, J[lev][0]->nGrowVect());
+        J[lev][1]->setVal(0.0, J[lev][1]->nGrowVect());
+        J[lev][2]->setVal(0.0, J[lev][2]->nGrowVect());
+    }
+
+    // Call the deposition kernel for each species
+    for (int ispecies = 0; ispecies < nSpecies(); ispecies++)
+    {
+        WarpXParticleContainer& species = GetParticleContainer(ispecies);
+        species.DepositCurrent(J, dt, relative_t);
+    }
+
+#ifdef WARPX_DIM_RZ
+    for (int lev = 0; lev < J.size(); ++lev)
+    {
+        WarpX::GetInstance().ApplyInverseVolumeScalingToCurrentDensity(J[lev][0].get(), J[lev][1].get(), J[lev][2].get(), lev);
+    }
+#endif
+}
+
+void
+MultiParticleContainer::DepositCharge (
+    amrex::Vector<std::unique_ptr<amrex::MultiFab> >& rho,
+    const amrex::Real relative_t, const int icomp)
+{
+    // Reset the rho array
+    for (int lev = 0; lev < rho.size(); ++lev)
+    {
+        int const nc = WarpX::ncomps;
+        rho[lev]->setVal(0.0, icomp*nc, nc, rho[lev]->nGrowVect());
+    }
+
+    // Push the particles in time, if needed
+    if (relative_t != 0.) PushX(relative_t);
+
+    // Call the deposition kernel for each species
+    for (int ispecies = 0; ispecies < nSpecies(); ispecies++)
+    {
+        WarpXParticleContainer& species = GetParticleContainer(ispecies);
+        bool const local = true;
+        bool const reset = false;
+        bool const do_rz_volume_scaling = false;
+        bool const interpolate_across_levels = false;
+        species.DepositCharge(rho, local, reset, do_rz_volume_scaling,
+                              interpolate_across_levels, icomp);
+    }
+
+    // Push the particles back in time
+    if (relative_t != 0.) PushX(-relative_t);
+
+#ifdef WARPX_DIM_RZ
+    for (int lev = 0; lev < rho.size(); ++lev)
+    {
+        WarpX::GetInstance().ApplyInverseVolumeScalingToChargeDensity(rho[lev].get(), lev);
+    }
+#endif
+}
+
 std::unique_ptr<MultiFab>
 MultiParticleContainer::GetChargeDensity (int lev, bool local)
 {