pressureplot

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"]
 

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<std::string, std::unique_ptr<HighFiveFile>> fileData_;
 };

src/diagnostic/detail/types/electromag.hpp

@@ -74,7 +74,7 @@ void ElectromagDiagnosticWriter<H5Writer>::getDataSetInfo(DiagnosticProperties&
             // highfive doesn't accept uint32 which ndarray.shape() is
             auto const& array_shape = vecF.getComponent(type).shape();
             attr[name][id]          = std::vector<std::size_t>(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<std::size_t>(ghosts[0]);
             if constexpr (GridLayout::dimension > 1)
@@ -100,7 +100,7 @@ void ElectromagDiagnosticWriter<H5Writer>::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<H5Writer>::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<H5Writer>::writeAttributes(
         patchAttributes,
     std::size_t maxLevel)
 {
-    writeAttributes_(diagnostic, *fileData_.at(diagnostic.quantity), fileAttributes,
+    writeAttributes_(diagnostic, Super::h5FileForQuantity(diagnostic), fileAttributes,
                      patchAttributes, maxLevel);
 }
 

src/diagnostic/detail/types/fluid.hpp

@@ -230,7 +230,7 @@ void FluidDiagnosticWriter<H5Writer>::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<H5Writer>::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<GridLayout::dimension>(path, field.data());
@@ -341,7 +341,7 @@ void FluidDiagnosticWriter<H5Writer>::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)

src/diagnostic/detail/types/info.hpp

@@ -113,7 +113,7 @@ void InfoDiagnosticWriter<H5Writer>::writeAttributes(
         defaultPatchAttributes["particle_count"] = std::size_t{0};
     }
 
-    writeAttributes_(diagnostic, *fileData_.at(diagnostic.quantity), fileAttributes,
+    writeAttributes_(diagnostic, Super::h5FileForQuantity(diagnostic), fileAttributes,
                      patchAttributes, maxLevel, defaultPatchAttributes);
 }
 

src/diagnostic/detail/types/meta.hpp

@@ -98,7 +98,7 @@ void MetaDiagnosticWriter<H5Writer>::initDataSets(
 
         if (diagnostic.quantity == "/tags")
             h5Writer.template createDataSet<int>(
-                *fileData_.at(diagnostic.quantity), path + "/tags",
+                Super::h5FileForQuantity(diagnostic), path + "/tags",
                 null or tags.count(path) == 0
                     ? std::vector<std::size_t>(GridLayout::dimension, 0)
                     : attr["tags"].template to<std::vector<std::size_t>>());
@@ -120,7 +120,7 @@ void MetaDiagnosticWriter<H5Writer>::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<GridLayout::dimension>(path + "/tags", tags[path]);
             tags.erase(path);
         }
@@ -135,7 +135,7 @@ void MetaDiagnosticWriter<H5Writer>::writeAttributes(
         patchAttributes,
     std::size_t maxLevel)
 {
-    writeAttributes_(diagnostic, *fileData_.at(diagnostic.quantity), fileAttributes,
+    writeAttributes_(diagnostic, Super::h5FileForQuantity(diagnostic), fileAttributes,
                      patchAttributes, maxLevel);
 }
 

src/diagnostic/detail/types/particle.hpp

@@ -114,7 +114,7 @@ void ParticlesDiagnosticWriter<H5Writer>::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<decltype(value)>;
@@ -171,7 +171,7 @@ void ParticlesDiagnosticWriter<H5Writer>::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<H5Writer>::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)

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):

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):