From 48f5fcbcf8d11c2496dd4ed373992dc002c72e2b Mon Sep 17 00:00:00 2001 From: pawel Date: Fri, 6 Sep 2024 12:05:42 +0200 Subject: [PATCH 01/10] move jinja templates to docs folder --- .github/workflows/tests+pypi.yml | 2 +- {pdoc_templates => docs/templates}/index.html.jinja2 | 0 2 files changed, 1 insertion(+), 1 deletion(-) rename {pdoc_templates => docs/templates}/index.html.jinja2 (100%) diff --git a/.github/workflows/tests+pypi.yml b/.github/workflows/tests+pypi.yml index afa57efe..74db009a 100644 --- a/.github/workflows/tests+pypi.yml +++ b/.github/workflows/tests+pypi.yml @@ -98,7 +98,7 @@ jobs: pip3 install pdoc pip install -e . -e ./examples export PDOC_ALLOW_EXEC=1 - python -We -m pdoc -o html PyMPDATA examples/PyMPDATA_examples -t pdoc_templates + python -We -m pdoc -o html PyMPDATA examples/PyMPDATA_examples -t docs/templates - if: ${{ github.ref == 'refs/heads/main' && matrix.platform == 'ubuntu-latest' }} uses: JamesIves/github-pages-deploy-action@4.1.1 with: diff --git a/pdoc_templates/index.html.jinja2 b/docs/templates/index.html.jinja2 similarity index 100% rename from pdoc_templates/index.html.jinja2 rename to docs/templates/index.html.jinja2 From c2c28154935b5b2bd6901ab119ce9464d624898c Mon Sep 17 00:00:00 2001 From: pawel Date: Tue, 10 Sep 2024 23:58:20 +0200 Subject: [PATCH 02/10] move some info from readme to index.html --- PyMPDATA/__init__.py | 3 ++ docs/templates/index.html.jinja2 | 84 ++++++++++++++++++++++++++++++-- 2 files changed, 83 insertions(+), 4 deletions(-) diff --git a/PyMPDATA/__init__.py b/PyMPDATA/__init__.py index b737f921..2d2600b1 100644 --- a/PyMPDATA/__init__.py +++ b/PyMPDATA/__init__.py @@ -6,6 +6,9 @@ `PyMPDATA.scalar_field.ScalarField` and `PyMPDATA.vector_field.VectorField` elements: ![](https://github.com/atmos-cloud-sim-uj/PyMPDATA/releases/download/tip/readme_grid.png) + +.. include:: ../README.md + """ # pylint: disable=invalid-name diff --git a/docs/templates/index.html.jinja2 b/docs/templates/index.html.jinja2 index 20118138..9d8c9ad4 100644 --- a/docs/templates/index.html.jinja2 +++ b/docs/templates/index.html.jinja2 @@ -16,10 +16,86 @@ pympdata logo - {% filter to_html %} -# PyMPDATA documentation home page -Work in progress... - {% endfilter %} +
+

What is PyMPDATA?

+

+ PyMPDATA is a Numba-accelerated Pythonic implementation of the MPDATA algorithm + of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond for + numerically solving generalised convection-diffusion PDEs in 1D, 2D and 3D structured meshes + with coordinate transformations. +

+
+
+

What is the difference between PyMPDATA and PyMPDATA-examples?

+

+ PyMPDATA is a Python package that provides the MPDATA algorithm implementation. + It is a library that can be used in your own projects. +

+

+ PyMPDATA-examples is a Python package that provides examples of how to use PyMPDATA. + It includes common Python modules used in PyMPDATA smoke tests and in example notebooks, + but the package wheels do not include the notebooks. +

+

They exist separately on PyPI, but their development and issue tracking is in the same GitHub repository.

+

Dependencies

+

+ PyMPDATA depends only on NumPy and Numba. +

+

+ PyMPDATA-examples requires additional packages listed in install_requires in + setup.py. + Amongst them is PyMPDATA. +

+
+
+

Installation

+

+ PyMPDATA is available on PyPI and can be installed using pip: +

+ 
pip install PyMPDATA
+

Note: the way above will not install PyMPDATA-examples, to install them both you can run:

+ 
pip install PyMPDATA[tests]
+
+

+ PyMPDATA-examples is available on PyPI and can be installed using pip: +

+ 
pip install PyMPDATA-examples
+

Note: this will also install PyMPDATA

+
+
+

Contributing, reporting issues and seeking support

+

+ Submitting new code to both packages is done through the same GitHub repository via + Pull requests. +

+

+ Issues regarding any incorrect, unintuitive or undocumented behaviour of PyMPDATA or PyMPDATA-examples + are best to be reported on the + GitHub issue tracker. +

+

+ We encourage to use the GitHub Discussions + feature (rather than the issue tracker) for seeking support in understanding, + using and extending PyMPDATA code. +

+
+
+

Credits and License

+

+ Development of PyMPDATA was supported by the EU through a grant of the + Foundation for Polish Science (POIR.04.04.00-00-5E1C/18). +

+

+ licence: GPL v3 +

+

+ As of now, the copyright to the entire PyMPDATA codebase is with the Jagiellonian + University, and code contributions are assumed to imply transfer of copyright. + Should there be a need to make an exception, please indicate it when creating + a pull request or contributing code in any other way. In any case, + the license of the contributed code must be compatible with GPL v3. +

+
{% include "search.html.jinja2" %} From e5bfc2315f58f2fc52e35464ad0208138b2230c8 Mon Sep 17 00:00:00 2001 From: Pawel Magnuszewski Date: Wed, 11 Sep 2024 00:58:42 +0200 Subject: [PATCH 03/10] remove image from init --- .pre-commit-config.yaml | 4 ++-- PyMPDATA/__init__.py | 1 - 2 files changed, 2 insertions(+), 3 deletions(-) diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index b7727740..1b295b5d 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -4,7 +4,7 @@ default_stages: [commit] repos: - repo: https://github.com/psf/black - rev: 24.1.1 + rev: 24.8.0 hooks: - id: black @@ -15,7 +15,7 @@ repos: args: ["--profile", "black"] - repo: https://github.com/pre-commit/pre-commit-hooks - rev: v4.5.0 + rev: v4.6.0 hooks: - id: trailing-whitespace - id: end-of-file-fixer diff --git a/PyMPDATA/__init__.py b/PyMPDATA/__init__.py index 2d2600b1..8a00957b 100644 --- a/PyMPDATA/__init__.py +++ b/PyMPDATA/__init__.py @@ -5,7 +5,6 @@ PyMPDATA uses staggered grid with the following node placement for `PyMPDATA.scalar_field.ScalarField` and `PyMPDATA.vector_field.VectorField` elements: -![](https://github.com/atmos-cloud-sim-uj/PyMPDATA/releases/download/tip/readme_grid.png) .. include:: ../README.md From daa78072b99c05a2d188bc8b6788dd3bf7bfa8c8 Mon Sep 17 00:00:00 2001 From: Pawel Magnuszewski Date: Wed, 18 Sep 2024 14:41:37 +0200 Subject: [PATCH 04/10] move some docs outside readme --- PyMPDATA/__init__.py | 15 ++- README.md | 3 +- docs/markdown/debugging_contribution.md | 89 +++++++++++++++ docs/markdown/dependencies.md | 9 ++ docs/markdown/examples.md | 14 +++ docs/markdown/grid.md | 14 +++ docs/markdown/grid_2.md | 19 ++++ docs/markdown/grid_3.md | 10 ++ docs/markdown/intro.md | 77 +++++++++++++ docs/markdown/julia_grid.md | 27 +++++ docs/markdown/matlab_grid.md | 35 ++++++ docs/markdown/options.md | 55 +++++++++ docs/markdown/package_structure.md | 16 +++ docs/markdown/python_grid.md | 38 +++++++ docs/markdown/python_grid_2.md | 30 +++++ docs/markdown/python_solver.md | 35 ++++++ docs/markdown/stepper_solver.md | 143 ++++++++++++++++++++++++ 17 files changed, 627 insertions(+), 2 deletions(-) create mode 100644 docs/markdown/debugging_contribution.md create mode 100644 docs/markdown/dependencies.md create mode 100644 docs/markdown/examples.md create mode 100644 docs/markdown/grid.md create mode 100644 docs/markdown/grid_2.md create mode 100644 docs/markdown/grid_3.md create mode 100644 docs/markdown/intro.md create mode 100644 docs/markdown/julia_grid.md create mode 100644 docs/markdown/matlab_grid.md create mode 100644 docs/markdown/options.md create mode 100644 docs/markdown/package_structure.md create mode 100644 docs/markdown/python_grid.md create mode 100644 docs/markdown/python_grid_2.md create mode 100644 docs/markdown/python_solver.md create mode 100644 docs/markdown/stepper_solver.md diff --git a/PyMPDATA/__init__.py b/PyMPDATA/__init__.py index 8a00957b..888223e2 100644 --- a/PyMPDATA/__init__.py +++ b/PyMPDATA/__init__.py @@ -6,7 +6,20 @@ `PyMPDATA.scalar_field.ScalarField` and `PyMPDATA.vector_field.VectorField` elements: -.. include:: ../README.md +.. include:: ../docs/markdown/intro.md +.. include:: ../docs/markdown/dependencies.md +.. include:: ../docs/markdown/examples.md +.. include:: ../docs/markdown/options.md +.. include:: ../docs/markdown/grid.md +.. include:: ../docs/markdown/python_grid.md +.. include:: ../docs/markdown/grid_2.md +.. include:: ../docs/markdown/julia_grid.md +.. include:: ../docs/markdown/matlab_grid.md +.. include:: ../docs/markdown/python_grid_2.md +.. include:: ../docs/markdown/grid_3.md +.. include:: ../docs/markdown/stepper_solver.md +.. include:: ../docs/markdown/python_solver.md +.. include:: ../docs/markdown/debugging_contribution.md """ diff --git a/README.md b/README.md index 9c70a715..251b1f9d 100644 --- a/README.md +++ b/README.md @@ -1,5 +1,4 @@ # pympdata logo - # PyMPDATA [![Python 3](https://img.shields.io/static/v1?label=Python&logo=Python&color=3776AB&message=3)](https://www.python.org/) @@ -154,6 +153,7 @@ Options = pyimport("PyMPDATA").Options options = Options(n_iters=2) ``` +
Matlab code (click to expand) @@ -162,6 +162,7 @@ Options = py.importlib.import_module('PyMPDATA').Options; options = Options(pyargs('n_iters', 2)); ```
+
Python code (click to expand) diff --git a/docs/markdown/debugging_contribution.md b/docs/markdown/debugging_contribution.md new file mode 100644 index 00000000..fb0f2f6a --- /dev/null +++ b/docs/markdown/debugging_contribution.md @@ -0,0 +1,89 @@ +#### Debugging + +PyMPDATA relies heavily on Numba to provide high-performance +number crunching operations. Arguably, one of the key advantage +of embracing Numba is that it can be easily switched off. This +brings multiple-order-of-magnitude drop in performance, yet +it also make the entire code of the library susceptible to +interactive debugging, one way of enabling it is by setting the +following environment variable before importing PyMPDATA: +
+Julia code (click to expand) + +```Julia +ENV["NUMBA_DISABLE_JIT"] = "1" +``` +
+
+Matlab code (click to expand) + +```Matlab +setenv('NUMBA_DISABLE_JIT', '1'); +``` +
+
+Python code (click to expand) + +```Python +import os +os.environ["NUMBA_DISABLE_JIT"] = "1" +``` +
+ +## Contributing, reporting issues, seeking support + +Submitting new code to the project, please preferably use [GitHub pull requests](https://github.com/open-atmos/PyMPDATA/pulls) +(or the [PyMPDATA-examples PR site](https://github.com/open-atmos/PyMPDATA-examples/pulls) if working on examples) - it helps to keep record of code authorship, +track and archive the code review workflow and allows to benefit +from the continuous integration setup which automates execution of tests +with the newly added code. + +As of now, the copyright to the entire PyMPDATA codebase is with the Jagiellonian +University, and code contributions are assumed to imply transfer of copyright. +Should there be a need to make an exception, please indicate it when creating +a pull request or contributing code in any other way. In any case, +the license of the contributed code must be compatible with GPL v3. + +Developing the code, we follow [The Way of Python](https://www.python.org/dev/peps/pep-0020/) and +the [KISS principle](https://en.wikipedia.org/wiki/KISS_principle). +The codebase has greatly benefited from [PyCharm code inspections](https://www.jetbrains.com/help/pycharm/code-inspection.html) +and [Pylint](https://pylint.org) code analysis (Pylint checks are part of the +CI workflows). + +Issues regarding any incorrect, unintuitive or undocumented bahaviour of +PyMPDATA are best to be reported on the [GitHub issue tracker](https://github.com/open-atmos/PyMPDATA/issues/new). +Feature requests are recorded in the "Ideas..." [PyMPDATA wiki page](https://github.com/open-atmos/PyMPDATA/wiki/Ideas-for-new-features-and-examples). + +We encourage to use the [GitHub Discussions](https://github.com/open-atmos/PyMPDATA/discussions) feature +(rather than the issue tracker) for seeking support in understanding, using and extending PyMPDATA code. + +Please use the PyMPDATA issue-tracking and dicsussion infrastructure for `PyMPDATA-examples` as well. +We look forward to your contributions and feedback. + +## Credits: +Development of PyMPDATA was supported by the EU through a grant of the [Foundation for Polish Science](http://fnp.org.pl) (POIR.04.04.00-00-5E1C/18). + +copyright: Jagiellonian University +licence: GPL v3 + +## Other open-source MPDATA implementations: +- mpdat_2d in babyEULAG (FORTRAN) + https://github.com/igfuw/bE_SDs/blob/master/babyEULAG.SDs.for#L741 +- mpdata-oop (C++, Fortran, Python) + https://github.com/igfuw/mpdata-oop +- apc-llc/mpdata (C++) + https://github.com/apc-llc/mpdata +- libmpdata++ (C++): + https://github.com/igfuw/libmpdataxx +- AtmosFOAM: + https://github.com/AtmosFOAM/AtmosFOAM/tree/947b192f69d973ea4a7cfab077eb5c6c6fa8b0cf/applications/solvers/advection/MPDATAadvectionFoam + +## Other Python packages for solving hyperbolic transport equations + +- PyPDE: https://pypi.org/project/PyPDE/ +- FiPy: https://pypi.org/project/FiPy/ +- ader: https://pypi.org/project/ader/ +- centpy: https://pypi.org/project/centpy/ +- mattflow: https://pypi.org/project/mattflow/ +- FastFD: https://pypi.org/project/FastFD/ +- Pyclaw: https://www.clawpack.org/pyclaw/ \ No newline at end of file diff --git a/docs/markdown/dependencies.md b/docs/markdown/dependencies.md new file mode 100644 index 00000000..983d02aa --- /dev/null +++ b/docs/markdown/dependencies.md @@ -0,0 +1,9 @@ +## Dependencies and installation + +To install PyMPDATA, one may use: ``pip install PyMPDATA`` (or +``pip install git+https://github.com/open-atmos/PyMPDATA.git`` to get updates beyond the latest release). +PyMPDATA depends on ``NumPy`` and ``Numba``. + +Running the tests shipped with the package requires additional packages listed in the +[test-time-requirements.txt](https://github.com/open-atmos/PyMPDATA/blob/main/test-time-requirements.txt) file +(which include ``PyMPDATA-examples``, see below). \ No newline at end of file diff --git a/docs/markdown/examples.md b/docs/markdown/examples.md new file mode 100644 index 00000000..bf481719 --- /dev/null +++ b/docs/markdown/examples.md @@ -0,0 +1,14 @@ +## Examples (Jupyter notebooks reproducing results from literature): + +PyMPDATA examples are bundled with PyMPDATA and located in the `examples` subfolder. +They constitute a separate [``PyMPDATA_examples``](https://pypi.org/p/PyMPDATA-examples) Python package which is also available at PyPI. +The examples have additional dependencies listed in [``PyMPDATA_examples`` package ``setup.py``](https://github.com/open-atmos/PyMPDATA/blob/main/examples/setup.py) file. +Running the examples requires the ``PyMPDATA_examples`` package to be installed. +Since the examples package includes Jupyter notebooks (and their execution requires write access), the suggested install and launch steps are: +``` +git clone https://github.com/open-atmos/PyMPDATA-examples.git +cd PyMPDATA-examples +pip install -e . +jupyter-notebook +``` +Alternatively, one can also install the examples package from pypi.org by using ``pip install PyMPDATA-examples``. diff --git a/docs/markdown/grid.md b/docs/markdown/grid.md new file mode 100644 index 00000000..925281cb --- /dev/null +++ b/docs/markdown/grid.md @@ -0,0 +1,14 @@ +#### Arakawa-C grid layer and boundary conditions + +In PyMPDATA, the solution domain is assumed to extend from the +first cell's boundary to the last cell's boundary (thus the +first scalar field value is at $\[\Delta x/2, \Delta y/2\]$. +The [``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) +and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) classes implement the +[Arakawa-C staggered grid](https://en.wikipedia.org/wiki/Arakawa_grids#Arakawa_C-grid) logic +in which: +- scalar fields are discretised onto cell centres (one value per cell), +- vector field components are discretised onto cell walls. + +The schematic of the employed grid/domain layout in two dimensions is given below +(with the Python code snippet generating the figure as a part of CI workflow): \ No newline at end of file diff --git a/docs/markdown/grid_2.md b/docs/markdown/grid_2.md new file mode 100644 index 00000000..c7a7bd72 --- /dev/null +++ b/docs/markdown/grid_2.md @@ -0,0 +1,19 @@ +The ``__init__`` methods of +[``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) +and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) +have the following signatures: +- [``ScalarField(data: np.ndarray, halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/scalar_field.py) +- [``VectorField(data: Tuple[np.ndarray, ...], halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/vector_field.py) +The ``data`` parameters are expected to be Numpy arrays or tuples of Numpy arrays, respectively. +The ``halo`` parameter is the extent of ghost-cell region that will surround the +data and will be used to implement boundary conditions. +Its value (in practice 1 or 2) is +dependent on maximal stencil extent for the MPDATA variant used and +can be easily obtained using the ``Options.n_halo`` property. + +As an example, the code below shows how to instantiate a scalar +and a vector field given a 2D constant-velocity problem, +using a grid of 24x24 points, Courant numbers of -0.5 and -0.25 +in "x" and "y" directions, respectively, with periodic boundary +conditions and with an initial Gaussian signal in the scalar field +(settings as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379)): \ No newline at end of file diff --git a/docs/markdown/grid_3.md b/docs/markdown/grid_3.md new file mode 100644 index 00000000..fd5ca36b --- /dev/null +++ b/docs/markdown/grid_3.md @@ -0,0 +1,10 @@ +Note that the shapes of arrays representing components +of the velocity field are different than the shape of +the scalar field array due to employment of the staggered grid. + +Besides the exemplified [``Periodic``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/periodic.html) class representing +periodic boundary conditions, PyMPDATA supports +[``Extrapolated``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/extrapolated.html), +[``Constant``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/constant.html) and +[``Polar``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/polar.html) +boundary conditions. \ No newline at end of file diff --git a/docs/markdown/intro.md b/docs/markdown/intro.md new file mode 100644 index 00000000..132a8a7f --- /dev/null +++ b/docs/markdown/intro.md @@ -0,0 +1,77 @@ +pympdata logo + +# PyMPDATA + +[![Python 3](https://img.shields.io/static/v1?label=Python&logo=Python&color=3776AB&message=3)](https://www.python.org/) +[![LLVM](https://img.shields.io/static/v1?label=LLVM&logo=LLVM&color=gold&message=Numba)](https://www.numba.org) +[![Linux OK](https://img.shields.io/static/v1?label=Linux&logo=Linux&color=yellow&message=%E2%9C%93)](https://en.wikipedia.org/wiki/Linux) +[![macOS OK](https://img.shields.io/static/v1?label=macOS&logo=Apple&color=silver&message=%E2%9C%93)](https://en.wikipedia.org/wiki/macOS) +[![Windows OK](https://img.shields.io/static/v1?label=Windows&logo=Windows&color=white&message=%E2%9C%93)](https://en.wikipedia.org/wiki/Windows) +[![Jupyter](https://img.shields.io/static/v1?label=Jupyter&logo=Jupyter&color=f37626&message=%E2%9C%93)](https://jupyter.org/) +[![Maintenance](https://img.shields.io/badge/Maintained%3F-yes-green.svg)](https://github.com/open-atmos/PyMPDATA/graphs/commit-activity) +[![OpenHub](https://www.openhub.net/p/atmos-cloud-sim-uj-PyMPDATA/widgets/project_thin_badge?format=gif)](https://www.openhub.net/p/atmos-cloud-sim-uj-PyMPDATA) +[![status](https://joss.theoj.org/papers/10e7361e43785dbb1b3d659c5b01757a/status.svg)](https://joss.theoj.org/papers/10e7361e43785dbb1b3d659c5b01757a) +[![DOI](https://zenodo.org/badge/225610671.svg)](https://zenodo.org/badge/latestdoi/225610671) +[![EU Funding](https://img.shields.io/static/v1?label=EU%20Funding%20by&color=103069&message=FNP&logoWidth=25&logo=image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAC4AAAAeCAYAAABTwyyaAAAEzklEQVRYw9WYS2yUVRiGn3P5ZzozpZ3aUsrNgoKlKBINmkhpCCwwxIAhsDCpBBIWhmCMMYTEhSJ4i9EgBnSBEm81MRrFBhNXEuUSMCopiRWLQqEGLNgr085M5//POS46NNYFzHQ6qGc1i5nzP/P973m/9ztCrf7A8T9csiibCocUbvTzfxLcAcaM3cY3imXz25lT3Y34G7gQYAKV3+bFAHcATlBTPogJNADG92iY28FHW97kyPbnuW/W7xgzAhukQ9xe04PJeOT0HkQRwK0TlEeGWb/kOO9v3kdD3a8YK9GhDMfa6mg9fxunOm/lWPtcpDI4K7n/jnN8+uQbrFrUSiwU/DtSEUB/MsKKBT+zYslJqiYNgVE4JwhHkzy86wlWvrKVWDSZ/YFjZlU39yw4y/rGoyQGowWB67zl4QQue+jssMdXrQvZ/00jyeHwqCgDKwnsiJjSvkYAxsG5K9WsenYbJdqAtAjhCIxCSZt/4fK1w5A2WCvxrUAKCHwNVoA2aGmvq11jJQQapEXrgMBKqmJJugejKGWLIxXrBPFoigfv/omd675gRkU/xgqUDlAhH3UDaAAlLSqUQekAYyVTyhLs3tDMsvntlIYzOFcEcOcEGd9jx9oDbGs6QO0t/Tijxi9S4bhzxiWaVh5m94Zm0n7oui4ybo0raUlcncQnxx+g+WgDF/vLoYDmoqSl/dJUnt7XRCoTZjij0Z6Pc2LiNS4EBBkNvoeOJXN+yPWWSZeANOhwJq/98nKVwNdoL8B5AROxBKBL0gjh8DMhdCh3eJnrA0yqhLpplwmyup6IajvAOIGfKGVx3VmCRGnOMpe5QAdG0bT8CAeeep0d6z6nqjSJnQiZWEllLMWrmz6k+fE9rGk8MVqYgsGv5ZH2i1Opr+9kajzB5d74hKQ+KS3d/WVMLhtgdu1lriRiOR/4nDVunaR24x7qp3UV5Cb/fJvC83nv26W81LIK58SYNFmwq4hsGx/5BwKlzYRma2NUthgOJSew4i7ru9nJYCQF5tApb2yvjiDQKJV/IfJKh0o6qssSLKv/jcAoRKHQQzE2Lj2OMV5OkWFc4MZIpsev8uXWXRx6ZicbGk8QZLxxgwe+x/rlR3h3816+f2E7lbEU+ZDn3vKVpePCdFovzCISHqbl5EIoQOteKMPB1rto65zNyfOz+KOrGl06lHPQyi/WOohH0/T0l1MZH6A3GUEKl7Pmr2la6wBrBWWRDP2DUcqjKVKBGom9RZmABAykwnglafpSJSPQvsfiOR0EQ7ExVmazA8cY6N4K1iw6RdAXRwi4mgrheT5Dvs4LeuS81a15Ll/3dQisFVSVpnj7sf1sX/sZvhAc+6UOrQyBVUQ8gx/orFmDsZqtaw/y1qZ9zKjp5vDpenyjcNe+cLNmTiUdf/bEOddVQ0VpgsOn54ET+EYxvWKALSu+5tGG76it7MNaiZKGQ23zCIcMfUMxBnrjN3fmHHvCAlp+vJcXWx6itqoXpAEnUNLx8iMfo5Xh1i17R3PJYCpC2cZ3qK3sQ8WGEDDuXlAQuFKGHzpmopXhTNfk0bmxs7uC1w6uJul79AxFkMIiBJy5UoUWjrZLU5DCFdTARDHuDqVw+OkSwI0MCEW4gtNF2BPrBCo8fKNbtILWX9aUDqFqHnn7AAAAAElFTkSuQmCC)](https://www.fnp.org.pl/en/) +[![PL Funding](https://img.shields.io/static/v1?label=PL%20Funding%20by&color=d21132&message=NCN&logoWidth=25&logo=image/png;base64,iVBORw0KGgoAAAANSUhEUgAAABQAAAANCAYAAACpUE5eAAAABmJLR0QA/wD/AP+gvaeTAAAAKUlEQVQ4jWP8////fwYqAiZqGjZqIHUAy4dJS6lqIOMdEZvRZDPcDQQAb3cIaY1Sbi4AAAAASUVORK5CYII=)](https://www.ncn.gov.pl/?language=en) +[![License: GPL v3](https://img.shields.io/badge/License-GPL%20v3-blue.svg)](https://www.gnu.org/licenses/gpl-3.0.html) + +[![Github Actions Build Status](https://github.com/open-atmos/PyMPDATA/actions/workflows/tests+pypi.yml/badge.svg?branch=main)](https://github.com/open-atmos/PyMPDATA/actions) +[![Appveyor Build status](http://ci.appveyor.com/api/projects/status/github/open-atmos/PyMPDATA?branch=main&svg=true)](https://ci.appveyor.com/project/slayoo/pympdata/branch/main) +[![Coverage Status](https://codecov.io/gh/open-atmos/PyMPDATA/branch/main/graph/badge.svg)](https://app.codecov.io/gh/open-atmos/PyMPDATA) + +[![PyPI version](https://badge.fury.io/py/PyMPDATA.svg)](https://pypi.org/project/PyMPDATA) +[![API docs](https://shields.mitmproxy.org/badge/docs-pdoc.dev-brightgreen.svg)](https://open-atmos.github.io/PyMPDATA/index.html) + +PyMPDATA is a high-performance Numba-accelerated Pythonic implementation of the MPDATA + algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond. +MPDATA numerically solves generalised transport equations - + partial differential equations used to model conservation/balance laws, scalar-transport problems, + convection-diffusion phenomena. +As of the current version, PyMPDATA supports homogeneous transport + in 1D, 2D and 3D using structured meshes, optionally + generalised by employment of a Jacobian of coordinate transformation. +PyMPDATA includes implementation of a set of MPDATA variants including + the non-oscillatory option, infinite-gauge, divergent-flow, double-pass donor cell (DPDC) and + third-order-terms options. +It also features support for integration of Fickian-terms in advection-diffusion + problems using the pseudo-transport velocity approach. +In 2D and 3D simulations, domain-decomposition is used for multi-threaded parallelism. + +PyMPDATA is engineered purely in Python targeting both performance and usability, + the latter encompassing research users', developers' and maintainers' perspectives. +From researcher's perspective, PyMPDATA offers hassle-free installation on multitude + of platforms including Linux, OSX and Windows, and eliminates compilation stage + from the perspective of the user. +From developers' and maintainers' perspective, PyMPDATA offers a suite of unit tests, + multi-platform continuous integration setup, + seamless integration with Python development aids including debuggers and profilers. + +PyMPDATA design features + a custom-built multi-dimensional Arakawa-C grid layer allowing + to concisely represent multi-dimensional stencil operations on both + scalar and vector fields. +The grid layer is built on top of NumPy's ndarrays (using "C" ordering) + using the Numba's ``@njit`` functionality for high-performance array traversals. +It enables one to code once for multiple dimensions, and automatically + handles (and hides from the user) any halo-filling logic related with boundary conditions. +Numba ``prange()`` functionality is used for implementing multi-threading + (it offers analogous functionality to OpenMP parallel loop execution directives). +The Numba's deviation from Python semantics rendering [closure variables + as compile-time constants](https://numba.pydata.org/numba-doc/dev/reference/pysemantics.html#global-and-closure-variables) + is extensively exploited within ``PyMPDATA`` + code base enabling the just-in-time compilation to benefit from + information on domain extents, algorithm variant used and problem + characteristics (e.g., coordinate transformation used, or lack thereof). + +A separate project called [``PyMPDATA-MPI``](https://github.com/open-atmos/PyMPDATA-MPI) + depicts how [``numba-mpi``](https://pypi.org/project/numba-mpi) can be used + to enable distributed memory parallelism in PyMPDATA. + +The [``PyMPDATA-examples``](https://pypi.org/project/PyMPDATA-examples/) + package covers a set of examples presented in the form of Jupyer notebooks + offering single-click deployment in the cloud using [mybinder.org](https://mybinder.org) + or using [colab.research.google.com](https://colab.research.google.com/). +The examples reproduce results from several published + works on MPDATA and its applications, and provide a validation of the implementation + and its performance. \ No newline at end of file diff --git a/docs/markdown/julia_grid.md b/docs/markdown/julia_grid.md new file mode 100644 index 00000000..a45aefed --- /dev/null +++ b/docs/markdown/julia_grid.md @@ -0,0 +1,27 @@ +
+Julia code (click to expand) + +```Julia +ScalarField = pyimport("PyMPDATA").ScalarField +VectorField = pyimport("PyMPDATA").VectorField +Periodic = pyimport("PyMPDATA.boundary_conditions").Periodic + +nx, ny = 24, 24 +Cx, Cy = -.5, -.25 +idx = CartesianIndices((nx, ny)) +halo = options.n_halo +advectee = ScalarField( + data=exp.( + -(getindex.(idx, 1) .- .5 .- nx/2).^2 / (2*(nx/10)^2) + -(getindex.(idx, 2) .- .5 .- ny/2).^2 / (2*(ny/10)^2) + ), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +advector = VectorField( + data=(fill(Cx, (nx+1, ny)), fill(Cy, (nx, ny+1))), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +``` +
\ No newline at end of file diff --git a/docs/markdown/matlab_grid.md b/docs/markdown/matlab_grid.md new file mode 100644 index 00000000..1447250f --- /dev/null +++ b/docs/markdown/matlab_grid.md @@ -0,0 +1,35 @@ +
+Matlab code (click to expand) + +```Matlab +ScalarField = py.importlib.import_module('PyMPDATA').ScalarField; +VectorField = py.importlib.import_module('PyMPDATA').VectorField; +Periodic = py.importlib.import_module('PyMPDATA.boundary_conditions').Periodic; + +nx = int32(24); +ny = int32(24); + +Cx = -.5; +Cy = -.25; + +[xi, yi] = meshgrid(double(0:1:nx-1), double(0:1:ny-1)); + +halo = options.n_halo; +advectee = ScalarField(pyargs(... + 'data', py.numpy.array(exp( ... + -(xi+.5-double(nx)/2).^2 / (2*(double(nx)/10)^2) ... + -(yi+.5-double(ny)/2).^2 / (2*(double(ny)/10)^2) ... + )), ... + 'halo', halo, ... + 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... +)); +advector = VectorField(pyargs(... + 'data', py.tuple({ ... + Cx * py.numpy.ones(int32([nx+1 ny])), ... + Cy * py.numpy.ones(int32([nx ny+1])) ... + }), ... + 'halo', halo, ... + 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... +)); +``` +
\ No newline at end of file diff --git a/docs/markdown/options.md b/docs/markdown/options.md new file mode 100644 index 00000000..4c991179 --- /dev/null +++ b/docs/markdown/options.md @@ -0,0 +1,55 @@ +## Options class + +The [``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class +groups both algorithm variant options as well as some implementation-related +flags that need to be set at the first place. All are set at the time +of instantiation using the following keyword arguments of the constructor +(all having default values indicated below): +- ``n_iters: int = 2``: number of iterations (2 means upwind + one corrective iteration) +- ``infinite_gauge: bool = False``: flag enabling the infinite-gauge option (does not maintain sign of the advected field, thus in practice implies switching flux corrected transport on) +- ``divergent_flow: bool = False``: flag enabling divergent-flow terms when calculating antidiffusive velocity +- ``nonoscillatory: bool = False``: flag enabling the non-oscillatory or monotone variant (a.k.a flux-corrected transport option, FCT) +- ``third_order_terms: bool = False``: flag enabling third-order terms +- ``epsilon: float = 1e-15``: value added to potentially zero-valued denominators +- ``non_zero_mu_coeff: bool = False``: flag indicating if code for handling the Fickian term is to be optimised out +- ``DPDC: bool = False``: flag enabling double-pass donor cell option (recursive pseudovelocities) +- ``dimensionally_split: bool = False``: flag disabling cross-dimensional terms in antidiffusive velocity +- ``dtype: np.floating = np.float64``: floating point precision + +For a discussion of the above options, see e.g., [Smolarkiewicz & Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), +[Jaruga, Arabas et al. 2015](https://doi.org/10.5194/gmd-8-1005-2015) and [Olesik, Arabas et al. 2020](https://arxiv.org/abs/2011.14726) +(the last with examples using PyMPDATA). + +In most use cases of PyMPDATA, the first thing to do is to instantiate the +[``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class +with arguments suiting the problem at hand, e.g.: + +
+Julia code (click to expand) + +```Julia +using Pkg +Pkg.add("PyCall") +using PyCall +Options = pyimport("PyMPDATA").Options +options = Options(n_iters=2) +``` +
+ +
+Matlab code (click to expand) + +```Matlab +Options = py.importlib.import_module('PyMPDATA').Options; +options = Options(pyargs('n_iters', 2)); +``` +
+ +
+Python code (click to expand) + +```Python +from PyMPDATA import Options +options = Options(n_iters=2) +``` +
\ No newline at end of file diff --git a/docs/markdown/package_structure.md b/docs/markdown/package_structure.md new file mode 100644 index 00000000..d8c04e4f --- /dev/null +++ b/docs/markdown/package_structure.md @@ -0,0 +1,16 @@ +## Package structure and API: + +In short, PyMPDATA numerically solves the following equation: + +$$ \partial_t (G \psi) + \nabla \cdot (Gu \psi) + \mu \Delta (G \psi) = 0 $$ + +where scalar field $\psi$ is referred to as the advectee, +vector field u is referred to as advector, and the G factor corresponds to optional coordinate transformation. +The inclusion of the Fickian diffusion term is optional and is realised through modification of the +advective velocity field with MPDATA handling both the advection and diffusion (for discussion +see, e.g. [Smolarkiewicz and Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), sec. 3.5, par. 4). + +The key classes constituting the PyMPDATA interface are summarised below with code +snippets exemplifying usage of PyMPDATA from Python, Julia and Matlab. + +A [pdoc-generated](https://pdoc.dev/) documentation of PyMPDATA public API is maintained at: [https://open-atmos.github.io/PyMPDATA]( https://open-atmos.github.io/PyMPDATA/index.html) diff --git a/docs/markdown/python_grid.md b/docs/markdown/python_grid.md new file mode 100644 index 00000000..1a3c08b0 --- /dev/null +++ b/docs/markdown/python_grid.md @@ -0,0 +1,38 @@ +
+Python code (click to expand) + +```Python +import numpy as np +from matplotlib import pyplot + +dx, dy = .2, .3 +grid = (10, 5) + +pyplot.scatter(*np.mgrid[ + dx / 2 : grid[0] * dx : dx, + dy / 2 : grid[1] * dy : dy + ], color='red', + label='scalar-field values at cell centres' +) +pyplot.quiver(*np.mgrid[ + 0 : (grid[0]+1) * dx : dx, + dy / 2 : grid[1] * dy : dy + ], 1, 0, pivot='mid', color='green', width=.005, + label='vector-field x-component values at cell walls' +) +pyplot.quiver(*np.mgrid[ + dx / 2 : grid[0] * dx : dx, + 0: (grid[1] + 1) * dy : dy + ], 0, 1, pivot='mid', color='blue', width=.005, + label='vector-field y-component values at cell walls' +) +pyplot.xticks(np.linspace(0, grid[0]*dx, grid[0]+1)) +pyplot.yticks(np.linspace(0, grid[1]*dy, grid[1]+1)) +pyplot.title(f'staggered grid layout (grid={grid}, dx={dx}, dy={dy})') +pyplot.xlabel('x') +pyplot.ylabel('y') +pyplot.legend(bbox_to_anchor=(.1, -.1), loc='upper left', ncol=1) +pyplot.grid() +pyplot.savefig('readme_grid.png') +``` +
diff --git a/docs/markdown/python_grid_2.md b/docs/markdown/python_grid_2.md new file mode 100644 index 00000000..87aed5c4 --- /dev/null +++ b/docs/markdown/python_grid_2.md @@ -0,0 +1,30 @@ +
+
+Python code (click to expand) + +```Python +from PyMPDATA import ScalarField +from PyMPDATA import VectorField +from PyMPDATA.boundary_conditions import Periodic +import numpy as np + +nx, ny = 24, 24 +Cx, Cy = -.5, -.25 +halo = options.n_halo + +xi, yi = np.indices((nx, ny), dtype=float) +advectee = ScalarField( + data=np.exp( + -(xi+.5-nx/2)**2 / (2*(nx/10)**2) + -(yi+.5-ny/2)**2 / (2*(ny/10)**2) + ), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +advector = VectorField( + data=(np.full((nx + 1, ny), Cx), np.full((nx, ny + 1), Cy)), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +``` +
diff --git a/docs/markdown/python_solver.md b/docs/markdown/python_solver.md new file mode 100644 index 00000000..cfe49941 --- /dev/null +++ b/docs/markdown/python_solver.md @@ -0,0 +1,35 @@ +
+Python code (click to expand) + +```Python +def plot(psi, zlim, norm=None): + xi, yi = np.indices(psi.shape) + fig, ax = pyplot.subplots(subplot_kw={"projection": "3d"}) + pyplot.gca().plot_wireframe( + xi+.5, yi+.5, + psi, color='red', linewidth=.5 + ) + ax.set_zlim(zlim) + for axis in (ax.xaxis, ax.yaxis, ax.zaxis): + axis.pane.fill = False + axis.pane.set_edgecolor('black') + axis.pane.set_alpha(1) + ax.grid(False) + ax.set_zticks([]) + ax.set_xlabel('x/dx') + ax.set_ylabel('y/dy') + ax.set_proj_type('ortho') + cnt = ax.contourf(xi+.5, yi+.5, psi, zdir='z', offset=-1, norm=norm) + cbar = pyplot.colorbar(cnt, pad=.1, aspect=10, fraction=.04) + return cbar.norm + +zlim = (-1, 1) +norm = plot(state_0, zlim) +pyplot.savefig('readme_gauss_0.png') +plot(state, zlim, norm) +pyplot.savefig('readme_gauss.png') +``` +
+ +![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss_0.png) +![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss.png) diff --git a/docs/markdown/stepper_solver.md b/docs/markdown/stepper_solver.md new file mode 100644 index 00000000..7e93972e --- /dev/null +++ b/docs/markdown/stepper_solver.md @@ -0,0 +1,143 @@ +#### Stepper + +The logic of the MPDATA iterative solver is represented +in PyMPDATA by the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class. + +When instantiating the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html), the user has a choice +of either supplying just the number of dimensions or specialising the stepper for a given grid: +
+Julia code (click to expand) + +```Julia +Stepper = pyimport("PyMPDATA").Stepper + +stepper = Stepper(options=options, n_dims=2) +``` +
+
+Matlab code (click to expand) + +```Matlab +Stepper = py.importlib.import_module('PyMPDATA').Stepper; + +stepper = Stepper(pyargs(... + 'options', options, ... + 'n_dims', int32(2) ... +)); +``` +
+
+Python code (click to expand) + +```Python +from PyMPDATA import Stepper + +stepper = Stepper(options=options, n_dims=2) +``` +
+or +
+Julia code (click to expand) + +```Julia +stepper = Stepper(options=options, grid=(nx, ny)) +``` +
+
+Matlab code (click to expand) + +```Matlab +stepper = Stepper(pyargs(... + 'options', options, ... + 'grid', py.tuple({nx, ny}) ... +)); +``` +
+
+Python code (click to expand) + +```Python +stepper = Stepper(options=options, grid=(nx, ny)) +``` +
+ +In the latter case, noticeably +faster execution can be expected, however the resultant +stepper is less versatile as bound to the given grid size. +If number of dimensions is supplied only, the integration +might take longer, yet same instance of the +stepper can be used for different grids. + +Since creating an instance of the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class +involves time-consuming compilation of the algorithm code, +the class is equipped with a cache logic - subsequent +calls with same arguments return references to previously +instantiated objects. Instances of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) contain no +mutable data and are (thread-)safe to be reused. + +The init method of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) has an optional +``non_unit_g_factor`` argument which is a Boolean flag +enabling handling of the G factor term which can be used to +represent coordinate transformations and/or variable fluid density. + +Optionally, the number of threads to use for domain decomposition +in the first (non-contiguous) dimension during 2D and 3D calculations +may be specified using the optional ``n_threads`` argument with a +default value of ``numba.get_num_threads()``. The multi-threaded +logic of PyMPDATA depends thus on settings of numba, namely on the +selected threading layer (either via ``NUMBA_THREADING_LAYER`` env +var or via ``numba.config.THREADING_LAYER``) and the selected size of the +thread pool (``NUMBA_NUM_THREADS`` env var or ``numba.config.NUMBA_NUM_THREADS``). + + +#### Solver + +Instances of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class are used to control +the integration and access solution data. During instantiation, +additional memory required by the solver is +allocated according to the options provided. + +The only method of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class besides the +init is [``advance(n_steps, mu_coeff, ...)``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html#Solver.advance) +which advances the solution by ``n_steps`` timesteps, optionally +taking into account a given diffusion coefficient ``mu_coeff``. + +Solution state is accessible through the ``Solver.advectee`` property. +Multiple solver[s] can share a single stepper, e.g., as exemplified in the shallow-water system solution in the examples package. + +Continuing with the above code snippets, instantiating +a solver and making 75 integration steps looks as follows: +
+Julia code (click to expand) + +```Julia +Solver = pyimport("PyMPDATA").Solver +solver = Solver(stepper=stepper, advectee=advectee, advector=advector) +solver.advance(n_steps=75) +state = solver.advectee.get() +``` +
+
+Matlab code (click to expand) + +```Matlab +Solver = py.importlib.import_module('PyMPDATA').Solver; +solver = Solver(pyargs('stepper', stepper, 'advectee', advectee, 'advector', advector)); +solver.advance(pyargs('n_steps', 75)); +state = solver.advectee.get(); +``` +
+
+Python code (click to expand) + +```Python +from PyMPDATA import Solver + +solver = Solver(stepper=stepper, advectee=advectee, advector=advector) +state_0 = solver.advectee.get().copy() +solver.advance(n_steps=75) +state = solver.advectee.get() +``` +
+ +Now let's plot the results using `matplotlib` roughly as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379): \ No newline at end of file From 67ba2fa984f3c55c37f4aa2845e15e3926023944 Mon Sep 17 00:00:00 2001 From: pawel Date: Tue, 24 Sep 2024 12:13:20 +0200 Subject: [PATCH 05/10] touch ups to readme --- README.md | 507 +++--------------------------------------------------- 1 file changed, 26 insertions(+), 481 deletions(-) diff --git a/README.md b/README.md index 251b1f9d..a72762f9 100644 --- a/README.md +++ b/README.md @@ -22,59 +22,29 @@ [![PyPI version](https://badge.fury.io/py/PyMPDATA.svg)](https://pypi.org/project/PyMPDATA) [![API docs](https://shields.mitmproxy.org/badge/docs-pdoc.dev-brightgreen.svg)](https://open-atmos.github.io/PyMPDATA/index.html) + PyMPDATA is a high-performance Numba-accelerated Pythonic implementation of the MPDATA - algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond. -MPDATA numerically solves generalised transport equations - - partial differential equations used to model conservation/balance laws, scalar-transport problems, - convection-diffusion phenomena. -As of the current version, PyMPDATA supports homogeneous transport - in 1D, 2D and 3D using structured meshes, optionally - generalised by employment of a Jacobian of coordinate transformation. -PyMPDATA includes implementation of a set of MPDATA variants including - the non-oscillatory option, infinite-gauge, divergent-flow, double-pass donor cell (DPDC) and - third-order-terms options. -It also features support for integration of Fickian-terms in advection-diffusion - problems using the pseudo-transport velocity approach. -In 2D and 3D simulations, domain-decomposition is used for multi-threaded parallelism. - -PyMPDATA is engineered purely in Python targeting both performance and usability, - the latter encompassing research users', developers' and maintainers' perspectives. -From researcher's perspective, PyMPDATA offers hassle-free installation on multitude - of platforms including Linux, OSX and Windows, and eliminates compilation stage - from the perspective of the user. -From developers' and maintainers' perspective, PyMPDATA offers a suite of unit tests, - multi-platform continuous integration setup, - seamless integration with Python development aids including debuggers and profilers. - -PyMPDATA design features - a custom-built multi-dimensional Arakawa-C grid layer allowing - to concisely represent multi-dimensional stencil operations on both - scalar and vector fields. -The grid layer is built on top of NumPy's ndarrays (using "C" ordering) - using the Numba's ``@njit`` functionality for high-performance array traversals. -It enables one to code once for multiple dimensions, and automatically - handles (and hides from the user) any halo-filling logic related with boundary conditions. -Numba ``prange()`` functionality is used for implementing multi-threading - (it offers analogous functionality to OpenMP parallel loop execution directives). -The Numba's deviation from Python semantics rendering [closure variables - as compile-time constants](https://numba.pydata.org/numba-doc/dev/reference/pysemantics.html#global-and-closure-variables) - is extensively exploited within ``PyMPDATA`` - code base enabling the just-in-time compilation to benefit from - information on domain extents, algorithm variant used and problem - characteristics (e.g., coordinate transformation used, or lack thereof). +algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond for +numerically solving generalised convection-diffusion PDEs in 1D, 2D and 3D structured meshes +with coordinate transformations. + +In short, PyMPDATA numerically solves the following equation: + +$$ \partial_t (G \psi) + \nabla \cdot (Gu \psi) + \mu \Delta (G \psi) = 0 $$ + +where scalar field $\psi$ is referred to as the advectee, +vector field u is referred to as advector, and the G factor corresponds to optional coordinate transformation. +The inclusion of the Fickian diffusion term is optional and is realised through modification of the +advective velocity field with MPDATA handling both the advection and diffusion (for discussion +see, e.g. [Smolarkiewicz and Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), sec. 3.5, par. 4). + +PyMPDATA [documentation](https://open-atmos.github.io/PyMPDATA/index.html) is generated via [``pdoc``](https://pdoc.dev/). A separate project called [``PyMPDATA-MPI``](https://github.com/open-atmos/PyMPDATA-MPI) depicts how [``numba-mpi``](https://pypi.org/project/numba-mpi) can be used - to enable distributed memory parallelism in PyMPDATA. - -The [``PyMPDATA-examples``](https://pypi.org/project/PyMPDATA-examples/) - package covers a set of examples presented in the form of Jupyer notebooks - offering single-click deployment in the cloud using [mybinder.org](https://mybinder.org) - or using [colab.research.google.com](https://colab.research.google.com/). -The examples reproduce results from several published - works on MPDATA and its applications, and provide a validation of the implementation + to enable distributed memory parallelism in PyMPDATA.applications, and provide a validation of the implementation and its performance. - + ## Dependencies and installation To install PyMPDATA, one may use: ``pip install PyMPDATA`` (or @@ -102,464 +72,39 @@ Alternatively, one can also install the examples package from pypi.org by using ## Package structure and API: -In short, PyMPDATA numerically solves the following equation: - -$$ \partial_t (G \psi) + \nabla \cdot (Gu \psi) + \mu \Delta (G \psi) = 0 $$ - -where scalar field $\psi$ is referred to as the advectee, -vector field u is referred to as advector, and the G factor corresponds to optional coordinate transformation. -The inclusion of the Fickian diffusion term is optional and is realised through modification of the -advective velocity field with MPDATA handling both the advection and diffusion (for discussion -see, e.g. [Smolarkiewicz and Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), sec. 3.5, par. 4). - -The key classes constituting the PyMPDATA interface are summarised below with code -snippets exemplifying usage of PyMPDATA from Python, Julia and Matlab. - -A [pdoc-generated](https://pdoc.dev/) documentation of PyMPDATA public API is maintained at: [https://open-atmos.github.io/PyMPDATA]( https://open-atmos.github.io/PyMPDATA/index.html) +The key classes constituting the PyMPDATA interface are summarised below. #### Options class The [``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class groups both algorithm variant options as well as some implementation-related -flags that need to be set at the first place. All are set at the time -of instantiation using the following keyword arguments of the constructor -(all having default values indicated below): -- ``n_iters: int = 2``: number of iterations (2 means upwind + one corrective iteration) -- ``infinite_gauge: bool = False``: flag enabling the infinite-gauge option (does not maintain sign of the advected field, thus in practice implies switching flux corrected transport on) -- ``divergent_flow: bool = False``: flag enabling divergent-flow terms when calculating antidiffusive velocity -- ``nonoscillatory: bool = False``: flag enabling the non-oscillatory or monotone variant (a.k.a flux-corrected transport option, FCT) -- ``third_order_terms: bool = False``: flag enabling third-order terms -- ``epsilon: float = 1e-15``: value added to potentially zero-valued denominators -- ``non_zero_mu_coeff: bool = False``: flag indicating if code for handling the Fickian term is to be optimised out -- ``DPDC: bool = False``: flag enabling double-pass donor cell option (recursive pseudovelocities) -- ``dimensionally_split: bool = False``: flag disabling cross-dimensional terms in antidiffusive velocity -- ``dtype: np.floating = np.float64``: floating point precision - -For a discussion of the above options, see e.g., [Smolarkiewicz & Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), -[Jaruga, Arabas et al. 2015](https://doi.org/10.5194/gmd-8-1005-2015) and [Olesik, Arabas et al. 2020](https://arxiv.org/abs/2011.14726) -(the last with examples using PyMPDATA). - -In most use cases of PyMPDATA, the first thing to do is to instantiate the -[``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class -with arguments suiting the problem at hand, e.g.: -
-Julia code (click to expand) - -```Julia -using Pkg -Pkg.add("PyCall") -using PyCall -Options = pyimport("PyMPDATA").Options -options = Options(n_iters=2) -``` -
+flags that need to be set at the first place. -
-Matlab code (click to expand) - -```Matlab -Options = py.importlib.import_module('PyMPDATA').Options; -options = Options(pyargs('n_iters', 2)); -``` -
- -
-Python code (click to expand) - -```Python -from PyMPDATA import Options -options = Options(n_iters=2) -``` -
- -#### Arakawa-C grid layer and boundary conditions +#### Arakawa-C grid layer In PyMPDATA, the solution domain is assumed to extend from the first cell's boundary to the last cell's boundary (thus the first scalar field value is at $\[\Delta x/2, \Delta y/2\]$. The [``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) classes implement the -[Arakawa-C staggered grid](https://en.wikipedia.org/wiki/Arakawa_grids#Arakawa_C-grid) logic -in which: -- scalar fields are discretised onto cell centres (one value per cell), -- vector field components are discretised onto cell walls. - -The schematic of the employed grid/domain layout in two dimensions is given below -(with the Python code snippet generating the figure as a part of CI workflow): -
-Python code (click to expand) - -```Python -import numpy as np -from matplotlib import pyplot - -dx, dy = .2, .3 -grid = (10, 5) - -pyplot.scatter(*np.mgrid[ - dx / 2 : grid[0] * dx : dx, - dy / 2 : grid[1] * dy : dy - ], color='red', - label='scalar-field values at cell centres' -) -pyplot.quiver(*np.mgrid[ - 0 : (grid[0]+1) * dx : dx, - dy / 2 : grid[1] * dy : dy - ], 1, 0, pivot='mid', color='green', width=.005, - label='vector-field x-component values at cell walls' -) -pyplot.quiver(*np.mgrid[ - dx / 2 : grid[0] * dx : dx, - 0: (grid[1] + 1) * dy : dy - ], 0, 1, pivot='mid', color='blue', width=.005, - label='vector-field y-component values at cell walls' -) -pyplot.xticks(np.linspace(0, grid[0]*dx, grid[0]+1)) -pyplot.yticks(np.linspace(0, grid[1]*dy, grid[1]+1)) -pyplot.title(f'staggered grid layout (grid={grid}, dx={dx}, dy={dy})') -pyplot.xlabel('x') -pyplot.ylabel('y') -pyplot.legend(bbox_to_anchor=(.1, -.1), loc='upper left', ncol=1) -pyplot.grid() -pyplot.savefig('readme_grid.png') -``` -
- -![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_grid.png) - -The ``__init__`` methods of -[``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) -and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) -have the following signatures: -- [``ScalarField(data: np.ndarray, halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/scalar_field.py) -- [``VectorField(data: Tuple[np.ndarray, ...], halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/vector_field.py) -The ``data`` parameters are expected to be Numpy arrays or tuples of Numpy arrays, respectively. -The ``halo`` parameter is the extent of ghost-cell region that will surround the -data and will be used to implement boundary conditions. -Its value (in practice 1 or 2) is -dependent on maximal stencil extent for the MPDATA variant used and -can be easily obtained using the ``Options.n_halo`` property. - -As an example, the code below shows how to instantiate a scalar -and a vector field given a 2D constant-velocity problem, -using a grid of 24x24 points, Courant numbers of -0.5 and -0.25 -in "x" and "y" directions, respectively, with periodic boundary -conditions and with an initial Gaussian signal in the scalar field -(settings as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379)): -
-Julia code (click to expand) - -```Julia -ScalarField = pyimport("PyMPDATA").ScalarField -VectorField = pyimport("PyMPDATA").VectorField -Periodic = pyimport("PyMPDATA.boundary_conditions").Periodic - -nx, ny = 24, 24 -Cx, Cy = -.5, -.25 -idx = CartesianIndices((nx, ny)) -halo = options.n_halo -advectee = ScalarField( - data=exp.( - -(getindex.(idx, 1) .- .5 .- nx/2).^2 / (2*(nx/10)^2) - -(getindex.(idx, 2) .- .5 .- ny/2).^2 / (2*(ny/10)^2) - ), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -advector = VectorField( - data=(fill(Cx, (nx+1, ny)), fill(Cy, (nx, ny+1))), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -``` -
-
-Matlab code (click to expand) +[Arakawa-C staggered grid](https://en.wikipedia.org/wiki/Arakawa_grids#Arakawa_C-grid) logic. -```Matlab -ScalarField = py.importlib.import_module('PyMPDATA').ScalarField; -VectorField = py.importlib.import_module('PyMPDATA').VectorField; -Periodic = py.importlib.import_module('PyMPDATA.boundary_conditions').Periodic; +#### Boundary conditions -nx = int32(24); -ny = int32(24); - -Cx = -.5; -Cy = -.25; - -[xi, yi] = meshgrid(double(0:1:nx-1), double(0:1:ny-1)); - -halo = options.n_halo; -advectee = ScalarField(pyargs(... - 'data', py.numpy.array(exp( ... - -(xi+.5-double(nx)/2).^2 / (2*(double(nx)/10)^2) ... - -(yi+.5-double(ny)/2).^2 / (2*(double(ny)/10)^2) ... - )), ... - 'halo', halo, ... - 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... -)); -advector = VectorField(pyargs(... - 'data', py.tuple({ ... - Cx * py.numpy.ones(int32([nx+1 ny])), ... - Cy * py.numpy.ones(int32([nx ny+1])) ... - }), ... - 'halo', halo, ... - 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... -)); -``` -
-
-Python code (click to expand) - -```Python -from PyMPDATA import ScalarField -from PyMPDATA import VectorField -from PyMPDATA.boundary_conditions import Periodic -import numpy as np - -nx, ny = 24, 24 -Cx, Cy = -.5, -.25 -halo = options.n_halo - -xi, yi = np.indices((nx, ny), dtype=float) -advectee = ScalarField( - data=np.exp( - -(xi+.5-nx/2)**2 / (2*(nx/10)**2) - -(yi+.5-ny/2)**2 / (2*(ny/10)**2) - ), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -advector = VectorField( - data=(np.full((nx + 1, ny), Cx), np.full((nx, ny + 1), Cy)), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -``` -
- -Note that the shapes of arrays representing components -of the velocity field are different than the shape of -the scalar field array due to employment of the staggered grid. - -Besides the exemplified [``Periodic``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/periodic.html) class representing -periodic boundary conditions, PyMPDATA supports -[``Extrapolated``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/extrapolated.html), -[``Constant``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/constant.html) and -[``Polar``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/polar.html) -boundary conditions. +Boundary conditions are implemented as classes defined in +[``BoundaryCondition``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions.html). #### Stepper The logic of the MPDATA iterative solver is represented in PyMPDATA by the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class. -When instantiating the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html), the user has a choice -of either supplying just the number of dimensions or specialising the stepper for a given grid: -
-Julia code (click to expand) - -```Julia -Stepper = pyimport("PyMPDATA").Stepper - -stepper = Stepper(options=options, n_dims=2) -``` -
-
-Matlab code (click to expand) - -```Matlab -Stepper = py.importlib.import_module('PyMPDATA').Stepper; - -stepper = Stepper(pyargs(... - 'options', options, ... - 'n_dims', int32(2) ... -)); -``` -
-
-Python code (click to expand) - -```Python -from PyMPDATA import Stepper - -stepper = Stepper(options=options, n_dims=2) -``` -
-or -
-Julia code (click to expand) - -```Julia -stepper = Stepper(options=options, grid=(nx, ny)) -``` -
-
-Matlab code (click to expand) - -```Matlab -stepper = Stepper(pyargs(... - 'options', options, ... - 'grid', py.tuple({nx, ny}) ... -)); -``` -
-
-Python code (click to expand) - -```Python -stepper = Stepper(options=options, grid=(nx, ny)) -``` -
- -In the latter case, noticeably -faster execution can be expected, however the resultant -stepper is less versatile as bound to the given grid size. -If number of dimensions is supplied only, the integration -might take longer, yet same instance of the -stepper can be used for different grids. - -Since creating an instance of the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class -involves time-consuming compilation of the algorithm code, -the class is equipped with a cache logic - subsequent -calls with same arguments return references to previously -instantiated objects. Instances of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) contain no -mutable data and are (thread-)safe to be reused. - -The init method of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) has an optional -``non_unit_g_factor`` argument which is a Boolean flag -enabling handling of the G factor term which can be used to -represent coordinate transformations and/or variable fluid density. - -Optionally, the number of threads to use for domain decomposition -in the first (non-contiguous) dimension during 2D and 3D calculations -may be specified using the optional ``n_threads`` argument with a -default value of ``numba.get_num_threads()``. The multi-threaded -logic of PyMPDATA depends thus on settings of numba, namely on the -selected threading layer (either via ``NUMBA_THREADING_LAYER`` env -var or via ``numba.config.THREADING_LAYER``) and the selected size of the -thread pool (``NUMBA_NUM_THREADS`` env var or ``numba.config.NUMBA_NUM_THREADS``). - - #### Solver Instances of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class are used to control the integration and access solution data. During instantiation, additional memory required by the solver is -allocated according to the options provided. - -The only method of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class besides the -init is [``advance(n_steps, mu_coeff, ...)``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html#Solver.advance) -which advances the solution by ``n_steps`` timesteps, optionally -taking into account a given diffusion coefficient ``mu_coeff``. - -Solution state is accessible through the ``Solver.advectee`` property. -Multiple solver[s] can share a single stepper, e.g., as exemplified in the shallow-water system solution in the examples package. - -Continuing with the above code snippets, instantiating -a solver and making 75 integration steps looks as follows: -
-Julia code (click to expand) - -```Julia -Solver = pyimport("PyMPDATA").Solver -solver = Solver(stepper=stepper, advectee=advectee, advector=advector) -solver.advance(n_steps=75) -state = solver.advectee.get() -``` -
-
-Matlab code (click to expand) - -```Matlab -Solver = py.importlib.import_module('PyMPDATA').Solver; -solver = Solver(pyargs('stepper', stepper, 'advectee', advectee, 'advector', advector)); -solver.advance(pyargs('n_steps', 75)); -state = solver.advectee.get(); -``` -
-
-Python code (click to expand) - -```Python -from PyMPDATA import Solver - -solver = Solver(stepper=stepper, advectee=advectee, advector=advector) -state_0 = solver.advectee.get().copy() -solver.advance(n_steps=75) -state = solver.advectee.get() -``` -
- -Now let's plot the results using `matplotlib` roughly as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379): - -
-Python code (click to expand) - -```Python -def plot(psi, zlim, norm=None): - xi, yi = np.indices(psi.shape) - fig, ax = pyplot.subplots(subplot_kw={"projection": "3d"}) - pyplot.gca().plot_wireframe( - xi+.5, yi+.5, - psi, color='red', linewidth=.5 - ) - ax.set_zlim(zlim) - for axis in (ax.xaxis, ax.yaxis, ax.zaxis): - axis.pane.fill = False - axis.pane.set_edgecolor('black') - axis.pane.set_alpha(1) - ax.grid(False) - ax.set_zticks([]) - ax.set_xlabel('x/dx') - ax.set_ylabel('y/dy') - ax.set_proj_type('ortho') - cnt = ax.contourf(xi+.5, yi+.5, psi, zdir='z', offset=-1, norm=norm) - cbar = pyplot.colorbar(cnt, pad=.1, aspect=10, fraction=.04) - return cbar.norm - -zlim = (-1, 1) -norm = plot(state_0, zlim) -pyplot.savefig('readme_gauss_0.png') -plot(state, zlim, norm) -pyplot.savefig('readme_gauss.png') -``` -
- -![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss_0.png) -![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss.png) - -#### Debugging - -PyMPDATA relies heavily on Numba to provide high-performance -number crunching operations. Arguably, one of the key advantage -of embracing Numba is that it can be easily switched off. This -brings multiple-order-of-magnitude drop in performance, yet -it also make the entire code of the library susceptible to -interactive debugging, one way of enabling it is by setting the -following environment variable before importing PyMPDATA: -
-Julia code (click to expand) - -```Julia -ENV["NUMBA_DISABLE_JIT"] = "1" -``` -
-
-Matlab code (click to expand) - -```Matlab -setenv('NUMBA_DISABLE_JIT', '1'); -``` -
-
-Python code (click to expand) - -```Python -import os -os.environ["NUMBA_DISABLE_JIT"] = "1" -``` -
+allocated according to the options provided. ## Contributing, reporting issues, seeking support From 4cf0a92caf906cb3d481be9662e404c89eef6d9b Mon Sep 17 00:00:00 2001 From: Sylwester Arabas Date: Thu, 26 Sep 2024 11:12:12 +0200 Subject: [PATCH 06/10] pair-programming on the docs --- PyMPDATA/__init__.py | 1 - README.md | 34 +++----------- docs/markdown/debugging_contribution.md | 36 +-------------- docs/markdown/dependencies.md | 9 ---- docs/markdown/grid.md | 4 +- docs/markdown/intro.md | 23 +--------- docs/markdown/options.md | 6 ++- docs/markdown/stepper_solver.md | 4 +- docs/templates/index.html.jinja2 | 60 ++++++++++++++----------- 9 files changed, 51 insertions(+), 126 deletions(-) delete mode 100644 docs/markdown/dependencies.md diff --git a/PyMPDATA/__init__.py b/PyMPDATA/__init__.py index 888223e2..ad33444a 100644 --- a/PyMPDATA/__init__.py +++ b/PyMPDATA/__init__.py @@ -7,7 +7,6 @@ `PyMPDATA.vector_field.VectorField` elements: .. include:: ../docs/markdown/intro.md -.. include:: ../docs/markdown/dependencies.md .. include:: ../docs/markdown/examples.md .. include:: ../docs/markdown/options.md .. include:: ../docs/markdown/grid.md diff --git a/README.md b/README.md index a72762f9..3fe3f39e 100644 --- a/README.md +++ b/README.md @@ -51,9 +51,8 @@ To install PyMPDATA, one may use: ``pip install PyMPDATA`` (or ``pip install git+https://github.com/open-atmos/PyMPDATA.git`` to get updates beyond the latest release). PyMPDATA depends on ``NumPy`` and ``Numba``. -Running the tests shipped with the package requires additional packages listed in the -[test-time-requirements.txt](https://github.com/open-atmos/PyMPDATA/blob/main/test-time-requirements.txt) file -(which include ``PyMPDATA-examples``, see below). +Running the tests shipped with the package requires additional packages that are installed +if pip is invoked with: ``pip install PyMPDATA[tests]``. ## Examples (Jupyter notebooks reproducing results from literature): @@ -115,7 +114,8 @@ from the continuous integration setup which automates execution of tests with the newly added code. As of now, the copyright to the entire PyMPDATA codebase is with the Jagiellonian -University, and code contributions are assumed to imply transfer of copyright. +University (2019-2023) and AGH University of Krakow (2023 onwards) - work places of the main maintainer. +Code contributions are assumed to imply transfer of copyright. Should there be a need to make an exception, please indicate it when creating a pull request or contributing code in any other way. In any case, the license of the contributed code must be compatible with GPL v3. @@ -137,29 +137,9 @@ Please use the PyMPDATA issue-tracking and dicsussion infrastructure for `PyMPDA We look forward to your contributions and feedback. ## Credits: -Development of PyMPDATA was supported by the EU through a grant of the [Foundation for Polish Science](http://fnp.org.pl) (POIR.04.04.00-00-5E1C/18). +Development of PyMPDATA was supported by the EU through a grant of the [Foundation for Polish Science](http://fnp.org.pl) (POIR.04.04.00-00-5E1C/18) +and by the [Polish National Science Centre](https://ncn.gov.pl/en) (grant no. 2020/39/D/ST10/01220) -copyright: Jagiellonian University +copyright: Jagiellonian University (2019-2023) & AGH University of Krakow (2023 onwards) licence: GPL v3 -## Other open-source MPDATA implementations: -- mpdat_2d in babyEULAG (FORTRAN) - https://github.com/igfuw/bE_SDs/blob/master/babyEULAG.SDs.for#L741 -- mpdata-oop (C++, Fortran, Python) - https://github.com/igfuw/mpdata-oop -- apc-llc/mpdata (C++) - https://github.com/apc-llc/mpdata -- libmpdata++ (C++): - https://github.com/igfuw/libmpdataxx -- AtmosFOAM: - https://github.com/AtmosFOAM/AtmosFOAM/tree/947b192f69d973ea4a7cfab077eb5c6c6fa8b0cf/applications/solvers/advection/MPDATAadvectionFoam - -## Other Python packages for solving hyperbolic transport equations - -- PyPDE: https://pypi.org/project/PyPDE/ -- FiPy: https://pypi.org/project/FiPy/ -- ader: https://pypi.org/project/ader/ -- centpy: https://pypi.org/project/centpy/ -- mattflow: https://pypi.org/project/mattflow/ -- FastFD: https://pypi.org/project/FastFD/ -- Pyclaw: https://www.clawpack.org/pyclaw/ diff --git a/docs/markdown/debugging_contribution.md b/docs/markdown/debugging_contribution.md index fb0f2f6a..f15b392a 100644 --- a/docs/markdown/debugging_contribution.md +++ b/docs/markdown/debugging_contribution.md @@ -32,39 +32,7 @@ os.environ["NUMBA_DISABLE_JIT"] = "1" ## Contributing, reporting issues, seeking support -Submitting new code to the project, please preferably use [GitHub pull requests](https://github.com/open-atmos/PyMPDATA/pulls) -(or the [PyMPDATA-examples PR site](https://github.com/open-atmos/PyMPDATA-examples/pulls) if working on examples) - it helps to keep record of code authorship, -track and archive the code review workflow and allows to benefit -from the continuous integration setup which automates execution of tests -with the newly added code. - -As of now, the copyright to the entire PyMPDATA codebase is with the Jagiellonian -University, and code contributions are assumed to imply transfer of copyright. -Should there be a need to make an exception, please indicate it when creating -a pull request or contributing code in any other way. In any case, -the license of the contributed code must be compatible with GPL v3. - -Developing the code, we follow [The Way of Python](https://www.python.org/dev/peps/pep-0020/) and -the [KISS principle](https://en.wikipedia.org/wiki/KISS_principle). -The codebase has greatly benefited from [PyCharm code inspections](https://www.jetbrains.com/help/pycharm/code-inspection.html) -and [Pylint](https://pylint.org) code analysis (Pylint checks are part of the -CI workflows). - -Issues regarding any incorrect, unintuitive or undocumented bahaviour of -PyMPDATA are best to be reported on the [GitHub issue tracker](https://github.com/open-atmos/PyMPDATA/issues/new). -Feature requests are recorded in the "Ideas..." [PyMPDATA wiki page](https://github.com/open-atmos/PyMPDATA/wiki/Ideas-for-new-features-and-examples). - -We encourage to use the [GitHub Discussions](https://github.com/open-atmos/PyMPDATA/discussions) feature -(rather than the issue tracker) for seeking support in understanding, using and extending PyMPDATA code. - -Please use the PyMPDATA issue-tracking and dicsussion infrastructure for `PyMPDATA-examples` as well. -We look forward to your contributions and feedback. - -## Credits: -Development of PyMPDATA was supported by the EU through a grant of the [Foundation for Polish Science](http://fnp.org.pl) (POIR.04.04.00-00-5E1C/18). - -copyright: Jagiellonian University -licence: GPL v3 +See [https://github.com/open-atmos/PyMPDATA/tree/main/README.md](README.md). ## Other open-source MPDATA implementations: - mpdat_2d in babyEULAG (FORTRAN) @@ -86,4 +54,4 @@ licence: GPL v3 - centpy: https://pypi.org/project/centpy/ - mattflow: https://pypi.org/project/mattflow/ - FastFD: https://pypi.org/project/FastFD/ -- Pyclaw: https://www.clawpack.org/pyclaw/ \ No newline at end of file +- Pyclaw: https://www.clawpack.org/pyclaw/ diff --git a/docs/markdown/dependencies.md b/docs/markdown/dependencies.md deleted file mode 100644 index 983d02aa..00000000 --- a/docs/markdown/dependencies.md +++ /dev/null @@ -1,9 +0,0 @@ -## Dependencies and installation - -To install PyMPDATA, one may use: ``pip install PyMPDATA`` (or -``pip install git+https://github.com/open-atmos/PyMPDATA.git`` to get updates beyond the latest release). -PyMPDATA depends on ``NumPy`` and ``Numba``. - -Running the tests shipped with the package requires additional packages listed in the -[test-time-requirements.txt](https://github.com/open-atmos/PyMPDATA/blob/main/test-time-requirements.txt) file -(which include ``PyMPDATA-examples``, see below). \ No newline at end of file diff --git a/docs/markdown/grid.md b/docs/markdown/grid.md index 925281cb..f06b555d 100644 --- a/docs/markdown/grid.md +++ b/docs/markdown/grid.md @@ -1,4 +1,4 @@ -#### Arakawa-C grid layer and boundary conditions +### Arakawa-C grid layer and boundary conditions In PyMPDATA, the solution domain is assumed to extend from the first cell's boundary to the last cell's boundary (thus the @@ -11,4 +11,4 @@ in which: - vector field components are discretised onto cell walls. The schematic of the employed grid/domain layout in two dimensions is given below -(with the Python code snippet generating the figure as a part of CI workflow): \ No newline at end of file +(with the Python code snippet generating the figure as a part of CI workflow): diff --git a/docs/markdown/intro.md b/docs/markdown/intro.md index 132a8a7f..d792b801 100644 --- a/docs/markdown/intro.md +++ b/docs/markdown/intro.md @@ -2,27 +2,6 @@ # PyMPDATA -[![Python 3](https://img.shields.io/static/v1?label=Python&logo=Python&color=3776AB&message=3)](https://www.python.org/) -[![LLVM](https://img.shields.io/static/v1?label=LLVM&logo=LLVM&color=gold&message=Numba)](https://www.numba.org) -[![Linux OK](https://img.shields.io/static/v1?label=Linux&logo=Linux&color=yellow&message=%E2%9C%93)](https://en.wikipedia.org/wiki/Linux) -[![macOS OK](https://img.shields.io/static/v1?label=macOS&logo=Apple&color=silver&message=%E2%9C%93)](https://en.wikipedia.org/wiki/macOS) -[![Windows OK](https://img.shields.io/static/v1?label=Windows&logo=Windows&color=white&message=%E2%9C%93)](https://en.wikipedia.org/wiki/Windows) -[![Jupyter](https://img.shields.io/static/v1?label=Jupyter&logo=Jupyter&color=f37626&message=%E2%9C%93)](https://jupyter.org/) -[![Maintenance](https://img.shields.io/badge/Maintained%3F-yes-green.svg)](https://github.com/open-atmos/PyMPDATA/graphs/commit-activity) -[![OpenHub](https://www.openhub.net/p/atmos-cloud-sim-uj-PyMPDATA/widgets/project_thin_badge?format=gif)](https://www.openhub.net/p/atmos-cloud-sim-uj-PyMPDATA) -[![status](https://joss.theoj.org/papers/10e7361e43785dbb1b3d659c5b01757a/status.svg)](https://joss.theoj.org/papers/10e7361e43785dbb1b3d659c5b01757a) -[![DOI](https://zenodo.org/badge/225610671.svg)](https://zenodo.org/badge/latestdoi/225610671) -[![EU Funding](https://img.shields.io/static/v1?label=EU%20Funding%20by&color=103069&message=FNP&logoWidth=25&logo=image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAC4AAAAeCAYAAABTwyyaAAAEzklEQVRYw9WYS2yUVRiGn3P5ZzozpZ3aUsrNgoKlKBINmkhpCCwwxIAhsDCpBBIWhmCMMYTEhSJ4i9EgBnSBEm81MRrFBhNXEuUSMCopiRWLQqEGLNgr085M5//POS46NNYFzHQ6qGc1i5nzP/P973m/9ztCrf7A8T9csiibCocUbvTzfxLcAcaM3cY3imXz25lT3Y34G7gQYAKV3+bFAHcATlBTPogJNADG92iY28FHW97kyPbnuW/W7xgzAhukQ9xe04PJeOT0HkQRwK0TlEeGWb/kOO9v3kdD3a8YK9GhDMfa6mg9fxunOm/lWPtcpDI4K7n/jnN8+uQbrFrUSiwU/DtSEUB/MsKKBT+zYslJqiYNgVE4JwhHkzy86wlWvrKVWDSZ/YFjZlU39yw4y/rGoyQGowWB67zl4QQue+jssMdXrQvZ/00jyeHwqCgDKwnsiJjSvkYAxsG5K9WsenYbJdqAtAjhCIxCSZt/4fK1w5A2WCvxrUAKCHwNVoA2aGmvq11jJQQapEXrgMBKqmJJugejKGWLIxXrBPFoigfv/omd675gRkU/xgqUDlAhH3UDaAAlLSqUQekAYyVTyhLs3tDMsvntlIYzOFcEcOcEGd9jx9oDbGs6QO0t/Tijxi9S4bhzxiWaVh5m94Zm0n7oui4ybo0raUlcncQnxx+g+WgDF/vLoYDmoqSl/dJUnt7XRCoTZjij0Z6Pc2LiNS4EBBkNvoeOJXN+yPWWSZeANOhwJq/98nKVwNdoL8B5AROxBKBL0gjh8DMhdCh3eJnrA0yqhLpplwmyup6IajvAOIGfKGVx3VmCRGnOMpe5QAdG0bT8CAeeep0d6z6nqjSJnQiZWEllLMWrmz6k+fE9rGk8MVqYgsGv5ZH2i1Opr+9kajzB5d74hKQ+KS3d/WVMLhtgdu1lriRiOR/4nDVunaR24x7qp3UV5Cb/fJvC83nv26W81LIK58SYNFmwq4hsGx/5BwKlzYRma2NUthgOJSew4i7ru9nJYCQF5tApb2yvjiDQKJV/IfJKh0o6qssSLKv/jcAoRKHQQzE2Lj2OMV5OkWFc4MZIpsev8uXWXRx6ZicbGk8QZLxxgwe+x/rlR3h3816+f2E7lbEU+ZDn3vKVpePCdFovzCISHqbl5EIoQOteKMPB1rto65zNyfOz+KOrGl06lHPQyi/WOohH0/T0l1MZH6A3GUEKl7Pmr2la6wBrBWWRDP2DUcqjKVKBGom9RZmABAykwnglafpSJSPQvsfiOR0EQ7ExVmazA8cY6N4K1iw6RdAXRwi4mgrheT5Dvs4LeuS81a15Ll/3dQisFVSVpnj7sf1sX/sZvhAc+6UOrQyBVUQ8gx/orFmDsZqtaw/y1qZ9zKjp5vDpenyjcNe+cLNmTiUdf/bEOddVQ0VpgsOn54ET+EYxvWKALSu+5tGG76it7MNaiZKGQ23zCIcMfUMxBnrjN3fmHHvCAlp+vJcXWx6itqoXpAEnUNLx8iMfo5Xh1i17R3PJYCpC2cZ3qK3sQ8WGEDDuXlAQuFKGHzpmopXhTNfk0bmxs7uC1w6uJul79AxFkMIiBJy5UoUWjrZLU5DCFdTARDHuDqVw+OkSwI0MCEW4gtNF2BPrBCo8fKNbtILWX9aUDqFqHnn7AAAAAElFTkSuQmCC)](https://www.fnp.org.pl/en/) -[![PL Funding](https://img.shields.io/static/v1?label=PL%20Funding%20by&color=d21132&message=NCN&logoWidth=25&logo=image/png;base64,iVBORw0KGgoAAAANSUhEUgAAABQAAAANCAYAAACpUE5eAAAABmJLR0QA/wD/AP+gvaeTAAAAKUlEQVQ4jWP8////fwYqAiZqGjZqIHUAy4dJS6lqIOMdEZvRZDPcDQQAb3cIaY1Sbi4AAAAASUVORK5CYII=)](https://www.ncn.gov.pl/?language=en) -[![License: GPL v3](https://img.shields.io/badge/License-GPL%20v3-blue.svg)](https://www.gnu.org/licenses/gpl-3.0.html) - -[![Github Actions Build Status](https://github.com/open-atmos/PyMPDATA/actions/workflows/tests+pypi.yml/badge.svg?branch=main)](https://github.com/open-atmos/PyMPDATA/actions) -[![Appveyor Build status](http://ci.appveyor.com/api/projects/status/github/open-atmos/PyMPDATA?branch=main&svg=true)](https://ci.appveyor.com/project/slayoo/pympdata/branch/main) -[![Coverage Status](https://codecov.io/gh/open-atmos/PyMPDATA/branch/main/graph/badge.svg)](https://app.codecov.io/gh/open-atmos/PyMPDATA) - -[![PyPI version](https://badge.fury.io/py/PyMPDATA.svg)](https://pypi.org/project/PyMPDATA) -[![API docs](https://shields.mitmproxy.org/badge/docs-pdoc.dev-brightgreen.svg)](https://open-atmos.github.io/PyMPDATA/index.html) - PyMPDATA is a high-performance Numba-accelerated Pythonic implementation of the MPDATA algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond. MPDATA numerically solves generalised transport equations - @@ -74,4 +53,4 @@ The [``PyMPDATA-examples``](https://pypi.org/project/PyMPDATA-examples/) or using [colab.research.google.com](https://colab.research.google.com/). The examples reproduce results from several published works on MPDATA and its applications, and provide a validation of the implementation - and its performance. \ No newline at end of file + and its performance. diff --git a/docs/markdown/options.md b/docs/markdown/options.md index 4c991179..31d8526b 100644 --- a/docs/markdown/options.md +++ b/docs/markdown/options.md @@ -1,4 +1,6 @@ -## Options class +## Tutorial (in Python, Julia and Matlab) + +### Options class The [``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class groups both algorithm variant options as well as some implementation-related @@ -52,4 +54,4 @@ options = Options(pyargs('n_iters', 2)); from PyMPDATA import Options options = Options(n_iters=2) ``` - \ No newline at end of file + diff --git a/docs/markdown/stepper_solver.md b/docs/markdown/stepper_solver.md index 7e93972e..fc9a9551 100644 --- a/docs/markdown/stepper_solver.md +++ b/docs/markdown/stepper_solver.md @@ -1,4 +1,4 @@ -#### Stepper +### Stepper The logic of the MPDATA iterative solver is represented in PyMPDATA by the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class. @@ -140,4 +140,4 @@ state = solver.advectee.get() ``` -Now let's plot the results using `matplotlib` roughly as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379): \ No newline at end of file +Now let's plot the results using `matplotlib` roughly as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379): diff --git a/docs/templates/index.html.jinja2 b/docs/templates/index.html.jinja2 index 9d8c9ad4..e6e03856 100644 --- a/docs/templates/index.html.jinja2 +++ b/docs/templates/index.html.jinja2 @@ -3,28 +3,51 @@ {% block title %}PyMPDATA module list{% endblock %} {% block nav %} -

