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

Replace .at() interpolation in 2D callbacks #373

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

Conversation

cpjordan
Copy link
Contributor

.at() raises an error when the number of processors used is large, use a VertexOnlyMesh instead for interpolation

- .at() raises an error when the number of processors used is large, use a VertexOnlyMesh instead for interpolation
@cpjordan
Copy link
Contributor Author

Prevents the Firedrake issue #2790 whilst Firedrake issue #3080 is waiting to be resolved. This method is slower than .at() so an alternative is to add a try except with .at() before moving onto using a VertexOnlyMesh.

There are some remaining direct applications of .at() but they are for 3D examples where either the callbacks are only being used in serial or where the problem is so large that .at() is functioning fine.

@cpjordan cpjordan marked this pull request as ready for review July 31, 2024 08:38
@thangel thangel requested review from tkarna and swarder October 9, 2024 14:52
Copy link
Contributor

@tkarna tkarna left a comment

Choose a reason for hiding this comment

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

Thanks, looks good. Only minor comments

thetis/callback.py Outdated Show resolved Hide resolved
thetis/callback.py Outdated Show resolved Hide resolved
thetis/callback.py Outdated Show resolved Hide resolved
thetis/callback.py Outdated Show resolved Hide resolved
Comment on lines 534 to 544
self.vom = VertexOnlyMesh(self.mesh, detector_locations, redundant=True)
self.function_spaces = {}
for field_name in field_names:
field = solver_obj.fields[field_name]
if isinstance(field.function_space().ufl_element(), VectorElement):
P0DG = VectorFunctionSpace(self.vom, "DG", 0)
P0DG_input_ordering = VectorFunctionSpace(self.vom.input_ordering, "DG", 0)
else:
P0DG = FunctionSpace(self.vom, "DG", 0)
P0DG_input_ordering = FunctionSpace(self.vom.input_ordering, "DG", 0)
self.function_spaces[field_name] = (P0DG, P0DG_input_ordering)
Copy link
Contributor

Choose a reason for hiding this comment

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

minor: there's quite a lot of duplicated code in the vom creation and evaluation. Could plausibly be refactored in to helper functions for example

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I've created a new function in utility.py L1159: vom_interpolator_function. It generates and returns the Functions only as suggested so we remove all the unnecessary Function duplication I had. We could probably go a step further and create a vom_interpolator class which also performs the interpolation as well rather than repeating these lines:

f_at_points, f_at_input_points = self.interp_functions[field_name]
f_at_points.interpolate(field)
f_at_input_points.interpolate(f_at_points)
return f_at_input_points.dat.data_ro[:]

Though maybe it's fine how it is for now. When Firedrake issue 3080 is dealt with we can probably revert to .at() anyway.

thetis/callback.py Outdated Show resolved Hide resolved
Comment on lines +756 to +757
f_at_points.interpolate(field)
f_at_input_points.interpolate(f_at_points)
Copy link
Contributor

Choose a reason for hiding this comment

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

Do you know if this vom evaluation works for 2D mesh on the sphere (in which case the point coordinates are 3d x,y,z)? The whole eval_func option was added to enable at()-like evaluation in this case (eval_func was a modified version of Firedrake's at)

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Unsure as I only usually do boring plain 2D meshes! I mainly fixed because of the parallelisation issues for big simulations. Happy to look into it in the future or perhaps someone in the Firedrake group might have answer for this already.

@cpjordan
Copy link
Contributor Author

The current version of Firedrake is failing, which seems to be causing all the Thetis tests to fail. Waiting for this to be updated.

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.

2 participants