Contact surface needs to carry triangle face normals

[sherm1]

Per Drake hydrostatics notes document and many f2f discussions, dissipation depends on face normals (unit vectors). These have to be consistently oriented so that all the normals on one side point "outward" and all the normals on the other side point "inward". Note that this is easy for rigid/soft contact since the contact surface inherits outward normals from the rigid object. It is much harder for soft/soft but we don't need to resolve that yet.

Currently the contact surface carries triangle area but not triangle normals (though they do get calculated during the area computation). The "outward" normal should be stored per-triangle. Also, there is no need for a gradient field to be stored at the contact surface vertices. That should be removed. For soft/soft we will likely need one more scalar field (local stiffness) but still no vector field.

To be clear: I am only talking about the contact surface. The tet meshes still need to carry separately-sampled fields for extent ε(𝒙) and extent gradient ∇ε(𝒙).

[SeanCurtis-TRI]

We'll discuss this face to face, but I wanted to officially register my position on this.

I feel that the statement:

The "outward" normal should be stored per-triangle.

risks being limiting and not generally correct, but only if that concept bleeds through to the API.

The intent I infer is to guarantee the following points:

• In soft-rigid contact, the contact surface (if one exists) is always a subset of the rigid surface.
• The normal on every point on the contact surface should be the same as the normal on the rigid surface at the coincident point.

I wholly embrace those two points. If that's what is being addressed by the original post, then I fully agree that in some cases, the contact surface's normals are constant over the domain of a triangle and the ContactSurface should report it as such. But that is not generally true (obviously not for soft-soft contact and potentially not for rigid-soft depending on how we implement it -- think rigid sphere vs soft volume mesh).

Definitely any consumer of of ContactSurface should not be concerned with how the normals are defined and should simply ask for the normal at point.

[sherm1]

My point is stronger than that. Clarified with Ryan this morning: when we intersect two tets and solve for the contact plane, that plane has a normal n. We then clip that plane and triangulate it. The plane normal n is the face normal for those triangles. We need that normal in order to properly calculate tractions. That normal cannot be reliably recalculated later and needs to be stored with the contact surface triangles.

That is the basis for the hydrostatically-motivated physical model we are using. It is conceivable that there is some other model that would involve a normal field that has to be interpolated across triangular faces. We do not currently have that model. Right now we need those face normals.

[edrumwri]

The one question I had yesterday and which I think Sean broached this morning:

We know that we'll need a tetrahedralized sphere for the hydroelastic model performing, e.g., sphere/half-space intersection. The contact surface output will then be polyhedral. But the hydrostatic model would use the normal to the sphere, not the normal to the polyhedral approximation to the sphere, right?

[sherm1]

It's true that a smooth ball in water gets pressure normal to the surface everywhere. But the input to our hydroelastic method is not a smooth ball -- it is a crude (by design) collection of tets. We use those to find planes from which we can construct a contact surface. Our publications go to great lengths to tout (a) the crudeness of the allowed tesselations, and (b) the linearity of the plane-finding method. As a result we can't make arguments depending on limits of fine tesselation as we could if we were doing FEM. A different (and maybe better) model could be designed that worked with approximations of smooth surfaces, but our current hydroelastic model is not that -- it is a model whose derivation assumes discrete shapes.

Perhaps the "locus of pressure-equilibrium points" description of the contact surface has caused confusion. The correct description would be more like "set of triangles whose vertices are pressure-equilibrium points".

[SeanCurtis-TRI]

I believe I had an epiphany this morning that finally allows to reconcile what has been our otherwise irreconcilable differences.

The continuous math

It starts with the functions ε(x⃗) (the extent function) and p₀(ε(x⃗)) (the static pressure function).

• ε(x⃗)
  • Defined in a volume in ℝ³ between a "soft" boundary and a "rigid" boundary.
  • The boundaries don't need to be smooth, but ε and ∇ε are smooth between the boundaries.
• p₀(ε(x⃗))
  • scalar function mapping extent to static pressure.
  • It must be a non-decreasing function.
  • p₀ and ∇p₀ are likewise smooth between the boundaries.
• hᴬᴮ = p₀ᴬ - p₀ᴮ (and it and its gradient are smooth between the boundaries).
• The contact surface is where hᴬᴮ = 0. It is a smooth surface and ∇hᴬᴮ is normal to the surface at all points on the surface.

Discretizing the math

For lack of better notation, I'm using L(⋅) to mean a linear approximation of the given function. (Is there standard notation for that?)

It's how we discretize the problem that I believe we've had an unrecognized difference which gives rise to different conclusions.

In my mental model we were doing the following:

• Define two fields: L(p₀) and L(∇p₀).
• Define hᴬᴮ = L(p₀ᴬ) - L(p₀ᴮ) and ∇hᴬᴮ = L(∇p₀ᴬ) - L(∇p₀ᴮ)
• As a result, ∇hᴬᴮ is piecewise linear and, in some sense, is supposed to be a piecewise linear approximation of the continuous ∇hᴬᴮ. Therefore, the normals of the contact surface are piecewise linear.

If, however, I shift my mental model to this, I can draw your conclusions:

• Define one field: L(p₀)
• hᴬᴮ = L(p₀ᴬ) - L(p₀ᴮ) is still the same, but ∇hᴬᴮ = ∇L(p₀ᴬ) - ∇L(p₀ᴮ)
• In other words, I don't independently sample ∇p₀. I literally take the gradient of L(p₀). Because L(p₀) is piecewise linear ∇L(p₀) is piecewise constant and, therefore, hᴬᴮ is also piecewise constant. Therefore, the triangle normals (vis a vis plane normals) are the normals to the contact surface.

Currently, our code doesn't reflect this. We are storing L(∇p₀) instead of ∇L(p₀).

[sherm1]

Thanks -- nicely articulated and I like the clear notation (I don't know of a standard notation). The same issue came up in yesterday's meeting with @ryanelandt -- the in-progress paper is currently fuzzy about the distinction you drew. It needs to be clear about L(∇p) vs. ∇L(p) and you are right that Ryan's assumptions about the contact surface are based on ∇L(p), i.e. that the pressure field gradient is constant over each triangle.

BTW I'm not skeptical about the possibility of improving hydroelastics with smarter normals in cases where we have reliable information about curvature. But developing that theory is a "future work" topic we should not try to squeeze in to the current implementation.

[DamrongGuoy]

Fixed by #12101.