Inelastic constraints very sensive to mesh sizing/spacing on curved surface

I am simulating stress/strain on an assembly that contains 8 solid bodies and 12 contact joints. I have implemented the contacts using dolfinx_mpc and it sometimes works end to end.

The whole assembly:

An example contact joint:

The mesh:

The solution:

Where it starts to go wrong is that if I regenerate my mesh to be slightly more fine or more coarse, some of the mpc constraints fail to add. They raise:

mpc_data = dolfinx_mpc.cpp.mpc.create_contact_inelastic_condition(
RuntimeError: Newton method failed to converge for non-affine geometry

The mesh shown above succeeds and it has a mesh size set to 5% of the length of the overall model. If I regenerate the mesh size set to 4.5% instead, I get this mesh:

Which, when I run the code, succeeds at adding 7 of the contacts, then fails on the 8th with:

Traceback (most recent call last):
  File "/Users/matthewferraro/Documents/git/fenicsx_demo/backend/test_solver.py", line 479, in <module>
    test_dino(0.45, 0.01)
  File "/Users/matthewferraro/Documents/git/fenicsx_demo/backend/test_solver.py", line 366, in test_dino
    domain, uh, stresses = solver.solve_static(
  File "/Users/matthewferraro/Documents/git/fenicsx_demo/backend/solver.py", line 188, in solve_static
    mpc.create_contact_inelastic_condition(
  File "/Users/matthewferraro/miniconda3/envs/fenicsx-demo-backend/lib/python3.10/site-packages/dolfinx_mpc/multipointconstraint.py", line 461, in create_contact_inelastic_condition
    mpc_data = dolfinx_mpc.cpp.mpc.create_contact_inelastic_condition(
RuntimeError: Newton method failed to converge for non-affine geometry

I can trigger the same failure by regenerating the mesh with a cell size of 5.5%:

In fact I can generate many meshes at many different cell sizes and there seems to be no pattern between what works and what doesn’t work. At each cell size the result is repeatable.

What causes the exception to be raised? What is happening inside dolfinx_mpc that involves Newton’s method at the time the constraints are being created? I have watched this talk but when I look at the c++ code I get a bit lost.

If a minimum viable example is required I can create one, but I would also be happy with just a conceptual explanation.

What happens is that when you create the inelastic contact conditions, the “coordinates” of the degrees of freedom that are constrained, are mapped to the degrees of freedom on the other surface. In this process, there is a pull-back routine, which takes a physical coordinate, and pulls it back to the reference element of the other cell. For non-affine geometries, this is done through a Newton iteration.

Quite recently (december 2024), I added some tolerances that you can adjust for both collision detection and pull-backs: Move tolerances out of hardcoded part of code by jorgensd · Pull Request #140 · jorgensd/dolfinx_mpc · GitHub
which is used in the construction:
dolfinx_mpc/cpp/utils.h at 817cd4c385e1c61761e22c6fedc0c6b31094cd2a · jorgensd/dolfinx_mpc · GitHub
Could you try updating to this version of DOLFINx_MPC, and adjust the tolerance?