diff --git a/.github/actions/install_h5py/action.yml b/.github/actions/install_h5py/action.yml index caa37a512..e0a8fb30c 100644 --- a/.github/actions/install_h5py/action.yml +++ b/.github/actions/install_h5py/action.yml @@ -21,6 +21,7 @@ runs: run: | export CC="mpicc" export HDF5_MPI="ON" + pip install 'mpi4py>=4.0.0' pip install h5py --no-cache-dir --no-binary h5py pip list diff --git a/.github/workflows/documentation.yml b/.github/workflows/documentation.yml index 5eb1b20a7..60885bc19 100644 --- a/.github/workflows/documentation.yml +++ b/.github/workflows/documentation.yml @@ -23,9 +23,11 @@ permissions: jobs: build_docs: runs-on: ubuntu-latest + env: GITHUB_PAT: ${{ secrets.GITHUB_TOKEN}} OMP_NUM_THREADS: 2 + steps: - name: Checkout uses: actions/checkout@v6 @@ -38,17 +40,7 @@ jobs: python-version: '3.10' - name: Install non-Python dependencies on Ubuntu - uses: awalsh128/cache-apt-pkgs-action@latest - with: - packages: gfortran openmpi-bin libopenmpi-dev libhdf5-openmpi-dev - version: 1.0 - execute_install_scripts: true - - - name: Reconfigure non-Python dependencies on Ubuntu - run: | - sudo apt-get update - sudo apt-get install --reinstall openmpi-bin libhdf5-openmpi-dev liblapack-dev libblas-dev - sudo apt install graphviz pandoc + uses: ./.github/actions/ubuntu_installations - name: Print information on MPI and HDF5 libraries run: | @@ -70,6 +62,10 @@ jobs: pip install .[test] pip freeze + - name: Install non-Python dependencies for Documentation + run: | + sudo apt install graphviz pandoc + - name: Install Python dependencies for Documentation run: | pip install -r docs/requirements.txt diff --git a/.github/workflows/testing.yml b/.github/workflows/testing.yml index 1b4c2ff7b..9a05f63a5 100644 --- a/.github/workflows/testing.yml +++ b/.github/workflows/testing.yml @@ -40,7 +40,8 @@ jobs: OMP_NUM_THREADS: 2 steps: - - uses: actions/checkout@v6 + - name: Checkout + uses: actions/checkout@v6 with: submodules: true @@ -125,17 +126,17 @@ jobs: if: matrix.os == 'macos-14' working-directory: ./scratch run: | - coverage report --ignore-errors --show-missing --sort=cover + psydac test --mod psydac.mapping -x -v - name: Run MPI tests with Pytest working-directory: ./scratch run: | - psydac test --mpi + psydac test --mod psydac.fem -x -v - name: Run single-process PETSc tests with Pytest working-directory: ./scratch run: | - psydac test --petsc + psydac test --mod psydac.feec -x -v - name: Run MPI PETSc tests with Pytest working-directory: ./scratch diff --git a/CHANGELOG.md b/CHANGELOG.md index d8fc460e8..f9dc0e1c3 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,19 +11,23 @@ All notable changes to this project will be documented in this file. - #576 : Add module `psydac.utilities.parallel_utils` for parallel execution and gathering variable-length arrays - #576 : Add module `psydac.utilities.operators` with class `Laplacian` used in some old examples - #577 : Add an installation configuration option to choose the backend language +- #567 : Improve `psydac test` command (with several new features) +- #565 : Expand editable install info in `README.md` - [DEVELOPER] Create action `install_petsc4py` to install PETSc & `petsc4py` w/ complex support +- [DEVELOPER] Configure Pylint in `pyproject.toml` +- [DEVELOPER] Add `CLAUDE.md` with developer rules ### Fixed -- #576 : Fix bug in `TensorFemSpace.eval_field` caused by round-off at MPI subdomain boundaries +- #527 : Require `sympde==0.20.0` which changes how multipatch interfaces are defined - #576 : Require `sympde==0.19.3` which fixes a bug in the linearity checks - #579 : Require `h5py>=3.16` which installs correctly with `setuptools>=81.0` +- #576 : Fix bug in `TensorFemSpace.eval_field` caused by round-off at MPI subdomain boundaries - #579 : Don't run postprocessing unit tests with `pytest-xdist` because `h5py` is not thread-safe - #579 : Return error code on failure of the `psydac test` and `psydac compile` commands - #577 : Fix installation following release of Pyccel 2.2 - #571 : Fix correct application of the sum factorization algorithm -- #570 : Optimize PSYDAC logo -- #565 : Expand editable install info in `README.md` +- #567 : Fix parallel creation of folder `__psydac__` in `psydac.api.fem_bilinear_form` - #566 : Fix command `psydac test --mpi` on Ubuntu machines - [DEVELOPER] Add missing 'description' properties (required!) to our GitHub actions - [DEVELOPER] Update CI installation of `petsc4py` after release of `setuptools` 81.0 @@ -32,6 +36,7 @@ All notable changes to this project will be documented in this file. ### Changed +- #527 : Improve `Geometry` class in module `psydac.cad.geometry` - #576 : Use latest version of Igakit (commit dalcinl/igakit@92ee097 of 2026/07/24) which supports NumPy >= 2.4 - #595 : Use PETSc 3.25.5 whose Python bindings `petsc4py` are built correctly with `cython>=3` - #580 : Use PETSc 3.25.0 whose Python bindings `petsc4py` install correctly with `setuptools>=81.0` @@ -39,6 +44,7 @@ All notable changes to this project will be documented in this file. - #579 : Require `numpy>=2.1` to support Python >= 3.10 - #579 : Require `pytest>=9.0` and use `pytest.toml` instead of `pytest.ini` for Pytest configuration - #579 : Move coverage configuration from `pyproject.toml` to `psydac/pytest.toml` +- #570 : Optimize PSYDAC logo - [DEVELOPER] Update GitHub Actions for repository checkout and Python setup - [DEVELOPER] Rename actions: `macos/ubuntu_install` -> `macos/ubuntu_installations` - [DEVELOPER] Do not check file changes to trigger testing workflow on PRs diff --git a/CLAUDE.md b/CLAUDE.md new file mode 100644 index 000000000..d25c1bb03 --- /dev/null +++ b/CLAUDE.md @@ -0,0 +1,74 @@ +# Project: PSYDAC + +## Overview +Python 3 library for isogeometric analysis (IGA). Solve general systems of partial +differential equations (PDEs) in weak form, defined using the domain-specific +language provided by SymPDE. Supports finite element exterior calculus (FEEC) +with tensor-product spline spaces. Handles multi-patch geometries in various ways; +usually broken-FEEC a.k.a. CONGA (conforming/non-conforming Galerkin). +Python code is automatically generated for the assembly of user-defined functionals, +linear forms, and bilinear forms. This Python code is then accelerated to C/Fortran +speed using Pyccel. The library enables large parallel computations on distributed- +memory supercomputers using MPI and OpenMP. + +## Tech stack +- Python 3.10+, type hints required everywhere +- meson-python + mesonpy for building +- pip with submodules (for igakit) +- pytest + pytest-cov + pytest-mpi + pytest-xdist for testing +- sympde for symbolic definition of weak formulations +- pyccel for translation of Python kernels to Fortran or C +- mpi4py for MPI parallelization +- h5py for parallel I/O +- petsc4py for direct linear solvers (optional install) + +## Code standards +- Follow PEP 8 style guide whenever possible +- Docstrings of public functions and classes follow Numpydoc conventions +- pylint + black + isort for linting/typing +- Computational kernels to be accelerated with Pyccel named as `_kernels.py` +- Minimize code duplication. Never add code that already exists +- Strive for concise, clean, and human-readable code. Avoid useless verbosity +- Self-explanatory names for variables, functions, and classes +- Do not reassign value to an existing variable, especially if the type changes +- Use short comments (one-liners or inline) to explain "why" rather than "what" or "how" + +## Testing conventions +- Each library subpackage contains a `tests/` folder with an `__init__.py` file +- Unit tests are part of the library and shipped with it +- Unit tests can be run with `psydac test` CLI command (see `README.md`) +- Test file names follow `test_.py` +- Test function names follow `test_` or `test__` +- Parametrize unit tests with `@pytest.parametrize` to minimize code duplication +- Keep run time of tests at a minimum +- Aim for 100% coverage on newly committed code + +## File structure +- psydac/api - high-level Python interface +- psydac/cad - computer-aided design (CAD) functionality +- psydac/cmd - command-line interface (CLI) commands +- psydac/core - splines functionality +- psydac/ddm - MPI decomposition of domain, vectors of spline coefficients, and matrices +- psydac/feec - finite element exterior calculus (FEEC) +- psydac/fem - "finite element method" middle-level Python interface (spaces and fields) +- psydac/linalg - "linear algebra" low-level Python interface (vectors, linear operators, iterative solvers) +- psydac/polar - implementation of polar splines for H1 spaces +- psydac/pyccel - distillation of an old version of Pyccel for Python code generation (to be removed) +- psydac/utilities - various generic functions used across the library +- docs - Sphinx documentation +- examples - Various examples of library usage (mostly Jupyter notebooks) +- scripts - Various scripts for library maintenance +- subprojects - Git submodules + +## Always do +- Before claiming done: + . Run `python -m black --check $(git ls-files "*.py")` + . Run `python -m isort --check $(git ls-files "*.py" | grep -v "FOUND_DUPLICATED_IMPORT.py")` + . Run `pylint --disable=all --enable=unused-import ./psydac --ignore-paths="__pyccel*" --ignore-paths="__epyccel*"` +- Unless specified differently, assume process with `rank=0` in MPI communicator is root +- Use root MPI process to perform serial operations (e.g. open/save image) + +## Never do +- Commit secrets or `.env` files +- Add new dependencies without team discussion +- Use `eval` function together with `input` diff --git a/psydac/api/ast/tests/test_poisson.py b/psydac/api/ast/tests/test_poisson.py index d6470a32f..c6025cd32 100644 --- a/psydac/api/ast/tests/test_poisson.py +++ b/psydac/api/ast/tests/test_poisson.py @@ -12,8 +12,7 @@ from sympde.calculus import grad, dot from sympde.topology import ScalarFunctionSpace from sympde.topology import elements_of, LogicalExpr -from sympde.topology import Square -from sympde.topology import Mapping, IdentityMapping, PolarMapping +from sympde.topology import Domain from sympde.expr import integral from sympde.expr import LinearForm from sympde.expr import BilinearForm @@ -38,9 +37,8 @@ #============================================================================== def test_codegen(): - domain = Square() - M = Mapping('M',2) - domain = M(domain) + domain = Domain.from_file(filename) + M = domain.mapping V = ScalarFunctionSpace('V', domain) u,v = elements_of(V, names='u,v') diff --git a/psydac/api/ast/tests/test_system_1.py b/psydac/api/ast/tests/test_system_1.py index 95fc81744..3694b5503 100644 --- a/psydac/api/ast/tests/test_system_1.py +++ b/psydac/api/ast/tests/test_system_1.py @@ -12,7 +12,7 @@ from sympde.topology import VectorFunctionSpace from sympde.topology import element_of from sympde.topology import Square -from sympde.topology import Mapping#, IdentityMapping, PolarMapping +from sympde.topology import Domain, Mapping from sympde.expr import integral from sympde.expr import LinearForm from sympde.expr import BilinearForm @@ -61,7 +61,9 @@ def test_codegen(): h1norm_F = SemiNorm(error, domain, kind='h1') # Create computational domain from topological domain - domain_h = discretize(domain, filename=filename) + # The spaces are defined on the logical domain, but the geometry needs + # the mapped domain defined in the file + domain_h = discretize(Domain.from_file(filename), filename=filename) # Discrete spaces Vh = discretize(V, domain_h) diff --git a/psydac/api/ast/tests/test_system_2.py b/psydac/api/ast/tests/test_system_2.py index 72ea2dbcd..3aecbddeb 100644 --- a/psydac/api/ast/tests/test_system_2.py +++ b/psydac/api/ast/tests/test_system_2.py @@ -12,7 +12,7 @@ from sympde.topology import ScalarFunctionSpace, VectorFunctionSpace from sympde.topology import element_of from sympde.topology import Square -from sympde.topology import Mapping +from sympde.topology import Domain, Mapping from sympde.expr import integral from sympde.expr import LinearForm from sympde.expr import BilinearForm @@ -59,7 +59,9 @@ def test_codegen(): l = LinearForm((q,v), int_0(f1*q[0]+f2*q[1]+v)) # Create computational domain from topological domain - domain_h = discretize(domain, filename=filename) + # The spaces are defined on the logical domain, but the geometry needs + # the mapped domain defined in the file + domain_h = discretize(Domain.from_file(filename), filename=filename) # Discrete spaces Vh = discretize(V1*V2, domain_h) diff --git a/psydac/api/ast/tests/test_system_3.py b/psydac/api/ast/tests/test_system_3.py index d5bc276a4..576df45b4 100644 --- a/psydac/api/ast/tests/test_system_3.py +++ b/psydac/api/ast/tests/test_system_3.py @@ -12,7 +12,7 @@ from sympde.topology import VectorFunctionSpace from sympde.topology import element_of from sympde.topology import Square -from sympde.topology import Mapping#, IdentityMapping, PolarMapping +from sympde.topology import Domain, Mapping from sympde.expr import integral from sympde.expr import BilinearForm from sympde.expr.evaluation import TerminalExpr @@ -54,7 +54,9 @@ def test_codegen(): b = BilinearForm((u,v), int_0(div(v)*div(u))) # Create computational domain from topological domain - domain_h = discretize(domain, filename=filename) + # The spaces are defined on the logical domain, but the geometry needs + # the mapped domain defined in the file + domain_h = discretize(Domain.from_file(filename), filename=filename) # Discrete spaces Vh = discretize(V, domain_h) diff --git a/psydac/api/discretization.py b/psydac/api/discretization.py index a7a1a5e31..00e645461 100644 --- a/psydac/api/discretization.py +++ b/psydac/api/discretization.py @@ -578,7 +578,7 @@ def discretize_domain(domain, *, filename=None, ncells=None, periodic=None, comm raise ValueError("Cannot provide both 'filename' and 'ncells'") elif filename: - return Geometry(filename=filename, comm=comm, mpi_dims_mask=mpi_dims_mask) + return Geometry.from_file(filename, domain=domain, comm=comm, mpi_dims_mask=mpi_dims_mask) elif ncells: return Geometry.from_topological_domain(domain, ncells, periodic=periodic, comm=comm, mpi_dims_mask=mpi_dims_mask) diff --git a/psydac/api/tests/test_2d_complex.py b/psydac/api/tests/test_2d_complex.py index 85c2479dc..014ef355f 100644 --- a/psydac/api/tests/test_2d_complex.py +++ b/psydac/api/tests/test_2d_complex.py @@ -363,7 +363,7 @@ def test_complex_poisson_2d_multipatch(): A = Square('A',bounds1=(0, 0.5), bounds2=(0, 1)) B = Square('B',bounds1=(0.5, 1.), bounds2=(0, 1)) - domain = Domain.join([A, B], [((0, 0, 1), (1, 0, -1))], 'domain') + domain = Domain.join([A, B], [((0, 0, 1), (1, 0, -1), 1)], 'domain') x, y = domain.coordinates diff --git a/psydac/api/tests/test_2d_multipatch_mapping_maxwell.py b/psydac/api/tests/test_2d_multipatch_mapping_maxwell.py index 336066c71..2ae7ea7b5 100644 --- a/psydac/api/tests/test_2d_multipatch_mapping_maxwell.py +++ b/psydac/api/tests/test_2d_multipatch_mapping_maxwell.py @@ -7,25 +7,22 @@ from pathlib import Path import pytest -import numpy as np from mpi4py import MPI from sympy import pi, sin, cos, Tuple, Matrix -from sympde.calculus import grad, dot, curl, cross +from sympde.calculus import dot, curl, cross from sympde.calculus import minus, plus from sympde.topology import VectorFunctionSpace from sympde.topology import elements_of from sympde.topology import NormalVector -from sympde.topology import Square, Domain -from sympde.topology import IdentityMapping, PolarMapping +from sympde.topology import Domain from sympde.expr.expr import LinearForm, BilinearForm from sympde.expr.expr import integral from sympde.expr.expr import Norm -from sympde.expr.equation import find, EssentialBC +from sympde.expr.equation import find from psydac.api.discretization import discretize from psydac.api.tests.build_domain import build_11_patch_pretzel, build_2_patch_annulus -from psydac.fem.basic import FemField from psydac.api.settings import PSYDAC_BACKEND_GPYCCEL from psydac.feec.pull_push import pull_2d_hcurl @@ -204,10 +201,10 @@ def teardown_function(): if __name__ == '__main__': - from collections import OrderedDict - from sympy import lambdify + from collections import OrderedDict + from sympy import lambdify from psydac.fem.plotting_utilities import get_plotting_grid, get_grid_vals - from psydac.fem.plotting_utilities import get_patch_knots_gridlines, my_small_plot + from psydac.fem.plotting_utilities import my_small_plot domain = build_11_patch_pretzel() x,y = domain.coordinates diff --git a/psydac/api/tests/test_2d_multipatch_mapping_poisson.py b/psydac/api/tests/test_2d_multipatch_mapping_poisson.py index f730f6bc6..3f7a0de19 100644 --- a/psydac/api/tests/test_2d_multipatch_mapping_poisson.py +++ b/psydac/api/tests/test_2d_multipatch_mapping_poisson.py @@ -15,7 +15,7 @@ from sympde.calculus import minus, plus from sympde.topology import ScalarFunctionSpace from sympde.topology import elements_of -from sympde.topology import NormalVector, Union +from sympde.topology import NormalVector from sympde.topology import Square, Domain from sympde.topology import IdentityMapping, PolarMapping, AffineMapping from sympde.expr.expr import LinearForm, BilinearForm @@ -131,7 +131,7 @@ def test_poisson_2d_2_patches_dirichlet_1(): assert ( abs(h1_error - expected_h1_error) < 1e-7 ) #------------------------------------------------------------------------------ -def test_poisson_2d_3_patches_dirichlet_2(): +def test_poisson_2d_3_patches_dirichlet_2(tmp_path): mapping_1 = IdentityMapping('M1', 2) mapping_2 = PolarMapping ('M2', 2, c1 = 0., c2 = 0.5, rmin = 0., rmax=1.) @@ -141,12 +141,12 @@ def test_poisson_2d_3_patches_dirichlet_2(): B = Square('B',bounds1=(0.5, 1.), bounds2=(0, np.pi)) C = Square('C',bounds1=(0.5, 1.), bounds2=(np.pi-0.5, np.pi + 1)) - D1 = mapping_1(A) - D2 = mapping_2(B) - D3 = mapping_3(C) + D1 = mapping_1(A) + D2 = mapping_2(B) + D3 = mapping_3(C) - connectivity = [((0,1,1),(1,1,-1)), ((1,1,1),(2,1,-1))] patches = [D1, D2, D3] + connectivity = [((0, 1, 1), (1, 1,-1), 1), ((1, 1, 1), (2, 1,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates @@ -155,7 +155,7 @@ def test_poisson_2d_3_patches_dirichlet_2(): l2_error, h1_error, uh = run_poisson_2d(solution, f, domain, ncells=[2**2,2**2], degree=[2,2]) - plot_fn=f'uh_multipatch_poisson.pdf' + plot_fn = str(tmp_path / 'uh_multipatch_poisson.pdf') plot_field(fem_field=uh, Vh=uh.space, domain=domain, title='uh', filename=plot_fn, hide_plot=True) expected_l2_error = 0.0019402242901236006 @@ -252,8 +252,11 @@ def test_poisson_2d_4_patch_dirichlet_0(): D3 = mapping_3(C) D4 = mapping_4(D) - connectivity = [((0,1,1),(1,1,-1)), ((2,1,1),(3,1,-1)), ((0,0,1),(2,0,-1)),((1,0,1),(3,0,-1))] patches = [D1, D2, D3, D4] + connectivity = [((0, 1, 1), (1, 1,-1), 1), + ((2, 1, 1), (3, 1,-1), 1), + ((0, 0, 1), (2, 0,-1), 1), + ((1, 0, 1), (3, 0,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates diff --git a/psydac/api/tests/test_2d_multipatch_poisson.py b/psydac/api/tests/test_2d_multipatch_poisson.py index 297975eba..b4f230ebb 100644 --- a/psydac/api/tests/test_2d_multipatch_poisson.py +++ b/psydac/api/tests/test_2d_multipatch_poisson.py @@ -78,8 +78,8 @@ def test_poisson_2d_2_patch_dirichlet_0(): A = Square('A', bounds1=(0, 0.5), bounds2=(0, 1)) B = Square('B', bounds1=(0.5, 1), bounds2=(0, 1)) - connectivity = [((0,0,1), (1,0,-1), 1)] patches = [A, B] + connectivity = [((0, 0, 1), (1, 0,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates @@ -100,8 +100,8 @@ def test_poisson_2d_2_patch_dirichlet_1(): A = Square('A', bounds1=(0, 0.5), bounds2=(0, 1)) B = Square('B', bounds1=(0.5, 1), bounds2=(0, 1)) - connectivity = [((0,0,1), (1,0,-1), 1)] patches = [A, B] + connectivity = [((0, 0, 1), (1, 0,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates @@ -121,8 +121,8 @@ def test_poisson_2d_2_patch_dirichlet_2(): A = Square('A', bounds1=(0, 0.5), bounds2=(0, 1)) B = Square('B', bounds1=(0.5, 1), bounds2=(0, 1)) - connectivity = [((0,0,1), (1,0,-1), 1)] patches = [A, B] + connectivity = [((0, 0, 1), (1, 0,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates @@ -150,8 +150,8 @@ def test_poisson_2d_2_patch_dirichlet_3(): D1 = M1(A) D2 = M2(B) - connectivity = [((0,0,1), (1,0,1), 1)] patches = [D1, D2] + connectivity = [((0, 0, 1), (1, 0, 1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates @@ -182,8 +182,8 @@ def test_poisson_2d_2_patch_dirichlet_4(): D1 = M1(A) D2 = M2(B) - connectivity = [((0,0,-1), (1,0,-1), 1)] - patches = [D1,D2] + patches = [D1, D2] + connectivity = [((0, 0, -1), (1, 0, -1), 1)] domain = Domain.join(patches, connectivity, 'domain') x,y = domain.coordinates diff --git a/psydac/api/tests/test_postprocessing.py b/psydac/api/tests/test_postprocessing.py index 73cc941f9..b592f2070 100644 --- a/psydac/api/tests/test_postprocessing.py +++ b/psydac/api/tests/test_postprocessing.py @@ -57,7 +57,7 @@ def build_2_mapped_squares(): D2 = mapping_2(B) patches = [D1, D2] - connectivity = [((0,1,1),(1,1,-1))] + connectivity = [((0, 1, 1), (1, 1,-1), 1)] return Domain.join(patches, connectivity, 'domain') @@ -66,7 +66,7 @@ def build_2_squares(): B = Square('B',bounds1=(0.5, 1.), bounds2=(np.pi/2, np.pi)) patches = [A, B] - connectivity = [((0,1,1),(1,1,-1))] + connectivity = [((0, 1, 1), (1, 1,-1), 1)] return Domain.join(patches, connectivity, 'domain') @@ -75,10 +75,9 @@ def build_2_cubes(): B = Cube('B',bounds1=(0.5, 1.), bounds2=(np.pi/2, np.pi), bounds3=(0, 1)) patches = [A, B] - connectivity = [((0,1,1),(1,1,-1))] + connectivity = [((0, 1, 1), (1, 1,-1), (1, 1, 1))] return Domain.join(patches, connectivity, 'domain') - ############################################################################### # Output Manager tests # ############################################################################### @@ -538,8 +537,8 @@ def test_reconstruct_multipatch(dtype): A = Square('A',bounds1=bounds1, bounds2=bounds2_A) B = Square('B',bounds1=bounds1, bounds2=bounds2_B) - connectivity = [((0,1,1),(1,1,-1))] - patches = [A,B] + patches = [A, B] + connectivity = [((0, 1, 1), (1, 1,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') Va = ScalarFunctionSpace('Va', A) diff --git a/psydac/cad/__init__.py b/psydac/cad/__init__.py index 6f39d9a03..419109b64 100644 --- a/psydac/cad/__init__.py +++ b/psydac/cad/__init__.py @@ -3,10 +3,3 @@ # LICENSE file or go to https://github.com/pyccel/psydac/blob/devel/LICENSE # # for full license details. # #---------------------------------------------------------------------------# - -__all__ = ['geometry'] - -from psydac.cad import geometry -from psydac.cad import cad -from psydac.cad import gallery -from psydac.cad import utils diff --git a/psydac/cad/cad.py b/psydac/cad/cad.py index abc39baab..cd18ff4db 100644 --- a/psydac/cad/cad.py +++ b/psydac/cad/cad.py @@ -3,6 +3,8 @@ # LICENSE file or go to https://github.com/pyccel/psydac/blob/devel/LICENSE # # for full license details. # #---------------------------------------------------------------------------# +from typing import Iterable + import numpy as np from psydac.fem.splines import SplineSpace @@ -12,26 +14,48 @@ from psydac.ddm.cart import DomainDecomposition #============================================================================== -def translate(mapping, displ): +def translate(mapping : SplineMapping, displ : Iterable[float]): + """ + Translate a CAD geometry by a given vector displacement. + + Translate a single-patch CAD geometry (given as a spline or NURBS mapping) + by adding a given displacement vector to the control points of the mapping. + + Parameters + ---------- + mapping : SplineMapping + The discrete mapping to be translated, which represents a CAD geometry. + + displ : Iterable[float] + The vector displacement by which to translate the geometry. + + Returns + ------- + SplineMapping + A new discrete mapping representing the translated geometry. + """ + assert isinstance(mapping, SplineMapping) + assert isinstance(displ, Iterable) + displ = np.array(displ) - assert( mapping.pdim == len(displ) ) + assert mapping.pdim == len(displ) pdim = mapping.pdim space = mapping.space control_points = mapping.control_points - fields = [FemField( space ) for d in range( pdim )] + fields = [FemField(space) for d in range(pdim)] # Get spline coefficients for each coordinate X_i starts = space.coeff_space.starts ends = space.coeff_space.ends - idx_to = tuple( slice( s, e+1 ) for s,e in zip( starts, ends ) ) - for i,field in enumerate( fields ): + idx_to = tuple(slice(s, e+1) for s, e in zip(starts, ends)) + for i, field in enumerate(fields): idx_from = tuple(list(idx_to)+[i]) field.coeffs[idx_to] = control_points[idx_from] + displ[i] field.coeffs.update_ghost_regions() - return SplineMapping( *fields ) + return SplineMapping(*fields) #============================================================================== def elevate(mapping, axis, times): @@ -179,16 +203,3 @@ def refine(mapping, axis, values): return NurbsMapping( *fields ) return SplineMapping( *fields ) - - - -###################################### -if __name__ == '__main__': - from psydac.cad.geometry import Geometry - - geo = Geometry('square_0.h5') - mapping = geo.patches[0] - new = translate(mapping, [1., 0., 0.]) - - geo = Geometry(patches=[mapping, new]) - geo.export('square_mp.h5') diff --git a/psydac/cad/geometry.py b/psydac/cad/geometry.py index cf89c85d1..53f78d61e 100644 --- a/psydac/cad/geometry.py +++ b/psydac/cad/geometry.py @@ -8,20 +8,19 @@ # the topology i.e. connectivity, boundaries # For the moment, it is used as a container, that can be loaded from a file # (hdf5) -from itertools import product -from collections import abc -import string -import random -import yaml import os -import string -import random -import warnings +from typing import Iterable +from itertools import chain import numpy as np import h5py +import yaml from mpi4py import MPI +from sympde.topology import Domain, Interface, Line, Square, Cube, NCubeInterior, Mapping, NCube +from sympde.topology.basic import Union +from sympde.topology.callable_mapping import BasicCallableMapping + from psydac.fem.splines import SplineSpace from psydac.fem.tensor import TensorFemSpace from psydac.fem.partitioning import create_cart, construct_connectivity, construct_interface_spaces @@ -29,39 +28,48 @@ from psydac.linalg.block import BlockVectorSpace, BlockVector from psydac.ddm.cart import DomainDecomposition, MultiPatchDomainDecomposition +__all__ = ( + 'Geometry', + 'export_nurbs_to_hdf5', + 'import_geopdes_to_nurbs', + 'refine_knots', + 'refine_nurbs', +) -from sympde.topology import Domain, Interface, Line, Square, Cube, NCubeInterior, Mapping, NCube -from sympde.topology.basic import Union +NoneType = type(None) #============================================================================== class Geometry: """ Distributed discrete geometry that works for single and multiple patches. - The Geometry object can be created in two ways: - - case 1 : through a geometry file whos name can be given to the constructor - - case 2 : provide the ncells, the periodicity and the mapping objects of each patch. + + The Geometry object can be created in four ways: + - case 0 : providing a `Domain` to `__init__` with detailed parameters for each patch. + - case 1 : passing the path to a geometry file to `from_file`. + - case 2 : passing a `SplineMapping` to `from_discrete_mapping` (single patch). + - case 3 : passing a `Domain`, ncells, and periodicity to `from_topological_domain` (single or multi-patch). Parameters ---------- domain : Sympde.topology.Domain - The symbolic domain to be discretized. + The symbolic topological domain to be discretized. - ncells : list | tuple | dict - The number of cells of the discretized topological domain in each direction. + pdim : int + Number of physical dimensions of the Geometry object (pdim >= ldim). - periodic : list | tuple | dict - The periodicity of the topological domain in each direction. + ncells : dict[str, Iterable[int]] + The number of cells of the discretized domain in each direction. - mappings : dict - The Mapping of each patch. + periodic : dict[str, Iterable[bool]], optional + The periodicity of the topological domain in each direction. - filename: str - The path to the geometry file. + mappings : dict[str, BasicCallableMapping], optional + The discrete mappings of each patch. - comm: MPI.Comm + comm: MPI.Intracomm, optional MPI intra-communicator. - - mpi_dims_mask: list of bool + + mpi_dims_mask: Iterable[bool], optional True if the dimension is to be used in the domain decomposition (=default for each dimension). If mpi_dims_mask[i]=False, the i-th dimension will not be decomposed. @@ -71,68 +79,154 @@ class Geometry: _patches = [] _topology = None + def __init__(self, + domain : Domain, + *, + pdim : int, + ncells : dict[str, Iterable[int]], + mappings : dict[str, SplineMapping | None] = None, + periodic : dict[str, Iterable[bool]] = None, + comm : MPI.Intracomm = None, + mpi_dims_mask : Iterable[bool] = None): + + # Type checks + assert isinstance(pdim, int) + assert isinstance(domain, Domain) + assert isinstance(ncells, dict) + assert isinstance(mappings, (NoneType, dict)) + assert isinstance(periodic, (NoneType, dict)) + assert isinstance(comm, (NoneType, MPI.Intracomm)) + assert isinstance(mpi_dims_mask, (NoneType, Iterable)) + + # Extract info from domain + ldim : int = domain.dim + interior_names : list = domain.interior_names + set_interior_names = set(interior_names) + + # Check sanity of pdim + assert pdim >= ldim + + # Check sanity of ncells + assert set(ncells.keys()) == set_interior_names + assert all(len(n) == ldim for n in ncells.values()) + assert all(isinstance(ni, (int, np.integer)) for ni in chain(*ncells.values())) + assert all(ni > 0 for ni in chain(*ncells.values())) + + # Although we allow the iterable values in ncells to contain NumPy + # integers, we convert them to lists of Python integers for consistency + ncells = {patch: [int(ni) for ni in n] for patch, n in ncells.items()} + + # Check sanity of periodic + if periodic is None: + periodic = {patch: [False] * len(n) for patch, n in ncells.items()} + else: + assert set(periodic.keys()) == set_interior_names + assert all(len(p) == ldim for p in periodic.values()) + assert all(isinstance(pi, bool) for pi in chain(*periodic.values())) + + # Check sanity of mappings + if mappings is None: + mappings = {name : None for name in interior_names} + else: + assert set(mappings.keys()) == set_interior_names + assert all(isinstance(m, (BasicCallableMapping, NoneType)) for m in mappings.values()) + assert all(m.pdim == pdim for m in mappings.values() if m is not None) + + # Check sanity of mpi_dims_mask + if mpi_dims_mask is not None: + assert len(mpi_dims_mask) == ldim + assert all(isinstance(mask, bool) for mask in mpi_dims_mask) + + # Create a (multi-patch) domain decomposition + if len(domain) == 1: + #name = domain.name + name = interior_names[0] + ddm = DomainDecomposition( + ncells = ncells[name], + periods = periodic[name], + comm = comm, + mpi_dims_mask = mpi_dims_mask, + ) + else: + ddm = MultiPatchDomainDecomposition( + ncells = [ ncells[itr] for itr in interior_names], + periods = [periodic[itr] for itr in interior_names], + comm = comm, + ) + + # Add attributes to the new object + self._domain = domain + self._ldim = domain.dim + self._pdim = pdim + self._ncells = ncells + self._mappings = mappings + self._periodic = periodic + self._comm = comm + self._ddm = ddm + self._cart = None + #-------------------------------------------------------------------------- - # Option [1]: from a (domain, mappings) or a file + # Option [1]: from a file #-------------------------------------------------------------------------- - def __init__(self, domain=None, ncells=None, periodic=None, mappings=None, - filename=None, comm=None, mpi_dims_mask=None): - - # ... read the geometry if the filename is given - if filename is not None: - self.read(filename, comm=comm, mpi_dims_mask=mpi_dims_mask) - - elif domain is not None: - assert isinstance(domain, Domain) - assert isinstance(ncells, dict) - assert isinstance(mappings, dict) - if periodic is not None: - assert isinstance(periodic, dict) - - # ... check sanity - interior_names = domain.interior_names - mappings_keys = sorted(list(mappings.keys())) - - assert sorted(interior_names) == mappings_keys - # ... - - if periodic is None: - periodic = {patch: [False]*len(ncells_i) for patch, ncells_i in ncells.items()} - - self._domain = domain - self._ldim = domain.dim - self._pdim = domain.dim # TODO must be given => only dim is defined for a Domain - self._ncells = ncells - self._periodic = periodic - self._mappings = mappings - self._cart = None - self._is_parallel = comm is not None - - if len(domain) == 1: - #name = domain.name - name = interior_names[0] - self._ddm = DomainDecomposition(ncells[name], periodic[name], comm=comm, mpi_dims_mask=mpi_dims_mask) - else: - ncells = [ncells[itr] for itr in interior_names] - periodic = [periodic[itr] for itr in interior_names] - self._ddm = MultiPatchDomainDecomposition(ncells, periodic, comm=comm) + @classmethod + def from_file(cls, + filename : str, + *, + domain : Domain = None, + comm : MPI.Intracomm = None, + mpi_dims_mask : Iterable[bool] = None): - else: - raise ValueError('Wrong input') - # ... + """ + Create a Geometry instance from an HDF5 input file in Psydac format. - self._comm = comm + Parameters + ---------- + filename: str + The path to the geometry file. + + domain: sympde.topology.Domain, optional + The topological domain of the geometry. It must match the domain + defined in the file: same dimension, patch names, patches and + interfaces. Each patch must have a non-analytical mapping, to which + the corresponding spline mapping is attached; hence the same domain + should not be used with different geometry files. If not given, the + domain is created from the file. + + comm: MPI.Intracomm, optional + The MPI intra-communicator. + + mpi_dims_mask: Iterable[bool], optional + True if the dimension is to be used in the domain decomposition + (=default for each dimension). If mpi_dims_mask[i]=False, the i-th + dimension will not be decomposed. + + Returns + ------- + Geometry + The new instance. + + Raises + ------ + ValueError + If the given domain does not match the domain defined in the file, + or if one of its patches has no mapping or an analytical mapping. + """ + geo = super().__new__(cls) + geo.read(filename, comm=comm, mpi_dims_mask=mpi_dims_mask, domain=domain) + return geo #-------------------------------------------------------------------------- # Option [2]: from a discrete mapping #-------------------------------------------------------------------------- @classmethod def from_discrete_mapping(cls, mapping, *, comm=None, mpi_dims_mask=None, name=None): - """Create a geometry from one discrete mapping. + """ + Create a single-patch Geometry instance from one discrete mapping. Parameters ---------- - mapping : SplineMapping - The Mapping from the unit square to the physical domain. + mapping : BasicCallableMapping + The mapping from the unit square to the physical domain. comm : MPI.Comm MPI intra-communicator. @@ -141,62 +235,80 @@ def from_discrete_mapping(cls, mapping, *, comm=None, mpi_dims_mask=None, name=N True if the dimension is to be used in the domain decomposition (=default for each dimension). If mpi_dims_mask[i]=False, the i-th dimension will not be decomposed. - name : string - Optional name for the Mapping that will be created. - Needed to avoid conflicts in case several mappings are created + name : str + Optional name for the symbolic Mapping that will be created. + Needed to avoid conflicts in case several mappings are created. + + Returns + ------- + Geometry + The new instance. """ mapping_name = name if name else 'mapping' - dim = mapping.ldim - M = Mapping(mapping_name, dim = dim) + dim = mapping.ldim + M = Mapping(mapping_name, dim = dim) # this is a symbolic mapping domain = M(NCube(name = 'Omega', dim = dim, min_coords = [0.] * dim, max_coords = [1.] * dim)) M.set_callable_mapping(mapping) + pdim = mapping.pdim mappings = {domain.name: mapping} ncells = {domain.name: mapping.space.domain_decomposition.ncells} periodic = {domain.name: mapping.space.domain_decomposition.periods} - return Geometry(domain=domain, ncells=ncells, periodic=periodic, mappings=mappings, comm=comm, mpi_dims_mask=mpi_dims_mask) - + return Geometry(domain = domain, + pdim = pdim, + ncells = ncells, + periodic = periodic, + mappings = mappings, + comm = comm, + mpi_dims_mask = mpi_dims_mask) #-------------------------------------------------------------------------- # Option [3]: discrete topological line/square/cube #-------------------------------------------------------------------------- @classmethod def from_topological_domain(cls, domain, ncells, *, periodic=None, comm=None, mpi_dims_mask=None): + assert isinstance(domain, Domain) + interior = domain.interior if not isinstance(interior, Union): interior = [interior] for itr in interior: if not isinstance(itr, NCubeInterior): - msg = "Topological domain must be an NCube;"\ + msg = "The topological domain of each patch must be an NCube;"\ " got {} instead.".format(type(itr)) raise TypeError(msg) - mappings = {itr.name:None for itr in interior} + mappings = {itr.name : None for itr in interior} + pdim = next(iter(interior)).dim if isinstance(ncells, (list, tuple)): - ncells = {itr.name:ncells for itr in interior} + ncells = {itr.name : ncells for itr in interior} if periodic is None: - periodic = [False]*domain.dim + periodic = [False] * domain.dim else: if len(interior) > 1 and True in periodic: + import warnings msg = "Discretizing a multipatch domain with a periodic flag is not advised -- continue at your own risk." # [MCP 18.12.2025] the following line may be causing a strange error in the CI (MPI tests for macos-14/Python 3.10) # warnings.warn(msg, Warning) warnings.warn(msg, UserWarning) - if isinstance(periodic, (list, tuple)): - periodic = {itr.name:periodic for itr in interior} - - geo = Geometry(domain=domain, mappings=mappings, ncells=ncells, periodic=periodic, comm=comm, mpi_dims_mask=mpi_dims_mask) + periodic = {itr.name : periodic for itr in interior} - return geo + return Geometry(domain = domain, + pdim = pdim, + ncells = ncells, + periodic = periodic, + mappings = mappings, + comm = comm, + mpi_dims_mask = mpi_dims_mask) #-------------------------------------------------------------------------- @property @@ -227,10 +339,6 @@ def domain(self): def ddm(self): return self._ddm - @property - def is_parallel(self): - return self._is_parallel - @property def mappings(self): return self._mappings @@ -238,49 +346,48 @@ def mappings(self): def __len__(self): return len(self.domain) - def read(self, filename, comm=None, mpi_dims_mask=None): + def read(self, filename, comm=None, mpi_dims_mask=None, domain=None): # ... check extension of the file - basename, ext = os.path.splitext(filename) - if not(ext == '.h5'): + _, ext = os.path.splitext(filename) + if ext != '.h5': raise ValueError('> Only h5 files are supported') # ... - # read the topological domain - domain = Domain.from_file(filename) - connectivity = construct_connectivity(domain) - - if len(domain)==1: - interiors = [domain.interior] + # read the topological domain, or check the given one against it + file_domain = Domain.from_file(filename) + if domain is None: + domain = file_domain else: - interiors = list(domain.interior.args) + _check_domain_matches_file(domain, file_domain, filename) - if not(comm is None): - kwargs = dict( driver='mpio', comm=comm ) if comm.size > 1 else {} + connectivity = construct_connectivity(domain) + interiors = _get_interiors(domain) + if comm is not None: + kwargs = dict(driver='mpio', comm=comm) if comm.size > 1 else {} else: kwargs = {} - h5 = h5py.File( filename, mode='r', **kwargs ) - yml = yaml.load( h5['geometry.yml'][()], Loader=yaml.SafeLoader ) + h5 = h5py.File(filename, mode='r', **kwargs) + yml = yaml.load(h5['geometry.yml'][()], Loader=yaml.SafeLoader) ldim = yml['ldim'] pdim = yml['pdim'] - n_patches = len( yml['patches'] ) + n_patches = len(yml['patches']) # ... if n_patches == 0: - h5.close() - raise ValueError( "Input file contains no patches." ) + raise ValueError("Input file contains no patches.") # ... - # ... read patchs + # ... read patches mappings = {} ncells = {} periodic = {} - spaces = [None]*n_patches - for i_patch in range( n_patches ): + spaces = [None] * n_patches + for i_patch in range(n_patches): item = yml['patches'][i_patch] patch_name = item['name'] @@ -291,38 +398,36 @@ def read(self, filename, comm=None, mpi_dims_mask=None): degree = [int (p) for p in patch.attrs['degree' ]] periodic_i = [bool(b) for b in patch.attrs['periodic']] - knots = [patch['knots_{}'.format(d)][:] for d in range( ldim )] - space_i = [SplineSpace( degree=p, knots=k, periodic=P ) - for p,k,P in zip( degree, knots, periodic_i )] + knots = [patch['knots_{}'.format(d)][:] for d in range(ldim)] + space_i = [SplineSpace(degree=p, knots=k, periodic=P) + for p, k, P in zip(degree, knots, periodic_i)] spaces[i_patch] = space_i ncells [interiors[i_patch].name] = [sp.ncells for sp in space_i] periodic[interiors[i_patch].name] = periodic_i - self._cart = None if n_patches == 1: - self._ddm = DomainDecomposition(ncells[domain.name], periodic[domain.name], comm=comm, mpi_dims_mask=mpi_dims_mask) - ddms = [self._ddm] + ddm = DomainDecomposition(ncells[domain.name], periodic[domain.name], comm=comm, mpi_dims_mask=mpi_dims_mask) + ddms = [ddm] else: - ncells_ = [ncells[itr.name] for itr in interiors] - periodic = [periodic[itr.name] for itr in interiors] - self._ddm = MultiPatchDomainDecomposition(ncells_, periodic, comm=comm) - ddms = self._ddm.domains + ncells_ = [ncells[itr.name] for itr in interiors] + periodic = [periodic[itr.name] for itr in interiors] + ddm = MultiPatchDomainDecomposition(ncells_, periodic, comm=comm) + ddms = ddm.domains carts = create_cart(ddms, spaces) - g_spaces = {inter:TensorFemSpace( ddms[i], *spaces[i], cart=carts[i]) for i,inter in enumerate(interiors)} + g_spaces = {inter:TensorFemSpace(ddms[i], *spaces[i], cart=carts[i]) for i,inter in enumerate(interiors)} - for i,j in connectivity: - ((axis_i, ext_i), (axis_j , ext_j)) = connectivity[i, j] + for i, j in connectivity: minus = interiors[i] plus = interiors[j] - max_ncells = [max(ni,nj) for ni,nj in zip(ncells[minus.name],ncells[plus.name])] + max_ncells = [max(ni, nj) for ni, nj in zip(ncells[minus.name], ncells[plus.name])] g_spaces[minus].add_refined_space(ncells=max_ncells) - g_spaces[plus].add_refined_space(ncells=max_ncells) + g_spaces[plus ].add_refined_space(ncells=max_ncells) # ... construct interface spaces - construct_interface_spaces(self._ddm, g_spaces, carts, interiors, connectivity) + construct_interface_spaces(ddm, g_spaces, carts, interiors, connectivity) for i_patch in range( n_patches ): @@ -336,26 +441,26 @@ def read(self, filename, comm=None, mpi_dims_mask=None): tensor_space = g_spaces[interiors[i_patch]] if dtype == 'SplineMapping': - mapping = SplineMapping.from_control_points( tensor_space, - patch['points'][..., :pdim] ) + mapping = SplineMapping.from_control_points(tensor_space, + patch['points'][..., :pdim]) elif dtype == 'NurbsMapping': - mapping = NurbsMapping.from_control_points_weights( tensor_space, - patch['points'][..., :pdim], - patch['weights'] ) + mapping = NurbsMapping.from_control_points_weights(tensor_space, + patch['points'][..., :pdim], + patch['weights']) - mapping.set_name( item['name'] ) + mapping.set_name(item['name']) mappings[patch_name] = mapping - if n_patches>1: - coeffs = [[e._coeffs for e in mapping._fields] for mapping in mappings.values()] - spaces = [[coeffs_ij.space for coeffs_ij in coeffs_i] for coeffs_i in coeffs] - spaces = [BlockVectorSpace(*space) for space in spaces] - w_spaces = [sp.spaces[0] for sp in spaces] - space = BlockVectorSpace(*spaces, connectivity=connectivity) - w_space = BlockVectorSpace(*w_spaces, connectivity=connectivity) - v = BlockVector(space) - w = BlockVector(w_space) + # ... Update ghost regions within each patch and across interfaces + if n_patches > 1: + coeffs = [[e.coeffs for e in mapping.fields] for mapping in mappings.values()] + patch_spaces = [BlockVectorSpace(*[c_ij.space for c_ij in c_i]) for c_i in coeffs] + patch_spaces_w = [c_i[0].space for c_i in coeffs] + space = BlockVectorSpace(*patch_spaces , connectivity=connectivity) + space_w = BlockVectorSpace(*patch_spaces_w, connectivity=connectivity) + v = BlockVector(space) + w = BlockVector(space_w) mapping_list = list(mappings.values()) for i in range(n_patches): for j in range(len(coeffs[i])): @@ -377,7 +482,7 @@ def read(self, filename, comm=None, mpi_dims_mask=None): if isinstance(mapping, NurbsMapping): mapping.weights_field.coeffs.update_ghost_regions() - + # ... # ... close the h5 file h5.close() @@ -389,14 +494,15 @@ def read(self, filename, comm=None, mpi_dims_mask=None): patch.mapping.set_callable_mapping(F) # ... + self._domain = domain self._ldim = ldim self._pdim = pdim - self._mappings = mappings - self._domain = domain - self._comm = comm self._ncells = ncells + self._mappings = mappings self._periodic = periodic - self._is_parallel = comm is not None + self._comm = comm + self._ddm = ddm + self._cart = None # ... def export( self, filename ): @@ -437,7 +543,6 @@ def export( self, filename ): yml['patches'] = patches_info # ... - # ... topology topo_yml = self.domain.todict() # ... @@ -883,3 +988,49 @@ def _read_patch(lines, i_patch, n_lines_per_patch, list_begin_line): nrb = NURBS(knots, control=points, weights=W) return nrb + +#============================================================================== +def _get_interiors(domain : Domain) -> list: + """Return the list of interior patches of a (multi-patch) domain.""" + return [domain.interior] if len(domain) == 1 else list(domain.interior.args) + +def _get_interface_keys(domain : Domain) -> set: + """Return the interfaces of a domain as hashable tuples, with orientation.""" + interfaces = domain.interfaces + if interfaces is None: + return set() + if not isinstance(interfaces, Union): + interfaces = [interfaces] + return {(i.minus.domain.name, i.minus.axis, i.minus.ext, + i.plus .domain.name, i.plus .axis, i.plus .ext, i.ornt) for i in interfaces} + +def _check_domain_matches_file(domain : Domain, file_domain : Domain, filename : str): + """Raise a ValueError if the domain cannot be used with the geometry file.""" + + # Properties are compared one by one, because domain == file_domain raises + # an error for multi-patch domains (see pyccel/sympde#197) + + # Analytical mappings would be used in assembly instead of the spline ones + for patch in _get_interiors(domain): + if patch.mapping is None or patch.mapping.is_analytical: + raise ValueError(f"Patch {patch.name} must have a non-analytical " + f"mapping to be used with geometry file {filename}") + + if domain.dim != file_domain.dim: + raise ValueError(f"Domain dimension {domain.dim} does not match " + f"dimension {file_domain.dim} in file {filename}") + + if domain.interior_names != file_domain.interior_names: + raise ValueError(f"Patch names {domain.interior_names} do not match " + f"{file_domain.interior_names} in file {filename}") + + # Patch equality ignores the parametric bounds (see pyccel/sympde#198) + for patch, file_patch in zip(_get_interiors(domain), _get_interiors(file_domain)): + bounds = (patch .logical_domain.min_coords, patch .logical_domain.max_coords) + file_bounds = (file_patch.logical_domain.min_coords, file_patch.logical_domain.max_coords) + if bounds != file_bounds: + raise ValueError(f"Parametric bounds {bounds} of patch {patch.name} " + f"do not match {file_bounds} in file {filename}") + + if _get_interface_keys(domain) != _get_interface_keys(file_domain): + raise ValueError(f"Interfaces of the domain do not match those in file {filename}") diff --git a/psydac/cad/mesh/multipatch/magnet.h5 b/psydac/cad/mesh/multipatch/magnet.h5 index c18c688ec..9e497c2a0 100644 Binary files a/psydac/cad/mesh/multipatch/magnet.h5 and b/psydac/cad/mesh/multipatch/magnet.h5 differ diff --git a/psydac/cad/mesh/multipatch/square.h5 b/psydac/cad/mesh/multipatch/square.h5 index ae3f87a27..781c05feb 100644 Binary files a/psydac/cad/mesh/multipatch/square.h5 and b/psydac/cad/mesh/multipatch/square.h5 differ diff --git a/psydac/cad/mesh/multipatch/square_repeated_knots.h5 b/psydac/cad/mesh/multipatch/square_repeated_knots.h5 index 471b157f3..a17343290 100644 Binary files a/psydac/cad/mesh/multipatch/square_repeated_knots.h5 and b/psydac/cad/mesh/multipatch/square_repeated_knots.h5 differ diff --git a/psydac/cad/tests/test_geometry.py b/psydac/cad/tests/test_geometry.py index 6671d1870..1c9625488 100644 --- a/psydac/cad/tests/test_geometry.py +++ b/psydac/cad/tests/test_geometry.py @@ -9,11 +9,11 @@ import numpy as np from mpi4py import MPI -from sympde.topology import Domain, Line, Square, Cube, Mapping +from sympde.topology import Domain, Line, Square, Cube, Mapping, IdentityMapping from psydac.cad.geometry import Geometry, export_nurbs_to_hdf5, refine_nurbs from psydac.cad.geometry import import_geopdes_to_nurbs -from psydac.cad.cad import elevate, refine +from psydac.cad.cad import elevate, refine, translate from psydac.cad.gallery import quart_circle, circle from psydac.mapping.discrete import SplineMapping, NurbsMapping from psydac.mapping.discrete_gallery import discrete_mapping @@ -25,6 +25,7 @@ base_dir = os.path.dirname(os.path.realpath(__file__)) #============================================================================== +@pytest.mark.xdist_group('h5py') def test_geometry_2d_1(): ncells = [1,1] @@ -43,13 +44,13 @@ def test_geometry_2d_1(): ncells = {domain.name:ncells} # create a geometry from a topological domain and the dict of mappings - geo = Geometry(domain=domain, ncells=ncells, mappings=mappings) + geo = Geometry(domain=domain, pdim=2, ncells=ncells, mappings=mappings) # export the geometry geo.export('geo.h5') # read it again - geo_0 = Geometry(filename='geo.h5') + geo_0 = Geometry.from_file('geo.h5') # export it again geo_0.export('geo_0.h5') @@ -61,6 +62,7 @@ def test_geometry_2d_1(): geo_1.export('geo_1.h5') #============================================================================== +@pytest.mark.xdist_group('h5py') def test_geometry_2d_2(): # create a nurbs mapping @@ -92,13 +94,13 @@ def test_geometry_2d_2(): periodic = {domain.name:[space.periodic for space in mapping.space.spaces]} # create a geometry from a topological domain and the dict of mappings - geo = Geometry(domain=domain, ncells=ncells, periodic=periodic, mappings=mappings) + geo = Geometry(domain=domain, pdim=2, ncells=ncells, periodic=periodic, mappings=mappings) # export the geometry geo.export('quart_circle.h5') # read it again - geo_0 = Geometry(filename='quart_circle.h5') + geo_0 = Geometry.from_file('quart_circle.h5') # export it again geo_0.export('quart_circle_0.h5') @@ -111,6 +113,7 @@ def test_geometry_2d_2(): #============================================================================== # TODO to be removed +@pytest.mark.xdist_group('h5py') def test_geometry_2d_3(): # create a nurbs mapping @@ -144,6 +147,7 @@ def test_geometry_2d_3(): #============================================================================== # TODO to be removed +@pytest.mark.xdist_group('h5py') def test_geometry_2d_4(): # create a nurbs mapping @@ -204,11 +208,11 @@ def test_geometry_with_mpi_dims_mask(): # Create a geometry from a topological domain and the dict of mappings # Here we allow for any distribution of the domain: mpi_dims_mask is not passed - geo = Geometry(domain=domain, ncells=d_ncells, mappings=mappings, comm=comm) + geo = Geometry(domain=domain, pdim=3, ncells=d_ncells, mappings=mappings, comm=comm) geo.export('geo_mpi_dims.h5') # Read geometry file in parallel, but using mpi_dims_mask - geo_from_file = Geometry(filename='geo_mpi_dims.h5', comm=comm, mpi_dims_mask=mpi_dims_mask) + geo_from_file = Geometry.from_file(filename='geo_mpi_dims.h5', comm=comm, mpi_dims_mask=mpi_dims_mask) # Verify that the domain is distributed as expected assert geo_from_file.ddm.starts == expected_starts @@ -219,7 +223,6 @@ def test_geometry_with_mpi_dims_mask(): if rank == 0: os.remove('geo_mpi_dims.h5') - # ============================================================================== @pytest.mark.mpi def test_from_discrete_mapping(): @@ -227,7 +230,7 @@ def test_from_discrete_mapping(): comm = MPI.COMM_WORLD rank = comm.rank size = comm.size - mpi_dims_mask = [False, False, True] # We swill verify that this has an effect + mpi_dims_mask = [False, False, True] # We will verify that this has an effect ncells = [4, 8, 2 * size] # Each process should have two cells along x3 degree = [3, 3, 3] @@ -268,9 +271,81 @@ def test_from_topological_domain(): assert geo_from_domain.ddm.starts == expected_starts assert geo_from_domain.ddm.ends == expected_ends +# ============================================================================== +@pytest.mark.parametrize('npatches', [1, 2]) +def test_geometry_init_without_mappings(npatches: int) -> None: + + if npatches == 1: + domain = Square(name='A') + else: + A = Square('A', bounds1=(0, 1), bounds2=(0, 1)) + B = Square('B', bounds1=(1, 2), bounds2=(0, 1)) + domain = Domain.join(patches=[A, B], + connectivity=[((0, 0, 1), (1, 0, -1), 1)], + name='Omega') + + ncells = {name: [4, 4] for name in domain.interior_names} + geo = Geometry(domain, pdim=2, ncells=ncells) + + assert geo.mappings == {name: None for name in domain.interior_names} + +#============================================================================== +def make_domain(npatches: int, *, ornt: int = 1) -> Domain: + """Create a domain made of one or two unit squares, with generic mappings.""" + if npatches == 1: + return Mapping('G', dim=2)(Square('P')) + + A = Mapping('GA', dim=2)(Square('PA')) + B = Mapping('GB', dim=2)(Square('PB')) + return Domain.join([A, B], [((0, 0, 1), (1, 0, -1), ornt)], 'Omega') + +def export_geometry(domain: Domain, filename: str) -> None: + """Export a spline geometry on the domain, with patches side by side along x.""" + F = discrete_mapping('identity', ncells=[2, 2], degree=[2, 2]) + names = domain.interior_names + mappings = {name: translate(F, [float(i), 0.0]) for i, name in enumerate(names)} + ncells = {name: [2, 2] for name in names} + Geometry(domain, pdim=2, ncells=ncells, mappings=mappings).export(filename) + +#============================================================================== +@pytest.mark.parametrize('npatches', [1, 2]) +@pytest.mark.xdist_group('h5py') +def test_geometry_from_file_with_domain(npatches: int, tmp_path) -> None: + + domain = make_domain(npatches) + filename = str(tmp_path / 'geo.h5') + export_geometry(domain, filename) + + geo = Geometry.from_file(filename, domain=domain) + + # The given domain is kept, and the spline mappings are attached to it + assert geo.domain is domain + patches = [domain.interior] if npatches == 1 else domain.interior.args + for patch, F in zip(patches, geo.mappings.values()): + assert patch.mapping.get_callable_mapping() is F + +#============================================================================== +@pytest.mark.parametrize(('npatches', 'make_wrong_domain', 'message'), [ + (1, lambda: IdentityMapping('G', dim=2)(Square('P')), 'non-analytical'), + (1, lambda: Square('P'), 'non-analytical'), + (1, lambda: Mapping('G', dim=3)(Cube('P')), 'dimension'), + (1, lambda: Mapping('H', dim=2)(Square('P')), 'Patch names'), + (1, lambda: Mapping('G', dim=2)(Square('P', bounds1=(0, 2))), 'Parametric bounds'), + (2, lambda: make_domain(2, ornt=-1), 'Interfaces'), +], ids=['analytical', 'no-mapping', 'dimension', 'names', 'bounds', 'orientation']) +@pytest.mark.xdist_group('h5py') +def test_geometry_from_file_with_wrong_domain(npatches: int, make_wrong_domain, message: str, tmp_path) -> None: + + filename = str(tmp_path / 'geo.h5') + export_geometry(make_domain(npatches), filename) + + with pytest.raises(ValueError, match=message): + Geometry.from_file(filename, domain=make_wrong_domain()) + #============================================================================== @pytest.mark.parametrize( 'ncells', [[8,8], [12,12], [14,14]] ) @pytest.mark.parametrize( 'degree', [[2,2], [3,2], [2,3], [3,3], [4,4]] ) +@pytest.mark.xdist_group('h5py') def test_export_nurbs_to_hdf5(ncells, degree): # create pipe geometry @@ -288,7 +363,7 @@ def test_export_nurbs_to_hdf5(ncells, degree): export_nurbs_to_hdf5(filename, new_pipe) # read the geometry - geo = Geometry(filename=filename) + geo = Geometry.from_file(filename) domain = geo.domain min_coords = domain.logical_domain.min_coords @@ -324,9 +399,9 @@ def test_export_nurbs_to_hdf5(ncells, degree): #============================================================================== @pytest.mark.parametrize( 'ncells', [[8,8], [12,12], [14,14]] ) @pytest.mark.parametrize( 'degree', [[2,2], [3,2], [2,3], [3,3], [4,4]] ) +@pytest.mark.xdist_group('h5py') def test_import_geopdes_to_nurbs(ncells, degree): - filename = os.path.join(base_dir, "geo_Lshaped_C1.txt") L_shaped = import_geopdes_to_nurbs(filename) @@ -337,7 +412,7 @@ def test_import_geopdes_to_nurbs(ncells, degree): export_nurbs_to_hdf5(filename, L_shaped) # read the geometry - geo = Geometry(filename=filename) + geo = Geometry.from_file(filename) domain = geo.domain min_coords = domain.logical_domain.min_coords @@ -361,14 +436,6 @@ def test_import_geopdes_to_nurbs(ncells, degree): if isinstance(mapping, NurbsMapping): assert np.allclose(L_shaped.weights.flatten(), mapping._weights_field.coeffs.toarray(), 1e-15, 1e-15) -#============================================================================== -@pytest.mark.xfail -def test_geometry_1(): - - line = Geometry.as_line(ncells=[10]) - square = Geometry.as_square(ncells=[10, 10]) - cube = Geometry.as_cube(ncells=[10, 10, 10]) - #============================================================================== # CLEAN UP SYMPY NAMESPACE #============================================================================== diff --git a/psydac/feec/polar/examples/poisson_2d.py b/psydac/feec/polar/examples/poisson_2d.py index 688f821d1..1f1558d62 100644 --- a/psydac/feec/polar/examples/poisson_2d.py +++ b/psydac/feec/polar/examples/poisson_2d.py @@ -19,18 +19,14 @@ `mpirun -n 6 python poisson_2d.py -S -d 3 3 -t disk -D 0.2 -m 'C0conga'` """ -from dataclasses import dataclass -from time import sleep, time +# Imports are local to keep the CLI fast (e.g. --help) +# pylint: disable=import-outside-toplevel -import numpy as np -from mpi4py import MPI -from pyccel import lambdify -from sympy import Rational, cos, pi, sin, symbols +from dataclasses import dataclass +# Base classes are needed at class definition from psydac.feec.polar.examples.polar_model_2d import PolarModel2D from psydac.linalg.basic import LinearOperator -from psydac.linalg.stencil import StencilMatrix, StencilVector -from psydac.utilities.operators import Laplacian # ============================================================================== @@ -56,6 +52,8 @@ class Poisson2D(PolarModel2D): """ def __init__(self, domain_log, mapping, phi_log, rho_log): + from pyccel import lambdify + super().__init__(domain_log, mapping) self._domain_log = domain_log @@ -93,8 +91,10 @@ def disk(R, shift_D): and source term. """ + import numpy as np from sympde.topology.analytical_mapping import TargetMapping from sympde.topology.domain import Square + from sympy import cos, pi, sin, symbols domain_log = Square("Omega", bounds1=(0, R), bounds2=(0, 2 * np.pi)) params = dict(c1=shift_D * R * R, c2=0, k=0, D=shift_D) @@ -134,8 +134,12 @@ def target(): $\phi(x,y) = (1 - s^8)\sin(k_x(x - 0.5))\cos(k_y y)$. """ + import numpy as np from sympde.topology.analytical_mapping import TargetMapping from sympde.topology.domain import Square + from sympy import Rational, cos, pi, sin + + from psydac.utilities.operators import Laplacian domain_log = Square("Omega", bounds1=(0, 1), bounds2=(0, 2 * np.pi)) params = dict(c1=0, c2=0, k=Rational(3, 10), D=Rational(2, 10)) @@ -167,8 +171,12 @@ def czarny(): $\phi(x,y) = (1 - s^8)\sin(\pi x)\cos(\pi y)$. """ + import numpy as np from sympde.topology.analytical_mapping import CzarnyMapping from sympde.topology.domain import Square + from sympy import Rational, cos, pi, sin + + from psydac.utilities.operators import Laplacian domain_log = Square("Omega", bounds1=(0, 1), bounds2=(0, 2 * np.pi)) params = dict(c1=0, c2=0, eps=Rational(1, 5), b=Rational(7, 5)) @@ -240,6 +248,7 @@ def __init__(self, S, M, P, alpha): C0PolarProjection_V0, C1PolarProjection_U0, ) + from psydac.linalg.stencil import StencilMatrix assert isinstance(S, StencilMatrix) assert isinstance(M, StencilMatrix) @@ -257,6 +266,8 @@ def __init__(self, S, M, P, alpha): self.W0 = W0 def dot(self, x, out=None): + from psydac.linalg.stencil import StencilVector + if out is None: y = self.M._domain.zeros() else: @@ -339,7 +350,7 @@ def compute_errors(phi, phi_ref, M, S): ErrorDiagnostics Reference L2 and H1 norms and relative L2 and H1 errors. """ - + import numpy as np # L2 and H1 norms ref_l2_2 = M.dot_inner(phi_ref.coeffs, phi_ref.coeffs) @@ -386,6 +397,8 @@ def plot_solution(use_spline_mapping, model, ncells, periodic, V0_h, refine=10): Refinement factor for plotting. """ import matplotlib.pyplot as plt + import numpy as np + from mpi4py import MPI from psydac.cad.geometry import Geometry from psydac.ddm.cart import DomainDecomposition @@ -394,7 +407,7 @@ def plot_solution(use_spline_mapping, model, ncells, periodic, V0_h, refine=10): from psydac.utilities.utils import refine_array_1d if use_spline_mapping: - geometry = Geometry(filename="geo.h5", comm=MPI.COMM_SELF) + geometry = Geometry.from_file("geo.h5", comm=MPI.COMM_SELF) map_discrete = [*geometry.mappings.values()].pop() Vnew = map_discrete.space map_plot = map_discrete @@ -488,7 +501,10 @@ def run_poisson_2d( mpi_comm, # given by function 'parallel_run_from_cli' ): import os + from time import sleep, time + import numpy as np + from mpi4py import MPI from sympde.calculus import dot, grad from sympde.expr import BilinearForm, LinearForm, integral from sympde.topology import ScalarFunctionSpace, elements_of diff --git a/psydac/feec/tests/test_feec_conf_projectors_cart_2d.py b/psydac/feec/tests/test_feec_conf_projectors_cart_2d.py index cf64ecf27..1c45c4b33 100644 --- a/psydac/feec/tests/test_feec_conf_projectors_cart_2d.py +++ b/psydac/feec/tests/test_feec_conf_projectors_cart_2d.py @@ -43,7 +43,7 @@ def get_polynomial_function(degree, hom_bc_axes, domain): g0_y = (y - 0.75)**degree[1] expr = g0_x * g0_y - callable_function = lambdify(domain.coordinates, expr) + callable_function = lambdify(domain.coordinates, expr, modules=['numpy']) return expr, callable_function diff --git a/psydac/fem/tests/test_plot_field_2d.py b/psydac/fem/tests/test_plot_field_2d.py index a099aee02..bf870bdf7 100644 --- a/psydac/fem/tests/test_plot_field_2d.py +++ b/psydac/fem/tests/test_plot_field_2d.py @@ -46,16 +46,16 @@ def test_plot_field(use_scalar_field, use_multipatch): degree = [2, 2] A = Square('A',bounds1=(0.5, 1.), bounds2=(0, np.pi/2)) - mapping_1 = PolarMapping('M1',2, c1= 0., c2= 0., rmin = 0., rmax=1.) + mapping_1 = PolarMapping('M1', 2, c1= 0., c2= 0., rmin = 0., rmax=1.) D1 = mapping_1(A) if use_multipatch: B = Square('B',bounds1=(0.5, 1.), bounds2=(np.pi/2, np.pi)) - mapping_2 = PolarMapping('M2',2, c1= 0., c2= 0., rmin = 0., rmax=1.) - D2 = mapping_2(B) + mapping_2 = PolarMapping('M2', 2, c1= 0., c2= 0., rmin = 0., rmax=1.) + D2 = mapping_2(B) - connectivity = [((0,1,1),(1,1,-1))] - patches = [D1,D2] + patches = [D1, D2] + connectivity = [((0, 1, 1), (1, 1,-1), 1)] domain = Domain.join(patches, connectivity, 'domain') else: domain = D1 diff --git a/psydac/mapping/discrete_gallery.py b/psydac/mapping/discrete_gallery.py index 5075d8c85..e7a9f2820 100644 --- a/psydac/mapping/discrete_gallery.py +++ b/psydac/mapping/discrete_gallery.py @@ -3,6 +3,8 @@ # LICENSE file or go to https://github.com/pyccel/psydac/blob/devel/LICENSE # # for full license details. # #---------------------------------------------------------------------------# +from typing import Iterable + import numpy as np from mpi4py import MPI @@ -33,6 +35,11 @@ 'spherical_shell', ) +__all__ = ( + 'get_available_mappings', + 'discrete_mapping', +) + class Collela3D( Mapping ): _expressions = {'x':'2.*(x1 + 0.1*sin(2.*pi*x1)*sin(2.*pi*x2)) - 1.', @@ -40,78 +47,143 @@ class Collela3D( Mapping ): 'z':'2.*x3 - 1.'} #============================================================================== -def discrete_mapping(mapping, ncells, degree, **kwargs): +def get_available_mappings(ldim): + """ + Get a list of `mapping` values accepted as argument to `discrete_mapping`. + + Parameters + ---------- + ldim : int + The number of logical dimensions of the topological domain. + + Returns + ------- + tuple + All the accepted values for the `mapping` parameter. + """ + assert isinstance(ldim, int), f'ldim must be int, got {type(ldim).__name__} instead' + assert ldim > 0, f'ldim must be > 0, got {ldim} instead' + + if ldim == 2: + return ('identity', 'collela', 'circle', 'annulus', 'quarter_annulus', + 'target', 'czarny') + elif ldim == 3: + return ('identity', 'collela', 'spherical_shell') + else: + return () - comm = kwargs.pop('comm', MPI.COMM_WORLD) - return_space = kwargs.pop('return_space', False) +#============================================================================== +def discrete_mapping(mapping, ncells, degree, *, + comm = MPI.COMM_WORLD, + return_space = False): + """ + Create a SplineMapping by interpolating one of the available analytical mappings. + + Parameters + ---------- + mapping : str + The name of the mapping. See `available_mappings` to get the options. + + ncells : Iterable[int] + The number of cells along each logical dimension. + + degree : Iterable[int] + The spline degree along each logical dimension. + + comm : MPI.Intracomm, optional + The MPI intracommunicator. + + return_space : bool, optional + Whether this function should also return the discrete space it creates. + + Returns + ------- + map_discrete : SplineMapping + The spline mapping created. + + space : TensorFemSpace + The space of the components of the spline mapping. + Only returned if `return_space` is True. + """ + # Check types + assert isinstance(mapping, str) + assert isinstance(ncells, Iterable) + assert isinstance(degree, Iterable) + assert isinstance(comm, MPI.Intracomm) or comm is None + assert isinstance(return_space, bool) + + # Check consistency of ncells and degree + assert all(isinstance(n, int) and n >= 1 for n in ncells) + assert all(isinstance(d, int) and d >= 0 for d in degree) + assert len(ncells) == len(degree) mapping = mapping.lower() - dim = len(ncells) - if dim not in [2, 3]: + ldim = len(ncells) + if ldim not in [2, 3]: raise NotImplementedError('Only 2D and 3D mappings are available') # ... - if dim == 2: + if ldim == 2: # Input parameters if mapping == 'identity': - map_symbolic = IdentityMapping('M', dim=dim) + map_symbolic = IdentityMapping('M', dim=ldim) limits = ((0, 1), (0, 1)) periodic = (False, False) elif mapping == 'collela': default_params = dict(k1=1.0, k2=1.0, eps=0.1) - map_symbolic = CollelaMapping2D('M', dim=dim, **default_params) + map_symbolic = CollelaMapping2D('M', dim=ldim, **default_params) limits = ((0, 1), (0, 1)) periodic = (False, False) elif mapping == 'circle': default_params = dict(rmin=0.0, rmax=1.0, c1=0.0, c2=0.0) - map_symbolic = PolarMapping('M', dim=dim, **default_params) + map_symbolic = PolarMapping('M', dim=ldim, **default_params) limits = ((0, 1), (0, 2*np.pi)) periodic = (False, True) elif mapping == 'annulus': default_params = dict(rmin=0.0, rmax=1.0, c1=0.0, c2=0.0) - map_symbolic = PolarMapping('M', dim=dim, **default_params) + map_symbolic = PolarMapping('M', dim=ldim, **default_params) limits = ((1, 4), (0, 2*np.pi)) periodic = (False, True) elif mapping == 'quarter_annulus': default_params = dict(rmin=0.0, rmax=1.0, c1=0.0, c2=0.0) - map_symbolic = PolarMapping('M', dim=dim, **default_params) + map_symbolic = PolarMapping('M', dim=ldim, **default_params) limits = ((1, 4), (0, np.pi/2)) periodic = (False, False) elif mapping == 'target': default_params = dict(c1=0, c2=0, k=0.3, D=0.2) - map_symbolic = TargetMapping('M', dim=dim, **default_params) + map_symbolic = TargetMapping('M', dim=ldim, **default_params) limits = ((0, 1), (0, 2*np.pi)) periodic = (False, True) elif mapping == 'czarny': default_params = dict(c2=0, b=1.4, eps=0.3) - map_symbolic = CzarnyMapping('M', dim=dim, **default_params) + map_symbolic = CzarnyMapping('M', dim=ldim, **default_params) limits = ((0, 1), (0, 2*np.pi)) periodic = (False, True) else: raise ValueError("Required 2D mapping not available") - elif dim == 3: + elif ldim == 3: # Input parameters if mapping == 'identity': - map_symbolic = IdentityMapping('M', dim=dim) + map_symbolic = IdentityMapping('M', dim=ldim) limits = ((0, 1), (0, 1), (0, 1)) periodic = ( False, False, False) elif mapping == 'collela': - map_symbolic = Collela3D('M', dim=dim) + map_symbolic = Collela3D('M', dim=ldim) limits = ((0, 1), (0, 1), (0, 1)) periodic = ( False, False, False) elif mapping == 'spherical_shell': - map_symbolic = SphericalMapping('M', dim=dim) + map_symbolic = SphericalMapping('M', dim=ldim) limits = ((1, 4), (0, np.pi), (0, np.pi/2)) periodic = ( False, False, False) diff --git a/pyproject.toml b/pyproject.toml index 9492defb3..36c290a3e 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -31,7 +31,7 @@ dependencies = [ 'termcolor', # Our packages from PyPi - 'sympde == 0.19.3', + 'sympde == 0.20.0', 'pyccel >= 2.2.3', 'gelato == 0.12', @@ -76,3 +76,37 @@ psydac = "psydac.cmd.main:psydac_command" setup = ['--default-library=static'] dist = ['--include-subprojects'] install = ['--only=python-modules'] + +[tool.pylint.main] +py-version = '3.10' +ignore-paths = ['psydac/pyccel/.*', '.*/__e?pyccel__/.*', '.*/__psydac__/.*'] +ignored-modules = ['mpi4py'] # C extension, members not visible to Pylint + +[tool.pylint.basic] +no-docstring-rgx = '^(_|test_|Test)' + +[tool.pylint.'messages control'] +disable = [ + # formatting owned by black / isort + 'line-too-long', + 'trailing-whitespace', + 'bad-indentation', + 'multiple-statements', + 'missing-final-newline', + 'trailing-newlines', + 'wrong-import-order', + # style choices for mathematical code + 'invalid-name', + 'too-many-arguments', + 'too-many-positional-arguments', + 'too-many-locals', + 'too-many-branches', + 'too-many-statements', + 'too-many-instance-attributes', + 'too-few-public-methods', + 'too-many-nested-blocks', + # pytest fixtures shadow outer names by design + 'redefined-outer-name', + 'duplicate-code', + 'fixme', +]