(Arbitrary-) Polygon elements #4074

roystgnr · 2025-02-10T18:59:30Z

This ended up being a larger diff than I'd hoped for - we need to change a whole swath of Quadrature and FE APIs, to take an Elem * or Elem & instead of just an ElemType, if we want to be able to handle elements where a single type might have differing subelement decompositions and/or different values for things like n_dofs from element to element.

With just Polygon subclasses, it might been conceivable to just add PENTAGON, HEXAGON, etc. enumerated types until we reached n_sides()=17 or some such too-crazy-to-exceed value, but I'm hoping that the refactoring here also lets us handle an upcoming Polyhedron class hierarchy.

Since I'm deprecating a bunch of functions, I also uncommented the libmesh_deprecated() on a bunch of others while I was at it, and wrapped everything I could in #ifdef LIBMESH_ENABLE_DEPRECATED to try and catch use cases at compile time.

Some InfFE code was using FE/FEInterface calls that were deprecated by @lindsayad as part of his work on p-refinement flexibility, and I'm not sure if I updated them correctly for p_level() > 0 cases. I don't think anybody has yet tried combining infinite elements with p refinement (@BalticPinguin?), but if somebody did and I just broke it then on the one hand I apologize but on the other I consider that acceptable punishment for not getting an example with feature coverage into CI.

moosebuild · 2025-02-12T00:22:56Z

Job Coverage, step Generate coverage on 8d05c27 wanted to post the following:

Coverage

	159f51	#4074 8d05c2
	Total	Total	+/-	New
Rate	63.22%	63.29%	+0.06%	75.23%
Hits	74214	74643	+429	498
Misses	43168	43301	+133	164

Diff coverage report

Full coverage report

Warnings

New new line coverage rate 75.23% is less than the suggested 90.0%

This comment will be updated on new commits.

moosebuild · 2025-03-04T02:56:58Z

Job Coverage, step Verify coverage on cf92f3f wanted to post the following:

The following coverage requirement(s) failed:

Failed to generate coverage rate (required: 55.0%)

roystgnr · 2025-03-04T13:31:13Z

Rebasing to fix a Makefile.in conflict, but the last CI only shows the Test MOOSE ARM mac MOOSE bug and a multiple-shas issue in Coverage, so this should be safe to merge finally.

Some subclasses like Polygon might have to do this for themselves, because they can't heap allocate until after initializing Elem.

Throwing in other method tests here while I'm at it, just to get a quick check before I try to instantiate a Polygon1 in our usual Elem tests

If we have a POLYGON1 or a future POLYHEDRON1 or anything else where the quadrature rule will need to depend on more than just a single element type, we'll need to initialize quadrature objects with the physical element.

Obviously we can't rely on static maps here in general anymore, since we're not going to have a different ElemType and map entry for every possible polyhedron. So, add some invalid_uint entries for the general types, and handle invalid_uint where we might run into it, and query the Elem itself instead of these maps whereever we can.

This lets us do a lot in Polygon that might otherwise have to be pushed down to subclasses, which will be good if we ever have any higher-order or internal-node-having subclasses.

I'll finish Lagrange+Polygon1 support soon, but this will be faster even then, and it'll get Polygon elem unit tests working earlier.

Not sure how I missed this earlier

Some of the tests here don't apply to a non-affine pentagon, or apply with different numbers to a distorted pentagon, or won't apply until we can do a Lagrange basis on a pentagon.

This should be useful in distinguishing Polygon/Polyhedron from more static types.

This should keep caching from breaking on polygon/polyhedra of mixed runtime topology

This should be more forwards compatible than hardcoding a POLYGON1 test.

Putting this in FEInterface was kind of an unintuitive decision in the first place, and now that we're adding runtime-dependent element topology like Polygons we're going to *need* a physical Elem to work with in those cases.

We make this one an equilateral pentagon because a non-skewed polygon won't give us an exact representation of linears for p>1 MONOMIAL

We can't build the old Side objects for Polygon1; the build code there relies on `side_nodes_map`.

This is a bit more descriptive of an important aspect (triangulated mapping), and we can still add terms to indicate higher order in later higher-order versions.

QComposite was hitting undefined behavior here, creating temporary subelements for an internal Lagrange FE, and then when that FE would reinit afterward it would test `get_elem_type()` and hit a dangling pointer.

Shouldn't have left it in the header after seeing how big it needed to be.

jwpeterson

My main comment is that we should prefer libmesh_error_msg() to libmesh_error() even for unlikely errors, since it's just better documentation even when we don't hit the errors.

jwpeterson · 2025-03-04T16:31:22Z

src/fe/fe_abstract.C

+  if (Elem::type_to_n_nodes_map[itemType] == invalid_uint)
+    libmesh_error();
+
  nodes.resize(Elem::type_to_n_nodes_map[itemType]);


It seems a little clunky to do the Elem::type_to_n_nodes_map[itemType] lookup twice (once for the error check, once for the resize) and I think you could further improve this diff by using a libmesh_error_msg() rather than giving no information if there's an error.

jwpeterson · 2025-03-04T16:33:50Z

