Attaching 2 rods to a platform with different objects

Hello everyone,

Maybe my post was not as clear as I thought so I created something easier. I created 2 beams. What I want is to connect the first indices of both beams.

I tried two things:

– I tried to use AttachConstraint but it did not work. The two beams were linked but when I moved one the other one did not follow even if there was the link between them.
I tried to take some inspiration from

https://www.sofa-framework.org/community/forum/topic/how-to-do-attachconstraint/

I tried to update the position and rest position of the beam but it did not work.

– I also tried to use BilateralInteractionConstraint as propose in the following post:

https://www.sofa-framework.org/community/forum/topic/best-way-to-implement-ligaments-attaching-finger-bones-together/

But error saying that the determinant of a matrix is too small but I do not see how I can solve that

You can find the BeamMultiCompMoving.py file where both options are implemented at the following link.

Can you help me on this connection between two points of different objects, please ?

Cheers,
Benjamin

Hi @Hugo,

Thank you for your answer. Understandably, you take time to answer due to the number of people asking for your help.

That’s true that my SOFA version starts being old. You gave me that version (containing the soft robot plugin and the inverse plugin) on a GDR day. If I want to update it, I don’t know if it will be easy to get the same plugins.

– For the AttachConstraint and the BilaterallC, I’m using the same data fields (two objects (beams) with only one index).
– If I commented on the BilaterallC, there is no error, just a warning about the denominator threshold of GC that is reached at the first iteration.

Sorry for the link to GitHub. So here is the code that I’m using.

import Sofa
from splib.constants import Key
# mm and kg units
# BEAM_LENGTH = 100  # in mm
# BEAM_RADIUS = 0.5  # in mm
# BEAM_DENSITY = 6.45e-6  # in kg/mm^3
# YOUNG_MODULUS = 160e6  # 160 GPa = 160e9 Pa = 160e9 N/(m.s^2) = 160e6 N/(mm.s^2)
# ELEMENT_NUMBER = 10
# GRAVITATIONAL_ACCELERATION = 9806.65  # in mm/s^2
# STEP_TIME = 0.001  # in s
# LOAD_MASS = 0.5  # in kg

# SI Units
BEAM_LENGTH = 0.1  # in m
BEAM_RADIUS = 0.5e-3  # in m
BEAM_DENSITY = 6450  # in kg/m^3
YOUNG_MODULUS = 160e9  # in Pa
ELEMENT_COUNT = 5
GRAVITATIONAL_ACCELERATION = 9.80665  # in m/s^2
STEP_TIME = 0.01  # in s
LOAD_MASS = 0.1  # in kg

# OpenGL
SLICE_COUNT = 100
SECTOR_COUNT = 5

