Improve Geometry class and fix CAD functions in psydac.cad

.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

.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
 

AGENTS.md

@@ -0,0 +1,104 @@
+# 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 in new code
+- meson-python (Meson build backend) 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)
+- sphinx + numpydoc + nbsphinx + myst-parser for documentation
+  (see `docs/requirements.txt`; build with `make -C docs html`)
+
+## Code standards
+- Follow PEP 8 style guide whenever possible
+- Docstrings of public functions and classes follow Numpydoc conventions
+- pylint + black + isort for linting and formatting
+- Computational kernels to be accelerated with Pyccel named as `<module>_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"
+- In docstrings, use single backticks for names and short code,
+  and double backticks only for code expressions
+- Write new code in clean style (e.g. no superfluous parentheses, f-strings),
+  even if the surrounding legacy code does not
+- Do not add author credits in code or comments;
+  credit goes in commit messages and PR descriptions
+
+## 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`)
+- Parallel tests are marked with `@pytest.mark.mpi` (from pytest-mpi) and run with
+  `psydac test --mpi`; serial tests run with `psydac test`
+- Test file names follow `test_<module>.py`
+- Test function names follow `test_<function>` or `test_<class>_<method>`
+- Parametrize unit tests with `@pytest.mark.parametrize` to minimize code duplication
+- Keep run time of tests at a minimum
+- Aim for 100% coverage on newly committed code
+- Every bugfix comes with a regression test that fails before the fix
+- Use distinct sizes in each direction (e.g. 3 x 5 cells, not 4 x 4),
+  so that swapped indices are detected
+- Write files only to the `tmp_path` fixture (`mpi_tmp_path` in MPI tests),
+  never to the working directory
+- After changing a module, also run the tests of the subpackages which use it
+  (e.g. `psydac/api/tests` after changing `psydac/cad`)
+
+## 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/mapping - spline and NURBS mappings, and a gallery of discrete mappings
+- 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
+
+## Workflow
+- PRs are squash-merged and the PR description becomes the commit message:
+  keep it complete, and write "(fixes #N)" next to each fix
+- Every PR adds entries to `CHANGELOG.md` under Added/Changed/Fixed/Removed
+- Format of the CHANGELOG entries:
+  . `-   #<PR> : <change>`, followed by `, see #<issue>` if there is one
+  . `-   [DEVELOPER] <change>` for changes which only concern developers
+- Issues describe bugs as they appear on `devel`, with a reproducer that runs on `devel`
+- Do not commit changes to submodule pointers (`subprojects/`) unless intended
+
+## Always do
+- Before claiming done, run these checks on the Python files you changed:
+  . Run `python -m black --check $(git diff --name-only --diff-filter=d devel -- "*.py")`
+  . Run `python -m isort --check $(git diff --name-only --diff-filter=d devel -- "*.py")`
+  . Run `pylint --disable=all --enable=unused-import $(git diff --name-only --diff-filter=d devel -- "*.py")`
+- Until the code base is reformatted, black and isort fail on most existing files:
+  do not reformat code outside your change. The pylint check must pass.
+- 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`

CHANGELOG.md

@@ -5,28 +5,53 @@ All notable changes to this project will be documented in this file.
 ## Unreleased
 
 ### Added
+
+-   #527 : Add constructor `Geometry.from_file` with option `domain` to reuse an existing SymPDE domain, after checking that it matches the file
+-   #527 : Add function `get_available_mappings` and option `mpi_dims_mask` to `discrete_mapping` in `psydac.mapping.discrete_gallery`
+-   #527 : Add property `mpi_dims_mask` to `DomainDecomposition`
 -   #596 : Add `topetsc` method to `psydac.linalg.basic.ZeroOperator`
 -   #596 : Add `out` argument to `array_to_psydac` and `vec_topetsc`
 -   #576 : Add module `psydac.feec.polar.conga_projections` with broken FEEC polar projections in 2D
 -   #576 : Add 2D Poisson and TE Maxwell examples on polar mapped domains in `psydac/feec/polar/examples`
 -   #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] Configure isort in `pyproject.toml`, with the profile of black
+-   [DEVELOPER] Add `AGENTS.md` with developer rules
 
 ### Fixed
 
+-   #527 : Require `sympde==0.20.0` which changes how multipatch interfaces are defined
+-   #527 : Fix `psydac.cad.elevate` and `refine` for NURBS mappings, see #608
+-   #527 : Keep the NURBS weights in `psydac.cad.translate`, see #609
+-   #527 : Fix type checks of input mappings in module `psydac.cad.cad`, see #610
+-   #527 : Fix `psydac.cad.elevate` and `refine` for distributed mappings, see #611
+-   #527 : Check the weights field in `NurbsMapping.__init__`, see #616
+-   #527 : Fix reading multipatch geometries with `pdim=1`, see #617
+-   #527 : Fix `psydac.cad.elevate` and `refine` for mappings with `pdim=1`, see #618
+-   #527 : Fix `psydac.cad.refine` when inserting a value which is already a knot, see #619
+-   #527 : Keep the `mpi_dims_mask` of a `DomainDecomposition` in `TensorFemSpace.add_refined_space` and `psydac.cad.cad.refine`, see #622
+-   #527 : Fix `local_ncells` of the decomposition returned by `DomainDecomposition.refine`, see #623
+-   #527 : Raise a `ValueError` in `Geometry.export` if a patch has no spline mapping, instead of an `AttributeError` after creating the file, see #631
+-   #527 : Use the physical dimension of the mappings in `Geometry.from_topological_domain` and `from_discrete_mapping`
+-   #527 : Pair the patches of a geometry file with the domain interiors by name, not by position
+-   #527 : Avoid setting the `shape` of NumPy arrays (deprecated since NumPy 2.5) in `KroneckerLinearSolver`
+-   #527 : Fix error message in `psydac.api.discretization.get_max_degree_of_one_space`
 -   #596 : Allow export to VTK with a prescribed grid containing exact domain bounds
 -   #596 : Fix global geometric projectors with complex `dtype`
 -   #576 : Fix bug in `TensorFemSpace.eval_field` caused by round-off at MPI subdomain boundaries
 -   #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
@@ -35,13 +60,22 @@ All notable changes to this project will be documented in this file.
 
 ### Changed
 
+-   #527 : Pass a domain decomposition `ddm` and physical dimension `pdim` to `Geometry.__init__`; remove arguments `ncells`, `periodic`, `filename`, `comm` and `mpi_dims_mask`
+-   #527 : Use the domain decomposition of the mapping's space in `Geometry.from_discrete_mapping`; remove arguments `comm` and `mpi_dims_mask`
+-   #527 : Key `Geometry.mappings` by the names of the domain interiors
+-   #527 : Return Cartesian (not homogeneous) control points from the NURBS functions in `psydac.cad.gallery`, see #608
+-   #527 : Keep the given SymPDE domain in `discretize(domain, filename=...)`, after checking that it matches the geometry file
+-   #527 : Do not import the submodules of `psydac.cad` in its `__init__.py`
+-   #527 : Use latest version of Igakit (commit dalcinl/igakit@2a0c9ee of 2026/10/14) which supports NumPy >= 2.5
+-   #527 : Raise a `ValueError` with an explanation, instead of an `AssertionError`, in `psydac.fem.partitioning.partition_coefficients` when an MPI process owns too few coefficients
 -   #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`
 -   #579 : Require `pyccel>=2.2.3` which can compile all kernels with C
 -   #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
