Skip to content
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

Parallel Compatible Cross-Mesh Interpolation #3067

Merged
merged 3 commits into from
Sep 1, 2023

Commits on Sep 1, 2023

  1. Configuration menu
    Copy the full SHA
    bb04a73 View commit details
    Browse the repository at this point in the history
  2. Fix bugged input ordering interpolation action

    I had managed to not test the whole of the interpolation API when I
    originally implemented this and it turned out that I would return None
    when, for example, calling interpolate(f_input_ordering, P0DG).
    ReubenHill committed Sep 1, 2023
    Configuration menu
    Copy the full SHA
    2ddd133 View commit details
    Browse the repository at this point in the history
  3. Implement cross mesh interpolation

    This implements cross-mesh interpolation for interpolating into function
    spaces which have point evaluation nodes. Full documentation is added to
    the manual.
    
    Two new classes are added to interpolation.py: CrossMeshInterpolator and
    SameMeshInterpolator, whilst Interpolator is made an abstract base
    class. To maintain the API of interpolate and Interpolator, the __new__
    method of Interpolator is  overridden to return an instance of the
    appropriate subclass.
    
    Two new keyword arguments are added to interpolate and Interpolator to
    allow for target meshes which extend outside the source mesh domain::
    see their docstrings for details.
    
    Docstrings are also added for some undocumented keyword arguments of
    Interpolator and interpolate.
    
    A full suite of tests is found in
    tests/regression/test_interpolate_cross_mesh.py
    
    Note that as part of this, I have changed the error by VertexOnlyMesh
    when points are outside the domain from a ValueError to a
    VertexOnlyMeshMissingPointsError
    
    Details of how this works from the PR description:
    
    We use VertexOnlyMesh as an intermediary for the global locations of the
    point evaluation nodes of the target function space:
    
    1. We get the point evaluation node locations using the method described
    in the interpolation from external data section of the manual. This will
    have the parallel domain decomposition of the target mesh.
    2. Next we create a VertexOnlyMesh (A) at those locations within the
    source mesh such that we inherit the source mesh's parallel domain
    decomposition.
    3. We interpolate our expression in our source function space onto a
    P0DG function space on VertexOnlyMesh (A), which has the effect of point
    evaluating at the target function space node locations.
    4. This VertexOnlyMesh (A) has an input_ordering VertexOnlyMesh (B)
    whose vertices have the ordering and parallel domain decomposition of
    the target function space global node locations. We interpolate from
    P0DG on (A) onto P0DG on (B). Under the hood, this is an SF reduce
    operation which moves the point evaluations from (A) to (B).
    5. We now have a Function on the input_ordering VertexOnlyMesh (B) which
    has point evaluations from our source mesh function space at the target
    mesh function space node locations. These are in the correct order and
    have the correct domain decomposition. We can therefore safely set the
    dat of a function in our target function space to the values of this
    function.
    
    For this to work for the general case, we would need to create a
    VertexOnlyMesh at the global quadrature points of the target function
    space, which is rather more complicated than the work I've done here.
    
    Some important notes:
    
     - This does not require one mesh to be a structured refinement of the
     other. This, for example, should allow you to solve a PDE on two
     different unstructured meshes, one of which is finer than the other,
     and directly check the difference in solutions by interpolating from
     one mesh to the other within Firedrake.
     - Crucially, this is entirely parallel compatible!
     - Since we can interpolate onto immersed manifolds we can perform line,
     surface and volume integrals by interpolating onto a mesh which has the
     domain of integration as an immersed manifold. This is demonstrated in
     the test_interpolate_cross_mesh.py test suite.
    
    Other notes:
    
    - The VertexOnlyMesh required is stored in the Interpolator, as are the
     underlying Interpolators
    - For interpolation between mixed spaces, I create sub_interpolators for
    each space and evaluate them as necessary when calling interpolate
    - Interpolation into a mixed space therefore requires the function space
     being interpolated from to be another mixed space. I throw a
     NotImplementedError if not.
    
    Regarding the manual:
    Note that not all the comments in the manual file are included in the
    literalinclude text of the manual, instead they are approximately
    rewritten as prose in the manual.
    ReubenHill committed Sep 1, 2023
    Configuration menu
    Copy the full SHA
    396cd76 View commit details
    Browse the repository at this point in the history