From b97d06489c156ff0213d8b6ad425ac202559497b Mon Sep 17 00:00:00 2001 From: Sylwester Arabas Date: Fri, 1 Dec 2023 01:05:39 +0100 Subject: [PATCH 01/11] add smoke test checking error measures for the Shima_et_al_2009 box example. Closes #327 --- .../Shima_et_al_2009/spectrum_plotter.py | 32 ++++---- .../box/shima_et_al_2009/test_convergence.py | 76 +++++++++++++++++++ 2 files changed, 94 insertions(+), 14 deletions(-) create mode 100644 tests/smoke_tests/box/shima_et_al_2009/test_convergence.py diff --git a/examples/PySDM_examples/Shima_et_al_2009/spectrum_plotter.py b/examples/PySDM_examples/Shima_et_al_2009/spectrum_plotter.py index 20285f0a8a..f989cec287 100644 --- a/examples/PySDM_examples/Shima_et_al_2009/spectrum_plotter.py +++ b/examples/PySDM_examples/Shima_et_al_2009/spectrum_plotter.py @@ -81,12 +81,16 @@ def save(self, file): self.finish() pyplot.savefig(file, format=self.format) - def plot(self, spectrum, t): - error = self.plot_analytic_solution(self.settings, t, spectrum) - self.plot_data(self.settings, t, spectrum) + def plot( + self, spectrum, t, label=None, color=None, title=None, add_error_to_label=False + ): + error = self.plot_analytic_solution(self.settings, t, spectrum, title) + if label is not None and add_error_to_label: + label += f" error={error:.4g}" + self.plot_data(self.settings, t, spectrum, label, color) return error - def plot_analytic_solution(self, settings, t, spectrum=None): + def plot_analytic_solution(self, settings, t, spectrum, title): if t == 0: analytic_solution = settings.spectrum.size_distribution else: @@ -123,11 +127,13 @@ def analytic_solution(x): if spectrum is not None: y = spectrum * si.kilograms / si.grams error = error_measure(y, y_true, x) - self.title = f"error measure: {error:.2f}" # TODO #327 relative error + self.title = ( + title or f"error measure: {error:.2f}" + ) # TODO #327 relative error return error return None - def plot_data(self, settings, t, spectrum): + def plot_data(self, settings, t, spectrum, label, color): if self.smooth: scope = self.smooth_scope if t != 0: @@ -144,18 +150,16 @@ def plot_data(self, settings, t, spectrum): self.ax.plot( (x[:-1] + dx / 2) * si.metres / si.micrometres, spectrum[:-scope] * si.kilograms / si.grams, - label=f"t = {t}s", - color=self.colors( - t / (self.settings.output_steps[-1] * self.settings.dt) - ), + label=label or f"t = {t}s", + color=color + or self.colors(t / (self.settings.output_steps[-1] * self.settings.dt)), ) else: self.ax.step( settings.radius_bins_edges[:-1] * si.metres / si.micrometres, spectrum * si.kilograms / si.grams, where="post", - label=f"t = {t}s", - color=self.colors( - t / (self.settings.output_steps[-1] * self.settings.dt) - ), + label=label or f"t = {t}s", + color=color + or self.colors(t / (self.settings.output_steps[-1] * self.settings.dt)), ) diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py new file mode 100644 index 0000000000..dffa230244 --- /dev/null +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -0,0 +1,76 @@ +""" +Tests checking if, for the example case from Fig 4 in +[Shima et al. 2009](https://doi.org/10.1002/qj.441), +the simulations converge towards analytic solution. +""" +from inspect import signature +from itertools import islice + +import pytest +from matplotlib import pyplot +from PySDM_examples.Shima_et_al_2009.example import run +from PySDM_examples.Shima_et_al_2009.settings import Settings +from PySDM_examples.Shima_et_al_2009.spectrum_plotter import SpectrumPlotter + +from PySDM.physics import si + +COLORS = ("red", "green", "blue") + + +class TestConvergence: # pylint: disable=missing-class-docstring + @staticmethod + @pytest.mark.parametrize( + "adaptive, dt", + ( + pytest.param(False, 100 * si.s, marks=pytest.mark.xfail(strict=True)), + (True, 100 * si.s), + pytest.param(False, 50 * si.s, marks=pytest.mark.xfail(strict=True)), + (False, 5 * si.s), + ), + ) + def test_convergence_with_sd_count(dt, adaptive, plot=False): + """check if increasing the number of super particles indeed + reduces the error of the simulation (vs. analytic solution)""" + # arrange + settings = Settings(steps=[3600]) + settings.adaptive = adaptive + plotter = SpectrumPlotter(settings) + errors = {} + + # act + for i, ln2_nsd in enumerate((13, 15, 17)): + settings.dt = dt + settings.n_sd = 2**ln2_nsd + values, _ = run(settings) + + errors[ln2_nsd] = plotter.plot( + **dict( + islice( + { # supporting older versions of PySDM-examples + "t": settings.steps[-1], + "spectrum": values[tuple(values.keys())[-1]], + "label": f"{ln2_nsd=}", + "color": COLORS[i], + "title": f"{settings.dt=} settings.times={settings.steps} {settings.adaptive=}", + "add_error_to_label": True, + }.items(), + len(signature(plotter.plot).parameters), + ) + ) + ) + + # plot + if plot: + plotter.show() + else: + pyplot.clf() + + # assert monotonicity (i.e., the larger the sd count, the smaller the error) + print(errors) + assert tuple(errors.keys()) == tuple(sorted(errors.keys())) + assert tuple(errors.values()) == tuple(reversed(sorted(errors.values()))) + + @staticmethod + def test_convergence_with_timestep(): + """ditto for timestep""" + pytest.skip("# TODO #1189") From ef160c5192e43cd32f1b7f4e87ab3dca8750b42a Mon Sep 17 00:00:00 2001 From: Sylwester Arabas Date: Fri, 1 Dec 2023 01:18:50 +0100 Subject: [PATCH 02/11] make pylint happier --- .../smoke_tests/box/shima_et_al_2009/test_convergence.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py index dffa230244..89d430ca3c 100644 --- a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -43,6 +43,13 @@ def test_convergence_with_sd_count(dt, adaptive, plot=False): settings.n_sd = 2**ln2_nsd values, _ = run(settings) + title = ( + "" + if i != 0 + else ( + f"{settings.dt=} settings.times={settings.steps} {settings.adaptive=}" + ) + ) errors[ln2_nsd] = plotter.plot( **dict( islice( @@ -51,7 +58,7 @@ def test_convergence_with_sd_count(dt, adaptive, plot=False): "spectrum": values[tuple(values.keys())[-1]], "label": f"{ln2_nsd=}", "color": COLORS[i], - "title": f"{settings.dt=} settings.times={settings.steps} {settings.adaptive=}", + "title": title, "add_error_to_label": True, }.items(), len(signature(plotter.plot).parameters), From fa45e3f7500593693717a78c5185d175ae416617 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Thu, 29 May 2025 14:44:58 +0200 Subject: [PATCH 03/11] update bibliography --- docs/bibliography.json | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/docs/bibliography.json b/docs/bibliography.json index 69ee7458e5..5ea579922c 100644 --- a/docs/bibliography.json +++ b/docs/bibliography.json @@ -121,7 +121,8 @@ "examples/PySDM_examples/Shima_et_al_2009/__init__.py", "PySDM/dynamics/collisions/collision.py", "tutorials/collisions/collisions_playground.ipynb", - "tutorials/wikipedia/sdm.ipynb" + "tutorials/wikipedia/sdm.ipynb", + "tests/smoke_tests/box/shima_et_al_2009/test_convergence.py" ], "label": "Shima et al. 2009 (Q. J. R. Meteorol. Soc. 135)", "title": "The super-droplet method for the numerical simulation of clouds and precipitation: a particle-based and probabilistic microphysics model coupled with a non-hydrostatic model" From adb86873abeeb004a5b65b007505082bcb9e2c78 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Thu, 29 May 2025 14:55:39 +0200 Subject: [PATCH 04/11] remove plot from Berry spectrumplotter as it is in parent class --- examples/PySDM_examples/Berry_1967/spectrum_plotter.py | 4 ---- tests/smoke_tests/box/shima_et_al_2009/test_convergence.py | 1 + 2 files changed, 1 insertion(+), 4 deletions(-) diff --git a/examples/PySDM_examples/Berry_1967/spectrum_plotter.py b/examples/PySDM_examples/Berry_1967/spectrum_plotter.py index 3b92adc5dd..6d461e0528 100644 --- a/examples/PySDM_examples/Berry_1967/spectrum_plotter.py +++ b/examples/PySDM_examples/Berry_1967/spectrum_plotter.py @@ -25,7 +25,3 @@ def show(self): self.finish() self.ticks() show_plot() - - def plot(self, spectrum, t): - settings = self.settings - self.plot_data(settings, t, spectrum) diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py index 89d430ca3c..b7fd81c29e 100644 --- a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -3,6 +3,7 @@ [Shima et al. 2009](https://doi.org/10.1002/qj.441), the simulations converge towards analytic solution. """ + from inspect import signature from itertools import islice From a2196ad18ba1fa2d28989783e413eed1c5b21df0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Sat, 31 May 2025 02:24:12 +0200 Subject: [PATCH 05/11] change parameterisation and hope it is right --- tests/smoke_tests/box/shima_et_al_2009/test_convergence.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py index b7fd81c29e..08738e5939 100644 --- a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -26,7 +26,7 @@ class TestConvergence: # pylint: disable=missing-class-docstring pytest.param(False, 100 * si.s, marks=pytest.mark.xfail(strict=True)), (True, 100 * si.s), pytest.param(False, 50 * si.s, marks=pytest.mark.xfail(strict=True)), - (False, 5 * si.s), + (True, 50 * si.s), ), ) def test_convergence_with_sd_count(dt, adaptive, plot=False): @@ -74,7 +74,6 @@ def test_convergence_with_sd_count(dt, adaptive, plot=False): pyplot.clf() # assert monotonicity (i.e., the larger the sd count, the smaller the error) - print(errors) assert tuple(errors.keys()) == tuple(sorted(errors.keys())) assert tuple(errors.values()) == tuple(reversed(sorted(errors.values()))) From 70523697ddd918b80f5d67a0c025338e51ddbdc0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Sat, 31 May 2025 02:40:38 +0200 Subject: [PATCH 06/11] add plot function to Berry (deleted one) with changed signature --- examples/PySDM_examples/Berry_1967/spectrum_plotter.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/examples/PySDM_examples/Berry_1967/spectrum_plotter.py b/examples/PySDM_examples/Berry_1967/spectrum_plotter.py index 6d461e0528..032611f31e 100644 --- a/examples/PySDM_examples/Berry_1967/spectrum_plotter.py +++ b/examples/PySDM_examples/Berry_1967/spectrum_plotter.py @@ -25,3 +25,9 @@ def show(self): self.finish() self.ticks() show_plot() + + def plot( + self, spectrum, t, label=None, color=None, title=None, add_error_to_label=False + ): + settings = self.settings + self.plot_data(settings, t, spectrum, label, color) From 18e0e32b6be35e50000c720229af7718db0e8cfd Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Sun, 1 Jun 2025 01:19:16 +0200 Subject: [PATCH 07/11] change the list of considered particles number to exemplify reduced error --- tests/smoke_tests/box/shima_et_al_2009/test_convergence.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py index 08738e5939..6ca85ccf65 100644 --- a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -39,7 +39,7 @@ def test_convergence_with_sd_count(dt, adaptive, plot=False): errors = {} # act - for i, ln2_nsd in enumerate((13, 15, 17)): + for i, ln2_nsd in enumerate((11, 15, 19)): settings.dt = dt settings.n_sd = 2**ln2_nsd values, _ = run(settings) From a3f476bfbdaca1a32f5c339d6a98b40e4f818870 Mon Sep 17 00:00:00 2001 From: AgnieszkaZaba <56157996+AgnieszkaZaba@users.noreply.github.com> Date: Tue, 3 Jun 2025 13:10:06 +0200 Subject: [PATCH 08/11] bump matplotlib CI requirement --- examples/setup.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/examples/setup.py b/examples/setup.py index 5e4930198f..9f554c78ea 100644 --- a/examples/setup.py +++ b/examples/setup.py @@ -28,8 +28,7 @@ def get_long_description(): "PyMPDATA" + (">=1.0.15" if CI else ""), "open-atmos-jupyter-utils", "pystrict", - # https://github.com/matplotlib/matplotlib/issues/28551 - "matplotlib" + ("!=3.9.1" if CI else ""), + "matplotlib" + (">=3.10.3" if CI else ""), "joblib", "ipywidgets", "seaborn", From 90181157b449e210b26b537a69045b5f909683b8 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Tue, 3 Jun 2025 19:49:48 +0200 Subject: [PATCH 09/11] del axes before clf --- tests/smoke_tests/box/shima_et_al_2009/test_convergence.py | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py index 6ca85ccf65..603d2c2a40 100644 --- a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -71,6 +71,7 @@ def test_convergence_with_sd_count(dt, adaptive, plot=False): if plot: plotter.show() else: + pyplot.delaxes() pyplot.clf() # assert monotonicity (i.e., the larger the sd count, the smaller the error) From 955e1aca3db665428bb4be61b239c7468fdbf5f7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Tue, 3 Jun 2025 19:57:20 +0200 Subject: [PATCH 10/11] update setup --- examples/setup.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/setup.py b/examples/setup.py index 9f554c78ea..77f4e9f84d 100644 --- a/examples/setup.py +++ b/examples/setup.py @@ -28,7 +28,7 @@ def get_long_description(): "PyMPDATA" + (">=1.0.15" if CI else ""), "open-atmos-jupyter-utils", "pystrict", - "matplotlib" + (">=3.10.3" if CI else ""), + "matplotlib" + (">=3.9.4" if CI else ""), "joblib", "ipywidgets", "seaborn", From 7ef56b3a0e19856e15c6f4f6d7da7f262ef20c2c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Agnieszka=20=C5=BBaba?= Date: Thu, 5 Jun 2025 00:17:24 +0200 Subject: [PATCH 11/11] use previous matplotlib setup; add workaround matplotlib bug --- examples/setup.py | 3 ++- tests/smoke_tests/box/shima_et_al_2009/test_convergence.py | 3 ++- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/examples/setup.py b/examples/setup.py index 77f4e9f84d..5e4930198f 100644 --- a/examples/setup.py +++ b/examples/setup.py @@ -28,7 +28,8 @@ def get_long_description(): "PyMPDATA" + (">=1.0.15" if CI else ""), "open-atmos-jupyter-utils", "pystrict", - "matplotlib" + (">=3.9.4" if CI else ""), + # https://github.com/matplotlib/matplotlib/issues/28551 + "matplotlib" + ("!=3.9.1" if CI else ""), "joblib", "ipywidgets", "seaborn", diff --git a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py index 603d2c2a40..f16ec2d552 100644 --- a/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py +++ b/tests/smoke_tests/box/shima_et_al_2009/test_convergence.py @@ -71,7 +71,8 @@ def test_convergence_with_sd_count(dt, adaptive, plot=False): if plot: plotter.show() else: - pyplot.delaxes() + # https://github.com/matplotlib/matplotlib/issues/9970 + pyplot.xscale("linear") pyplot.clf() # assert monotonicity (i.e., the larger the sd count, the smaller the error)