diff --git a/src/Hipace.H b/src/Hipace.H index 7b8a73db7a..7595ec8432 100644 --- a/src/Hipace.H +++ b/src/Hipace.H @@ -313,7 +313,7 @@ private: Diagnostic m_diags; void ResizeFDiagFAB (const amrex::Box box, const int lev) { m_diags.ResizeFDiagFAB(box, lev); }; - void FillDiagnostics (int lev, int i_slice); + void FillDiagnostics (int i_slice); /** \brief get diagnostics Component names of Fields to output */ amrex::Vector& getDiagComps () { return m_diags.getComps(); }; /** \brief get diagnostics Component names of Fields to output */ diff --git a/src/Hipace.cpp b/src/Hipace.cpp index da4bc1f314..8445810fcd 100644 --- a/src/Hipace.cpp +++ b/src/Hipace.cpp @@ -246,7 +246,6 @@ void Hipace::MakeNewLevelFromScratch ( int lev, amrex::Real /*time*/, const amrex::BoxArray& ba, const amrex::DistributionMapping&) { - AMREX_ALWAYS_ASSERT(lev == 0); // We are going to ignore the DistributionMapping argument and build our own. amrex::DistributionMapping dm; @@ -360,7 +359,7 @@ Hipace::Evolve () for (int step = m_numprocs_z - 1 - m_rank_z; step <= m_max_step; step += m_numprocs_z) { #ifdef HIPACE_USE_OPENPMD - m_openpmd_writer.InitDiagnostics(step, m_output_period, m_max_step); + m_openpmd_writer.InitDiagnostics(step, m_output_period, m_max_step, finestLevel()+1); #endif if (m_verbose>=1) std::cout<<"Rank "< beamnames = getDiagBeamNames(); #ifdef HIPACE_USE_OPENPMD - constexpr int lev = 0; m_openpmd_writer.WriteDiagnostics(getDiagF(), m_multi_beam, getDiagGeom(), - m_physical_time, output_step, lev, getDiagSliceDir(), varnames, beamnames, - it, m_box_sorters, geom[lev], call_type); + m_physical_time, output_step, finestLevel()+1, getDiagSliceDir(), varnames, beamnames, + it, m_box_sorters, geom, call_type); #else amrex::Print()<<"WARNING: HiPACE++ compiled without openPMD support, the simulation has no I/O.\n"; #endif diff --git a/src/diagnostics/OpenPMDWriter.H b/src/diagnostics/OpenPMDWriter.H index f033240031..eb371ed00d 100644 --- a/src/diagnostics/OpenPMDWriter.H +++ b/src/diagnostics/OpenPMDWriter.H @@ -67,12 +67,13 @@ private: * \param[in] a_box_sorter_vec Vector (over species) of particles sorted by box * \param[in] geom Geometry of the simulation, to get the cell size etc. * \param[in] beamnames list of the names of the beam to be written to file + * \param[in] lev MR level */ void WriteBeamParticleData (MultiBeam& beams, openPMD::Iteration iteration, const int output_step, const int it, const amrex::Vector& a_box_sorter_vec, const amrex::Geometry& geom, - const amrex::Vector< std::string > beamnames); + const amrex::Vector< std::string > beamnames, const int lev); /** \brief writing openPMD field data * @@ -82,24 +83,25 @@ private: * \param[in] varnames list of variable names for the fields (ExmBy, EypBx, Ey, ...) * \param[in,out] iteration openPMD iteration to which the data is written * \param[in] output_step current time step to dump + * \param[in] lev MR level */ void WriteFieldData (amrex::FArrayBox const& fab, amrex::Geometry const& geom, const int slice_dir, const amrex::Vector< std::string > varnames, - openPMD::Iteration iteration, const int output_step); + openPMD::Iteration iteration, const int output_step, const int lev); /** Named Beam SoA attributes per particle as defined in BeamIdx */ amrex::Vector m_real_names {"weighting", "momentum_x", "momentum_y", "momentum_z"}; - /** Parallel openPMD-api Series object for output */ - std::unique_ptr< openPMD::Series > m_outputSeries; + /** vector over levels of openPMD-api Series object for output */ + amrex::Vector> m_outputSeries; /** openPMD backend: h5, bp, or json. Default depends on what is available */ std::string m_openpmd_backend = "default"; /** Last iteration that was written to file. * This is stored to make sure we don't write the last iteration multiple times. */ - int m_last_output_dumped = -1; + amrex::Vector m_last_output_dumped; /** vector of length nbeams with the numbers of particles already written to file */ amrex::Vector m_offset; @@ -114,8 +116,10 @@ public: * \param[in] output_step current iteration * \param[in] output_period output period * \param[in] max_step maximum time step of the simulation + * \param[in] nlev number of mesh refinement levels */ - void InitDiagnostics (const int output_step, const int output_period, const int max_step); + void InitDiagnostics (const int output_step, const int output_period, const int max_step, + const int nlev); /** \brief writing openPMD data * @@ -124,7 +128,7 @@ public: * \param[in] geom Geometry of the simulation, to get the cell size etc. * \param[in] physical_time Physical time of the currenerationt it. * \param[in] output_step current iteration to be written to file - * \param[in] lev MR level + * \param[in] nlev number of MR levels * \param[in] slice_dir direction of slicing. 0=x, 1=y, -1=no slicing (3D array) * \param[in] varnames list of variable names for the fields (ExmBy, EypBx, Ey, ...) * \param[in] beamnames list of the names of the beam to be written to file @@ -136,13 +140,14 @@ public: void WriteDiagnostics( amrex::Vector const& a_mf, MultiBeam& a_multi_beam, amrex::Vector const& geom, - const amrex::Real physical_time, const int output_step, const int lev, + const amrex::Real physical_time, const int output_step, const int nlev, const int slice_dir, const amrex::Vector< std::string > varnames, const amrex::Vector< std::string > beamnames, const int it, - const amrex::Vector& a_box_sorter_vec, const amrex::Geometry& geom3D, + const amrex::Vector& a_box_sorter_vec, amrex::Vector const& geom3D, const OpenPMDWriterCallType call_type); - void reset () {m_outputSeries.reset();}; + /** \brief Resets the openPMD series of all levels */ + void reset (); /** Prefix/path for the output files */ std::string m_file_prefix = "diags/hdf5"; diff --git a/src/diagnostics/OpenPMDWriter.cpp b/src/diagnostics/OpenPMDWriter.cpp index 2697894f6c..7c5c1063df 100644 --- a/src/diagnostics/OpenPMDWriter.cpp +++ b/src/diagnostics/OpenPMDWriter.cpp @@ -20,7 +20,8 @@ OpenPMDWriter::OpenPMDWriter () } void -OpenPMDWriter::InitDiagnostics (const int output_step, const int output_period, const int max_step) +OpenPMDWriter::InitDiagnostics (const int output_step, const int output_period, const int max_step, + const int nlev) { HIPACE_PROFILE("OpenPMDWriter::InitDiagnostics()"); @@ -28,8 +29,8 @@ OpenPMDWriter::InitDiagnostics (const int output_step, const int output_period, if (output_period < 0 || (!(output_step == max_step) && output_step % output_period != 0)) return; - // pick first available backend if default is chosen - if( m_openpmd_backend == "default" ) { + // pick first available backend if default is chosen + if( m_openpmd_backend == "default" ) { #if openPMD_HAVE_HDF5==1 m_openpmd_backend = "h5"; #elif openPMD_HAVE_ADIOS2==1 @@ -37,12 +38,26 @@ OpenPMDWriter::InitDiagnostics (const int output_step, const int output_period, #else m_openpmd_backend = "json"; #endif - } + } + + if (nlev > 1) { + for (int lev=0; lev( + filename, openPMD::Access::CREATE) ); + m_last_output_dumped.push_back(-1); + } + } else { + std::string filename = m_file_prefix + "/openpmd_%06T." + m_openpmd_backend; + + m_outputSeries.push_back(std::make_unique< openPMD::Series >( + filename, openPMD::Access::CREATE) ); + m_last_output_dumped.push_back(-1); + } - std::string filename = m_file_prefix + "/openpmd_%06T." + m_openpmd_backend; - m_outputSeries = std::make_unique< openPMD::Series >( - filename, openPMD::Access::CREATE); // TODO: meta-data: author, mesh path, extensions, software } @@ -51,32 +66,35 @@ void OpenPMDWriter::WriteDiagnostics ( amrex::Vector const& a_mf, MultiBeam& a_multi_beam, amrex::Vector const& geom, - const amrex::Real physical_time, const int output_step, const int lev, + const amrex::Real physical_time, const int output_step, const int nlev, const int slice_dir, const amrex::Vector< std::string > varnames, const amrex::Vector< std::string > beamnames, const int it, - const amrex::Vector& a_box_sorter_vec, const amrex::Geometry& geom3D, + const amrex::Vector& a_box_sorter_vec, amrex::Vector const& geom3D, const OpenPMDWriterCallType call_type) { - openPMD::Iteration iteration = m_outputSeries->iterations[output_step]; + for (int lev=0; leviterations[output_step]; - if (call_type == OpenPMDWriterCallType::beams ) { - iteration.setTime(physical_time); - WriteBeamParticleData(a_multi_beam, iteration, output_step, it, a_box_sorter_vec, geom3D, beamnames); - m_outputSeries->flush(); + if (call_type == OpenPMDWriterCallType::beams ) { + iteration.setTime(physical_time); + if (lev == 0) { + WriteBeamParticleData(a_multi_beam, iteration, output_step, it, a_box_sorter_vec, geom3D[lev], beamnames, lev); + } + m_outputSeries[lev]->flush(); - } else if (call_type == OpenPMDWriterCallType::fields ) { - WriteFieldData(a_mf[lev], geom[lev], slice_dir, varnames, iteration, output_step); - m_outputSeries->flush(); - m_last_output_dumped = output_step; + } else if (call_type == OpenPMDWriterCallType::fields ) { + WriteFieldData(a_mf[lev], geom[lev], slice_dir, varnames, iteration, output_step, lev); + m_outputSeries[lev]->flush(); + m_last_output_dumped[lev] = output_step; + } } - } void OpenPMDWriter::WriteFieldData ( amrex::FArrayBox const& fab, amrex::Geometry const& geom, const int slice_dir, const amrex::Vector< std::string > varnames, - openPMD::Iteration iteration, const int output_step) + openPMD::Iteration iteration, const int output_step, const int lev) { // todo: periodicity/boundary, field solver, particle pusher, etc. auto meshes = iteration.meshes; @@ -118,7 +136,7 @@ OpenPMDWriter::WriteFieldData ( // If slicing requested, remove number of points for the slicing direction if (slice_dir >= 0) global_size.erase(global_size.begin() + 2-slice_dir); - if (m_last_output_dumped != output_step) { + if (m_last_output_dumped[lev] != output_step) { openPMD::Dataset dataset(datatype, global_size); field_comp.resetDataset(dataset); } @@ -129,7 +147,6 @@ OpenPMDWriter::WriteFieldData ( data = openPMD::shareRaw( fab.dataPtr( icomp ) ); // non-owning view until flush() - // Determine the offset and size of this data chunk in the global output amrex::IntVect const box_offset = data_box.smallEnd(); openPMD::Offset chunk_offset = utils::getReversedVec(box_offset); @@ -148,7 +165,7 @@ OpenPMDWriter::WriteBeamParticleData (MultiBeam& beams, openPMD::Iteration itera const int output_step, const int it, const amrex::Vector& a_box_sorter_vec, const amrex::Geometry& geom, - const amrex::Vector< std::string > beamnames) + const amrex::Vector< std::string > beamnames, const int lev) { HIPACE_PROFILE("WriteBeamParticleData()"); @@ -165,7 +182,7 @@ OpenPMDWriter::WriteBeamParticleData (MultiBeam& beams, openPMD::Iteration itera auto& beam = beams.getBeam(ibeam); const unsigned long long np = beams.get_total_num_particles(ibeam); - if (m_last_output_dumped != output_step) { + if (m_last_output_dumped[lev] != output_step) { SetupPos(beam_species, beam, np, geom); SetupRealProperties(beam_species, m_real_names, np); } @@ -365,4 +382,11 @@ OpenPMDWriter::SaveRealProperty (BeamParticleContainer& pc, } } +void OpenPMDWriter::reset () +{ + for (int lev = 0; lev0) return; // The Poisson solver operates on transverse slices only. // The constructor takes the BoxArray and the DistributionMap of a slice, // so the FFTPlans are built on a slice.