Available Modules

-
    - {% for submodule in all_modules if "." not in submodule %} -
  • {{ submodule }}
    • - {% endfor %} -
{% endblock %} + {% block content %}
- - pympdata logo - +
+ + pympdata logo + +

What is PyMPDATA?

PyMPDATA is a Numba-accelerated Pythonic implementation of the MPDATA algorithm - of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond for + of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond for numerically solving generalised convection-diffusion PDEs in 1D, 2D and 3D structured meshes with coordinate transformations.

+ +

PyMPDATA resources

+ + + + + + + + + + + +
PyMPDATAPyMPDATA-examplesPyMPDATA-MPI
+ + + + +
+

What is the difference between PyMPDATA and PyMPDATA-examples?

@@ -79,23 +102,6 @@ using and extending PyMPDATA code.

-
-

Credits and License

-

- Development of PyMPDATA was supported by the EU through a grant of the - Foundation for Polish Science (POIR.04.04.00-00-5E1C/18). -

-

- licence: GPL v3 -

-

- As of now, the copyright to the entire PyMPDATA codebase is with the Jagiellonian - University, and code contributions are assumed to imply transfer of copyright. - Should there be a need to make an exception, please indicate it when creating - a pull request or contributing code in any other way. In any case, - the license of the contributed code must be compatible with GPL v3. -

