Surface Pressure Constraint causes unrealistic sliding

I have been evolving the geometries and layouts of simple soft robots which are made up of soft parts, rigid parts and actuating parts. I have successfully managed to do this using Sofa and SofaPython3 where I use and control the pressure value in a Surface Pressure Forcefield component to actuate the actuating parts. Unfortunately, I was not getting the desired result. As the pressures used would be sufficient when the actuating part was large but would barely do anything when the actuating part was small. I then attempted to use the Surface pressure constraint from the Soft Robots plugin instead which seemed to solve the issue of smaller volumes not actuating as much as larger volumes, and also allowed me to specifically actuate the volume growth iteself. Unfortunately, it created a new issue – when the soft robot’s actuator inflates the soft robot begins to slide along the surface, despite friction contact being implemented. It appears that when the actuator inflates, a reaction force is created causing the robot to slide.

1. Please could someone explain how both the Surface Pressure forcefield and the surface pressure constraint works and why they don’t seem to produce similar results (for example a pressure value of 150 for the forcefield component seems to equate to a pressure value of 1.5 for the contraint component).

2. How can I prevent the unrealistic sliding caused by using the surface constraint component?

The code illustrating an example is provided below:

import sys, os
import math

sys.path.append(os.path.abspath("./bindings/Sofa/package"))
sys.path.append(os.path.abspath("./bindings/SofaRuntime/package"))
sys.path.append(os.path.abspath("./bindings/SofaTypes/package"))

import Sofa
# Choose in your script to activate or not the GUI
USE_GUI = True

class FingerController(Sofa.Core.Controller):
    """ This is a custom controller to perform actions when events are triggered """

    def __init__(self, *args, **kwargs):
        # These are needed (and the normal way to override from a python class)
        Sofa.Core.Controller.__init__(self, *args, **kwargs)
        self.node = kwargs["node"]
        self.rootNode = kwargs["rootNode"]
        self.act1_size = kwargs["act1_size"]
        return

    def onAnimateBeginEvent(self, event):
        if self.rootNode.time.value > 1.5:
            self.node.Surface.surfaceConstraint.value.value = [0.6*math.sin(self.rootNode.time.value-1.5-math.pi/2)+0.6]

class FingerController2(Sofa.Core.Controller):
    """ This is a custom controller to perform actions when events are triggered """

    def __init__(self, *args, **kwargs):
        # These are needed (and the normal way to override from a python class)
        Sofa.Core.Controller.__init__(self, *args, **kwargs)
        self.node = kwargs["node"]
        self.rootNode = kwargs["rootNode"]
        self.act2_size = kwargs["act2_size"]
        return

    def onAnimateBeginEvent(self, event):
        if self.rootNode.time.value > 1.5+math.pi:
            self.node.Surface.surfaceConstraint.value.value = [-0.6*math.sin(self.rootNode.time.value-1.5-math.pi/2)+0.6]

