aboutsummaryrefslogtreecommitdiff
path: root/Source
diff options
context:
space:
mode:
Diffstat (limited to 'Source')
-rw-r--r--Source/Diagnostics/ReducedDiags/BeamRelevant.cpp9
-rw-r--r--Source/Diagnostics/ReducedDiags/CMakeLists.txt1
-rw-r--r--Source/Diagnostics/ReducedDiags/FieldEnergy.cpp7
-rw-r--r--Source/Diagnostics/ReducedDiags/FieldMaximum.cpp7
-rw-r--r--Source/Diagnostics/ReducedDiags/LoadBalanceCosts.cpp8
-rw-r--r--Source/Diagnostics/ReducedDiags/Make.package1
-rw-r--r--Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp6
-rw-r--r--Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp7
-rw-r--r--Source/Diagnostics/ReducedDiags/ParticleHistogram.cpp5
-rw-r--r--Source/Diagnostics/ReducedDiags/ParticleNumber.cpp13
-rw-r--r--Source/Diagnostics/ReducedDiags/RhoMaximum.H40
-rw-r--r--Source/Diagnostics/ReducedDiags/RhoMaximum.cpp170
12 files changed, 244 insertions, 30 deletions
diff --git a/Source/Diagnostics/ReducedDiags/BeamRelevant.cpp b/Source/Diagnostics/ReducedDiags/BeamRelevant.cpp
index 57ddc64c5..fa7ff377f 100644
--- a/Source/Diagnostics/ReducedDiags/BeamRelevant.cpp
+++ b/Source/Diagnostics/ReducedDiags/BeamRelevant.cpp
@@ -37,7 +37,7 @@ BeamRelevant::BeamRelevant (std::string rd_name)
// 13: rms gamma
// 14,15,16: emittance x,y,z
// 17: charge
- m_data.resize(18,0.0);
+ m_data.resize(18, 0.0_rt);
#elif (defined WARPX_DIM_XZ)
// 0, 1: mean x,z
// 2, 3, 4: mean px,py,pz
@@ -47,7 +47,7 @@ BeamRelevant::BeamRelevant (std::string rd_name)
// 11: rms gamma
// 12,13: emittance x,z
// 14: charge
- m_data.resize(15,0.0);
+ m_data.resize(15, 0.0_rt);
#endif
if (ParallelDescriptor::IOProcessor())
@@ -55,8 +55,7 @@ BeamRelevant::BeamRelevant (std::string rd_name)
if ( m_IsNotRestart )
{
// open file
- std::ofstream ofs{m_path + m_rd_name + "." + m_extension,
- std::ofstream::out | std::ofstream::app};
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
// write header row
#if (defined WARPX_DIM_3D || defined WARPX_DIM_RZ)
ofs << "#";
@@ -161,7 +160,7 @@ void BeamRelevant::ComputeDiags (int step)
if (w_sum < std::numeric_limits<Real>::min() )
{
for (int i = 0; i < static_cast<int>(m_data.size()); ++i){
- m_data[i] = 0.0;
+ m_data[i] = 0.0_rt;
}
return;
}
diff --git a/Source/Diagnostics/ReducedDiags/CMakeLists.txt b/Source/Diagnostics/ReducedDiags/CMakeLists.txt
index 08276f254..070b6815b 100644
--- a/Source/Diagnostics/ReducedDiags/CMakeLists.txt
+++ b/Source/Diagnostics/ReducedDiags/CMakeLists.txt
@@ -9,5 +9,6 @@ target_sources(WarpX
ReducedDiags.cpp
FieldMaximum.cpp
ParticleExtrema.cpp
+ RhoMaximum.cpp
ParticleNumber.cpp
)
diff --git a/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp b/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp
index 892ea9c7c..cf02e08fe 100644
--- a/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp
+++ b/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp
@@ -37,15 +37,14 @@ FieldEnergy::FieldEnergy (std::string rd_name)
constexpr int noutputs = 3; // total energy, E-field energy and B-field energy
// resize data array
- m_data.resize(noutputs*nLevel,0.0);
+ m_data.resize(noutputs*nLevel, 0.0_rt);
if (ParallelDescriptor::IOProcessor())
{
if ( m_IsNotRestart )
{
// open file
- std::ofstream ofs{m_path + m_rd_name + "." + m_extension,
- std::ofstream::out | std::ofstream::app};
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
// write header row
ofs << "#";
ofs << "[1]step()";
@@ -82,7 +81,7 @@ void FieldEnergy::ComputeDiags (int step)
// Judge if the diags should be done
if (!m_intervals.contains(step+1)) { return; }
- // get WarpX class object
+ // get a reference to WarpX instance
auto & warpx = WarpX::GetInstance();
// get number of level
diff --git a/Source/Diagnostics/ReducedDiags/FieldMaximum.cpp b/Source/Diagnostics/ReducedDiags/FieldMaximum.cpp
index a7250c5aa..ad2ee1a05 100644
--- a/Source/Diagnostics/ReducedDiags/FieldMaximum.cpp
+++ b/Source/Diagnostics/ReducedDiags/FieldMaximum.cpp
@@ -30,15 +30,14 @@ FieldMaximum::FieldMaximum (std::string rd_name)
constexpr int noutputs = 8; // total energy, E-field energy and B-field energy
// resize data array
- m_data.resize(noutputs*nLevel,0.0); // max of Ex,Ey,Ez,|E|,Bx,By,Bz and |B|
+ m_data.resize(noutputs*nLevel, 0.0_rt); // max of Ex,Ey,Ez,|E|,Bx,By,Bz and |B|
if (ParallelDescriptor::IOProcessor())
{
if ( m_IsNotRestart )
{
// open file
- std::ofstream ofs{m_path + m_rd_name + "." + m_extension,
- std::ofstream::out | std::ofstream::app};
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
// write header row
ofs << "#";
ofs << "[1]step()";
@@ -94,7 +93,7 @@ void FieldMaximum::ComputeDiags (int step)
// Judge if the diags should be done
if (!m_intervals.contains(step+1)) { return; }
- // get WarpX class object
+ // get a reference to WarpX instance
auto & warpx = WarpX::GetInstance();
// get number of level
diff --git a/Source/Diagnostics/ReducedDiags/LoadBalanceCosts.cpp b/Source/Diagnostics/ReducedDiags/LoadBalanceCosts.cpp
index 048aa7dab..3ec0f0d45 100644
--- a/Source/Diagnostics/ReducedDiags/LoadBalanceCosts.cpp
+++ b/Source/Diagnostics/ReducedDiags/LoadBalanceCosts.cpp
@@ -23,7 +23,7 @@ LoadBalanceCosts::LoadBalanceCosts (std::string rd_name)
// function that gathers costs
void LoadBalanceCosts::ComputeDiags (int step)
{
- // get WarpX class object
+ // get a reference to WarpX instance
auto& warpx = WarpX::GetInstance();
const amrex::LayoutData<amrex::Real>* cost = warpx.getCosts(0);
@@ -50,8 +50,8 @@ void LoadBalanceCosts::ComputeDiags (int step)
const size_t dataSize =
static_cast<size_t>(m_nDataFields)*
static_cast<size_t>(nBoxes);
- m_data.resize(dataSize, 0.0);
- m_data.assign(dataSize, 0.0);
+ m_data.resize(dataSize, 0.0_rt);
+ m_data.assign(dataSize, 0.0_rt);
// read in WarpX costs to local copy; compute if using `Heuristic` update
amrex::Vector<std::unique_ptr<amrex::LayoutData<amrex::Real> > > costs;
@@ -218,7 +218,7 @@ void LoadBalanceCosts::WriteToFile (int step) const
ofs.close();
- // get WarpX class object
+ // get a reference to WarpX instance
auto& warpx = WarpX::GetInstance();
if (!ParallelDescriptor::IOProcessor()) return;
diff --git a/Source/Diagnostics/ReducedDiags/Make.package b/Source/Diagnostics/ReducedDiags/Make.package
index 3a8361693..22ba3dd34 100644
--- a/Source/Diagnostics/ReducedDiags/Make.package
+++ b/Source/Diagnostics/ReducedDiags/Make.package
@@ -7,6 +7,7 @@ CEXE_sources += LoadBalanceCosts.cpp
CEXE_sources += ParticleHistogram.cpp
CEXE_sources += FieldMaximum.cpp
CEXE_sources += ParticleExtrema.cpp
+CEXE_sources += RhoMaximum.cpp
CEXE_sources += ParticleNumber.cpp
VPATH_LOCATIONS += $(WARPX_HOME)/Source/Diagnostics/ReducedDiags
diff --git a/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp b/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp
index 9c67cea9e..c7ce7b18c 100644
--- a/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp
+++ b/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp
@@ -12,6 +12,7 @@
#include "ParticleExtrema.H"
#include "FieldEnergy.H"
#include "FieldMaximum.H"
+#include "RhoMaximum.H"
#include "ParticleNumber.H"
#include "MultiReducedDiags.H"
@@ -62,6 +63,11 @@ MultiReducedDiags::MultiReducedDiags ()
m_multi_rd[i_rd] =
std::make_unique<FieldMaximum>(m_rd_names[i_rd]);
}
+ else if (rd_type.compare("RhoMaximum") == 0)
+ {
+ m_multi_rd[i_rd] =
+ std::make_unique<RhoMaximum>(m_rd_names[i_rd]);
+ }
else if (rd_type.compare("BeamRelevant") == 0)
{
m_multi_rd[i_rd] =
diff --git a/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp b/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp
index 9d2b3eac2..72395cc84 100644
--- a/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp
+++ b/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp
@@ -23,7 +23,7 @@ using namespace amrex;
ParticleEnergy::ParticleEnergy (std::string rd_name)
: ReducedDiags{rd_name}
{
- // get WarpX class object
+ // get a reference to WarpX instance
auto & warpx = WarpX::GetInstance();
// get MultiParticleContainer class object
@@ -33,7 +33,7 @@ ParticleEnergy::ParticleEnergy (std::string rd_name)
const auto nSpecies = mypc.nSpecies();
// resize data array
- m_data.resize(2*nSpecies+2,0.0);
+ m_data.resize(2*nSpecies+2, 0.0_rt);
// get species names (std::vector<std::string>)
const auto species_names = mypc.GetSpeciesNames();
@@ -43,8 +43,7 @@ ParticleEnergy::ParticleEnergy (std::string rd_name)
if ( m_IsNotRestart )
{
// open file
- std::ofstream ofs{m_path + m_rd_name + "." + m_extension,
- std::ofstream::out | std::ofstream::app};
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
// write header row
ofs << "#";
ofs << "[1]step()";
diff --git a/Source/Diagnostics/ReducedDiags/ParticleHistogram.cpp b/Source/Diagnostics/ReducedDiags/ParticleHistogram.cpp
index f0377ca32..442f248f5 100644
--- a/Source/Diagnostics/ReducedDiags/ParticleHistogram.cpp
+++ b/Source/Diagnostics/ReducedDiags/ParticleHistogram.cpp
@@ -91,8 +91,7 @@ ParticleHistogram::ParticleHistogram (std::string rd_name)
if ( m_IsNotRestart )
{
// open file
- std::ofstream ofs{m_path + m_rd_name + "." + m_extension,
- std::ofstream::out | std::ofstream::app};
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
// write header row
ofs << "#";
ofs << "[1]step()";
@@ -122,7 +121,7 @@ void ParticleHistogram::ComputeDiags (int step)
// Judge if the diags should be done
if (!m_intervals.contains(step+1)) return;
- // get WarpX class object
+ // get a reference to WarpX instance
auto & warpx = WarpX::GetInstance();
// get time at level 0
diff --git a/Source/Diagnostics/ReducedDiags/ParticleNumber.cpp b/Source/Diagnostics/ReducedDiags/ParticleNumber.cpp
index 701698039..111ab07c0 100644
--- a/Source/Diagnostics/ReducedDiags/ParticleNumber.cpp
+++ b/Source/Diagnostics/ReducedDiags/ParticleNumber.cpp
@@ -8,11 +8,13 @@
#include "ParticleNumber.H"
#include "WarpX.H"
+using namespace amrex::literals;
+
// constructor
ParticleNumber::ParticleNumber (std::string rd_name)
: ReducedDiags{rd_name}
{
- // get WarpX class object
+ // get a reference to WarpX instance
auto & warpx = WarpX::GetInstance();
// get MultiParticleContainer class object
@@ -23,7 +25,7 @@ ParticleNumber::ParticleNumber (std::string rd_name)
// resize data array to 2*(nSpecies+1) (each species + sum over all species
// for both number of macroparticles and of physical particles)
- m_data.resize(2*(nSpecies+1), amrex::Real(0.));
+ m_data.resize(2*(nSpecies+1), 0.0_rt);
// get species names (std::vector<std::string>)
const auto species_names = mypc.GetSpeciesNames();
@@ -33,8 +35,7 @@ ParticleNumber::ParticleNumber (std::string rd_name)
if ( m_IsNotRestart )
{
// open file
- std::ofstream ofs{m_path + m_rd_name + "." + m_extension,
- std::ofstream::out | std::ofstream::app};
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
// write header row
ofs << "#";
ofs << "[1]step()";
@@ -95,8 +96,8 @@ void ParticleNumber::ComputeDiags (int step)
const int idx_first_species_sum_weight = idx_total_sum_weight + 1;
// Initialize total number of macroparticles and total weight (all species) to 0
- m_data[idx_total_macroparticles] = amrex::Real(0.);
- m_data[idx_total_sum_weight] = amrex::Real(0.);
+ m_data[idx_total_macroparticles] = 0.0_rt;
+ m_data[idx_total_sum_weight] = 0.0_rt;
// loop over species
for (int i_s = 0; i_s < nSpecies; ++i_s)
diff --git a/Source/Diagnostics/ReducedDiags/RhoMaximum.H b/Source/Diagnostics/ReducedDiags/RhoMaximum.H
new file mode 100644
index 000000000..0638fae16
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/RhoMaximum.H
@@ -0,0 +1,40 @@
+/* Copyright 2020 Neil Zaim
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#ifndef WARPX_DIAGNOSTICS_REDUCEDDIAGS_RHOMAXIMUM_H_
+#define WARPX_DIAGNOSTICS_REDUCEDDIAGS_RHOMAXIMUM_H_
+
+#include "ReducedDiags.H"
+#include "Diagnostics/ComputeDiagFunctors/ComputeDiagFunctor.H"
+
+/**
+ * This class mainly contains a function that computes the extrema of the total charge density
+ * and of the charge density of each charged species.
+ */
+class RhoMaximum : public ReducedDiags
+{
+public:
+
+ /** constructor
+ * @param[in] rd_name reduced diags names
+ */
+ RhoMaximum(std::string rd_name);
+
+ /** This function computes the maximum and minimum values of rho (summed over all species) and
+ * the maximum absolute value of rho for each species.
+ */
+ virtual void ComputeDiags(int step) override final;
+
+private:
+ /** Vector of (pointers to) functors to compute rho, per level, per species. We reuse here the
+ * same functors as those used for regular diagnostics.
+ * */
+ amrex::Vector< amrex::Vector <std::unique_ptr<ComputeDiagFunctor > > > m_rho_functors;
+
+};
+
+#endif // WARPX_DIAGNOSTICS_REDUCEDDIAGS_RHOMAXIMUM_H_
diff --git a/Source/Diagnostics/ReducedDiags/RhoMaximum.cpp b/Source/Diagnostics/ReducedDiags/RhoMaximum.cpp
new file mode 100644
index 000000000..e804d4f3d
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/RhoMaximum.cpp
@@ -0,0 +1,170 @@
+/* Copyright 2020 Neil Zaim
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#include "RhoMaximum.H"
+#include "Diagnostics/ComputeDiagFunctors/RhoFunctor.H"
+#include "WarpX.H"
+
+using namespace amrex::literals;
+
+// constructor
+RhoMaximum::RhoMaximum (std::string rd_name)
+: ReducedDiags{rd_name}
+{
+
+ // RZ coordinate is not working
+#if (defined WARPX_DIM_RZ)
+ AMREX_ALWAYS_ASSERT_WITH_MESSAGE(false,
+ "RhoMaximum reduced diagnostics does not work for RZ coordinate.");
+#endif
+
+ // read number of levels
+ int nLevel = 0;
+ amrex::ParmParse pp("amr");
+ pp.query("max_level", nLevel);
+ nLevel += 1;
+ m_rho_functors.resize(nLevel);
+
+ // We do not use coarsening in the functors for this diag.
+ const amrex::IntVect crse_ratio = amrex::IntVect(1);
+
+ // Initialize functors for the total charge density
+ for (int lev = 0; lev < nLevel; ++lev)
+ {
+ m_rho_functors[lev].push_back(std::make_unique<RhoFunctor>(lev, crse_ratio));
+ }
+
+ // get MultiParticleContainer class object
+ const auto & mypc = WarpX::GetInstance().GetPartContainer();
+
+ // get number of species (int)
+ const auto nSpecies = mypc.nSpecies();
+
+ // A vector to store the indices of charges species in mypc
+ amrex::Vector<int> indices_charged_species;
+
+ // number of charged species
+ int n_charged_species = 0;
+
+ for (int i = 0; i < nSpecies; ++i)
+ {
+ // Only charged species are relevant for this diag
+ if (mypc.GetParticleContainer(i).getCharge() != 0.0_rt)
+ {
+ indices_charged_species.push_back(i);
+ n_charged_species += 1;
+ for (int lev = 0; lev < nLevel; ++lev)
+ {
+ // Initialize functors for the charge density of each charged species
+ m_rho_functors[lev].push_back(std::make_unique<RhoFunctor>(lev, crse_ratio, i));
+ }
+ }
+ }
+
+ // get species names (std::vector<std::string>)
+ const auto species_names = mypc.GetSpeciesNames();
+
+ // Min and max of total rho + max of |rho| for each species
+ const int noutputs_per_level = 2+n_charged_species;
+ m_data.resize(static_cast<std::size_t>(nLevel*noutputs_per_level), 0.0_rt);
+
+ if (amrex::ParallelDescriptor::IOProcessor())
+ {
+ if ( m_IsNotRestart )
+ {
+ // open file
+ std::ofstream ofs{m_path + m_rd_name + "." + m_extension, std::ofstream::out};
+ // write header row
+ ofs << "#";
+ ofs << "[1]step()";
+ ofs << m_sep;
+ ofs << "[2]time(s)";
+ constexpr int shift_max_rho = 3;
+ constexpr int shift_min_rho = 4;
+ constexpr int shift_first_species = 5;
+ for (int lev = 0; lev < nLevel; ++lev)
+ {
+ ofs << m_sep;
+ ofs << "[" + std::to_string(shift_max_rho+noutputs_per_level*lev) + "]";
+ ofs << "max_rho_lev"+std::to_string(lev)+" (C/m^3)";
+ ofs << m_sep;
+ ofs << "[" + std::to_string(shift_min_rho+noutputs_per_level*lev) + "]";
+ ofs << "min_rho_lev"+std::to_string(lev)+" (C/m^3)";
+ for (int i = 0; i < n_charged_species; ++i)
+ {
+ ofs << m_sep;
+ ofs << "[" + std::to_string(shift_first_species+i+noutputs_per_level*lev) + "]";
+ ofs << "max_" + species_names[indices_charged_species[i]]
+ + "_|rho|_lev"+std::to_string(lev)+" (C/m^3)";
+ }
+ }
+ ofs << std::endl;
+ // close file
+ ofs.close();
+ }
+ }
+
+}
+// end constructor
+
+// function that computes maximum charge density values
+void RhoMaximum::ComputeDiags (int step)
+{
+ // Judge if the diags should be done
+ if (!m_intervals.contains(step+1)) { return; }
+
+ // get a reference to WarpX instance
+ auto & warpx = WarpX::GetInstance();
+
+ // get number of levels
+ const auto nLevel = warpx.finestLevel() + 1;
+
+ const int n_charged_species = m_rho_functors[0].size() - 1;
+ // Min and max of total rho + max of |rho| for each species
+ const int noutputs_per_level = 2+n_charged_species;
+
+ // loop over refinement levels
+ for (int lev = 0; lev < nLevel; ++lev)
+ {
+ // Declare a temporary MultiFAB to store the charge densities.
+ amrex::BoxArray ba = warpx.boxArray(lev);
+ amrex::DistributionMapping dmap = warpx.DistributionMap(lev);
+ constexpr int ncomp = 1;
+ constexpr int ngrow = 0;
+ amrex::MultiFab mf_temp(ba, dmap, ncomp, ngrow);
+
+ // Fill temporary MultiFAB with total charge density
+ constexpr int idx_total_rho_functor = 0;
+ constexpr int idx_first_species_functor = 1;
+ constexpr int icomp = 0;
+ constexpr int i_buffer = 0;
+ m_rho_functors[lev][idx_total_rho_functor]->operator()(mf_temp, icomp, i_buffer);
+
+ constexpr int idx_max_rho_data = 0;
+ constexpr int idx_min_rho_data = 1;
+ constexpr int idx_first_species_data = 2;
+
+ // Fill output array with min and max of total rho
+ m_data[lev*noutputs_per_level + idx_max_rho_data] = mf_temp.max(icomp);
+ m_data[lev*noutputs_per_level + idx_min_rho_data] = mf_temp.min(icomp);
+
+ // Loop over all charged species
+ for (int i = 0; i < n_charged_species; ++i)
+ {
+ // Fill temporary MultiFAB with the species charge density
+ m_rho_functors[lev][idx_first_species_functor+i]->operator()(mf_temp, icomp, i_buffer);
+ // Fill output array with max |rho| of species
+ m_data[lev*noutputs_per_level + idx_first_species_data + i] = mf_temp.norm0();
+ }
+ }
+ // end loop over refinement levels
+
+ /* m_data now contains up-to-date values for:
+ * [max(rho), min(rho), max(|rho_charged_species1|), max(|rho_charged_species2|), ...] */
+
+}
+// end void RhoMaximum::ComputeDiags