vtkpp

.github/workflows/cmake_ubuntu.yml

@@ -57,7 +57,7 @@ jobs:
               -DCMAKE_C_COMPILER_LAUNCHER=ccache \
               -DCMAKE_CXX_COMPILER_LAUNCHER=ccache \
               -DlowResourceTests=ON -DdevMode=ON -Dbench=ON \
-              -DCMAKE_CXX_FLAGS="-O3 -DPHARE_DIAG_DOUBLES=1 " -Dphare_configurator=ON
+              -DCMAKE_CXX_FLAGS="-O3 -DPHARE_DIAG_DOUBLES=1" -Dphare_configurator=ON
 
     - name: Build
       working-directory: ${{github.workspace}}/build

pyphare/pyphare/pharein/__init__.py

@@ -328,6 +328,9 @@ def as_paths(rb):
             )
 
     if len(simulation.diagnostics) > 0:
+        if simulation.diag_options is not None and "format" in simulation.diag_options:
+            add_string(diag_path + "format", simulation.diag_options["format"])
+
         if simulation.diag_options is not None and "options" in simulation.diag_options:
             add_string(
                 diag_path + "filePath", simulation.diag_options["options"]["dir"]

pyphare/pyphare/pharein/diagnostics.py

@@ -62,9 +62,9 @@ def validate_timestamps(clazz, key, **kwargs):
     timestamps = phare_utilities.np_array_ify(kwargs.get(key, []))
 
     if np.any(timestamps < init_time):
-        raise RuntimeError(
-            f"Error: timestamp({sim.time_step_nbr}) cannot be less than simulation.init_time({init_time}))"
-        )
+        timestamps = timestamps[timestamps >= init_time]
+        print(f"Warning: some timestamps below ({init_time}) are filtered")
+
     if np.any(timestamps > sim.final_time):
         raise RuntimeError(
             f"Error: timestamp({sim.time_step_nbr}) cannot be greater than simulation.final_time({sim.final_time}))"

pyphare/pyphare/pharein/simulation.py

@@ -160,31 +160,39 @@ def check_time(**kwargs):
         and "time_step" not in kwargs
     )
 
+    if sum([final_and_dt, final_and_nsteps, nsteps_and_dt]) != 1:
+        raise ValueError(
+            "Error: Specify either 'final_time' and 'time_step' or 'time_step_nbr' and 'time_step'"
+            + " or 'final_time' and 'time_step_nbr'"
+        )
+
     start_time = kwargs.get("restart_options", {}).get("restart_time", 0)
 
+    def _final_time():
+        if "final_time" in kwargs:
+            return kwargs["final_time"]
+        return start_time + kwargs["time_step"] * kwargs["time_step_nbr"]
+
+    final_time = _final_time()
+    total_time = final_time - start_time
+    if total_time < 0:
+        raise RuntimeError("Simulation time cannot be negative - review inputs")
+    if phare_utilities.fp_equal(total_time, 0):
+        return 0, 0, final_time
+
     if final_and_dt:
-        time_step_nbr = int(kwargs["final_time"] / kwargs["time_step"])
-        time_step = kwargs["final_time"] / time_step_nbr
+        time_step_nbr = int(total_time / kwargs["time_step"])
+        time_step = total_time / time_step_nbr
 
     elif final_and_nsteps:
-        time_step = kwargs["final_time"] / kwargs["time_step_nbr"]
+        time_step = total_time / kwargs["time_step_nbr"]
         time_step_nbr = kwargs["time_step_nbr"]
 
-    elif nsteps_and_dt:
+    else:  # must be nsteps_and_dt
         time_step = kwargs["time_step"]
         time_step_nbr = kwargs["time_step_nbr"]
 
-    else:
-        raise ValueError(
-            "Error: Specify either 'final_time' and 'time_step' or 'time_step_nbr' and 'time_step'"
-            + " or 'final_time' and 'time_step_nbr'"
-        )
-
-    return (
-        time_step_nbr,
-        time_step,
-        kwargs.get("final_time", start_time + time_step * time_step_nbr),
-    )
+    return time_step_nbr, time_step, final_time
 
 
 # ------------------------------------------------------------------------------
@@ -485,7 +493,7 @@ def check_directory(directory, key):
 # diag_options = {"format":"phareh5", "options": {"dir": "phare_ouputs/"}}
 def check_diag_options(**kwargs):
     diag_options = kwargs.get("diag_options", None)
-    formats = ["phareh5"]
+    formats = ["phareh5", "pharevtkhdf"]
     if diag_options is not None and "format" in diag_options:
         if diag_options["format"] not in formats:
             raise ValueError("Error - diag_options format is invalid")
@@ -731,8 +739,10 @@ def wrapper(simulation_object, **kwargs):
                 kwargs["max_nbr_levels"],
             ) = check_refinement_boxes(ndim, **kwargs)
         else:
-            kwargs["max_nbr_levels"] = kwargs.get("max_nbr_levels", None)
-            assert kwargs["max_nbr_levels"] is not None  # this needs setting otherwise
+            if "max_nbr_levels" not in kwargs:
+                print("WARNING, 'max_nbr_levels' is not set, defaulting to 1")
+            kwargs["max_nbr_levels"] = kwargs.get("max_nbr_levels", 1)
+
             kwargs["refinement_boxes"] = None
             kwargs["tagging_threshold"] = kwargs.get("tagging_threshold", 0.1)
 

pyphare/pyphare/pharesee/phare_vtk/__init__.py

@@ -0,0 +1,7 @@
+#
+#
+#
+
+from .plot import plot
+
+__all__ = ["plot"]

pyphare/pyphare/pharesee/phare_vtk/base.py