-
{% include "search.html.jinja2" %} From 518b09d9844f152c683595583a7d9344f2797ccf Mon Sep 17 00:00:00 2001 From: pawel Date: Fri, 27 Sep 2024 23:26:54 +0200 Subject: [PATCH 07/10] merge markdown files for pympdata landing page --- PyMPDATA/__init__.py | 22 +- docs/markdown/debugging_contribution.md | 57 --- docs/markdown/examples.md | 14 - docs/markdown/grid.md | 14 - docs/markdown/grid_2.md | 19 - docs/markdown/grid_3.md | 10 - docs/markdown/intro.md | 56 --- docs/markdown/julia_grid.md | 27 -- docs/markdown/matlab_grid.md | 35 -- docs/markdown/options.md | 57 --- docs/markdown/package_structure.md | 16 - docs/markdown/pympdata_landing.md | 514 ++++++++++++++++++++++++ docs/markdown/python_grid.md | 38 -- docs/markdown/python_grid_2.md | 30 -- docs/markdown/python_solver.md | 35 -- docs/markdown/stepper_solver.md | 143 ------- 16 files changed, 515 insertions(+), 572 deletions(-) delete mode 100644 docs/markdown/debugging_contribution.md delete mode 100644 docs/markdown/examples.md delete mode 100644 docs/markdown/grid.md delete mode 100644 docs/markdown/grid_2.md delete mode 100644 docs/markdown/grid_3.md delete mode 100644 docs/markdown/intro.md delete mode 100644 docs/markdown/julia_grid.md delete mode 100644 docs/markdown/matlab_grid.md delete mode 100644 docs/markdown/options.md delete mode 100644 docs/markdown/package_structure.md create mode 100644 docs/markdown/pympdata_landing.md delete mode 100644 docs/markdown/python_grid.md delete mode 100644 docs/markdown/python_grid_2.md delete mode 100644 docs/markdown/python_solver.md delete mode 100644 docs/markdown/stepper_solver.md diff --git a/PyMPDATA/__init__.py b/PyMPDATA/__init__.py index ad33444a..d88924a8 100644 --- a/PyMPDATA/__init__.py +++ b/PyMPDATA/__init__.py @@ -1,25 +1,5 @@ """ -Numba-accelerated Pythonic implementation of Multidimensional Positive Definite -Advection Transport Algorithm (MPDATA) with examples in Python, Julia and Matlab - -PyMPDATA uses staggered grid with the following node placement for -`PyMPDATA.scalar_field.ScalarField` and -`PyMPDATA.vector_field.VectorField` elements: - -.. include:: ../docs/markdown/intro.md -.. include:: ../docs/markdown/examples.md -.. include:: ../docs/markdown/options.md -.. include:: ../docs/markdown/grid.md -.. include:: ../docs/markdown/python_grid.md -.. include:: ../docs/markdown/grid_2.md -.. include:: ../docs/markdown/julia_grid.md -.. include:: ../docs/markdown/matlab_grid.md -.. include:: ../docs/markdown/python_grid_2.md -.. include:: ../docs/markdown/grid_3.md -.. include:: ../docs/markdown/stepper_solver.md -.. include:: ../docs/markdown/python_solver.md -.. include:: ../docs/markdown/debugging_contribution.md - +.. include:: ../docs/markdown/pympdata_landing.md """ # pylint: disable=invalid-name diff --git a/docs/markdown/debugging_contribution.md b/docs/markdown/debugging_contribution.md deleted file mode 100644 index f15b392a..00000000 --- a/docs/markdown/debugging_contribution.md +++ /dev/null @@ -1,57 +0,0 @@ -#### Debugging - -PyMPDATA relies heavily on Numba to provide high-performance -number crunching operations. Arguably, one of the key advantage -of embracing Numba is that it can be easily switched off. This -brings multiple-order-of-magnitude drop in performance, yet -it also make the entire code of the library susceptible to -interactive debugging, one way of enabling it is by setting the -following environment variable before importing PyMPDATA: -
-Julia code (click to expand) - -```Julia -ENV["NUMBA_DISABLE_JIT"] = "1" -``` -
-
-Matlab code (click to expand) - -```Matlab -setenv('NUMBA_DISABLE_JIT', '1'); -``` -
-
-Python code (click to expand) - -```Python -import os -os.environ["NUMBA_DISABLE_JIT"] = "1" -``` -
- -## Contributing, reporting issues, seeking support - -See [https://github.com/open-atmos/PyMPDATA/tree/main/README.md](README.md). - -## Other open-source MPDATA implementations: -- mpdat_2d in babyEULAG (FORTRAN) - https://github.com/igfuw/bE_SDs/blob/master/babyEULAG.SDs.for#L741 -- mpdata-oop (C++, Fortran, Python) - https://github.com/igfuw/mpdata-oop -- apc-llc/mpdata (C++) - https://github.com/apc-llc/mpdata -- libmpdata++ (C++): - https://github.com/igfuw/libmpdataxx -- AtmosFOAM: - https://github.com/AtmosFOAM/AtmosFOAM/tree/947b192f69d973ea4a7cfab077eb5c6c6fa8b0cf/applications/solvers/advection/MPDATAadvectionFoam - -## Other Python packages for solving hyperbolic transport equations - -- PyPDE: https://pypi.org/project/PyPDE/ -- FiPy: https://pypi.org/project/FiPy/ -- ader: https://pypi.org/project/ader/ -- centpy: https://pypi.org/project/centpy/ -- mattflow: https://pypi.org/project/mattflow/ -- FastFD: https://pypi.org/project/FastFD/ -- Pyclaw: https://www.clawpack.org/pyclaw/ diff --git a/docs/markdown/examples.md b/docs/markdown/examples.md deleted file mode 100644 index bf481719..00000000 --- a/docs/markdown/examples.md +++ /dev/null @@ -1,14 +0,0 @@ -## Examples (Jupyter notebooks reproducing results from literature): - -PyMPDATA examples are bundled with PyMPDATA and located in the `examples` subfolder. -They constitute a separate [``PyMPDATA_examples``](https://pypi.org/p/PyMPDATA-examples) Python package which is also available at PyPI. -The examples have additional dependencies listed in [``PyMPDATA_examples`` package ``setup.py``](https://github.com/open-atmos/PyMPDATA/blob/main/examples/setup.py) file. -Running the examples requires the ``PyMPDATA_examples`` package to be installed. -Since the examples package includes Jupyter notebooks (and their execution requires write access), the suggested install and launch steps are: -``` -git clone https://github.com/open-atmos/PyMPDATA-examples.git -cd PyMPDATA-examples -pip install -e . -jupyter-notebook -``` -Alternatively, one can also install the examples package from pypi.org by using ``pip install PyMPDATA-examples``. diff --git a/docs/markdown/grid.md b/docs/markdown/grid.md deleted file mode 100644 index f06b555d..00000000 --- a/docs/markdown/grid.md +++ /dev/null @@ -1,14 +0,0 @@ -### Arakawa-C grid layer and boundary conditions - -In PyMPDATA, the solution domain is assumed to extend from the -first cell's boundary to the last cell's boundary (thus the -first scalar field value is at $\[\Delta x/2, \Delta y/2\]$. -The [``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) -and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) classes implement the -[Arakawa-C staggered grid](https://en.wikipedia.org/wiki/Arakawa_grids#Arakawa_C-grid) logic -in which: -- scalar fields are discretised onto cell centres (one value per cell), -- vector field components are discretised onto cell walls. - -The schematic of the employed grid/domain layout in two dimensions is given below -(with the Python code snippet generating the figure as a part of CI workflow): diff --git a/docs/markdown/grid_2.md b/docs/markdown/grid_2.md deleted file mode 100644 index c7a7bd72..00000000 --- a/docs/markdown/grid_2.md +++ /dev/null @@ -1,19 +0,0 @@ -The ``__init__`` methods of -[``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) -and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) -have the following signatures: -- [``ScalarField(data: np.ndarray, halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/scalar_field.py) -- [``VectorField(data: Tuple[np.ndarray, ...], halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/vector_field.py) -The ``data`` parameters are expected to be Numpy arrays or tuples of Numpy arrays, respectively. -The ``halo`` parameter is the extent of ghost-cell region that will surround the -data and will be used to implement boundary conditions. -Its value (in practice 1 or 2) is -dependent on maximal stencil extent for the MPDATA variant used and -can be easily obtained using the ``Options.n_halo`` property. - -As an example, the code below shows how to instantiate a scalar -and a vector field given a 2D constant-velocity problem, -using a grid of 24x24 points, Courant numbers of -0.5 and -0.25 -in "x" and "y" directions, respectively, with periodic boundary -conditions and with an initial Gaussian signal in the scalar field -(settings as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379)): \ No newline at end of file diff --git a/docs/markdown/grid_3.md b/docs/markdown/grid_3.md deleted file mode 100644 index fd5ca36b..00000000 --- a/docs/markdown/grid_3.md +++ /dev/null @@ -1,10 +0,0 @@ -Note that the shapes of arrays representing components -of the velocity field are different than the shape of -the scalar field array due to employment of the staggered grid. - -Besides the exemplified [``Periodic``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/periodic.html) class representing -periodic boundary conditions, PyMPDATA supports -[``Extrapolated``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/extrapolated.html), -[``Constant``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/constant.html) and -[``Polar``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/polar.html) -boundary conditions. \ No newline at end of file diff --git a/docs/markdown/intro.md b/docs/markdown/intro.md deleted file mode 100644 index d792b801..00000000 --- a/docs/markdown/intro.md +++ /dev/null @@ -1,56 +0,0 @@ -pympdata logo - -# PyMPDATA - -PyMPDATA is a high-performance Numba-accelerated Pythonic implementation of the MPDATA - algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond. -MPDATA numerically solves generalised transport equations - - partial differential equations used to model conservation/balance laws, scalar-transport problems, - convection-diffusion phenomena. -As of the current version, PyMPDATA supports homogeneous transport - in 1D, 2D and 3D using structured meshes, optionally - generalised by employment of a Jacobian of coordinate transformation. -PyMPDATA includes implementation of a set of MPDATA variants including - the non-oscillatory option, infinite-gauge, divergent-flow, double-pass donor cell (DPDC) and - third-order-terms options. -It also features support for integration of Fickian-terms in advection-diffusion - problems using the pseudo-transport velocity approach. -In 2D and 3D simulations, domain-decomposition is used for multi-threaded parallelism. - -PyMPDATA is engineered purely in Python targeting both performance and usability, - the latter encompassing research users', developers' and maintainers' perspectives. -From researcher's perspective, PyMPDATA offers hassle-free installation on multitude - of platforms including Linux, OSX and Windows, and eliminates compilation stage - from the perspective of the user. -From developers' and maintainers' perspective, PyMPDATA offers a suite of unit tests, - multi-platform continuous integration setup, - seamless integration with Python development aids including debuggers and profilers. - -PyMPDATA design features - a custom-built multi-dimensional Arakawa-C grid layer allowing - to concisely represent multi-dimensional stencil operations on both - scalar and vector fields. -The grid layer is built on top of NumPy's ndarrays (using "C" ordering) - using the Numba's ``@njit`` functionality for high-performance array traversals. -It enables one to code once for multiple dimensions, and automatically - handles (and hides from the user) any halo-filling logic related with boundary conditions. -Numba ``prange()`` functionality is used for implementing multi-threading - (it offers analogous functionality to OpenMP parallel loop execution directives). -The Numba's deviation from Python semantics rendering [closure variables - as compile-time constants](https://numba.pydata.org/numba-doc/dev/reference/pysemantics.html#global-and-closure-variables) - is extensively exploited within ``PyMPDATA`` - code base enabling the just-in-time compilation to benefit from - information on domain extents, algorithm variant used and problem - characteristics (e.g., coordinate transformation used, or lack thereof). - -A separate project called [``PyMPDATA-MPI``](https://github.com/open-atmos/PyMPDATA-MPI) - depicts how [``numba-mpi``](https://pypi.org/project/numba-mpi) can be used - to enable distributed memory parallelism in PyMPDATA. - -The [``PyMPDATA-examples``](https://pypi.org/project/PyMPDATA-examples/) - package covers a set of examples presented in the form of Jupyer notebooks - offering single-click deployment in the cloud using [mybinder.org](https://mybinder.org) - or using [colab.research.google.com](https://colab.research.google.com/). -The examples reproduce results from several published - works on MPDATA and its applications, and provide a validation of the implementation - and its performance. diff --git a/docs/markdown/julia_grid.md b/docs/markdown/julia_grid.md deleted file mode 100644 index a45aefed..00000000 --- a/docs/markdown/julia_grid.md +++ /dev/null @@ -1,27 +0,0 @@ -
-Julia code (click to expand) - -```Julia -ScalarField = pyimport("PyMPDATA").ScalarField -VectorField = pyimport("PyMPDATA").VectorField -Periodic = pyimport("PyMPDATA.boundary_conditions").Periodic - -nx, ny = 24, 24 -Cx, Cy = -.5, -.25 -idx = CartesianIndices((nx, ny)) -halo = options.n_halo -advectee = ScalarField( - data=exp.( - -(getindex.(idx, 1) .- .5 .- nx/2).^2 / (2*(nx/10)^2) - -(getindex.(idx, 2) .- .5 .- ny/2).^2 / (2*(ny/10)^2) - ), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -advector = VectorField( - data=(fill(Cx, (nx+1, ny)), fill(Cy, (nx, ny+1))), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -``` -
\ No newline at end of file diff --git a/docs/markdown/matlab_grid.md b/docs/markdown/matlab_grid.md deleted file mode 100644 index 1447250f..00000000 --- a/docs/markdown/matlab_grid.md +++ /dev/null @@ -1,35 +0,0 @@ -
-Matlab code (click to expand) - -```Matlab -ScalarField = py.importlib.import_module('PyMPDATA').ScalarField; -VectorField = py.importlib.import_module('PyMPDATA').VectorField; -Periodic = py.importlib.import_module('PyMPDATA.boundary_conditions').Periodic; - -nx = int32(24); -ny = int32(24); - -Cx = -.5; -Cy = -.25; - -[xi, yi] = meshgrid(double(0:1:nx-1), double(0:1:ny-1)); - -halo = options.n_halo; -advectee = ScalarField(pyargs(... - 'data', py.numpy.array(exp( ... - -(xi+.5-double(nx)/2).^2 / (2*(double(nx)/10)^2) ... - -(yi+.5-double(ny)/2).^2 / (2*(double(ny)/10)^2) ... - )), ... - 'halo', halo, ... - 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... -)); -advector = VectorField(pyargs(... - 'data', py.tuple({ ... - Cx * py.numpy.ones(int32([nx+1 ny])), ... - Cy * py.numpy.ones(int32([nx ny+1])) ... - }), ... - 'halo', halo, ... - 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... -)); -``` -
\ No newline at end of file diff --git a/docs/markdown/options.md b/docs/markdown/options.md deleted file mode 100644 index 31d8526b..00000000 --- a/docs/markdown/options.md +++ /dev/null @@ -1,57 +0,0 @@ -## Tutorial (in Python, Julia and Matlab) - -### Options class - -The [``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class -groups both algorithm variant options as well as some implementation-related -flags that need to be set at the first place. All are set at the time -of instantiation using the following keyword arguments of the constructor -(all having default values indicated below): -- ``n_iters: int = 2``: number of iterations (2 means upwind + one corrective iteration) -- ``infinite_gauge: bool = False``: flag enabling the infinite-gauge option (does not maintain sign of the advected field, thus in practice implies switching flux corrected transport on) -- ``divergent_flow: bool = False``: flag enabling divergent-flow terms when calculating antidiffusive velocity -- ``nonoscillatory: bool = False``: flag enabling the non-oscillatory or monotone variant (a.k.a flux-corrected transport option, FCT) -- ``third_order_terms: bool = False``: flag enabling third-order terms -- ``epsilon: float = 1e-15``: value added to potentially zero-valued denominators -- ``non_zero_mu_coeff: bool = False``: flag indicating if code for handling the Fickian term is to be optimised out -- ``DPDC: bool = False``: flag enabling double-pass donor cell option (recursive pseudovelocities) -- ``dimensionally_split: bool = False``: flag disabling cross-dimensional terms in antidiffusive velocity -- ``dtype: np.floating = np.float64``: floating point precision - -For a discussion of the above options, see e.g., [Smolarkiewicz & Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), -[Jaruga, Arabas et al. 2015](https://doi.org/10.5194/gmd-8-1005-2015) and [Olesik, Arabas et al. 2020](https://arxiv.org/abs/2011.14726) -(the last with examples using PyMPDATA). - -In most use cases of PyMPDATA, the first thing to do is to instantiate the -[``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class -with arguments suiting the problem at hand, e.g.: - -
-Julia code (click to expand) - -```Julia -using Pkg -Pkg.add("PyCall") -using PyCall -Options = pyimport("PyMPDATA").Options -options = Options(n_iters=2) -``` -
- -
-Matlab code (click to expand) - -```Matlab -Options = py.importlib.import_module('PyMPDATA').Options; -options = Options(pyargs('n_iters', 2)); -``` -
- -
-Python code (click to expand) - -```Python -from PyMPDATA import Options -options = Options(n_iters=2) -``` -
diff --git a/docs/markdown/package_structure.md b/docs/markdown/package_structure.md deleted file mode 100644 index d8c04e4f..00000000 --- a/docs/markdown/package_structure.md +++ /dev/null @@ -1,16 +0,0 @@ -## Package structure and API: - -In short, PyMPDATA numerically solves the following equation: - -$$ \partial_t (G \psi) + \nabla \cdot (Gu \psi) + \mu \Delta (G \psi) = 0 $$ - -where scalar field $\psi$ is referred to as the advectee, -vector field u is referred to as advector, and the G factor corresponds to optional coordinate transformation. -The inclusion of the Fickian diffusion term is optional and is realised through modification of the -advective velocity field with MPDATA handling both the advection and diffusion (for discussion -see, e.g. [Smolarkiewicz and Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), sec. 3.5, par. 4). - -The key classes constituting the PyMPDATA interface are summarised below with code -snippets exemplifying usage of PyMPDATA from Python, Julia and Matlab. - -A [pdoc-generated](https://pdoc.dev/) documentation of PyMPDATA public API is maintained at: [https://open-atmos.github.io/PyMPDATA]( https://open-atmos.github.io/PyMPDATA/index.html) diff --git a/docs/markdown/pympdata_landing.md b/docs/markdown/pympdata_landing.md new file mode 100644 index 00000000..d48ee473 --- /dev/null +++ b/docs/markdown/pympdata_landing.md @@ -0,0 +1,514 @@ +pympdata logo +# Introduction +PyMPDATA is a high-performance Numba-accelerated Pythonic implementation of the MPDATA + algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond. +MPDATA numerically solves generalised transport equations - + partial differential equations used to model conservation/balance laws, scalar-transport problems, + convection-diffusion phenomena. +As of the current version, PyMPDATA supports homogeneous transport + in 1D, 2D and 3D using structured meshes, optionally + generalised by employment of a Jacobian of coordinate transformation. +PyMPDATA includes implementation of a set of MPDATA variants including + the non-oscillatory option, infinite-gauge, divergent-flow, double-pass donor cell (DPDC) and + third-order-terms options. +It also features support for integration of Fickian-terms in advection-diffusion + problems using the pseudo-transport velocity approach. +In 2D and 3D simulations, domain-decomposition is used for multi-threaded parallelism. + +PyMPDATA is engineered purely in Python targeting both performance and usability, + the latter encompassing research users', developers' and maintainers' perspectives. +From researcher's perspective, PyMPDATA offers hassle-free installation on multitude + of platforms including Linux, OSX and Windows, and eliminates compilation stage + from the perspective of the user. +From developers' and maintainers' perspective, PyMPDATA offers a suite of unit tests, + multi-platform continuous integration setup, + seamless integration with Python development aids including debuggers and profilers. + +PyMPDATA design features + a custom-built multi-dimensional Arakawa-C grid layer allowing + to concisely represent multi-dimensional stencil operations on both + scalar and vector fields. +The grid layer is built on top of NumPy's ndarrays (using "C" ordering) + using the Numba's ``@njit`` functionality for high-performance array traversals. +It enables one to code once for multiple dimensions, and automatically + handles (and hides from the user) any halo-filling logic related with boundary conditions. +Numba ``prange()`` functionality is used for implementing multi-threading + (it offers analogous functionality to OpenMP parallel loop execution directives). +The Numba's deviation from Python semantics rendering [closure variables + as compile-time constants](https://numba.pydata.org/numba-doc/dev/reference/pysemantics.html#global-and-closure-variables) + is extensively exploited within ``PyMPDATA`` + code base enabling the just-in-time compilation to benefit from + information on domain extents, algorithm variant used and problem + characteristics (e.g., coordinate transformation used, or lack thereof). + +# Tutorial (in Python, Julia and Matlab) +## Options class + +The [``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class +groups both algorithm variant options as well as some implementation-related +flags that need to be set at the first place. All are set at the time +of instantiation using the following keyword arguments of the constructor +(all having default values indicated below): +- ``n_iters: int = 2``: number of iterations (2 means upwind + one corrective iteration) +- ``infinite_gauge: bool = False``: flag enabling the infinite-gauge option (does not maintain sign of the advected field, thus in practice implies switching flux corrected transport on) +- ``divergent_flow: bool = False``: flag enabling divergent-flow terms when calculating antidiffusive velocity +- ``nonoscillatory: bool = False``: flag enabling the non-oscillatory or monotone variant (a.k.a flux-corrected transport option, FCT) +- ``third_order_terms: bool = False``: flag enabling third-order terms +- ``epsilon: float = 1e-15``: value added to potentially zero-valued denominators +- ``non_zero_mu_coeff: bool = False``: flag indicating if code for handling the Fickian term is to be optimised out +- ``DPDC: bool = False``: flag enabling double-pass donor cell option (recursive pseudovelocities) +- ``dimensionally_split: bool = False``: flag disabling cross-dimensional terms in antidiffusive velocity +- ``dtype: np.floating = np.float64``: floating point precision + +For a discussion of the above options, see e.g., [Smolarkiewicz & Margolin 1998](https://doi.org/10.1006/jcph.1998.5901), +[Jaruga, Arabas et al. 2015](https://doi.org/10.5194/gmd-8-1005-2015) and [Olesik, Arabas et al. 2020](https://arxiv.org/abs/2011.14726) +(the last with examples using PyMPDATA). + +In most use cases of PyMPDATA, the first thing to do is to instantiate the +[``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class +with arguments suiting the problem at hand, e.g.: + +
+Julia code (click to expand) + +```Julia +using Pkg +Pkg.add("PyCall") +using PyCall +Options = pyimport("PyMPDATA").Options +options = Options(n_iters=2) +``` +
+ +
+Matlab code (click to expand) + +```Matlab +Options = py.importlib.import_module('PyMPDATA').Options; +options = Options(pyargs('n_iters', 2)); +``` +
+ +
+Python code (click to expand) + +```Python +from PyMPDATA import Options +options = Options(n_iters=2) +``` +
+ +## Arakawa-C grid layer and boundary conditions + +In PyMPDATA, the solution domain is assumed to extend from the +first cell's boundary to the last cell's boundary (thus the +first scalar field value is at $\[\Delta x/2, \Delta y/2\]$. +The [``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) +and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) classes implement the +[Arakawa-C staggered grid](https://en.wikipedia.org/wiki/Arakawa_grids#Arakawa_C-grid) logic +in which: +- scalar fields are discretised onto cell centres (one value per cell), +- vector field components are discretised onto cell walls. + +The schematic of the employed grid/domain layout in two dimensions is given below +(with the Python code snippet generating the figure as a part of CI workflow): +![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_grid.png) + +
+Python code (click to expand) + +```Python +import numpy as np +from matplotlib import pyplot + +dx, dy = .2, .3 +grid = (10, 5) + +pyplot.scatter(*np.mgrid[ + dx / 2 : grid[0] * dx : dx, + dy / 2 : grid[1] * dy : dy + ], color='red', + label='scalar-field values at cell centres' +) +pyplot.quiver(*np.mgrid[ + 0 : (grid[0]+1) * dx : dx, + dy / 2 : grid[1] * dy : dy + ], 1, 0, pivot='mid', color='green', width=.005, + label='vector-field x-component values at cell walls' +) +pyplot.quiver(*np.mgrid[ + dx / 2 : grid[0] * dx : dx, + 0: (grid[1] + 1) * dy : dy + ], 0, 1, pivot='mid', color='blue', width=.005, + label='vector-field y-component values at cell walls' +) +pyplot.xticks(np.linspace(0, grid[0]*dx, grid[0]+1)) +pyplot.yticks(np.linspace(0, grid[1]*dy, grid[1]+1)) +pyplot.title(f'staggered grid layout (grid={grid}, dx={dx}, dy={dy})') +pyplot.xlabel('x') +pyplot.ylabel('y') +pyplot.legend(bbox_to_anchor=(.1, -.1), loc='upper left', ncol=1) +pyplot.grid() +pyplot.savefig('readme_grid.png') +``` +
+ +The ``__init__`` methods of +[``ScalarField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/scalar_field.html) +and [``VectorField``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/vector_field.html) +have the following signatures: +- [``ScalarField(data: np.ndarray, halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/scalar_field.py) +- [``VectorField(data: Tuple[np.ndarray, ...], halo: int, boundary_conditions)``](https://github.com/open-atmos/PyMPDATA/blob/main/PyMPDATA/vector_field.py) +The ``data`` parameters are expected to be Numpy arrays or tuples of Numpy arrays, respectively. +The ``halo`` parameter is the extent of ghost-cell region that will surround the +data and will be used to implement boundary conditions. +Its value (in practice 1 or 2) is +dependent on maximal stencil extent for the MPDATA variant used and +can be easily obtained using the ``Options.n_halo`` property. + +As an example, the code below shows how to instantiate a scalar +and a vector field given a 2D constant-velocity problem, +using a grid of 24x24 points, Courant numbers of -0.5 and -0.25 +in "x" and "y" directions, respectively, with periodic boundary +conditions and with an initial Gaussian signal in the scalar field +(settings as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379)): + +
+Julia code (click to expand) + +```Julia +ScalarField = pyimport("PyMPDATA").ScalarField +VectorField = pyimport("PyMPDATA").VectorField +Periodic = pyimport("PyMPDATA.boundary_conditions").Periodic + +nx, ny = 24, 24 +Cx, Cy = -.5, -.25 +idx = CartesianIndices((nx, ny)) +halo = options.n_halo +advectee = ScalarField( + data=exp.( + -(getindex.(idx, 1) .- .5 .- nx/2).^2 / (2*(nx/10)^2) + -(getindex.(idx, 2) .- .5 .- ny/2).^2 / (2*(ny/10)^2) + ), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +advector = VectorField( + data=(fill(Cx, (nx+1, ny)), fill(Cy, (nx, ny+1))), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +``` +
+
+Matlab code (click to expand) + +```Matlab +ScalarField = py.importlib.import_module('PyMPDATA').ScalarField; +VectorField = py.importlib.import_module('PyMPDATA').VectorField; +Periodic = py.importlib.import_module('PyMPDATA.boundary_conditions').Periodic; + +nx = int32(24); +ny = int32(24); + +Cx = -.5; +Cy = -.25; + +[xi, yi] = meshgrid(double(0:1:nx-1), double(0:1:ny-1)); + +halo = options.n_halo; +advectee = ScalarField(pyargs(... + 'data', py.numpy.array(exp( ... + -(xi+.5-double(nx)/2).^2 / (2*(double(nx)/10)^2) ... + -(yi+.5-double(ny)/2).^2 / (2*(double(ny)/10)^2) ... + )), ... + 'halo', halo, ... + 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... +)); +advector = VectorField(pyargs(... + 'data', py.tuple({ ... + Cx * py.numpy.ones(int32([nx+1 ny])), ... + Cy * py.numpy.ones(int32([nx ny+1])) ... + }), ... + 'halo', halo, ... + 'boundary_conditions', py.tuple({Periodic(), Periodic()}) ... +)); +``` +
+
+Python code (click to expand) + +```Python +from PyMPDATA import ScalarField +from PyMPDATA import VectorField +from PyMPDATA.boundary_conditions import Periodic +import numpy as np + +nx, ny = 24, 24 +Cx, Cy = -.5, -.25 +halo = options.n_halo + +xi, yi = np.indices((nx, ny), dtype=float) +advectee = ScalarField( + data=np.exp( + -(xi+.5-nx/2)**2 / (2*(nx/10)**2) + -(yi+.5-ny/2)**2 / (2*(ny/10)**2) + ), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +advector = VectorField( + data=(np.full((nx + 1, ny), Cx), np.full((nx, ny + 1), Cy)), + halo=halo, + boundary_conditions=(Periodic(), Periodic()) +) +``` +
+ +Note that the shapes of arrays representing components +of the velocity field are different than the shape of +the scalar field array due to employment of the staggered grid. + +Besides the exemplified [``Periodic``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/periodic.html) class representing +periodic boundary conditions, PyMPDATA supports +[``Extrapolated``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/extrapolated.html), +[``Constant``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/constant.html) and +[``Polar``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/boundary_conditions/polar.html) +boundary conditions. + +## Stepper + +The logic of the MPDATA iterative solver is represented +in PyMPDATA by the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class. + +When instantiating the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html), the user has a choice +of either supplying just the number of dimensions or specialising the stepper for a given grid: +
+Julia code (click to expand) + +```Julia +Stepper = pyimport("PyMPDATA").Stepper + +stepper = Stepper(options=options, n_dims=2) +``` +
+
+Matlab code (click to expand) + +```Matlab +Stepper = py.importlib.import_module('PyMPDATA').Stepper; + +stepper = Stepper(pyargs(... + 'options', options, ... + 'n_dims', int32(2) ... +)); +``` +
+
+Python code (click to expand) + +```Python +from PyMPDATA import Stepper + +stepper = Stepper(options=options, n_dims=2) +``` +
+or +
+Julia code (click to expand) + +```Julia +stepper = Stepper(options=options, grid=(nx, ny)) +``` +
+
+Matlab code (click to expand) + +```Matlab +stepper = Stepper(pyargs(... + 'options', options, ... + 'grid', py.tuple({nx, ny}) ... +)); +``` +
+
+Python code (click to expand) + +```Python +stepper = Stepper(options=options, grid=(nx, ny)) +``` +
+ +In the latter case, noticeably +faster execution can be expected, however the resultant +stepper is less versatile as bound to the given grid size. +If number of dimensions is supplied only, the integration +might take longer, yet same instance of the +stepper can be used for different grids. + +Since creating an instance of the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class +involves time-consuming compilation of the algorithm code, +the class is equipped with a cache logic - subsequent +calls with same arguments return references to previously +instantiated objects. Instances of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) contain no +mutable data and are (thread-)safe to be reused. + +The init method of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) has an optional +``non_unit_g_factor`` argument which is a Boolean flag +enabling handling of the G factor term which can be used to +represent coordinate transformations and/or variable fluid density. + +Optionally, the number of threads to use for domain decomposition +in the first (non-contiguous) dimension during 2D and 3D calculations +may be specified using the optional ``n_threads`` argument with a +default value of ``numba.get_num_threads()``. The multi-threaded +logic of PyMPDATA depends thus on settings of numba, namely on the +selected threading layer (either via ``NUMBA_THREADING_LAYER`` env +var or via ``numba.config.THREADING_LAYER``) and the selected size of the +thread pool (``NUMBA_NUM_THREADS`` env var or ``numba.config.NUMBA_NUM_THREADS``). + + +## Solver + +Instances of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class are used to control +the integration and access solution data. During instantiation, +additional memory required by the solver is +allocated according to the options provided. + +The only method of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class besides the +init is [``advance(n_steps, mu_coeff, ...)``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html#Solver.advance) +which advances the solution by ``n_steps`` timesteps, optionally +taking into account a given diffusion coefficient ``mu_coeff``. + +Solution state is accessible through the ``Solver.advectee`` property. +Multiple solver[s] can share a single stepper, e.g., as exemplified in the shallow-water system solution in the examples package. + +Continuing with the above code snippets, instantiating +a solver and making 75 integration steps looks as follows: +
+Julia code (click to expand) + +```Julia +Solver = pyimport("PyMPDATA").Solver +solver = Solver(stepper=stepper, advectee=advectee, advector=advector) +solver.advance(n_steps=75) +state = solver.advectee.get() +``` +
+
+Matlab code (click to expand) + +```Matlab +Solver = py.importlib.import_module('PyMPDATA').Solver; +solver = Solver(pyargs('stepper', stepper, 'advectee', advectee, 'advector', advector)); +solver.advance(pyargs('n_steps', 75)); +state = solver.advectee.get(); +``` +
+
+Python code (click to expand) + +```Python +from PyMPDATA import Solver + +solver = Solver(stepper=stepper, advectee=advectee, advector=advector) +state_0 = solver.advectee.get().copy() +solver.advance(n_steps=75) +state = solver.advectee.get() +``` +
+ +Now let's plot the results using `matplotlib` roughly as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379): +
+Python code (click to expand) + +```Python +def plot(psi, zlim, norm=None): + xi, yi = np.indices(psi.shape) + fig, ax = pyplot.subplots(subplot_kw={"projection": "3d"}) + pyplot.gca().plot_wireframe( + xi+.5, yi+.5, + psi, color='red', linewidth=.5 + ) + ax.set_zlim(zlim) + for axis in (ax.xaxis, ax.yaxis, ax.zaxis): + axis.pane.fill = False + axis.pane.set_edgecolor('black') + axis.pane.set_alpha(1) + ax.grid(False) + ax.set_zticks([]) + ax.set_xlabel('x/dx') + ax.set_ylabel('y/dy') + ax.set_proj_type('ortho') + cnt = ax.contourf(xi+.5, yi+.5, psi, zdir='z', offset=-1, norm=norm) + cbar = pyplot.colorbar(cnt, pad=.1, aspect=10, fraction=.04) + return cbar.norm + +zlim = (-1, 1) +norm = plot(state_0, zlim) +pyplot.savefig('readme_gauss_0.png') +plot(state, zlim, norm) +pyplot.savefig('readme_gauss.png') +``` +
+ +![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss_0.png) +![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss.png) + +# Debugging + +PyMPDATA relies heavily on Numba to provide high-performance +number crunching operations. Arguably, one of the key advantage +of embracing Numba is that it can be easily switched off. This +brings multiple-order-of-magnitude drop in performance, yet +it also make the entire code of the library susceptible to +interactive debugging, one way of enabling it is by setting the +following environment variable before importing PyMPDATA: +
+Julia code (click to expand) + +```Julia +ENV["NUMBA_DISABLE_JIT"] = "1" +``` +
+
+Matlab code (click to expand) + +```Matlab +setenv('NUMBA_DISABLE_JIT', '1'); +``` +
+
+Python code (click to expand) + +```Python +import os +os.environ["NUMBA_DISABLE_JIT"] = "1" +``` +
+ +# Contributing, reporting issues, seeking support + +See [https://github.com/open-atmos/PyMPDATA/tree/main/README.md](README.md). + +# Other open-source MPDATA implementations: +- mpdat_2d in babyEULAG (FORTRAN) + https://github.com/igfuw/bE_SDs/blob/master/babyEULAG.SDs.for#L741 +- mpdata-oop (C++, Fortran, Python) + https://github.com/igfuw/mpdata-oop +- apc-llc/mpdata (C++) + https://github.com/apc-llc/mpdata +- libmpdata++ (C++): + https://github.com/igfuw/libmpdataxx +- AtmosFOAM: + https://github.com/AtmosFOAM/AtmosFOAM/tree/947b192f69d973ea4a7cfab077eb5c6c6fa8b0cf/applications/solvers/advection/MPDATAadvectionFoam + +# Other Python packages for solving hyperbolic transport equations + +- PyPDE: https://pypi.org/project/PyPDE/ +- FiPy: https://pypi.org/project/FiPy/ +- ader: https://pypi.org/project/ader/ +- centpy: https://pypi.org/project/centpy/ +- mattflow: https://pypi.org/project/mattflow/ +- FastFD: https://pypi.org/project/FastFD/ +- Pyclaw: https://www.clawpack.org/pyclaw/ diff --git a/docs/markdown/python_grid.md b/docs/markdown/python_grid.md deleted file mode 100644 index 1a3c08b0..00000000 --- a/docs/markdown/python_grid.md +++ /dev/null @@ -1,38 +0,0 @@ -
-Python code (click to expand) - -```Python -import numpy as np -from matplotlib import pyplot - -dx, dy = .2, .3 -grid = (10, 5) - -pyplot.scatter(*np.mgrid[ - dx / 2 : grid[0] * dx : dx, - dy / 2 : grid[1] * dy : dy - ], color='red', - label='scalar-field values at cell centres' -) -pyplot.quiver(*np.mgrid[ - 0 : (grid[0]+1) * dx : dx, - dy / 2 : grid[1] * dy : dy - ], 1, 0, pivot='mid', color='green', width=.005, - label='vector-field x-component values at cell walls' -) -pyplot.quiver(*np.mgrid[ - dx / 2 : grid[0] * dx : dx, - 0: (grid[1] + 1) * dy : dy - ], 0, 1, pivot='mid', color='blue', width=.005, - label='vector-field y-component values at cell walls' -) -pyplot.xticks(np.linspace(0, grid[0]*dx, grid[0]+1)) -pyplot.yticks(np.linspace(0, grid[1]*dy, grid[1]+1)) -pyplot.title(f'staggered grid layout (grid={grid}, dx={dx}, dy={dy})') -pyplot.xlabel('x') -pyplot.ylabel('y') -pyplot.legend(bbox_to_anchor=(.1, -.1), loc='upper left', ncol=1) -pyplot.grid() -pyplot.savefig('readme_grid.png') -``` -
diff --git a/docs/markdown/python_grid_2.md b/docs/markdown/python_grid_2.md deleted file mode 100644 index 87aed5c4..00000000 --- a/docs/markdown/python_grid_2.md +++ /dev/null @@ -1,30 +0,0 @@ - -
-Python code (click to expand) - -```Python -from PyMPDATA import ScalarField -from PyMPDATA import VectorField -from PyMPDATA.boundary_conditions import Periodic -import numpy as np - -nx, ny = 24, 24 -Cx, Cy = -.5, -.25 -halo = options.n_halo - -xi, yi = np.indices((nx, ny), dtype=float) -advectee = ScalarField( - data=np.exp( - -(xi+.5-nx/2)**2 / (2*(nx/10)**2) - -(yi+.5-ny/2)**2 / (2*(ny/10)**2) - ), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -advector = VectorField( - data=(np.full((nx + 1, ny), Cx), np.full((nx, ny + 1), Cy)), - halo=halo, - boundary_conditions=(Periodic(), Periodic()) -) -``` -
diff --git a/docs/markdown/python_solver.md b/docs/markdown/python_solver.md deleted file mode 100644 index cfe49941..00000000 --- a/docs/markdown/python_solver.md +++ /dev/null @@ -1,35 +0,0 @@ -
-Python code (click to expand) - -```Python -def plot(psi, zlim, norm=None): - xi, yi = np.indices(psi.shape) - fig, ax = pyplot.subplots(subplot_kw={"projection": "3d"}) - pyplot.gca().plot_wireframe( - xi+.5, yi+.5, - psi, color='red', linewidth=.5 - ) - ax.set_zlim(zlim) - for axis in (ax.xaxis, ax.yaxis, ax.zaxis): - axis.pane.fill = False - axis.pane.set_edgecolor('black') - axis.pane.set_alpha(1) - ax.grid(False) - ax.set_zticks([]) - ax.set_xlabel('x/dx') - ax.set_ylabel('y/dy') - ax.set_proj_type('ortho') - cnt = ax.contourf(xi+.5, yi+.5, psi, zdir='z', offset=-1, norm=norm) - cbar = pyplot.colorbar(cnt, pad=.1, aspect=10, fraction=.04) - return cbar.norm - -zlim = (-1, 1) -norm = plot(state_0, zlim) -pyplot.savefig('readme_gauss_0.png') -plot(state, zlim, norm) -pyplot.savefig('readme_gauss.png') -``` -
- -![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss_0.png) -![plot](https://github.com/open-atmos/PyMPDATA/releases/download/tip/readme_gauss.png) diff --git a/docs/markdown/stepper_solver.md b/docs/markdown/stepper_solver.md deleted file mode 100644 index fc9a9551..00000000 --- a/docs/markdown/stepper_solver.md +++ /dev/null @@ -1,143 +0,0 @@ -### Stepper - -The logic of the MPDATA iterative solver is represented -in PyMPDATA by the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class. - -When instantiating the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html), the user has a choice -of either supplying just the number of dimensions or specialising the stepper for a given grid: -
-Julia code (click to expand) - -```Julia -Stepper = pyimport("PyMPDATA").Stepper - -stepper = Stepper(options=options, n_dims=2) -``` -
-
-Matlab code (click to expand) - -```Matlab -Stepper = py.importlib.import_module('PyMPDATA').Stepper; - -stepper = Stepper(pyargs(... - 'options', options, ... - 'n_dims', int32(2) ... -)); -``` -
-
-Python code (click to expand) - -```Python -from PyMPDATA import Stepper - -stepper = Stepper(options=options, n_dims=2) -``` -
-or -
-Julia code (click to expand) - -```Julia -stepper = Stepper(options=options, grid=(nx, ny)) -``` -
-
-Matlab code (click to expand) - -```Matlab -stepper = Stepper(pyargs(... - 'options', options, ... - 'grid', py.tuple({nx, ny}) ... -)); -``` -
-
-Python code (click to expand) - -```Python -stepper = Stepper(options=options, grid=(nx, ny)) -``` -
- -In the latter case, noticeably -faster execution can be expected, however the resultant -stepper is less versatile as bound to the given grid size. -If number of dimensions is supplied only, the integration -might take longer, yet same instance of the -stepper can be used for different grids. - -Since creating an instance of the [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) class -involves time-consuming compilation of the algorithm code, -the class is equipped with a cache logic - subsequent -calls with same arguments return references to previously -instantiated objects. Instances of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) contain no -mutable data and are (thread-)safe to be reused. - -The init method of [``Stepper``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/stepper.html) has an optional -``non_unit_g_factor`` argument which is a Boolean flag -enabling handling of the G factor term which can be used to -represent coordinate transformations and/or variable fluid density. - -Optionally, the number of threads to use for domain decomposition -in the first (non-contiguous) dimension during 2D and 3D calculations -may be specified using the optional ``n_threads`` argument with a -default value of ``numba.get_num_threads()``. The multi-threaded -logic of PyMPDATA depends thus on settings of numba, namely on the -selected threading layer (either via ``NUMBA_THREADING_LAYER`` env -var or via ``numba.config.THREADING_LAYER``) and the selected size of the -thread pool (``NUMBA_NUM_THREADS`` env var or ``numba.config.NUMBA_NUM_THREADS``). - - -#### Solver - -Instances of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class are used to control -the integration and access solution data. During instantiation, -additional memory required by the solver is -allocated according to the options provided. - -The only method of the [``Solver``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html) class besides the -init is [``advance(n_steps, mu_coeff, ...)``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/solver.html#Solver.advance) -which advances the solution by ``n_steps`` timesteps, optionally -taking into account a given diffusion coefficient ``mu_coeff``. - -Solution state is accessible through the ``Solver.advectee`` property. -Multiple solver[s] can share a single stepper, e.g., as exemplified in the shallow-water system solution in the examples package. - -Continuing with the above code snippets, instantiating -a solver and making 75 integration steps looks as follows: -
-Julia code (click to expand) - -```Julia -Solver = pyimport("PyMPDATA").Solver -solver = Solver(stepper=stepper, advectee=advectee, advector=advector) -solver.advance(n_steps=75) -state = solver.advectee.get() -``` -
-
-Matlab code (click to expand) - -```Matlab -Solver = py.importlib.import_module('PyMPDATA').Solver; -solver = Solver(pyargs('stepper', stepper, 'advectee', advectee, 'advector', advector)); -solver.advance(pyargs('n_steps', 75)); -state = solver.advectee.get(); -``` -
-
-Python code (click to expand) - -```Python -from PyMPDATA import Solver - -solver = Solver(stepper=stepper, advectee=advectee, advector=advector) -state_0 = solver.advectee.get().copy() -solver.advance(n_steps=75) -state = solver.advectee.get() -``` -
- -Now let's plot the results using `matplotlib` roughly as in Fig. 5 in [Arabas et al. 2014](https://doi.org/10.3233/SPR-140379): From d05246abdff6abe5197d0a14669b762b8b0e6de2 Mon Sep 17 00:00:00 2001 From: pawel Date: Sat, 28 Sep 2024 00:34:37 +0200 Subject: [PATCH 08/10] add simple landing page for pympdata-examples --- docs/markdown/pympdata_examples_landing.md | 24 ++++++++++++++++++++++ examples/PyMPDATA_examples/__init__.py | 5 +---- 2 files changed, 25 insertions(+), 4 deletions(-) create mode 100644 docs/markdown/pympdata_examples_landing.md diff --git a/docs/markdown/pympdata_examples_landing.md b/docs/markdown/pympdata_examples_landing.md new file mode 100644 index 00000000..473620d4 --- /dev/null +++ b/docs/markdown/pympdata_examples_landing.md @@ -0,0 +1,24 @@ +# Introduction +PyMPDATA examples are bundled with PyMPDATA and located in the examples subfolder. +They constitute a separate PyMPDATA_examples Python package which is also available at PyPI. +The examples have additional dependencies listed in PyMPDATA_examples package setup.py file. +Running the examples requires the PyMPDATA_examples package to be installed. + +Below is an example of how to use the PyMPDATA_examples package to run a simple advection-diffusion in 2D +`PyMPDATA_examples.advection_diffusion_2d` +![adv_diff](https://github.com/open-atmos/PyMPDATA/releases/download/tip/advection_diffusion.gif) + +# Installation +Since the examples package includes Jupyter notebooks (and their execution requires write access), the suggested install and launch steps are: + +``` +git clone https://github.com/open-atmos/PyMPDATA-examples.git +cd PyMPDATA-examples +pip install -e . +jupyter-notebook +``` + +Alternatively, one can also install the examples package from pypi.org by using +``` +pip install PyMPDATA-examples. +``` diff --git a/examples/PyMPDATA_examples/__init__.py b/examples/PyMPDATA_examples/__init__.py index 4e88a725..7cc7d8cd 100644 --- a/examples/PyMPDATA_examples/__init__.py +++ b/examples/PyMPDATA_examples/__init__.py @@ -1,8 +1,5 @@ """ -PyMPDATA_examples package includes common Python modules used in PyMPDATA smoke tests -and in example notebooks (but the package wheels do not include the notebooks). -![adv_diff](https://github.com/open-atmos/PyMPDATA/releases/download/tip/advection_diffusion.gif) - +.. include:: ../../docs/markdown/pympdata_examples_landing.md """ from importlib.metadata import PackageNotFoundError, version From 62cebd35fb9df2860c840f0e541deb20856700f2 Mon Sep 17 00:00:00 2001 From: pawel Date: Sat, 28 Sep 2024 01:22:56 +0200 Subject: [PATCH 09/10] updates to css on main landing page --- docs/templates/index.html.jinja2 | 113 +++++++++++++++++++------------ 1 file changed, 68 insertions(+), 45 deletions(-) diff --git a/docs/templates/index.html.jinja2 b/docs/templates/index.html.jinja2 index e6e03856..3199adce 100644 --- a/docs/templates/index.html.jinja2 +++ b/docs/templates/index.html.jinja2 @@ -8,12 +8,13 @@ {% block content %}
-
- - pympdata logo - -
-
+
+ + pympdata logo + +

