diff --git a/.github/workflows/pdoc.yml b/.github/workflows/pdoc.yml
index af23934bd7..7c4261fa8d 100644
--- a/.github/workflows/pdoc.yml
+++ b/.github/workflows/pdoc.yml
@@ -28,7 +28,7 @@ jobs:
python-version: "3.12"
- name: exclude pvpython scripts
- run: rm examples/PySDM_examples/utils/pvanim.py
+ run: grep -rl '^#!/usr/bin/env pvpython' examples/PySDM_examples | xargs rm -f
- run: pip install pdoc nbformat gitpython
- run: pip install -e .
diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml
index 52e4d9ee5a..ea7e1ca990 100644
--- a/.github/workflows/tests.yml
+++ b/.github/workflows/tests.yml
@@ -300,8 +300,10 @@ jobs:
rm /tmp/pytest/test_run_notebooks_*current
ls /tmp/pytest/test_run_notebooks_*/
ls /tmp/pytest/test_run_notebooks_*/output/
- ffmpeg -i /tmp/pytest/test_run_notebooks_*/output/docs_intro_animation.ogv -loop 0 -vf "fps=10,split[s0][s1];[s0]palettegen[p];[s1][p]paletteuse" $HOME/work/_temp/_github_home/figures/docs_intro_animation_${{ matrix.platform }}.gif
mv /tmp/pytest/test_run_notebooks_*/output/last_animation_frame.pdf $HOME/work/_temp/_github_home/figures/last_animation_frame_${{ matrix.platform }}.pdf
+ for anim_name in docs_intro_animation parcel_animation; do
+ ffmpeg -i /tmp/pytest/test_run_notebooks_*/output/${anim_name}.ogv -loop 0 -vf "fps=10,split[s0][s1];[s0]palettegen[p];[s1][p]paletteuse" $HOME/work/_temp/_github_home/figures/${anim_name}_${{ matrix.platform }}.gif
+ done
ls $HOME/work/_temp/_github_home/figures
- name: animation movie upload
diff --git a/PySDM/exporters/parcel_vtk_exporter.py b/PySDM/exporters/parcel_vtk_exporter.py
new file mode 100644
index 0000000000..96557b45d5
--- /dev/null
+++ b/PySDM/exporters/parcel_vtk_exporter.py
@@ -0,0 +1,148 @@
+"""
+VTK Exporter for Pyrcel PySDM simulations.
+
+This module defines `VTKExporterPyrcel`, a subclass of `PySDM.exporters.VTKExporter`,
+that writes simulation outputs to VTK format using `pyevtk`. It exports product
+profiles (e.g., relative humidity) as unstructured grids and particle attributes
+as point clouds, along with `.pvd` collection files for time-series visualization
+in ParaView.
+"""
+
+from pyevtk.hl import unstructuredGridToVTK, pointsToVTK
+from pyevtk.vtk import VtkHexahedron, VtkGroup
+import numpy as np
+
+from PySDM.exporters import VTKExporter
+
+
+class VTKExporterPyrcel(VTKExporter):
+ """
+ Custom VTK exporter for Pyrcel PySDM, exporting products as grids
+ and attributes as point clouds for ParaView visualization.
+ """
+
+ def __init__(self, n_sd, output, mass_of_dry_air):
+ super().__init__()
+ self.output = output
+ self.x_coords = np.random.random(n_sd)
+ self.y_coords = np.random.random(n_sd)
+ self.z_coords = np.random.random(n_sd)
+ self.half_diagonal = []
+ self.n_levels = len(self.output["products"]["z"])
+
+ volume = mass_of_dry_air / output["products"]["rhod"][0]
+ delta_z = output["products"]["z"][1] - output["products"]["z"][0]
+ for level in range(self.n_levels):
+ if level > 0:
+ prev_to_curr_density_ratio = (
+ output["products"]["rhod"][level - 1]
+ / output["products"]["rhod"][level]
+ )
+ volume *= prev_to_curr_density_ratio
+ delta_z = (
+ output["products"]["z"][level] - output["products"]["z"][level - 1]
+ )
+ area = volume / delta_z
+ self.half_diagonal.append((2 * area) ** 0.5)
+
+ def write_pvd(self):
+ pvd_attributes = VtkGroup(self.attributes_file_path)
+ for key, value in self.exported_times["attributes"].items():
+ pvd_attributes.addFile(key + ".vtu", sim_time=value)
+ pvd_attributes.save()
+
+ pvd_products = VtkGroup(self.products_file_path)
+ for key, value in self.exported_times["products"].items():
+ pvd_products.addFile(key + ".vtu", sim_time=value)
+ pvd_products.save()
+
+ def export_products(
+ self, step, simulation
+ ): # pylint: disable=arguments-differ, too-many-locals
+ path = self.products_file_path + "_num" + self.add_leading_zeros(step)
+ self.exported_times["products"][path] = step * simulation.particulator.dt
+
+ n_vertices = 4 * self.n_levels
+ x_vertices = np.zeros(n_vertices)
+ y_vertices = np.zeros(n_vertices)
+ z_vertices = np.zeros(n_vertices)
+ connectivity = [0, 1, 2, 3]
+ cell_type = np.zeros(self.n_levels - 1)
+ cell_type[:] = VtkHexahedron.tid
+
+ for level in range(self.n_levels):
+ half_diag = self.half_diagonal[level]
+ current_z = self.output["products"]["z"][level]
+ idx = level * 4
+ x_vertices[idx], y_vertices[idx], z_vertices[idx] = (
+ -half_diag,
+ -half_diag,
+ current_z,
+ )
+ idx += 1
+ x_vertices[idx], y_vertices[idx], z_vertices[idx] = (
+ -half_diag,
+ half_diag,
+ current_z,
+ )
+ idx += 1
+ x_vertices[idx], y_vertices[idx], z_vertices[idx] = (
+ half_diag,
+ half_diag,
+ current_z,
+ )
+ idx += 1
+ x_vertices[idx], y_vertices[idx], z_vertices[idx] = (
+ half_diag,
+ -half_diag,
+ current_z,
+ )
+ connectivity += [*range(4 * (level + 1), 4 * (level + 2))] * 2
+ connectivity = np.asarray(connectivity[:-4])
+ offset = np.asarray(range(8, 8 * self.n_levels, 8))
+
+ _ = {"test_pd": np.array([44] * n_vertices)}
+
+ relative_humidity = self.output["products"]["S_max_percent"]
+ cell_data = {
+ "RH": np.full(shape=(len(relative_humidity) - 1,), fill_value=np.nan)
+ }
+ cell_data["RH"][:step] = (
+ np.array(relative_humidity[:-1] + np.diff(relative_humidity) / 2)
+ )[:step]
+ unstructuredGridToVTK(
+ path,
+ x_vertices,
+ y_vertices,
+ z_vertices,
+ connectivity=connectivity,
+ offsets=offset,
+ cell_types=cell_type,
+ cellData=cell_data,
+ )
+
+ def export_attributes(self, step, simulation): # pylint: disable=arguments-differ
+ path = self.attributes_file_path + "_num" + self.add_leading_zeros(step)
+ self.exported_times["attributes"][path] = step * simulation.particulator.dt
+
+ payload = {}
+
+ for key in self.output["attributes"].keys():
+ payload[key] = np.asarray(self.output["attributes"][key])[:, step].copy()
+ if step != 0:
+ delta_z = (
+ self.output["products"]["z"][step]
+ - self.output["products"]["z"][step - 1]
+ )
+ if step == 1:
+ self.z_coords *= delta_z
+ else:
+ self.z_coords += delta_z
+
+ pointsToVTK(
+ path,
+ 2 * (self.x_coords - 0.5) * self.half_diagonal[step],
+ 2 * (self.y_coords - 0.5) * self.half_diagonal[step],
+ self.z_coords,
+ data=payload,
+ )
diff --git a/examples/PySDM_examples/Pyrcel/simulation.py b/examples/PySDM_examples/Pyrcel/simulation.py
index aa0f92b709..a65d22d958 100644
--- a/examples/PySDM_examples/Pyrcel/simulation.py
+++ b/examples/PySDM_examples/Pyrcel/simulation.py
@@ -13,7 +13,13 @@
class Simulation(BasicSimulation):
def __init__(
- self, settings, products=None, scipy_solver=False, rtol_thd=1e-10, rtol_x=1e-10
+ self,
+ settings,
+ products=None,
+ scipy_solver=False,
+ rtol_thd=1e-10,
+ rtol_x=1e-10,
+ mass_of_dry_air=44 * si.kg,
):
n_sd = sum(settings.n_sd_per_mode)
builder = Builder(
@@ -27,7 +33,7 @@ def __init__(
initial_water_vapour_mixing_ratio=settings.initial_vapour_mixing_ratio,
T0=settings.initial_temperature,
w=settings.vertical_velocity,
- mass_of_dry_air=44 * si.kg,
+ mass_of_dry_air=mass_of_dry_air,
),
)
builder.add_dynamic(AmbientThermodynamics())
diff --git a/examples/PySDM_examples/Strzabala_2025_BEng/__init__.py b/examples/PySDM_examples/Strzabala_2025_BEng/__init__.py
new file mode 100644
index 0000000000..224822c350
--- /dev/null
+++ b/examples/PySDM_examples/Strzabala_2025_BEng/__init__.py
@@ -0,0 +1,8 @@
+# pylint: disable=invalid-name
+"""
+ParaView visualisation of a parcel-model example based on the test case from
+[Pyrcel package docs](https://pyrcel.readthedocs.io/)
+
+paraview_parcel_model.ipynb:
+.. include:: ./paraview_parcel_model.ipynb.badges.md
+"""
diff --git a/examples/PySDM_examples/Strzabala_2025_BEng/paraview.ipynb b/examples/PySDM_examples/Strzabala_2025_BEng/paraview.ipynb
new file mode 100644
index 0000000000..2ec8f17e03
--- /dev/null
+++ b/examples/PySDM_examples/Strzabala_2025_BEng/paraview.ipynb
@@ -0,0 +1,164 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "id": "6dcbcc6d",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "\"\"\"Run Pyrcel simulation and export Paraview files.\"\"\"\n",
+ "\n",
+ "import numpy as np\n",
+ "import subprocess\n",
+ "import platform\n",
+ "import pathlib\n",
+ "from PySDM import Formulae\n",
+ "from PySDM.physics import si\n",
+ "from PySDM.initialisation.spectra import Lognormal\n",
+ "from PySDM.products import (\n",
+ " ParcelDisplacement,\n",
+ " AmbientTemperature,\n",
+ " AmbientDryAirDensity,\n",
+ " AmbientRelativeHumidity,\n",
+ " ParticleSizeSpectrumPerVolume,\n",
+ " ParticleVolumeVersusRadiusLogarithmSpectrum,\n",
+ ")\n",
+ "\n",
+ "from PySDM_examples.Pyrcel import Settings, Simulation\n",
+ "\n",
+ "from PySDM.exporters.parcel_vtk_exporter import VTKExporterPyrcel"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "cd6288b5",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "settings = Settings(\n",
+ " dz=10 * si.m,\n",
+ " n_sd_per_mode=(25, 25),\n",
+ " aerosol_modes_by_kappa={\n",
+ " 0.54: Lognormal(\n",
+ " norm_factor=850 / si.cm**3,\n",
+ " m_mode=15 * si.nm,\n",
+ " s_geom=1.6,\n",
+ " ),\n",
+ " 1.2: Lognormal(\n",
+ " norm_factor=10 / si.cm**3,\n",
+ " m_mode=850 * si.nm,\n",
+ " s_geom=1.2,\n",
+ " ),\n",
+ " },\n",
+ " vertical_velocity=1.0 * si.m / si.s,\n",
+ " initial_pressure=775 * si.mbar,\n",
+ " initial_temperature=274 * si.K,\n",
+ " initial_relative_humidity=0.90,\n",
+ " displacement=1000 * si.m,\n",
+ " formulae=Formulae(constants={\"MAC\": 0.3}),\n",
+ ")\n",
+ "\n",
+ "dry_radius_bin_edges = np.logspace(\n",
+ " np.log10(1e-3 * si.um), np.log10(5e0 * si.um), 33, endpoint=False\n",
+ ")\n",
+ "\n",
+ "simulation = Simulation(\n",
+ " settings,\n",
+ " products=(\n",
+ " ParcelDisplacement(name=\"z\"),\n",
+ " AmbientRelativeHumidity(name=\"S_max_percent\", unit=\"%\", var=\"RH\"),\n",
+ " AmbientTemperature(name=\"T\"),\n",
+ " ParticleSizeSpectrumPerVolume(\n",
+ " name=\"dry:dN/dR\", radius_bins_edges=dry_radius_bin_edges, dry=True\n",
+ " ),\n",
+ " ParticleVolumeVersusRadiusLogarithmSpectrum(\n",
+ " name=\"dry:dV/dlnR\", radius_bins_edges=dry_radius_bin_edges, dry=True\n",
+ " ),\n",
+ " AmbientDryAirDensity(),\n",
+ " ),\n",
+ " mass_of_dry_air = 66666 * si.kg,\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "id": "3ff50a25",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "output = simulation.run()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "id": "15d59657",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "e = VTKExporterPyrcel(n_sd=simulation.particulator.n_sd, output=output, mass_of_dry_air=simulation.particulator.environment.mass_of_dry_air)\n",
+ "for step in settings.output_steps:\n",
+ " e.export_products(step, simulation)\n",
+ " e.export_attributes(step, simulation)\n",
+ "e.write_pvd()"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "id": "4a6ea0e6",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "product = pathlib.Path(\"./output/sd_products.pvd\").absolute()\n",
+ "attributes = pathlib.Path(\"./output/sd_attributes.pvd\").absolute()\n",
+ "\n",
+ "paraview_script = pathlib.Path(\"./paraview_parcel_model.py\").absolute()\n",
+ "\n",
+ "args = [\n",
+ " \"pvpython\",\n",
+ " \"--force-offscreen-rendering\",\n",
+ " str(paraview_script),\n",
+ " \"--sd-products-pvd\",\n",
+ " str(product),\n",
+ " \"--sd-attributes-pvd\",\n",
+ " str(attributes),\n",
+ " \"--output-animation-path\",\n",
+ " str(pathlib.Path(\"./output/parcel_animation.ogv\").absolute()),\n",
+ " \"--output-screenshot-path\",\n",
+ " str(pathlib.Path(\"./output/last_frame.pdf\").absolute()),\n",
+ "]\n",
+ "result = subprocess.run(\n",
+ " args,\n",
+ " check=platform.system() != \"Windows\",\n",
+ " capture_output=True,\n",
+ " text=True,\n",
+ ")"
+ ]
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "base",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.12.3"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/examples/PySDM_examples/Strzabala_2025_BEng/paraview_parcel_model.py b/examples/PySDM_examples/Strzabala_2025_BEng/paraview_parcel_model.py
new file mode 100644
index 0000000000..953e9eab09
--- /dev/null
+++ b/examples/PySDM_examples/Strzabala_2025_BEng/paraview_parcel_model.py
@@ -0,0 +1,244 @@
+#!/usr/bin/env pvpython
+"""
+ParaView pvpython script for visualizing sd_products and sd_attributes.
+"""
+
+import argparse
+from collections import namedtuple
+import pathlib
+from paraview import simple as pvs # pylint: disable=import-error
+
+pvs._DisableFirstRenderCameraReset() # pylint: disable=protected-access
+
+
+def cli_using_argparse(argparse_parser):
+ """
+ Command line interface using argparse.
+ """
+ argparse_parser.add_argument(
+ "--sd-products-pvd",
+ dest="sd_products_pvd",
+ help="Path to sd_products.pvd",
+ )
+ argparse_parser.add_argument(
+ "--sd-attributes-pvd",
+ dest="sd_attributes_pvd",
+ help="Path to sd_attributes.pvd",
+ )
+ argparse_parser.add_argument(
+ "--output-animation-path",
+ help="Output path for the animation file.",
+ )
+ argparse_parser.add_argument(
+ "--output-screenshot-path",
+ help="Output path for the screenshot file.",
+ )
+
+ argparse_parser.add_argument(
+ "--scalarbar-title-size",
+ type=int,
+ default=30,
+ help="Font size for scalar bar titles.",
+ )
+ argparse_parser.add_argument(
+ "--scalarbar-label-size",
+ type=int,
+ default=30,
+ help="Font size for scalar bar labels.",
+ )
+ argparse_parser.add_argument(
+ "--axes-title-size",
+ type=int,
+ default=30,
+ help="Font size for axes titles.",
+ )
+ argparse_parser.add_argument(
+ "--axes-label-size",
+ type=int,
+ default=30,
+ help="Font size for axes labels.",
+ )
+
+
+parser = argparse.ArgumentParser(
+ description="ParaView pvpython script for visualizing sd_products and sd_attributes PVD files."
+)
+cli_using_argparse(parser)
+args = parser.parse_args()
+
+sd_productspvd = pvs.PVDReader(
+ registrationName="sd_products.pvd", FileName=args.sd_products_pvd
+)
+sd_attributespvd = pvs.PVDReader(
+ registrationName="sd_attributes.pvd", FileName=args.sd_attributes_pvd
+)
+
+setup = {
+ "renderView1": pvs.GetActiveViewOrCreate("RenderView"),
+ "sd_productspvdDisplay": pvs.Show(
+ sd_productspvd,
+ pvs.GetActiveViewOrCreate("RenderView"),
+ "UnstructuredGridRepresentation",
+ ),
+ "sd_attributespvdDisplay": pvs.Show(
+ sd_attributespvd,
+ pvs.GetActiveViewOrCreate("RenderView"),
+ "UnstructuredGridRepresentation",
+ ),
+ "rHlookup_table": pvs.GetColorTransferFunction("RH"),
+ "volumelookup_table": pvs.GetColorTransferFunction("volume"),
+}
+
+animation_setup = namedtuple("setup", setup.keys())(**setup)
+
+
+def configure_color_bar_and_display(
+ display, lookup_table, kind, *, anim_setup=animation_setup
+):
+ """
+ Configure color bar and display settings.
+ """
+ display.Representation = "Surface"
+ display.SetScalarBarVisibility(anim_setup.renderView1, True)
+ color_bar = pvs.GetScalarBar(lookup_table, anim_setup.renderView1)
+ color_bar.LabelColor = [0.0, 0.0, 0.0]
+ color_bar.DrawScalarBarOutline = 1
+ color_bar.ScalarBarOutlineColor = [0.0, 0.0, 0.0]
+ color_bar.TitleColor = [0.0, 0.0, 0.0]
+
+ color_bar.TitleFontSize = args.scalarbar_title_size
+ color_bar.LabelFontSize = args.scalarbar_label_size
+
+ if kind == "prod":
+ display.Opacity = 0.4
+ display.DisableLighting = 1
+ display.Diffuse = 0.76
+ lookup_table.RescaleTransferFunction(90.0, 101.0)
+ lookup_table.ApplyPreset("Black, Blue and White", True)
+ lookup_table.NanColor = [0.67, 1.0, 1.0]
+ else:
+ display.PointSize = 13.0
+ display.RenderPointsAsSpheres = 1
+ display.Interpolation = "PBR"
+ lookup_table.RescaleTransferFunction(1e-18, 1e-13)
+ lookup_table.ApplyPreset("Cold and Hot", True)
+ lookup_table.MapControlPointsToLogSpace()
+ lookup_table.UseLogScale = 1
+ lookup_table.NumberOfTableValues = 16
+ lookup_table.InvertTransferFunction()
+
+
+def configure_data_axes_grid(display, kind):
+ """
+ Configure data axes grid settings.
+ """
+ display.DataAxesGrid.GridAxesVisibility = 1
+ display.DataAxesGrid.XTitle = ""
+ display.DataAxesGrid.YTitle = ""
+ display.DataAxesGrid.XAxisUseCustomLabels = 1
+ display.DataAxesGrid.YAxisUseCustomLabels = 1
+
+ display.DataAxesGrid.XTitleFontSize = args.axes_title_size
+ display.DataAxesGrid.XLabelFontSize = args.axes_label_size
+ display.DataAxesGrid.YTitleFontSize = args.axes_title_size
+ display.DataAxesGrid.YLabelFontSize = args.axes_label_size
+
+ if kind == "prod":
+ display.DataAxesGrid.ZTitle = ""
+ display.DataAxesGrid.GridColor = [0.0, 0.0, 0.0]
+ display.DataAxesGrid.ShowGrid = 1
+ display.DataAxesGrid.LabelUniqueEdgesOnly = 0
+ display.DataAxesGrid.ZAxisUseCustomLabels = 1
+ display.DataAxesGrid.ZTitleFontSize = args.axes_title_size
+ display.DataAxesGrid.ZLabelFontSize = args.axes_label_size
+ else:
+ display.DataAxesGrid.ZTitle = "Z [m]"
+ display.DataAxesGrid.ZLabelColor = [0.0, 0.0, 0.0]
+ display.DataAxesGrid.ZTitleColor = [0.0, 0.0, 0.0]
+ display.DataAxesGrid.ZTitleFontSize = args.axes_title_size
+ display.DataAxesGrid.ZLabelFontSize = args.axes_label_size
+
+
+def configure_view_appearance(render_view):
+ """
+ Configure view appearance settings."""
+ render_view.OrientationAxesLabelColor = [0.0, 0.0, 0.0]
+ render_view.OrientationAxesOutlineColor = [0.0, 0.0, 0.0]
+ render_view.OrientationAxesXVisibility = 0
+ render_view.OrientationAxesYVisibility = 0
+ render_view.OrientationAxesZColor = [0.0, 0.0, 0.0]
+ render_view.UseColorPaletteForBackground = 0
+ render_view.Background = [1.0, 1.0, 1.0]
+
+
+def set_camera_view(render_view):
+ """
+ Set camera view settings.
+ """
+ layout1 = pvs.GetLayout()
+ layout1.SetSize(1592, 1128)
+ render_view.CameraPosition = [
+ 1548.95,
+ -1349.49,
+ 699.27,
+ ]
+ render_view.CameraFocalPoint = [
+ -1.37e-13,
+ 3.18e-13,
+ 505.00,
+ ]
+ render_view.CameraViewUp = [
+ -0.07,
+ 0.06,
+ 0.99,
+ ]
+ render_view.CameraParallelScale = 534.08
+
+
+def save_animation_and_screenshot(render_view, animation_path, screenshot_path):
+ """
+ Save animation and screenshot.
+ """
+ animation_scene = pvs.GetAnimationScene()
+ animation_scene.UpdateAnimationUsingDataTimeSteps()
+
+ pvs.SaveAnimation(animation_path, render_view, FrameRate=15)
+
+ if not screenshot_path:
+ return
+
+ time_steps = sd_productspvd.TimestepValues
+ if not time_steps:
+ return
+
+ last_time = time_steps[-1]
+ render_view.ViewTime = last_time
+
+ for reader in (sd_productspvd, sd_attributespvd):
+ reader.UpdatePipeline(last_time)
+
+ pvs.ExportView(
+ filename=str(pathlib.Path(screenshot_path)),
+ view=render_view,
+ Rasterize3Dgeometry=False,
+ GL2PSdepthsortmethod="BSP sorting (slow, best)",
+ )
+
+ pvs.RenderAllViews()
+
+
+configure_color_bar_and_display(
+ animation_setup.sd_productspvdDisplay, animation_setup.rHlookup_table, kind="prod"
+)
+configure_data_axes_grid(animation_setup.sd_productspvdDisplay, kind="prod")
+configure_view_appearance(animation_setup.renderView1)
+configure_color_bar_and_display(
+ animation_setup.sd_attributespvdDisplay,
+ animation_setup.volumelookup_table,
+ kind="attr",
+)
+configure_data_axes_grid(animation_setup.sd_attributespvdDisplay, kind="attr")
+set_camera_view(animation_setup.renderView1)
+save_animation_and_screenshot(
+ animation_setup.renderView1, args.output_animation_path, args.output_screenshot_path
+)
diff --git a/examples/docs/pysdm_examples_landing.md b/examples/docs/pysdm_examples_landing.md
index 8479131c91..f5a8109371 100644
--- a/examples/docs/pysdm_examples_landing.md
+++ b/examples/docs/pysdm_examples_landing.md
@@ -86,6 +86,7 @@ Example notebooks include:
- `PySDM_examples.Yang_et_al_2018`: polydisperse particle spectrum, activation/deactivation cycles
- `PySDM_examples.Abdul_Razzak_Ghan_2000`: polydisperse activation, comparison against GCM parameterisation
- `PySDM_examples.Pyrcel`: polydisperse activation, mimicking example test case from Pyrcel documentation
+ - `PySDM_examples.Strzabala_2025_BEng`: ParaView visualisation example
- `PySDM_examples.Lowe_et_al_2019`: externally mixed polydisperse size spectrum with surface-active organics case
- `PySDM_examples.Grabowski_and_Pawlowska_2023`: polydisperse activation, focus on ripening
- `PySDM_examples.Jensen_and_Nugent_2017`: polydisperse activation featuring giant CCN
@@ -97,6 +98,10 @@ Example notebooks include:
The parcel environment is also featured in the PySDM tutorials.
+
+ 
+
+
## 0D box environment
The box environment is void of any spatial or thermodynamic context, it constitutes the most basic framework.
diff --git a/tests/unit_tests/exporters/test_parcel_vtk_exporter.py b/tests/unit_tests/exporters/test_parcel_vtk_exporter.py
new file mode 100644
index 0000000000..156f0b460c
--- /dev/null
+++ b/tests/unit_tests/exporters/test_parcel_vtk_exporter.py
@@ -0,0 +1,205 @@
+"""Test module for VTKExporterPyrcel."""
+
+from unittest.mock import Mock
+import numpy as np
+import pytest
+
+from PySDM.exporters.parcel_vtk_exporter import VTKExporterPyrcel
+
+
+class TestVTKExporterPyrcel:
+ """Test class for VTKExporterPyrcel functionality."""
+
+ @staticmethod
+ @pytest.fixture
+ def mock_simulation():
+ """Mock simulation object with necessary attributes."""
+ simulation = Mock()
+ simulation.particulator = Mock()
+ simulation.particulator.dt = 1.0
+ simulation.particulator.n_sd = 100
+ simulation.particulator.environment = Mock()
+ simulation.particulator.environment.mass_of_dry_air = 66666
+
+ return simulation
+
+ @staticmethod
+ @pytest.fixture
+ def mock_output():
+ """Mock output data structure."""
+ n_levels = 5
+ output = {
+ "products": {
+ "z": np.linspace(0, 100, n_levels),
+ "rhod": np.ones(n_levels) * 1.2,
+ "S_max_percent": np.linspace(100, 110, n_levels),
+ },
+ "attributes": {
+ "radius": np.random.random((100, 10)),
+ "n": np.random.random((100, 10)),
+ },
+ }
+ return output
+
+ def test_vtk_exporter_pyrcel_initialization(self, mock_output):
+ """Test that VTKExporterPyrcel initializes correctly."""
+ # Arrange
+ n_sd = 100
+ mass_of_dry_air = 66666
+
+ # Act
+ exporter = VTKExporterPyrcel(n_sd, mock_output, mass_of_dry_air)
+
+ # Assert
+ assert exporter.output == mock_output
+ assert len(exporter.x_coords) == n_sd
+ assert len(exporter.y_coords) == n_sd
+ assert len(exporter.z_coords) == n_sd
+ assert len(exporter.half_diagonal) == exporter.n_levels
+ assert exporter.n_levels == len(mock_output["products"]["z"])
+
+ def test_vtk_exporter_pyrcel_export_products(
+ self, mock_output, mock_simulation, tmp_path
+ ):
+ """Test exporting products to VTK format."""
+ # Arrange
+ n_sd = 100
+ mass_of_dry_air = 66666
+
+ exporter = VTKExporterPyrcel(n_sd, mock_output, mass_of_dry_air)
+ exporter.path = str(tmp_path)
+ exporter.attributes_file_path = str(tmp_path / "sd_attributes")
+ exporter.products_file_path = str(tmp_path / "sd_products")
+ exporter.exported_times = {"products": {}, "attributes": {}}
+
+ step = 1
+
+ # Act
+ exporter.export_products(step, mock_simulation)
+
+ # Assert
+ expected_file = (
+ tmp_path / f"sd_products_num{exporter.add_leading_zeros(step)}.vtu"
+ )
+ assert expected_file.exists()
+
+ assert len(exporter.exported_times["products"]) == 1
+ expected_time = step * mock_simulation.particulator.dt
+ assert list(exporter.exported_times["products"].values())[0] == expected_time
+
+ def test_vtk_exporter_pyrcel_export_attributes(
+ self, mock_output, mock_simulation, tmp_path
+ ):
+ """Test exporting attributes to VTK format."""
+ # Arrange
+ n_sd = 100
+ mass_of_dry_air = 66666
+
+ exporter = VTKExporterPyrcel(n_sd, mock_output, mass_of_dry_air)
+ exporter.path = str(tmp_path)
+ exporter.attributes_file_path = str(tmp_path / "sd_attributes")
+ exporter.products_file_path = str(tmp_path / "sd_products")
+ exporter.exported_times = {"products": {}, "attributes": {}}
+
+ step = 1
+
+ # Act
+ exporter.export_attributes(step, mock_simulation)
+
+ # Assert
+ expected_file = (
+ tmp_path / f"sd_attributes_num{exporter.add_leading_zeros(step)}.vtu"
+ )
+ assert expected_file.exists()
+
+ assert len(exporter.exported_times["attributes"]) == 1
+ expected_time = step * mock_simulation.particulator.dt
+ assert list(exporter.exported_times["attributes"].values())[0] == expected_time
+
+ def test_vtk_exporter_pyrcel_write_pvd(
+ self, mock_output, mock_simulation, tmp_path
+ ):
+ """Test writing PVD collection files."""
+ # Arrange
+ n_sd = 100
+ mass_of_dry_air = 66666
+
+ exporter = VTKExporterPyrcel(n_sd, mock_output, mass_of_dry_air)
+ exporter.path = str(tmp_path)
+ exporter.attributes_file_path = str(tmp_path / "sd_attributes")
+ exporter.products_file_path = str(tmp_path / "sd_products")
+ exporter.exported_times = {"products": {}, "attributes": {}}
+
+ steps = [1, 2]
+ for step in steps:
+ exporter.export_products(step, mock_simulation)
+ exporter.export_attributes(step, mock_simulation)
+
+ # Act
+ exporter.write_pvd()
+
+ # Assert
+ attributes_pvd = tmp_path / "sd_attributes.pvd"
+ products_pvd = tmp_path / "sd_products.pvd"
+
+ assert attributes_pvd.exists()
+ assert products_pvd.exists()
+
+ def test_vtk_exporter_pyrcel_coordinate_calculation(
+ self, mock_output, mock_simulation, tmp_path
+ ):
+ """Test that coordinate calculations work correctly across multiple steps."""
+ # Arrange
+ exporter = VTKExporterPyrcel(100, mock_output, 666666)
+ exporter.path = str(tmp_path)
+ exporter.attributes_file_path = str(tmp_path / "sd_attributes")
+ exporter.products_file_path = str(tmp_path / "sd_products")
+ exporter.exported_times = {"products": {}, "attributes": {}}
+
+ initial_z_coords = exporter.z_coords.copy()
+
+ # Act
+ step1 = 1
+ exporter.export_attributes(step1, mock_simulation)
+ z_after_step1 = exporter.z_coords.copy()
+
+ # Act
+ step2 = 2
+ exporter.export_attributes(step2, mock_simulation)
+ z_after_step2 = exporter.z_coords.copy()
+
+ # Assert
+ if step1 != 0:
+ delta_z_step1 = (
+ mock_output["products"]["z"][step1]
+ - mock_output["products"]["z"][step1 - 1]
+ )
+ expected_z_step1 = initial_z_coords * delta_z_step1
+ np.testing.assert_array_equal(z_after_step1, expected_z_step1)
+
+ if step2 != 0:
+ delta_z_step2 = (
+ mock_output["products"]["z"][step2]
+ - mock_output["products"]["z"][step2 - 1]
+ )
+ expected_z_step2 = z_after_step1 + delta_z_step2
+ np.testing.assert_array_equal(z_after_step2, expected_z_step2)
+
+ def test_vtk_exporter_pyrcel_half_diagonal_calculation(self, mock_output):
+ """Test the half_diagonal calculation logic."""
+ # Arrange
+ n_sd = 100
+ mass_of_dry_air = 66666
+
+ # Act
+ exporter = VTKExporterPyrcel(n_sd, mock_output, mass_of_dry_air)
+
+ # Assert
+ assert len(exporter.half_diagonal) == exporter.n_levels
+ assert all(hd > 0 for hd in exporter.half_diagonal)
+
+ volume_0 = mass_of_dry_air / mock_output["products"]["rhod"][0]
+ delta_z_0 = mock_output["products"]["z"][1] - mock_output["products"]["z"][0]
+ area_0 = volume_0 / delta_z_0
+ expected_hd_0 = (2 * area_0) ** 0.5
+ assert abs(exporter.half_diagonal[0] - expected_hd_0) < 1e-10