@@ -0,0 +1,110 @@
+#
+#
+#
+
+import vtk
+
+
+class PhaseOutput:
+    def __init__(self, **kwargs):
+        self.kwargs = {**kwargs}
+
+    def __iter__(self):
+        return self.kwargs.items().__iter__()
+
+
+def surface_filter(output, **kwargs):
+    surface = vtk.vtkDataSetSurfaceFilter()
+    surface.SetInputConnection(output.GetOutputPort())
+    return PhaseOutput(surface=surface)
+
+
+def composite_data_geometry_filter(output, **kwargs):
+    geom = vtk.vtkCompositeDataGeometryFilter()
+    geom.SetInputConnection(output.GetOutputPort())
+    return PhaseOutput(geom=geom)
+
+
+def poly_data_mapper(output, **kwargs):
+    array_name = kwargs.get("array_name", "data")
+    mapper = vtk.vtkPolyDataMapper()
+    mapper.SetInputConnection(output.GetOutputPort())
+    mapper.SelectColorArray(array_name)
+    mapper.SetScalarModeToUsePointData()
+    mapper.SetScalarRange(
+        output.GetOutput().GetPointData().GetArray(array_name).GetRange(0)
+    )
+    return PhaseOutput(mapper=mapper)
+
+
+def all_times_in(reader):
+    info = reader.GetOutputInformation(0)
+
+    time_key = vtk.vtkStreamingDemandDrivenPipeline.TIME_STEPS()
+
+    if not info.Has(time_key):
+        raise RuntimeError("VTK Error: cannot get times from file")
+
+    times = info.Get(time_key)
+    return times
+
+
+def bounds_in(reader):
+    amr = reader.GetOutput()
+    bbox = vtk.vtkBoundingBox()
+
+    for level in range(amr.GetNumberOfLevels()):
+        for idx in range(amr.GetNumberOfDataSets(level)):
+            ds = amr.GetDataSet(level, idx)
+            if ds is not None:
+                b = [0.0] * 6
+                ds.GetBounds(b)
+                bbox.AddBounds(b)
+
+    global_bounds = [0.0] * 6
+    bbox.GetBounds(global_bounds)
+    return global_bounds
+
+
+def _default_phases():
+    # Typical use case for PHARE fields
+    return [surface_filter, composite_data_geometry_filter, poly_data_mapper]
+
+
+class VtkFile:
+    def __init__(self, filename, time=None, array_name="data", phases=None):
+        phases = phases if phases is not None else _default_phases()
+        if len(phases) == 0:
+            raise RuntimeError("Error: Zero phases!")
+
+        self.array_name = array_name
+        self.reader = vtk.vtkHDFReader()
+        self.reader.SetFileName(filename)
+        self.reader.Update()
+        self.reader.UpdateInformation()
+
+        self.times = all_times_in(self.reader)
+        self.reader.UpdateTimeStep(time if time is not None else self.times[-1])
+
+        _in = self.reader
+        for i in range(len(phases)):
+            ret = phases[i](_in, array_name=array_name)
+            for key, val in ret:
+                if hasattr(val, "Update"):
+                    val.Update()
+                setattr(self, key, val)
+                _in = val
+
+        self.mapper = _in
+        self.bounds = bounds_in(self.reader)
+        self.spacing = self.reader.GetOutput().GetDataSet(0, 0).GetSpacing()
+
+
+class VtkFieldFile(VtkFile):
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)
+
+
+class VtkTensorFieldFile(VtkFile):
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)

pyphare/pyphare/pharesee/phare_vtk/plot.py

@@ -0,0 +1,38 @@
+#
+#
+#
+
+
+import vtk
+from .base import VtkTensorFieldFile
+
+
+def plot(vtk_file, out_file="vtk.png"):
+    if isinstance(vtk_file, str):
+        vtk_file = VtkTensorFieldFile(vtk_file)
+    vtk_file.geom.GetOutput().GetPointData().SetActiveScalars("data")
+
+    colors = vtk.vtkNamedColors()
+
+    renderer = vtk.vtkRenderer()
+    renderer.AddActor(vtk.vtkActor(mapper=vtk_file.mapper))
+    renderer.SetBackground(colors.GetColor3d("Silver"))
+
+    ren_win = vtk.vtkRenderWindow()
+    ren_win.AddRenderer(renderer)
+    ren_win.SetSize(800, 600)
+    ren_win.OffScreenRenderingOn()  # do not flash image on screen
+
+    iren = vtk.vtkRenderWindowInteractor()
+    iren.SetRenderWindow(ren_win)
+
+    ren_win.Render()
+
+    w2if = vtk.vtkWindowToImageFilter()
+    w2if.SetInput(ren_win)
+    w2if.Update()
+
+    writer = vtk.vtkPNGWriter()
+    writer.SetFileName(out_file)
+    writer.SetInputConnection(w2if.GetOutputPort())
+    writer.Write()

pyphare/pyphare/pharesee/plotting.py

@@ -325,7 +325,7 @@ def finest_field_plot(run_path, qty, **kwargs):
         ax.set_aspect("equal")
         divider = make_axes_locatable(ax)
         cax = divider.append_axes("right", size="5%", pad=0.08)
-        cb = plt.colorbar(im, ax=ax, cax=cax)
+        plt.colorbar(im, ax=ax, cax=cax)
     else:
         raise ValueError("finest_field_plot not yet ready for 3d")
 

pyphare/pyphare_tests/test_pharesee/test_geometry.py

@@ -116,7 +116,6 @@ def test_particle_ghost_area_boxes(self):
         }
 
         gaboxes = ghost_area_boxes(hierarchy, "particles")
-        particles = "particles"
 
         # same number of levels
         self.assertEqual(len(expected), len(gaboxes))