Introduce source and target spaces for operator domains by keileg · Pull Request #1747 · pmgbergen/porepy

keileg · 2026-08-14T13:13:27Z

UPDATE: REVIEW GUIDANCE

An initial QC on the overall design was provided by @Yuriyzabegaev, leading to significant improvements in the code. All these comments should now have been taken care of. I suggest that the next review stage focus (at least) on the following points:

The syntax for operator spaces and when they must be declared. I have left some GH review comments in constitutive_laws.py and test_constitutive_laws.py to highlight what I think are the main points.
Higher level documentation: Should the OperatorSpace be documented beyond usage in code, tutorials, and examples? The natural place to put it would be equations.ipynb, but that would be more natural with a slight reconfiguration of that file. Let me know what you think.
Also consider the other open GH comments. Most of these point to changes that are more detailed but still worth highlighting.

Overview

This PR ensures that all Ad operators have a source and target (a better terminology would be 'domain' and 'range', but both these words are already overloaded in PorePy and Python in general). This streamlines the existing code, where some but far from all operators know where they are defined. Also, tests of operator size compatibility are added.

Technical changes

This PR introduces the following changes (in increasing level of complexity - this is also a suggested order of review):

A GridEntity enum is now used to encode grid quantity. This replaces strings "cell" etc.
A DomainType enum is encodes different types of domain. This replaces strings "subdomain" etc. It was necessary to construct some special domain types to cover unclear or flexible cases, see code documentation.
An OperatorSpace class, which wraps dof information, a domain (list of grids, mortar- og boundary grids) and domain type.
All operators must be constructed with a source and target OperatorSpace. Rolling this out required quite a few changes in the models folder, but these are mainly shallow edits.
The Ad Operations class (which now resides in its own module) infers its source and range from the source and range of its children operands. That is, if A and B are two operands and x is an operation, the resulting operator A x B will have its source and target automatically identified. If the triplet {A, B, x} are incompatible, an error will be raised at the time of construction (in contrast, this error will currently be raised during parsing, often with a less then intelligible error message). Also, for the classes DenseAdArray and SparseAdArray additional checks of array sizes are added to enhance the chances of catching errors.

Tests have been updated and new tests (in particular unit tests for the OperatorSpace class) have been added.

Resolves #1601

Types of changes

What types of changes does this PR introduce to PorePy?
Put an x in the boxes that apply.

Minor change (e.g., dependency bumps, broken links).
Bugfix (non-breaking change which fixes an issue).
New feature (non-breaking change which adds functionality).
Breaking change (fix or feature that would cause existing functionality to not work as expected).
Testing (contribution related to testing of existing or new functionality).
Documentation (contribution related to adding, improving, or fixing documentation).
Maintenance (e.g., improve logic and performance, remove obsolete code).
Other:

Checklist

Put an x in the boxes that apply or explain briefly why the box is not relevant.

The documentation is up-to-date.
Static typing is included in the update.
This PR does not duplicate existing functionality.
The update is covered by the test suite (including tests added in the PR).
If new skipped tests have been introduced in this PR, pytest was run with the --run-skipped flag.

review-notebook-app · 2026-08-14T13:13:31Z

Check out this pull request on

See visual diffs & provide feedback on Jupyter Notebooks.

Powered by ReviewNB

Yuriyzabegaev

I have read through everything, great job!

I only looked at the tests less thoroughly. Is there a test that tries to combine (add / multiply / apply / etc) incompatible operators and catches the error early? (before the actual AD parsing). Sorry if I didn't see it.

See my comment about self-mapping operators, that's the only thing worth thinking about (and not sure that it's worth doing anything about).

Other than that, just minor comments, everything looks good!

keileg · 2026-08-26T10:23:30Z

-        operator and its inverse;
-    test_trace and test_divergence: Operators for discrete traces and divergences;
-    test_ad_discretization_class: test for AD discretizations.
+"""Tests for the AD representations of grid-related operators: projections between


I ended up refactoring this file as part of augmenting some of the tests contained. Summarized, the changes are:

Tests for subdomain, interface and boundary projections are grouped into separate classes.

Each of these sets of tests use parametrization and fixtures/common setups to a much larger extent.

The tests also cover the use of OperatorSpaces.
Only the last change is considers the core changes in this PR.

I have verified that the de facto test coverage did not decrease (in some places in was improved). In opinion, this file need not be a focal point for the review.

keileg · 2026-08-27T07:42:15Z

+        # here.
+        t_hf_inv = one / (
+            pp.ad.SparseArray(
+                d_n_by_dist, source=hf_source, target=pp.ad.OperatorSpace.waived()


Note: This is one of few places where the waived space is used (the other cases are similarly related to the differentiable tpfa functions). A cost-benefit analysis made it very clear it is not worth the effort to extend the machinery.

Would you say this is in principle doable in the future?

My feeling is, not unless we extend the OperatorSpace concept. With half-transmissibilities we get two items for each internal face, one for each external. Altertatively, we get one item per cell-face combination. Both can be implemented, but both will require processing grid topology is some sense.

…ction

Old dict format is kept for backwards compatibility

…operator space Should be a cosmetic change, the code is currently not in use

Removed dict-like structure

The return annotations of get_row_dof_info and get_col_dof_info reference GridEntities, which is imported only under TYPE_CHECKING. Neither module had "from __future__ import annotations", so on Python 3.13 and earlier the annotation is evaluated when the class body runs, and importing porepy fails with a NameError. Python 3.14 defers annotation evaluation (PEP 649), which is why this was only caught by the Python 3.13 leg of the test workflow. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The shortened GridEntities return annotations let the signatures of get_row_dof_info and get_col_dof_info fit on a single line. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Co-authored-by: Ivar Stefansson <ivar.stefansson@uib.no>

Use entities instead of dof terminology, as the latter is only relevant for columns

…ified

In infer_source_target, the case of two scalar operands tested whether either operand carried a domain-bearing space. Since the branch is only entered when both sources have the scalar domain type, those tests were always false and the domain-preserving code was unreachable. A Scalar constructed with domains does not have the scalar space, and is already handled by the general elementwise branches, which is now covered by tests. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

keileg requested review from IvarStefansson and jwboth as code owners August 14, 2026 13:13