def createScene(root):
    # Set time step and gravity of the Scene
    root.gravity=[0, -9.81, 0]
    root.dt=0.01
    root.addObject('RequiredPlugin', pluginName="SoftRobots SofaGeneralObjectInteraction SofaConstraint SofaGeneralSimpleFem SofaOpenglVisual SofaGeneralLoader SofaImplicitOdeSolver SofaLoader SofaMiscCollision SofaSimpleFem")
    # Scene must now include a VisualLoop
    root.addObject('DefaultVisualManagerLoop')
    # Scene must now include a AnimationLoop
    root.addObject('FreeMotionAnimationLoop')
    root.addObject('GenericConstraintSolver', tolerance=1e-3, maxIt=1000)

    root.addObject('DefaultPipeline', verbose=0)
    root.addObject('BruteForceDetection', name="N2")
    root.addObject('DefaultContactManager', name="Response", response="FrictionContact", responseParams="mu=1")
    root.addObject('NewProximityIntersection', name="Proximity", alarmDistance=0.8, contactDistance=0.5)
    root.addObject('DefaultCollisionGroupManager', name="Group")

    floor = root.addChild('Floor')
    floor.addObject('MeshGmshLoader', name="floorLoader", filename="/home/mike/SOFA_v20.12.01_Linux/share/sofa/mesh/test_meshes/floor.msh")
    floor.addObject('QuadSetTopologyContainer', name="quadContainer", src="@floorLoader")
    floor.addObject('TriangleSetTopologyContainer', name='triangleContainer')
    floor.addObject('TriangleSetTopologyModifier', name='triangleModifier')
    floor.addObject('Quad2TriangleTopologicalMapping', input="@quadContainer", output="@triangleContainer")
    floor.addObject('MechanicalObject', name="FloorDOF", src="@triangleContainer", scale=50, translation=[-50.0, 0.0, -50.0], rotation=[0.0, 0.0, 0.0])
    floor.addObject('TriangleCollisionModel', src="@triangleContainer", simulated=0, moving=0, bothSide=1, group=1)
    floor.addObject('LineCollisionModel', src="@triangleContainer", simulated=0, moving=0, group=1)
    floor.addObject('PointCollisionModel', src="@triangleContainer", simulated=0, moving=0, group=1)
    floor.addObject('OglModel', name="FloorVis", color="white")

    soft = root.addChild('soft')
    soft.addObject('EulerImplicitSolver', rayleighStiffness=0)
    soft.addObject('CGLinearSolver', iterations=200, tolerance=1e-06, threshold=1e-06)
    soft.addObject('MeshGmshLoader', name="meshLoader", filename="/home/mike/meshes/FourPart/3*3*3/5softpart.msh")
    soft.addObject('HexahedronSetTopologyContainer', name="topo", src="@meshLoader")
    soft.addObject('MechanicalObject', name="softDOF", template="Vec3d", translation=[0.0, 1.0, 0.0])
    soft.addObject('UniformMass', vertexMass=0.3)
    soft.addObject('PrecomputedConstraintCorrection')
    soft.addObject('HexahedralFEMForceField', name="FEM", youngModulus=500, poissonRatio=0.3, method="large")
    surface = soft.addChild('Surface')
    surface.addObject('QuadSetTopologyContainer', name="QuadContainer")
    surface.addObject('QuadSetTopologyModifier', name="QuadModifier")
    surface.addObject('Hexa2QuadTopologicalMapping', input="@../topo", output="@QuadContainer")
    collision = surface.addChild('Collsion')
    collision.addObject('TriangleSetTopologyContainer', name="TriContainer")
    collision.addObject('TriangleSetTopologyModifier', name="TriModifier")
    collision.addObject('Quad2TriangleTopologicalMapping', input="@../QuadContainer", output="@TriContainer")
    collision.addObject('TriangleCollisionModel', name="CollisionMO", template="Vec3d", group=2)
    collision.addObject('LineCollisionModel', group=2)
    collision.addObject('PointCollisionModel', group=2)
    collision.addObject('OglModel', name="softVis", color="red")
    collision.addObject('BarycentricMapping', input="@../../softDOF", output="@softVis")

    actuator1 = root.addChild('actuator1')
    actuator1.addObject('EulerImplicitSolver', rayleighStiffness=0)
    actuator1.addObject('CGLinearSolver', iterations=200, tolerance=1e-06, threshold=1e-06)
    actuator1.addObject('MeshGmshLoader', name="meshLoader", filename="/home/mike/meshes/FourPart/3*3*3/5actuator1part.msh")
    actuator1.addObject('HexahedronSetTopologyContainer', name="topo", src="@meshLoader")
    actuator1.addObject('MechanicalObject', name="actuator1DOF", translation=[0.0, 1.0, 0.0])
    actuator1.addObject('UniformMass', vertexMass=0.3)
    actuator1.addObject('PrecomputedConstraintCorrection')
    surface = actuator1.addChild('Surface')
    surface.addObject('QuadSetTopologyContainer', name="QuadContainer")
    surface.addObject('QuadSetTopologyModifier', name="QuadModifier")
    surface.addObject('Hexa2QuadTopologicalMapping', input="@../topo", output="@QuadContainer")
    surface.addObject('MeshSpringForceField', name="Springs", stiffness=200, damping=1)
    # surface.addObject('SurfacePressureForceField', name="Pressure1", pressure=0, pulseMode=0)
    surface.addObject('SurfacePressureConstraint', name='surfaceConstraint', template='Vec3d', value=0, valueType='pressure')
    collision = surface.addChild('Collision')
    collision.addObject('TriangleSetTopologyContainer', name="TriContainer")
    collision.addObject('TriangleSetTopologyModifier', name="TriModifier")
    collision.addObject('Quad2TriangleTopologicalMapping', input="@../QuadContainer", output="@TriContainer")
    collision.addObject('TriangleCollisionModel', name="CollisionMO", template="Vec3d", group=2)
    collision.addObject('LineCollisionModel', group=2)
    collision.addObject('PointCollisionModel', group=2)
    collision.addObject('OglModel', name="actuator1Vis", color="green")
    collision.addObject('BarycentricMapping', input="@../../actuator1DOF", output="@actuator1Vis")

    actuator2 = root.addChild('actuator2')
    actuator2.addObject('EulerImplicitSolver', rayleighStiffness=0)
    actuator2.addObject('CGLinearSolver', iterations=200, tolerance=1e-06, threshold=1e-06)
    actuator2.addObject('MeshGmshLoader', name="meshLoader", filename="/home/mike/meshes/FourPart/3*3*3/5actuator2part.msh")
    actuator2.addObject('HexahedronSetTopologyContainer', name="topo", src="@meshLoader")
    actuator2.addObject('MechanicalObject', name="actuator2DOF", translation=[0.0, 1.0, 0.0])
    actuator2.addObject('UniformMass', vertexMass=0.3)
    actuator2.addObject('PrecomputedConstraintCorrection')
    surface = actuator2.addChild('Surface')
    surface.addObject('QuadSetTopologyContainer', name="QuadContainer")
    surface.addObject('QuadSetTopologyModifier', name="QuadModifier")
    surface.addObject('Hexa2QuadTopologicalMapping', input="@../topo", output="@QuadContainer")
    surface.addObject('MeshSpringForceField', name="Springs", stiffness=200, damping=1)
    # surface.addObject('SurfacePressureForceField', name="Pressure2", pressure=2, pulseMode=0)
    surface.addObject('SurfacePressureConstraint', name='surfaceConstraint', template='Vec3d', value=0, valueType='pressure')
    collision = surface.addChild('Collision')
    collision.addObject('TriangleSetTopologyContainer', name="TriContainer")
    collision.addObject('TriangleSetTopologyModifier', name="TriModifier")
    collision.addObject('Quad2TriangleTopologicalMapping', input="@../QuadContainer", output="@TriContainer")
    collision.addObject('TriangleCollisionModel', name="CollisionMO", template="Vec3d", group=2)
    collision.addObject('LineCollisionModel', group=2)
    collision.addObject('PointCollisionModel', group=2)
    collision.addObject('OglModel', name="actuator2Vis", color="yellow")
    collision.addObject('BarycentricMapping', input="@../../actuator2DOF", output="@actuator2Vis")

    root.addObject('AttachConstraint', object1='@soft', object2='@actuator1', twoWay=True, indices1=[0, 1, 2, 3, 4, 5], indices2=[6, 7, 8, 9, 10, 11])
    root.addObject('AttachConstraint', object1='@soft', object2='@actuator2', twoWay=True, indices1=[2, 5, 8, 11], indices2=[4, 6, 8, 10])
    root.addObject('AttachConstraint', object1='@actuator1', object2='@actuator2', twoWay=True, indices1=[2, 5, 8, 11], indices2=[0, 2, 4, 6])

    root.addObject(FingerController(name="ActuationController", rootNode=root, node=actuator1, act1_size=2))
    root.addObject(FingerController2(name="ActuationController2", rootNode=root, node=actuator2, act2_size=2))