Documentation

+
+

What is PyMPDATA?

PyMPDATA is a Numba-accelerated Pythonic implementation of the MPDATA algorithm @@ -21,34 +22,11 @@ numerically solving generalised convection-diffusion PDEs in 1D, 2D and 3D structured meshes with coordinate transformations.

+

+ A separate project called PyMPDATA-MPI depicts how numba-mpi can be used to enable distributed memory parallelism in PyMPDATA. +

- -

PyMPDATA resources

- - - - - - - - - - - -
PyMPDATAPyMPDATA-examplesPyMPDATA-MPI
- - - - -
- -
+

What is the difference between PyMPDATA and PyMPDATA-examples?

PyMPDATA is a Python package that provides the MPDATA algorithm implementation. @@ -60,17 +38,52 @@ but the package wheels do not include the notebooks.

They exist separately on PyPI, but their development and issue tracking is in the same GitHub repository.

-

Dependencies

-

- PyMPDATA depends only on NumPy and Numba. -

-

- PyMPDATA-examples requires additional packages listed in install_requires in - setup.py. - Amongst them is PyMPDATA. -

-
+
+

Important links

+ + + + + + + + + + + + + + + +
PyMPDATAPyMPDATA-examplesPyMPDATA-MPI
+ + + + + +
+ + + +
+
+

