Repair the dual basis on QUAD8 and TRI6 faces - #4536
Conversation
The weight \int N_k used to build the biorthogonal dual basis is not positive on the serendipity quadratic faces. It is exactly 0 at a TRI6 vertex, which leaves those vertex dual shape functions identically zero, and -1/3 at a QUAD8 corner, which flips the sign of any per-node normalization built from it. Neither case is currently diagnosed. Biorthogonalize instead against the locally quadratic transformed basis of Popp et al., SIAM J. Sci. Comput. 34(4):B421-B446, 2012, Sec. 4.4.1, on exactly those two element types. Each vertex absorbs a fraction alpha of its adjacent mid-edge shapes and each mid-edge shape is scaled by 1 - 2*alpha; making the dual biorthogonal to that basis replaces D by T^-1 diag(T d), which for this T reduces to a sparse update of D needing no extra linear solve. With alpha = 1/5 the transformed weights are strictly positive (QUAD8 -> 1/5 and 4/5, TRI6 -> 1/15 and 1/10) and the partition of unity is preserved, so the dual basis still reproduces constants. The transform is the identity on every other element type, and the alpha = 0 limit reproduces diag(d) exactly, so all other dual bases are unchanged. dual_shape_verification_test gains QUAD8 and TRI6 cases covering positivity, reproduction of constants, and biorthogonality against the transformed basis. Refs idaholab/moose#33558 Co-authored-by: Claude <noreply@anthropic.com>
|
I've started up CI. Walk me through the math on this, and make sure I'm not misunderstanding? The standard dual basis on these elements can't be made to satisfy both biorthogonality and the partition-of-unity property at once without going degenerate, so we do a slight change-of-basis on the primal bases and then we find a dual to that changed basis instead? How does that work at the application level, when we no longer have a dual to the original basis but we aren't exposing the altered basis? I think I'm okay with this, under the circumstances (though it pains me to see a non-serendipity element behaving as badly and requiring as much of a hack as a serendipity!), but we've got to document what's going on in the headers, not just the implementation. Somehow I was an idiot and approved the original addition of |
compute_dual_shape_coeffsbuilds the dual basis fromD = diag(\int N_k), but that weight isexactly 0 at a TRI6 vertex -- leaving those vertex dual shape functions identically zero -- and -1/3
at a QUAD8 corner. Neither is diagnosed.
This biorthogonalizes against the locally quadratic transformed basis of Popp et al. (2012, Sec.
4.4.1) on those two element types, which reduces to a sparse update of
Dwith no extra solve. Thetransform is the identity everywhere else and
alpha = 0reproducesdiag(d)exactly, so nothingelse changes, there is nothing to regold, and no API changes.
New QUAD8/TRI6 cases in
dual_shape_verification_testfail without the fix. Full unit suite passes;MOOSE first-order dual mortar output is bit-identical before and after.
Refs idaholab/moose#33558