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, 35 insertions, 6 deletions

diff --git a/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H b/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H
index 7b917284b..b7f349bb0 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H
+++ b/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H
@@ -41,6 +41,8 @@ class PsatdAlgorithm : public SpectralBaseAlgorithm
          * \param[in] dt time step of the simulation
          * \param[in] update_with_rho whether the update equation for E uses rho or not
          * \param[in] time_averaging whether to use time averaging for large time steps
+         * \param[in] J_linear_in_time whether to use two currents computed at the beginning and the end
+         *            of the time interval (instead of using one current computed at half time)
          */
         PsatdAlgorithm (
             const SpectralKSpace& spectral_kspace,
@@ -52,7 +54,8 @@ class PsatdAlgorithm : public SpectralBaseAlgorithm
             const amrex::Array<amrex::Real,3>& v_galilean,
             const amrex::Real dt,
             const bool update_with_rho,
-            const bool time_averaging);
+            const bool time_averaging,
+            const bool J_linear_in_time);
 
         /**
          * \brief Updates the E and B fields in spectral space, according to the relevant PSATD equations
@@ -66,12 +69,22 @@ class PsatdAlgorithm : public SpectralBaseAlgorithm
          */
         virtual int getRequiredNumberOfFields () const override final
         {
-            if (m_time_averaging) {
-                return SpectralAvgFieldIndex::n_fields;
-            } else {
-                return SpectralFieldIndex::n_fields;
+            if (m_J_linear_in_time)
+            {
+                return SpectralFieldIndexJLinearInTime::n_fields;
             }
-        }
+            else
+            {
+                if (m_time_averaging)
+                {
+                    return SpectralFieldIndexTimeAveraging::n_fields;
+                }
+                else
+                {
+                    return SpectralFieldIndex::n_fields;
+                }
+            }
+        };
 
         /**
          * \brief Initializes the coefficients used in \c pushSpectralFields to update the E and B fields
@@ -86,6 +99,19 @@ class PsatdAlgorithm : public SpectralBaseAlgorithm
             const amrex::Real dt);
 
         /**
+         * \brief Initialize additional coefficients used in \c pushSpectralFields to update E,B,
+         *        required only when using time averaging with the assumption that J is linear in time
+         *
+         * \param[in] spectral_kspace spectral space
+         * \param[in] dm distribution mapping
+         * \param[in] dt time step of the simulation
+         */
+        void InitializeSpectralCoefficientsAvgLin (
+            const SpectralKSpace& spectral_kspace,
+            const amrex::DistributionMapping& dm,
+            const amrex::Real dt);
+
+        /**
          * \brief Initializes additional coefficients used in \c pushSpectralFields to update the E and B fields,
          *        required only when using time averaging with large time steps
          *
@@ -138,6 +164,8 @@ class PsatdAlgorithm : public SpectralBaseAlgorithm
         SpectralRealCoefficients C_coef, S_ck_coef;
         SpectralComplexCoefficients T2_coef, X1_coef, X2_coef, X3_coef, X4_coef;
 
+        SpectralComplexCoefficients X5_coef, X6_coef;
+
         // These real and complex coefficients are allocated only with averaged Galilean PSATD
         SpectralComplexCoefficients Psi1_coef, Psi2_coef, Y1_coef, Y2_coef, Y3_coef, Y4_coef;
 
@@ -153,6 +181,7 @@ class PsatdAlgorithm : public SpectralBaseAlgorithm
         amrex::Real m_dt;
         bool m_update_with_rho;
         bool m_time_averaging;
+        bool m_J_linear_in_time;
         bool m_is_galilean;
 };
 #endif // WARPX_USE_PSATD