@@ -52,6 +86,9 @@ All notable changes to this project will be documented in this file.
 
 ### Removed
 
+-   #527 : Remove method `Geometry.read`, use constructor `Geometry.from_file` instead
+-   #527 : Remove property `Geometry.is_parallel`
+
 ## [1.0.0] - 2026-01-19
 
 The first official release on PyPI. A complete overhaul since version 0.1.

examples/old_examples/poisson_2d_mapping.py

@@ -632,7 +632,7 @@ def main(*, test_case, ncells, degree, nquads,
         map_discrete = SplineMapping.from_mapping(V, map_analytic)
         # Write discrete geometry to HDF5 file
         t0 = time()
-        geometry = Geometry.from_discrete_mapping(map_discrete, comm=mpi_comm)
+        geometry = Geometry.from_discrete_mapping(map_discrete)
         geometry.export('geo.h5')
         t1 = time()
         timing['export'] += t1-t0

examples/performance/compare_3d_matrix_assembly_speed.py

@@ -80,7 +80,7 @@ def make_square_torus_geometry_3d(ncells, degree, comm=None):
     mapping = SquareTorus('S')
     map_discrete = SplineMapping.from_mapping(V_h, mapping.get_callable_mapping())
 
-    geometry = Geometry.from_discrete_mapping(map_discrete, comm=comm)
+    geometry = Geometry.from_discrete_mapping(map_discrete)
 
     if mpi_rank == 0:
         if not os.path.isdir('geometry'):

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')
 

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)

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)

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)

psydac/api/discretization.py

@@ -110,7 +110,7 @@ def get_max_degree_of_one_space(Vh):
     Returns
     -------
     list[int]
-        The maximum polynomial degre of Vh with respect to each coordinate.
+        The maximum polynomial degree of Vh with respect to each coordinate.
 
     """
 
@@ -125,7 +125,7 @@ def get_max_degree_of_one_space(Vh):
         return [max(p) for p in zip(*degree)]
 
     else:
-        raise TypeError(f'Type({V}) not understood')
+        raise TypeError(f'Type {type(Vh).__name__} not understood')
 
 
 def get_max_degree(*spaces):
@@ -454,11 +454,8 @@ def discretize_space(V, domain_h, *, degree=None, multiplicity=None, knots=None,
 
     connectivity = construct_connectivity(domain)
     if isinstance(domain_h, Geometry) and all(domain_h.mappings.values()):
-        # from a discrete geoemtry
-        if interiors[0].name in domain_h.mappings:
-            mappings  = [domain_h.mappings[inter.name] for inter in interiors]
-        else:
-            mappings  = [domain_h.mappings[inter.logical_domain.name] for inter in interiors]
+        # from a discrete geometry
+        mappings  = [domain_h.mappings[inter.name] for inter in interiors]
 
         # Get all the FEM spaces from the mapping and convert their coeff_space at the dtype needed
         spaces    = [change_dtype(m.space, dtype) for m in mappings]
@@ -578,7 +575,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)