src/fe/fe_lagrange.C

+          case C0POLYGON:
+            // C0Polygon requires using newer FE APIs
+            if (!e)
+              libmesh_error();


Prefer a libmesh_error_msg() with some details here as well.

jwpeterson · 2025-03-04T16:34:51Z

src/fe/fe_lagrange_shape_2D.C

+            {
+              // C0Polygon requires using newer FE APIs
+              if (!elem)
+                libmesh_error();


libmesh_error_msg()

jwpeterson · 2025-03-04T17:09:49Z

src/fe/fe_lagrange_shape_2D.C

+            {
+              // C0Polygon requires using newer FE APIs
+              if (!elem)
+                libmesh_error();


libmesh_error_msg()

jwpeterson · 2025-03-04T17:15:00Z

src/geom/elem.C

+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+    invalid_uint,  // unused
+


This is because of the enumeration value you chose for C0POLYGON = 50, yes? What is the reason for having the gap, exactly? It doesn't look like we've left gaps for other types of Elems in the past so my first instinct would have just been to keep the numbering sequential for the new type...

Well, this is embarrassing... honestly, the gaps are there because I have another branch that was going to use them, for higher-order nodal (Lagrange + Rational-Bezier-Bernstein mappings) geometric element support. That branch is probably less than half finished and it's turned out to be less of a priority than I'd expected, though, so I'll reshuffle this now, and whenever I get back to that I'll move its numbers up when I hit the rebase conflicts.

jwpeterson · 2025-03-04T17:17:17Z

src/geom/elem.C

+
+  const auto n_elem_ints = this->n_extra_integers();
+  returnval->add_extra_integers(n_elem_ints);
+  for (unsigned int i = 0; i != n_elem_ints; ++i)


Could be a make_range() loop but not a big preference for me.

jwpeterson · 2025-03-04T17:18:18Z

src/geom/face_c0polygon.C

+      swap2nodes(s,ns-1-s);
+      swap2neighbors(s,ns-2-s);
+      swap2boundarysides(s,ns-2-s,boundary_info);
+      swap2boundaryedges(s,ns-2-s,boundary_info);


Need spaces after commas, at least that is how we usually format comma-separated parameter lists?

Okay, but I'm not changing any of the 200ish spots in older files where I forgot the spaces. ;-)

jwpeterson · 2025-03-04T17:21:04Z

src/geom/face_polygon.C

+    case MIN_ANGLE:
+      bounds.first  = 90. - (180./this->n_sides());
+      bounds.second = 180. - (360./this->n_sides());
+      break;


I like that you added these bounds for completeness, but there is no Elem::quality() implementation yet for Polygons is there?

There's no specialization, but the base class implementation actually works fine for polygons, and we now test with it in elem_test.C - I had to add a special case for them because a squashed pentagon I construct for one instantiation has an angle that's even worse than our standard pyramids.

jwpeterson · 2025-03-04T17:21:53Z

src/quadrature/quadrature.C

+  // Some element types require data from a specific element, so can
+  // only be used with newer APIs.
+  if (t == C0POLYGON)
+    libmesh_error();


libmesh_error_msg()

The base class `Elem::volume()` call caught bad quadrature weights even when something (cancellation?) prevented them from being caught by our FE tests.

This avoids a bunch of unused values in between, and the only point of it was to reduce conflicts with another feature branch that may be delayed a while regardless.

This was referenced Feb 11, 2025

Add transient_ex3 example demonstrating explicit DG/FV formulation of 2D advection equation #4073

Merged

Quadrature API changes #4076

Merged

FE API deprecation/changes #4077

Merged

roystgnr force-pushed the polygon_elements branch from 48c594b to cf92f3f Compare March 3, 2025 20:16

roystgnr added 23 commits March 4, 2025 07:32

Enum entry for a polygon element type

80c4f0d

Add POLYGON1 string-to-enum

4582501

Make Elem::Elem() links initialization optional

ce29acd

Some subclasses like Polygon might have to do this for themselves, because they can't heap allocate until after initializing Elem.

Add Polygon1 Elem subclass

2d23848

Build new Polygon classes

cbfbcc9

Add some Polygon volume tests

ea39119

Throwing in other method tests here while I'm at it, just to get a quick check before I try to instantiate a Polygon1 in our usual Elem tests

Add Elem::disconnected_clone()

c821836

Use disconnected_clone() in create_submesh()

4480a1f

Use disconnected_clone() in true_centroid()

ca1054f

disconnected_clone() in copy_nodes_and_elements()

363f621

Detect POLYGON1 in older QBase::init() API

22a84c6

If we have a POLYGON1 or a future POLYHEDRON1 or anything else where the quadrature rule will need to depend on more than just a single element type, we'll need to initialize quadrature objects with the physical element.

QGauss support for Polygon quadrature

f0b91c5

Polygon::n_nodes_per_side()

041f486

This lets us do a lot in Polygon that might otherwise have to be pushed down to subclasses, which will be good if we ever have any higher-order or internal-node-having subclasses.

Overload Polygon1::true_centroid()

a0f9a5d

I'll finish Lagrange+Polygon1 support soon, but this will be faster even then, and it'll get Polygon elem unit tests working earlier.

Expect unsigned_uint in POLYGON1 maps

04d70c2

Don't try to refine Polygon1 in tests

d722189

One more QBase::init flag to add

aa9940e

Not sure how I missed this earlier

Assert we're within Polygon::master_point() limits

88bc446

Nodal quadrature for POLYGON1

35826df

Polygon1::has_affine_map() implementation

7a25010

Prepare ElemTest unit tests for a pentagon elem

467567d

Some of the tests here don't apply to a non-affine pentagon, or apply with different numbers to a distorted pentagon, or won't apply until we can do a Lagrange basis on a pentagon.

Instantiate POLYGON1 in more unit tests

327c1ee

roystgnr added 20 commits March 4, 2025 07:32

Add Elem::runtime_topology() bool

6834fcc

This should be useful in distinguishing Polygon/Polyhedron from more static types.

Consider elem, not just elem type, for FE caching

7942194

This should keep caching from breaking on polygon/polyhedra of mixed runtime topology

Use Elem::runtime_topology() where appropriate

42f064e

This should be more forwards compatible than hardcoding a POLYGON1 test.

Comment on future test expansion

0bf0af9

Generally useful Polygon methods

0ac5591

runtime-elem-topology-compatible n_shape_functions

ba2108b

Clarify an old comment

297330a

Added Elem::on_reference_element()

b9408aa

Putting this in FEInterface was kind of an unintuitive decision in the first place, and now that we're adding runtime-dependent element topology like Polygons we're going to *need* a physical Elem to work with in those cases.

Lagrange FE support for Polygon1

dcc51c4

Enable the rest of the Polygon1 unit tests

d368696

FEInterface::max_order() entries for Polygon1

8833b48

Add Monomial support to Polygon1

a3154fa

Add a Polygon1 setup option to FE unit tests

6768eab

We make this one an equilateral pentagon because a non-skewed polygon won't give us an exact representation of linears for p>1 MONOMIAL

Instantiate Polygon1 FE tests

44692f8

Fix Polygon1 with --enable-deprecated

377c460

We can't build the old Side objects for Polygon1; the build code there relies on `side_nodes_map`.

Clearer comment about Polygon1 master points

660e368

Remember to mark overriding functions as override

5c6253b

Polygon1 -> C0Polygon rename

460c29f

This is a bit more descriptive of an important aspect (triangulated mapping), and we can still add terms to indicate higher order in later higher-order versions.

Avoid dereferencing dangling FEAbstract::_elem

b3f478d

QComposite was hitting undefined behavior here, creating temporary subelements for an internal Lagrange FE, and then when that FE would reinit afterward it would test `get_elem_type()` and hit a dangling pointer.

Re-bootstrap

faed578

roystgnr force-pushed the polygon_elements branch from cf92f3f to faed578 Compare March 4, 2025 13:32

roystgnr added 3 commits March 4, 2025 09:14

Update copyright dates

acfa020

Fix doc comment

cfd6e10

Move Polygon constructor to .C class

e1cc433

Shouldn't have left it in the header after seeing how big it needed to be.

jwpeterson requested changes Mar 4, 2025

View reviewed changes

roystgnr added 5 commits March 4, 2025 14:28

Fix C0Polygon quadrature weights

97b3780

Refactor and increase Polygon volume testing

ab768ab

The base class `Elem::volume()` call caught bad quadrature weights even when something (cancellation?) prevented them from being caught by our FE tests.

libmesh_error() -> libmesh_error_msg()

e726013

Spaces after commas

2764ed2

Bump up C0POLYGON enum value

8d05c27

This avoids a bunch of unused values in between, and the only point of it was to reduce conflicts with another feature branch that may be delayed a while regardless.

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

(Arbitrary-) Polygon elements #4074

(Arbitrary-) Polygon elements #4074

roystgnr commented Feb 10, 2025

moosebuild commented Feb 12, 2025 •

edited

Loading

moosebuild commented Mar 4, 2025

roystgnr commented Mar 4, 2025

jwpeterson left a comment

jwpeterson Mar 4, 2025

jwpeterson Mar 4, 2025

jwpeterson Mar 4, 2025

jwpeterson Mar 4, 2025

jwpeterson Mar 4, 2025

roystgnr Mar 4, 2025

jwpeterson Mar 4, 2025

jwpeterson Mar 4, 2025 •

edited

Loading

roystgnr Mar 4, 2025

jwpeterson Mar 4, 2025

roystgnr Mar 4, 2025

jwpeterson Mar 4, 2025

(Arbitrary-) Polygon elements #4074

Are you sure you want to change the base?

(Arbitrary-) Polygon elements #4074

Conversation

roystgnr commented Feb 10, 2025

moosebuild commented Feb 12, 2025 • edited Loading

Coverage

Warnings

moosebuild commented Mar 4, 2025

roystgnr commented Mar 4, 2025

jwpeterson left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jwpeterson Mar 4, 2025 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

moosebuild commented Feb 12, 2025 •

edited

Loading

jwpeterson Mar 4, 2025 •

edited

Loading