def createScene(rootNode):
    print('---------- Entering createScene ----------')

    rootNode.createObject('RequiredPlugin',pluginName='SoftRobots BeamAdapter SofaPython SofaSparseSolver')

    rootNode.createObject('VisualStyle',
                          displayFlags='showVisualModels showBehaviorModels showCollisionModels hideBoundingCollisionModels showForceFields showInteractionForceFields hideWireframe')

    rootNode.createObject('BackgroundSetting', color=[0.0, 0.0, 0.0, 1.0])

    rootNode.dt = STEP_TIME
    rootNode.gravity = [0.0, -GRAVITATIONAL_ACCELERATION, 0.0]
    rootNode.createObject('GenericConstraintSolver', name="genericConstraintSolver1", maxIterations="10000", tolerance="1e-12")
    rootNode.createObject('FreeMotionAnimationLoop')
    # rootNode.createObject('DefaultPipeline', verbose="0")
    # rootNode.createObject('BruteForceDetection', name="N2")
    # rootNode.createObject('DefaultContactManager', response="FrictionContact", responseParams="mu=0.6")
    # rootNode.createObject('LocalMinDistance', name="Proximity", alarmDistance="0.05", contactDistance="0.001", angleCone="0.00")
    rootNode.createObject('OglSceneFrame', style="Arrows", alignment="TopRight")
    # rootNode.createObject('DefaultPipeline', verbose="0")
    # rootNode.createObject('BruteForceDetection', name="N2")
    # rootNode.createObject('DefaultPipeline', verbose="0")
    rootNode.createObject('PythonScriptController', classname="Controller")
    
    ##### Compute the beam topology using python
    edgesList = []
    for i in range(0, ELEMENT_COUNT):
        edgesList.append([i, i + 1])

    print("edgesList = ", edgesList)

    positionsList = []
    for i in range(0, ELEMENT_COUNT + 1):
        dx = BEAM_LENGTH / ELEMENT_COUNT
        positionsList.append([dx * i, 0.0, 0.0, 0.0, 0.0, 0.0, 1])

    print("positionsList = ", positionsList)

    wall = rootNode.createChild('wall')
    wall.createObject('EulerImplicitSolver', name="odesolverwall", rayleighStiffness="0.1", rayleighMass="0.1")
    # wall.createObject('SparseLDLSolver')
    wall.createObject('CGLinearSolver', iterations="100", tolerance="1.0e-9", threshold="1.0e-9" )
    # wall.createObject('BTDLinearSolver', printLog="false", verbose="false")
    wall.createObject('MechanicalObject', name='wallFrame', template="Rigid3", position=['0.0 -0.05 0.0  0.0 0.0 0.707 0.707','0.0 0.0 0.0  0.0 0.0 0.707 0.707','0.0 0.05 0.0  0.0 0.0 0.707 0.707'])
    wall.createObject('MeshTopology', edges=['0 1','1 2'])
    # beam.createObject('BeamInterpolation', name='interpolation', crossSectionShape='rectangular', lengthY=2*BEAM_RADIUS, lengthZ=2*BEAM_RADIUS, defaultYoungModulus=YOUNG_MODULUS)
    wall.createObject('BeamInterpolation', name='wallinter', crossSectionShape='circular', radius=BEAM_RADIUS,
                      innerRadius=0.0, defaultYoungModulus=YOUNG_MODULUS, dofsAndBeamsAligned=1)
    wall.createObject('AdaptiveBeamForceFieldAndMass', name='wallForceField', computeMass=1, massDensity=BEAM_DENSITY)
    wall.createObject('PartialFixedConstraint', name='clampingProxEnd', indices=['0','1','2'], fixedDirections='0 1 1 1 1 1')
    wall.createObject('UniformMass', totalMass="0.1", showAxisSizeFactor="0.01")
    # wall.createObject('LinearSolverConstraintCorrection')
    wall.createObject('UncoupledConstraintCorrection')
    
    ##### Beam model
    beam = rootNode.createChild('beam')
    beam.createObject('EulerImplicitSolver', name="odesolverbeam", rayleighStiffness="0.1", rayleighMass="0.1")
    beam.createObject('CGLinearSolver', iterations="100", tolerance="1.0e-9", threshold="1.0e-9" )
    # beam.createObject('SparseLDLSolver')
    # beam.createObject('BTDLinearSolver', printLog="false", verbose="false")
    beam.createObject('MechanicalObject', name='frame', template="Rigid3", position=positionsList)
    beam.createObject('MeshTopology', edges=edgesList)
    # beam.createObject('BeamInterpolation', name='interpolation', crossSectionShape='rectangular', lengthY=2*BEAM_RADIUS, lengthZ=2*BEAM_RADIUS, defaultYoungModulus=YOUNG_MODULUS)
    beam.createObject('BeamInterpolation', name='interpolation', crossSectionShape='circular', radius=BEAM_RADIUS,
                      innerRadius=0.0, defaultYoungModulus=YOUNG_MODULUS, DOF0TransformNode0=['0.0 0.0 0.0  0.0 0.0 0.0 0.0','0.0 0.0 0.0  0.0 0.0 0.0 0.0','0.0 0.0 0.0  0.0 0.0 0.0 0.0'],
                      DOF1TransformNode1=['0.0 0.0 0.0  0.0 0.0 0.0 0.0','0.0 0.0 0.0  0.0 0.0 0.0 0.0','0.0 0.0 0.0  0.0 0.0 0.0 0.0'], dofsAndBeamsAligned=1.0)
    beam.createObject('AdaptiveBeamForceFieldAndMass', name='BeamForceField', computeMass=1, massDensity=BEAM_DENSITY)
    beam.createObject('UniformMass', totalMass="0.1", showAxisSizeFactor="0.01" )
    # beam.createObject('LinearSolverConstraintCorrection')
    beam.createObject('UncoupledConstraintCorrection')
    ##### Clamping
    # beam.createObject('FixedConstraint', name='clamping', indices='0')
    beam.createObject('PartialFixedConstraint', name='clampingProxEnd', indices='0', fixedDirections='0 1 1 1 1 1')
    beam.createObject('FixedConstraint', name='clampingDisEnd', indices=ELEMENT_COUNT)

    # rootNode.createObject('AttachConstraint', name='constforPal', object1="@wall", object2="@beam", indices1="0", indices2="0", constraintFactor="1",twoWays='1')
    # rootNode.createObject('BilateralInteractionConstraint', template="Rigid3", object1="@wall", object2="@beam", first_point="1", second_point="0")
    
    print('---------- Exiting createScene ----------')

    return rootNode

class Controller(Sofa.PythonScriptController):

    def createGraph(self, node):
        print('---------- Entering createGraph ----------')
        self.node = node
        print('---------- Exiting createGraph ----------')

    def reset(self):
        print('---------- Entering reset ----------')
        self.time = 0.0
        print('---------- Exiting reset ----------')
        return 0

    def onBeginAnimationStep(self, dt):
        return 0

    def onEndAnimationStep(self, dt):
        return 0

    def onKeyPressed(self, key):
        print("Key Pressed")
        if key == Key.U: # mouvements suivant l'axe X
            self.node.wall.wallFrame.rest_position=self.node.wall.wallFrame.position
        return 0