Attempting to model a bonding rigid-flexible interface - help?

[NonchalantWombat]

I have a sim where I am compressing a rubber torus (MPM elastic entity) under high load from a rigid cylinder, and it is deforming in a non-physical way, which I believe is due to both the collision/penetration of the rubber entity into the cylinder from the large forces and the fact that the rubber is not "bonded" to the surface of the cylinder. How do I achieve this sort of a "chemical bond," where on circular face of the torus is bonded to the cylinder and the rest of it is allowed to compress and bulge outward while constrained by this contact? Also, how do I prevent excessive penetration due to the large weight? I would like it to be as close to real as possible. See the attached image for what I mean, where the green layer should be bonded between the cylinder and torus.

Code:

import numpy as np
import pandas as pd
import os

# List of the geometry files to simulate
stl_files = ["Bubbles.stl", "Snake.stl", "Spine.stl", "Wavy.stl", "Bulk.stl"]

gs.init(backend=gs.gpu)

for stl_file in stl_files:
    print(f"--- Starting simulation for: {stl_file} ---")

    if not os.path.exists(stl_file):
        print(f"Warning: {stl_file} not found. Skipping.")
        continue

    scene = gs.Scene(show_viewer=True,
        sim_options=gs.options.SimOptions(
            dt       = 1.0e-3,
            substeps = 10,
        ),
        mpm_options=gs.options.MPMOptions(
            lower_bound   = (-0.5, -0.5, -0.1),
            upper_bound   = (0.5, 0.5, 1.0),
            particle_size = 0.0005,
            #enable_CPIC = True,
        ),
        #pbd_options=gs.options.PBDOptions(
        #    lower_bound   = (-0.5, -0.5, -0.1),
        #    upper_bound   = (0.5, 0.5, 1.0),
        #    particle_size = 0.0005,
        #    #enable_CPIC = True,
        #),
        viewer_options=gs.options.ViewerOptions(
            camera_pos=(0.3, 0.3, 0.05),
            camera_lookat=(0.0, 0.0, 0.0),
            camera_fov=50,
            max_FPS=120
        ),
        vis_options = gs.options.VisOptions(
            ambient_light    = (0.25, 0.25, 0.25),
        ),
    )
    # rigid ground plane
    plane = scene.add_entity(
        morph=gs.morphs.Plane(),
    )
    
    # Neoprene footpad with the current STL file
    pancake = scene.add_entity(
        material=gs.materials.MPM.Elastic(
            E=3430000.0,
            nu=0.465,
            rho=1230.0,
            model="neohooken"
        ),
        #material=gs.materials.PBD.Elastic(
        #    rho=1230.0,
        #    kinetic_friction=1.0,
        #    static_friction=3.0,
        #    volume_compliance=1.0,
        #    bending_compliance=1.0
        #),
        morph=gs.morphs.Mesh(
            file=stl_file,
            pos=(0.0, 0.0, 0.0),
            scale=0.001,
        ),    
        surface=gs.surfaces.Default(
            color=(1.0, 0.4, 0.4),
            vis_mode="particle",
        )
    )
    
    # rigid heavy plate pressing down
    plate = scene.add_entity(
        material=gs.materials.Rigid(
            #rho = 1500,
            rho = 115000,
            #coup_restitution = 0.02,
            coup_friction = 3.0
        ),
        morph=gs.morphs.Cylinder(
            pos  = (0.0, 0.0, 0.05),
            radius = 0.075,
            height = 0.05
        ),
        surface=gs.surfaces.Default(
            color    = (0.3, 0.3, 0.3),
            opacity   = 0.25,
            vis_mode = 'visual',
        )
    )
    scene.build()
    horizon = 500  # Run for 0.5 seconds (200 steps * 0.001s/step)
    dt = scene.sim_options.dt
    times = []
    thicknesses = []
    
    for i in range(horizon):
        scene.step()
        pts = pancake.get_particles()
        if pts.shape[0] > 0 and pts.shape[1] > 0:
            pts = pts[0]
            z_min = np.min(pts[:, 2])
            z_max = np.max(pts[:, 2])
            thickness = (z_max - z_min) * 1000  # Convert to mm
            times.append(i * dt)
            thicknesses.append(thickness)
    
    df = pd.DataFrame({
        'Time [s]': times,
        'Thickness [mm]': thicknesses
    })
    
    csv_filename = f"results_{stl_file.split('.')[0]}.csv"
    df.to_csv(csv_filename, index=False)
    print(f"--- Simulation for {stl_file} complete. Data saved to {csv_filename} ---\n")

print("All simulations finished.")```

[jgillick]

I'm guessing it'll be easier for them to help if you can submit the basic code of what you're doing today so it can be reproduced.

[NonchalantWombat]

Added into my original post, thanks