def main():
    import SofaRuntime
    import Sofa.Gui
    # Make sure to load all SOFA libraries
    SofaRuntime.importPlugin("SofaComponentAll")
    #Create the root node
    root = Sofa.Core.Node("root")
    # Call the below 'createScene' function to create the scene graph
    createScene(root)
    Sofa.Simulation.init(root)

    if not USE_GUI:
        for iteration in range(1000):
            Sofa.Simulation.animate(root, root.dt.value)
    else:
        # Find out the supported GUIs
        print ("Supported GUIs are: " + Sofa.Gui.GUIManager.ListSupportedGUI(","))
        # Launch the GUI (qt or qglviewer)
        Sofa.Gui.GUIManager.Init("myscene", "qglviewer")
        Sofa.Gui.GUIManager.createGUI(root, __file__)
        Sofa.Gui.GUIManager.SetDimension(1080, 1080)
        # Initialization of the scene will be done here
        Sofa.Gui.GUIManager.MainLoop(root)
        Sofa.Gui.GUIManager.closeGUI()
        print("GUI was closed")

    print("Simulation is done.")
    softTotal = [0.0, 0.0, 0.0]
    softPos = root.soft.softDOF.position.value
    for vertex in softPos:
        softTotal += vertex
    actuator1Pos = root.actuator1.actuator1DOF.position.value
    for vertex in actuator1Pos:
        softTotal += vertex
    softAverage = softTotal/(len(softPos)+len(actuator1Pos))
    print(softAverage)
    fitness = ((softAverage[0])**2 + (softAverage[1])**2 + (softAverage[2])**2)**0.5
    print(fitness)

# Function used only if this script is called from a python environment
if __name__ == '__main__':
    main()

The mesh files can be found here