Installation

PyMPDATA is available on PyPI and can be installed using pip: @@ -85,7 +98,18 @@ 

pip install PyMPDATA-examples

Note: this will also install PyMPDATA

-
+
+

Dependencies

+

+ PyMPDATA depends only on NumPy and Numba. +

+

+ PyMPDATA-examples requires additional packages listed in install_requires in + setup.py. + Amongst them is PyMPDATA. +

+
+

Contributing, reporting issues and seeking support

Submitting new code to both packages is done through the same GitHub repository via @@ -102,7 +126,6 @@ using and extending PyMPDATA code.

-
{% include "search.html.jinja2" %} {% endblock %} From be86b3b07c9e44de0c247668f944e06e5d87f8fb Mon Sep 17 00:00:00 2001 From: Sylwester Arabas Date: Tue, 1 Oct 2024 01:38:02 +0200 Subject: [PATCH 10/10] small changes to the text, links, etc --- README.md | 2 +- docs/markdown/pympdata_landing.md | 7 ++++--- docs/templates/index.html.jinja2 | 12 ++++++------ 3 files changed, 11 insertions(+), 10 deletions(-) diff --git a/README.md b/README.md index 3fe3f39e..2d9a3784 100644 --- a/README.md +++ b/README.md @@ -77,7 +77,7 @@ The key classes constituting the PyMPDATA interface are summarised below. The [``Options``](https://open-atmos.github.io/PyMPDATA/PyMPDATA/options.html) class groups both algorithm variant options as well as some implementation-related -flags that need to be set at the first place. +flags. #### Arakawa-C grid layer diff --git a/docs/markdown/pympdata_landing.md b/docs/markdown/pympdata_landing.md index d48ee473..c3313be0 100644 --- a/docs/markdown/pympdata_landing.md +++ b/docs/markdown/pympdata_landing.md @@ -1,13 +1,14 @@ pympdata logo # Introduction PyMPDATA is a high-performance Numba-accelerated Pythonic implementation of the MPDATA - algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond. + algorithm of [Smolarkiewicz et al.](https://doi.org/10.1016%2F0021-9991%2884%2990121-9) + used in geophysical fluid dynamics and beyond. MPDATA numerically solves generalised transport equations - partial differential equations used to model conservation/balance laws, scalar-transport problems, convection-diffusion phenomena. As of the current version, PyMPDATA supports homogeneous transport in 1D, 2D and 3D using structured meshes, optionally - generalised by employment of a Jacobian of coordinate transformation. + generalised by employment of a Jacobian of coordinate transformation or a fluid density profile. PyMPDATA includes implementation of a set of MPDATA variants including the non-oscillatory option, infinite-gauge, divergent-flow, double-pass donor cell (DPDC) and third-order-terms options. @@ -29,7 +30,7 @@ PyMPDATA design features to concisely represent multi-dimensional stencil operations on both scalar and vector fields. The grid layer is built on top of NumPy's ndarrays (using "C" ordering) - using the Numba's ``@njit`` functionality for high-performance array traversals. + using the Numba's ``@jit(nopython=True)`` functionality for high-performance array traversals. It enables one to code once for multiple dimensions, and automatically handles (and hides from the user) any halo-filling logic related with boundary conditions. Numba ``prange()`` functionality is used for implementing multi-threading diff --git a/docs/templates/index.html.jinja2 b/docs/templates/index.html.jinja2 index 3199adce..722e99b3 100644 --- a/docs/templates/index.html.jinja2 +++ b/docs/templates/index.html.jinja2 @@ -19,8 +19,8 @@

PyMPDATA is a Numba-accelerated Pythonic implementation of the MPDATA algorithm of Smolarkiewicz et al. used in geophysical fluid dynamics and beyond for - numerically solving generalised convection-diffusion PDEs in 1D, 2D and 3D structured meshes - with coordinate transformations. + numerically solving generalised convection-diffusion PDEs. PyMPDATA supports integration in 1D, 2D and 3D structured meshes + with optional coordinate transformations.

A separate project called PyMPDATA-MPI depicts how numba-mpi can be used to enable distributed memory parallelism in PyMPDATA. @@ -34,10 +34,10 @@

PyMPDATA-examples is a Python package that provides examples of how to use PyMPDATA. - It includes common Python modules used in PyMPDATA smoke tests and in example notebooks, - but the package wheels do not include the notebooks. + It includes common Python modules used in PyMPDATA smoke tests and in example Jupyter notebooks + (but the package wheels do not include the notebooks, only .py files imported from the notebooks and PyMPDATA tests).

-

They exist separately on PyPI, but their development and issue tracking is in the same GitHub repository.

+

The two projects exist separately on PyPI, but their development and issue tracking is hosted at the same GitHub repository.

Important links

@@ -101,7 +101,7 @@

Dependencies

- PyMPDATA depends only on NumPy and Numba. + PyMPDATA depends on NumPy, Numba and pystrict.

PyMPDATA-examples requires additional packages listed in install_requires in