dealii
dealii copied to clipboard
dealii
[Post 9.3] Triangulation::get_manifold(): assert that manifold has been attached
Currently the function Triangulation::get_manifold() returns for a manifold id for which no manifold has been attached a flat manifold. In most cases this is not what you want is rather a sign of a mistake. Possible case: when you create a copy of a triangulation.
This is an incompatible change but would simplify life (debugging!) of users much!
This is actually a pretty big breaking change. In particular, if someone is misusing manifold ids to stand in for some intermediate value, their code is now going to break.
We might want to consider this for 10.0 (there really should be a manifold for each manifold id) but I'd be surprised if this passes the test suite.
I agree that this is a change in behavior, but I also agree with @peterrum that this will usually point to a bug. The current behavior dates back to a time when we had geometry information only attached to boundaries, and used the boundary_id for this purpose. In other words, the boundary_id was used for both indicating what kinds of boundary conditions might hold on a given part of the boundary, and what its geometry was. To make this work, one had to assume a straight boundary (=flat manifold) for boundary ids for which no specific boundary object was attached.
We no longer have this kind of ambiguity, and I would support a change in behavior (even in 9.3) because I agree that the current behavior could be considered a bug. It's also easy to fix in backward compatible ways for everyone who now runs into the new assertion.
Let's see what the test suite has to say about this: /rebuild
People abusing manifold ids may just explicitly attach a flat manifold, so it makes perfect sense for me to have this change.
Quite many tests have failed. I went through all of them. Let me share some insight here what the reasons for the failing tests was:
- Some functions in the namespace
GridGeneratordid not attach some manifolds due to implicit assumptions (easy fix). - Many tests call
tr.reset_manifold(0);. The consequence of this function call is that the manifold with the given id is removed from the triangulation (and not replaced by anything new (flat) one). For the time being, I have changed the behavior so that the manifold is replaced by a flat manifold. Personally, I would not like to change the behavior of this function (I did it only because it was the quickest fix). I'll do the changes once we have agreed how to proceed. - There are some tests where
GridTools::copy_boundary_to_manifold_id(triangulation);is called without attaching new manifolds. - The function
GridTools::get_coarse_mesh_description()does not attach any manifolds so that the new mesh created based on this description does not have any manifolds. The user needs to attach these manually. - The function
GridGenerator::extract_boundary_mesh()does not attach any manifolds (as documented in the documentation). However, since the function might refine the mesh, manifolds have to be attached before calling the function even if they are only flat!
Overall, I think the implicit assumption is really error prone: we had one instance at the workshop (when the convert_hypercube_to_simplex_mesh() was used); I know at least two different persons having problems due to this; and I had at least two times to debug quite a long time to get the reason why the results did not match the expectation. I am not sure what the right way is to proceed here...
Many tests call tr.reset_manifold(0);. The consequence of this function call is that the manifold with the given id is removed from the triangulation (and not replaced by anything new (flat) one). For the time being, I have changed the behavior so that the manifold is replaced by a flat manifold. Personally, I would not like to change the behavior of this function (I did it only because it was the quickest fix). I'll do the changes once we have agreed how to proceed.
I think this is the right approach.
There are some tests where GridTools::copy_boundary_to_manifold_id(triangulation); is called without attaching new manifolds.
That's a mistake, if the manifolds are not set.
The function GridTools::get_coarse_mesh_description() does not attach any manifolds so that the new mesh created based on this description does not have any manifolds. The user needs to attach these manually.
Is there a reason not to have manifolds on the coarse mesh? I'd copy the manifolds over.
The function GridGenerator::extract_boundary_mesh() does not attach any manifolds (as documented in the documentation). However, since the function might refine the mesh, manifolds have to be attached before calling the function even if they are only flat!
This should be fixed by creating a SubManifold class that takes a Manifold and interprets it as a Manifold on a grid with dimension dim-1. Let me open an issue with this. For the moment we can create flat manifolds, and attach them to the extracted mesh (if a manifold is not attached yet), so that users won't see a difference in the behaviour: if you know how to attach a submanifold to the codim one grid, attach it (and the extract function does nothing), if not you leave it as it is, and the extract function attaches a submanifold.
That's a mistake, if the manifolds are not set.
Is there a reason not to have manifolds on the coarse mesh? I'd copy the manifolds over.
Copying manifolds from one triangulation to another one is an ongoing issue. For example, TransfiniteInterpolationManifold stores the reference of the underlying triangulation and currently we don't have a way to update that reference. There might be a solution by @tamiko and @kronbichler but that is not in place yet. Let's copy/clone the manifolds in a follow-up PR.
For the moment we can create flat manifolds, and attach them to the extracted mesh (if a manifold is not attached yet), so that users won't see a difference in the behaviour:
Done.
Can you document this change of behavior, and same for the reset_all_manifolds() function below? I think it would be useful to state that the function doesn't actually remove the entry from the array, but just overwrite them with a flat manifold.
Done.
Of course, an even better solution would be to only allow this function to be called if no entity in the triangulation still uses this manifold id. I think that was the intent of the function: To un-register the manifold_id.
Let's postpone this to a follow-up PR.
This has the potential to be quite disruptive. Let's do this right after the release.
I think we are too far into the release cycle to make a large change like this one. Lets do it in the near future but not in time for 9.4.
I also think we should postpone this to after the release. Let me adjust the milestone.
I have updated the PR. We are still long before release so that we could attempt to merge this now.