Skip to content

Mesh-mesh GJK/EPA contact is not invariant to common world translation #1546

Description

@mmichelis

Follow-up to #1527 after #1530. For the five-vertex convex mesh in the repro below, MuJoCo Warp creates a valid constraint at y=0 but misses it after both bodies are translated to y=1. A common translation should not change the result.

Tested with MuJoCo Warp 3.10.0.2 at 081c74e276cfc7f2484e98b80b1377aa4b3c64bd, MuJoCo 3.10.0, Warp 1.16.0.dev20260716, and CUDA 12.9 on an RTX 5090.

y=0.0: native_contacts=3, native_distance=-0.0816492254
y=0.0: mjwarp_contacts=1, constraint=True,  distance=-0.08164963871, addresses=[[0, 1, 2, -1]]
y=1.0: native_contacts=3, native_distance=-0.0816492254
y=1.0: mjwarp_contacts=1, constraint=False, distance=1.19209238e-08, addresses=[[-1, -1, -1, -1]]

Root cause

Both meshes receive the same translation, so their relative transform is unchanged: (A + t) - (B + t) = A - B. MJWarp instead transforms each support point to absolute float32 world coordinates before subtracting the pair, effectively computing float32(A + t) - float32(B + t).

Near y=1, float32 spacing is about 1.19e-7, so relevant simplex components of a few 1e-9 round away. GJK consequently produces a three-point triangle instead of a four-point tetrahedron containing the origin. EPA cannot initialize a full-dimensional polytope and leaves the tiny positive GJK distance unchanged.

The zero contact margin requires a negative distance to create a constraint. The per-geom gap=0.01 values sum to 0.02, so MJWarp retains the positive-distance candidate as a raw contact but leaves its constraint addresses at -1. The gap exposes the failure; it does not cause it.

The five-vertex mesh is a minimized, near-degenerate trigger for the precision loss; none of its four-vertex subsets reproduces the same failure.

Proposed fix

Run rigid non-heightfield GJK, EPA, and multicontact in a pair-local frame, then restore the origin when writing the contact position:

origin = geom1.pos
geom1.pos = wp.vec3(0.0)
geom2.pos -= origin

# GJK, EPA, and multicontact

contact_pos = 0.5 * (witness1 + witness2) + origin

This cancels the common translation before support-point rounding. It restores constraints at both translations with distance -0.08164963871; only the final contact position needs the world origin restored.

Repro

import mujoco
import mujoco_warp as mjw
import numpy as np
import warp as wp


DEVICE = "cuda:0"
VERTICES = np.array(
  [
    [0.0, 0.0, -0.3],
    [0.2, 0.0, -0.2],
    [-0.2, 0.0, -0.2],
    [-0.1, -0.2, -0.2],
    [0.0, -0.2, -0.2],
  ],
  dtype=np.float32,
)


def make_model(y_offset):
  vertices = " ".join(str(value) for value in VERTICES.reshape(-1))
  xml = f"""
  <mujoco>
    <option cone="elliptic"/>
    <asset>
      <mesh name="a" vertex="{vertices}"/>
      <mesh name="b" vertex="{vertices}"/>
    </asset>
    <worldbody>
      <body pos="0 {y_offset} 0.5">
        <freejoint/>
        <geom type="mesh" mesh="a" gap="0.01"/>
      </body>
      <body pos="0.2 {y_offset} 0.5">
        <freejoint/>
        <geom type="mesh" mesh="b" gap="0.01"/>
      </body>
    </worldbody>
  </mujoco>
  """
  return mujoco.MjModel.from_xml_string(xml)


def run(y_offset):
  host_model = make_model(y_offset)
  native_data = mujoco.MjData(host_model)
  mujoco.mj_forward(host_model, native_data)
  native_distance = min(
    (native_data.contact[i].dist for i in range(native_data.ncon)), default=float("inf")
  )

  with wp.ScopedDevice(DEVICE):
    model = mjw.put_model(host_model)
    data = mjw.put_data(host_model, mujoco.MjData(host_model))
    mjw.forward(model, data)
    nacon = int(data.nacon.numpy()[0])
    distances = data.contact.dist.numpy().reshape(-1)[:nacon]
    addresses = data.contact.efc_address.numpy().reshape(-1, 4)[:nacon]

  distance = float(np.min(distances)) if nacon else float("inf")
  constrained = bool(np.any(addresses >= 0))
  print(f"y={y_offset}: native_contacts={native_data.ncon}, native_distance={native_distance:.10g}")
  print(
    f"y={y_offset}: mjwarp_contacts={nacon}, constraint={constrained}, "
    f"distance={distance:.10g}, addresses={addresses.tolist()}"
  )
  return native_data.ncon, native_distance, nacon, distance, constrained


def main():
  y0 = run(0.0)
  y1 = run(1.0)
  assert y0[0] > 0 and y0[1] < 0.0 and y0[2] > 0 and y0[3] < 0.0 and y0[4]
  assert y1[0] > 0 and y1[1] < 0.0 and y1[2] > 0 and y1[3] >= 0.0 and not y1[4]


if __name__ == "__main__":
  main()

The y=1 assertion documents the expected unpatched failure.

Metadata

Metadata

Assignees

Labels

No labels
No labels

Type

No type

Projects

No projects

Milestone

No milestone

Relationships

None yet

Development

No branches or pull requests

Issue actions