From 6f0fbb9a685717070ffbf363d96a81343890526c Mon Sep 17 00:00:00 2001 From: David Grote Date: Thu, 24 Sep 2020 21:10:05 -0700 Subject: RZ spectral current correction and Galilean (#1216) * Added stub for current correction in RZ spectral solver * Fixed comments in RZ spectral for current correction stub * Modified automated test for Galilean PSATD (#1033) * Impemented current correction in RZ spectral * Implementation Galilean version of RZ spectral solver * For RZ spectral, do forward and backward transform with current correction * Big fix in DivEFunctor.cpp for RZ spectral * Added RZ rho diagnostic for saving the modes * Implemented fft_periodic_single_box for RZ spectral * Moved routines from SpectralSolverRZ.H to .cpp * Added hook for VayDeposition in GalileanPsatdAlgorithmRZ * Bug fix in DivEFunctor * Fixes and cleanup for GalileanPsatdAlgorithmRZ * Fix line spacing in SpectralSolverRZ.H * Fix factor 1/2 in update of Ep_m * Fix factor 1/2 in update of Em_m * Fix sign error in current correction in GalileanPsatdAlgorithmRZ.cpp Co-authored-by: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> * Add Langmuir RZ PSATD test with current correction * Add Galilean tests with/without current correction * For RZ psatd, simplified copy for forward transform * Added GalileanPsatdAlgorithmRZ.cpp to CMakeLists * Minor cleanup in RZ spectral solver * In GalileanPsatdAlgorithmRZ.cpp use member initialization for m_v_galilean Co-authored-by: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> * Added some _rt to GalileanPsatdAlgorithmRZ.cpp Co-authored-by: Olga Shapoval <30510597+oshapoval@users.noreply.github.com> Co-authored-by: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> Co-authored-by: Edoardo Zoni --- .../SpectralSolver/SpectralSolverRZ.cpp | 45 ++++++++++++++++++---- 1 file changed, 38 insertions(+), 7 deletions(-) (limited to 'Source/FieldSolver/SpectralSolver/SpectralSolverRZ.cpp') diff --git a/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.cpp b/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.cpp index c6457a036..38851ef50 100644 --- a/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.cpp +++ b/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.cpp @@ -7,6 +7,7 @@ #include "SpectralKSpaceRZ.H" #include "SpectralSolverRZ.H" #include "SpectralAlgorithms/PsatdAlgorithmRZ.H" +#include "SpectralAlgorithms/GalileanPsatdAlgorithmRZ.H" #include "WarpX.H" #include "Utils/WarpXProfilerWrapper.H" @@ -28,6 +29,7 @@ SpectralSolverRZ::SpectralSolverRZ (amrex::BoxArray const & realspace_ba, amrex::DistributionMapping const & dm, int const n_rz_azimuthal_modes, int const norder_z, bool const nodal, + const amrex::Array& v_galilean, amrex::RealVect const dx, amrex::Real const dt, int const lev) : k_space(realspace_ba, dm, dx) @@ -41,8 +43,15 @@ SpectralSolverRZ::SpectralSolverRZ (amrex::BoxArray const & realspace_ba, // - Select the algorithm depending on the input parameters // Initialize the corresponding coefficients over k space // PML is not supported. - algorithm = std::unique_ptr( - new PsatdAlgorithmRZ(k_space, dm, n_rz_azimuthal_modes, norder_z, nodal, dt)); + if (v_galilean[2] == 0) { + // v_galilean is 0: use standard PSATD algorithm + algorithm = std::unique_ptr( + new PsatdAlgorithmRZ(k_space, dm, n_rz_azimuthal_modes, norder_z, nodal, dt)); + } else { + // Otherwise: use the Galilean algorithm + algorithm = std::unique_ptr( + new GalileanPsatdAlgorithmRZ(k_space, dm, n_rz_azimuthal_modes, norder_z, nodal, v_galilean, dt)); + } // - Initialize arrays for fields in spectral space + FFT plans field_data = SpectralFieldDataRZ(realspace_ba, k_space, dm, @@ -50,15 +59,19 @@ SpectralSolverRZ::SpectralSolverRZ (amrex::BoxArray const & realspace_ba, n_rz_azimuthal_modes, lev); } - +/* \brief Transform the component `i_comp` of MultiFab `field_mf` + * to spectral space, and store the corresponding result internally + * (in the spectral field specified by `field_index`) */ void -SpectralSolverRZ::ForwardTransform (amrex::MultiFab const & field_mf, - int const field_index, +SpectralSolverRZ::ForwardTransform (amrex::MultiFab const & field_mf, int const field_index, int const i_comp) { WARPX_PROFILE("SpectralSolverRZ::ForwardTransform"); field_data.ForwardTransform(field_mf, field_index, i_comp); } +/* \brief Transform the two MultiFabs `field_mf1` and `field_mf2` + * to spectral space, and store the corresponding results internally + * (in the spectral field specified by `field_index1` and `field_index2`) */ void SpectralSolverRZ::ForwardTransform (amrex::MultiFab const & field_mf1, int const field_index1, amrex::MultiFab const & field_mf2, int const field_index2) { @@ -67,14 +80,17 @@ SpectralSolverRZ::ForwardTransform (amrex::MultiFab const & field_mf1, int const field_mf2, field_index2); } +/* \brief Transform spectral field specified by `field_index` back to + * real space, and store it in the component `i_comp` of `field_mf` */ void -SpectralSolverRZ::BackwardTransform (amrex::MultiFab& field_mf, - int const field_index, +SpectralSolverRZ::BackwardTransform (amrex::MultiFab& field_mf, int const field_index, int const i_comp) { WARPX_PROFILE("SpectralSolverRZ::BackwardTransform"); field_data.BackwardTransform(field_mf, field_index, i_comp); } +/* \brief Transform spectral fields specified by `field_index1` and `field_index2` + * back to real space, and store it in `field_mf1` and `field_mf2`*/ void SpectralSolverRZ::BackwardTransform (amrex::MultiFab& field_mf1, int const field_index1, amrex::MultiFab& field_mf2, int const field_index2) { @@ -83,6 +99,7 @@ SpectralSolverRZ::BackwardTransform (amrex::MultiFab& field_mf1, int const field field_mf2, field_index2); } +/* \brief Update the fields in spectral space, over one timestep */ void SpectralSolverRZ::pushSpectralFields () { WARPX_PROFILE("SpectralSolverRZ::pushSpectralFields"); @@ -92,12 +109,26 @@ SpectralSolverRZ::pushSpectralFields () { algorithm->pushSpectralFields(field_data); } +/** + * \brief Public interface to call the member function ComputeSpectralDivE + * of the base class SpectralBaseAlgorithmRZ from objects of class SpectralSolverRZ + */ void SpectralSolverRZ::ComputeSpectralDivE (const std::array,3>& Efield, amrex::MultiFab& divE) { algorithm->ComputeSpectralDivE(field_data, Efield, divE); } +/** + * \brief Public interface to call the virtual function \c CurrentCorrection, + * defined in the base class SpectralBaseAlgorithmRZ and possibly overridden + * by its derived classes (e.g. PsatdAlgorithmRZ), from + * objects of class SpectralSolverRZ through the private unique pointer \c algorithm + * + * \param[in,out] current two-dimensional array of unique pointers to MultiFab + * storing the three components of the current density + * \param[in] rho unique pointer to MultiFab storing the charge density + */ void SpectralSolverRZ::CurrentCorrection (std::array,3>& current, const std::unique_ptr& rho) { -- cgit v1.2.3