aboutsummaryrefslogtreecommitdiff
path: root/Source
diff options
context:
space:
mode:
Diffstat (limited to 'Source')
-rw-r--r--Source/Diagnostics/Make.package2
-rw-r--r--Source/Diagnostics/ReducedDiags/FieldEnergy.H36
-rw-r--r--Source/Diagnostics/ReducedDiags/FieldEnergy.cpp139
-rw-r--r--Source/Diagnostics/ReducedDiags/Make.package14
-rw-r--r--Source/Diagnostics/ReducedDiags/MultiReducedDiags.H47
-rw-r--r--Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp94
-rw-r--r--Source/Diagnostics/ReducedDiags/ParticleEnergy.H37
-rw-r--r--Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp166
-rw-r--r--Source/Diagnostics/ReducedDiags/ReducedDiags.H59
-rw-r--r--Source/Diagnostics/ReducedDiags/ReducedDiags.cpp92
-rw-r--r--Source/Evolve/WarpXEvolveEM.cpp8
-rw-r--r--Source/Initialization/WarpXInitData.cpp7
-rw-r--r--Source/WarpX.H4
-rw-r--r--Source/WarpX.cpp5
14 files changed, 709 insertions, 1 deletions
diff --git a/Source/Diagnostics/Make.package b/Source/Diagnostics/Make.package
index 710e4399b..12560b49e 100644
--- a/Source/Diagnostics/Make.package
+++ b/Source/Diagnostics/Make.package
@@ -15,5 +15,7 @@ ifeq ($(USE_OPENPMD), TRUE)
CEXE_sources += WarpXOpenPMD.cpp
endif
+include $(WARPX_HOME)/Source/Diagnostics/ReducedDiags/Make.package
+
INCLUDE_LOCATIONS += $(WARPX_HOME)/Source/Diagnostics
VPATH_LOCATIONS += $(WARPX_HOME)/Source/Diagnostics
diff --git a/Source/Diagnostics/ReducedDiags/FieldEnergy.H b/Source/Diagnostics/ReducedDiags/FieldEnergy.H
new file mode 100644
index 000000000..82fa4b6c4
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/FieldEnergy.H
@@ -0,0 +1,36 @@
+/* Copyright 2019-2020 Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#ifndef WARPX_DIAGNOSTICS_REDUCEDDIAGS_FIELDENERGY_H_
+#define WARPX_DIAGNOSTICS_REDUCEDDIAGS_FIELDENERGY_H_
+
+#include "ReducedDiags.H"
+#include <fstream>
+
+/**
+ * This class mainly contains a function that
+ * computes the field energy.
+ */
+class FieldEnergy : public ReducedDiags
+{
+public:
+
+ /** constructor
+ * @param[in] rd_name reduced diags names */
+ FieldEnergy(std::string rd_name);
+
+ /** This funciton computes the field energy (EF).
+ * EF = E eps / 2 + B / mu / 2,
+ * where E is the electric field,
+ * B is the magnetic field,
+ * eps is the vacuum permittivity,
+ * mu is the vacuum permeability. */
+ virtual void ComputeDiags(int step) override final;
+
+};
+
+#endif
diff --git a/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp b/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp
new file mode 100644
index 000000000..73e6a1c9a
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/FieldEnergy.cpp
@@ -0,0 +1,139 @@
+/* Copyright 2019-2020 Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#include "FieldEnergy.H"
+#include "WarpX.H"
+#include "WarpXConst.H"
+#include "AMReX_REAL.H"
+#include "AMReX_ParticleReduce.H"
+#include <iostream>
+#include <cmath>
+
+using namespace amrex;
+
+// constructor
+FieldEnergy::FieldEnergy (std::string rd_name)
+: ReducedDiags{rd_name}
+{
+
+ // RZ coordinate is not working
+ #if (defined WARPX_DIM_RZ)
+ AMREX_ALWAYS_ASSERT_WITH_MESSAGE(false,
+ "FieldEnergy reduced diagnostics does not work for RZ coordinate.");
+ #endif
+
+ // get WarpX class object
+ auto & warpx = WarpX::GetInstance();
+
+ // read number of levels
+ int nLevel = 0;
+ ParmParse pp("amr");
+ pp.query("max_level", nLevel);
+ nLevel += 1;
+
+ // resize data array
+ m_data.resize(3*nLevel,0.0);
+
+ if (ParallelDescriptor::IOProcessor())
+ {
+ if ( m_IsNotRestart )
+ {
+ // open file
+ std::ofstream ofs;
+ ofs.open(m_path + m_rd_name + "." + m_extension,
+ std::ofstream::out | std::ofstream::app);
+ // write header row
+ ofs << "#";
+ ofs << "[1]step";
+ ofs << m_sep;
+ ofs << "[2]time(s)";
+ for (int lev = 0; lev < nLevel; ++lev)
+ {
+ ofs << m_sep;
+ ofs << "[" + std::to_string(3+3*lev) + "]";
+ ofs << "total(J)lev"+std::to_string(lev);
+ ofs << m_sep;
+ ofs << "[" + std::to_string(4+3*lev) + "]";
+ ofs << "E(J)lev"+std::to_string(lev);
+ ofs << m_sep;
+ ofs << "[" + std::to_string(5+3*lev) + "]";
+ ofs << "B(J)lev"+std::to_string(lev);
+ }
+ ofs << std::endl;
+ // close file
+ ofs.close();
+ }
+ }
+
+}
+// end constructor
+
+// function that computes field energy
+void FieldEnergy::ComputeDiags (int step)
+{
+
+ // Judge if the diags should be done
+ if ( (step+1) % m_freq != 0 ) { return; }
+
+ // get WarpX class object
+ auto & warpx = WarpX::GetInstance();
+
+ // get number of level
+ auto nLevel = warpx.finestLevel() + 1;
+
+ // loop over refinement levels
+ for (int lev = 0; lev < nLevel; ++lev)
+ {
+
+ // get MultiFab data at lev
+ const MultiFab & Ex = warpx.getEfield(lev,0);
+ const MultiFab & Ey = warpx.getEfield(lev,1);
+ const MultiFab & Ez = warpx.getEfield(lev,2);
+ const MultiFab & Bx = warpx.getBfield(lev,0);
+ const MultiFab & By = warpx.getBfield(lev,1);
+ const MultiFab & Bz = warpx.getBfield(lev,2);
+
+ // get cell size
+ Geometry const & geom = warpx.Geom(lev);
+ auto domain_box = geom.Domain();
+ #if (AMREX_SPACEDIM == 2)
+ auto dV = geom.CellSize(0) * geom.CellSize(1);
+ #elif (AMREX_SPACEDIM == 3)
+ auto dV = geom.CellSize(0) * geom.CellSize(1) * geom.CellSize(2);
+ #endif
+
+ // compute E squared
+ Real tmpx = Ex.norm2(0,geom.periodicity());
+ Real tmpy = Ey.norm2(0,geom.periodicity());
+ Real tmpz = Ez.norm2(0,geom.periodicity());
+ Real Es = tmpx*tmpx + tmpy*tmpy + tmpz*tmpz;
+
+ // compute B squared
+ tmpx = Bx.norm2(0,geom.periodicity());
+ tmpy = By.norm2(0,geom.periodicity());
+ tmpz = Bz.norm2(0,geom.periodicity());
+ Real Bs = tmpx*tmpx + tmpy*tmpy + tmpz*tmpz;
+
+ // save data
+ m_data[lev*3+1] = 0.5 * Es * PhysConst::ep0 * dV;
+ m_data[lev*3+2] = 0.5 * Bs / PhysConst::mu0 * dV;
+ m_data[lev*3+0] = m_data[lev*3+1] + m_data[lev*3+2];
+
+ }
+ // end loop over refinement levels
+
+ /* m_data now contains up-to-date values for:
+ * [total field energy at level 0,
+ * electric field energy at level 0,
+ * magnetic field energy at level 0,
+ * total field energy at level 1,
+ * electric field energy at level 1,
+ * magnetic field energy at level 1,
+ * ......] */
+
+}
+// end void FieldEnergy::ComputeDiags
diff --git a/Source/Diagnostics/ReducedDiags/Make.package b/Source/Diagnostics/ReducedDiags/Make.package
new file mode 100644
index 000000000..37f76d3d5
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/Make.package
@@ -0,0 +1,14 @@
+CEXE_headers += MultiReducedDiags.H
+CEXE_sources += MultiReducedDiags.cpp
+
+CEXE_headers += ReducedDiags.H
+CEXE_sources += ReducedDiags.cpp
+
+CEXE_headers += ParticleEnergy.H
+CEXE_sources += ParticleEnergy.cpp
+
+CEXE_headers += FieldEnergy.H
+CEXE_sources += FieldEnergy.cpp
+
+INCLUDE_LOCATIONS += $(WARPX_HOME)/Source/Diagnostics/ReducedDiags
+VPATH_LOCATIONS += $(WARPX_HOME)/Source/Diagnostics/ReducedDiags
diff --git a/Source/Diagnostics/ReducedDiags/MultiReducedDiags.H b/Source/Diagnostics/ReducedDiags/MultiReducedDiags.H
new file mode 100644
index 000000000..79f487d81
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/MultiReducedDiags.H
@@ -0,0 +1,47 @@
+/* Copyright 2019-2020 Maxence Thevenet, Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#ifndef WARPX_DIAGNOSTICS_REDUCEDDIAGS_MULTIREDUCEDDIAGS_H_
+#define WARPX_DIAGNOSTICS_REDUCEDDIAGS_MULTIREDUCEDDIAGS_H_
+
+#include "ReducedDiags.H"
+#include <vector>
+#include <string>
+#include <memory>
+
+/**
+ * This class holds multiple instances of ReducedDiagnostics, and contains
+ * general functions to initialize, compute, and write these diagnostics
+ * to file.
+ */
+class MultiReducedDiags
+{
+public:
+
+ /// Bool: whether or not reduced diagnostics are activated
+ int m_plot_rd = 0;
+
+ /// names of reduced diagnostics
+ std::vector<std::string> m_rd_names;
+
+ /// m_multi_rd stores a pointer to each reduced diagnostics
+ std::vector<std::unique_ptr<ReducedDiags>> m_multi_rd;
+
+ /// constructor
+ MultiReducedDiags();
+
+ /** Loop over all ReducedDiags and call their ComputeDiags
+ * @param[in] step current iteration time */
+ void ComputeDiags(int step);
+
+ /** Loop over all ReducedDiags and call their WriteToFile
+ * @param[in] step current iteration time */
+ void WriteToFile(int step);
+
+};
+
+#endif
diff --git a/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp b/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp
new file mode 100644
index 000000000..71a33924b
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/MultiReducedDiags.cpp
@@ -0,0 +1,94 @@
+/* Copyright 2019-2020 Maxence Thevenet, Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#include "ParticleEnergy.H"
+#include "FieldEnergy.H"
+#include "MultiReducedDiags.H"
+#include "AMReX_ParmParse.H"
+#include "AMReX_ParallelDescriptor.H"
+#include <fstream>
+
+using namespace amrex;
+
+// constructor
+MultiReducedDiags::MultiReducedDiags ()
+{
+
+ // read reduced diags names
+ ParmParse pp("warpx");
+ m_plot_rd = pp.queryarr("reduced_diags_names", m_rd_names);
+
+ // if names are not given, reduced diags will not be done
+ if ( m_plot_rd == 0 ) { return; }
+
+ // resize
+ m_multi_rd.resize(m_rd_names.size());
+
+ // loop over all reduced diags
+ for (int i_rd = 0; i_rd < m_rd_names.size(); ++i_rd)
+ {
+
+ ParmParse pp(m_rd_names[i_rd]);
+
+ // read reduced diags type
+ std::string rd_type;
+ pp.query("type", rd_type);
+
+ // match diags
+ if (rd_type.compare("ParticleEnergy") == 0)
+ {
+ m_multi_rd[i_rd].reset
+ ( new ParticleEnergy(m_rd_names[i_rd]));
+ }
+ else if (rd_type.compare("FieldEnergy") == 0)
+ {
+ m_multi_rd[i_rd].reset
+ ( new FieldEnergy(m_rd_names[i_rd]));
+ }
+ else
+ { Abort("No matching reduced diagnostics type found."); }
+ // end if match diags
+
+ }
+ // end loop over all reduced diags
+
+}
+// end constructor
+
+// call functions to compute diags
+void MultiReducedDiags::ComputeDiags (int step)
+{
+ // loop over all reduced diags
+ for (int i_rd = 0; i_rd < m_rd_names.size(); ++i_rd)
+ {
+ m_multi_rd[i_rd] -> ComputeDiags(step);
+ }
+ // end loop over all reduced diags
+}
+// end void MultiReducedDiags::ComputeDiags
+
+// funciton to write data
+void MultiReducedDiags::WriteToFile (int step)
+{
+
+ // Only the I/O rank does
+ if ( !ParallelDescriptor::IOProcessor() ) { return; }
+
+ // loop over all reduced diags
+ for (int i_rd = 0; i_rd < m_rd_names.size(); ++i_rd)
+ {
+
+ // Judge if the diags should be done
+ if ( (step+1) % m_multi_rd[i_rd]->m_freq != 0 ) { return; }
+
+ // call the write to file function
+ m_multi_rd[i_rd]->WriteToFile(step);
+
+ }
+ // end loop over all reduced diags
+}
+// end void MultiReducedDiags::WriteToFile
diff --git a/Source/Diagnostics/ReducedDiags/ParticleEnergy.H b/Source/Diagnostics/ReducedDiags/ParticleEnergy.H
new file mode 100644
index 000000000..d7c60a24e
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/ParticleEnergy.H
@@ -0,0 +1,37 @@
+/* Copyright 2019-2020 Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#ifndef WARPX_DIAGNOSTICS_REDUCEDDIAGS_PARTICLEENERGY_H_
+#define WARPX_DIAGNOSTICS_REDUCEDDIAGS_PARTICLEENERGY_H_
+
+#include "ReducedDiags.H"
+#include <fstream>
+
+/**
+ * This class mainly contains a function that
+ * computes the particle relativistic kinetic energy
+ * of each species.
+ */
+class ParticleEnergy : public ReducedDiags
+{
+public:
+
+ /** constructor
+ * @param[in] rd_name reduced diags names */
+ ParticleEnergy(std::string rd_name);
+
+ /** This funciton computes the particle relativistic
+ * kinetic energy (EP).
+ * \param [in] step current time step
+ * EP = sqrt( p^2 c^2 + m^2 c^4 ) - m c^2,
+ * where p is the relativistic momentum,
+ * m is the particle rest mass. */
+ virtual void ComputeDiags(int step) override final;
+
+};
+
+#endif
diff --git a/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp b/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp
new file mode 100644
index 000000000..132ad2165
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/ParticleEnergy.cpp
@@ -0,0 +1,166 @@
+/* Copyright 2019-2020 Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#include "ParticleEnergy.H"
+#include "WarpX.H"
+#include "WarpXConst.H"
+#include "AMReX_REAL.H"
+#include "AMReX_ParticleReduce.H"
+#include <iostream>
+#include <cmath>
+#include <limits>
+
+using namespace amrex;
+
+// constructor
+ParticleEnergy::ParticleEnergy (std::string rd_name)
+: ReducedDiags{rd_name}
+{
+ // get WarpX class object
+ auto & warpx = WarpX::GetInstance();
+
+ // get MultiParticleContainer class object
+ auto & mypc = warpx.GetPartContainer();
+
+ // get number of species (int)
+ auto nSpecies = mypc.nSpecies();
+
+ // resize data array
+ m_data.resize(2*nSpecies+2,0.0);
+
+ // get species names (std::vector<std::string>)
+ auto species_names = mypc.GetSpeciesNames();
+
+ if (ParallelDescriptor::IOProcessor())
+ {
+ if ( m_IsNotRestart )
+ {
+ // open file
+ std::ofstream ofs;
+ ofs.open(m_path + m_rd_name + "." + m_extension,
+ std::ofstream::out | std::ofstream::app);
+ // write header row
+ ofs << "#";
+ ofs << "[1]step";
+ ofs << m_sep;
+ ofs << "[2]time(s)";
+ ofs << m_sep;
+ ofs << "[3]total(J)";
+ for (int i = 0; i < nSpecies; ++i)
+ {
+ ofs << m_sep;
+ ofs << "[" + std::to_string(4+i) + "]";
+ ofs << species_names[i]+"(J)";
+ }
+ ofs << m_sep;
+ ofs << "[" + std::to_string(4+nSpecies) + "]";
+ ofs << "total.mean(J)";
+ for (int i = 0; i < nSpecies; ++i)
+ {
+ ofs << m_sep;
+ ofs << "[" + std::to_string(5+nSpecies+i) + "]";
+ ofs << species_names[i]+".mean(J)";
+ }
+ ofs << std::endl;
+ // close file
+ ofs.close();
+ }
+ }
+
+}
+// end constructor
+
+// function that computes kinetic energy
+void ParticleEnergy::ComputeDiags (int step)
+{
+
+ // Judge if the diags should be done
+ if ( (step+1) % m_freq != 0 ) { return; }
+
+ // get MultiParticleContainer class object
+ auto & mypc = WarpX::GetInstance().GetPartContainer();
+
+ // get number of species (int)
+ auto nSpecies = mypc.nSpecies();
+
+ // get species names (std::vector<std::string>)
+ auto species_names = mypc.GetSpeciesNames();
+
+ // speed of light squared
+ auto c2 = PhysConst::c * PhysConst::c;
+
+ // loop over species
+ for (int i_s = 0; i_s < nSpecies; ++i_s)
+ {
+ // get WarpXParticleContainer class object
+ auto & myspc = mypc.GetParticleContainer(i_s);
+
+ // get mass (Real)
+ auto m = myspc.getMass();
+
+ using PType = typename WarpXParticleContainer::SuperParticleType;
+
+ // Use amex::ReduceSum to compute the sum of energies of all particles
+ // held by the current MPI rank, for this species. This involves a loop over all
+ // boxes held by this MPI rank.
+ auto Etot = ReduceSum( myspc,
+ [=] AMREX_GPU_HOST_DEVICE (const PType& p) -> Real
+ {
+ auto w = p.rdata(PIdx::w);
+ auto ux = p.rdata(PIdx::ux);
+ auto uy = p.rdata(PIdx::uy);
+ auto uz = p.rdata(PIdx::uz);
+ auto us = (ux*ux + uy*uy + uz*uz);
+ return ( std::sqrt(us*c2 + c2*c2) - c2 ) * m * w;
+ });
+
+ // Same thing for the particles weights.
+ auto Wtot = ReduceSum( myspc,
+ [=] AMREX_GPU_HOST_DEVICE (const PType& p) -> Real
+ {
+ return p.rdata(PIdx::w);
+ });
+
+ // reduced sum over mpi ranks
+ ParallelDescriptor::ReduceRealSum
+ (Etot, ParallelDescriptor::IOProcessorNumber());
+ ParallelDescriptor::ReduceRealSum
+ (Wtot, ParallelDescriptor::IOProcessorNumber());
+
+ // save results for this species i_s into m_data
+ m_data[i_s+1] = Etot;
+ if ( Wtot > std::numeric_limits<Real>::min() )
+ { m_data[nSpecies+2+i_s] = Etot / Wtot; }
+ else
+ { m_data[nSpecies+2+i_s] = 0.0; }
+
+ }
+ // end loop over species
+
+ // save total energy
+ // loop over species
+ m_data[0] = 0.0; // total energy
+ m_data[nSpecies+1] = 0.0; // total mean energy
+ for (int i_s = 0; i_s < nSpecies; ++i_s)
+ {
+ m_data[0] += m_data[i_s+1];
+ m_data[nSpecies+1] += m_data[nSpecies+2+i_s];
+ }
+ // end loop over species
+
+ /* m_data now contains up-to-date values for:
+ * [total energy (all species),
+ * total energy (species 1),
+ * ...,
+ * total energy (species n),
+ * mean energy (all species),
+ * mean energy (species 1),
+ * ...,
+ * mean energy (species n)] */
+
+}
+// end void ParticleEnergy::ComputeDiags
diff --git a/Source/Diagnostics/ReducedDiags/ReducedDiags.H b/Source/Diagnostics/ReducedDiags/ReducedDiags.H
new file mode 100644
index 000000000..7ff065f49
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/ReducedDiags.H
@@ -0,0 +1,59 @@
+/* Copyright 2019-2020 Maxence Thevenet, Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#ifndef WARPX_DIAGNOSTICS_REDUCEDDIAGS_REDUCEDDIAGS_H_
+#define WARPX_DIAGNOSTICS_REDUCEDDIAGS_REDUCEDDIAGS_H_
+
+#include "AMReX_REAL.H"
+#include <string>
+#include <vector>
+#include <fstream>
+
+/**
+ * Base class for reduced diagnostics. Each type of reduced diagnostics is
+ * implemented in a derived class, and must override the (pure virtual)
+ * function ComputeDiags.
+ */
+class ReducedDiags
+{
+public:
+
+ /// output path (default)
+ std::string m_path = "./diags/reducedfiles/";
+
+ /// output extension (default)
+ std::string m_extension = "txt";
+
+ /// diags name
+ std::string m_rd_name;
+
+ /// output frequency
+ int m_freq = 1;
+
+ /// check if it is a restart run
+ int m_IsNotRestart = 1;
+
+ /// separator in the output file
+ std::string m_sep = ",";
+
+ /// output data
+ std::vector<amrex::Real> m_data;
+
+ /** constructor
+ * @param[in] rd_name reduced diags name */
+ ReducedDiags(std::string rd_name);
+
+ /// function to compute diags
+ virtual void ComputeDiags(int step) = 0;
+
+ /** write to file function
+ * @param[in] step time step */
+ virtual void WriteToFile(int step) const;
+
+};
+
+#endif
diff --git a/Source/Diagnostics/ReducedDiags/ReducedDiags.cpp b/Source/Diagnostics/ReducedDiags/ReducedDiags.cpp
new file mode 100644
index 000000000..81831aa79
--- /dev/null
+++ b/Source/Diagnostics/ReducedDiags/ReducedDiags.cpp
@@ -0,0 +1,92 @@
+/* Copyright 2019-2020 Maxence Thevenet, Yinjian Zhao
+ *
+ * This file is part of WarpX.
+ *
+ * License: BSD-3-Clause-LBNL
+ */
+
+#include "ReducedDiags.H"
+#include "WarpX.H"
+#include "AMReX_ParmParse.H"
+#include "AMReX_Utility.H"
+#include <iomanip>
+
+using namespace amrex;
+
+// constructor
+ReducedDiags::ReducedDiags (std::string rd_name)
+{
+
+ m_rd_name = rd_name;
+
+ ParmParse pp(m_rd_name);
+
+ // read path
+ pp.query("path", m_path);
+
+ // read extension
+ pp.query("extension", m_extension);
+
+ // creater folder
+ if (!UtilCreateDirectory(m_path, 0755))
+ { CreateDirectoryFailed(m_path); }
+
+ // check if it is a restart run
+ std::string restart_chkfile = "";
+ ParmParse pp_amr("amr");
+ pp_amr.query("restart", restart_chkfile);
+ m_IsNotRestart = restart_chkfile.empty();
+
+ // replace / create output file
+ if ( m_IsNotRestart ) // not a restart
+ {
+ std::ofstream ofs;
+ ofs.open(m_path+m_rd_name+"."+m_extension, std::ios::trunc);
+ ofs.close();
+ }
+
+ // read reduced diags frequency
+ pp.query("frequency", m_freq);
+
+ // read separator
+ pp.query("separator", m_sep);
+
+}
+// end constructor
+
+// write to file function
+void ReducedDiags::WriteToFile (int step) const
+{
+
+ // open file
+ std::ofstream ofs;
+ ofs.open(m_path + m_rd_name + "." + m_extension,
+ std::ofstream::out | std::ofstream::app);
+
+ // write step
+ ofs << step+1;
+
+ ofs << m_sep;
+
+ // set precision
+ ofs << std::fixed << std::setprecision(14) << std::scientific;
+
+ // write time
+ ofs << WarpX::GetInstance().gett_new(0);
+
+ // loop over data size and write
+ for (int i = 0; i < m_data.size(); ++i)
+ {
+ ofs << m_sep;
+ ofs << m_data[i];
+ }
+ // end loop over data size
+
+ // end line
+ ofs << std::endl;
+
+ // close file
+ ofs.close();
+
+}
+// end ReducedDiags::WriteToFile
diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp
index 84b1491c1..bb1300562 100644
--- a/Source/Evolve/WarpXEvolveEM.cpp
+++ b/Source/Evolve/WarpXEvolveEM.cpp
@@ -24,7 +24,6 @@
#include <AMReX_AmrMeshInSituBridge.H>
#endif
-
using namespace amrex;
void
@@ -206,6 +205,13 @@ WarpX::EvolveEM (int numsteps)
t_new[i] = cur_time;
}
+ /// reduced diags
+ if (reduced_diags->m_plot_rd != 0)
+ {
+ reduced_diags->ComputeDiags(step);
+ reduced_diags->WriteToFile(step);
+ }
+
// slice gen //
if (to_make_plot || to_write_openPMD || do_insitu || to_make_slice_plot)
{
diff --git a/Source/Initialization/WarpXInitData.cpp b/Source/Initialization/WarpXInitData.cpp
index b046c667e..8b2fe1831 100644
--- a/Source/Initialization/WarpXInitData.cpp
+++ b/Source/Initialization/WarpXInitData.cpp
@@ -91,6 +91,13 @@ WarpX::InitData ()
if ((insitu_int > 0) && (insitu_start == 0))
UpdateInSitu();
+
+ // Write reduced diagnostics before the first iteration.
+ if (reduced_diags->m_plot_rd != 0)
+ {
+ reduced_diags->ComputeDiags(-1);
+ reduced_diags->WriteToFile(-1);
+ }
}
}
diff --git a/Source/WarpX.H b/Source/WarpX.H
index c368e2be9..ae25a8168 100644
--- a/Source/WarpX.H
+++ b/Source/WarpX.H
@@ -18,6 +18,7 @@
#include <BackTransformedDiagnostic.H>
#include <BilinearFilter.H>
#include <NCIGodfreyFilter.H>
+#include "MultiReducedDiags.H"
#ifdef WARPX_USE_PSATD
# include <SpectralSolver.H>
@@ -212,6 +213,9 @@ public:
amrex::Vector<amrex::Real> mirror_z_width;
amrex::Vector<int> mirror_z_npoints;
+ /// object with all reduced diagnotics, similar to MultiParticleContainer for species.
+ MultiReducedDiags* reduced_diags;
+
void applyMirrors(amrex::Real time);
void ComputeDt ();
diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp
index 71dd422ee..d93fab7df 100644
--- a/Source/WarpX.cpp
+++ b/Source/WarpX.cpp
@@ -201,6 +201,9 @@ WarpX::WarpX ()
}
do_back_transformed_particles = mypc->doBackTransformedDiagnostics();
+ /** create object for reduced diagnostics */
+ reduced_diags = new MultiReducedDiags();
+
Efield_aux.resize(nlevs_max);
Bfield_aux.resize(nlevs_max);
@@ -298,6 +301,8 @@ WarpX::~WarpX ()
ClearLevel(lev);
}
+ delete reduced_diags;
+
#ifdef BL_USE_SENSEI_INSITU
delete insitu_bridge;
#endif