Skip to content

RestrictedFunctionSpace: support Fieldsplit, multigrid, and python PC #4169

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 10 commits into
base: master
Choose a base branch
from

Conversation

pbrubeck
Copy link
Contributor

@pbrubeck pbrubeck commented Mar 28, 2025

Description

This PR adds support for pc_type fieldsplit, mg, and python with RestrictedFunctionSpace.

Apparently setting name of the RFS also changes the name of the split in the PETSc solver option.
With the default name, the options prefix becomes fieldsplit_Restricted_on_boundary, instead of fieldsplit_x.
This PR ensures that the restricted spaces that are created within solvers preserve the name of the unrestricted parent space so that the same set of fieldsplit options works for both restrict=False and restrict=True.

@pbrubeck pbrubeck requested a review from ksagiyam March 28, 2025 15:05
@pbrubeck pbrubeck force-pushed the pbrubeck/fix/restricted-dmhooks branch from 4193431 to 1b20212 Compare March 28, 2025 18:22
Comment on lines 419 to 420
V = type(self).make_function_space(mesh, element, name=name,
boundary_set=V_parent.boundary_set)
Copy link
Contributor Author

@pbrubeck pbrubeck Mar 31, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ksagiyam I want V.reconstruct() to preserve the boundary_set of the original V. I have added the code to construct a WithGeometry(RestrictedFunctionSpace) above in make_function_space, but maybe this is not the right place.

Comment on lines +244 to +249
def get_coordinates(V):
coords = V.mesh().coordinates
if V.boundary_set:
W = V.reconstruct(element=coords.function_space().ufl_element())
coords = firedrake.Function(W).interpolate(coords)
return coords
Copy link
Contributor Author

@pbrubeck pbrubeck Mar 31, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ksagiyam I need to restrict V.mesh().coordinates, is there a better way to do this?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we really need to do that here? Normally, we shouldn't need to do that as restricted function spaces are just regular function spaces with different DoF ordering.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

prolong and restrict need the coordinates in the same dof ordering as the restricted spaces.

@pbrubeck pbrubeck changed the title dmhooks: support RestrictedFunctionSpace RestrictedFunctionSpace: support Fieldsplit, multigrid, and python PC Mar 31, 2025
Comment on lines +244 to +249
def get_coordinates(V):
coords = V.mesh().coordinates
if V.boundary_set:
W = V.reconstruct(element=coords.function_space().ufl_element())
coords = firedrake.Function(W).interpolate(coords)
return coords
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we really need to do that here? Normally, we shouldn't need to do that as restricted function spaces are just regular function spaces with different DoF ordering.

@@ -376,12 +376,18 @@ def make_function_space(cls, mesh, element, name=None):
if mesh is not topology:
# Create a concrete WithGeometry or FiredrakeDualSpace on this mesh
new = cls.create(new, mesh)

if boundary_set:
new = RestrictedFunctionSpace(new, boundary_set=boundary_set)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

RestrictedFunctionSpace here is topological, so this should be done around line 374, I think.
At some point, I think we need to tidy up reconstruct and make_function_space.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The way things are now, the geometric resticted space is wrapping a geometric unrestricted space:

from firedrake import *
mesh = UnitSquareMesh(1, 1)
V = FunctionSpace(mesh, "CG", 1)
problem = LinearVariationalProblem(inner(TrialFunction(V), TestFunction(V))*dx, 0, Function(V), bcs=DirichletBC(V, 0, [1]), restrict=True)

problem.u_restrict.function_space()

WithGeometry(RestrictedFunctionSpace("WithGeometry(FunctionSpace(<firedrake.mesh.MeshTopology object at 0x7ee09b678f0>, FiniteElement('Lagrange', triangle, 1), name=None), <Mesh #0>)", name=None, boundary_set=frozenset({1})), Mesh(VectorElement(FiniteElement('Lagrange', triangle, 1), dim=2), 0))

@pbrubeck pbrubeck force-pushed the pbrubeck/fix/restricted-dmhooks branch from 9df7a04 to 7a07dff Compare March 31, 2025 16:35
@pbrubeck pbrubeck force-pushed the pbrubeck/fix/restricted-dmhooks branch from 9a449c5 to f165b4b Compare March 31, 2025 22:12
try:
return cache[key]
except KeyError:
Vc = firedrake.VectorFunctionSpace(mesh, element)
Vc = V.collapse().reconstruct(element=finat.ufl.VectorElement(element, dim=mesh.geometric_dimension()))
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
Vc = V.collapse().reconstruct(element=finat.ufl.VectorElement(element, dim=mesh.geometric_dimension()))
Vc = V.collapse().reconstruct(element=finat.ufl.VectorElement(element))

@@ -876,6 +886,14 @@ class RestrictedFunctionSpace(FunctionSpace):
If using this class to solve or similar, a list of DirichletBCs will still
need to be specified on this space and passed into the function.
"""
def __new__(cls, function_space, boundary_set=frozenset(), name=None):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Use RestrictedFunctionSpace from functionspace.py (i.e. the function) instead

cell = mesh.topology.ufl_cell()
if len(kwargs) > 0 or element.cell != cell:
element = element.reconstruct(cell=cell, **kwargs)

V = type(self).make_function_space(mesh, element, name=name)
for i in reversed(indices):
V = V.sub(i)

if boundary_set:
V = RestrictedFunctionSpace(V, boundary_set=boundary_set, name=V.name)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For separate discussion: should RestrictedFunctionSpace be hidden from the user?

r"""Reconstruct this :class:`.WithGeometryBase` .

:kwarg mesh: the new :func:`~.Mesh` (defaults to same mesh)
:kwarg element: the new :class:`finat.ufl.FiniteElement` (defaults to same element)
:kwarg boundary_set: boundary subdomain labels defining a new
:func:`~.RestrictedFunctionSpace` (defaults to same boundary_set)
:kwarg name: the new name (defaults to None)
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It seems that we can reconstruct named spaces without too much worry, as the new name defaults to None

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants