diff --git a/.typos.toml b/.typos.toml index 432385eef0..595090c271 100644 --- a/.typos.toml +++ b/.typos.toml @@ -28,6 +28,7 @@ choises = "choises" # appears in comment explaining validation purpose ordr = "ordr" # typo for "order" in "weno_ordr" - tests param suggestions unknwn = "unknwn" # typo for "unknown" - tests unknown param detection tru = "tru" # typo for "true" in "when_tru" - tests dependency keys +PNGs = "PNGs" [files] -extend-exclude = ["docs/documentation/references*", "docs/references.bib", "tests/", "toolchain/cce_simulation_workgroup_256.sh", "build-docs/"] +extend-exclude = ["docs/documentation/references*", "docs/references.bib", "tests/", "toolchain/cce_simulation_workgroup_256.sh", "build-docs/", "build/", "build_test/"] diff --git a/docs/documentation/case.md b/docs/documentation/case.md index 356f6d90fc..2f884ab94b 100644 --- a/docs/documentation/case.md +++ b/docs/documentation/case.md @@ -651,7 +651,7 @@ To restart the simulation from $k$-th time step, see @ref running "Restarting Ca The table lists formatted database output parameters. The parameters define variables that are outputted from simulation and file types and formats of data as well as options for post-processing. -- `format` specifies the choice of the file format of data file outputted by MFC by an integer of 1 and 2. `format = 1` and `2` correspond to Silo-HDF5 format and binary format, respectively. +- `format` specifies the choice of the file format of data file outputted by MFC by an integer of 1 and 2. `format = 1` and `2` correspond to Silo-HDF5 format and binary format, respectively. Both formats are supported by `./mfc.sh viz` (see @ref visualization "Flow Visualization"). Silo-HDF5 requires the h5py Python package; binary has no extra dependencies. - `precision` specifies the choice of the floating-point format of the data file outputted by MFC by an integer of 1 and 2. `precision = 1` and `2` correspond to single-precision and double-precision formats, respectively. diff --git a/docs/documentation/getting-started.md b/docs/documentation/getting-started.md index abb4c7d55b..101a569116 100644 --- a/docs/documentation/getting-started.md +++ b/docs/documentation/getting-started.md @@ -204,6 +204,24 @@ MFC is **unit-agnostic**: the solver performs no internal unit conversions. What The only requirement is **consistency** — all inputs must use the same unit system. Note that some parameters use **transformed stored forms** rather than standard physical values (e.g., `gamma` expects \f$1/(\gamma-1)\f$, not \f$\gamma\f$ itself). See @ref sec-stored-forms for details. +## Visualizing Results + +After running post_process, visualize the output directly from the command line: + +```shell +# List available variables +./mfc.sh viz examples/2D_shockbubble/ --list-vars --step 0 + +# Render a pressure snapshot +./mfc.sh viz examples/2D_shockbubble/ --var pres --step 1000 + +# Generate a video +./mfc.sh viz examples/2D_shockbubble/ --var pres --step all --mp4 +``` + +Output images and videos are saved to the `viz/` subdirectory of the case. +For more options, see @ref visualization "Flow Visualization" or run `./mfc.sh viz -h`. + ## Helpful Tools ### Parameter Lookup diff --git a/docs/documentation/running.md b/docs/documentation/running.md index 9039ab465f..cd663634fd 100644 --- a/docs/documentation/running.md +++ b/docs/documentation/running.md @@ -73,6 +73,14 @@ using 4 cores: ./mfc.sh run examples/2D_shockbubble/case.py -t simulation post_process -n 4 ``` +- Visualizing post-processed output: + +```shell +./mfc.sh viz examples/2D_shockbubble/ --var pres --step 1000 +``` + +See @ref visualization "Flow Visualization" for the full set of visualization options. + --- ## Running on GPUs diff --git a/docs/documentation/troubleshooting.md b/docs/documentation/troubleshooting.md index e272897986..9e7dd2c468 100644 --- a/docs/documentation/troubleshooting.md +++ b/docs/documentation/troubleshooting.md @@ -37,6 +37,7 @@ This guide covers debugging tools, common issues, and troubleshooting workflows ./mfc.sh run case.py -v # Run with verbose output ./mfc.sh test --only # Run a specific test ./mfc.sh clean # Clean and start fresh +./mfc.sh viz case_dir/ --list-vars --step 0 # Inspect post-processed data ``` --- @@ -457,6 +458,47 @@ Common issues: --- +## Visualization Issues + +### "No 'binary/' or 'silo_hdf5/' directory found" + +**Cause:** Post-processing has not been run, or the case directory path is wrong. + +**Fix:** +1. Run post_process first: + ```bash + ./mfc.sh run case.py -t post_process + ``` +2. Verify the path points to the case directory (containing `binary/` or `silo_hdf5/`) + +### "Variable 'X' not found" + +**Cause:** The requested variable was not written during post-processing. + +**Fix:** +1. List available variables: + ```bash + ./mfc.sh viz case_dir/ --list-vars --step 0 + ``` +2. Ensure your case file enables the desired output (e.g., ``prim_vars_wrt = 'T'``, ``cons_vars_wrt = 'T'``) + +### "h5py is required to read Silo-HDF5 files" + +**Cause:** The case was post-processed with `format=1` (Silo-HDF5) but `h5py` is not installed. + +**Fix:** +- Install h5py: `pip install h5py` +- Or re-run post_process with `format=2` in your case file to produce binary output + +### Visualization looks wrong or has artifacts + +**Possible causes and fixes:** +1. **Color range:** Try setting explicit `--vmin` and `--vmax` values +2. **Wrong variable:** Use `--list-vars` to check available variables +3. **3D slice position:** Adjust `--slice-axis` and `--slice-value` to view the correct plane + +--- + ## Getting Help If you can't resolve an issue: diff --git a/docs/documentation/visualization.md b/docs/documentation/visualization.md index 2f5ed308d0..6a77f47eff 100644 --- a/docs/documentation/visualization.md +++ b/docs/documentation/visualization.md @@ -2,13 +2,209 @@ # Flow visualization -A post-processed database in Silo-HDF5 format can be visualized and analyzed using Paraview and VisIt. -After the post-processing of simulation data (see section @ref running "Running"), a directory named `silo_hdf5` contains a silo-HDF5 database. -Here, `silo_hdf5/` includes a directory named `root/` that contains index files for flow field data at each saved time step. +After running `post_process` on a simulation (see @ref running "Running"), MFC produces output in either Silo-HDF5 format (`format=1`) or binary format (`format=2`). +These can be visualized using MFC's built-in CLI tool or external tools like ParaView and VisIt. -### Visualizing with Paraview +--- -Paraview is an open-source interactive parallel visualization and graphical analysis tool for viewing scientific data. +## Quick visualization with `./mfc.sh viz` + +MFC includes a built-in visualization command that renders images and videos directly from post-processed output — no external GUI tools needed. + +### Basic usage + +```bash +# Plot pressure at the last available timestep (--step defaults to 'last') +./mfc.sh viz case_dir/ --var pres + +# Plot density at all available timesteps +./mfc.sh viz case_dir/ --var rho --step all +``` + +The command auto-detects the output format (binary or Silo-HDF5) and dimensionality (1D, 2D, or 3D). +Output images are saved to `case_dir/viz/` by default. +The default colormap is `viridis`, default DPI is 150, and `--step` defaults to `last`. + +### Exploring available data + +Before plotting, you can inspect what data is available: + +```bash +# List all available timesteps +./mfc.sh viz case_dir/ --list-steps + +# List all available variables at a given timestep +./mfc.sh viz case_dir/ --list-vars --step 0 +``` + +### Timestep selection + +The `--step` argument accepts several formats: + +| Format | Example | Description | +|--------|---------|-------------| +| Single | `--step 1000` | One timestep | +| Range | `--step 0:10000:500` | Start:end:stride (inclusive) | +| Last | `--step last` | Most recent available timestep | +| All | `--step all` | Every available timestep | + +### Rendering options + +Customize the appearance of plots: + +```bash +# Custom colormap and color range +./mfc.sh viz case_dir/ --var rho --step 1000 --cmap RdBu --vmin 0.5 --vmax 2.0 + +# Higher resolution +./mfc.sh viz case_dir/ --var pres --step 500 --dpi 300 + +# Logarithmic color scale +./mfc.sh viz case_dir/ --var schlieren --step 1000 --log-scale +``` + +| Option | Description | Default | +|--------|-------------|---------| +| `--cmap` | Matplotlib colormap name | `viridis` | +| `--vmin` | Minimum color scale value | auto | +| `--vmax` | Maximum color scale value | auto | +| `--dpi` | Image resolution (dots per inch) | 150 | +| `--log-scale` | Use logarithmic color scale | off | +| `--output` | Output directory for images | `case_dir/viz/` | + +### 3D slicing + +For 3D simulations, `viz` extracts a 2D slice for plotting. +By default, it slices at the midplane along the z-axis. + +> [!NOTE] +> To limit memory use, 3D batch rendering is capped at 500 timesteps and +> `--interactive` mode at 50. Use `--step start:end:stride` to stay within +> these limits when processing many steps. + +```bash +# Default z-midplane slice +./mfc.sh viz case_dir/ --var pres --step 500 + +# Slice along the x-axis at x=0.25 +./mfc.sh viz case_dir/ --var pres --step 500 --slice-axis x --slice-value 0.25 + +# Slice by array index +./mfc.sh viz case_dir/ --var pres --step 500 --slice-axis y --slice-index 50 +``` + +### Video generation + +Generate MP4 videos from a range of timesteps: + +```bash +# Basic video (10 fps) +./mfc.sh viz case_dir/ --var pres --step 0:10000:100 --mp4 + +# Custom frame rate +./mfc.sh viz case_dir/ --var schlieren --step all --mp4 --fps 24 + +# Video with fixed color range +./mfc.sh viz case_dir/ --var rho --step 0:5000:50 --mp4 --vmin 0.1 --vmax 1.0 +``` + +Videos are saved as `case_dir/viz/.mp4`. +The color range is automatically computed from the first, middle, and last frames unless `--vmin`/`--vmax` are specified. + +### Tiled 1D rendering + +For 1D cases, omitting `--var` (or passing `--var all`) renders all variables in a single tiled figure: + +```bash +# Tiled plot of all variables at the last timestep +./mfc.sh viz case_dir/ --step last + +# Equivalent explicit form +./mfc.sh viz case_dir/ --var all --step last +``` + +Each variable gets its own subplot with automatic LaTeX-style axis labels. +Tiled mode is available for 1D and 2D data. For 3D data, omitting `--var` auto-selects the first variable. + +### Interactive mode + +Launch a browser-based interactive viewer with `--interactive`: + +```bash +./mfc.sh viz case_dir/ --interactive + +# Custom port +./mfc.sh viz case_dir/ --interactive --port 9000 +``` + +The interactive viewer provides a Dash web UI with: +- Variable and timestep selection +- Live plot updates +- Pan, zoom, and hover inspection + +> [!NOTE] +> Interactive mode requires the `dash` Python package. + +### Terminal UI (TUI) + +For environments without a browser — such as SSH sessions or HPC login nodes — use `--tui` to launch a live terminal UI: + +```bash +./mfc.sh viz case_dir/ --tui + +# Start with a specific variable pre-selected +./mfc.sh viz case_dir/ --var pres --tui +``` + +The TUI loads all timesteps and renders plots directly in the terminal using Unicode block characters. +It supports 1D and 2D data only (use `--interactive` for 3D). + +**Keyboard shortcuts:** + +| Key | Action | +|-----|--------| +| `.` / `→` | Next timestep | +| `,` / `←` | Previous timestep | +| `Space` | Toggle autoplay | +| `l` | Toggle logarithmic scale | +| `f` | Freeze / unfreeze color range | +| `↑` / `↓` | Select variable (in sidebar) | +| `c` | Cycle colormap | +| `q` | Quit | + +> [!NOTE] +> The TUI requires the `textual` and `textual-plotext` Python packages, which are part of the optional `[viz]` extras and are auto-installed on the first `./mfc.sh viz --tui` run. + +### Plot styling + +Axis labels use LaTeX-style math notation — for example, `pres` is labeled as \f$p\f$, `vel1` as \f$u\f$, and `alpha1` as \f$\alpha_1\f$. +Plots use serif fonts and the Computer Modern math font for consistency with publication figures. + +### Format selection + +The output format is auto-detected from the case directory. +To override: + +```bash +./mfc.sh viz case_dir/ --var pres --step 0 --format binary +./mfc.sh viz case_dir/ --var pres --step 0 --format silo +``` + +> [!NOTE] +> Reading Silo-HDF5 files requires the `h5py` Python package. +> If it is not installed, you will see a clear error message with installation instructions. +> Alternatively, use `format=2` (binary) in your case file to produce binary output, which has no extra dependencies. + +### Complete option reference + +Run `./mfc.sh viz -h` for a full list of options. + +--- + +## Visualizing with ParaView + +ParaView is an open-source interactive parallel visualization and graphical analysis tool for viewing scientific data. +Post-processed data in Silo-HDF5 format (`format=1`) can be opened directly in ParaView. Paraview 5.11.0 has been confirmed to work with the MFC databases for some parallel environments. Nevertheless, the installation and configuration of Paraview can be environment-dependent and are left to the user. diff --git a/examples/1D_inert_shocktube/viz.py b/examples/1D_inert_shocktube/viz.py deleted file mode 100644 index 6d96f4a8bb..0000000000 --- a/examples/1D_inert_shocktube/viz.py +++ /dev/null @@ -1,80 +0,0 @@ -import mfc.viz -import os - -import subprocess -import seaborn as sns -import matplotlib.pyplot as plt -from tqdm import tqdm - -from case import sol_L as sol - -case = mfc.viz.Case(".") - -os.makedirs("viz", exist_ok=True) - -# sns.set_theme(style=mfc.viz.generate_cpg_style()) - -Y_VARS = ["H2", "O2", "H2O", "N2"] - -variables = [ - ("rho", "prim.1"), - ("u_x", "prim.2"), - ("p", "prim.3"), - ("E", "cons.3"), - *[(f"Y_{name}", f"prim.{5 + sol.species_index(name)}") for name in Y_VARS], - ("T", "prim.15"), -] - -for variable in tqdm(variables, desc="Loading Variables"): - case.load_variable(*variable) - -for step in tqdm(case.get_timesteps(), desc="Rendering Frames"): - fig, axes = plt.subplots(2, 3, figsize=(16, 9)) - - def pad_ylim(ylim, pad=0.1): - return ( - ylim[0] - pad * (ylim[1] - ylim[0]), - ylim[1] + pad * (ylim[1] - ylim[0]), - ) - - case.plot_step(step, "rho", ax=axes[0, 0]) - axes[0, 0].set_ylim(*pad_ylim(case.get_minmax_time("rho"))) - axes[0, 0].set_ylabel("$\\rho$") - case.plot_step(step, "u_x", ax=axes[0, 1]) - axes[0, 1].set_ylim(*pad_ylim(case.get_minmax_time("u_x"))) - axes[0, 1].set_ylabel("$u_x$") - case.plot_step(step, "p", ax=axes[1, 0]) - axes[1, 0].set_ylim(*pad_ylim(case.get_minmax_time("p"))) - axes[1, 0].set_ylabel("$p$") - for y in Y_VARS: - case.plot_step(step, f"Y_{y}", ax=axes[1, 1], label=y) - axes[1, 1].set_ylim(0, 1.1 * max(case.get_minmax_time(f"Y_{y}")[1] for y in Y_VARS)) - axes[1, 1].set_ylabel("$Y_k$") - case.plot_step(step, "T", ax=axes[1, 2]) - axes[1, 2].set_ylim(*pad_ylim(case.get_minmax_time("T"))) - axes[1, 2].set_ylabel("$T$") - case.plot_step(step, "E", ax=axes[0, 2]) - axes[0, 2].set_ylim(*pad_ylim(case.get_minmax_time("E"))) - axes[0, 2].set_ylabel("$E$") - - plt.tight_layout() - plt.savefig(f"viz/{step:06d}.png") - plt.close() - -subprocess.run( - [ - "ffmpeg", - "-y", - "-framerate", - "60", - "-pattern_type", - "glob", - "-i", - "viz/*.png", - "-c:v", - "libx264", - "-pix_fmt", - "yuv420p", - "viz.mp4", - ] -) diff --git a/examples/1D_reactive_shocktube/viz.py b/examples/1D_reactive_shocktube/viz.py deleted file mode 100644 index 63c769fca9..0000000000 --- a/examples/1D_reactive_shocktube/viz.py +++ /dev/null @@ -1,80 +0,0 @@ -import mfc.viz -import os - -import subprocess -import seaborn as sns -import matplotlib.pyplot as plt -from tqdm import tqdm - -from case import sol_L as sol - -case = mfc.viz.Case(".") - -os.makedirs("viz", exist_ok=True) - -sns.set_theme(style=mfc.viz.generate_cpg_style()) - -Y_VARS = ["H2", "O2", "H2O", "N2"] - -variables = [ - ("rho", "prim.1"), - ("u_x", "prim.2"), - ("p", "prim.3"), - ("E", "cons.3"), - *[(f"Y_{name}", f"prim.{5 + sol.species_index(name)}") for name in Y_VARS], - ("T", "prim.15"), -] - -for variable in tqdm(variables, desc="Loading Variables"): - case.load_variable(*variable) - -for step in tqdm(case.get_timesteps(), desc="Rendering Frames"): - fig, axes = plt.subplots(2, 3, figsize=(16, 9)) - - def pad_ylim(ylim, pad=0.1): - return ( - ylim[0] - pad * (ylim[1] - ylim[0]), - ylim[1] + pad * (ylim[1] - ylim[0]), - ) - - case.plot_step(step, "rho", ax=axes[0, 0]) - axes[0, 0].set_ylim(*pad_ylim(case.get_minmax_time("rho"))) - axes[0, 0].set_ylabel("$\\rho$") - case.plot_step(step, "u_x", ax=axes[0, 1]) - axes[0, 1].set_ylim(*pad_ylim(case.get_minmax_time("u_x"))) - axes[0, 1].set_ylabel("$u_x$") - case.plot_step(step, "p", ax=axes[1, 0]) - axes[1, 0].set_ylim(*pad_ylim(case.get_minmax_time("p"))) - axes[1, 0].set_ylabel("$p$") - for y in Y_VARS: - case.plot_step(step, f"Y_{y}", ax=axes[1, 1], label=y) - axes[1, 1].set_ylim(0, 1.1 * max(case.get_minmax_time(f"Y_{y}")[1] for y in Y_VARS)) - axes[1, 1].set_ylabel("$Y_k$") - case.plot_step(step, "T", ax=axes[1, 2]) - axes[1, 2].set_ylim(*pad_ylim(case.get_minmax_time("T"))) - axes[1, 2].set_ylabel("$T$") - case.plot_step(step, "E", ax=axes[0, 2]) - axes[0, 2].set_ylim(*pad_ylim(case.get_minmax_time("E"))) - axes[0, 2].set_ylabel("$E$") - - plt.tight_layout() - plt.savefig(f"viz/{step:06d}.png") - plt.close() - -subprocess.run( - [ - "ffmpeg", - "-y", - "-framerate", - "60", - "-pattern_type", - "glob", - "-i", - "viz/*.png", - "-c:v", - "libx264", - "-pix_fmt", - "yuv420p", - "viz.mp4", - ] -) diff --git a/examples/nD_perfect_reactor/README.md b/examples/nD_perfect_reactor/README.md index 00c2cec69a..869a411f5a 100644 --- a/examples/nD_perfect_reactor/README.md +++ b/examples/nD_perfect_reactor/README.md @@ -3,12 +3,4 @@ Reference: > G. B. Skinner and G. H. Ringrose, “Ignition Delays of a Hydrogen—Oxygen—Argon Mixture at Relatively Low Temperatures”, J. Chem. Phys., vol. 42, no. 6, pp. 2190–2192, Mar. 1965. Accessed: Oct. 13, 2024. -```bash -$ python3 analyze.py -Induction Times ([OH] >= 1e-6): - + Skinner et al.: 5.200e-05 s - + Cantera: 5.130e-05 s - + (Che)MFC: 5.130e-05 s -``` - diff --git a/examples/nD_perfect_reactor/analyze.py b/examples/nD_perfect_reactor/analyze.py deleted file mode 100644 index b2bd4e1fa8..0000000000 --- a/examples/nD_perfect_reactor/analyze.py +++ /dev/null @@ -1,107 +0,0 @@ -import cantera as ct -import seaborn as sns -from tqdm import tqdm -import matplotlib.pyplot as plt - -import mfc.viz -from case import dt, Tend, SAVE_COUNT, sol - -case = mfc.viz.Case(".", dt) - -sns.set_theme(style=mfc.viz.generate_cpg_style()) - -Y_MAJORS = set(["H", "O", "OH", "HO2"]) -Y_MINORS = set(["H2O", "H2O2"]) -Y_VARS = Y_MAJORS | Y_MINORS - -for name in tqdm(Y_VARS, desc="Loading Variables"): - case.load_variable(f"Y_{name}", f"prim.{5 + sol.species_index(name)}") -case.load_variable("rho", "prim.1") - -fig, axes = plt.subplots(1, 2, figsize=(12, 6)) - -mfc_plots = [[], []] -for y in Y_MAJORS: - mfc_plots[0].append(case.plot_time(f"Y_{y}", ax=axes[0], label=f"${y}$")) -for y in Y_MINORS: - mfc_plots[1].append(case.plot_time(f"Y_{y}", ax=axes[1], label=f"${y}$")) - -time_save = Tend / SAVE_COUNT - -oh_idx = sol.species_index("OH") - - -def generate_ct_saves() -> tuple: - reactor = ct.IdealGasReactor(sol) - reactor_network = ct.ReactorNet([reactor]) - - ct_time = 0.0 - ct_ts, ct_Ys, ct_rhos = [0.0], [reactor.thermo.Y], [reactor.thermo.density] - - while ct_time < Tend: - reactor_network.advance(ct_time + time_save) - ct_time += time_save - ct_ts.append(ct_time) - ct_Ys.append(reactor.thermo.Y) - ct_rhos.append(reactor.thermo.density) - - return ct_ts, ct_Ys, ct_rhos - - -ct_ts, ct_Ys, ct_rhos = generate_ct_saves() -for y in Y_VARS: - sns.lineplot( - x=ct_ts, - y=[yt[sol.species_index(y)] for yt in ct_Ys], - linestyle=":", - ax=axes[0 if y in Y_MAJORS else 1], - color="white", - alpha=0.5, - label=f"{y} (Cantera)", - ) - - -def find_induction_time(ts: list, Ys: list, rhos: list) -> float: - for t, y, rho in zip(ts, Ys, rhos): - if (y * rho / sol.molecular_weights[oh_idx]) >= 1e-6: - return t - - return None - - -skinner_induction_time = 0.052e-3 -ct_induction_time = find_induction_time(ct_ts, [y[oh_idx] for y in ct_Ys], [rho for rho in ct_rhos]) -mfc_induction_time = find_induction_time( - sorted(case.get_timestamps()), - [case.get_data()[step]["Y_OH"][0] for step in sorted(case.get_timesteps())], - [case.get_data()[step]["rho"][0] for step in sorted(case.get_timesteps())], -) - -print("Induction Times ([OH] >= 1e-6):") -print(f" + Skinner et al.: {skinner_induction_time:.3e} s") -print(f" + Cantera: {ct_induction_time:.3e} s") -print(f" + (Che)MFC: {mfc_induction_time:.3e} s") - -axes[0].add_artist( - axes[0].legend( - [ - axes[0].axvline(x=skinner_induction_time, color="r", linestyle="-"), - axes[0].axvline(x=mfc_induction_time, color="b", linestyle="-."), - axes[0].axvline(x=ct_induction_time, color="g", linestyle=":"), - ], - ["Skinner et al.", "(Che)MFC", "Cantera"], - title="Induction Times", - loc="lower right", - ) -) - -for i in range(2): - axes[i].legend(title="Species", ncol=2) - axes[i].set_ylabel("$Y_k$") - axes[i].set_xscale("log") - axes[i].set_yscale("log") - axes[i].set_xlabel("Time") - -plt.tight_layout() -plt.savefig(f"plots.png", dpi=300) -plt.close() diff --git a/examples/nD_perfect_reactor/export.py b/examples/nD_perfect_reactor/export.py deleted file mode 100644 index 112622a32f..0000000000 --- a/examples/nD_perfect_reactor/export.py +++ /dev/null @@ -1,76 +0,0 @@ -import csv -import cantera as ct -from tqdm import tqdm - -import mfc.viz -from case import dt, NS, Tend, SAVE_COUNT, sol - -case = mfc.viz.Case(".", dt) - -for name in tqdm(sol.species_names, desc="Loading Variables"): - case.load_variable(f"Y_{name}", f"prim.{5 + sol.species_index(name)}") -case.load_variable("rho", "prim.1") - -time_save = Tend / SAVE_COUNT - -oh_idx = sol.species_index("OH") - - -def generate_ct_saves() -> tuple: - reactor = ct.IdealGasReactor(sol) - reactor_network = ct.ReactorNet([reactor]) - - ct_time = 0.0 - ct_ts, ct_Ys, ct_rhos = [0.0], [reactor.thermo.Y], [reactor.thermo.density] - - while ct_time < Tend: - reactor_network.advance(ct_time + time_save) - ct_time += time_save - ct_ts.append(ct_time) - ct_Ys.append(reactor.thermo.Y) - ct_rhos.append(reactor.thermo.density) - - return ct_ts, ct_Ys, ct_rhos - - -ct_ts, ct_Ys, ct_rhos = generate_ct_saves() - -with open("mfc.csv", "w") as f: - writer = csv.writer(f) - keys = ["t"] + list(set(case.get_data()[0].keys()) - set(["x"])) - writer.writerow(keys) - for i, t_step in enumerate(sorted(case.get_timesteps())): - t = t_step * dt - row = [t] + [case.get_data()[t_step][key][0] for key in keys[1:]] - writer.writerow(row) - -with open("cantera.csv", "w") as f: - writer = csv.writer(f) - keys = ["t"] + [f"Y_{_}" for _ in list(sol.species_names)] + ["rho"] - writer.writerow(keys) - for step in range(len(ct_ts)): - row = [ct_ts[step]] + [ct_Ys[step][i] for i in range(len(sol.species_names))] + [ct_rhos[step]] - print([ct_ts[step]], row) - writer.writerow(row) - - -def find_induction_time(ts: list, Ys: list, rhos: list) -> float: - for t, y, rho in zip(ts, Ys, rhos): - if (y * rho / sol.molecular_weights[oh_idx]) >= 1e-6: - return t - - return None - - -skinner_induction_time = 0.052e-3 -ct_induction_time = find_induction_time(ct_ts, [y[oh_idx] for y in ct_Ys], [rho for rho in ct_rhos]) -mfc_induction_time = find_induction_time( - sorted(case.get_timestamps()), - [case.get_data()[step]["Y_OH"][0] for step in sorted(case.get_timesteps())], - [case.get_data()[step]["rho"][0] for step in sorted(case.get_timesteps())], -) - -print("Induction Times ([OH] >= 1e-6):") -print(f" + Skinner et al.: {skinner_induction_time} s") -print(f" + Cantera: {ct_induction_time} s") -print(f" + (Che)MFC: {mfc_induction_time} s") diff --git a/toolchain/bootstrap/lint.sh b/toolchain/bootstrap/lint.sh index 2a28bd18ef..1fb641d1c9 100644 --- a/toolchain/bootstrap/lint.sh +++ b/toolchain/bootstrap/lint.sh @@ -12,9 +12,27 @@ for arg in "$@"; do esac done +# Install viz optional deps only if not already present. On air-gapped systems or +# networks where PyPI is unreachable the install may fail; in that case we skip the +# viz-specific lint and tests rather than aborting the entire precheck. +VIZ_LINT=true +if ! python3 -c "import matplotlib, dash, textual, imageio, h5py" 2>/dev/null; then + log "(venv) Installing$MAGENTA viz$COLOR_RESET optional dependencies for linting..." + if ! { uv pip install -q "$(pwd)/toolchain[viz]" 2>/dev/null \ + || python3 -m pip install -q "$(pwd)/toolchain[viz]" 2>/dev/null; }; then + log "${YELLOW}Warning:${COLOR_RESET} viz optional dependencies could not be installed (no network?). Skipping viz lint/tests." + VIZ_LINT=false + fi +fi + log "(venv) Running$MAGENTA pylint$COLOR_RESET on$MAGENTA MFC$COLOR_RESET's $MAGENTA""toolchain$COLOR_RESET." -pylint -d R1722,W0718,C0301,C0116,C0115,C0114,C0410,W0622,W0640,C0103,W1309,C0411,W1514,R0401,W0511,C0321,C3001,R0801,R0911,R0912 "$(pwd)/toolchain/" +# Exclude the viz subpackage from pylint when its optional deps are unavailable, +# since pylint needs to import matplotlib/dash/etc. to analyse those modules. +PYLINT_VIZ_OPT="" +[ "$VIZ_LINT" = false ] && PYLINT_VIZ_OPT="--ignore-paths=.*/mfc/viz/.*" +# shellcheck disable=SC2086 +pylint -d R1722,W0718,C0301,C0116,C0115,C0114,C0410,W0622,W0640,C0103,W1309,C0411,W1514,R0401,W0511,C0321,C3001,R0801,R0911,R0912 $PYLINT_VIZ_OPT "$(pwd)/toolchain/" log "(venv) Running$MAGENTA pylint$COLOR_RESET on$MAGENTA MFC$COLOR_RESET's $MAGENTA""examples$COLOR_RESET." @@ -32,6 +50,9 @@ if [ "$RUN_TESTS" = true ]; then cd "$(pwd)/toolchain" python3 -m unittest mfc.params_tests.test_registry mfc.params_tests.test_definitions mfc.params_tests.test_validate mfc.params_tests.test_integration -v python3 -m unittest mfc.cli.test_cli -v + if [ "$VIZ_LINT" = true ]; then + python3 -m unittest mfc.viz.test_viz -v + fi cd - > /dev/null fi diff --git a/toolchain/main.py b/toolchain/main.py index 568db1db67..d024fb46d9 100644 --- a/toolchain/main.py +++ b/toolchain/main.py @@ -125,6 +125,8 @@ def __print_greeting(): def __checks(): + if ARG("command") == "viz": + return if not does_command_exist("cmake"): raise MFCException("CMake is required to build MFC but couldn't be located on your system. Please ensure it installed and discoverable (e.g in your system's $PATH).") @@ -175,6 +177,9 @@ def __run(): # pylint: disable=too-many-branches elif cmd == "generate": from mfc import generate # pylint: disable=import-outside-toplevel generate.generate() + elif cmd == "viz": + from mfc.viz import viz # pylint: disable=import-outside-toplevel + viz.viz() elif cmd == "params": from mfc import params_cmd # pylint: disable=import-outside-toplevel params_cmd.params() diff --git a/toolchain/mfc/args.py b/toolchain/mfc/args.py index a05fdd5d72..80e6e0888c 100644 --- a/toolchain/mfc/args.py +++ b/toolchain/mfc/args.py @@ -109,7 +109,7 @@ def custom_error(message): # Add default arguments of other subparsers # This ensures all argument keys exist even for commands that don't define them # Only process subparsers that have common arguments we need - relevant_subparsers = ["run", "test", "build", "clean", "count", "count_diff", "validate"] + relevant_subparsers = ["run", "test", "build", "clean", "count", "count_diff", "validate", "viz"] for name in relevant_subparsers: if args["command"] == name: continue @@ -120,8 +120,8 @@ def custom_error(message): # Parse with dummy input to get defaults (suppress errors for required positionals) try: # Commands with required positional input need a dummy value - if name in ["run", "validate"]: - vals, _ = subparser.parse_known_args(["dummy_input.py"]) + if name in ["run", "validate", "viz"]: + vals, _ = subparser.parse_known_args(["dummy_dir/"]) elif name == "build": vals, _ = subparser.parse_known_args([]) else: diff --git a/toolchain/mfc/cli/commands.py b/toolchain/mfc/cli/commands.py index 8ad8c4bd07..2a41e45aa3 100644 --- a/toolchain/mfc/cli/commands.py +++ b/toolchain/mfc/cli/commands.py @@ -861,6 +861,260 @@ include_common=["targets", "mfc_config", "jobs", "verbose", "debug_log"], ) +VIZ_COMMAND = Command( + name="viz", + help="Visualize post-processed MFC output.", + description=( + "Render post-processed MFC output as PNG images or MP4 video, or explore " + "interactively. Supports 1D line plots, 2D colormaps, 3D midplane slices, " + "and tiled all-variable views. PNG files are saved to case_dir/viz/ by default.\n\n" + "Output modes:\n" + " (default) Save PNG image(s) to case_dir/viz/\n" + " --mp4 Encode frames into an MP4 video\n" + " --interactive Launch a Dash web UI in your browser\n" + " --tui Launch a terminal UI (works over SSH, no browser needed)\n\n" + "Variable selection:\n" + " --var NAME Plot a single variable\n" + " (omit --var) 1D/2D: tiled layout of all variables; 3D: first variable\n\n" + "Quick-start workflow:\n" + " 1. ./mfc.sh viz case_dir/ --list-steps\n" + " 2. ./mfc.sh viz case_dir/ --list-vars --step 0\n" + " 3. ./mfc.sh viz case_dir/ --var pres --step 1000" + ), + positionals=[ + Positional( + name="input", + help="Path to the case directory containing binary/ or silo_hdf5/ output.", + completion=Completion(type=CompletionType.DIRECTORIES), + ), + ], + arguments=[ + Argument( + name="var", + help=( + "Variable to visualize (e.g. pres, rho, vel1, schlieren). " + "Omit (or pass 'all') for a tiled layout of all variables " + "(1D and 2D data) or the first variable (3D data). " + "Use --list-vars to see available names." + ), + type=str, + default=None, + metavar="VAR", + ), + Argument( + name="step", + help=( + "Timestep(s) to render. Formats: a single integer (e.g. 1000), " + "a range start:end:stride (e.g. 0:5000:500), " + "a comma list (e.g. 0,100,200), " + "an ellipsis list (e.g. 0,100,...,1000), " + "'last' (default — renders the final step only), or 'all'. " + "Use --list-steps to see available timesteps." + ), + type=str, + default='last', + metavar="STEP", + ), + Argument( + name="format", + short="f", + help="Input data format: binary or silo (auto-detected from directory structure if omitted).", + type=str, + default=None, + choices=["binary", "silo"], + completion=Completion(type=CompletionType.CHOICES, choices=["binary", "silo"]), + ), + Argument( + name="output", + short="o", + help="Directory for saved PNG images or MP4 video (default: case_dir/viz/).", + type=str, + default=None, + metavar="DIR", + completion=Completion(type=CompletionType.DIRECTORIES), + ), + Argument( + name="cmap", + help="Matplotlib colormap name (default: viridis).", + type=str, + default='viridis', + metavar="CMAP", + completion=Completion(type=CompletionType.CHOICES, choices=[ + # Perceptually uniform sequential + "viridis", "plasma", "inferno", "magma", "cividis", + # Diverging + "RdBu", "RdYlBu", "RdYlGn", "RdGy", "coolwarm", "bwr", "seismic", + "PiYG", "PRGn", "BrBG", "PuOr", "Spectral", + # Sequential + "Blues", "Greens", "Oranges", "Reds", "Purples", "Greys", + "YlOrRd", "YlOrBr", "YlGn", "YlGnBu", "GnBu", "BuGn", + "BuPu", "PuBu", "PuBuGn", "PuRd", "RdPu", + # Sequential 2 + "hot", "afmhot", "gist_heat", "copper", + "bone", "gray", "pink", "spring", "summer", "autumn", "winter", "cool", + "binary", "gist_yarg", "gist_gray", + # Cyclic + "twilight", "twilight_shifted", "hsv", + # Qualitative + "tab10", "tab20", "tab20b", "tab20c", + "Set1", "Set2", "Set3", "Paired", "Accent", "Dark2", "Pastel1", "Pastel2", + # Miscellaneous + "turbo", "jet", "rainbow", "nipy_spectral", "gist_ncar", + "gist_rainbow", "gist_stern", "gist_earth", "ocean", "terrain", + "gnuplot", "gnuplot2", "CMRmap", "cubehelix", "brg", "flag", "prism", + "Wistia", + ]), + ), + Argument( + name="vmin", + help="Minimum value for color scale.", + type=float, + default=None, + metavar="VMIN", + ), + Argument( + name="vmax", + help="Maximum value for color scale.", + type=float, + default=None, + metavar="VMAX", + ), + Argument( + name="dpi", + help="Image resolution in DPI (default: 150).", + type=int, + default=150, + metavar="DPI", + ), + Argument( + name="slice-axis", + help="Axis for 3D slice: x, y, or z (default: z).", + type=str, + default='z', + choices=["x", "y", "z"], + dest="slice_axis", + completion=Completion(type=CompletionType.CHOICES, choices=["x", "y", "z"]), + ), + Argument( + name="slice-value", + help="Coordinate value at which to take the 3D slice.", + type=float, + default=None, + dest="slice_value", + metavar="VAL", + ), + Argument( + name="slice-index", + help="Array index at which to take the 3D slice.", + type=int, + default=None, + dest="slice_index", + metavar="IDX", + ), + Argument( + name="mp4", + help="Encode all rendered frames into an MP4 video (requires --step with multiple timesteps).", + action=ArgAction.STORE_TRUE, + default=False, + ), + Argument( + name="fps", + help="Frames per second for MP4 output (default: 10).", + type=int, + default=10, + metavar="FPS", + ), + Argument( + name="list-vars", + help="Print the variable names available at the given timestep and exit.", + action=ArgAction.STORE_TRUE, + default=False, + dest="list_vars", + ), + Argument( + name="list-steps", + help="Print all available timesteps and exit.", + action=ArgAction.STORE_TRUE, + default=False, + dest="list_steps", + ), + Argument( + name="log-scale", + help="Use a logarithmic color/y scale (skips non-positive values).", + action=ArgAction.STORE_TRUE, + default=False, + dest="log_scale", + ), + Argument( + name="interactive", + short="i", + help=( + "Launch an interactive Dash web UI in your browser. " + "Loads all timesteps (or the set given by --step) and lets you " + "scrub through them and switch variables live." + ), + action=ArgAction.STORE_TRUE, + default=False, + ), + Argument( + name="port", + help="Port for the interactive web server (default: 8050).", + type=int, + default=8050, + metavar="PORT", + ), + Argument( + name="host", + help="Host/bind address for the interactive web server (default: 127.0.0.1).", + default="127.0.0.1", + metavar="HOST", + ), + Argument( + name="tui", + help=( + "Launch an interactive terminal UI (1D/2D only). " + "Works over SSH with no browser required. " + "Use arrow keys to step through timesteps." + ), + action=ArgAction.STORE_TRUE, + default=False, + ), + ], + examples=[ + Example("./mfc.sh viz case_dir/ --list-steps", "Discover available timesteps"), + Example("./mfc.sh viz case_dir/ --list-vars --step 0", "Discover available variables at step 0"), + Example("./mfc.sh viz case_dir/ --var pres --step 1000", "Save pressure PNG at step 1000 → case_dir/viz/"), + Example("./mfc.sh viz case_dir/ --step 1000", "Save tiled PNG of all variables (1D/2D) at step 1000"), + Example("./mfc.sh viz case_dir/ --var schlieren --step 0:10000:500 --mp4", "Encode schlieren MP4 from range"), + Example("./mfc.sh viz case_dir/ --step 0,100,200,...,1000", "Render all steps 0–1000 (stride inferred)"), + Example("./mfc.sh viz case_dir/ --var pres --step 500 --slice-axis x", "3D: x-plane slice of pressure"), + Example("./mfc.sh viz case_dir/ --var pres --interactive", "Browser UI — scrub timesteps and switch vars"), + Example("./mfc.sh viz case_dir/ --var pres --tui", "Terminal UI over SSH (1D/2D, no browser)"), + ], + key_options=[ + ("-- Discovery --", ""), + ("--list-steps", "Print available timesteps and exit"), + ("--list-vars", "Print available variable names and exit"), + ("-- Variable / step selection --", ""), + ("--var NAME", "Variable to plot (omit for tiled all-vars layout)"), + ("--step STEP", "last (default), int, start:stop:stride, list, or 'all'"), + ("-- Output modes --", ""), + ("(default)", "Save PNG to case_dir/viz/; use -o DIR to change"), + ("--mp4", "Encode frames into an MP4 video"), + ("--interactive / -i", "Dash web UI in browser (supports 1D/2D/3D)"), + ("--tui", "Terminal UI over SSH — no browser needed (1D/2D)"), + ("-- Appearance --", ""), + ("--cmap NAME", "Matplotlib colormap (default: viridis)"), + ("--vmin / --vmax", "Fix color-scale limits"), + ("--log-scale", "Logarithmic color/y axis"), + ("--dpi N", "Image resolution (default: 150)"), + ("-- 3D options --", ""), + ("--slice-axis x|y|z", "Plane to slice (default: z midplane)"), + ("--slice-value VAL", "Slice at coordinate value"), + ("--slice-index IDX", "Slice at array index"), + ], +) + PARAMS_COMMAND = Command( name="params", help="Search and explore MFC case parameters.", @@ -1002,6 +1256,7 @@ CLEAN_COMMAND, VALIDATE_COMMAND, NEW_COMMAND, + VIZ_COMMAND, PARAMS_COMMAND, PACKER_COMMAND, COMPLETION_COMMAND, diff --git a/toolchain/mfc/cli/docs_gen.py b/toolchain/mfc/cli/docs_gen.py index a91213a7be..989ea6a612 100644 --- a/toolchain/mfc/cli/docs_gen.py +++ b/toolchain/mfc/cli/docs_gen.py @@ -243,7 +243,7 @@ def generate_cli_reference(schema: CLISchema) -> str: # Command categories core_commands = ["build", "run", "test", "clean", "validate"] - utility_commands = ["new", "params", "packer", "completion", "generate", "help"] + utility_commands = ["new", "viz", "params", "packer", "completion", "generate", "help"] dev_commands = ["lint", "format", "spelling", "precheck", "count", "count_diff"] ci_commands = ["bench", "bench_diff"] other_commands = ["load", "interactive"] diff --git a/toolchain/mfc/viz.py b/toolchain/mfc/viz.py deleted file mode 100644 index 1c97862db3..0000000000 --- a/toolchain/mfc/viz.py +++ /dev/null @@ -1,150 +0,0 @@ -import os -import glob -import typing - -import pandas as pd -import seaborn as sns - - -def generate_cpg_style() -> dict: - BG_COLOR = '#1a1a1a' - TX_COLOR = '#FFFFFF' - - return { - 'axes.facecolor': '#121212', - 'axes.edgecolor': BG_COLOR, - 'axes.labelcolor': TX_COLOR, - 'text.color': TX_COLOR, - 'xtick.color': TX_COLOR, - 'ytick.color': TX_COLOR, - 'grid.color': BG_COLOR, - 'figure.facecolor': BG_COLOR, - 'figure.edgecolor': BG_COLOR, - 'savefig.facecolor': BG_COLOR, - 'savefig.edgecolor': BG_COLOR, - } - - -# pylint: disable=too-many-instance-attributes -class Case: - def __init__(self, dirpath: str, dt = None, parallel_io: bool = False): - assert not parallel_io, "Parallel I/O is not supported yet." - - self._dirpath = dirpath - self._data = {} - self._procs = set() - self._timesteps = set() - self._timestamps = set() - self._ndims = 0 - self._axes = [] - self._coords = [set(), set(), set()] - self._dt = dt - - self._minmax_time = {} - self._minmax_step = {} - - for f in glob.glob(os.path.join(self._dirpath, 'D', f'cons.1.*.*.dat')): - self._procs.add(int(f.split('.')[-3])) - step = int(f.split('.')[-2]) - self._timesteps.add(step) - self._timestamps.add(step * (self._dt or 1)) - - df_t0_p0 = self._read_csv('cons.1', 0, 0) - self._ndims = len(df_t0_p0.columns) - 1 - for dim in range(self._ndims): - self._coords[dim] = set(df_t0_p0.iloc[:, dim]) - self._axes.append(['x', 'y', 'z'][dim]) - - for t_step in self._timesteps: - df = pd.DataFrame() - for proc in self._procs: - df = pd.concat([ - df, - self._read_csv( - 'cons.1', proc, t_step, - names=self._axes, usecols=self._axes - ) - ]) - - self._data[t_step] = df - self._minmax_step[t_step] = {} - - for axis in self._axes: - self._compute_minmax(axis) - - def _read_csv(self, path: str, proc: int, t_step: int, **kwargs) -> pd.DataFrame: - return pd.read_csv( - os.path.join(self._dirpath, 'D', f'{path}.{proc:02d}.{t_step:06d}.dat'), - sep=r'\s+', header=None, **kwargs - ) - - def get_ndims(self) -> int: return self._ndims - def get_coords(self) -> set: return self._coords - def get_timesteps(self) -> set: return self._timesteps - def get_timestamps(self) -> set: return self._timestamps - def get_procs(self) -> set: return self._procs - def get_data(self) -> dict: return self._data - - def define_variable(self, name: str, func: typing.Callable): - for t_step, data in self._data.items(): - data[name] = data.apply( - lambda row: func(t_step, row[self._axes]), - axis=1 - ) - - self._compute_minmax(name) - - def load_variable(self, name: str, path: str): - for t_step in self._timesteps: - dfs = [] - for proc in self._procs: - dfs.append(self._read_csv( - path, proc, t_step, - names=[*self._axes, name] - )) - - self._data[t_step] = pd.merge(self._data[t_step], pd.concat(dfs)) - - self._compute_minmax(name) - - def _compute_minmax(self, varname: str): - lmins, lmaxs = set(), set() - for t_step in self._timesteps: - lmin, lmax = self._data[t_step][varname].min(), self._data[t_step][varname].max() - self._minmax_step[t_step][varname] = (lmin, lmax) - lmins.add(lmin); lmaxs.add(lmax) - - self._minmax_time[varname] = (min(lmins), max(lmaxs)) - - def get_minmax_time(self, varname: str) -> tuple: - return self._minmax_time[varname] - - def get_minmax_step(self, varname: str, t_step: int) -> tuple: - return self._minmax_step[t_step][varname] - - def plot_time(self, varname: str, aggregator: typing.Callable = None, **kwargs): - if aggregator is None: - aggregator = lambda x: x.mean() - - return sns.lineplot( - x=list((self._dt or 1) * t for t in self._timesteps), - y=[ - aggregator(self._data[t_step][varname]) - for t_step in self._timesteps - ], **kwargs) - - def plot_step(self, t_step: int, varname: str, axes: str = None, **kwargs): - axes = axes or self._axes - - if len(axes) == 1: - return sns.lineplot(self._data[t_step], x=axes[0], y=varname, **kwargs) - - if len(axes) == 2: - return sns.heatmap( - self._data[t_step].pivot( - index=axes[0], columns=axes[1], values=varname - ), - **kwargs - ) - - assert False, "3D plotting is not supported yet." diff --git a/toolchain/mfc/viz/__init__.py b/toolchain/mfc/viz/__init__.py new file mode 100644 index 0000000000..e69de29bb2 diff --git a/toolchain/mfc/viz/_step_cache.py b/toolchain/mfc/viz/_step_cache.py new file mode 100644 index 0000000000..ec4780565a --- /dev/null +++ b/toolchain/mfc/viz/_step_cache.py @@ -0,0 +1,57 @@ +"""Bounded FIFO step cache shared by tui.py and interactive.py. + +Keeps up to CACHE_MAX assembled timesteps in memory, evicting the oldest +entry when the cap is reached. The module-level state is intentional: +both the TUI and the interactive server are single-instance; a per-session +cache avoids redundant disk reads while bounding peak memory usage. + +Thread safety: Dash runs Flask with threading enabled. All mutations are +protected by _lock so concurrent callbacks from multiple browser tabs cannot +corrupt _cache_order or exceed CACHE_MAX. +""" + +import threading +from typing import Callable + +CACHE_MAX: int = 50 +_cache: dict = {} +_cache_order: list = [] +_lock = threading.Lock() + + +def load(step: int, read_func: Callable) -> object: + """Return cached data for *step*, calling *read_func* on a miss. + + read_func is called *before* eviction so that a failed read (e.g. a + missing or corrupt file) does not discard a valid cache entry. + """ + with _lock: + if step in _cache: + return _cache[step] + # Read outside-the-lock would allow concurrent loads of the same + # step; keeping it inside is simpler and safe since read_func is + # plain file I/O that never calls load() recursively. + data = read_func(step) + # Evict only after a successful read. + if len(_cache) >= CACHE_MAX: + evict = _cache_order.pop(0) + _cache.pop(evict, None) + _cache[step] = data + _cache_order.append(step) + return data + + +def seed(step: int, data: object) -> None: + """Clear the cache and pre-populate it with already-loaded data.""" + with _lock: + _cache.clear() + _cache_order.clear() + _cache[step] = data + _cache_order.append(step) + + +def clear() -> None: + """Reset the cache to empty (useful for test teardown).""" + with _lock: + _cache.clear() + _cache_order.clear() diff --git a/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/0.dat b/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/0.dat new file mode 100644 index 0000000000..eb7725ff71 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/0.dat differ diff --git a/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/1.dat b/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/1.dat new file mode 100644 index 0000000000..9c77b3a488 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/1.dat differ diff --git a/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/2.dat b/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/2.dat new file mode 100644 index 0000000000..c1880e9625 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_binary/binary/p0/2.dat differ diff --git a/toolchain/mfc/viz/fixtures/1d_binary/binary/root/0.dat b/toolchain/mfc/viz/fixtures/1d_binary/binary/root/0.dat new file mode 100644 index 0000000000..eb7725ff71 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_binary/binary/root/0.dat differ diff --git a/toolchain/mfc/viz/fixtures/1d_binary/binary/root/1.dat b/toolchain/mfc/viz/fixtures/1d_binary/binary/root/1.dat new file mode 100644 index 0000000000..9c77b3a488 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_binary/binary/root/1.dat differ diff --git a/toolchain/mfc/viz/fixtures/1d_binary/binary/root/2.dat b/toolchain/mfc/viz/fixtures/1d_binary/binary/root/2.dat new file mode 100644 index 0000000000..c1880e9625 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_binary/binary/root/2.dat differ diff --git a/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/0.silo b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/0.silo new file mode 100644 index 0000000000..360cdeab89 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/0.silo differ diff --git a/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/1.silo b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/1.silo new file mode 100644 index 0000000000..14e09a31f9 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/1.silo differ diff --git a/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/2.silo b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/2.silo new file mode 100644 index 0000000000..2f5a4e468c Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/p0/2.silo differ diff --git a/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_0.silo b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_0.silo new file mode 100644 index 0000000000..af6de87b83 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_0.silo differ diff --git a/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_1.silo b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_1.silo new file mode 100644 index 0000000000..4a9c85e880 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_1.silo differ diff --git a/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_2.silo b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_2.silo new file mode 100644 index 0000000000..8bf7d5eb88 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/1d_silo/silo_hdf5/root/collection_2.silo differ diff --git a/toolchain/mfc/viz/fixtures/2d_binary/binary/p0/0.dat b/toolchain/mfc/viz/fixtures/2d_binary/binary/p0/0.dat new file mode 100644 index 0000000000..7bbb7fc340 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/2d_binary/binary/p0/0.dat differ diff --git a/toolchain/mfc/viz/fixtures/2d_binary/binary/p0/1.dat b/toolchain/mfc/viz/fixtures/2d_binary/binary/p0/1.dat new file mode 100644 index 0000000000..9cf99145e1 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/2d_binary/binary/p0/1.dat differ diff --git a/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/p0/0.silo b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/p0/0.silo new file mode 100644 index 0000000000..23901b7fca Binary files /dev/null and b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/p0/0.silo differ diff --git a/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/p0/1.silo b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/p0/1.silo new file mode 100644 index 0000000000..15ed802b45 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/p0/1.silo differ diff --git a/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/root/collection_0.silo b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/root/collection_0.silo new file mode 100644 index 0000000000..16abbd65ca Binary files /dev/null and b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/root/collection_0.silo differ diff --git a/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/root/collection_1.silo b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/root/collection_1.silo new file mode 100644 index 0000000000..61bbdf1396 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/2d_silo/silo_hdf5/root/collection_1.silo differ diff --git a/toolchain/mfc/viz/fixtures/3d_binary/binary/p0/0.dat b/toolchain/mfc/viz/fixtures/3d_binary/binary/p0/0.dat new file mode 100644 index 0000000000..a717f5bf70 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/3d_binary/binary/p0/0.dat differ diff --git a/toolchain/mfc/viz/fixtures/3d_binary/binary/p0/1.dat b/toolchain/mfc/viz/fixtures/3d_binary/binary/p0/1.dat new file mode 100644 index 0000000000..8269df27cc Binary files /dev/null and b/toolchain/mfc/viz/fixtures/3d_binary/binary/p0/1.dat differ diff --git a/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/p0/0.silo b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/p0/0.silo new file mode 100644 index 0000000000..d90063d6f2 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/p0/0.silo differ diff --git a/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/p0/1.silo b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/p0/1.silo new file mode 100644 index 0000000000..ec3f35f6d7 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/p0/1.silo differ diff --git a/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/root/collection_0.silo b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/root/collection_0.silo new file mode 100644 index 0000000000..65dd70f157 Binary files /dev/null and b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/root/collection_0.silo differ diff --git a/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/root/collection_1.silo b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/root/collection_1.silo new file mode 100644 index 0000000000..8366ba186f Binary files /dev/null and b/toolchain/mfc/viz/fixtures/3d_silo/silo_hdf5/root/collection_1.silo differ diff --git a/toolchain/mfc/viz/interactive.py b/toolchain/mfc/viz/interactive.py new file mode 100644 index 0000000000..6d0e0189e0 --- /dev/null +++ b/toolchain/mfc/viz/interactive.py @@ -0,0 +1,664 @@ +""" +Interactive Dash-based visualization for MFC post-processed data. + +Launched via ``./mfc.sh viz --var --step all --interactive``. +Opens a dark-themed web UI in your browser (or via SSH tunnel) with live +controls for slice position, isosurface thresholds, volume opacity, +colormap, log scale, vmin/vmax, and timestep playback. +""" +# pylint: disable=use-dict-literal + +from typing import List, Callable, Optional + +import numpy as np +import plotly.graph_objects as go +from dash import Dash, dcc, html, Input, Output, State, callback_context, no_update + +from mfc.printer import cons +from . import _step_cache + +# --------------------------------------------------------------------------- +# Colormaps available in the picker +# --------------------------------------------------------------------------- +_CMAPS = [ + "viridis", "plasma", "inferno", "magma", "cividis", + "turbo", "jet", "rainbow", "nipy_spectral", + "RdBu", "RdYlBu", "RdYlGn", "coolwarm", "bwr", "seismic", "Spectral", + "hot", "afmhot", "gist_heat", "copper", + "bone", "gray", "spring", "summer", "autumn", "winter", "cool", "pink", + "Blues", "Greens", "Oranges", "Reds", "Purples", "Greys", + "twilight", "twilight_shifted", "hsv", + "tab10", "tab20", "terrain", "ocean", "gist_earth", + "gnuplot", "gnuplot2", "CMRmap", "cubehelix", "Wistia", +] + +# --------------------------------------------------------------------------- +# Catppuccin Mocha palette +# --------------------------------------------------------------------------- +_BG = '#181825' +_SURF = '#1e1e2e' +_OVER = '#313244' +_BORDER = '#45475a' +_TEXT = '#cdd6f4' +_SUB = '#a6adc8' +_MUTED = '#6c7086' +_ACCENT = '#cba6f7' +_GREEN = '#a6e3a1' +_RED = '#f38ba8' +_BLUE = '#89b4fa' +_TEAL = '#94e2d5' +_YELLOW = '#f9e2af' + +# --------------------------------------------------------------------------- +# Server-side data cache {step -> AssembledData} (bounded to avoid OOM) +# --------------------------------------------------------------------------- +_load = _step_cache.load +_CACHE_MAX = _step_cache.CACHE_MAX + + +# --------------------------------------------------------------------------- +# Layout helpers +# --------------------------------------------------------------------------- + +def _section(title, *children): + return html.Div([ + html.Div(title, style={ + 'fontSize': '10px', 'fontWeight': 'bold', + 'textTransform': 'uppercase', 'letterSpacing': '0.08em', + 'color': _MUTED, 'borderBottom': f'1px solid {_OVER}', + 'paddingBottom': '4px', 'marginTop': '16px', 'marginBottom': '6px', + }), + *children, + ]) + + +def _lbl(text): + return html.Div(text, style={ + 'fontSize': '11px', 'color': _SUB, + 'marginBottom': '2px', 'marginTop': '6px', + }) + + +def _slider(sid, lo, hi, step, val, marks=None): # pylint: disable=too-many-arguments,too-many-positional-arguments + return dcc.Slider( + id=sid, min=lo, max=hi, step=step, value=val, + marks=marks or {}, updatemode='mouseup', + tooltip={'placement': 'bottom', 'always_visible': True}, + ) + + +def _btn(bid, label, color=_TEXT): + return html.Button(label, id=bid, n_clicks=0, style={ + 'flex': '1', 'padding': '5px 8px', 'fontSize': '12px', + 'backgroundColor': _OVER, 'color': color, + 'border': f'1px solid {_BORDER}', 'borderRadius': '4px', + 'cursor': 'pointer', 'fontFamily': 'monospace', + }) + + +def _num(sid, placeholder='auto'): + return dcc.Input( + id=sid, type='number', placeholder=placeholder, debounce=True, + style={ + 'width': '100%', 'backgroundColor': _OVER, 'color': _TEXT, + 'border': f'1px solid {_BORDER}', 'borderRadius': '4px', + 'padding': '4px 6px', 'fontSize': '12px', 'fontFamily': 'monospace', + 'boxSizing': 'border-box', + }, + ) + + +# --------------------------------------------------------------------------- +# 3D figure builder +# --------------------------------------------------------------------------- + +def _build_3d(ad, raw, varname, step, mode, cmap, # pylint: disable=too-many-arguments,too-many-positional-arguments,too-many-locals,too-many-branches,too-many-statements + log_fn, cmin, cmax, cbar_title, + slice_axis, slice_pos, + iso_min_frac, iso_max_frac, iso_n, iso_caps, + vol_opacity, vol_nsurf, vol_min_frac, vol_max_frac): + """Return (trace, title) for a 3D assembled dataset.""" + cbar = dict( + title=dict(text=cbar_title, font=dict(color=_TEXT)), + tickfont=dict(color=_TEXT), + ) + rng = cmax - cmin if cmax > cmin else 1.0 + + if mode == 'slice': + axis_coords = {'x': ad.x_cc, 'y': ad.y_cc, 'z': ad.z_cc} + coords = axis_coords[slice_axis] + coord_val = coords[0] + (coords[-1] - coords[0]) * slice_pos + idx = int(np.clip(np.argmin(np.abs(coords - coord_val)), 0, len(coords) - 1)) + actual = float(coords[idx]) + + if slice_axis == 'x': + sliced = log_fn(raw[idx, :, :]) # (ny, nz) + YY, ZZ = np.meshgrid(ad.y_cc, ad.z_cc, indexing='ij') + trace = go.Surface( + x=np.full_like(YY, actual), y=YY, z=ZZ, + surfacecolor=sliced, cmin=cmin, cmax=cmax, + colorscale=cmap, colorbar=cbar, showscale=True, + ) + elif slice_axis == 'y': + sliced = log_fn(raw[:, idx, :]) # (nx, nz) + XX, ZZ = np.meshgrid(ad.x_cc, ad.z_cc, indexing='ij') + trace = go.Surface( + x=XX, y=np.full_like(XX, actual), z=ZZ, + surfacecolor=sliced, cmin=cmin, cmax=cmax, + colorscale=cmap, colorbar=cbar, showscale=True, + ) + else: # z + sliced = log_fn(raw[:, :, idx]) # (nx, ny) + XX, YY = np.meshgrid(ad.x_cc, ad.y_cc, indexing='ij') + trace = go.Surface( + x=XX, y=YY, z=np.full_like(XX, actual), + surfacecolor=sliced, cmin=cmin, cmax=cmax, + colorscale=cmap, colorbar=cbar, showscale=True, + ) + title = f'{varname} · {slice_axis} = {actual:.4g} · step {step}' + + elif mode == 'isosurface': + X3, Y3, Z3 = np.meshgrid(ad.x_cc, ad.y_cc, ad.z_cc, indexing='ij') + vf = log_fn(raw.ravel()) + ilo = cmin + rng * iso_min_frac + ihi = cmin + rng * max(iso_max_frac, iso_min_frac + 0.01) + caps = dict(x_show=iso_caps, y_show=iso_caps, z_show=iso_caps) + trace = go.Isosurface( + x=X3.ravel(), y=Y3.ravel(), z=Z3.ravel(), value=vf, + isomin=ilo, isomax=ihi, surface_count=int(iso_n), + colorscale=cmap, cmin=cmin, cmax=cmax, + caps=caps, colorbar=cbar, + ) + title = f'{varname} · {int(iso_n)} isosurfaces · step {step}' + + else: # volume + X3, Y3, Z3 = np.meshgrid(ad.x_cc, ad.y_cc, ad.z_cc, indexing='ij') + vf = log_fn(raw.ravel()) + vlo = cmin + rng * vol_min_frac + vhi = cmin + rng * max(vol_max_frac, vol_min_frac + 0.01) + trace = go.Volume( + x=X3.ravel(), y=Y3.ravel(), z=Z3.ravel(), value=vf, + isomin=vlo, isomax=vhi, + opacity=float(vol_opacity), surface_count=int(vol_nsurf), + colorscale=cmap, colorbar=cbar, + ) + title = f'{varname} · volume · step {step}' + + return trace, title + + +# --------------------------------------------------------------------------- +# Main entry point +# --------------------------------------------------------------------------- + +def run_interactive( # pylint: disable=too-many-locals,too-many-statements,too-many-arguments,too-many-positional-arguments + varname: str, + steps: List[int], + read_func: Callable, + port: int = 8050, + host: str = '127.0.0.1', + bubble_func: Optional[Callable] = None, +): + """Launch the interactive Dash visualization server.""" + app = Dash( + __name__, + title=f'MFC viz · {varname}', + suppress_callback_exceptions=True, + ) + + # Load first step to know dimensionality and available variables + init = _load(steps[0], read_func) + ndim = init.ndim + all_varnames = sorted(init.variables.keys()) + if varname not in all_varnames: + varname = all_varnames[0] if all_varnames else varname + + step_opts = [{'label': str(s), 'value': s} for s in steps] + var_opts = [{'label': v, 'value': v} for v in all_varnames] + cmap_opts = [{'label': c, 'value': c} for c in _CMAPS] + + if ndim == 3: + mode_opts = [ + {'label': ' Slice', 'value': 'slice'}, + {'label': ' Isosurface', 'value': 'isosurface'}, + {'label': ' Volume', 'value': 'volume'}, + ] + elif ndim == 2: + mode_opts = [{'label': ' Heatmap', 'value': 'heatmap'}] + else: + mode_opts = [{'label': ' Line', 'value': 'line'}] + default_mode = mode_opts[0]['value'] + + # ------------------------------------------------------------------ + # Sidebar layout + # ------------------------------------------------------------------ + sidebar = html.Div([ + + # Header + html.Div('MFC viz', style={ + 'fontSize': '16px', 'fontWeight': 'bold', 'color': _ACCENT, + }), + html.Div( + f'{ndim}D · {len(steps)} step{"s" if len(steps) != 1 else ""}', + style={'fontSize': '11px', 'color': _MUTED}, + ), + + # ── Variable ────────────────────────────────────────────────── + _section('Variable', + dcc.Dropdown( + id='var-sel', options=var_opts, value=varname, clearable=False, + style={'fontSize': '12px', 'backgroundColor': _OVER, + 'border': f'1px solid {_BORDER}'}, + ), + ), + + # ── Timestep ────────────────────────────────────────────────── + _section('Timestep', + dcc.Dropdown( + id='step-sel', options=step_opts, value=steps[0], clearable=False, + style={'fontSize': '12px', 'backgroundColor': _OVER, + 'border': f'1px solid {_BORDER}'}, + ), + html.Div([ + _btn('play-btn', '▶ Play', _GREEN), + html.Div(style={'width': '6px'}), + _btn('stop-btn', '■ Stop', _RED), + ], style={'display': 'flex', 'marginTop': '6px'}), + _lbl('Playback speed (fps)'), + _slider('fps-sl', 0.5, 10, 0.5, 2, + marks={0.5: '0.5', 2: '2', 5: '5', 10: '10'}), + dcc.Checklist( + id='loop-chk', + options=[{'label': ' Loop', 'value': 'loop'}], value=['loop'], + style={'fontSize': '12px', 'color': _SUB, 'marginTop': '4px'}, + ), + ), + + # ── Viz mode ────────────────────────────────────────────────── + _section('Viz Mode', + dcc.RadioItems( + id='mode-sel', options=mode_opts, value=default_mode, + style={'fontSize': '12px', 'color': _TEXT}, + inputStyle={'marginRight': '6px', 'cursor': 'pointer'}, + labelStyle={'display': 'block', 'marginBottom': '5px', 'cursor': 'pointer'}, + ), + ), + + # ── Slice ───────────────────────────────────────────────────── + html.Div(id='ctrl-slice', children=[ + _section('Slice', + _lbl('Axis'), + dcc.RadioItems( + id='slice-axis', options=['x', 'y', 'z'], value='z', + inline=True, style={'fontSize': '12px', 'color': _TEXT}, + inputStyle={'marginRight': '4px'}, + labelStyle={'marginRight': '14px'}, + ), + _lbl('Position (0 = start, 1 = end)'), + _slider('slice-pos', 0.0, 1.0, 0.01, 0.5, + marks={0: '0', 0.5: '½', 1: '1'}), + ), + ]), + + # ── Isosurface ──────────────────────────────────────────────── + html.Div(id='ctrl-iso', style={'display': 'none'}, children=[ + _section('Isosurface', + _lbl('Min threshold (fraction of color range)'), + _slider('iso-min', 0.0, 1.0, 0.01, 0.2, + marks={0: '0', 0.5: '0.5', 1: '1'}), + _lbl('Max threshold (fraction of color range)'), + _slider('iso-max', 0.0, 1.0, 0.01, 0.8, + marks={0: '0', 0.5: '0.5', 1: '1'}), + _lbl('Number of isosurfaces'), + _slider('iso-n', 1, 10, 1, 3, + marks={1: '1', 3: '3', 5: '5', 10: '10'}), + dcc.Checklist( + id='iso-caps', + options=[{'label': ' Show end-caps', 'value': 'caps'}], value=[], + style={'fontSize': '12px', 'color': _SUB, 'marginTop': '6px'}, + ), + ), + ]), + + # ── Volume ──────────────────────────────────────────────────── + html.Div(id='ctrl-vol', style={'display': 'none'}, children=[ + _section('Volume', + _lbl('Opacity per shell'), + _slider('vol-opacity', 0.01, 0.5, 0.01, 0.1, + marks={0.01: '0', 0.25: '.25', 0.5: '.5'}), + _lbl('Number of shells'), + _slider('vol-nsurf', 3, 30, 1, 15, + marks={3: '3', 15: '15', 30: '30'}), + _lbl('Isomin (fraction of color range)'), + _slider('vol-min', 0.0, 1.0, 0.01, 0.0, + marks={0: '0', 0.5: '0.5', 1: '1'}), + _lbl('Isomax (fraction of color range)'), + _slider('vol-max', 0.0, 1.0, 0.01, 1.0, + marks={0: '0', 0.5: '0.5', 1: '1'}), + ), + ]), + + # ── Color ───────────────────────────────────────────────────── + _section('Color', + _lbl('Colormap'), + dcc.Dropdown( + id='cmap-sel', options=cmap_opts, value='viridis', clearable=False, + style={'fontSize': '12px', 'backgroundColor': _OVER, + 'border': f'1px solid {_BORDER}'}, + ), + dcc.Checklist( + id='log-chk', + options=[{'label': ' Log scale', 'value': 'log'}], value=[], + style={'fontSize': '12px', 'color': _SUB, 'marginTop': '6px'}, + ), + html.Div([ + html.Div([_lbl('vmin'), _num('vmin-inp')], + style={'flex': 1, 'marginRight': '6px'}), + html.Div([_lbl('vmax'), _num('vmax-inp')], + style={'flex': 1}), + ], style={'display': 'flex'}), + html.Button('↺ Auto range', id='reset-btn', n_clicks=0, style={ + 'marginTop': '8px', 'padding': '4px 8px', 'fontSize': '11px', + 'width': '100%', 'backgroundColor': _OVER, 'color': _TEAL, + 'border': f'1px solid {_BORDER}', 'borderRadius': '4px', + 'cursor': 'pointer', 'fontFamily': 'monospace', + }), + ), + + # ── Status ──────────────────────────────────────────────────── + html.Div(id='status-bar', style={ + 'marginTop': 'auto', 'paddingTop': '12px', + 'fontSize': '11px', 'color': _MUTED, + 'borderTop': f'1px solid {_OVER}', 'lineHeight': '1.7', + }), + + ], style={ + 'width': '265px', 'minWidth': '265px', + 'backgroundColor': _SURF, 'padding': '14px', + 'height': '100vh', 'overflowY': 'auto', + 'display': 'flex', 'flexDirection': 'column', + 'fontFamily': 'monospace', 'color': _TEXT, + 'boxSizing': 'border-box', + }) + + app.layout = html.Div([ + sidebar, + html.Div([ + dcc.Graph( + id='viz-graph', style={'height': '100vh'}, + config={ + 'displayModeBar': True, 'scrollZoom': True, + 'modeBarButtonsToRemove': ['select2d', 'lasso2d'], + 'toImageButtonOptions': {'format': 'png', 'scale': 2}, + }, + ), + ], style={'flex': '1', 'overflow': 'hidden', 'backgroundColor': _BG}), + + dcc.Interval(id='play-iv', interval=500, n_intervals=0, disabled=True), + dcc.Store(id='playing-st', data=False), + ], style={ + 'display': 'flex', 'height': '100vh', 'overflow': 'hidden', + 'backgroundColor': _BG, 'fontFamily': 'monospace', + }) + + # ------------------------------------------------------------------ + # Callbacks + # ------------------------------------------------------------------ + + @app.callback( + Output('play-iv', 'disabled'), + Output('play-iv', 'interval'), + Output('playing-st', 'data'), + Output('play-btn', 'children'), + Input('play-btn', 'n_clicks'), + Input('stop-btn', 'n_clicks'), + Input('fps-sl', 'value'), + State('playing-st', 'data'), + prevent_initial_call=True, + ) + def _toggle_play(_, __, fps, is_playing): # pylint: disable=unused-argument + iv = max(int(1000 / max(float(fps or 2), 0.1)), 50) + trig = (callback_context.triggered or [{}])[0].get('prop_id', '') + if 'stop-btn' in trig: + return True, iv, False, '▶ Play' + if 'play-btn' in trig: + playing = not is_playing + return not playing, iv, playing, ('⏸ Pause' if playing else '▶ Play') + return not is_playing, iv, is_playing, no_update # fps-only change + + @app.callback( + Output('step-sel', 'value'), + Input('play-iv', 'n_intervals'), + State('step-sel', 'value'), + State('loop-chk', 'value'), + prevent_initial_call=True, + ) + def _advance_step(_, current, loop_val): + try: + idx = steps.index(current) + except ValueError: + idx = 0 + nxt = idx + 1 + if nxt >= len(steps): + return steps[0] if ('loop' in (loop_val or [])) else no_update + return steps[nxt] + + @app.callback( + Output('ctrl-slice', 'style'), + Output('ctrl-iso', 'style'), + Output('ctrl-vol', 'style'), + Input('mode-sel', 'value'), + ) + def _toggle_controls(mode): + show, hide = {'display': 'block'}, {'display': 'none'} + return ( + show if mode == 'slice' else hide, + show if mode == 'isosurface' else hide, + show if mode == 'volume' else hide, + ) + + @app.callback( + Output('vmin-inp', 'value'), + Output('vmax-inp', 'value'), + Input('reset-btn', 'n_clicks'), + Input('var-sel', 'value'), + prevent_initial_call=True, + ) + def _reset_range(_reset, _var): + return None, None + + @app.callback( + Output('viz-graph', 'figure'), + Output('status-bar', 'children'), + Input('var-sel', 'value'), + Input('step-sel', 'value'), + Input('mode-sel', 'value'), + Input('slice-axis', 'value'), + Input('slice-pos', 'value'), + Input('iso-min', 'value'), + Input('iso-max', 'value'), + Input('iso-n', 'value'), + Input('iso-caps', 'value'), + Input('vol-opacity', 'value'), + Input('vol-nsurf', 'value'), + Input('vol-min', 'value'), + Input('vol-max', 'value'), + Input('cmap-sel', 'value'), + Input('log-chk', 'value'), + Input('vmin-inp', 'value'), + Input('vmax-inp', 'value'), + ) + def _update(var_sel, step, mode, # pylint: disable=too-many-arguments,too-many-positional-arguments,too-many-locals,too-many-branches,too-many-statements + slice_axis, slice_pos, + iso_min_frac, iso_max_frac, iso_n, iso_caps, + vol_opacity, vol_nsurf, vol_min_frac, vol_max_frac, + cmap, log_chk, vmin_in, vmax_in): + + selected_var = var_sel or varname + try: + ad = _load(step, read_func) + except Exception as exc: # pylint: disable=broad-except + return no_update, [html.Span(f' Error loading step {step}: {exc}', + style={'color': _RED})] + if selected_var not in ad.variables: + avail = ', '.join(sorted(ad.variables)) + return no_update, [html.Span( + f' Variable {selected_var!r} not in step {step} ' + f'(available: {avail})', style={'color': _RED})] + raw = ad.variables[selected_var] + log = bool(log_chk and 'log' in log_chk) + cmap = cmap or 'viridis' + + # Color range + if vmin_in is not None: + vmin = float(vmin_in) + else: + safe = raw[raw > 0] if log and np.any(raw > 0) else raw + vmin = float(np.nanmin(safe)) + if vmax_in is not None: + vmax = float(vmax_in) + else: + vmax = float(np.nanmax(raw)) + if vmax <= vmin: + vmax = vmin + 1e-10 + + if log: + def _tf(arr): + return np.where(arr > 0, np.log10(np.maximum(arr, 1e-300)), np.nan) + cmin = float(np.log10(max(vmin, 1e-300))) + cmax = float(np.log10(max(vmax, 1e-300))) + cbar_title = f'log\u2081\u2080({selected_var})' + else: + def _tf(arr): return arr + cmin, cmax = vmin, vmax + cbar_title = selected_var + + fig = go.Figure() + title = '' + + if ad.ndim == 3: + trace, title = _build_3d( + ad, raw, selected_var, step, mode, cmap, _tf, cmin, cmax, cbar_title, + slice_axis or 'z', float(slice_pos or 0.5), + float(iso_min_frac or 0.2), float(iso_max_frac or 0.8), + int(iso_n or 3), bool(iso_caps and 'caps' in iso_caps), + float(vol_opacity or 0.1), int(vol_nsurf or 15), + float(vol_min_frac or 0.0), float(vol_max_frac or 1.0), + ) + fig.add_trace(trace) + # Bubble overlay for 3D + if bubble_func is not None: + try: + bubbles = bubble_func(step) + if bubbles is not None and len(bubbles) > 0: + if mode == 'slice': + s_axis = slice_axis or 'z' + s_col = {'x': 0, 'y': 1, 'z': 2}[s_axis] + ax_coords = {'x': ad.x_cc, 'y': ad.y_cc, 'z': ad.z_cc}[s_axis] + s_coord = ax_coords[0] + (ax_coords[-1] - ax_coords[0]) * float(slice_pos or 0.5) + near = np.abs(bubbles[:, s_col] - s_coord) <= bubbles[:, 3] + vis = bubbles[near] if near.any() else None + else: + vis = bubbles + if vis is not None and len(vis) > 0: + fig.add_trace(go.Scatter3d( + x=vis[:, 0], y=vis[:, 1], z=vis[:, 2], + mode='markers', + marker=dict(size=4, color='white', opacity=0.6, symbol='circle'), + showlegend=False, + hovertemplate='x=%{x:.3g}
y=%{y:.3g}
z=%{z:.3g}bubble', + )) + except Exception: # pylint: disable=broad-except + pass + # Compute aspect ratio from domain extents so slices (which + # have a constant coordinate on one axis) don't collapse that axis. + dx = float(ad.x_cc[-1] - ad.x_cc[0]) if len(ad.x_cc) > 1 else 1.0 + dy = float(ad.y_cc[-1] - ad.y_cc[0]) if len(ad.y_cc) > 1 else 1.0 + dz = float(ad.z_cc[-1] - ad.z_cc[0]) if len(ad.z_cc) > 1 else 1.0 + max_d = max(dx, dy, dz, 1e-30) + fig.update_layout(scene=dict( + xaxis=dict(title='x', backgroundcolor=_SURF, + gridcolor=_OVER, color=_TEXT), + yaxis=dict(title='y', backgroundcolor=_SURF, + gridcolor=_OVER, color=_TEXT), + zaxis=dict(title='z', backgroundcolor=_SURF, + gridcolor=_OVER, color=_TEXT), + bgcolor=_BG, + aspectmode='manual', + aspectratio=dict(x=dx/max_d, y=dy/max_d, z=dz/max_d), + )) + + elif ad.ndim == 2: + cbar = dict(title=dict(text=cbar_title, font=dict(color=_TEXT)), + tickfont=dict(color=_TEXT)) + fig.add_trace(go.Heatmap( + x=ad.x_cc, y=ad.y_cc, z=_tf(raw).T, + zmin=cmin, zmax=cmax, colorscale=cmap, colorbar=cbar, + )) + fig.update_layout( + xaxis=dict(title='x', color=_TEXT, gridcolor=_OVER, scaleanchor='y'), + yaxis=dict(title='y', color=_TEXT, gridcolor=_OVER), + plot_bgcolor=_BG, + ) + # Bubble overlay for 2D + if bubble_func is not None: + try: + bubbles = bubble_func(step) + if bubbles is not None and len(bubbles) > 0: + shapes = [ + dict( + type='circle', + xref='x', yref='y', + x0=float(b[0] - b[3]), y0=float(b[1] - b[3]), + x1=float(b[0] + b[3]), y1=float(b[1] + b[3]), + line=dict(color='white', width=0.8), + fillcolor='rgba(0,0,0,0)', + ) + for b in bubbles + ] + fig.update_layout(shapes=shapes) + except Exception: # pylint: disable=broad-except + pass + title = f'{selected_var} · step {step}' + + else: # 1D + plot_y = _tf(raw) if log else raw + fig.add_trace(go.Scatter( + x=ad.x_cc, y=plot_y, mode='lines', + line=dict(color=_ACCENT, width=2), name=selected_var, + )) + fig.update_layout( + xaxis=dict(title='x', color=_TEXT, gridcolor=_OVER), + yaxis=dict(title=cbar_title, color=_TEXT, gridcolor=_OVER, + range=[cmin, cmax] if (vmin_in or vmax_in) else None), + plot_bgcolor=_BG, + ) + title = f'{selected_var} · step {step}' + + fig.update_layout( + title=dict(text=title, font=dict(color=_TEXT, size=13, family='monospace')), + paper_bgcolor=_BG, + font=dict(color=_TEXT, family='monospace'), + margin=dict(l=0, r=0, t=36, b=0), + uirevision=mode, # preserve camera angle within a mode + ) + + dmin, dmax = float(np.nanmin(raw)), float(np.nanmax(raw)) + status = html.Div([ + html.Span(f'step {step}', style={'color': _YELLOW}), + html.Span(f' · shape {raw.shape}', style={'color': _MUTED}), + html.Br(), + html.Span('min ', style={'color': _MUTED}), + html.Span(f'{dmin:.4g}', style={'color': _BLUE}), + html.Span(' max ', style={'color': _MUTED}), + html.Span(f'{dmax:.4g}', style={'color': _RED}), + ]) + return fig, status + + # ------------------------------------------------------------------ + cons.print(f'\n[bold green]Interactive viz server:[/bold green] ' + f'[bold]http://{host}:{port}[/bold]') + if host in ('127.0.0.1', 'localhost'): + cons.print(f'[dim]SSH tunnel: ssh -L {port}:localhost:{port} [/dim]') + cons.print('[dim]Ctrl+C to stop.[/dim]\n') + app.run(debug=False, port=port, host=host) diff --git a/toolchain/mfc/viz/reader.py b/toolchain/mfc/viz/reader.py new file mode 100644 index 0000000000..8f60d57d00 --- /dev/null +++ b/toolchain/mfc/viz/reader.py @@ -0,0 +1,518 @@ +""" +Binary format reader for MFC post-processed output. + +Reads Fortran unformatted binary files produced by MFC's post_process +with format=2. Each Fortran `write` produces one record: + [4-byte record-marker][payload][4-byte record-marker] + +File layout per processor: + Record 1 (header): m(int32), n(int32), p(int32), dbvars(int32) + Record 2 (grid): x_cb [, y_cb [, z_cb]] (float32 or float64) + Records 3..N (vars): varname(50-char) + data((m+1)*(n+1)*(p+1) floats) +""" + +import glob +import itertools +import os +import struct +import warnings +from dataclasses import dataclass, field +from functools import lru_cache +from typing import Dict, List, Optional, Tuple + +import numpy as np + + +NAME_LEN = 50 # Fortran character length for variable names + + +@dataclass +class ProcessorData: + """Data from a single processor file. + + m, n, p store dimension metadata but their exact semantics differ + between readers: + - Binary: m = Fortran header value, x_cb has m+2 elements, + cell count is m+1. + - Silo: m = cell count = len(x_cb) - 1. + Assembly code intentionally avoids using m/n/p for array sizing — + it derives everything from x_cb/y_cb/z_cb lengths. If future code + needs m directly, this discrepancy must be resolved. + """ + m: int + n: int + p: int + x_cb: np.ndarray + y_cb: np.ndarray + z_cb: np.ndarray + variables: Dict[str, np.ndarray] = field(default_factory=dict) + + +@dataclass +class AssembledData: + """Assembled multi-processor data on a global grid.""" + ndim: int + x_cc: np.ndarray + y_cc: np.ndarray + z_cc: np.ndarray + variables: Dict[str, np.ndarray] = field(default_factory=dict) + + + +def _detect_endianness(path: str) -> str: + """Detect endianness from the first record marker. + + The header record contains 4 int32s (m, n, p, dbvars) = 16 bytes, + so the leading Fortran record marker must be 16. + """ + with open(path, 'rb') as f: + raw = f.read(4) + if len(raw) < 4: + raise EOFError(f"File too short to detect endianness: {path}") + le = struct.unpack('i', raw)[0] + if be == 16: + return '>' + raise ValueError( + f"Cannot detect endianness: first record marker is {le} (LE) / {be} (BE), expected 16" + ) + + +def _read_record_endian(f, endian: str) -> bytes: + """Read one Fortran unformatted record with known endianness.""" + raw = f.read(4) + if len(raw) < 4: + raise EOFError("Unexpected end of file reading record marker") + rec_len = struct.unpack(f'{endian}i', raw)[0] + if rec_len < 0: + raise ValueError(f"Invalid Fortran record length: {rec_len}") + payload = f.read(rec_len) + if len(payload) < rec_len: + raise EOFError("Unexpected end of file reading record payload") + trail = f.read(4) + if len(trail) < 4: + raise EOFError("Unexpected end of file reading trailing record marker") + trail_len = struct.unpack(f'{endian}i', trail)[0] + if trail_len != rec_len: + raise ValueError( + f"Fortran record marker mismatch: leading={rec_len}, trailing={trail_len}. " + "File may be corrupted." + ) + return payload + + +def read_binary_file(path: str, var_filter: Optional[str] = None) -> ProcessorData: # pylint: disable=too-many-locals,too-many-statements + """ + Read a single MFC binary post-process file. + + Args: + path: Path to the .dat file. + var_filter: If given, only load this variable (skip others). + + Returns: + ProcessorData with grid and variable data. + """ + endian = _detect_endianness(path) + + with open(path, 'rb') as f: + # Record 1: header [m, n, p, dbvars] — 4 int32 + hdr = _read_record_endian(f, endian) + m, n, p, dbvars = struct.unpack(f'{endian}4i', hdr) + if m < 0 or n < 0 or p < 0 or dbvars < 0: + raise ValueError( + f"Invalid header in {path}: m={m}, n={n}, p={p}, dbvars={dbvars}" + ) + + # Record 2: grid coordinates — all in one record + grid_raw = _read_record_endian(f, endian) + grid_bytes = len(grid_raw) + + # Determine number of grid values + if p > 0: + n_vals = (m + 2) + (n + 2) + (p + 2) + elif n > 0: + n_vals = (m + 2) + (n + 2) + else: + n_vals = m + 2 + + # Auto-detect grid precision from record size + if grid_bytes == n_vals * 8: + grid_dtype = np.dtype(f'{endian}f8') + elif grid_bytes == n_vals * 4: + grid_dtype = np.dtype(f'{endian}f4') + else: + bytes_per_val = grid_bytes / n_vals if n_vals else 0 + raise ValueError( + f"Cannot determine grid precision: {grid_bytes} bytes for {n_vals} values " + f"({bytes_per_val:.1f} bytes/value)" + ) + + grid_arr = np.frombuffer(grid_raw, dtype=grid_dtype) + + # Split into x_cb, y_cb, z_cb + offset = 0 + x_cb = grid_arr[offset:offset + m + 2].astype(np.float64) + offset += m + 2 + if n > 0: + y_cb = grid_arr[offset:offset + n + 2].astype(np.float64) + offset += n + 2 + else: + y_cb = np.array([0.0]) + if p > 0: + z_cb = grid_arr[offset:offset + p + 2].astype(np.float64) + else: + z_cb = np.array([0.0]) + + # Records 3..N: variables + variables: Dict[str, np.ndarray] = {} + data_size = (m + 1) * max(n + 1, 1) * max(p + 1, 1) + + for _ in range(dbvars): + var_raw = _read_record_endian(f, endian) + varname = var_raw[:NAME_LEN].decode('ascii', errors='replace').strip() + + if var_filter is not None and varname != var_filter: + continue + + # Auto-detect variable data precision from record size + data_bytes = len(var_raw) - NAME_LEN + if data_bytes == data_size * 8: + var_dtype = np.dtype(f'{endian}f8') + elif data_bytes == data_size * 4: + var_dtype = np.dtype(f'{endian}f4') + else: + var_bpv = data_bytes / data_size if data_size else 0 + raise ValueError( + f"Cannot determine variable precision for '{varname}': " + f"{data_bytes} bytes for {data_size} values ({var_bpv:.1f} bytes/value)" + ) + + data = np.frombuffer(var_raw[NAME_LEN:], dtype=var_dtype).astype(np.float64) + + # Reshape for multi-dimensional data (Fortran column-major order) + if p > 0: + data = data.reshape((m + 1, n + 1, p + 1), order='F') + elif n > 0: + data = data.reshape((m + 1, n + 1), order='F') + + variables[varname] = data + + return ProcessorData(m=m, n=n, p=p, x_cb=x_cb, y_cb=y_cb, z_cb=z_cb, variables=variables) + + +def discover_format(case_dir: str) -> str: + """Detect whether case has binary or silo_hdf5 output.""" + if os.path.isdir(os.path.join(case_dir, 'binary')): + return 'binary' + if os.path.isdir(os.path.join(case_dir, 'silo_hdf5')): + return 'silo' + raise FileNotFoundError( + f"No 'binary/' or 'silo_hdf5/' directory found in {case_dir}. " + "Run post_process with format=1 (Silo) or format=2 (binary) first." + ) + + +def discover_timesteps(case_dir: str, fmt: str) -> List[int]: + """Return sorted list of available timesteps.""" + if fmt not in ('binary', 'silo'): + raise ValueError(f"Unknown format '{fmt}'. Supported: 'binary', 'silo'.") + + if fmt == 'binary': + # Check root/ first (1D), then p0/ + root_dir = os.path.join(case_dir, 'binary', 'root') + if os.path.isdir(root_dir): + steps = set() + for fname in os.listdir(root_dir): + if fname.endswith('.dat'): + try: + steps.add(int(fname[:-4])) + except ValueError: + pass + if steps: + return sorted(steps) + + # Multi-dimensional: look in p0/ + p0_dir = os.path.join(case_dir, 'binary', 'p0') + if os.path.isdir(p0_dir): + steps = set() + for fname in os.listdir(p0_dir): + if fname.endswith('.dat'): + try: + steps.add(int(fname[:-4])) + except ValueError: + pass + return sorted(steps) + + elif fmt == 'silo': + p0_dir = os.path.join(case_dir, 'silo_hdf5', 'p0') + if os.path.isdir(p0_dir): + steps = set() + for fname in os.listdir(p0_dir): + if fname.endswith('.silo') and not fname.startswith('collection'): + try: + steps.add(int(fname[:-5])) + except ValueError: + pass + return sorted(steps) + + return [] # no timestep files found in expected directories + + +def _discover_processors(case_dir: str, fmt: str) -> List[int]: + """Return sorted list of processor ranks.""" + if fmt == 'binary': + base = os.path.join(case_dir, 'binary') + else: + base = os.path.join(case_dir, 'silo_hdf5') + + ranks = [] + if not os.path.isdir(base): + return ranks + for entry in os.listdir(base): + if entry.startswith('p') and entry[1:].isdigit(): + ranks.append(int(entry[1:])) + return sorted(ranks) + + +def _is_1d(case_dir: str) -> bool: + """Check if the output is 1D (binary/root/ directory exists and contains .dat files).""" + root = os.path.join(case_dir, 'binary', 'root') + return os.path.isdir(root) and any(f.endswith('.dat') for f in os.listdir(root)) + + +def assemble_from_proc_data( # pylint: disable=too-many-locals,too-many-statements + proc_data: List[Tuple[int, ProcessorData]], +) -> AssembledData: + """ + Assemble multi-processor data into a single global grid. + + Shared helper used by both binary and silo assembly paths. + Handles ghost/buffer cell overlap between processors by using + per-cell coordinate lookup (np.unique + np.searchsorted + np.ix_). + """ + if not proc_data: + raise ValueError("No processor data to assemble") + + # Single processor — fast path + if len(proc_data) == 1: + _, pd = proc_data[0] + x_cc = (pd.x_cb[:-1] + pd.x_cb[1:]) / 2.0 + y_cc = (pd.y_cb[:-1] + pd.y_cb[1:]) / 2.0 if pd.n > 0 else np.array([0.0]) + z_cc = (pd.z_cb[:-1] + pd.z_cb[1:]) / 2.0 if pd.p > 0 else np.array([0.0]) + ndim = 1 + (pd.n > 0) + (pd.p > 0) + return AssembledData( + ndim=ndim, x_cc=x_cc, y_cc=y_cc, z_cc=z_cc, + variables=pd.variables, + ) + + # Multi-processor assembly + sample = proc_data[0][1] + ndim = 1 + (sample.n > 0) + (sample.p > 0) + + # Compute cell centers for each processor + proc_centers = [] + for rank, pd in proc_data: + x_cc = (pd.x_cb[:-1] + pd.x_cb[1:]) / 2.0 + y_cc = (pd.y_cb[:-1] + pd.y_cb[1:]) / 2.0 if pd.n > 0 else np.array([0.0]) + z_cc = (pd.z_cb[:-1] + pd.z_cb[1:]) / 2.0 if pd.p > 0 else np.array([0.0]) + proc_centers.append((rank, pd, x_cc, y_cc, z_cc)) + + # Build unique sorted global coordinate arrays (handles ghost overlap). + # Normalize each axis by its extent before rounding so that precision is + # always 12 significant digits *relative to the domain size*. This is + # correct for both micro-scale domains (extent ~ 1e-10) and large-scale + # domains (extent > 1e12) where the old formula (decimals = -log10(extent) + # + 12) could go negative, causing np.round to round to the nearest 10 or + # 100 and incorrectly merging distinct cell centers. + def _dedup(arr): + extent = arr.max() - arr.min() + if extent > 0: + origin = arr.min() + norm = np.round((arr - origin) / extent, 12) + return origin + np.unique(norm) * extent, origin, extent + return np.unique(arr), arr.min(), 0.0 + + def _norm_round(arr, origin, extent): + """Round *arr* with the same relative tolerance used by _dedup.""" + if extent > 0: + return origin + np.round((arr - origin) / extent, 12) * extent + return arr + + all_x = np.concatenate([xc for _, _, xc, _, _ in proc_centers]) + global_x, x_orig, x_ext = _dedup(all_x) + if ndim >= 2: + all_y = np.concatenate([yc for _, _, _, yc, _ in proc_centers]) + global_y, y_orig, y_ext = _dedup(all_y) + else: + global_y, y_orig, y_ext = np.array([0.0]), 0.0, 0.0 + if ndim >= 3: + all_z = np.concatenate([zc for _, _, _, _, zc in proc_centers]) + global_z, z_orig, z_ext = _dedup(all_z) + else: + global_z, z_orig, z_ext = np.array([0.0]), 0.0, 0.0 + + varnames = sorted({vn for _, pd in proc_data for vn in pd.variables}) + nx, ny, nz = len(global_x), len(global_y), len(global_z) + + global_vars: Dict[str, np.ndarray] = {} + for vn in varnames: + if ndim == 3: + global_vars[vn] = np.zeros((nx, ny, nz)) + elif ndim == 2: + global_vars[vn] = np.zeros((nx, ny)) + else: + global_vars[vn] = np.zeros(nx) + + # Place each processor's data using per-cell coordinate lookup. + # Apply the same normalized rounding used by _dedup so that lookup + # coordinates match the global grid entries exactly. + for _rank, pd, x_cc, y_cc, z_cc in proc_centers: + xi = np.clip(np.searchsorted(global_x, _norm_round(x_cc, x_orig, x_ext)), 0, nx - 1) + yi = np.clip(np.searchsorted(global_y, _norm_round(y_cc, y_orig, y_ext)), 0, ny - 1) if ndim >= 2 else np.array([0]) + zi = np.clip(np.searchsorted(global_z, _norm_round(z_cc, z_orig, z_ext)), 0, nz - 1) if ndim >= 3 else np.array([0]) + + for vn, data in pd.variables.items(): + if vn not in global_vars: + continue + if ndim == 3: + global_vars[vn][np.ix_(xi, yi, zi)] = data + elif ndim == 2: + global_vars[vn][np.ix_(xi, yi)] = data + else: + global_vars[vn][xi] = data + + return AssembledData( + ndim=ndim, x_cc=global_x, y_cc=global_y, z_cc=global_z, + variables=global_vars, + ) + + +def assemble(case_dir: str, step: int, fmt: str = 'binary', # pylint: disable=too-many-locals + var: Optional[str] = None) -> AssembledData: + """ + Read and assemble multi-processor data for a given timestep. + + For 1D, reads the root file directly. + For 2D/3D, reads all processor files and assembles into global arrays. + """ + if fmt != 'binary': + raise ValueError(f"Format '{fmt}' not supported by binary reader. Use silo_reader.") + + # 1D case: read root file directly + if _is_1d(case_dir): + root_path = os.path.join(case_dir, 'binary', 'root', f'{step}.dat') + if not os.path.isfile(root_path): + raise FileNotFoundError(f"Root file not found: {root_path}") + pdata = read_binary_file(root_path, var_filter=var) + x_cc = (pdata.x_cb[:-1] + pdata.x_cb[1:]) / 2.0 + return AssembledData( + ndim=1, x_cc=x_cc, + y_cc=np.array([0.0]), z_cc=np.array([0.0]), + variables=pdata.variables, + ) + + # Multi-dimensional: read all processor files + ranks = _discover_processors(case_dir, fmt) + if not ranks: + raise FileNotFoundError(f"No processor directories found in {case_dir}/binary/") + + proc_data: List[Tuple[int, ProcessorData]] = [] + for rank in ranks: + fpath = os.path.join(case_dir, 'binary', f'p{rank}', f'{step}.dat') + if not os.path.isfile(fpath): + raise FileNotFoundError( + f"Processor file not found: {fpath}. " + "Incomplete output (missing rank) would produce incorrect data." + ) + pdata = read_binary_file(fpath, var_filter=var) + if pdata.m == 0 and pdata.n == 0 and pdata.p == 0: + warnings.warn(f"Processor p{rank} has zero dimensions, skipping", stacklevel=2) + continue + proc_data.append((rank, pdata)) + + if not proc_data: + raise FileNotFoundError(f"No valid processor data found for step {step}") + + return assemble_from_proc_data(proc_data) + + +# --------------------------------------------------------------------------- +# Lagrange bubble position reader +# --------------------------------------------------------------------------- + +@lru_cache(maxsize=32) +def _nBubs_per_step(path: str) -> int: + """Count how many bubble rows share the first time value in *path*. + + The result is cached so repeated calls for the same file (across + different steps in an MP4 render) only scan the file once. + """ + with open(path) as f: + f.readline() # skip header + first = f.readline() + if not first.strip(): + return 0 + t0 = first.split()[0] + n = 1 + for line in f: + parts = line.split() + if parts and parts[0] == t0: + n += 1 + else: + break + return n + + +def read_lag_bubbles_at_step(case_dir: str, step: int) -> Optional[np.ndarray]: + """Read Lagrange bubble positions at a given save-step index. + + Reads ``D/lag_bubble_evol_.dat`` files written by MFC when + ``lag_params%write_bubbles = T``. Each file contains one row per + bubble per simulation time-step (appended after every completed RK + stage). This function seeks efficiently to the block for *step* by + counting the fixed number of bubbles per rank. + + Returns an ``(N, 4)`` float64 array of ``(x, y, z, r)`` in + simulation-normalized units across all MPI ranks, or ``None`` when + no bubble data is found. + """ + d_dir = os.path.join(case_dir, 'D') + if not os.path.isdir(d_dir): + return None + + files = sorted(glob.glob(os.path.join(d_dir, 'lag_bubble_evol_*.dat'))) + if not files: + return None + + chunks: List[np.ndarray] = [] + for fpath in files: + try: + nBubs = _nBubs_per_step(fpath) + if nBubs == 0: + continue + # Line layout: 1 header + step * nBubs prior data rows + skip = 1 + step * nBubs + rows = [] + with open(fpath) as f: + for _ in itertools.islice(f, skip): + pass + for line in itertools.islice(f, nBubs): + parts = line.split() + if len(parts) >= 8: + # cols: time id x y z mv conc r [rdot p] + rows.append((float(parts[2]), float(parts[3]), + float(parts[4]), float(parts[7]))) + if rows: + chunks.append(np.array(rows, dtype=np.float64)) + except (OSError, ValueError): + continue + + return np.concatenate(chunks, axis=0) if chunks else None + + +def has_lag_bubble_evol(case_dir: str) -> bool: + """Return True if ``D/lag_bubble_evol_*.dat`` files exist in *case_dir*.""" + d_dir = os.path.join(case_dir, 'D') + return bool(glob.glob(os.path.join(d_dir, 'lag_bubble_evol_*.dat'))) diff --git a/toolchain/mfc/viz/renderer.py b/toolchain/mfc/viz/renderer.py new file mode 100644 index 0000000000..ac69f960ed --- /dev/null +++ b/toolchain/mfc/viz/renderer.py @@ -0,0 +1,548 @@ +""" +Image and video rendering for MFC visualization. + +Produces PNG images (1D line plots, 2D colormaps) and MP4 videos +from assembled MFC data. Uses matplotlib with the Agg backend +for headless rendering. +""" + +import math +import os +import re +import tempfile + +import numpy as np + +import imageio + +import matplotlib +try: + matplotlib.use('Agg') +except Exception: # pylint: disable=broad-except + pass +import matplotlib.pyplot as plt # pylint: disable=wrong-import-position +from matplotlib.colors import LogNorm # pylint: disable=wrong-import-position + +matplotlib.rcParams.update({ + 'mathtext.fontset': 'cm', + 'font.family': 'serif', +}) + +# LaTeX-style labels for known MFC variable names +_LABEL_MAP = { + 'pres': r'$p$', + 'rho': r'$\rho$', + 'E': r'$E$', + 'T': r'$T$', + 'D': r'$D$', + 'c': r'$c$', + 'gamma': r'$\gamma$', + 'pi_inf': r'$\pi_\infty$', + 'pres_inf': r'$p_\infty$', + 'heat_ratio': r'$\gamma$', + 'schlieren': r'$|\nabla \rho|$', + 'psi': r'$\psi$', + 'n': r'$n$', + 'qm': r'$q_m$', + 'Bx': r'$B_x$', 'By': r'$B_y$', 'Bz': r'$B_z$', + 'voidFraction': r'void fraction', + 'liutex_mag': r'$|\lambda|$', + 'damage_state': r'damage', +} + +_INDEXED_PATTERNS = [ + (r'^vel(\d+)$', lambda m: [r'$u$', r'$v$', r'$w$'][int(m.group(1)) - 1] + if int(m.group(1)) <= 3 else rf'$v_{m.group(1)}$'), + (r'^mom(\d+)$', lambda m: rf'$\rho {["u", "v", "w"][int(m.group(1)) - 1]}$' + if int(m.group(1)) <= 3 else rf'$m_{m.group(1)}$'), + (r'^alpha(\d+)$', lambda m: rf'$\alpha_{m.group(1)}$'), + (r'^alpha_rho(\d+)$', lambda m: rf'$\alpha_{m.group(1)}\rho_{m.group(1)}$'), + (r'^alpha_rho_e(\d+)$', lambda m: rf'$\alpha_{m.group(1)}\rho_{m.group(1)}e_{m.group(1)}$'), + (r'^omega(\d+)$', lambda m: rf'$\omega_{m.group(1)}$'), + (r'^tau(\d+)$', lambda m: rf'$\tau_{m.group(1)}$'), + (r'^xi(\d+)$', lambda m: rf'$\xi_{m.group(1)}$'), + (r'^flux(\d+)$', lambda m: rf'$F_{m.group(1)}$'), + (r'^liutex_axis(\d+)$', lambda m: rf'$\lambda_{m.group(1)}$'), + (r'^rho(\d+)$', lambda m: rf'$\rho_{m.group(1)}$'), + (r'^Y_(.+)$', lambda m: rf'$Y_{{\mathrm{{{m.group(1)}}}}}$'), + (r'^nR(\d+)$', lambda m: rf'$nR_{{{m.group(1)}}}$'), + (r'^nV(\d+)$', lambda m: rf'$nV_{{{m.group(1)}}}$'), + (r'^nP(\d+)$', lambda m: rf'$nP_{{{m.group(1)}}}$'), + (r'^nM(\d+)$', lambda m: rf'$nM_{{{m.group(1)}}}$'), + (r'^color_function(\d+)$', lambda m: rf'color $f_{m.group(1)}$'), +] + + +def _overlay_bubbles(ax, bubbles, scale: float = 1.0) -> None: + """Overlay Lagrange bubble positions as circles on *ax*. + + Args: + ax: matplotlib Axes to draw on. + bubbles: (N, 4) array of (x, y, z, r) in simulation-normalized units. + scale: Multiply rendered radius by this factor for visibility. + """ + if bubbles is None or len(bubbles) == 0: + return + from matplotlib.patches import Circle # pylint: disable=import-outside-toplevel + from matplotlib.collections import PatchCollection # pylint: disable=import-outside-toplevel + circles = [Circle((b[0], b[1]), b[3] * scale) for b in bubbles] + pc = PatchCollection(circles, facecolors='none', edgecolors='white', + linewidths=0.5, alpha=0.8) + ax.add_collection(pc) + + +def pretty_label(varname): + """Map an MFC variable name to a LaTeX-style label for plots.""" + if varname in _LABEL_MAP: + return _LABEL_MAP[varname] + for pattern, formatter in _INDEXED_PATTERNS: + m = re.match(pattern, varname) + if m: + return formatter(m) + return varname + + +def render_1d(x_cc, data, varname, step, output, **opts): # pylint: disable=too-many-arguments,too-many-positional-arguments + """Render a 1D line plot and save as PNG.""" + fig, ax = plt.subplots(figsize=opts.get('figsize', (10, 6))) + label = pretty_label(varname) + ax.plot(x_cc, data, linewidth=1.5) + ax.set_xlabel(r'$x$') + ax.set_ylabel(label) + ax.set_title(f'{label} (step {step})') + ax.grid(True, alpha=0.3) + ax.ticklabel_format(axis='y', style='sci', scilimits=(-3, 4), useMathText=True) + + log_scale = opts.get('log_scale', False) + if log_scale: + ax.set_yscale('log') + vmin = opts.get('vmin') + vmax = opts.get('vmax') + if vmin is not None or vmax is not None: + ax.set_ylim(vmin, vmax) + + fig.tight_layout() + fig.savefig(output, dpi=opts.get('dpi', 150)) + plt.close(fig) + + +def render_1d_tiled(x_cc, variables, step, output, **opts): # pylint: disable=too-many-locals + """Render all 1D variables in a tiled subplot grid and save as PNG.""" + varnames = sorted(variables.keys()) + n = len(varnames) + if n == 0: + return + if n == 1: + render_1d(x_cc, variables[varnames[0]], varnames[0], step, output, **opts) + return + + ncols = 2 if n <= 8 else 3 + nrows = math.ceil(n / ncols) + fig_w = 5 * ncols + fig_h = 2.8 * nrows + fig, axes = plt.subplots(nrows, ncols, + figsize=opts.get('figsize', (fig_w, fig_h)), + sharex=True, squeeze=False) + + log_scale = opts.get('log_scale', False) + for idx, vn in enumerate(varnames): + row, col = divmod(idx, ncols) + ax = axes[row][col] + ax.plot(x_cc, variables[vn], linewidth=1.2) + ax.set_ylabel(pretty_label(vn), fontsize=9) + ax.tick_params(labelsize=8) + ax.grid(True, alpha=0.3) + if log_scale: + ax.set_yscale('log') + + # Hide unused subplots + for idx in range(n, nrows * ncols): + row, col = divmod(idx, ncols) + axes[row][col].set_visible(False) + + # X-label only on bottom row + for col in range(ncols): + bottom_row = min(nrows - 1, (n - 1) // ncols) if col < (n % ncols or ncols) else nrows - 2 + axes[bottom_row][col].set_xlabel(r'$x$', fontsize=9) + + fig.suptitle(f'step {step}', fontsize=11, y=0.99) + fig.tight_layout(rect=[0, 0, 1, 0.97]) + fig.savefig(output, dpi=opts.get('dpi', 150)) + plt.close(fig) + + +def _figsize_for_domain(x_cc, y_cc, base=10): + """Compute figure size that matches the physical domain aspect ratio.""" + dx = float(x_cc[-1] - x_cc[0]) if len(x_cc) > 1 else 1.0 + dy = float(y_cc[-1] - y_cc[0]) if len(y_cc) > 1 else 1.0 + aspect = dy / dx if dx > 0 else 1.0 + # Clamp to avoid extremely tall/wide figures + aspect = max(0.2, min(aspect, 5.0)) + # Extra width for colorbar + fig_w = base + 1.5 + fig_h = max(base * aspect, 3.0) + return (fig_w, fig_h) + + +def render_2d(x_cc, y_cc, data, varname, step, output, **opts): # pylint: disable=too-many-arguments,too-many-positional-arguments,too-many-locals + """Render a 2D colormap via pcolormesh and save as PNG.""" + default_size = _figsize_for_domain(x_cc, y_cc) + fig, ax = plt.subplots(figsize=opts.get('figsize', default_size)) + + cmap = opts.get('cmap', 'viridis') + vmin = opts.get('vmin') + vmax = opts.get('vmax') + log_scale = opts.get('log_scale', False) + + norm = None + if log_scale: + lo = vmin if vmin is not None else np.nanmin(data[data > 0]) if np.any(data > 0) else 1e-10 + hi = vmax if vmax is not None else np.nanmax(data) + if not np.isfinite(hi) or hi <= 0: + hi = 1.0 + if not np.isfinite(lo) or lo <= 0 or lo >= hi: + lo = hi * 1e-10 + norm = LogNorm(vmin=lo, vmax=hi) + vmin = None + vmax = None + + # data shape is (nx, ny), pcolormesh expects (ny, nx) when using x_cc, y_cc + pcm = ax.pcolormesh(x_cc, y_cc, data.T, cmap=cmap, vmin=vmin, vmax=vmax, + norm=norm, shading='auto') + label = pretty_label(varname) + fig.colorbar(pcm, ax=ax, label=label) + ax.set_xlabel(r'$x$') + ax.set_ylabel(r'$y$') + ax.set_title(f'{label} (step {step})') + ax.set_aspect('equal', adjustable='box') + + _overlay_bubbles(ax, opts.get('bubbles'), scale=opts.get('bubble_scale', 1.0)) + + fig.tight_layout() + fig.savefig(output, dpi=opts.get('dpi', 150)) + plt.close(fig) + + +def render_2d_tiled(assembled, step, output, **opts): # pylint: disable=too-many-locals + """Render all 2D variables in a tiled subplot grid and save as PNG.""" + varnames = sorted(assembled.variables.keys()) + n = len(varnames) + if n == 0: + return + if n == 1: + render_2d(assembled.x_cc, assembled.y_cc, + assembled.variables[varnames[0]], varnames[0], step, output, **opts) + return + + ncols = min(n, 3) + nrows = math.ceil(n / ncols) + cell_w, cell_h = _figsize_for_domain(assembled.x_cc, assembled.y_cc, base=4) + fig, axes = plt.subplots(nrows, ncols, + figsize=opts.get('figsize', (cell_w * ncols, cell_h * nrows)), + squeeze=False) + + cmap = opts.get('cmap', 'viridis') + log_scale = opts.get('log_scale', False) + for idx, vn in enumerate(varnames): + row, col = divmod(idx, ncols) + ax = axes[row][col] + data = assembled.variables[vn] + norm = None + vmin = opts.get('vmin') + vmax = opts.get('vmax') + if log_scale and np.any(data > 0): + lo = float(np.nanmin(data[data > 0])) + hi = float(np.nanmax(data)) + if np.isfinite(lo) and np.isfinite(hi) and lo < hi: + norm = LogNorm(vmin=lo, vmax=hi) + vmin = None + vmax = None + pcm = ax.pcolormesh(assembled.x_cc, assembled.y_cc, data.T, + cmap=cmap, vmin=vmin, vmax=vmax, + norm=norm, shading='auto') + label = pretty_label(vn) + fig.colorbar(pcm, ax=ax, label=label) + ax.set_title(label, fontsize=9) + ax.set_aspect('equal', adjustable='box') + ax.tick_params(labelsize=7) + _overlay_bubbles(ax, opts.get('bubbles'), scale=opts.get('bubble_scale', 1.0)) + + for idx in range(n, nrows * ncols): + row, col = divmod(idx, ncols) + axes[row][col].set_visible(False) + + fig.suptitle(f'step {step}', fontsize=11, y=1.01) + fig.tight_layout() + fig.savefig(output, dpi=opts.get('dpi', 150), bbox_inches='tight') + plt.close(fig) + + +def render_3d_slice(assembled, varname, step, output, slice_axis='z', # pylint: disable=too-many-arguments,too-many-positional-arguments,too-many-locals,too-many-statements,too-many-branches + slice_index=None, slice_value=None, **opts): + """Extract a 2D slice from 3D data and render as a colormap.""" + data_3d = assembled.variables[varname] + + axis_map = {'x': 0, 'y': 1, 'z': 2} + if slice_axis not in axis_map: + raise ValueError( + f"Invalid slice_axis '{slice_axis}'. Must be one of: 'x', 'y', 'z'." + ) + axis_idx = axis_map[slice_axis] + + coords = [assembled.x_cc, assembled.y_cc, assembled.z_cc] + coord_along = coords[axis_idx] + + if slice_index is not None: + idx = slice_index + elif slice_value is not None: + idx = int(np.argmin(np.abs(coord_along - slice_value))) + else: + idx = len(coord_along) // 2 + + idx = max(0, min(idx, len(coord_along) - 1)) + + if axis_idx == 0: + sliced = data_3d[idx, :, :] + x_plot, y_plot = assembled.y_cc, assembled.z_cc + xlabel, ylabel = r'$y$', r'$z$' + elif axis_idx == 1: + sliced = data_3d[:, idx, :] + x_plot, y_plot = assembled.x_cc, assembled.z_cc + xlabel, ylabel = r'$x$', r'$z$' + else: + sliced = data_3d[:, :, idx] + x_plot, y_plot = assembled.x_cc, assembled.y_cc + xlabel, ylabel = r'$x$', r'$y$' + + default_size = _figsize_for_domain(x_plot, y_plot) + fig, ax = plt.subplots(figsize=opts.get('figsize', default_size)) + + cmap = opts.get('cmap', 'viridis') + vmin = opts.get('vmin') + vmax = opts.get('vmax') + log_scale = opts.get('log_scale', False) + + norm = None + if log_scale: + pos = sliced[sliced > 0] + lo = vmin if vmin is not None else np.nanmin(pos) if pos.size > 0 else 1e-10 + hi = vmax if vmax is not None else np.nanmax(sliced) + if not np.isfinite(hi) or hi <= 0: + hi = 1.0 + if not np.isfinite(lo) or lo <= 0 or lo >= hi: + lo = hi * 1e-10 + norm = LogNorm(vmin=lo, vmax=hi) + vmin = None + vmax = None + + # sliced shape depends on axis: need to transpose appropriately + pcm = ax.pcolormesh(x_plot, y_plot, sliced.T, cmap=cmap, vmin=vmin, + vmax=vmax, norm=norm, shading='auto') + label = pretty_label(varname) + fig.colorbar(pcm, ax=ax, label=label) + ax.set_xlabel(xlabel) + ax.set_ylabel(ylabel) + slice_coord = coord_along[idx] + ax.set_title(f'{label} (step {step}, {slice_axis}={slice_coord:.4g})') + ax.set_aspect('equal', adjustable='box') + + # Overlay bubbles that lie within one radius of the slice plane + bubbles = opts.get('bubbles') + if bubbles is not None and len(bubbles) > 0: + slice_col = {'x': 0, 'y': 1, 'z': 2}[slice_axis] + plot_cols = [c for c in (0, 1, 2) if c != slice_col] + near = np.abs(bubbles[:, slice_col] - slice_coord) <= bubbles[:, 3] + _overlay_bubbles(ax, + bubbles[near][:, [plot_cols[0], plot_cols[1], slice_col, 3]] + if near.any() else None, + scale=opts.get('bubble_scale', 1.0)) + + fig.tight_layout() + fig.savefig(output, dpi=opts.get('dpi', 150)) + plt.close(fig) + + +def render_mp4(varname, steps, output, fps=10, # pylint: disable=too-many-arguments,too-many-positional-arguments,too-many-locals,too-many-statements,too-many-branches + read_func=None, tiled=False, bubble_func=None, **opts): + """ + Generate an MP4 video by iterating over timesteps. + + Args: + varname: Variable name to plot (ignored when tiled=True). + steps: List of timestep integers. + output: Output MP4 file path. + fps: Frames per second. + read_func: Callable(step) -> AssembledData for loading each frame. + tiled: If True, render all 1D variables in a tiled layout per frame. + bubble_func: Optional callable ``(step: int) -> ndarray`` returning + ``(N, 4)`` bubble positions ``(x, y, z, r)`` for each frame. + **opts: Rendering options (cmap, vmin, vmax, dpi, log_scale, figsize, + slice_axis, slice_index, slice_value). + + Returns: + True if the MP4 was successfully written, False on failure + (e.g., missing imageio dependency or encoding error). + """ + if read_func is None: + raise ValueError("read_func must be provided for MP4 rendering") + + if not steps: + raise ValueError("No timesteps provided for MP4 generation") + + opts = dict(opts) # avoid mutating the caller's dict + + # Pre-compute vmin/vmax from first, middle, and last frames if not provided + # (not needed for tiled mode — each subplot auto-scales independently) + auto_vmin = opts.get('vmin') + auto_vmax = opts.get('vmax') + + if not tiled and (auto_vmin is None or auto_vmax is None): + sample_steps = [steps[0]] + if len(steps) > 1: + sample_steps.append(steps[-1]) + if len(steps) > 2: + sample_steps.append(steps[len(steps) // 2]) + + all_mins, all_maxs = [], [] + log_scale = opts.get('log_scale', False) + for s in sample_steps: + ad = read_func(s) + d = ad.variables.get(varname) + if d is not None: + if log_scale: + pos = d[d > 0] + if pos.size > 0: + all_mins.append(np.nanmin(pos)) + all_maxs.append(np.nanmax(pos)) + else: + all_mins.append(np.nanmin(d)) + all_maxs.append(np.nanmax(d)) + + if auto_vmin is None and all_mins: + opts['vmin'] = min(all_mins) + if auto_vmax is None and all_maxs: + opts['vmax'] = max(all_maxs) + + # Write frames to a unique temp directory to avoid concurrent-run conflicts + output_dir = os.path.dirname(os.path.abspath(output)) + os.makedirs(output_dir, exist_ok=True) + viz_dir = tempfile.mkdtemp(dir=output_dir, prefix='_frames_') + + def _cleanup(): + for fname in sorted(f for f in os.listdir(viz_dir) if f.endswith('.png')): + try: + os.remove(os.path.join(viz_dir, fname)) + except OSError: + pass + try: + os.rmdir(viz_dir) + except OSError: + pass + + try: + from tqdm import tqdm # pylint: disable=import-outside-toplevel + step_iter = tqdm(steps, desc='Rendering frames') + except ImportError: + step_iter = steps + + try: + for i, step in enumerate(step_iter): + assembled = read_func(step) + frame_path = os.path.join(viz_dir, f'{i:06d}.png') + + # Inject per-step bubble positions into opts if bubble_func provided + frame_opts = opts + if bubble_func is not None: + try: + frame_opts = dict(opts, bubbles=bubble_func(step)) + except Exception: # pylint: disable=broad-except + pass + + if tiled and assembled.ndim == 1: + render_1d_tiled(assembled.x_cc, assembled.variables, + step, frame_path, **frame_opts) + elif tiled and assembled.ndim == 2: + render_2d_tiled(assembled, step, frame_path, **frame_opts) + elif assembled.ndim == 1: + var_data = assembled.variables.get(varname) + if var_data is None: + continue + render_1d(assembled.x_cc, var_data, + varname, step, frame_path, **frame_opts) + elif assembled.ndim == 2: + var_data = assembled.variables.get(varname) + if var_data is None: + continue + render_2d(assembled.x_cc, assembled.y_cc, + var_data, + varname, step, frame_path, **frame_opts) + elif assembled.ndim == 3: + var_data = assembled.variables.get(varname) + if var_data is None: + continue + render_3d_slice(assembled, varname, step, frame_path, **frame_opts) + else: + raise ValueError( + f"Unsupported dimensionality ndim={assembled.ndim} for step {step}. " + "Expected 1, 2, or 3." + ) + except BaseException: + _cleanup() + raise + + # Combine frames into MP4 using imageio + imageio-ffmpeg (bundled ffmpeg) + frame_files = sorted(f for f in os.listdir(viz_dir) if f.endswith('.png')) + + success = False + try: + if not frame_files: + return False + + def _to_rgb(arr): + """Normalise an image array to uint8 RGB (3-channel). + + imageio may return RGBA (4-ch) or even grayscale depending on the + PNG source. libx264/yuv420p requires consistent 3-channel input. + """ + if arr.ndim == 2: # grayscale → RGB + arr = np.stack([arr, arr, arr], axis=-1) + elif arr.shape[2] == 4: # RGBA → RGB (drop alpha) + arr = arr[:, :, :3] + return arr.astype(np.uint8) + + # First pass: find the maximum frame dimensions across all frames. + # Round up to the nearest even pixel: yuv420p requires even width and height. + max_h, max_w = 0, 0 + for fname in frame_files: + arr = imageio.imread(os.path.join(viz_dir, fname)) + max_h = max(max_h, arr.shape[0]) + max_w = max(max_w, arr.shape[1]) + ref_h = max_h + max_h % 2 + ref_w = max_w + max_w % 2 + + def _uniform_frame(arr): + """Convert to RGB and pad with white to (ref_h, ref_w).""" + arr = _to_rgb(arr) + h, w = arr.shape[:2] + if h == ref_h and w == ref_w: + return arr + out = np.full((ref_h, ref_w, 3), 255, dtype=np.uint8) + out[:h, :w] = arr + return out + + # Second pass: encode. macro_block_size=1 disables imageio's own resize + # since we already ensured even dimensions above. + with imageio.get_writer( + output, fps=fps, codec='libx264', pixelformat='yuv420p', + macro_block_size=1, ffmpeg_log_level='error', + ) as writer: + for fname in frame_files: + writer.append_data(_uniform_frame( + imageio.imread(os.path.join(viz_dir, fname)) + )) + success = True + except Exception: # pylint: disable=broad-except + pass + finally: + _cleanup() + return success diff --git a/toolchain/mfc/viz/silo_reader.py b/toolchain/mfc/viz/silo_reader.py new file mode 100644 index 0000000000..58e4c06b7b --- /dev/null +++ b/toolchain/mfc/viz/silo_reader.py @@ -0,0 +1,159 @@ +""" +Silo-HDF5 reader for MFC post-processed output. + +Silo files produced by MFC are valid HDF5 underneath. Each Silo object +is stored as an HDF5 Named Datatype whose ``silo`` compound attribute +carries the metadata (mesh name, data-array path, dimensions, etc.). +Actual data lives in numbered datasets under the ``.silo/`` group. + +This reader uses h5py to navigate that structure. +""" + +import os +import warnings +from typing import Dict, List, Optional, Tuple + +import h5py +import numpy as np + +from .reader import AssembledData, ProcessorData, assemble_from_proc_data + +# Silo type constants (from silo.h) +_DB_QUADMESH = 130 +_DB_QUADVAR = 501 + + + +def _resolve_path(h5file, path_bytes): + """Resolve a silo internal path (e.g. b'/.silo/#000003') and return its data as a numpy array.""" + path = path_bytes.decode() if isinstance(path_bytes, (bytes, np.bytes_)) else str(path_bytes) + return np.array(h5file[path]) + + +def read_silo_file( # pylint: disable=too-many-locals + path: str, var_filter: Optional[str] = None +) -> ProcessorData: + """ + Read a single Silo-HDF5 file produced by MFC post_process. + + Args: + path: Path to the ``.silo`` file. + var_filter: If given, only load this variable (case-sensitive). + + Returns: + ProcessorData with grid coordinates and variable arrays. + """ + + with h5py.File(path, "r") as f: + # --- locate the mesh ------------------------------------------------ + mesh_name = None + mesh_attr = None + for key, obj in f.items(): + if key in ("..", ".silo"): + continue + if not isinstance(obj, h5py.Datatype): + continue + silo_type = obj.attrs.get("silo_type") + if silo_type is not None and int(silo_type) == _DB_QUADMESH: + if "silo" not in obj.attrs: + continue + mesh_name = key + mesh_attr = obj.attrs["silo"] + break + + if mesh_attr is None: + raise ValueError(f"No rectilinear mesh found in {path}") + + ndims = int(mesh_attr["ndims"]) + coord_paths = [] + for i in range(ndims): + coord_paths.append(mesh_attr[f"coord{i}"]) + + coords = [_resolve_path(f, cp) for cp in coord_paths] + + x_cb = coords[0] + y_cb = coords[1] if ndims >= 2 else np.array([0.0]) + z_cb = coords[2] if ndims >= 3 else np.array([0.0]) + + # Grid dimensions: node counts minus 1 give cell counts + m = len(x_cb) - 1 + n = (len(y_cb) - 1) if ndims >= 2 else 0 + p = (len(z_cb) - 1) if ndims >= 3 else 0 + + # --- locate variables ------------------------------------------------ + variables: Dict[str, np.ndarray] = {} + for key, obj in f.items(): + if key in ("..", ".silo", mesh_name): + continue + if not isinstance(obj, h5py.Datatype): + continue + silo_type = obj.attrs.get("silo_type") + if silo_type is None or int(silo_type) != _DB_QUADVAR: + continue + + # Apply variable filter + if var_filter is not None and key != var_filter: + continue + + if "silo" not in obj.attrs: + continue + attr = obj.attrs["silo"] + try: + data_path = attr["value0"] + except (KeyError, ValueError): + warnings.warn(f"Variable '{key}' missing 'value0' in silo attr, skipping", stacklevel=2) + continue + data = _resolve_path(f, data_path).astype(np.float64) + + # MFC's DBPUTQV1 passes the Fortran column-major array as a + # flat buffer. HDF5 stores it row-major. Reinterpret the + # bytes in Fortran order so data[i,j,k] = value at (x_i,y_j,z_k), + # matching the binary reader convention. + # Assumption: Silo/HDF5 preserves the Fortran dimension ordering + # (nx, ny, nz) as the dataset shape. If a future Silo version + # reverses the shape, this reshape would silently transpose data. + if data.ndim >= 2: + data = np.ascontiguousarray(data).ravel().reshape(data.shape, order='F') + variables[key] = data + + return ProcessorData( + m=m, n=n, p=p, x_cb=x_cb, y_cb=y_cb, z_cb=z_cb, variables=variables + ) + + +def assemble_silo( + case_dir: str, + step: int, + var: Optional[str] = None, +) -> AssembledData: + """ + Read and assemble multi-processor Silo-HDF5 data for a given timestep. + """ + + base = os.path.join(case_dir, "silo_hdf5") + if not os.path.isdir(base): + raise FileNotFoundError(f"Silo-HDF5 directory not found: {base}") + ranks: List[int] = [] + for entry in os.listdir(base): + if entry.startswith("p") and entry[1:].isdigit(): + ranks.append(int(entry[1:])) + ranks.sort() + + if not ranks: + raise FileNotFoundError(f"No processor directories in {base}") + + proc_data: List[Tuple[int, ProcessorData]] = [] + for rank in ranks: + silo_file = os.path.join(base, f"p{rank}", f"{step}.silo") + if not os.path.isfile(silo_file): + raise FileNotFoundError( + f"Processor file not found: {silo_file}. " + "Incomplete output (missing rank) would produce incorrect data." + ) + pdata = read_silo_file(silo_file, var_filter=var) + proc_data.append((rank, pdata)) + + if not proc_data: + raise FileNotFoundError(f"No Silo data found for step {step}") + + return assemble_from_proc_data(proc_data) diff --git a/toolchain/mfc/viz/test_viz.py b/toolchain/mfc/viz/test_viz.py new file mode 100644 index 0000000000..949d9e1c5f --- /dev/null +++ b/toolchain/mfc/viz/test_viz.py @@ -0,0 +1,797 @@ +""" +Tests for the viz module. + +Covers step parsing, label formatting, format/timestep discovery, +data assembly (binary + silo, 1D/2D/3D), and 1D rendering. +Uses checked-in fixture data generated from minimal MFC runs. +""" +# pylint: disable=import-outside-toplevel,protected-access + +import os +import tempfile +import unittest + + +FIXTURES = os.path.join(os.path.dirname(__file__), 'fixtures') + +# Fixture paths for each dimension + format +FIX_1D_BIN = os.path.join(FIXTURES, '1d_binary') +FIX_1D_SILO = os.path.join(FIXTURES, '1d_silo') +FIX_2D_BIN = os.path.join(FIXTURES, '2d_binary') +FIX_2D_SILO = os.path.join(FIXTURES, '2d_silo') +FIX_3D_BIN = os.path.join(FIXTURES, '3d_binary') +FIX_3D_SILO = os.path.join(FIXTURES, '3d_silo') + + +# --------------------------------------------------------------------------- +# Tests: _parse_steps +# --------------------------------------------------------------------------- + +class TestParseSteps(unittest.TestCase): + """Test _parse_steps() step-argument parsing.""" + + def _parse(self, arg, available): + from .viz import _parse_steps + return _parse_steps(arg, available) + + def test_all_keyword(self): + """'all' returns every available step.""" + self.assertEqual(self._parse('all', [0, 100, 200]), [0, 100, 200]) + + def test_none_returns_all(self): + """None returns every available step.""" + self.assertEqual(self._parse(None, [0, 100, 200]), [0, 100, 200]) + + def test_last_keyword(self): + """'last' returns only the final step.""" + self.assertEqual(self._parse('last', [0, 100, 200]), [200]) + + def test_last_empty(self): + """'last' with no available steps returns empty list.""" + self.assertEqual(self._parse('last', []), []) + + def test_single_int(self): + """Single integer selects that step.""" + self.assertEqual(self._parse('100', [0, 100, 200]), [100]) + + def test_single_int_missing(self): + """Single integer not in available returns empty list.""" + self.assertEqual(self._parse('999', [0, 100, 200]), []) + + def test_range(self): + """Range format 'start:end:stride' selects matching steps.""" + result = self._parse('0:500:200', [0, 100, 200, 300, 400, 500]) + self.assertEqual(result, [0, 200, 400]) + + def test_range_no_stride(self): + """Range format 'start:end' defaults to stride 1.""" + result = self._parse('0:2', [0, 1, 2, 3]) + self.assertEqual(result, [0, 1, 2]) + + def test_comma_list(self): + """Comma-separated list selects the intersection with available steps.""" + result = self._parse('0,100,200,1000', [0, 100, 200, 300, 1000]) + self.assertEqual(result, [0, 100, 200, 1000]) + + def test_comma_list_filters_unavailable(self): + """Comma list silently drops steps not in available.""" + result = self._parse('0,999', [0, 100, 200]) + self.assertEqual(result, [0]) + + def test_ellipsis_expansion(self): + """Ellipsis infers stride and expands the range.""" + result = self._parse('0,100,200,...,1000', + list(range(0, 1001, 100))) + self.assertEqual(result, list(range(0, 1001, 100))) + + def test_ellipsis_partial_available(self): + """Ellipsis expansion filters to only available steps.""" + # only even-numbered hundreds available + avail = [0, 200, 400, 600, 800, 1000] + result = self._parse('0,100,...,1000', avail) + self.assertEqual(result, [0, 200, 400, 600, 800, 1000]) + + def test_ellipsis_requires_two_prefix_values(self): + """Ellipsis with only one prefix value raises MFCException.""" + from mfc.common import MFCException + with self.assertRaises(MFCException): + self._parse('0,...,1000', [0, 100, 1000]) + + def test_ellipsis_must_be_second_to_last(self): + """Ellipsis not in second-to-last position raises MFCException.""" + from mfc.common import MFCException + with self.assertRaises(MFCException): + self._parse('0,100,...,500,1000', [0, 100, 500, 1000]) + + def test_invalid_value(self): + """Non-numeric, non-keyword input raises MFCException.""" + from mfc.common import MFCException + with self.assertRaises(MFCException): + self._parse('bogus', [0, 100]) + + +# --------------------------------------------------------------------------- +# Tests: pretty_label +# --------------------------------------------------------------------------- + +class TestPrettyLabel(unittest.TestCase): + """Test pretty_label() LaTeX label generation.""" + + def _label(self, varname): + from .renderer import pretty_label + return pretty_label(varname) + + def test_known_scalar(self): + """Known scalars map to LaTeX.""" + self.assertEqual(self._label('pres'), r'$p$') + self.assertEqual(self._label('rho'), r'$\rho$') + + def test_vel_indexed(self): + """vel1/vel2/vel3 map to u/v/w.""" + self.assertEqual(self._label('vel1'), r'$u$') + self.assertEqual(self._label('vel2'), r'$v$') + self.assertEqual(self._label('vel3'), r'$w$') + + def test_alpha_indexed(self): + """alpha maps to LaTeX subscript.""" + self.assertIn('2', self._label('alpha2')) + + def test_unknown_passthrough(self): + """Unknown variable names pass through unchanged.""" + self.assertEqual(self._label('my_custom_var'), 'my_custom_var') + + +# --------------------------------------------------------------------------- +# Tests: discover_format +# --------------------------------------------------------------------------- + +class TestDiscoverFormat(unittest.TestCase): + """Test discover_format() binary/silo detection.""" + + def test_binary_detection(self): + """Detects binary format from binary/ directory.""" + from .reader import discover_format + self.assertEqual(discover_format(FIX_1D_BIN), 'binary') + + def test_silo_detection(self): + """Detects silo format from silo_hdf5/ directory.""" + from .reader import discover_format + self.assertEqual(discover_format(FIX_1D_SILO), 'silo') + + def test_missing_dir_raises(self): + """Missing directories raise FileNotFoundError.""" + from .reader import discover_format + d = tempfile.mkdtemp() + try: + with self.assertRaises(FileNotFoundError): + discover_format(d) + finally: + os.rmdir(d) + + +# --------------------------------------------------------------------------- +# Tests: discover_timesteps +# --------------------------------------------------------------------------- + +class TestDiscoverTimesteps(unittest.TestCase): + """Test discover_timesteps() against fixture data.""" + + def test_binary_1d(self): + """Finds sorted timesteps from 1D binary fixture.""" + from .reader import discover_timesteps + steps = discover_timesteps(FIX_1D_BIN, 'binary') + self.assertEqual(steps, sorted(steps)) + self.assertIn(0, steps) + self.assertGreater(len(steps), 1) + + def test_silo_1d(self): + """Finds sorted timesteps from 1D silo fixture.""" + from .reader import discover_timesteps + steps = discover_timesteps(FIX_1D_SILO, 'silo') + self.assertEqual(steps, sorted(steps)) + self.assertIn(0, steps) + self.assertGreater(len(steps), 1) + + +# --------------------------------------------------------------------------- +# Tests: binary read + assemble (1D, 2D, 3D) +# --------------------------------------------------------------------------- + +class TestAssembleBinary1D(unittest.TestCase): + """Test binary reader with 1D fixture data.""" + + def test_ndim(self): + """1D fixture assembles with ndim=1.""" + from .reader import assemble + data = assemble(FIX_1D_BIN, 0, 'binary') + self.assertEqual(data.ndim, 1) + + def test_grid_and_vars(self): + """1D fixture has non-empty grid and expected variables.""" + from .reader import assemble + data = assemble(FIX_1D_BIN, 0, 'binary') + self.assertGreater(len(data.x_cc), 0) + self.assertIn('pres', data.variables) + self.assertIn('vel1', data.variables) + self.assertEqual(data.variables['pres'].shape, data.x_cc.shape) + + def test_var_filter(self): + """Passing var= loads only that variable.""" + from .reader import assemble + data = assemble(FIX_1D_BIN, 0, 'binary', var='pres') + self.assertIn('pres', data.variables) + self.assertNotIn('vel1', data.variables) + + +class TestAssembleBinary2D(unittest.TestCase): + """Test binary reader with 2D fixture data.""" + + def test_ndim(self): + """2D fixture assembles with ndim=2.""" + from .reader import assemble + data = assemble(FIX_2D_BIN, 0, 'binary') + self.assertEqual(data.ndim, 2) + + def test_grid_shape(self): + """2D fixture has 2D variable arrays matching grid.""" + from .reader import assemble + data = assemble(FIX_2D_BIN, 0, 'binary') + self.assertGreater(len(data.x_cc), 0) + self.assertGreater(len(data.y_cc), 0) + pres = data.variables['pres'] + self.assertEqual(pres.shape, (len(data.x_cc), len(data.y_cc))) + + +class TestAssembleBinary3D(unittest.TestCase): + """Test binary reader with 3D fixture data.""" + + def test_ndim(self): + """3D fixture assembles with ndim=3.""" + from .reader import assemble + data = assemble(FIX_3D_BIN, 0, 'binary') + self.assertEqual(data.ndim, 3) + + def test_grid_shape(self): + """3D fixture has 3D variable arrays matching grid.""" + from .reader import assemble + data = assemble(FIX_3D_BIN, 0, 'binary') + self.assertGreater(len(data.x_cc), 0) + self.assertGreater(len(data.y_cc), 0) + self.assertGreater(len(data.z_cc), 0) + pres = data.variables['pres'] + self.assertEqual(pres.shape, + (len(data.x_cc), len(data.y_cc), len(data.z_cc))) + + +# --------------------------------------------------------------------------- +# Tests: silo read + assemble (1D, 2D, 3D) +# --------------------------------------------------------------------------- + +class TestAssembleSilo1D(unittest.TestCase): + """Test silo reader with 1D fixture data.""" + + def test_ndim(self): + """1D silo fixture assembles with ndim=1.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_1D_SILO, 0) + self.assertEqual(data.ndim, 1) + + def test_grid_and_vars(self): + """1D silo fixture has non-empty grid and expected variables.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_1D_SILO, 0) + self.assertGreater(len(data.x_cc), 0) + self.assertIn('pres', data.variables) + self.assertEqual(data.variables['pres'].shape, data.x_cc.shape) + + +class TestAssembleSilo2D(unittest.TestCase): + """Test silo reader with 2D fixture data.""" + + def test_ndim(self): + """2D silo fixture assembles with ndim=2.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_2D_SILO, 0) + self.assertEqual(data.ndim, 2) + + def test_grid_shape(self): + """2D silo fixture has 2D variable arrays matching grid.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_2D_SILO, 0) + pres = data.variables['pres'] + self.assertEqual(pres.shape, (len(data.x_cc), len(data.y_cc))) + + +class TestAssembleSilo3D(unittest.TestCase): + """Test silo reader with 3D fixture data.""" + + def test_ndim(self): + """3D silo fixture assembles with ndim=3.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_3D_SILO, 0) + self.assertEqual(data.ndim, 3) + + def test_grid_shape(self): + """3D silo fixture has 3D variable arrays matching grid.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_3D_SILO, 0) + pres = data.variables['pres'] + self.assertEqual(pres.shape, + (len(data.x_cc), len(data.y_cc), len(data.z_cc))) + + +# --------------------------------------------------------------------------- +# Tests: binary vs silo consistency +# --------------------------------------------------------------------------- + +class TestBinarySiloConsistency(unittest.TestCase): + """Verify binary and silo readers produce consistent results.""" + + def test_1d_same_grid(self): + """Binary and silo 1D fixtures have the same grid.""" + from .reader import assemble + from .silo_reader import assemble_silo + import numpy as np + bin_data = assemble(FIX_1D_BIN, 0, 'binary') + silo_data = assemble_silo(FIX_1D_SILO, 0) + np.testing.assert_allclose(bin_data.x_cc, silo_data.x_cc, atol=1e-10) + + def test_1d_same_vars(self): + """Binary and silo 1D fixtures have the same variable names.""" + from .reader import assemble + from .silo_reader import assemble_silo + bin_data = assemble(FIX_1D_BIN, 0, 'binary') + silo_data = assemble_silo(FIX_1D_SILO, 0) + self.assertEqual(sorted(bin_data.variables.keys()), + sorted(silo_data.variables.keys())) + + +# --------------------------------------------------------------------------- +# Tests: 1D rendering (requires matplotlib/imageio) +# --------------------------------------------------------------------------- + +class TestRender1D(unittest.TestCase): + """Smoke test: render 1D plots from fixture data.""" + + def test_render_png(self): + """Renders a single-variable PNG that is non-empty.""" + from .reader import assemble + from .renderer import render_1d + data = assemble(FIX_1D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_1d(data.x_cc, data.variables['pres'], 'pres', 0, out) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + def test_render_tiled_png(self): + """Tiled render of all variables produces a non-empty PNG.""" + from .reader import assemble + from .renderer import render_1d_tiled + data = assemble(FIX_1D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_1d_tiled(data.x_cc, data.variables, 0, out) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + +class TestRender2D(unittest.TestCase): + """Smoke test: render a 2D PNG from fixture data.""" + + def test_render_2d_png(self): + """Renders a 2D colormap PNG that is non-empty.""" + from .reader import assemble + from .renderer import render_2d + data = assemble(FIX_2D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_2d(data.x_cc, data.y_cc, data.variables['pres'], + 'pres', 0, out) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + +class TestRender3DSlice(unittest.TestCase): + """Smoke test: render a 3D slice PNG from fixture data.""" + + def test_render_3d_slice_png(self): + """Renders a 3D midplane-slice PNG that is non-empty.""" + from .reader import assemble + from .renderer import render_3d_slice + data = assemble(FIX_3D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_3d_slice(data, 'pres', 0, out) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + +# --------------------------------------------------------------------------- +# Tests: _steps_hint +# --------------------------------------------------------------------------- + +class TestStepsHint(unittest.TestCase): + """Test _steps_hint() step preview for error messages.""" + + def _hint(self, steps, n=8): + from .viz import _steps_hint + return _steps_hint(steps, n) + + def test_empty(self): + """Empty steps returns 'none found'.""" + self.assertEqual(self._hint([]), "none found") + + def test_short_list_shows_all(self): + """Short list shows all steps without truncation.""" + result = self._hint([0, 100, 200]) + self.assertIn('0', result) + self.assertIn('200', result) + self.assertNotIn('...', result) + + def test_long_list_truncated(self): + """Long list includes count and truncation marker.""" + steps = list(range(0, 2000, 100)) # 20 steps + result = self._hint(steps, n=8) + self.assertIn('...', result) + self.assertIn('[20 total]', result) + self.assertIn('0', result) # head present + self.assertIn('1900', result) # tail present + + +# --------------------------------------------------------------------------- +# Tests: _validate_cmap +# --------------------------------------------------------------------------- + +class TestValidateCmap(unittest.TestCase): + """Test _validate_cmap() colormap validation.""" + + def _validate(self, name): + from .viz import _validate_cmap + _validate_cmap(name) + + def test_known_cmaps_pass(self): + """Known colormaps do not raise.""" + for name in ('viridis', 'plasma', 'coolwarm', 'gray'): + with self.subTest(name=name): + self._validate(name) + + def test_unknown_cmap_raises(self): + """Unknown colormap raises MFCException.""" + from mfc.common import MFCException + with self.assertRaises(MFCException): + self._validate('notacolormap_xyz_1234') + + def test_typo_suggests_correct(self): + """Typo in colormap name raises MFCException suggesting the correct spelling.""" + from mfc.common import MFCException + with self.assertRaises(MFCException) as ctx: + self._validate('virids') # typo of viridis + self.assertIn('viridis', str(ctx.exception)) + + +# --------------------------------------------------------------------------- +# Tests: bounded TUI cache +# --------------------------------------------------------------------------- + +class TestTuiCache(unittest.TestCase): + """Test that the shared step cache respects CACHE_MAX.""" + + def setUp(self): + import mfc.viz._step_cache as cache_mod + self._mod = cache_mod + cache_mod.clear() + + def tearDown(self): + self._mod.clear() + + def _read(self, step): + return f"data_{step}" + + def test_cache_stores_entry(self): + """Loaded step is stored in cache.""" + self._mod.load(0, self._read) + self.assertIn(0, self._mod._cache) + + def test_cache_hit_avoids_reload(self): + """Second load of same step does not call read_func again.""" + calls = [0] + def counting(step): + calls[0] += 1 + return step + self._mod.load(5, counting) + self._mod.load(5, counting) + self.assertEqual(calls[0], 1) + + def test_cache_evicts_oldest_at_cap(self): + """Oldest entry is evicted when CACHE_MAX is exceeded.""" + cap = self._mod.CACHE_MAX + for i in range(cap + 3): + self._mod.load(i, self._read) + self.assertLessEqual(len(self._mod._cache), cap) + self.assertNotIn(0, self._mod._cache) # first evicted + self.assertIn(cap + 2, self._mod._cache) # most recent kept + + def test_seed_clears_and_populates(self): + """seed() clears existing cache and pre-loads one entry.""" + self._mod.load(99, self._read) # put something in first + self._mod.seed(0, "preloaded") + self.assertEqual(len(self._mod._cache), 1) + self.assertEqual(self._mod._cache[0], "preloaded") + + +# --------------------------------------------------------------------------- +# Tests: log scale rendering (new feature smoke tests) +# --------------------------------------------------------------------------- + +class TestRenderLogScale(unittest.TestCase): + """Smoke test: log scale option produces valid PNG output.""" + + def test_render_1d_log_scale(self): + """render_1d with log_scale=True produces a non-empty PNG.""" + from .reader import assemble + from .renderer import render_1d + data = assemble(FIX_1D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_1d(data.x_cc, data.variables['pres'], 'pres', 0, out, + log_scale=True) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + def test_render_2d_log_scale(self): + """render_2d with log_scale=True produces a non-empty PNG.""" + from .reader import assemble + from .renderer import render_2d + data = assemble(FIX_2D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_2d(data.x_cc, data.y_cc, data.variables['pres'], + 'pres', 0, out, log_scale=True) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + +# --------------------------------------------------------------------------- +# Tests: multi-rank assembly (ghost-cell deduplication) +# --------------------------------------------------------------------------- + +class TestMultiRankAssembly(unittest.TestCase): + """Test assemble_from_proc_data with synthetic multi-processor data.""" + + def _make_proc(self, x_cb, pres): + """Build a minimal 1D ProcessorData from boundary coordinates.""" + import numpy as np + from .reader import ProcessorData + return ProcessorData( + m=len(x_cb) - 1, + n=0, + p=0, + x_cb=np.array(x_cb, dtype=np.float64), + y_cb=np.array([0.0]), + z_cb=np.array([0.0]), + variables={'pres': np.array(pres, dtype=np.float64)}, + ) + + def test_two_rank_1d_dedup(self): + """Two processors with one overlapping ghost cell assemble correctly.""" + import numpy as np + from .reader import assemble_from_proc_data + # Domain: 4 cells with centers at 0.125, 0.375, 0.625, 0.875 + # Proc 0 sees cells 0-2 (center 0.625 is ghost from proc 1) + # Proc 1 sees cells 1-3 (center 0.375 is ghost from proc 0) + p0 = self._make_proc([0.00, 0.25, 0.50, 0.75], + [1.0, 2.0, 3.0]) # centers: 0.125, 0.375, 0.625 + p1 = self._make_proc([0.25, 0.50, 0.75, 1.00], + [2.0, 3.0, 4.0]) # centers: 0.375, 0.625, 0.875 + + result = assemble_from_proc_data([(0, p0), (1, p1)]) + + self.assertEqual(result.ndim, 1) + self.assertEqual(len(result.x_cc), 4) + np.testing.assert_allclose(result.x_cc, [0.125, 0.375, 0.625, 0.875]) + np.testing.assert_allclose(result.variables['pres'], [1.0, 2.0, 3.0, 4.0]) + + def test_large_extent_dedup(self): + """Deduplication works correctly for large-extent domains (>1e6).""" + import numpy as np + from .reader import assemble_from_proc_data + # Scale up by 1e7: extent=1e7, decimals = ceil(-log10(1e7)) + 12 = 5 + scale = 1e7 + p0 = self._make_proc( + [0.00 * scale, 0.25 * scale, 0.50 * scale, 0.75 * scale], + [1.0, 2.0, 3.0], + ) + p1 = self._make_proc( + [0.25 * scale, 0.50 * scale, 0.75 * scale, 1.00 * scale], + [2.0, 3.0, 4.0], + ) + result = assemble_from_proc_data([(0, p0), (1, p1)]) + self.assertEqual(len(result.x_cc), 4) + np.testing.assert_allclose( + result.variables['pres'], [1.0, 2.0, 3.0, 4.0] + ) + + def test_very_large_extent_dedup_negative_decimals(self): + """Deduplication works for extent ~1e13 where decimals becomes negative. + + At scale=1e13: extent = 1e13, decimals = ceil(-log10(1e13)) + 12 = -1. + np.round with negative decimals rounds to the nearest 10^|d|, so + np.round(x, -1) rounds to the nearest 10. Cell widths of 2.5e12 + are >> 10, so distinct cell-centers must not be collapsed. + """ + import numpy as np + from .reader import assemble_from_proc_data + scale = 1e13 + p0 = self._make_proc( + [0.00 * scale, 0.25 * scale, 0.50 * scale, 0.75 * scale], + [1.0, 2.0, 3.0], + ) + p1 = self._make_proc( + [0.25 * scale, 0.50 * scale, 0.75 * scale, 1.00 * scale], + [2.0, 3.0, 4.0], + ) + result = assemble_from_proc_data([(0, p0), (1, p1)]) + self.assertEqual(len(result.x_cc), 4) + np.testing.assert_allclose( + result.variables['pres'], [1.0, 2.0, 3.0, 4.0] + ) + + +# --------------------------------------------------------------------------- +# Tests: render_2d_tiled +# --------------------------------------------------------------------------- + +class TestRender2DTiled(unittest.TestCase): + """Smoke test: render_2d_tiled produces a valid PNG from 2D fixture data.""" + + def test_render_2d_tiled_png(self): + """Tiled render of all 2D variables produces a non-empty PNG.""" + from .reader import assemble + from .renderer import render_2d_tiled + data = assemble(FIX_2D_BIN, 0, 'binary') + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_2d_tiled(data, 0, out) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + +# --------------------------------------------------------------------------- +# Tests: render_3d_slice non-default axes and selectors +# --------------------------------------------------------------------------- + +class TestRender3DSliceAxes(unittest.TestCase): + """Test render_3d_slice with non-default slice axes and selectors.""" + + def setUp(self): + from .reader import assemble + self._data = assemble(FIX_3D_BIN, 0, 'binary') + + def _render(self, **kwargs): + from .renderer import render_3d_slice + with tempfile.NamedTemporaryFile(suffix='.png', delete=False) as f: + out = f.name + try: + render_3d_slice(self._data, 'pres', 0, out, **kwargs) + self.assertTrue(os.path.isfile(out)) + self.assertGreater(os.path.getsize(out), 0) + finally: + os.unlink(out) + + def test_x_axis_slice(self): + """X-axis midplane slice produces a non-empty PNG.""" + self._render(slice_axis='x') + + def test_y_axis_slice(self): + """Y-axis midplane slice produces a non-empty PNG.""" + self._render(slice_axis='y') + + def test_slice_by_index(self): + """slice_index=0 selects first plane along default z axis.""" + self._render(slice_index=0) + + def test_slice_by_value(self): + """slice_value selects the plane nearest the given coordinate.""" + z_mid = float(self._data.z_cc[len(self._data.z_cc) // 2]) + self._render(slice_value=z_mid) + + +# --------------------------------------------------------------------------- +# Tests: render_mp4 +# --------------------------------------------------------------------------- + +class TestRenderMp4(unittest.TestCase): + """Smoke test: render_mp4 exercises frame rendering and returns a bool.""" + + def _make_read_func(self, case_dir, fmt): + from .reader import assemble + def _read(step): + return assemble(case_dir, step, fmt) + return _read + + def test_mp4_1d_returns_bool(self): + """render_mp4 with 1D data returns True or False without raising.""" + from .reader import discover_timesteps + from .renderer import render_mp4 + steps = discover_timesteps(FIX_1D_BIN, 'binary')[:2] + read_func = self._make_read_func(FIX_1D_BIN, 'binary') + with tempfile.TemporaryDirectory() as tmpdir: + out = os.path.join(tmpdir, 'test.mp4') + result = render_mp4('pres', steps, out, fps=2, read_func=read_func) + self.assertIsInstance(result, bool) + + def test_mp4_tiled_1d_returns_bool(self): + """render_mp4 with tiled=True returns True or False without raising.""" + from .reader import discover_timesteps + from .renderer import render_mp4 + steps = discover_timesteps(FIX_1D_BIN, 'binary')[:2] + read_func = self._make_read_func(FIX_1D_BIN, 'binary') + with tempfile.TemporaryDirectory() as tmpdir: + out = os.path.join(tmpdir, 'test_tiled.mp4') + result = render_mp4('pres', steps, out, fps=2, + read_func=read_func, tiled=True) + self.assertIsInstance(result, bool) + + def test_mp4_no_read_func_raises(self): + """render_mp4 with read_func=None raises ValueError.""" + from .renderer import render_mp4 + with self.assertRaises(ValueError): + render_mp4('pres', [0], '/tmp/unused.mp4', read_func=None) + + def test_mp4_empty_steps_raises(self): + """render_mp4 with empty steps raises ValueError.""" + from .renderer import render_mp4 + with self.assertRaises(ValueError): + render_mp4('pres', [], '/tmp/unused.mp4', + read_func=lambda s: None) + + +# --------------------------------------------------------------------------- +# Tests: silo assemble_silo var_filter +# --------------------------------------------------------------------------- + +class TestAssembleSiloVarFilter(unittest.TestCase): + """Test assemble_silo with var= filter to cover the silo var_filter path.""" + + def test_1d_var_filter_includes_only_requested(self): + """Silo 1D: var='pres' loads pres and excludes vel1.""" + from .silo_reader import assemble_silo + data = assemble_silo(FIX_1D_SILO, 0, var='pres') + self.assertIn('pres', data.variables) + self.assertNotIn('vel1', data.variables) + + def test_2d_var_filter_includes_only_requested(self): + """Silo 2D: var='pres' loads pres and excludes other variables.""" + from .silo_reader import assemble_silo + filtered = assemble_silo(FIX_2D_SILO, 0, var='pres') + all_data = assemble_silo(FIX_2D_SILO, 0) + self.assertIn('pres', filtered.variables) + other_vars = [v for v in all_data.variables if v != 'pres'] + if other_vars: + self.assertNotIn(other_vars[0], filtered.variables) + + +if __name__ == "__main__": + unittest.main() diff --git a/toolchain/mfc/viz/tui.py b/toolchain/mfc/viz/tui.py new file mode 100644 index 0000000000..4f9da4df6a --- /dev/null +++ b/toolchain/mfc/viz/tui.py @@ -0,0 +1,540 @@ +""" +Terminal UI (TUI) for MFC visualization using Textual + plotext. + +Launched via ``./mfc.sh viz --tui [--var VAR] [--step STEP]``. +Opens a full-terminal interactive viewer that works over SSH with no +browser or port-forwarding required. + +Supports 1D line plots and 2D heatmaps only. + +Requires: textual, textual-plotext, plotext +""" +from __future__ import annotations + +from typing import Callable, List, Optional + +import numpy as np + +from rich.color import Color as RichColor +from rich.console import Group as RichGroup +from rich.style import Style +from rich.text import Text as RichText + +from textual import on +from textual.app import App, ComposeResult +from textual.binding import Binding +from textual.containers import Horizontal, Vertical +from textual.reactive import reactive +from textual.widgets import Footer, Header, Label, ListItem, ListView, Static + +from textual_plotext import PlotextPlot + +from mfc.common import MFCException +from mfc.printer import cons +from . import _step_cache + +# Colormaps available via [c] cycling +_CMAPS: List[str] = [ + 'viridis', 'plasma', 'inferno', 'magma', 'cividis', + 'coolwarm', 'RdBu_r', 'seismic', 'gray', +] + +_load = _step_cache.load +_CACHE_MAX = _step_cache.CACHE_MAX + +# Terminal character cells are approximately twice as tall as they are wide in +# pixels (e.g. 8 px wide × 16 px tall). A square physical domain should +# therefore occupy a ~2:1 (col:row) character grid to look correct. +_CELL_RATIO: float = 2.0 + +# Physical domain aspect ratio is clamped to this range so that very elongated +# domains don't produce unusable slivers. +_ASPECT_MIN: float = 0.2 +_ASPECT_MAX: float = 5.0 + + +# --------------------------------------------------------------------------- +# Plot widget +# --------------------------------------------------------------------------- + +class MFCPlot(PlotextPlot): # pylint: disable=too-many-instance-attributes,too-few-public-methods + """Plotext plot widget. Caller sets ._x_cc / ._y_cc / ._data / ._ndim / + ._varname / ._step before calling .refresh().""" + + DEFAULT_CSS = """ + MFCPlot { + border: solid $accent; + width: 1fr; + height: 1fr; + } + """ + + def __init__(self, **kwargs): + super().__init__(**kwargs) + self._x_cc: Optional[np.ndarray] = None + self._y_cc: Optional[np.ndarray] = None + self._data: Optional[np.ndarray] = None + self._ndim: int = 1 + self._varname: str = "" + self._step: int = 0 + self._cmap_name: str = _CMAPS[0] + self._log_scale: bool = False + self._vmin: Optional[float] = None + self._vmax: Optional[float] = None + self._last_vmin: float = 0.0 + self._last_vmax: float = 1.0 + self._bubbles: Optional[np.ndarray] = None # (N,4) x,y,z,r + + def render(self): # pylint: disable=too-many-branches,too-many-locals,too-many-statements + data = self._data + x_cc = self._x_cc + self.plt.clear_figure() + + # 1D: use normal plotext path — gives proper axes and title for free. + if data is None or x_cc is None or self._ndim == 1: + if data is not None and x_cc is not None: + if self._log_scale: + plot_y = np.where(data > 0, np.log10(np.maximum(data, 1e-300)), np.nan) + ylabel = f"log\u2081\u2080({self._varname})" + title_tag = " [log]" + else: + plot_y = data + ylabel = self._varname + title_tag = "" + if self._vmin is not None or self._vmax is not None: + title_tag += " [frozen]" + finite = plot_y[np.isfinite(plot_y)] + self._last_vmin = float(finite.min()) if finite.size else 0.0 + self._last_vmax = float(finite.max()) if finite.size else 1.0 + self.plt.plot(x_cc.tolist(), plot_y.tolist()) + self.plt.xlabel("x") + self.plt.ylabel(ylabel) + self.plt.title(f"{self._varname} (step {self._step}){title_tag}") + if self._vmin is not None or self._vmax is not None: + lo = self._vmin if self._vmin is not None else self._last_vmin + hi = self._vmax if self._vmax is not None else self._last_vmax + self.plt.ylim(lo, hi) + else: + self.plt.title("No data loaded") + return super().render() + + # 2D: pure-Rich heatmap with vertical colorbar. + import matplotlib # pylint: disable=import-outside-toplevel + import matplotlib.colors as mcolors # pylint: disable=import-outside-toplevel + + # Content area = widget size minus 1-char border on each side. + # Reserve 1 row each for header and footer → h_plot rows for the image. + w_plot = max(self.size.width - 2, 4) + h_plot_avail = max(self.size.height - 4, 4) # -2 border, -2 header+footer + + # Right side: gap + gradient strip + value labels + _CB_GAP, _CB_W, _CB_LBL = 1, 2, 9 + w_map_avail = max(w_plot - _CB_GAP - _CB_W - _CB_LBL, 4) + + # Preserve the physical x/y aspect ratio so the heatmap is not + # stretched to fill the terminal. The domain ratio is clamped to + # avoid extremely wide or tall slivers. + y_cc_2d = self._y_cc if self._y_cc is not None else np.array([0.0, 1.0]) + x_extent = max(abs(float(x_cc[-1]) - float(x_cc[0])), 1e-30) # pylint: disable=unsubscriptable-object + y_extent = max(abs(float(y_cc_2d[-1]) - float(y_cc_2d[0])), 1e-30) + domain_ratio = float(np.clip(x_extent / y_extent, _ASPECT_MIN, _ASPECT_MAX)) + # Convert to character-grid ratio: 1 row ≈ _CELL_RATIO columns wide. + char_ratio = domain_ratio * _CELL_RATIO # desired w_map / h_plot + + # Fit within the available character budget: try height-constrained first, + # fall back to width-constrained if the ideal width exceeds w_map_avail. + w_ideal = int(round(h_plot_avail * char_ratio)) + if w_ideal <= w_map_avail: + w_map = max(w_ideal, 4) + h_plot = h_plot_avail + else: + h_plot = max(int(round(w_map_avail / char_ratio)), 4) + w_map = w_map_avail + + ix = np.linspace(0, data.shape[0] - 1, w_map, dtype=int) + iy = np.linspace(0, data.shape[1] - 1, h_plot, dtype=int) + ds = data[np.ix_(ix, iy)] # pylint: disable=unsubscriptable-object + + # Compute which screen cells to stamp with an open-circle glyph. + # col → x_cc[ix[col]], row 0 = y_max (flipped), row h_plot-1 = y_min. + bubble_cells: set = set() + bubbles = self._bubbles + if bubbles is not None and len(bubbles) > 0: + x_phys = x_cc[ix] # type: ignore[index] # pylint: disable=unsubscriptable-object + y_phys = y_cc_2d[iy] + x_min, x_max = float(x_phys[0]), float(x_phys[-1]) + y_min, y_max = float(y_phys[0]), float(y_phys[-1]) + x_range = max(abs(x_max - x_min), 1e-30) + y_range = max(abs(y_max - y_min), 1e-30) + for b in bubbles: # pylint: disable=not-an-iterable + bx, by, br = float(b[0]), float(b[1]), float(b[3]) + if bx < x_min - br or bx > x_max + br: + continue + if by < y_min - br or by > y_max + br: + continue + # Screen centre (col increases right, row 0 = top = y_max) + col_c = (bx - x_min) / x_range * (w_map - 1) + row_c = (y_max - by) / y_range * (h_plot - 1) + # Screen radius in character units + col_r = br / x_range * (w_map - 1) + row_r = br / y_range * (h_plot - 1) + if col_r < 0.5 and row_r < 0.5: + # Sub-cell bubble — mark centre only + c, r = int(round(col_c)), int(round(row_c)) + if 0 <= r < h_plot and 0 <= c < w_map: + bubble_cells.add((r, c)) + else: + # Parametric circle outline + n_pts = min(max(12, int(2 * np.pi * max(col_r, row_r))), 72) + for ti in range(n_pts): + angle = 2 * np.pi * ti / n_pts + c = int(round(col_c + col_r * np.cos(angle))) + r = int(round(row_c + row_r * np.sin(angle))) + if 0 <= r < h_plot and 0 <= c < w_map: + bubble_cells.add((r, c)) + + vmin = self._vmin if self._vmin is not None else float(ds.min()) + vmax = self._vmax if self._vmax is not None else float(ds.max()) + if vmax <= vmin: + vmax = vmin + 1e-10 + cmap = matplotlib.colormaps[self._cmap_name] + log_active = False + if self._log_scale: + pos = ds[ds > 0] + lo = float(np.nanmin(pos)) if pos.size > 0 else None + if lo is not None and lo < vmax: + norm = mcolors.LogNorm(vmin=lo, vmax=vmax) + vmin = lo + log_active = True + else: + norm = mcolors.Normalize(vmin=vmin, vmax=vmax) + else: + norm = mcolors.Normalize(vmin=vmin, vmax=vmax) + self._last_vmin = vmin + self._last_vmax = vmax + # Transpose + flip so y=0 appears at the bottom of the display. + rgba = cmap(norm(ds.T[::-1])) # (h_plot, w_map, 4) + + lines = [] + for row in range(h_plot): + line = RichText() + # Heatmap cells — one terminal character per data point. + for col in range(w_map): + r = int(rgba[row, col, 0] * 255) + g = int(rgba[row, col, 1] * 255) + b = int(rgba[row, col, 2] * 255) + bg = RichColor.from_rgb(r, g, b) + if (row, col) in bubble_cells: + line.append("○", style=Style(bgcolor=bg, color="white", bold=True)) + else: + line.append(" ", style=Style(bgcolor=bg)) + # Gap + line.append(" " * _CB_GAP) + # Colorbar gradient strip (t=1 at top = vmax, t=0 at bottom = vmin) + t = 1.0 - row / max(h_plot - 1, 1) + cb = cmap(t) + cr, cg, cbb = int(cb[0] * 255), int(cb[1] * 255), int(cb[2] * 255) + for _ in range(_CB_W): + line.append(" ", style=Style(bgcolor=RichColor.from_rgb(cr, cg, cbb))) + # Value labels at top, middle, bottom + if row == 0: + lbl = f" {vmax:.3g}" + elif row == h_plot - 1: + lbl = f" {vmin:.3g}" + elif row == h_plot // 2: + mid = np.sqrt(vmin * vmax) if (log_active and vmin > 0) else (vmin + vmax) / 2 + lbl = f" {mid:.3g}" + else: + lbl = "" + line.append(lbl.ljust(_CB_LBL)[:_CB_LBL]) + lines.append(line) + + y_cc = self._y_cc if self._y_cc is not None else np.array([0.0, 1.0]) + log_tag = " [log]" if log_active else (" [log n/a]" if self._log_scale else "") + frozen_tag = " [frozen]" if self._vmin is not None else "" + header = RichText( + f" {self._varname} (step {self._step})" + f" [{vmin:.3g}, {vmax:.3g}]{log_tag}{frozen_tag}", + style="bold" + ) + footer = RichText( + f" x: [{x_cc[0]:.3f} \u2026 {x_cc[-1]:.3f}]" # pylint: disable=unsubscriptable-object + f" y: [{y_cc[0]:.3f} \u2026 {y_cc[-1]:.3f}]", + style="dim" + ) + return RichGroup(header, *lines, footer) + + +# --------------------------------------------------------------------------- +# Main TUI app +# --------------------------------------------------------------------------- + +class MFCTuiApp(App): # pylint: disable=too-many-instance-attributes + """Textual TUI for MFC post-processed data.""" + + CSS = """ + Screen { + layers: base; + } + + #content { + height: 1fr; + layout: horizontal; + } + + #sidebar { + width: 22; + border-right: solid $accent; + padding: 0 1; + } + + #var-title { + text-style: bold; + color: $accent; + padding: 0 0 1 0; + } + + #var-list { + height: 1fr; + } + + #status { + dock: bottom; + height: 1; + background: $panel; + color: $text-muted; + padding: 0 1; + } + """ + + BINDINGS = [ + Binding("q", "quit", "Quit"), + Binding("comma", "prev_step", "◀ step"), + Binding("period", "next_step", "step ▶"), + Binding("left", "prev_step", "◀ step", show=False), + Binding("right", "next_step", "step ▶", show=False), + Binding("space", "toggle_play", "▶/⏸"), + Binding("c", "cycle_cmap", "cmap"), + Binding("l", "toggle_log", "log"), + Binding("f", "toggle_freeze", "freeze"), + ] + + step_idx: reactive[int] = reactive(0, always_update=True) + var_name: reactive[str] = reactive("", always_update=True) + cmap_name: reactive[str] = reactive(_CMAPS[0], always_update=True) + log_scale: reactive[bool] = reactive(False, always_update=True) + playing: reactive[bool] = reactive(False, always_update=True) + + def __init__( # pylint: disable=too-many-arguments,too-many-positional-arguments + self, + steps: List[int], + varnames: List[str], + read_func: Callable, + ndim: int, + init_var: Optional[str] = None, + bubble_func: Optional[Callable] = None, + **kwargs, + ): + super().__init__(**kwargs) + self._steps = steps + self._varnames = varnames + self._read_func = read_func + self._ndim = ndim + self._bubble_func = bubble_func + # Store init_var but don't set the reactive yet — the DOM doesn't exist + # until on_mount, and the watcher calls query_one which needs the DOM. + self._init_var = init_var or (varnames[0] if varnames else "") + self._frozen_range: Optional[tuple] = None + self._play_timer = None + + def compose(self) -> ComposeResult: + yield Header(show_clock=False) + with Horizontal(id="content"): + with Vertical(id="sidebar"): + yield Label("Variables", id="var-title") + yield ListView( + *[ListItem(Label(v), id=f"var-{v}") for v in self._varnames], + id="var-list", + ) + yield MFCPlot(id="plot") + yield Static(self._status_text(), id="status") + yield Footer() + + def on_mount(self) -> None: + # DOM is ready — now safe to set the reactive (fires watcher → _push_data) + self.var_name = self._init_var + # Highlight the initial variable in the sidebar list + lv = self.query_one("#var-list", ListView) + for i, v in enumerate(self._varnames): + if v == self.var_name: + lv.index = i + break + + # ------------------------------------------------------------------ + # Reactive watchers + # ------------------------------------------------------------------ + + def watch_step_idx(self, _old: int, _new: int) -> None: + self._push_data() + + def watch_var_name(self, _old: str, _new: str) -> None: + self._push_data() + + def watch_cmap_name(self, _old: str, _new: str) -> None: + self._push_data() + + def watch_log_scale(self, _old: bool, _new: bool) -> None: + self._push_data() + + def watch_playing(self, _old: bool, new: bool) -> None: + if new: + self._play_timer = self.set_interval(0.5, self._auto_advance) + else: + if self._play_timer is not None: + self._play_timer.stop() + self._play_timer = None + + # ------------------------------------------------------------------ + # Helpers + # ------------------------------------------------------------------ + + def _status_text(self) -> str: + step = self._steps[self.step_idx] if self._steps else 0 + total = len(self._steps) + flags = [] + if self.log_scale: + flags.append("log") + if self._frozen_range is not None: + flags.append("frozen") + if self.playing: + flags.append("▶") + flag_str = (" " + " ".join(flags)) if flags else "" + return ( + f" step {step} [{self.step_idx + 1}/{total}]" + f" var: {self.var_name}" + f" cmap: {self.cmap_name}" + f"{flag_str}" + ) + + def _push_data(self) -> None: + """Load the current step/var and push data into the plot widget.""" + if not self._steps or not self.var_name: + return + step = self._steps[self.step_idx] + try: + assembled = _load(step, self._read_func) + except Exception as exc: # pylint: disable=broad-except + self.query_one("#status", Static).update( + f" [red]Error loading step {step}: {exc}[/red]" + ) + return + + data = assembled.variables.get(self.var_name) + plot = self.query_one("#plot", MFCPlot) + plot._x_cc = assembled.x_cc # pylint: disable=protected-access + plot._y_cc = assembled.y_cc # pylint: disable=protected-access + plot._data = data # pylint: disable=protected-access + plot._ndim = self._ndim # pylint: disable=protected-access + plot._varname = self.var_name # pylint: disable=protected-access + plot._step = step # pylint: disable=protected-access + plot._cmap_name = self.cmap_name # pylint: disable=protected-access + plot._log_scale = self.log_scale # pylint: disable=protected-access + if self._bubble_func is not None and self._ndim == 2: + try: + plot._bubbles = self._bubble_func(step) # pylint: disable=protected-access + except Exception: # pylint: disable=broad-except + plot._bubbles = None # pylint: disable=protected-access + else: + plot._bubbles = None # pylint: disable=protected-access + if self._frozen_range is not None: + plot._vmin, plot._vmax = self._frozen_range # pylint: disable=protected-access + else: + plot._vmin = None # pylint: disable=protected-access + plot._vmax = None # pylint: disable=protected-access + plot.refresh() + + self.query_one("#status", Static).update(self._status_text()) + + # ------------------------------------------------------------------ + # Actions + # ------------------------------------------------------------------ + + @on(ListView.Selected, "#var-list") + def on_var_selected(self, event: ListView.Selected) -> None: + item_id = event.item.id or "" + if item_id.startswith("var-"): + self.var_name = item_id[4:] + + def action_prev_step(self) -> None: + if self.step_idx > 0: + self.step_idx -= 1 + + def action_next_step(self) -> None: + if self.step_idx < len(self._steps) - 1: + self.step_idx += 1 + + def action_cycle_cmap(self) -> None: + idx = (_CMAPS.index(self.cmap_name) + 1) % len(_CMAPS) + self.cmap_name = _CMAPS[idx] + + def action_toggle_log(self) -> None: + self.log_scale = not self.log_scale + + def action_toggle_freeze(self) -> None: + if self._frozen_range is not None: + self._frozen_range = None + else: + plot = self.query_one("#plot", MFCPlot) + self._frozen_range = (plot._last_vmin, plot._last_vmax) # pylint: disable=protected-access + self._push_data() + + def action_toggle_play(self) -> None: + self.playing = not self.playing + + def _auto_advance(self) -> None: + self.step_idx = (self.step_idx + 1) % len(self._steps) + + +# --------------------------------------------------------------------------- +# Public entry point +# --------------------------------------------------------------------------- + +def run_tui( + init_var: Optional[str], + steps: List[int], + read_func: Callable, + ndim: int, + bubble_func: Optional[Callable] = None, +) -> None: + """Launch the Textual TUI for MFC visualization (1D/2D only).""" + if ndim not in (1, 2): + raise MFCException( + f"--tui only supports 1D and 2D data (got ndim={ndim}). " + "Use --interactive for 3D data." + ) + + # Preload first step to discover variables + first = _load(steps[0], read_func) + varnames = sorted(first.variables.keys()) + if not varnames: + raise MFCException("No variables found in data") + if init_var not in varnames: + init_var = varnames[0] + + cons.print( + f"[bold]Launching TUI[/bold] — {len(steps)} step(s), " + f"{len(varnames)} variable(s)" + ) + cons.print("[dim] ,/. or ←/→ prev/next step • space play • l log • f freeze • ↑↓ variable • q quit[/dim]") + + _step_cache.seed(steps[0], first) + + app = MFCTuiApp( + steps=steps, + varnames=varnames, + read_func=read_func, + ndim=ndim, + init_var=init_var, + bubble_func=bubble_func, + ) + app.run() diff --git a/toolchain/mfc/viz/viz.py b/toolchain/mfc/viz/viz.py new file mode 100644 index 0000000000..ed1c511e5e --- /dev/null +++ b/toolchain/mfc/viz/viz.py @@ -0,0 +1,494 @@ +""" +Main entry point for the ``./mfc.sh viz`` command. + +Dispatches to reader + renderer based on CLI arguments. +""" + +import os +import importlib +import importlib.util +import shutil +import subprocess +import sys + +from mfc.state import ARG +from mfc.common import MFC_ROOT_DIR, MFCException +from mfc.printer import cons + + +def _ensure_viz_deps() -> None: + """Install the [viz] optional extras on first use. + + Checks one sentinel per feature group so that a user who has matplotlib + pre-installed (e.g. via Anaconda) but lacks imageio, textual, or h5py + still triggers the install. + """ + _SENTINELS = ("matplotlib", "imageio", "h5py", "textual", "dash", "plotext", "plotly") + if all(importlib.util.find_spec(p) is not None for p in _SENTINELS): + return # all present + + toolchain_path = os.path.join(MFC_ROOT_DIR, "toolchain") + cons.print("[bold]Installing viz dependencies[/bold] " + "(first run — this may take a minute)...") + + env = {**os.environ, "UV_LINK_MODE": "copy"} + if shutil.which("uv"): + cmd = ["uv", "pip", "install", f"{toolchain_path}[viz]"] + else: + cmd = [sys.executable, "-m", "pip", "install", f"{toolchain_path}[viz]"] + + try: + result = subprocess.run(cmd, env=env, check=False, timeout=300) + except subprocess.TimeoutExpired as exc: + raise MFCException( + "Timed out installing viz dependencies (network may be restricted). " + f"Try manually: pip install '{toolchain_path}[viz]'" + ) from exc + if result.returncode != 0: + raise MFCException( + "Failed to install viz dependencies. " + f"Try manually: pip install '{toolchain_path}[viz]'" + ) + + importlib.invalidate_caches() + cons.print("[bold green]Viz dependencies installed.[/bold green]\n") + + +_CMAP_POPULAR = ( + 'viridis, plasma, inferno, magma, turbo, ' + 'coolwarm, RdBu_r, bwr, hot, jet, gray, seismic' +) + + +def _validate_cmap(name): + """Raise a helpful MFCException if *name* is not a known matplotlib colormap.""" + import matplotlib # pylint: disable=import-outside-toplevel + if name in matplotlib.colormaps: + return + try: + from rapidfuzz import process # pylint: disable=import-outside-toplevel + matches = process.extract(name, list(matplotlib.colormaps), limit=3) + suggestions = ', '.join(m[0] for m in matches) + hint = f" Did you mean: {suggestions}?" + except ImportError: + hint = f" Popular choices: {_CMAP_POPULAR}." + raise MFCException(f"Unknown colormap '{name}'.{hint}") + + +def _steps_hint(steps, n=8): + """Short inline preview of available steps for error messages.""" + if not steps: + return "none found" + if len(steps) <= n: + return ', '.join(str(s) for s in steps) + half = n // 2 + head = ', '.join(str(s) for s in steps[:half]) + tail = ', '.join(str(s) for s in steps[-half:]) + return f"{head}, ... [{len(steps)} total] ..., {tail}" + + +def _parse_step_list(s, available_steps): + """ + Parse a comma-separated step list, with optional '...' ellipsis expansion. + + Examples: + "0,100,200,1000" -> [0, 100, 200, 1000] (intersection with available) + "0,100,200,...,1000" -> range(0, 1001, 100) (infers stride=100 from last pair) + """ + parts = [p.strip() for p in s.split(',')] + avail_set = set(available_steps) + + if '...' in parts: + idx = parts.index('...') + if idx < 2: + raise MFCException( + f"Invalid --step value '{s}'. " + "Ellipsis '...' requires at least two values before it, " + "e.g. '0,100,...,1000'." + ) + if idx != len(parts) - 2: + raise MFCException( + f"Invalid --step value '{s}'. " + "Ellipsis '...' must be the second-to-last item, " + "e.g. '0,100,...,1000'." + ) + try: + prefix = [int(p) for p in parts[:idx]] + end = int(parts[idx + 1]) + except ValueError as exc: + raise MFCException( + f"Invalid --step value '{s}': all values must be integers." + ) from exc + + stride = prefix[-1] - prefix[-2] + if stride <= 0: + raise MFCException( + f"Invalid --step value '{s}': " + f"ellipsis stride must be positive (got {stride})." + ) + requested = list(range(prefix[0], end + 1, stride)) + else: + try: + requested = [int(p) for p in parts] + except ValueError as exc: + raise MFCException( + f"Invalid --step value '{s}': all values must be integers." + ) from exc + + return [t for t in requested if t in avail_set] + + +def _parse_steps(step_arg, available_steps): + """ + Parse the --step argument into a list of timestep integers. + + Formats: + - Single int: "1000" + - Range: "0:10000:500" (start:end:stride) + - Comma list: "0,100,200,1000" + - Ellipsis list: "0,100,200,...,1000" (stride inferred from last pair) + - "last": last available timestep + - "all": all available timesteps + """ + if step_arg is None or step_arg == 'all': + return available_steps + + if step_arg == 'last': + return [available_steps[-1]] if available_steps else [] + + s = str(step_arg) + + if ',' in s: + return _parse_step_list(s, available_steps) + + try: + if ':' in s: + parts = s.split(':') + start = int(parts[0]) + end = int(parts[1]) + stride = int(parts[2]) if len(parts) > 2 else 1 + requested = list(range(start, end + 1, stride)) + return [t for t in requested if t in set(available_steps)] + + single = int(s) + except ValueError as exc: + raise MFCException( + f"Invalid --step value '{step_arg}'. " + "Expected an integer, a range (start:end:stride), " + "a comma list (0,100,200), an ellipsis list (0,100,...,1000), " + "'last', or 'all'." + ) from exc + + if available_steps and single not in set(available_steps): + return [] + return [single] + + +def viz(): # pylint: disable=too-many-locals,too-many-statements,too-many-branches + """Main viz command dispatcher.""" + _ensure_viz_deps() + + from .reader import discover_format, discover_timesteps, assemble, has_lag_bubble_evol, read_lag_bubbles_at_step # pylint: disable=import-outside-toplevel + from .renderer import render_1d, render_1d_tiled, render_2d, render_2d_tiled, render_3d_slice, render_mp4 # pylint: disable=import-outside-toplevel + + case_dir = ARG('input') + if case_dir is None: + raise MFCException("Please specify a case directory.") + + # Resolve case directory + if not os.path.isdir(case_dir): + raise MFCException(f"Directory not found: {case_dir}") + + # Auto-detect or use specified format + fmt_arg = ARG('format') + if fmt_arg: + if fmt_arg not in ('binary', 'silo'): + raise MFCException(f"Unknown format '{fmt_arg}'. " + "Supported formats: 'binary', 'silo'.") + fmt = fmt_arg + else: + try: + fmt = discover_format(case_dir) + except FileNotFoundError as exc: + msg = str(exc) + if os.path.isfile(os.path.join(case_dir, 'case.py')): + msg += (" This looks like an MFC case directory. " + "Did you forget to run post_process first?") + raise MFCException(msg) from exc + + cons.print(f"[bold]Format:[/bold] {fmt}") + + # Handle --list-steps + if ARG('list_steps'): + steps = discover_timesteps(case_dir, fmt) + if not steps: + cons.print("[yellow]No timesteps found.[/yellow]") + else: + cons.print(f"[bold]Available timesteps ({len(steps)}):[/bold]") + # Print in columns + line = "" + for i, s in enumerate(steps): + line += f"{s:>8}" + if (i + 1) % 10 == 0: + cons.print(line) + line = "" + if line: + cons.print(line) + return + + # Handle --list-vars (requires --step) + if ARG('list_vars'): + step_arg = ARG('step') + steps = discover_timesteps(case_dir, fmt) + if not steps: + raise MFCException( + f"No timesteps found in '{case_dir}' ({fmt} format). " + "Ensure post_process has been run and produced output files.") + + if step_arg is None or step_arg == 'all': + step = steps[0] + cons.print(f"[dim]Using first available timestep: {step}[/dim]") + elif step_arg == 'last': + step = steps[-1] + cons.print(f"[dim]Using last available timestep: {step}[/dim]") + else: + try: + step = int(step_arg) + except ValueError as exc: + raise MFCException(f"Invalid --step value '{step_arg}'. " + "Expected an integer, 'last', or 'all'.") from exc + if step not in steps: + raise MFCException( + f"Timestep {step} not found. " + f"Available steps: {_steps_hint(steps)}") + + if fmt == 'silo': + from .silo_reader import assemble_silo # pylint: disable=import-outside-toplevel + assembled = assemble_silo(case_dir, step) + else: + assembled = assemble(case_dir, step, fmt) + + varnames = sorted(assembled.variables.keys()) + cons.print(f"[bold]Available variables ({len(varnames)}):[/bold]") + for vn in varnames: + data = assembled.variables[vn] + cons.print(f" {vn:<20s} min={data.min():.6g} max={data.max():.6g}") + return + + # For rendering, --step is required; --var is optional for 1D/2D (shows all in tiled layout) + varname = ARG('var') + step_arg = ARG('step') + tiled = varname is None or varname == 'all' + + if ARG('interactive') or ARG('tui') or ARG('mp4'): + # Load all steps by default; honour an explicit --step so users can + # reduce the set for large 3D cases before hitting the step limit. + if step_arg == 'last': + step_arg = 'all' + + steps = discover_timesteps(case_dir, fmt) + if not steps: + raise MFCException( + f"No timesteps found in '{case_dir}' ({fmt} format). " + "Ensure post_process has been run and produced output files.") + + requested_steps = _parse_steps(step_arg, steps) + if not requested_steps: + raise MFCException( + f"No matching timesteps for --step {step_arg!r}. " + f"Available steps: {_steps_hint(steps)}") + + # Collect rendering options + render_opts = { + 'cmap': ARG('cmap'), + 'dpi': ARG('dpi'), + 'slice_axis': ARG('slice_axis'), + } + if ARG('vmin') is not None: + render_opts['vmin'] = float(ARG('vmin')) + if ARG('vmax') is not None: + render_opts['vmax'] = float(ARG('vmax')) + if ARG('log_scale'): + render_opts['log_scale'] = True + if ARG('slice_index') is not None: + render_opts['slice_index'] = int(ARG('slice_index')) + if ARG('slice_value') is not None: + render_opts['slice_value'] = float(ARG('slice_value')) + + interactive = ARG('interactive') + + # Lagrange bubble overlay: auto-detect D/lag_bubble_evol_*.dat files + bubble_func = None + if has_lag_bubble_evol(case_dir): + bubble_func = lambda s: read_lag_bubbles_at_step(case_dir, s) # noqa: E731 + cons.print("[dim]Lagrange bubble positions detected — overlaying on plots.[/dim]") + + # Load all variables when tiled, interactive, or TUI; filter otherwise. + # TUI needs all vars loaded so the sidebar can switch between them. + load_all = tiled or interactive or ARG('tui') + + def read_step(step): + if fmt == 'silo': + from .silo_reader import assemble_silo # pylint: disable=import-outside-toplevel + return assemble_silo(case_dir, step, var=None if load_all else varname) + return assemble(case_dir, step, fmt, var=None if load_all else varname) + + # Validate variable name / discover available variables + test_assembled = read_step(requested_steps[0]) + avail = sorted(test_assembled.variables.keys()) + + # Guard against loading too many 3D timesteps (memory). + # Interactive mode caches all steps simultaneously, so use a tighter limit. + _3d_limit = 50 if interactive else 500 + if test_assembled.ndim == 3 and len(requested_steps) > _3d_limit: + raise MFCException( + f"Refusing to load {len(requested_steps)} timesteps for 3D data " + f"(limit is {_3d_limit}). Use --step with a range or stride to reduce.") + + # Tiled mode works for 1D and 2D. For 3D, auto-select the first variable. + if tiled and not interactive and not ARG('tui'): + if test_assembled.ndim == 3: + varname = avail[0] if avail else None + if varname is None: + raise MFCException("No variables found in output.") + tiled = False + cons.print(f"[dim]Auto-selected variable: [bold]{varname}[/bold]" + " (use --var to specify)[/dim]") + + if not tiled and not interactive and not ARG('tui') and varname not in test_assembled.variables: + # test_assembled was loaded with var_filter=varname so its variables dict + # may be empty. Re-read without filter (errors only, so extra I/O is fine) + # to build a useful "available variables" list for the error message. + if not avail: + if fmt == 'silo': + from .silo_reader import assemble_silo # pylint: disable=import-outside-toplevel + _full = assemble_silo(case_dir, requested_steps[0]) + else: + _full = assemble(case_dir, requested_steps[0], fmt) + avail = sorted(_full.variables.keys()) + avail_str = ', '.join(avail) if avail else '(none — check post_process output)' + raise MFCException( + f"Variable '{varname}' not found. " + f"Available: {avail_str}. " + f"Use --list-vars to see variables at a given step." + ) + + # TUI mode — launch Textual terminal UI (1D/2D only) + if ARG('tui'): + if test_assembled.ndim == 3: + raise MFCException( + "--tui only supports 1D and 2D data. " + "Use --interactive for 3D data." + ) + from .tui import run_tui # pylint: disable=import-outside-toplevel + init_var = varname if varname in avail else (avail[0] if avail else None) + run_tui(init_var, requested_steps, read_step, ndim=test_assembled.ndim, + bubble_func=bubble_func) + return + + # Interactive mode — launch Dash web server + if interactive: + from .interactive import run_interactive # pylint: disable=import-outside-toplevel + port = ARG('port') + host = ARG('host') + # Default to first available variable if --var was not specified + init_var = varname if varname in avail else (avail[0] if avail else None) + run_interactive(init_var, requested_steps, read_step, + port=int(port), host=str(host), + bubble_func=bubble_func) + return + + # Validate colormap before any rendering + cmap_name = ARG('cmap') + if cmap_name: + _validate_cmap(cmap_name) + + # Create output directory + output_base = ARG('output') + if output_base is None: + output_base = os.path.join(case_dir, 'viz') + os.makedirs(output_base, exist_ok=True) + + # MP4 mode + if ARG('mp4'): + fps = float(ARG('fps')) + _FPS_DEFAULT = 10.0 + _MIN_DURATION = 5.0 # seconds + n_frames = len(requested_steps) + if fps == _FPS_DEFAULT and n_frames / fps < _MIN_DURATION: + fps = max(1.0, n_frames / _MIN_DURATION) + cons.print( + f"[dim]Auto-adjusted fps to {fps:.2g} " + f"so video is at least {_MIN_DURATION:.0f}s " + f"(use --fps to override).[/dim]" + ) + label = 'all' if tiled else varname + mp4_path = os.path.join(output_base, f'{label}.mp4') + cons.print(f"[bold]Generating MP4:[/bold] {mp4_path} ({n_frames} frames @ {fps:.2g} fps = {n_frames/fps:.1f}s)") + success = render_mp4(varname, requested_steps, mp4_path, + fps=fps, read_func=read_step, + tiled=tiled, bubble_func=bubble_func, **render_opts) + if success: + cons.print(f"[bold green]Done:[/bold green] {mp4_path}") + else: + raise MFCException(f"Failed to generate {mp4_path}. " + "Ensure imageio and imageio-ffmpeg are installed.") + return + + # Single or multiple PNG frames + try: + from tqdm import tqdm # pylint: disable=import-outside-toplevel + step_iter = tqdm(requested_steps, desc='Rendering') + except ImportError: + step_iter = requested_steps + + failures = [] + label = 'all' if tiled else varname + for step in step_iter: + try: + # Reuse the already-loaded probe data for the first step + assembled = test_assembled if step == requested_steps[0] else read_step(step) + except (FileNotFoundError, EOFError, ValueError) as exc: + cons.print(f"[yellow]Warning:[/yellow] Skipping step {step}: {exc}") + failures.append(step) + continue + + output_path = os.path.join(output_base, f'{label}_{step}.png') + + # Inject bubble positions for this step + step_opts = render_opts + if bubble_func is not None: + try: + step_opts = dict(render_opts, bubbles=bubble_func(step)) + except Exception: # pylint: disable=broad-except + pass + + if tiled and assembled.ndim == 1: + render_1d_tiled(assembled.x_cc, assembled.variables, + step, output_path, **step_opts) + elif tiled and assembled.ndim == 2: + render_2d_tiled(assembled, step, output_path, **step_opts) + elif assembled.ndim == 1: + render_1d(assembled.x_cc, assembled.variables[varname], + varname, step, output_path, **step_opts) + elif assembled.ndim == 2: + render_2d(assembled.x_cc, assembled.y_cc, + assembled.variables[varname], + varname, step, output_path, **step_opts) + elif assembled.ndim == 3: + render_3d_slice(assembled, varname, step, output_path, **step_opts) + else: + cons.print(f"[yellow]Warning:[/yellow] Unsupported ndim={assembled.ndim} " + f"for step {step}, skipping.") + failures.append(step) + continue + + if len(requested_steps) == 1: + cons.print(f"[bold green]Saved:[/bold green] {output_path}") + + rendered = len(requested_steps) - len(failures) + if failures: + cons.print(f"[yellow]Warning:[/yellow] {len(failures)} step(s) failed: " + f"{failures[:10]}{'...' if len(failures) > 10 else ''}") + if rendered > 1: + cons.print(f"[bold green]Saved {rendered} frames to:[/bold green] {output_base}/") diff --git a/toolchain/pyproject.toml b/toolchain/pyproject.toml index 53e2140290..4788c1d56e 100644 --- a/toolchain/pyproject.toml +++ b/toolchain/pyproject.toml @@ -9,7 +9,6 @@ dependencies = [ # General "rich", "wheel", - "typos", "PyYAML", "argparse", "dataclasses", @@ -33,10 +32,6 @@ dependencies = [ "numpy", "pandas", - # Plotting - "seaborn", - "matplotlib", - # Chemistry "cantera>=3.1.0", #"pyrometheus == 1.0.5", @@ -44,15 +39,33 @@ dependencies = [ # Frontier Profiling "astunparse==1.6.2", - "colorlover", - "dash>=1.12.0", "pymongo", - "tabulate", + "tabulate" +] + +[project.optional-dependencies] +viz = [ + # 2D/3D plotting (renderer, TUI) + "matplotlib", + + # Silo-HDF5 reader + "h5py", + + # MP4 export + "imageio>=2.33", + "imageio-ffmpeg>=0.5.0", + + # Terminal UI (--tui) + "plotext>=5.2.0", + "textual>=0.47.0", + "textual-plotext>=0.2.0", + + # Interactive web UI (--interactive) + "dash>=2.0", + "plotly", + + # Progress bar (PNG/MP4 batch rendering) "tqdm", - "dash-svg", - "dash-bootstrap-components", - "kaleido", - "plotille" ] [tool.hatch.metadata]