2 files changed, 38 insertions, 17 deletions

diff --git a/Source/FieldSolver/SpectralSolver/SpectralSolver.H b/Source/FieldSolver/SpectralSolver/SpectralSolver.H
index 5266025b6..d4019a9a3 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralSolver.H
+++ b/Source/FieldSolver/SpectralSolver/SpectralSolver.H
@@ -1,8 +1,7 @@
 #ifndef WARPX_SPECTRAL_SOLVER_H_
 #define WARPX_SPECTRAL_SOLVER_H_
 
-#include <SpectralKSpace.H>
-#include <PsatdAlgorithm.H>
+#include <SpectralBaseAlgorithm.H>
 #include <SpectralFieldData.H>
 
 /* \brief Top-level class for the electromagnetic spectral solver
@@ -14,28 +13,17 @@
 class SpectralSolver
 {
     public:
-        // Inline definition of the member functions of `SpectralSolver`
+        // Inline definition of the member functions of `SpectralSolver`,
+        // except the constructor (see `SpectralSolver.cpp`)
         // The body of these functions is short, since the work is done in the
         // underlying classes `SpectralFieldData` and `PsatdAlgorithm`
 
-        /* \brief Initialize the spectral solver */
+        // Constructor
         SpectralSolver( const amrex::BoxArray& realspace_ba,
                         const amrex::DistributionMapping& dm,
                         const int norder_x, const int norder_y,
                         const int norder_z, const bool nodal,
-                        const amrex::RealVect dx, const amrex::Real dt ) {
-            // Initialize all structures using the same distribution mapping dm
-
-            // - Initialize k space object (Contains info about the size of
-            // the spectral space corresponding to each box in `realspace_ba`,
-            // as well as the value of the corresponding k coordinates)
-            const SpectralKSpace k_space= SpectralKSpace(realspace_ba, dm, dx);
-            // - Initialize the algorithm (coefficients) over this space
-            algorithm = std::unique_ptr<PsatdAlgorithm>( new PsatdAlgorithm(
-                    k_space, dm, norder_x, norder_y, norder_z, nodal, dt ) );
-            // - Initialize arrays for fields in Fourier space + FFT plans
-            field_data = SpectralFieldData( realspace_ba, k_space, dm );
-        };
+                        const amrex::RealVect dx, const amrex::Real dt );
 
         /* \brief Transform the component `i_comp` of MultiFab `mf`
          *  to spectral space, and store the corresponding result internally
@@ -59,6 +47,9 @@ class SpectralSolver
         /* \brief Update the fields in spectral space, over one timestep */
         void pushSpectralFields(){
             BL_PROFILE("SpectralSolver::pushSpectralFields");
+            // Virtual function: the actual function used here depends
+            // on the sub-class of `SpectralBaseAlgorithm` that was
+            // initialized in the constructor of `SpectralSolver`
             algorithm->pushSpectralFields( field_data );
         };
 
diff --git a/Source/FieldSolver/SpectralSolver/SpectralSolver.cpp b/Source/FieldSolver/SpectralSolver/SpectralSolver.cpp
new file mode 100644
index 000000000..a04a68d2c
--- /dev/null
+++ b/Source/FieldSolver/SpectralSolver/SpectralSolver.cpp
@@ -0,0 +1,30 @@
+#include <SpectralKSpace.H>
+#include <SpectralSolver.H>
+#include <PsatdAlgorithm.H>
+
+/* \brief TODO Description
+ */
+SpectralSolver( const amrex::BoxArray& realspace_ba,
+                const amrex::DistributionMapping& dm,
+                const int norder_x, const int norder_y,
+                const int norder_z, const bool nodal,
+                const amrex::RealVect dx, const amrex::Real dt ) {
+
+    // Initialize all structures using the same distribution mapping dm
+
+    // - Initialize k space object (Contains info about the size of
+    // the spectral space corresponding to each box in `realspace_ba`,
+    // as well as the value of the corresponding k coordinates)
+    const SpectralKSpace k_space= SpectralKSpace(realspace_ba, dm, dx);
+
+    // - Select the algorithm depending on the input parameters
+    //   Initialize the corresponding coefficients over k space
+    // TODO: Add more algorithms + selection depending on input parameters
+    //       For the moment, this only uses the standard PsatdAlgorithm
+    algorithm = std::unique_ptr<PsatdAlgorithm>( new PsatdAlgorithm(
+            k_space, dm, norder_x, norder_y, norder_z, nodal, dt ) );
+
+    // - Initialize arrays for fields in spectral space + FFT plans
+    field_data = SpectralFieldData( realspace_ba, k_space, dm );
+
+};