From 94978e9173a11088e0d839e9b70544e9208c3d84 Mon Sep 17 00:00:00 2001 From: p Date: Fri, 11 Jul 2025 13:17:27 +0200 Subject: [PATCH] pressureplot --- pyphare/pyphare/pharesee/run/utils.py | 10 ++--- src/diagnostic/detail/h5typewriter.hpp | 9 +++++ src/diagnostic/detail/types/electromag.hpp | 8 ++-- src/diagnostic/detail/types/fluid.hpp | 6 +-- src/diagnostic/detail/types/info.hpp | 2 +- src/diagnostic/detail/types/meta.hpp | 6 +-- src/diagnostic/detail/types/particle.hpp | 6 +-- tests/functional/harris/harris_2d.py | 21 ++++++++-- tests/simulator/test_run.py | 45 ++++++++++++++++------ 9 files changed, 79 insertions(+), 34 deletions(-) diff --git a/pyphare/pyphare/pharesee/run/utils.py b/pyphare/pyphare/pharesee/run/utils.py index e9cb5b57d..d6ffac24c 100644 --- a/pyphare/pyphare/pharesee/run/utils.py +++ b/pyphare/pyphare/pharesee/run/utils.py @@ -378,7 +378,7 @@ def _compute_pressure(patch_datas, **kwargs): Myy = patch_datas["Myy"].dataset[:] Myz = patch_datas["Myz"].dataset[:] Mzz = patch_datas["Mzz"].dataset[:] - massDensity = patch_datas["rho"].dataset[:] + massDensity = patch_datas["value"].dataset[:] Vix = patch_datas["Vx"].dataset[:] Viy = patch_datas["Vy"].dataset[:] Viz = patch_datas["Vz"].dataset[:] @@ -417,10 +417,10 @@ def _compute_pop_pressure(patch_datas, **kwargs): Myy = patch_datas[popname + "_Myy"].dataset[:] Myz = patch_datas[popname + "_Myz"].dataset[:] Mzz = patch_datas[popname + "_Mzz"].dataset[:] - Fx = patch_datas[popname + "_Fx"].dataset[:] - Fy = patch_datas[popname + "_Fy"].dataset[:] - Fz = patch_datas[popname + "_Fz"].dataset[:] - N = patch_datas[popname + "_rho"].dataset[:] + Fx = patch_datas["x"].dataset[:] + Fy = patch_datas["y"].dataset[:] + Fz = patch_datas["z"].dataset[:] + N = patch_datas["value"].dataset[:] mass = kwargs["mass"] diff --git a/src/diagnostic/detail/h5typewriter.hpp b/src/diagnostic/detail/h5typewriter.hpp index c08157626..57d450eb0 100644 --- a/src/diagnostic/detail/h5typewriter.hpp +++ b/src/diagnostic/detail/h5typewriter.hpp @@ -142,6 +142,15 @@ class H5TypeWriter : public PHARE::diagnostic::TypeWriter } + auto& h5FileForQuantity(DiagnosticProperties& diagnostic) + { + if (!fileData_.count(diagnostic.quantity)) + throw std::runtime_error("Unknown Diagnostic Quantity: " + diagnostic.quantity); + + return *fileData_.at(diagnostic.quantity); + } + + Writer& h5Writer_; std::unordered_map> fileData_; }; diff --git a/src/diagnostic/detail/types/electromag.hpp b/src/diagnostic/detail/types/electromag.hpp index a640bb91e..96f1c3048 100644 --- a/src/diagnostic/detail/types/electromag.hpp +++ b/src/diagnostic/detail/types/electromag.hpp @@ -74,7 +74,7 @@ void ElectromagDiagnosticWriter::getDataSetInfo(DiagnosticProperties& // highfive doesn't accept uint32 which ndarray.shape() is auto const& array_shape = vecF.getComponent(type).shape(); attr[name][id] = std::vector(array_shape.data(), - array_shape.data() + array_shape.size()); + array_shape.data() + array_shape.size()); auto ghosts = GridLayout::nDNbrGhosts(vecF.getComponent(type).physicalQuantity()); attr[name][id + "_ghosts_x"] = static_cast(ghosts[0]); if constexpr (GridLayout::dimension > 1) @@ -100,7 +100,7 @@ void ElectromagDiagnosticWriter::initDataSets( Attributes& patchAttributes, std::size_t maxLevel) { auto& h5Writer = this->h5Writer_; - auto& h5file = *fileData_.at(diagnostic.quantity); + auto& h5file = Super::h5FileForQuantity(diagnostic); auto vecFields = h5Writer.modelView().getElectromagFields(); auto initVF = [&](auto& path, auto& attr, std::string key, auto null) { @@ -151,7 +151,7 @@ void ElectromagDiagnosticWriter::write(DiagnosticProperties& diagnosti for (auto* vecField : h5Writer.modelView().getElectromagFields()) if (diagnostic.quantity == "/" + vecField->name()) - h5Writer.writeTensorFieldAsDataset(*fileData_.at(diagnostic.quantity), + h5Writer.writeTensorFieldAsDataset(Super::h5FileForQuantity(diagnostic), h5Writer.patchPath() + "/" + vecField->name(), *vecField); } @@ -165,7 +165,7 @@ void ElectromagDiagnosticWriter::writeAttributes( patchAttributes, std::size_t maxLevel) { - writeAttributes_(diagnostic, *fileData_.at(diagnostic.quantity), fileAttributes, + writeAttributes_(diagnostic, Super::h5FileForQuantity(diagnostic), fileAttributes, patchAttributes, maxLevel); } diff --git a/src/diagnostic/detail/types/fluid.hpp b/src/diagnostic/detail/types/fluid.hpp index 89839c5d7..7221faa82 100644 --- a/src/diagnostic/detail/types/fluid.hpp +++ b/src/diagnostic/detail/types/fluid.hpp @@ -230,7 +230,7 @@ void FluidDiagnosticWriter::initDataSets( { auto& h5Writer = this->h5Writer_; auto& ions = h5Writer.modelView().getIons(); - auto& h5file = *fileData_.at(diagnostic.quantity); + auto& h5file = Super::h5FileForQuantity(diagnostic); auto writeGhosts = [&](auto& path, auto& attr, std::string key, auto null) { this->writeGhostsAttr_(h5file, path, @@ -299,7 +299,7 @@ void FluidDiagnosticWriter::write(DiagnosticProperties& diagnostic) { auto& h5Writer = this->h5Writer_; auto& ions = h5Writer.modelView().getIons(); - auto& h5file = *fileData_.at(diagnostic.quantity); + auto& h5file = Super::h5FileForQuantity(diagnostic); auto writeDS = [&](auto path, auto& field) { h5file.template write_data_set_flat(path, field.data()); @@ -341,7 +341,7 @@ void FluidDiagnosticWriter::writeAttributes( std::size_t maxLevel) { auto& h5Writer = this->h5Writer_; - auto& h5file = *fileData_.at(diagnostic.quantity); + auto& h5file = Super::h5FileForQuantity(diagnostic); auto checkWrite = [&](auto& tree, std::string qty, auto const& pop) { if (diagnostic.quantity == tree + qty) diff --git a/src/diagnostic/detail/types/info.hpp b/src/diagnostic/detail/types/info.hpp index 9f177029e..09c94fdec 100644 --- a/src/diagnostic/detail/types/info.hpp +++ b/src/diagnostic/detail/types/info.hpp @@ -113,7 +113,7 @@ void InfoDiagnosticWriter::writeAttributes( defaultPatchAttributes["particle_count"] = std::size_t{0}; } - writeAttributes_(diagnostic, *fileData_.at(diagnostic.quantity), fileAttributes, + writeAttributes_(diagnostic, Super::h5FileForQuantity(diagnostic), fileAttributes, patchAttributes, maxLevel, defaultPatchAttributes); } diff --git a/src/diagnostic/detail/types/meta.hpp b/src/diagnostic/detail/types/meta.hpp index 8ca27a85f..45652c258 100644 --- a/src/diagnostic/detail/types/meta.hpp +++ b/src/diagnostic/detail/types/meta.hpp @@ -98,7 +98,7 @@ void MetaDiagnosticWriter::initDataSets( if (diagnostic.quantity == "/tags") h5Writer.template createDataSet( - *fileData_.at(diagnostic.quantity), path + "/tags", + Super::h5FileForQuantity(diagnostic), path + "/tags", null or tags.count(path) == 0 ? std::vector(GridLayout::dimension, 0) : attr["tags"].template to>()); @@ -120,7 +120,7 @@ void MetaDiagnosticWriter::write(DiagnosticProperties& diagnostic) if (tags.count(path) > 0) { - auto& h5 = *fileData_.at(diagnostic.quantity); + auto& h5 = Super::h5FileForQuantity(diagnostic); h5.template write_data_set_flat(path + "/tags", tags[path]); tags.erase(path); } @@ -135,7 +135,7 @@ void MetaDiagnosticWriter::writeAttributes( patchAttributes, std::size_t maxLevel) { - writeAttributes_(diagnostic, *fileData_.at(diagnostic.quantity), fileAttributes, + writeAttributes_(diagnostic, Super::h5FileForQuantity(diagnostic), fileAttributes, patchAttributes, maxLevel); } diff --git a/src/diagnostic/detail/types/particle.hpp b/src/diagnostic/detail/types/particle.hpp index a85402b2b..72555f7b6 100644 --- a/src/diagnostic/detail/types/particle.hpp +++ b/src/diagnostic/detail/types/particle.hpp @@ -114,7 +114,7 @@ void ParticlesDiagnosticWriter::initDataSets( Attributes& patchAttributes, std::size_t maxLevel) { auto& h5Writer = this->h5Writer_; - auto& h5file = *fileData_.at(diagnostic.quantity); + auto& h5file = Super::h5FileForQuantity(diagnostic); auto createDataSet = [&](auto&& path, auto& attr, auto& key, auto& value, auto null) { using ValueType = std::decay_t; @@ -171,7 +171,7 @@ void ParticlesDiagnosticWriter::write(DiagnosticProperties& diagnostic auto checkWrite = [&](auto& tree, auto pType, auto& ps) { std::string active{tree + pType}; if (diagnostic.quantity == active && ps.size() > 0) - hdf5::ParticleWriter::write(*fileData_.at(diagnostic.quantity), ps, + hdf5::ParticleWriter::write(this->h5FileForQuantity(diagnostic), ps, h5Writer.patchPath() + "/"); }; @@ -193,7 +193,7 @@ void ParticlesDiagnosticWriter::writeAttributes( std::size_t maxLevel) { auto& h5Writer = this->h5Writer_; - auto& h5file = *fileData_.at(diagnostic.quantity); + auto& h5file = Super::h5FileForQuantity(diagnostic); auto checkWrite = [&](auto& tree, std::string pType, auto const& pop) { if (diagnostic.quantity == tree + pType) diff --git a/tests/functional/harris/harris_2d.py b/tests/functional/harris/harris_2d.py index 21dad469e..6054b6e3c 100644 --- a/tests/functional/harris/harris_2d.py +++ b/tests/functional/harris/harris_2d.py @@ -139,9 +139,11 @@ def vthz(x, y): for quantity in ["mass_density", "bulkVelocity"]: ph.FluidDiagnostics(quantity=quantity, write_timestamps=timestamps) - ph.FluidDiagnostics( - quantity="density", write_timestamps=timestamps, population_name="protons" - ) + for quantity in ["density", "pressure_tensor"]: + ph.FluidDiagnostics( + quantity=quantity, write_timestamps=timestamps, population_name="protons" + ) + ph.InfoDiagnostics(quantity="particle_count") ph.LoadBalancer(active=True, auto=True, mode="nppc", tol=0.05) @@ -155,6 +157,7 @@ def plot_file_for_qty(plot_dir, qty, time): def plot(diag_dir, plot_dir): run = Run(diag_dir) + pop_name = "protons" for time in timestamps: run.GetDivB(time).plot( filename=plot_file_for_qty(plot_dir, "divb", time), @@ -165,7 +168,7 @@ def plot(diag_dir, plot_dir): run.GetRanks(time).plot( filename=plot_file_for_qty(plot_dir, "Ranks", time), plot_patches=True ) - run.GetN(time, pop_name="protons").plot( + run.GetN(time, pop_name=pop_name).plot( filename=plot_file_for_qty(plot_dir, "N", time), plot_patches=True ) for c in ["x", "y", "z"]: @@ -181,6 +184,16 @@ def plot(diag_dir, plot_dir): vmin=-2, vmax=2, ) + run.GetPressure(time, pop_name=pop_name).plot( + filename=plot_file_for_qty(plot_dir, "Pxx", time), + qty=pop_name + "_Pxx", + plot_patches=True, + ) + run.GetPressure(time, pop_name=pop_name).plot( + filename=plot_file_for_qty(plot_dir, "Pzz", time), + qty=pop_name + "_Pzz", + plot_patches=True, + ) class HarrisTest(SimulatorTest): diff --git a/tests/simulator/test_run.py b/tests/simulator/test_run.py index f6a0a785f..c7ba7b46d 100644 --- a/tests/simulator/test_run.py +++ b/tests/simulator/test_run.py @@ -142,22 +142,24 @@ def vthz(x, y): for quantity in ["E", "B"]: ph.ElectromagDiagnostics(quantity=quantity, write_timestamps=timestamps) - for quantity in ["mass_density", "charge_density", "bulkVelocity"]: + + for quantity in [ + "mass_density", + "charge_density", + "bulkVelocity", + "pressure_tensor", + ]: ph.FluidDiagnostics(quantity=quantity, write_timestamps=timestamps) pop = "protons" ph.ParticleDiagnostics( - quantity="domain", - write_timestamps=timestamps, - population_name=pop, - ) - ph.FluidDiagnostics( - quantity="density", write_timestamps=timestamps, population_name=pop - ) - ph.FluidDiagnostics( - quantity="charge_density", write_timestamps=timestamps, population_name=pop + quantity="domain", write_timestamps=timestamps, population_name=pop ) + for quantity in ["density", "charge_density", "pressure_tensor"]: + ph.FluidDiagnostics( + quantity=quantity, write_timestamps=timestamps, population_name=pop + ) return sim @@ -167,6 +169,7 @@ def plot_file_for_qty(qty, time): def plot(diag_dir): run = Run(diag_dir) + pop_name = "protons" for time in timestamps: run.GetDivB(time).plot( filename=plot_file_for_qty("divb", time), @@ -177,7 +180,7 @@ def plot(diag_dir): run.GetRanks(time).plot( filename=plot_file_for_qty("Ranks", time), plot_patches=True ) - run.GetN(time, pop_name="protons").plot( + run.GetN(time, pop_name=pop_name).plot( filename=plot_file_for_qty("N", time), plot_patches=True ) for c in ["x", "y", "z"]: @@ -194,6 +197,26 @@ def plot(diag_dir): vmin=-2, vmax=2, ) + run.GetPressure(time, pop_name=pop_name).plot( + filename=plot_file_for_qty(f"{pop_name}_Pxx", time), + qty=pop_name + "_Pxx", + plot_patches=True, + ) + run.GetPressure(time, pop_name=pop_name).plot( + filename=plot_file_for_qty(f"{pop_name}_Pzz", time), + qty=pop_name + "_Pzz", + plot_patches=True, + ) + run.GetPi(time).plot( + filename=plot_file_for_qty("Pxx", time), + qty="Pxx", + plot_patches=True, + ) + run.GetPi(time).plot( + filename=plot_file_for_qty("Pzz", time), + qty="Pzz", + plot_patches=True, + ) def assert_file_exists_with_size_at_least(file, size=10000):