p-multigrid: remove custom interpolation code - #5289
Open
pbrubeck wants to merge 21 commits into
Open
Conversation
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
July 29, 2026 10:25
a79a5c0 to
2c9e9a1
Compare
pbrubeck
marked this pull request as ready for review
July 29, 2026 10:31
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
July 29, 2026 10:33
2c9e9a1 to
8b84e62
Compare
pbrubeck
force-pushed
the
pbrubeck/fix/dual-enriched
branch
from
July 29, 2026 10:33
cf08365 to
a9ee8fa
Compare
This was referenced Jul 29, 2026
leo-collins
reviewed
Jul 29, 2026
pbrubeck
commented
Jul 29, 2026
pbrubeck
commented
Jul 29, 2026
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
July 29, 2026 18:30
8b84e62 to
41c33c8
Compare
pbrubeck
force-pushed
the
pbrubeck/fix/dual-enriched
branch
2 times, most recently
from
July 29, 2026 23:29
168ee7b to
6072eee
Compare
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
July 29, 2026 23:29
3d81cdf to
33be7ba
Compare
pbrubeck
changed the base branch from
pbrubeck/fix/dual-enriched
to
pbrubeck/hdivcurl-delta
July 29, 2026 23:33
pbrubeck
force-pushed
the
pbrubeck/hdivcurl-delta
branch
from
August 18, 2026 22:22
a992610 to
537d1b5
Compare
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
August 18, 2026 22:22
2cdcc16 to
59e617b
Compare
pbrubeck
commented
Aug 18, 2026
pbrubeck
commented
Aug 24, 2026
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
August 24, 2026 22:24
e9e3a1e to
59e617b
Compare
An element whose functionals do not all evaluate on the same points is a direct sum, and has a dual basis only on each of its sub-elements, not one weight tensor on one point set. Reconstructing it as a runtime Quadrature space only ever wanted the points, so ask for those directly. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Assembling a two-form Interpolate matrix-free passed the same bcs as both the row and the column bcs of the ImplicitMatrixContext, and left on_diag at its default of True. Neither holds for interpolation between two different spaces: applying the target space's bcs to the source side raised "Cofunction defined on incompatible FunctionSpace", and on_diag made multTranspose try to write the row values into the column vector. Split the bcs by the space they are defined on, as the assembled path already does, and set on_diag from whether the two spaces coincide. A bc on a component of a mixed space is defined on an IndexedProxyFunctionSpace rather than on the space itself, so compare against the top-level parent, which BCBase.function_space now returns on request. Coarsen the p-multigrid nullspace with a locally built interpolation matrix, rather than through getNestSubMatrix, which only the removed MixedInterpolationMatrix provided. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
pbrubeck
commented
Aug 26, 2026
pbrubeck
commented
Aug 26, 2026
The dual basis of a restricted element is a Concatenate, and a Cofunction dual argument sums over the index that concatenates it. Forming that contraction first leaves the Concatenate indexed by a reduction index that no assignment variable carries, so unconcatenate has nothing to split it against and it reaches argument factorisation, which cannot break it up. Split the pair of the dual coefficient and the tabulation first, as assembly already does for coefficient evaluation, and contract each block over the basis indices of its own block. Interpolating between facet-restricted spaces, as p-multigrid does under FacetSplitPC, raised FactorisationError before this. The error escaped a PETSc callback, which left the preconditioner half built, and destroying it later segfaulted. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Contributor
|
This is failing an awful lot of tests. Is it ready for review? |
Contributor
Author
I'll be fixing those today. I'll ask for review when it's passing. In principle it should be a small fix to make everything pass |
Contributor
|
Can you check if this fixes #5391? I just debugged that down to it being (I think) to do with |
Every worker of the `tests/tsfc` session segfaults at its first `compile_form`, in `petsctools.cite` -> `PETSc.Sys.registerCitation`, because `PetscCitationsList` is still NULL: PETSc was never initialised. `petsc4py.lib.ImportPETSc` puts the extension module straight into `sys.modules`, bypassing the `petsc4py/PETSc.py` shim that calls `PETSc._initialize`. Once that has happened no later `from petsc4py import PETSc` initialises anything. The same fingerprint shows in the firedrake sessions: on `main` the conftest's `import pyadjoint` initialises PETSc before `firedrake/__init__.py` does, and `petsctools.init` warns that it is too late; on this branch that warning never appears. This does not reproduce locally, so record who calls `ImportPETSc` and report the state of PETSc once collection has imported every test module. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
This reverts commit 30edc16.
Contributor
Author
It does fix it. |
Contributor
Author
This is ready for review |
fiat pbrubeck/fix/dual-enriched (this PR's FIAT dependency) landed review fixes at firedrakeproject/fiat@46e28722. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01BjZ6kYK1tmLfm6wUxXwVtG
connorjward
requested changes
Sep 2, 2026
connorjward
left a comment
connorjward
left a comment
Contributor
There was a problem hiding this comment.
In general I think this is a wonderful change.
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
September 7, 2026 14:47
bd96c37 to
10666d2
Compare
The Cofunction dual argument contracts over the nodes, but the basis indices it sums were not declared as quadrature indices, so spectral.flatten dropped them from the sum factorisation tail ordering and the contraction ran dense. Declare them alongside the points. This pays off with fiat#268, which cancels Deltas over the indices a monomial contracts, so the node sum against the dual basis of a point evaluation element collapses rather than being carried out. That refactorisation also folds a constant out of an IndexSum, which test_refactorise pins. Interpolation between facet-restricted elements on a hexahedron costs 6.8x fewer flops at degree 11, growing as p^2.89 rather than p^4.20, matching what the coefficient path already achieved. Set on_diag in ImplicitMatrixContext.__init__ rather than as a class attribute the caller overwrites, and give BCBase a parent_function_space property in place of the function_space(parent=True) flag. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_015Vnd4VHtANbfmiMxSPSXSG
pbrubeck
force-pushed
the
pbrubeck/remove-pmg-parloop
branch
from
September 7, 2026 14:48
10666d2 to
df69102
Compare
The exact solution of that Poisson problem lies in the space, so what the test measures is round-off, and 1e-15 left it no margin: the error came out at 4.8e-16 here and 1.05e-15 in CI, on either side of the bound. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01F1K5upgeZTKWj9m4tHt9Vn
The folded constant this expected came from collect_monomials rebuilding every monomial's rest from a flattened factor list, which costs flops on form assembly. fiat now reflattens only when a Delta cancels, so the monomial keeps the grouping the recursion chose. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0158RTiFnSyhgKbsbz11JMak
connorjward
previously approved these changes
Sep 8, 2026
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
3 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Removes ~1050 lines of hand-written interpolation machinery from
firedrake/preconditioners/pmg.py—StandaloneInterpolationMatrix,MixedInterpolationMatrix,kronmxv,get_permutation_to_nodal_elements,evaluate_dualand friends — and replaces the p-multigrid transfer operators withthe generic
assemble(interpolate(TrialFunction(cV), fV)). The helpers that are stillneeded move to
firedrake/preconditioners/fdm.py.Depends on firedrakeproject/fiat#268
Fixes #5391 by removing the transfer matrix cache
🤖 Generated with Claude Code