-
Notifications
You must be signed in to change notification settings - Fork 116
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 non extruded ocean boundaries #334
Conversation
This is done by resetting the Coordinate field to OriginalCoordinate (squashing the domain back to its original form) before going into the adapt. After the adapt the fs mesh movement is reapplied. This may result in conservation issues.
Forgets to divide \rho in ct_rhs. fs mesh movement with variable_density still needs to be fixed (or options check to say that it doesn't work).
In geostrophic pressure initialisation after an adapt: don't assume the number of vector fields in new_states is the same as in old. In enforce_discrete_properties: don't assume the new states already have all meshes mentioned in the options.
Does not work yet.
…ion of LinearMomentum+MeshMovement and LinearMomentum+MeshMovement+MeshAdaptivity
Instead of duplicating a surface mesh indicated by surface_elements write more generic create_parallel_redudant_mesh that takes in any mesh (which could of course be a surface mesh). Also move this to fefields and mostly restore vertical extrapolation module.
Picking apart horizontal_picker and add code that creates horizontal interpolation mesh in parallel, non-extruded case (actual interpolation not hooked up yet). Spherical case has been simplified by re-adding boundary condition such that its surface element list now maps from elements in the horizontal interpolation mesh to facets of the 3D mesh.
…o parallel-non-extruded-ocean-boundaries
* Since we transform the physical domain back to its original shape during an adapt, but didn't update the physical position of detectors they ended up outside of their element but more importantly can end up outside the transformed physical (either local or global) domain. In that case we aren't able to find back a containing element in the new mesh for all detectors after an adapt. Now, the detectors are also moved (transformed back) with the mesh before the adapt, and moved back afterwards. This is done by keeping the local coordinates fixed and reinterpolating to a position within the element. Once the element has been transformed back this should gives us back the original pre-adapt detector position. * Deleted the second half of search_detectors as it has an incorrect and potentially dangerous call to picker_inquire with global=.true. This combination only works if we provide picker_inquire with the same positions to search on all processes, whereas search_for_detectors naively assumed each process could just ask for its own missing detector locations and picker_inquire would somehow figure it out. If in the first (and now only) call to picker_inquire with global=.false. not all detectors have been found, i.e. some detectors have left the local part of the physical domain that we can see, this should be picked up in a subsequent call to distribute_detectors(). In the case where we are absolutely sure that this doesn't (shouldn't) happen (i.e. all detectors can be found locally), we should check that search_detectors() indeed has found all detectors. This happens in adapt, where the local physical domain does not change during the adapt. We now check this in all runs (not just with --enable-debugging) as it is cheap to do and we will always crash in zoltan if we have lost detectors.
The update vector and K matrix should only be exchanged in exchange_detectors _during_ the RK stages - otherwise these are note even allocated. Rather than guessing based on obscure logic, like "present(positions): oh you must be calling from adaptivity cause we currently provide a specific positions in the call to distribute_detectors which passes positions on to us optionally, so you probably don't need to" - explicitly tell the routine whether we want exchanging or not.
…-ocean-boundaries
Add additional has_surface_field() and has_scalar_surface_field() (latter for scalar surface fields stored under the bc of a vector field) versions where now both the bc and the surface field can be specified by name (previously only existed with bc specified by number). This avoids trying to extract it with a stat= argument - if the surface field does not exist a nonzero stat is returned, but using sfield = extract_surface_field() it would still be copying from a non-initialized pointer location into sfield.
cf74c05
to
d9801eb
Compare
This all looks good! My only question would be whether we should add the Also, the |
Thanks, @cianwilson - much appreciated. I've changed ownership of these tests. On my laptop beaker_blob runs in just 46s. I'd be interested to know if it's much slower for you cause that could indicate solver issues... |
Not horrendously slow admittedly. 1:30 for me on |
This PR includes the following:
Other fixes related to free surface mesh movement in this PR (originally in #130):
OriginalCoordinates
before the adapt so that the domain with the (assumed flat) initial free surface is adapted. After the adapt, the mesh is moved back according to the interpolated free surface fieldFixes to particles with a moving free surface: