aboutsummaryrefslogtreecommitdiff
path: root/Source/Diagnostics/FullDiagnostics.cpp
diff options
context:
space:
mode:
authorGravatar Hannah Klion <klion@lbl.gov> 2022-03-13 22:57:14 -0700
committerGravatar GitHub <noreply@github.com> 2022-03-13 22:57:14 -0700
commit7df47510e16b785ff70254ef27aaea6ea9770003 (patch)
tree957629844173373c0b33612c1b80783dd9199b8b /Source/Diagnostics/FullDiagnostics.cpp
parent98d65284abdc1330f86c9fe9a6563a0fba486860 (diff)
downloadWarpX-7df47510e16b785ff70254ef27aaea6ea9770003.tar.gz
WarpX-7df47510e16b785ff70254ef27aaea6ea9770003.tar.zst
WarpX-7df47510e16b785ff70254ef27aaea6ea9770003.zip
Add by-cell averages of particle properties (#2892)
* Add ParticleReductionFunctor * Adds a by-cell particle average diagnostic * Averaged quantity is given by an input parser in (x,y,z,ux,uy,uz) * Coordinates in these parsers are in the WarpX convention * Users specify which species the averages are calculated for. The averages are calculated for each given species individually. * If there are no particles in a cell, the average will be zero. * add documentation for ParticleReductionFunctor * add test for cell-averaged particle diagnostics * set benchmarks for praticle_fields_diags and particle_fields_diags_single_precision * small changes to parameter docs and test script (addressing David's comments) * Apply suggestions from code review * remove unnecessary comments in python test script * remove unused variable * fixed inaccurate comments * use static_cast<int> pattern * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * add comment in ParticleReductionFunctor.H * warpx.serialize_ics -> warpx.serialize_initial_conditions * AMREX_SPACEDIM to WARPX_DIM_3D and WARPX_DIM_XZ * fix compile errors from dimensionality conditionals * fix compile errors * Apply suggestions from code review * clarify docs * add include headers * update docs Co-authored-by: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
Diffstat (limited to 'Source/Diagnostics/FullDiagnostics.cpp')
-rw-r--r--Source/Diagnostics/FullDiagnostics.cpp18
1 files changed, 15 insertions, 3 deletions
diff --git a/Source/Diagnostics/FullDiagnostics.cpp b/Source/Diagnostics/FullDiagnostics.cpp
index fb7a055f6..79791ae8e 100644
--- a/Source/Diagnostics/FullDiagnostics.cpp
+++ b/Source/Diagnostics/FullDiagnostics.cpp
@@ -5,6 +5,7 @@
#include "ComputeDiagFunctors/DivEFunctor.H"
#include "ComputeDiagFunctors/PartPerCellFunctor.H"
#include "ComputeDiagFunctors/PartPerGridFunctor.H"
+#include "ComputeDiagFunctors/ParticleReductionFunctor.H"
#include "ComputeDiagFunctors/RhoFunctor.H"
#include "Diagnostics/Diagnostics.H"
#include "Diagnostics/ParticleDiag/ParticleDiag.H"
@@ -408,10 +409,13 @@ FullDiagnostics::InitializeFieldFunctors (int lev)
// Species index to loop over species that dump rho per species
int i = 0;
- m_all_field_functors[lev].resize( m_varnames.size() );
+ const auto nvar = static_cast<int>(m_varnames_fields.size());
+ const auto nspec = static_cast<int>(m_pfield_species.size());
+ const auto ntot = static_cast<int>(nvar + m_pfield_varnames.size() * nspec);
+
+ m_all_field_functors[lev].resize(ntot);
// Fill vector of functors for all components except individual cylindrical modes.
- const auto n = static_cast<int>(m_all_field_functors[lev].size());
- for (int comp=0; comp<n; comp++){
+ for (int comp=0; comp<nvar; comp++){
if ( m_varnames[comp] == "Ex" ){
m_all_field_functors[lev][comp] = std::make_unique<CellCenterFunctor>(warpx.get_pointer_Efield_aux(lev, 0), lev, m_crse_ratio);
} else if ( m_varnames[comp] == "Ey" ){
@@ -456,6 +460,14 @@ FullDiagnostics::InitializeFieldFunctors (int lev)
amrex::Abort("Error: " + m_varnames[comp] + " is not a known field output type");
}
}
+ // Add functors for average particle data for each species
+ for (int pcomp=0; pcomp<int(m_pfield_varnames.size()); pcomp++) {
+ std::string varname = m_pfield_varnames[pcomp];
+ for (int ispec=0; ispec<int(m_pfield_species.size()); ispec++) {
+ m_all_field_functors[lev][nvar + pcomp * nspec + ispec] = std::make_unique<ParticleReductionFunctor>(nullptr,
+ lev, m_crse_ratio, m_pfield_strings[varname], m_pfield_species_index[ispec]);
+ }
+ }
AddRZModesToDiags( lev );
}