Skip to content

Common plane higher order quadrature#208

Draft
ebchin wants to merge 7 commits into
developfrom
feature/ebchin/common-plane-quad-rules
Draft

Common plane higher order quadrature#208
ebchin wants to merge 7 commits into
developfrom
feature/ebchin/common-plane-quad-rules

Conversation

@ebchin

@ebchin ebchin commented May 28, 2026

Copy link
Copy Markdown
Member

No description provided.

@srwopschall srwopschall left a comment

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ebchin - a lot of comments, but only a handful of things that really need to be ironed out. Please review. In the meantime I need to look at the entirety of CommonPlane.cpp and verify the code and when/where average face coordinates are used vs. current configuration and possibly warped face coordinates. I want to be sure we don't have any inconsistencies.

{
constexpr int numVerts = 4;

RealT x1[numVerts] = { 0.0, 1.0, 1.0, 0.0 };

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Minor: considering adding a comment about the configuration of the two blocks including the gap.

tribol::ExecutionMode::Sequential );

tribol::setPenaltyOptions( 0, tribol::KINEMATIC, tribol::KINEMATIC_CONSTANT );
tribol::setCommonPlaneIntegrationOptions( 0, rule, 3 );

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Consider clarifying what 0 and 3 are since someone may come to the tests to see how this is used.

return result;
}

struct EdgeLocalContactResult {

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this struct intended to hold a force result. Consider consistent naming similar to the WarpedQuadForceResult.

{
constexpr int numVerts = 2;

RealT x1[numVerts] = { 1.0, 0.0 };

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Minor: consider adding comment about configuration of the two edges. Is this an interpen interaction or full overlap?

tribol::ExecutionMode::Sequential );

tribol::setPenaltyOptions( 0, tribol::KINEMATIC, tribol::KINEMATIC_CONSTANT );
tribol::setCommonPlaneIntegrationOptions( 0, rule, 4 );

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Clarify what 0 and 4 are.

const RealT length = magnitude( x1 - x0, y1 - y0 );

for ( int qp = 0; qp < num_qpts; ++qp ) {
const RealT s = segment_rule_coords[qp];

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Minor: consider for clarity adding comment about what this mapping is.

const bool mapped_face2 =
EvalLinearEdgeAtProjectedPoint( actual_xf2, x_q, overlapNormal, x_qf2, phi_q2, use_rate ? dim : 0,
use_rate ? &actual_vf2[0] : nullptr, use_rate ? vel_q2 : nullptr );
if ( !mapped_face1 || !mapped_face2 ) {

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same comment. Add debug statement here? Is this bad if we get in this if-block?


has_contact = true;

RealT local_pressure = local_gap * penalty_stiff_per_area;

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can you clarify whether this routine before computed the pressure? As a sanity check is this line inside ComputeGapRatePressure? Shouldn't this only compute the pressure as done on line 675? Or did the original version of this routine compute it all?

// contact overlap patch //
////////////////////////////////////////////////////////////////////////
EvalWeakFormIntegral<COMMON_PLANE, SINGLE_POINT>( cntctElem, phi1, phi2 );
EvalWeakFormIntegralCommonPlane( cntctElem, pen_enfrc_options.common_plane_rule,

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think here we need to clarify which set of face coordinates cntctElem has and how they are used in this routine and whether that conflicts with multipoint integration using current configuration physical coords.

// Mortar retains the historical triangle rule for now.
GaussPolyIntTri( elem, integ, 3, TRI_RULE_LEGACY );

// call Taylor-Wingate-Bos integation rule. NOTE: this is not

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Since we are touching integration rules, the TWB integration rule never worked and was experimental at the time. It is not consistent with the underlying finite element basis functions. Let's consider removing.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants