fix the quadratic element issue in TriangulatorInterface
closes #4144
The original approach uses node ID as the key to find the corresponding midpoint. The ID seems to have changed during triangulation.
Job Coverage, step Generate coverage on 33866f9 wanted to post the following:
Coverage
| e50c9b | #4145 33866f | ||||
|---|---|---|---|---|---|
| Total | Total | +/- | New | ||
| Rate | 63.42% | 63.42% | +0.00% | 100.00% | |
| Hits | 74838 | 74842 | +4 | 2 | |
| Misses | 43166 | 43163 | -3 | 0 | |
This comment will be updated on new commits.
This invalidates the comment explaining
all_midpoints:// [...] all_midpoints[{p, m}] will be the mth midpoint // location following after point p (when traversing a triangle // counter-clockwise)And it doesn't change uses of
all_midpointsthat clearly make that assumption too, explicitly queryingelem->point(n)forn<3to get the point for the key. None of them needed to be changed?I think we need a new unit test to demonstrate the bug and cover the fix.
The only change this PR makes is to store the node point coordinates instead of the node id of the first node when a segment (i.e. two nodes) is created. Later on, this point coordinates are used as p in all_midpoints[{p, m}] instead of using _mesh.point() and the stored node id. Thus, I do not think anything else needs to change; the only issue of the original code is that the node ids change during triangulation.
I agree that a unit test would help.
Ah, I get it now. I saw segment_midpoints_keys and jumped to the conclusion that we were storing midpoints, rather than non-midpoint-points-used-as-keys-to-find-midpoints. We've got to get a clearer name for that, at the very least, and clearer comments explaining what we're doing and why.
To do that, of course, we have to understand why. I don't like fixing the effects of a problem whose cause I don't understand. It's like putting a bucket to catch water under a drip from the ceiling. Your carpet is no longer getting wet, that's great, but the drip shouldn't have been there and it might be damaging your drywall if you consider the problem solved instead of investigating further.
We're creating this->segments in elems_to_segments(), at the same time we're giving them segment_midpoints. How are we losing the connection between them afterward? Is this an insert_any_extra_boundary_points() case? Are we calling a prepare_for_use() somewhere I don't see that's renumbering nodes in the middle? We use segments repeatedly in insert_refinement_points() too, so if something with it is getting screwed up on iteration N then I don't want to just go on to N+1 because we've worked around the use case in increase_triangle_order().
Could you come up with a unit test? Aside from making sure we don't regress on this, that should be enough for me to go through in the debugger and find where and why we have node ids changing during triangulation.
Sure, I'll aim to work on a unit test later this week, once I tackle some pressing deadlines.
Added a unit test. It should work with the this PR's fix but will fail if run with the current devel
Not sure why some civet tests failed (not familiar with libmesh's testing system).
Your CPPUNIT_TEST() invocation is in the #ifdef LIBMESH_HAVE_POLY2TRI block, but you defined it in the middle of the #ifdef LIBMESH_HAVE_TRIANGLE block. For builds that have Poly2Tri enabled but not Triangle (which is common and includes all MOOSE builds - Triangle is GPL-only), the unit tests want to include your new test but don't see a definition for it.
The ideal thing to do for tests that can be supported by both Poly2Tri and Triangle backends is to create a single implementation that works with a reference-to-base-class and two instantiations that each pass in one of the different subclasses. There's plenty of examples in there to look at.
I added the unit test as directed. thanks.
@roystgnr just a friendly reminder, thanks.
Hi @roystgnr , can you take a look at this again. The unit test has been added so that you can check when the node ids are changed. Thanks.