[feature] add trace_over to superimpose a particle trajectory to a heatmap

.github/workflows/docs.yaml

@@ -3,6 +3,9 @@ on:
   push:
     branches:
       - main
+  pull_request:
+    branches:
+      - main
 
 permissions:
   contents: write
@@ -30,15 +33,22 @@ jobs:
         run: |
           for f in examples/*.py; do 
             echo "Processing $f..."
-            python "$f" -j 2
+            python "$f" -j 2 || exit 1
           done
 
       - name: Build Sphinx
         run: |
-          sphinx-build -b html docs/source docs/_build
+          sphinx-build -W -b html docs/source docs/_build/html
+
+      - name: Upload Documentation Artifact
+        uses: actions/upload-artifact@v4
+        with:
+          name: docs-preview
+          path: ./docs/_build/html
 
       - name: Deploy to GitHub Pages
+        if: github.event_name == 'push' && github.ref == 'refs/heads/main'
         uses: peaceiris/actions-gh-pages@v3
         with:
           github_token: ${{ secrets.GITHUB_TOKEN }}
-          publish_dir: ./docs/_build
+          publish_dir: ./docs/_build/html

README.md

@@ -9,22 +9,22 @@ The module is built around 4 components:
 +  **`RunContext`**: Handles IO, directory creation, and data discovery. Detects simulation geometry, dimensions, and available fields.
 +  **`PhysicsProcessor`**: Performs mathematical transformations. Handles grid conversions (e.g., converting internal coordinates to Cartesian $x, z$ for plotting), applies zooms, and computes derived quantities.
 +  **`SliceRenderer`**: Matplotlib engine. Manages multi-panel layouts, colorbars (Log, Linear, TwoSlope), streamlines, contours.
-+  **`Pipeline`**: The coordinator:
-    *   Pre-calculates global bounds across all frames.
-    *   Distributes rendering tasks across multiple CPU cores.
-    *   Manages the temporal evolution of SpaceTime heatmaps and particle data.
++  **`Pipeline`**: Coordinates the detection, processing, and rendering of the simulation data. Computes the bounds and distributes with `multiprocessing`.
 
 ## Supported quantity types
 
 Users define what to plot by passing lists of "Quantity" objects to the Pipeline:
 
-*   **`MapMovie2D`**: For $f(x, z, t)$ fields, 2D pcolormesh plots generated for every frame. Supports:
+*   **`MapMovie2D`**: For $f(x, z, t)$ fields that will be rendered as a heatmap (pcolormesh) animation. Supports:
     *   `streamlines`: Vector overlays (e.g., velocity fields).
-    *   `contours`: Scalar overlays (e.g., density levels).
+    *   `contours`: draw contour lines over the pcolormesh. (e.g., density levels).
     *   `compute`: Custom functions to calculate new fields on the fly.
-*   **`LineMovie1D`**: For $f(x,t)$ fields, 1D plot generated for every frame (e.g., radial profile over time).
-*   **`SpaceTimeHeatmap`**: For $f(x,t)$ fields, generates a space-time heapmap. Supports `ref_function` to overlay analytical trajectories.
-*   **`PartQuantity`**: Tracks Lagrangian particle properties (like `PART_X1`) over time. Supports `ref_function` to overlay analytical functions.
+*   **`LineMovie1D`**: For $f(x,t)$ fields, renders a line plot that updates every frame. (e.g., radial profile over time).
+*   **`SpaceTimeHeatmap`**: For $f(x,t)$ fields, renders a space-time heatmap. Supports `ref_function` to overlay analytical trajectories.
+*   **`PartQuantity`**: Tracks Lagrangian particle properties (like `PART_X1`) over time. 
+
+`SpaceTimeHeatmap` and `PartQuantity` supports `ref_function` to overlay analytical functions.
+One can also plot a particle quantity over a spacetime heatmap with `trace_over`
 
 ## Usage
 

docs/source/api.rst

@@ -7,8 +7,6 @@ RunContext
 ----------
 .. autoclass:: RunContext
    :members:
-   :undoc-members:
-   :show-inheritance:
 
 Pipeline
 --------
@@ -19,6 +17,10 @@ Pipeline
 Quantity Types
 --------------
 .. autoclass:: MapMovie2D
-   :members:
+   :show-inheritance:
+.. autoclass:: LineMovie1D
+   :show-inheritance:
 .. autoclass:: SpaceTimeHeatmap
-   :members:
+   :show-inheritance:
+.. autoclass:: PartQuantity
+   :show-inheritance:

docs/source/examples.rst

@@ -1,7 +1,7 @@
 Examples Gallery
 ================
 
-Select an example below to explore the capabilities of the pipeline.
+Here are three most basic examples to use the Pipeline, depending on the nature of the quantity: 1D, 2D, or particles.
 
 .. grid:: 1 2 2 3
    :gutter: 3
@@ -32,4 +32,17 @@ Select an example below to explore the capabilities of the pipeline.
 
    examples/run_1D
    examples/run_2D
-   examples/run_particles
+   examples/run_particles
+   examples/run_particles_over
+
+Here are more complex examples
+
+.. grid:: 1 2 2 3
+   :gutter: 3
+
+   .. grid-item-card:: Particle over a heatmap
+      :link: examples/run_particles_over
+      :link-type: doc
+      :img-top: ../../examples/data_test/frames/particles_over_test/particles_over_test_spacetimeheatmap.png
+
+      Superimpose a particle trajectory over a :math:`(x,t)` heatmap

docs/source/examples/run_1D.rst

@@ -41,6 +41,8 @@ On the heatmap, one can add a reference line (e.g., an analytical trajectory) by
    :language: python
    :linenos:
 
-.. image:: ../../../examples/data_test/frames/1D_test_withref/1D_test_spacetimeheatmap.png
+.. figure:: ../../../examples/data_test/frames/1D_test_withref/1D_test_spacetimeheatmap.png
    :width: 80%
-   :align: center
+   :align: center
+
+   1D_test_spacetimeheatmap.png

docs/source/examples/run_particles.rst

@@ -5,7 +5,9 @@ Particle quantities
    :language: python
    :linenos:
 
-.. image:: ../../../examples/data_test/frames/particles_test/particles_test_timeseries.png
+.. figure:: ../../../examples/data_test/frames/particles_test/particles_test_timeseries.png
     :width: 80%
     :align: center
 
+    particles_test_timeseries.png
+

docs/source/examples/run_particles_over.rst

@@ -0,0 +1,13 @@
+Particle quantities over heatmap
+================================
+
+.. literalinclude:: ../../../examples/run_particles_over.py
+   :language: python
+   :linenos:
+
+.. figure:: ../../../examples/data_test/frames/particles_over_test/particles_over_test_spacetimeheatmap.png
+    :width: 80%
+    :align: center
+
+    particles_over_test_spacetimeheatmap.png
+

examples/data_test/inputs/particles_over_test.ini

@@ -0,0 +1,55 @@
+##
+[Grid]
+X1-grid    1  1.99      1  u  2.01
+X2-grid    1  -0.0125      1  u  0.0125
+X3-grid    1  -0.3  512    u  0.3
+
+[TimeIntegrator]
+CFL         0.5
+tstop       20
+first_dt    1.0e-5
+nstages     2
+
+[Hydro]
+solver    hllc
+csiso     userdef
+# viscosity    explicit  userdef
+
+[Dust]
+nSpecies         3
+drag             userdef  1   0.2   0.04    # St=1, 0.2, 0.04
+drag_feedback    no
+
+[Particles]
+count            per_proc  1
+stopping_time    constant  1.0
+ParticleMass     3e-3
+DustToGas        3e-3
+
+[Gravity]
+potential    central 
+Mcentral     1.0
+
+[Boundary]
+X1-beg    userdef
+X1-end    userdef
+X2-beg    periodic
+X2-end    periodic
+X3-beg    userdef
+X3-end    userdef
+
+[Setup]
+sigma0        0.125
+sigmaSlope    -0.5
+CsSlope       -0.5
+h0            0.05
+alpha         1.0e-4
+
+[Output]
+vtk    0.1
+dmp    1000.0
+log        1000
+dmp_dir    /home/dp316/dp316/dc-fang1/IdefixRuns/VerticalSettling/outputs/SettlingL_Tau
+vtk_dir    /home/dp316/dp316/dc-fang1/IdefixRuns/VerticalSettling/outputs/SettlingL_Tau/vtks
+# File produced automatically by a Perl script
+# Do not edit

examples/run_particles_over.py

@@ -0,0 +1,46 @@
+from idefix2python import RunContext, Pipeline, PartQuantity, SpaceTimeHeatmap
+import utilities
+from pathlib import Path
+
+projectPath = Path(__file__).parent / "data_test"
+task = "particles_over_test"
+# By default the vtks are expected to be in {projetPath}/{task}/outputs/vtks/
+# In this example, the vtks/ folder contains both part*.vtk and data*.vtk
+
+
+def analytical_trajectory(t):
+    z0 = 0.1
+    fluid = utilities.Fluid(0.05, -0.5, 0.125, -0.5, Stokes0=1, z0=z0)
+    return utilities.solve_2nd_order_ode(fluid.azSettling, z0, 0, t)
+
+
+analytical_trajectory.plot_kwargs = {"ls": "--", "color": "cyan", "lw": 2}
+
+z_part = PartQuantity(
+    "PART_X3",
+    r"$z^\mathrm{part}$",
+    plot_coords=[0, 0],
+)
+
+custom_partQuantities = [z_part]
+SpaceTimeHeatmaps = [
+    SpaceTimeHeatmap(
+        "Dust0_RHO",
+        r"$\rho^\mathrm{dust}$",
+        plot_coords=[0, 0],
+        trace_over=[z_part],
+        ref_function=analytical_trajectory,
+    )
+]
+runContext = RunContext(
+    task,
+    projectPath,
+    active_directions=[2],  # currently necessary for lagrangian particles.)
+)
+pipeline = Pipeline(
+    runContext,
+    partQuantities=custom_partQuantities,
+    spaceTimeHeatmaps=SpaceTimeHeatmaps,
+)
+
+pipeline.run()

examples/utilities.py

@@ -25,7 +25,7 @@ def vK(r):
 
 
 class Fluid:
-    def __init__(self, cs0, csSlope, sigma0, sigmaSlope, Stokes0):
+    def __init__(self, cs0, csSlope, sigma0, sigmaSlope, Stokes0, z0=0.1):
         self.cs0 = cs0
         self.csSlope = csSlope
         self.sigma0 = sigma0
@@ -39,6 +39,8 @@ def __init__(self, cs0, csSlope, sigma0, sigmaSlope, Stokes0):
         # self.rhoSlope = -1.5
         self.Stokes0 = Stokes0
 
+        self.z0 = z0
+
     def eta(self, r):
         cs0 = self.cs0
         csSlope = self.csSlope
@@ -58,6 +60,26 @@ def vrDrift(self, r):
         st = self.Stokes(r)
         return self.eta(r) * vK(r) / (st + 1 / st)
 
+    def vzSettling_approx(self, z):
+        r = 2.0
+        st = self.Stokes(r)
+        OmegaK = r ** (-1.5)
+        return -OmegaK * st * z
+
+    def azSettling(self, z, vz, t):
+        r = 2.0
+        tstop = self.Stokes0
+        OmegaK = r ** (-1.5)
+        return -vz / tstop - z * OmegaK**2
+
+    # def z_drift(self, t):
+    #     r = 2.0
+    #     tstop = self.Stokes0
+    #     OmegaK = r ** (-1.5)
+    #     return self.z0 * np.exp(OmegaK*st*z**2/2)
+    #     z = self.z0 * np.exp(-t / (2 * tstop))
+    #     st = self.Stokes(r)
+
     # def eta(self, r):
     #     cs0 = self.cs0
     #     csSlope = self.csSlope
@@ -72,3 +94,19 @@ def vrDrift(self, r):
 
     # def vrDrift(self, r):
     #     return -self.eta(r) / np.sqrt(r) / (self.Stokes(r) + 1 / self.Stokes(r))
+
+
+def solve_2nd_order_ode(f, u0, du0, times):
+    """
+    u" = f(u',u,t)
+    """
+    u = [u0]
+    du = [du0]
+    dt = np.diff(times)
+    for i, t in enumerate(times[:-1]):
+        u_prev = u[-1]
+        du_prev = du[-1]
+
+        du += [du_prev + f(u_prev, du_prev, t) * dt[i]]
+        u += [u_prev + du_prev * dt[i]]
+    return u