Particle Buffer for Backtransformed Diagnostics (#1898)

* move BTD call before Redistribute

* define particle buffer, BTD particle functor, particle filter

* select particles on the slice

* add particle functors

* add Lorentz Transform

* fix conflict

* add call to loretnz operator

* storing Particles in buffer

* This is WIP. Added particle buffers, filled them, sent them for flush with number already flushed. Does not work for multiple flushes. Crashes for OpenPMD.

* trailing endif after rebase

* adding print statements and not flushing particles in OPENPMD if numpart is 0

* last timestep flush is ensured at the end of evolve loop

* fix bug in declaring uy uy new and computing uzp

* set particle Geom, BA, and DMAP for particle flush with plotfile

* set Particle BA Geom DMAP for particle buffer and no BTD transform for force flush

* separate compute and pack from flush

* Fix Typo: resizeable -> resizable

Fixes HDF5 BTD particle output.

* new class for plotfile particles for BTD

* copy particle_H and DATA and Header. some WIP print statements

* Merge plotfile

* clean print statements

* fix warning message

* struct declaration in header, fix warning

* doxygen comments and copyright

* clean print statements

* fix eol and override function warning

* tile data

* fix output species array size bug

* fix access for particle buffer size

* clean and move time-update

* add cur_time update back

* remove cur time update which was called twice

* dont access particles flushed already for full diagnostics

* cur time must be updated for RigidInjection BTD CI test to pass

* temporarily move call to BTD

* updating time and calling BTD before movewindow

* cleanup

* reset benchmarks and analysis script

* clean and add comments

* fix particle box array, geom, dmap

* reset benchmarks for multi_J rz and ElectrostaticSphereEB_mixedBC

* wip commit

* wip commit

* add SI conversion

* abort for openpmd bp backend if species is selected. Also write particle output for BTD only if write_species is 1

* add documentation for aborting if adios is used with openpmd and add other BTD input parameters

* Apply suggestions from code review

commit Axel's suggestions from review

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>

* use bool instead of int

* fix doxygen format

* using h5 as backend in example test to ensure consistency with abort for particle output.

* fix doxygen comment

* reset benchmark again for comoving_2d_psatd galilean_2d_psatd multi_J_rz_psatd

* reset benchmark for background_mcc

* self-review suggestions

* reset benchmarks. Update with last snapshot full info

* Axel's PR suggestions

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>

* Axel's doxygen fix

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>

* add comments

* fix eol

* improved exception handling for stringsteam

* PR suggestions

* Axels' suggestions from code review

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>

* Axel's suggestions :)

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>

* simplify logic

* suggestions from review (Axel/Reva)

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>

* variable name change for clarity

* if num particles in tmp array is 0, return

* Use new BTD inputs to set up BTD for particles in the corresponding particle container

* unused var

* fix logic error

* speciesID undefined

* separate particle and field buffer calls and initialization for BTD. Data common to both are initialized separately

* rename variable so it does not ghost existing varname

* add more comments

* Assert that fields are on

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* redundant function definition

* unused variable

* unused variable zp

* 1D 2D 3D definition

* fix else

* Apply suggestions from code review

Add Prabhat's suggestion

Co-authored-by: Prabhat Kumar <89051199+prkkumar@users.noreply.github.com>

* missing semicolon and ignore xp yp for 1D

* resetting benchmarks for boosted sims and mcc sim

* temporarily changing tolerance since the relative difference for momentum_z is 3.68e-3 and the current tolerance is 2.5e-3

Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Co-authored-by: Prabhat Kumar <89051199+prkkumar@users.noreply.github.com>

1 files changed, 27 insertions, 10 deletions

diff --git a/Source/Diagnostics/FullDiagnostics.cpp b/Source/Diagnostics/FullDiagnostics.cpp
index 6b057f574..81f04a330 100644
--- a/Source/Diagnostics/FullDiagnostics.cpp
+++ b/Source/Diagnostics/FullDiagnostics.cpp
@@ -62,9 +62,12 @@ FullDiagnostics::InitializeParticleBuffer ()
         }
     }
     // Initialize one ParticleDiag per species requested
-    for (auto const& species : m_output_species_names){
-        const int idx = mpc.getSpeciesID(species);
-        m_output_species.push_back(ParticleDiag(m_diag_name, species, mpc.GetParticleContainerPtr(idx)));
+    m_output_species.resize(m_num_buffers);
+    for (int i_buffer = 0; i_buffer < m_num_buffers; ++i_buffer) {
+        for (auto const& species : m_output_species_names){
+            const int idx = mpc.getSpeciesID(species);
+            m_output_species[i_buffer].push_back(ParticleDiag(m_diag_name, species, mpc.GetParticleContainerPtr(idx)));
+        }
     }
 }
 
