Improve Python Analysis of Restart Tests (#2384)

* Check All Particle Species, All Directions

* Check also Particle Momenta

* Add Check on Particle Weight, Fields

* Smaller Field Tolerance (1e-12)

* Single Loop over Particles and Fields

* Cleaning

* Add ds.force_periodicity() for yt-4.0

1 files changed, 31 insertions, 34 deletions

diff --git a/Examples/Tests/restart/analysis_restart.py b/Examples/Tests/restart/analysis_restart.py
index 51efc95a8..ad79ca969 100755
--- a/Examples/Tests/restart/analysis_restart.py
+++ b/Examples/Tests/restart/analysis_restart.py
@@ -6,41 +6,38 @@ import numpy as np
 sys.path.insert(1, '../../../../warpx/Regression/Checksum/')
 import checksumAPI
 
-tolerance = sys.float_info.epsilon
-print('tolerance = ', tolerance)
-
+# Load output data generated after restart
 filename = sys.argv[1]
-ds  = yt.load( filename )
-ad  = ds.all_data()
-xb  = ad['beam',     'particle_position_x'].to_ndarray()
-xe  = ad['plasma_e', 'particle_position_x'].to_ndarray()
-zb  = ad['beam',     'particle_position_z'].to_ndarray()
-ze  = ad['plasma_e', 'particle_position_z'].to_ndarray()
-
-filename = 'orig_' + filename
-ds  = yt.load( filename )
-ad  = ds.all_data()
-xb0 = ad['beam',     'particle_position_x'].to_ndarray()
-xe0 = ad['plasma_e', 'particle_position_x'].to_ndarray()
-zb0 = ad['beam',     'particle_position_z'].to_ndarray()
-ze0 = ad['plasma_e', 'particle_position_z'].to_ndarray()
-
-xb.sort()
-xb0.sort()
-assert(np.max(abs(xb-xb0))<tolerance)
-
-xe.sort()
-xe0.sort()
-assert(np.max(abs(xe-xe0))<tolerance)
-
-zb.sort()
-zb0.sort()
-assert(np.max(abs(zb-zb0))<tolerance)
-
-ze.sort()
-ze0.sort()
-assert(np.max(abs(ze-ze0))<tolerance)
-
+ds_restart = yt.load(filename)
+ds_restart.force_periodicity()
+ad_restart = ds_restart.covering_grid(level = 0,
+                                      left_edge = ds_restart.domain_left_edge,
+                                      dims = ds_restart.domain_dimensions)
+
+# Load output data generated from initial run
+benchmark = 'orig_' + filename
+ds_benchmark = yt.load(benchmark)
+ds_benchmark.force_periodicity()
+ad_benchmark = ds_benchmark.covering_grid(level = 0,
+                                          left_edge = ds_benchmark.domain_left_edge,
+                                          dims = ds_benchmark.domain_dimensions)
+
+# Loop over all fields (all particle species, all particle attributes, all grid fields)
+# and compare output data generated from initial run with output data generated after restart
+tolerance = 1e-12
+print('\ntolerance = {:g}'.format(tolerance))
+print()
+for field in ds_benchmark.field_list:
+    dr = ad_restart[field].squeeze().v
+    db = ad_benchmark[field].squeeze().v
+    error = np.amax(np.abs(dr - db))
+    if (np.amax(np.abs(db)) != 0.):
+        error /= np.amax(np.abs(db))
+    print('field: {}; error = {:g}'.format(field, error))
+    assert(error < tolerance)
+print()
+
+# Check-sum analysis
 filename = sys.argv[1]
 test_name = filename[:-9] # Could also be os.path.split(os.getcwd())[1]
 checksumAPI.evaluate_checksum(test_name, filename)