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 <ezoni@lbl.gov>

1 files changed, 20 insertions, 1 deletions

diff --git a/Examples/Tests/Langmuir/analysis_langmuir_multi_rz.py b/Examples/Tests/Langmuir/analysis_langmuir_multi_rz.py
index e98c24d7e..9e3942803 100755
--- a/Examples/Tests/Langmuir/analysis_langmuir_multi_rz.py
+++ b/Examples/Tests/Langmuir/analysis_langmuir_multi_rz.py
@@ -13,6 +13,7 @@
 # $$ E_r = -\partial_r \phi = \epsilon \,\frac{mc^2}{e}\frac{2\,r}{w_0^2} \exp\left(-\frac{r^2}{w_0^2}\right) \sin(k_0 z) \sin(\omega_p t)
 # $$ E_z = -\partial_z \phi = - \epsilon \,\frac{mc^2}{e} k_0 \exp\left(-\frac{r^2}{w_0^2}\right) \cos(k_0 z) \sin(\omega_p t)
 import sys
+import re
 import matplotlib
 matplotlib.use('Agg')
 import matplotlib.pyplot as plt
@@ -28,6 +29,9 @@ fn = sys.argv[1]
 
 test_name = fn[:-9] # Could also be os.path.split(os.getcwd())[1]
 
+# Parse test name and check if current correction (psatd.current_correction) is applied
+current_correction = True if re.search('current_correction', fn) else False
+
 # Parameters (these parameters must match the parameters in `inputs.multi.rz.rt`)
 epsilon = 0.01
 n = 2.e24
@@ -102,11 +106,26 @@ plt.tight_layout()
 plt.savefig(test_name+'_analysis.png')
 
 error_rel = overall_max_error
-tolerance_rel = 0.04
+
+if current_correction:
+   tolerance_rel = 0.06
+else:
+   tolerance_rel = 0.04
 
 print("error_rel    : " + str(error_rel))
 print("tolerance_rel: " + str(tolerance_rel))
 
 assert( error_rel < tolerance_rel )
 
+# Check charge conservation (relative L-infinity norm of error) with current correction
+if current_correction:
+    divE = data['divE'].to_ndarray()
+    rho  = data['rho' ].to_ndarray() / epsilon_0
+    error_rel = np.amax(np.abs(divE - rho)) / max(np.amax(divE), np.amax(rho))
+    tolerance = 1.e-9
+    print("Check charge conservation:")
+    print("error_rel = {}".format(error_rel))
+    print("tolerance = {}".format(tolerance))
+    assert( error_rel < tolerance )
+
 checksumAPI.evaluate_checksum(test_name, fn)