@@ -123,8 +126,8 @@ FullDiagnostics::Flush ( int i_buffer )
 
     m_flush_format->WriteToFile(
         m_varnames, m_mf_output[i_buffer], m_geom_output[i_buffer], warpx.getistep(),
-        warpx.gett_new(0), m_output_species, nlev_output, m_file_prefix, m_file_min_digits,
-        m_plot_raw_fields, m_plot_raw_fields_guards);
+        warpx.gett_new(0), m_output_species[i_buffer], nlev_output, m_file_prefix,
+        m_file_min_digits, m_plot_raw_fields, m_plot_raw_fields_guards);
 
     FlushRaw();
 }
@@ -274,7 +277,7 @@ FullDiagnostics::AddRZModesToOutputNames (const std::string& field, int ncomp){
 
 
 void
-FullDiagnostics::InitializeFieldBufferData (int i_buffer, int lev ) {
+FullDiagnostics::InitializeBufferData (int i_buffer, int lev ) {
     auto & warpx = WarpX::GetInstance();
     amrex::RealBox diag_dom;
     bool use_warpxba = true;
@@ -468,8 +471,15 @@ FullDiagnostics::PrepareFieldDataForOutput ()
     warpx.FillBoundaryAux(warpx.getngUpdateAux());
 
     // Update the RealBox used for the geometry filter in particle diags
-    for (int i = 0; i < m_output_species.size(); ++i) {
-        m_output_species[i].m_diag_domain = m_geom_output[0][0].ProbDomain();
+    // Note that full diagnostics every diag has only one buffer. (m_num_buffers = 1).
+    // For m_geom_output[i_buffer][lev], the first element is the buffer index, and
+    // second is level=0
+    // The level is set to 0, because the whole physical domain of the simulation is used
+    // to set the domain dimensions for the output particle container.
+    for (int i_buffer = 0; i_buffer < m_num_buffers; ++i_buffer) {
+        for (int i = 0; i < m_output_species[i_buffer].size(); ++i) {
+            m_output_species[i_buffer][i].m_diag_domain = m_geom_output[i_buffer][0].ProbDomain();
+        }
     }
 }
 
@@ -481,6 +491,10 @@ FullDiagnostics::MovingWindowAndGalileanDomainShift (int step)
     // Account for galilean shift
     amrex::Real new_lo[AMREX_SPACEDIM];
     amrex::Real new_hi[AMREX_SPACEDIM];
+    // Note that Full diagnostics has only one snapshot, m_num_buffers = 1
+    // m_geom_output[i_buffer][lev] below have values 0 and 0, respectively, because
+    // we need the physical extent from mesh-refinement level = 0,
+    // and only for the 0th snapshot, since full diagnostics has only one snapshot.
     const amrex::Real* current_lo = m_geom_output[0][0].ProbLo();
     const amrex::Real* current_hi = m_geom_output[0][0].ProbHi();
 
@@ -502,11 +516,12 @@ FullDiagnostics::MovingWindowAndGalileanDomainShift (int step)
         new_hi[0] = current_hi[0] + warpx.m_galilean_shift[2];
     }
 #endif
-    // Update RealBox of geometry with galilean-shifted boundary
+    // Update RealBox of geometry with galilean-shifted boundary.
     for (int lev = 0; lev < nmax_lev; ++lev) {
+        // Note that Full diagnostics has only one snapshot, m_num_buffers = 1
+        // Thus here we set the prob domain for the 0th snapshot only.
         m_geom_output[0][lev].ProbDomain( amrex::RealBox(new_lo, new_hi) );
     }
-
     // For Moving Window Shift
     if (warpx.moving_window_active(step+1)) {
         int moving_dir = warpx.moving_window_dir;
@@ -527,6 +542,8 @@ FullDiagnostics::MovingWindowAndGalileanDomainShift (int step)
         new_hi[moving_dir] = cur_hi[moving_dir] + num_shift_base*geom_dx[moving_dir];
         // Update RealBox of geometry with shifted domain geometry for moving-window
         for (int lev = 0; lev < nmax_lev; ++lev) {
+            // Note that Full diagnostics has only one snapshot, m_num_buffers = 1
+            // Thus here we set the prob domain for the 0th snapshot only.
             m_geom_output[0][lev].ProbDomain( amrex::RealBox(new_lo, new_hi) );
         }
     }