25 February 2020 at 3 h 16 min #15235
Hi SOFA Community,
I am working on a project where I want to try to model the biomechanics of a moving musculo-skeletal system. Right now I am relying on using the Soft Robotics plugin. For now, I want to try to implement a human finger with the bones and ligaments in the sense of having soft ligament models attach to the rigid bones models at the attachment points as it is in human anatomy. I have the bone models, and some simplified models of ligaments (pretty much rectangular prisms); how should I approach this?
As a corollary, from preliminary testing in SOFA, it seems like running the SOFA simulation with even 1 soft ligament (about 300 tetrahedrons) contacting a floor plane causes the simulation to run much slower (about 8 FPS). In contrast, just rigid bodies contacting each other runs consistently about 4 times faster (32 FPS). Should I even try to model soft bodies when I may be trying to have 12 ligaments per finger? Or should I try to model the finger bones attached to each other with “RestShapeSpringsForceField”‘s in place of ligaments (although from my understanding, this has the issue of not modeling the ligaments wrapping around the bone)?2 March 2020 at 3 h 23 min #15260
Right now I am trying to implement “attachConstraint” between a rigid and deformable part. I am not sure how to select the indices of a rigid part. I am able to select the indices of a deformable part using “BoxROI”, but it seems “BoxROI” does not output any vertices, lines, triangles, or anything else for the rigid part. Could I get advice on that?3 March 2020 at 10 h 31 min #15261WongParticipant
I think you can take a look at this link to see if it is suitable for you.9 March 2020 at 12 h 14 min #15343NickHockParticipant
- CSIRO, Australia
One option would be to use the tools-for-multi-material-simulations, on Sofa Marketplace. (Disclosure, I’m the author.)
You would build a partitioned mesh for tetrahedral FEM, from the surface meshes of the bones, ligaments, tendons and muscles. The ligaments may not segment well from a CT/MRI scan. You can generate model for them using the swept spline tool. You can set the ligaments and tendons to a lower modulus than the bones. Provided there is a finite gap between the bones, then the surface mesh of the union can be used as the collision surface. You would also be able to actuate the muscles by changing their modulus with a python script controller.
The extensor tendon hood is particularly complex, as are the dermal ligaments if you are modelling the skin and subcutis. I would suggest working from a high resolution MRI scan plus specialist hand anatomy texts. You can create anisotropy by adding springs of the right orientation, to mimic collagen fibres in the layers of the tendon hood. With the exception of the synovial sheath of the flexor tendons, the soft tissues of the finger are a continuous hyper-elastic matrix, but the extreme distortion of the finite elements is likely to make even Sofa unstable. You would be more likely to achieve a stable simulation by introducing collision surfaces between the bones and the extensor tendon hood.10 March 2020 at 8 h 52 min #15382
Hi @Wong, thanks for the link. That video does look really good, but I’m not sure how to implement that research for myself very easily. If someone could help point me to implementing that in a reasonable matter in a couple of weeks or using some otherwise similar software, I would really appreciate it.
Hi @NickHock, thank you for your advice. Do you have any resources on where I can get these models of high resolution MRI scans? From my (perhaps limited) understanding, in these MRI scans, it’s difficult to differentiate ligaments, and that’s why I am currently using some polygonal shapes made in SolidWorks that approximate them. I have good models for the bones themselves independent of soft tissue. My approach is to have rounded contact surfaces created from my bone models for the joint between the proximal phalange and metacarpal (the metacarpophalangeal joint) and to attach two linear collateral ligaments and two linear volar ligaments to link these bones (instead of a joint capsule with ligaments that completely surround the bones)
At this point, I am trying to use “Rigidify” as implemented in the “tripod” tutorial in the Soft Robotics plugin. Unfortunately, I am getting a very unstable simulation when attaching two bones with two linear ligaments collateral. If I push a bone at all, at some point the bone starts to fly away, pulling the ligaments along with them. I am using “RestShapeSpringsForceField” to attach the ligaments to the bones.
# -*- coding: utf-8 -*- from stlib.physics.deformable import ElasticMaterialObject from stlib.scene import Node from stlib.physics.rigid import Floor, RigidObject import Sofa # from bones import MCfab from splib.objectmodel import SofaPrefab, SofaObject from stlib.physics.mixedmaterial import Rigidify from splib.numerics import vec3 from splib.numerics.quat import Quat from tutorial import * def rigidifyHelper(rootNode, mn, nodeLig, gind, f, n, attachment, origin): rigidifiedstruct = Rigidify(mn, nodeLig, groupIndices=gind, frames=f, name=n) rigidifiedstruct.DeformableParts.createObject("UncoupledConstraintCorrection") rigidifiedstruct.RigidParts.createObject("UncoupledConstraintCorrection") # Use this to activate some rendering on the rigidified object setData(rigidifiedstruct.RigidParts.dofs, showObject=True, showObjectScale=5, drawMode=2) setData(rigidifiedstruct.RigidParts.RigidifiedParticules.dofs, showObject=True, showObjectScale=0.1, drawMode=1, showColor=[1., 1., 0., 1.]) setData(rigidifiedstruct.DeformableParts.dofs, showObject=True, showObjectScale=0.1, drawMode=2) rigidifiedstruct.RigidParts.createObject("RestShapeSpringsForceField", name="s1", points=0, external_rest_shape=attachment.mstate, stiffness="1e100", angularStiffness="1e100", drawSpring=True) rigidifiedstruct.RigidParts.createObject("RestShapeSpringsForceField", name="s2", points=1, external_rest_shape=origin.mstate, stiffness="1e100", angularStiffness="1e100", drawSpring=True) rigidifiedstruct.RigidParts.addChild(attachment) rigidifiedstruct.RigidParts.addChild(origin) rootNode.Simulation.addChild(rigidifiedstruct) def createScene(rootNode): # -*- coding: utf-8 -*- from stlib.scene import MainHeader, ContactHeader from stlib.visuals import ShowGrid from stlib.physics.rigid import Floor from stlib.physics.rigid import Cube m=MainHeader(rootNode, plugins=["SoftRobots", "SofaSparseSolver"], gravity=[0, -9.810, 0], dt=0.001) voffset = 30.0 cd = 0.25 ContactHeader(rootNode, alarmDistance=2, contactDistance=cd, frictionCoef=0.08) groupIndices1 =  frames1 =  groupIndices2 =  modelNode = rootNode.createChild("Modeling") DIP_AL1 = ElasticMaterialObject(modelNode, name="DIP_AL1", volumeMeshFileName="mesh/DIP_AL1_V2_2_noholes.vtk", surfaceMeshFileName="mesh/DIP_AL1_V2_2_noholes.stl", rotation = [90,-7,70], collisionMesh="mesh/DIP_AL1_V2_2_noholes.stl", translation=[-2.0,63.0,10.8], youngModulus=1e100, poissonRatio=0.45, totalMass=0.32e10, surfaceColor=[219.0/255,112.0/255,147.0/255]) DIP_AL2 = ElasticMaterialObject(modelNode, name="DIP_AL2", volumeMeshFileName="mesh/DIP_AL1_V2_2_noholes.vtk", surfaceMeshFileName="mesh/DIP_AL1_V2_2_noholes.stl", rotation = [90,7,70], collisionMesh="mesh/DIP_AL1_V2_2_noholes.stl", translation=[-2.0,63.0,-8.8], youngModulus=1E100, poissonRatio=0.45, totalMass=0.32e10, surfaceColor=[219.0/255,112.0/255,147.0/255]) # ElasticMaterialObject(rootNode, name="DIP_AL", volumeMeshFileName="mesh/DIP_AL1_V2_2_noholes.vtk", surfaceMeshFileName="mesh/DIP_AL1_V2_2_noholes.stl", rotation = [90,5,73], collisionMesh="mesh/DIP_AL1_V2_2_noholes.stl", translation=[-1.5,63.0,-7.0], surfaceColor=[219.0/255,112.0/255,147.0/255]) # DP = RigidObject(modelNode, name="DP", surfaceMeshFileName="mesh/simhandV2_DP_fine.obj", # translation=[21.06,128.8+voffset+cd+2,3.72], color=[1,1,1], totalMass = 0.39e-3) # MP = RigidObject(modelNode, name="MP", surfaceMeshFileName="mesh/simhandV2_MP_fine.obj", # translation=[11.59,101.51+voffset+cd+1,3.09], color=[1,1,1], totalMass = 1.56e-3) PP = RigidObject(modelNode, name="PP", surfaceMeshFileName="mesh/simhandV2_PP_fine.obj", translation=[2.83,58.02+voffset+cd,2.40], color=[1,1,1], totalMass = 4.73e-3) MC = RigidObject(modelNode, name="MC", surfaceMeshFileName="mesh/simhandV1_MC_fine.obj", translation=[0.0,0.0+voffset,0.0], color=[1,1,1], isAStaticObject=True, totalMass = 9.12e-3) Floor(modelNode, isAStaticObject=True, color=[0.9,0.9,0.9]) dis1 = DIP_AL1.createObject("SphereROI", name="distal_roi", template="Vec3", centers="1.5 72 12", radii="2", drawSphere="true", drawSize=1) pro1 = DIP_AL1.createObject("SphereROI", name="proximal_roi", template="Vec3", centers="-5.2 54 9.6", radii="2", drawSphere="true", drawSize=1) dis2 = DIP_AL2.createObject("SphereROI", name="distal_roi", template="Vec3", centers="1.5 72 -10.5", radii="2", drawSphere="true", drawSize=1) pro2 = DIP_AL2.createObject("SphereROI", name="proximal_roi", template="Vec3", centers="-5.2 54 -9.6", radii="2", drawSphere="true", drawSize=1) dis1.init() pro1.init() dis2.init() pro2.init() DIP_AL1.init() DIP_AL2.init() groupIndices1.append([ind for ind in dis1.indices]) groupIndices1.append([ind for ind in pro1.indices]) groupIndices2.append([ind for ind in dis2.indices]) groupIndices2.append([ind for ind in pro2.indices]) frames1.append([2.83,58.02+voffset+cd,2.40] + list(Quat.createFromEuler([0, 0, 0], inDegree=True))) frames1.append([0.0,0.0+voffset,0.0] + list(Quat.createFromEuler([0, 0, 0], inDegree=True))) simulationNode = rootNode.createChild("Simulation") simulationNode.createObject("EulerImplicitSolver") simulationNode.createObject("CGLinearSolver") rigidifyHelper(rootNode, modelNode, DIP_AL1, groupIndices1, frames1, "RigidifiedStructure1", rootNode.Modeling.PP, rootNode.Modeling.MC) rigidifyHelper(rootNode, modelNode, DIP_AL2, groupIndices2, frames1, "RigidifiedStructure2", rootNode.Modeling.PP, rootNode.Modeling.MC) m.getObject("VisualStyle").displayFlags="showForceFields showBehaviorModels" #ShowGrid(rootNode) return rootNode
has my files, using SoftRobotics v19.04 Windows 64-bit binary as found here: https://project.inria.fr/softrobot/install-get-started-2/download/10 March 2020 at 12 h 01 min #15386NickHockParticipant
- CSIRO, Australia
I have a full set of hand bones that I segmented from an MRI scan, and the original scans. If you contact me via email, I’ll send you a DropBox link. Using these requires an acknowledgement of the centre that took the scans for me.
Segmenting tissues takes time and practice. The ligaments are always hard to distinguish, and only the thickest ones are visible. The aponeuroses (sheet tendons) of the extensor hood are too fine to show up on MRI. You will always have to generate most ligaments from anatomical knowledge. The swept spline generator in the multi-material tools will help to generate appropriate shapes.
The Soft Robotics plugin can do many things, but requires setting up correctly to be stable. (i.e. Follow their instructions exactly and check and test each step of your install.) You should also be sure to start with simulation parameters (stiffness, damping, time step, mesh sizes etc) that are known to be stable.
You may well be able to do all that you want with just the standard Sofa Python plugin. The example .pyscn files in the multi-material-tools only use the Python plugin.13 April 2020 at 13 h 36 min #15689NouraParticipant
I came through this topic and I worked on muskuloskeletal simulations before, so I can share some thoughts which might be of interest.
I’m not familiar with Sofa Soft Robotics plugin, so I can not advice in this direction. Though, I think that it could be a convenient option.
– You mentioned that you are trying to implement
attachConstraintto connect the bones to the ligaments, I think that they are not the suitable choice.
These constraints are pretty simple, and don’t allow the propagation of forces in both directions.
Instead, you may need to use the
BilateralInteractionConstraintto connect bones with ligaments. The latter are implemented in SOFA based on Lagrange multipliers. A great documentation about them could be found here.
Note that if you go for this option, you have to use
FreeMotionAnimationLoopinstead of the
DefaultAnimationLoopwhich is the default loop set in any simulation.
In addition, you have to include the the constraint correction to define the way the compliance matrix is computed. A component called
UncoupledConstraintCorrectioncould be used for that, and an example scene called “BilateralInteractionConstraint.scn” is provided in SOFA.
Following this approach assumes that you model the soft tissue as 3D meshes and not as simple springs, you can get something like this simulation.
– Further, I noticed that you are using
RestShapeSpringsForceFieldin your simulation, If you mean If you intend to use simple springs, I’m not sure why you need to set attachment constraints?
StiffSpringForceFieldfor example, you can get a simple simuation like this (with the springs in green and red).
Hope this could inspire you!
Thanks for the interesting details about your approach. I totally agree with you that the segmentation and modelling of the wrist and carpal soft tissue are extremely challenging. I’ll have a look to the “swept spline tool” you referred to. It looks very interesting.
- You must be logged in to reply to this topic.