How to import a geometry (centreline) as .msh and use beam elements for the simulation

[athishsm]

Hi,

I’m modelling nitinol wire springs using 1D beam elements and trying to replicate my Abaqus simulations in SOFA. I imported the spring’s centerline, meshed a .step file in GMSH, and brought the .msh into SOFA.

However, when I activate the visual model (either in the view tab or in the code), SOFA crashes. It works fine when I only activate the collision models without the visual model. The spring appears as just a series of points (white dots in the image (black background)), and only one point moves when interacted with. Ideally, I’d like to interact with the entire spring, but I’m not sure what I’m doing wrong.

Please let me know if there's any resource that tells the user how to import a line geometry (.stp or .iges or .msh or .obj) file into SOFA, mesh it using beam elements and finally interact with it.

Also, could you direct me to any relevant resources or provide insight into what might be causing the crash when visual model is turned on? I’ve attached the .SCN file for reference in case they help.

The spring (as volume) looks like this (attached below). If volume mesh was imported with low Young’s modulus, I could interact with it (SHIFT + LeftMouseButton). With Young’s modulus of metals, the FPS drops, and I need to rethink the time step size and solver typer.

.scn file:

<?xml version="1.0" ?>
<Node name="root" dt="0.1" gravity="0 -9.81 0">

    <Node name="RequiredPlugins">
        <RequiredPlugin name="Sofa.Component.AnimationLoop"/> <!-- Needed to use components [FreeMotionAnimationLoop] -->
        <RequiredPlugin name="Sofa.Component.Collision.Detection.Algorithm"/> <!-- Needed to use components [BVHNarrowPhase BruteForceBroadPhase CollisionPipeline] -->
        <RequiredPlugin name="Sofa.Component.Collision.Detection.Intersection"/> <!-- Needed to use components [MinProximityIntersection] -->
        <RequiredPlugin name="Sofa.Component.Collision.Geometry"/> <!-- Needed to use components [LineCollisionModel PointCollisionModel TriangleCollisionModel] -->
        <RequiredPlugin name="Sofa.Component.Collision.Response.Contact"/> <!-- Needed to use components [CollisionResponse] -->
        <RequiredPlugin name="Sofa.Component.Constraint.Lagrangian.Correction"/> <!-- Needed to use components [UncoupledConstraintCorrection] -->
        <RequiredPlugin name="Sofa.Component.Constraint.Lagrangian.Solver"/> <!-- Needed to use components [LCPConstraintSolver] -->
        <RequiredPlugin name="Sofa.Component.Constraint.Projective"/> <!-- Needed to use components [FixedConstraint] -->
		<RequiredPlugin name="Sofa.Component.IO.Mesh"/> <!-- Needed to use components [MeshOBJLoader] -->
        <RequiredPlugin name="Sofa.Component.LinearSolver.Iterative"/> <!-- Needed to use components [CGLinearSolver] -->
        <RequiredPlugin name="Sofa.Component.LinearSystem"/> <!-- Needed to use components [MatrixLinearSystem] -->
        <RequiredPlugin name="Sofa.Component.Mapping.Linear"/> <!-- Needed to use components [BarycentricMapping] -->
        <RequiredPlugin name="Sofa.Component.Mass"/> <!-- Needed to use components [UniformMass] -->
        <RequiredPlugin name="Sofa.Component.MechanicalLoad"/> <!-- Needed to use components [ConstantForceField] -->
		<RequiredPlugin name="Sofa.Component.ODESolver.Backward"/> <!-- Needed to use components [EulerImplicitSolver] -->
        <RequiredPlugin name="Sofa.Component.SolidMechanics.FEM.Elastic"/> <!-- Needed to use components [HexahedronFEMForceField] -->
        <RequiredPlugin name="Sofa.Component.StateContainer"/> <!-- Needed to use components [MechanicalObject] -->
        <RequiredPlugin name="Sofa.Component.Topology.Container.Constant"/> <!-- Needed to use components [MeshTopology] -->
        <RequiredPlugin name="Sofa.Component.Topology.Container.Dynamic"/> <!-- Needed to use components [TetrahedronSetGeometryAlgorithms, TetrahedronSetTopologyContainer] -->
		<RequiredPlugin name="Sofa.Component.Topology.Container.Grid"/> <!-- Needed to use components [SparseGridRamificationTopology] -->
        <RequiredPlugin name="Sofa.Component.Topology.Mapping"/> <!-- Needed to use components [Tetra2TriangleTopologicalMapping] -->
		<RequiredPlugin name="Sofa.Component.Visual"/> <!-- Needed to use components [InteractiveCamera VisualStyle] -->
        <RequiredPlugin name="Sofa.GL.Component.Rendering3D"/> <!-- Needed to use components [OglModel] -->
        <RequiredPlugin name="Sofa.GL.Component.Shader"/> <!-- Needed to use components [LightManager SpotLight] -->
    </Node>
    
	<!--<DefaultAnimationLoop/>-->
    <DefaultAnimationLoop computeBoundingBox="false"/>

	<!--<VisualStyle displayFlags="showVisual  " />-->
    <VisualStyle displayFlags="showCollisionModels hideVisualModels" />

    <CollisionPipeline />
    <BruteForceBroadPhase/>
    <BVHNarrowPhase/>
    <MinProximityIntersection name="Proximity" alarmDistance="0.5" contactDistance="0.25"/>
    <CollisionResponse name="Response" response="PenalityContactForceField"/>

	<!--<LightManager />
    <SpotLight name="light1" color="1 1 1" position="0 80 25" direction="0 -1 -0.8" cutoff="30" exponent="1" />
    <SpotLight name="light2" color="1 1 1" position="0 40 100" direction="0 0 -1" cutoff="30" exponent="1" />
	-->
    <MeshGmshLoader name="meshLoaderCoarse" filename="mesh/Spring_Beam_Centreline.msh" />
    <MeshOBJLoader name="meshLoaderCoarseOBJ" filename="mesh/Spring_Beam_Centreline.obj" />
	
	<MeshGmshLoader name="meshLoaderCoarse_crimper" filename="mesh/crimper_for_beam_stent.msh" />
	<MeshOBJLoader name="loader_crimper" filename="mesh/crimper_for_beam_stent.obj" />
	
    <Node name="stent_30D">
        <EulerImplicitSolver name="cg_odesolver"   rayleighStiffness="0.1" rayleighMass="0.1" />
        <CGLinearSolver name="linear solver" iterations="25" tolerance="1e-09" threshold="1e-09" />

        <EdgeSetTopologyContainer name="topo" src="@../meshLoaderCoarse" />
		<EdgeSetTopologyModifier name="colliseEdgeModifier" />
        <EdgeSetGeometryAlgorithms template="Vec3d" name="GeomAlgo"/>
        
        <MechanicalObject template="Vec3d" name="MechanicalModel" showObject="1"/>

        <BeamFEMForceField name="FEM" youngModulus="3000000" poissonRatio="0.3" radius="0.01" />

        <!--<MeshMatrixMass massDensity="1" topology="@topo"/>
        <ConstantForceField totalForce="100 0 0" />-->

        <Node name="Visual" >
            <OglModel name="VisualModel" src="@../../meshLoaderCoarse" />
            <BarycentricMapping name="VisualMapping" input="@../MechanicalModel" output="@VisualModel" />
        </Node>

        <Node name="Collision">
            <MeshTopology src="@../../meshLoaderCoarse"/>
            <MechanicalObject name="CollisionMO" scale="1.0"/>
            <TriangleCollisionModel name="CollisionModel" contactStiffness="3"/>
            <BarycentricMapping name="CollisionMapping" input="@../MechanicalModel" output="@CollisionMO" />
        </Node>
    </Node>
	
    <Node name="crimper">
        <EulerImplicitSolver rayleighStiffness="0"/>
        <CGLinearSolver iterations="200" tolerance="1e-06" threshold="1e-06"/>
		
		<!-- Keep tetrahedral mesh for crimper if it is a 3D structure -->
		<TetrahedronSetTopologyContainer name="topo" src="@../meshLoaderCoarse_crimper" />
        <TetrahedronSetGeometryAlgorithms template="Vec3d" name="GeomAlgo"/>
		<MechanicalObject template="Vec3d" name="MechanicalModel_crimper" showObject="1"/>
		<TetrahedronFEMForceField name="FEM" youngModulus="100000" poissonRatio="0.3" method="large" />
        <MeshMatrixMass massDensity="1" topology="@topo"/>
		<ConstantForceField totalForce="100 0 0" />
		
		<Node name="Visual_crimper" >
            <OglModel name="VisualModel_crimper" src="@../../loader_crimper" />
            <BarycentricMapping name="VisualMapping" input="@../MechanicalModel_crimper" output="@VisualModel_crimper" />
        </Node>
		
		<Node name="Collision_crimper">
            <MeshTopology src="@../../meshLoaderCoarse_crimper"/>
            <MechanicalObject name="CollisionMO_crimper" scale="1.0"/>
            <TriangleCollisionModel name="CollisionModel_crimper" contactStiffness="3"/>
            <BarycentricMapping name="CollisionMapping_crimper" input="@../MechanicalModel_crimper" output="@CollisionMO_crimper" />
        </Node>
        
    </Node>
</Node>

[athishsm]

Any help with this issue?

[hugtalbot]

Hi @athishsm

Could you please copy paste the output from the command window please when launching the scene?

[athishsm]

Hi @hugtalbot,
This is the snippet of command window when visual model was not activated.

The command window also crashes when I activate the visual model.

My main goal is to import the centreline of the 3D spring and just use beam elements and metal's material properties.

[hugtalbot]

Keep us posted @athishsm about your progress, we would happily assist you further
Do not hesitate to open a new topic if this import is solved

Cheers

[athishsm]

Dear @hugtalbot,
Thanks for the kind help.
Yes, I managed to import it now.

However, the edges (or nodes) of line objects on the geometry are not connected together (geometry is separated into multiple line objects or pieces) (marked by a red box). You suggested I use Bilateral constraints to connect the nodes of the two pieces.

I ticked on the viewIndices in the visualisation tab of the GeomAlgo to find the node numbers. Then, I tried joining two nodes (eg: 39,40) using Bilateral Constraints (as shown inside the red box). But the nodes get separated when I interact with the geometry. I attached Code 1 below.

With the DefaultAnimationLoop, I could at least interact with the geometry. But for constraints, I need to have a FreeMotionAnimationLoop and ConstraintSolver. So, when I use FreemotionAnimationLoop with the same settings (as Code 1), the geometry is not interactable. Then I changed the solver from CGLinearSolver to BTDLinearSolver and the geometry fell apart due to a numerical convergence error (Code 2 below). I'm unsure what I missed or what I need to know to implement this correctly. Could you please help me sort this out?

Code 1:

<Node name="RequiredPlugins">
    <RequiredPlugin name="Sofa.Component.AnimationLoop"/> <!-- Needed to use components [FreeMotionAnimationLoop] -->
    <RequiredPlugin name="Sofa.Component.Collision.Detection.Algorithm"/> <!-- Needed to use components [BVHNarrowPhase BruteForceBroadPhase CollisionPipeline] -->
    <RequiredPlugin name="Sofa.Component.Collision.Detection.Intersection"/> <!-- Needed to use components [MinProximityIntersection] -->
    <RequiredPlugin name="Sofa.Component.Collision.Geometry"/> <!-- Needed to use components [LineCollisionModel PointCollisionModel TriangleCollisionModel] -->
    <RequiredPlugin name="Sofa.Component.Collision.Response.Contact"/> <!-- Needed to use components [CollisionResponse] -->
    <RequiredPlugin name="Sofa.Component.Constraint.Lagrangian.Correction"/> <!-- Needed to use components [UncoupledConstraintCorrection] -->
    <RequiredPlugin name="Sofa.Component.Constraint.Lagrangian.Solver"/> <!-- Needed to use components [LCPConstraintSolver] -->
    <RequiredPlugin name="Sofa.Component.Constraint.Projective"/> <!-- Needed to use components [FixedConstraint] -->
	<RequiredPlugin name="Sofa.Component.IO.Mesh"/> <!-- Needed to use components [MeshOBJLoader] -->
	<RequiredPlugin name="Sofa.Component.LinearSolver.Direct"/> <!-- Needed to use components [EigenSparseLU] -->
    <RequiredPlugin name="Sofa.Component.LinearSolver.Iterative"/> <!-- Needed to use components [CGLinearSolver] -->
    <RequiredPlugin name="Sofa.Component.LinearSystem"/> <!-- Needed to use components [MatrixLinearSystem] -->
    <RequiredPlugin name="Sofa.Component.Mapping.Linear"/> <!-- Needed to use components [BarycentricMapping] -->
    <RequiredPlugin name="Sofa.Component.Mass"/> <!-- Needed to use components [UniformMass] -->
    <RequiredPlugin name="Sofa.Component.MechanicalLoad"/> <!-- Needed to use components [ConstantForceField] -->
	<RequiredPlugin name="Sofa.Component.ODESolver.Backward"/> <!-- Needed to use components [EulerImplicitSolver] -->
    <RequiredPlugin name="Sofa.Component.SolidMechanics.FEM.Elastic"/> <!-- Needed to use components [HexahedronFEMForceField] -->
    <RequiredPlugin name="Sofa.Component.StateContainer"/> <!-- Needed to use components [MechanicalObject] -->
    <RequiredPlugin name="Sofa.Component.Topology.Container.Constant"/> <!-- Needed to use components [MeshTopology] -->
    <RequiredPlugin name="Sofa.Component.Topology.Container.Dynamic"/> <!-- Needed to use components [TetrahedronSetGeometryAlgorithms, TetrahedronSetTopologyContainer] -->
	<RequiredPlugin name="Sofa.Component.Topology.Container.Grid"/> <!-- Needed to use components [SparseGridRamificationTopology] -->
    <RequiredPlugin name="Sofa.Component.Topology.Mapping"/> <!-- Needed to use components [Tetra2TriangleTopologicalMapping] -->
	<RequiredPlugin name="Sofa.Component.Visual"/> <!-- Needed to use components [InteractiveCamera VisualStyle] -->
    <RequiredPlugin name="Sofa.GL.Component.Rendering3D"/> <!-- Needed to use components [OglModel] -->
    <RequiredPlugin name="Sofa.GL.Component.Shader"/> <!-- Needed to use components [LightManager SpotLight] -->
	<RequiredPlugin name="BeamAdapter"/> <!-- Needed to use components [AdaptiveBeamController,AdaptiveBeamForceFieldAndMass,BeamInterpolation] -->
</Node>

<DefaultAnimationLoop/>
<!--<FreeMotionAnimationLoop/>-->
<VisualStyle displayFlags="showBehaviorModels hideForceFields showCollisionModels showVisual showInteractionForceFields" />
<CollisionPipeline depth="6" verbose="0" draw="0" />
<BruteForceBroadPhase/>
<BVHNarrowPhase/>
<MinProximityIntersection name="Proximity" alarmDistance="0.3" contactDistance="0.2" />
<CollisionResponse name="Response" response="PenalityContactForceField" />

<MeshGmshLoader name="meshLoaderCoarse" filename="mesh/Spring_Beam_Centreline.msh" />
<MeshOBJLoader name="meshLoaderCoarseOBJ" filename="mesh/Spring_Beam_Centreline.obj" />

<MeshGmshLoader name="meshLoaderCoarse_crimper" filename="mesh/crimper_for_beam_stent.msh" />
<MeshOBJLoader name="loader_crimper" filename="mesh/crimper_for_beam_stent.obj" />

<Node name="stent_30D">
	<!--<EulerImplicitSolver rayleighStiffness="0" rayleighMass="0" printLog="false" />
	<EigenSparseLU template="CompressedRowSparseMatrixMat3x3d"/>
	
	<EulerImplicitSolver rayleighStiffness="0" printLog="false"  rayleighMass="0.1" />
    <BTDLinearSolver template="BTDMatrix6d" verbose="true" />-->
	
	<EulerImplicitSolver name="cg_odesolver"   rayleighStiffness="0.1" rayleighMass="0.1" />
    <CGLinearSolver name="linear solver" iterations="25" tolerance="1e-09" threshold="1e-09" />
	
    <EdgeSetTopologyContainer name="topo" src="@../meshLoaderCoarse" />
	<EdgeSetTopologyModifier name="colliseEdgeModifier" />		
    <EdgeSetGeometryAlgorithms template="Rigid3d" name="GeomAlgo" drawEdges="true" edgeScale="1.0"/>
	
    <MechanicalObject template="Rigid3" name="MechanicalModel" showObject="1" />
	
	<!--<BeamInterpolation name="BeamInterpolation" radius="0.1"/> 
	<AdaptiveBeamForceFieldAndMass name="BeamForceField"   computeMass="1" massDensity="500"/> 
	<FixedProjectiveConstraint name="FixedProjectiveConstraint" indices="0" />-->
	<BeamFEMForceField name="FEM" radius="0.00025" radiusInner="0" youngModulus="2.03e11" poissonRatio="0.3"/>
	
	<Node name="Visual1" >
        <OglModel name="VisualModel" src="@../../meshLoaderCoarse" />
    </Node>
	
	<Node name="Collision1">
        <MechanicalObject />
		<IdentityMapping  />
    </Node>

</Node>
<BilateralLagrangianConstraint template="Rigid3" object1="@stent_30D/MechanicalModel" object2="@stent_30D/MechanicalModel" first_point="39" second_point="40" />

Code 2:

<Node name="RequiredPlugins">
    <RequiredPlugin name="Sofa.Component.AnimationLoop"/> <!-- Needed to use components [FreeMotionAnimationLoop] -->
    <RequiredPlugin name="Sofa.Component.Collision.Detection.Algorithm"/> <!-- Needed to use components [BVHNarrowPhase BruteForceBroadPhase CollisionPipeline] -->
    <RequiredPlugin name="Sofa.Component.Collision.Detection.Intersection"/> <!-- Needed to use components [MinProximityIntersection] -->
    <RequiredPlugin name="Sofa.Component.Collision.Geometry"/> <!-- Needed to use components [LineCollisionModel PointCollisionModel TriangleCollisionModel] -->
    <RequiredPlugin name="Sofa.Component.Collision.Response.Contact"/> <!-- Needed to use components [CollisionResponse] -->
    <RequiredPlugin name="Sofa.Component.Constraint.Lagrangian.Correction"/> <!-- Needed to use components [UncoupledConstraintCorrection] -->
    <RequiredPlugin name="Sofa.Component.Constraint.Lagrangian.Solver"/> <!-- Needed to use components [LCPConstraintSolver] -->
    <RequiredPlugin name="Sofa.Component.Constraint.Projective"/> <!-- Needed to use components [FixedConstraint] -->
	<RequiredPlugin name="Sofa.Component.IO.Mesh"/> <!-- Needed to use components [MeshOBJLoader] -->
	<RequiredPlugin name="Sofa.Component.LinearSolver.Direct"/> <!-- Needed to use components [EigenSparseLU] -->
    <RequiredPlugin name="Sofa.Component.LinearSolver.Iterative"/> <!-- Needed to use components [CGLinearSolver] -->
    <RequiredPlugin name="Sofa.Component.LinearSystem"/> <!-- Needed to use components [MatrixLinearSystem] -->
    <RequiredPlugin name="Sofa.Component.Mapping.Linear"/> <!-- Needed to use components [BarycentricMapping] -->
    <RequiredPlugin name="Sofa.Component.Mass"/> <!-- Needed to use components [UniformMass] -->
    <RequiredPlugin name="Sofa.Component.MechanicalLoad"/> <!-- Needed to use components [ConstantForceField] -->
	<RequiredPlugin name="Sofa.Component.ODESolver.Backward"/> <!-- Needed to use components [EulerImplicitSolver] -->
    <RequiredPlugin name="Sofa.Component.SolidMechanics.FEM.Elastic"/> <!-- Needed to use components [HexahedronFEMForceField] -->
    <RequiredPlugin name="Sofa.Component.StateContainer"/> <!-- Needed to use components [MechanicalObject] -->
    <RequiredPlugin name="Sofa.Component.Topology.Container.Constant"/> <!-- Needed to use components [MeshTopology] -->
    <RequiredPlugin name="Sofa.Component.Topology.Container.Dynamic"/> <!-- Needed to use components [TetrahedronSetGeometryAlgorithms, TetrahedronSetTopologyContainer] -->
	<RequiredPlugin name="Sofa.Component.Topology.Container.Grid"/> <!-- Needed to use components [SparseGridRamificationTopology] -->
    <RequiredPlugin name="Sofa.Component.Topology.Mapping"/> <!-- Needed to use components [Tetra2TriangleTopologicalMapping] -->
	<RequiredPlugin name="Sofa.Component.Visual"/> <!-- Needed to use components [InteractiveCamera VisualStyle] -->
    <RequiredPlugin name="Sofa.GL.Component.Rendering3D"/> <!-- Needed to use components [OglModel] -->
    <RequiredPlugin name="Sofa.GL.Component.Shader"/> <!-- Needed to use components [LightManager SpotLight] -->
	<RequiredPlugin name="BeamAdapter"/> <!-- Needed to use components [AdaptiveBeamController,AdaptiveBeamForceFieldAndMass,BeamInterpolation] -->
</Node>

<!--<DefaultAnimationLoop/>-->
<FreeMotionAnimationLoop/>
<GenericConstraintSolver tolerance="0.00001" maxIterations="1000"/>
<VisualStyle displayFlags="showBehaviorModels hideForceFields showCollisionModels showVisual showInteractionForceFields" />
<CollisionPipeline depth="6" verbose="0" draw="0" />
<BruteForceBroadPhase/>
<BVHNarrowPhase/>
<MinProximityIntersection name="Proximity" alarmDistance="0.3" contactDistance="0.2" />
<CollisionResponse name="Response" response="PenalityContactForceField" />

<MeshGmshLoader name="meshLoaderCoarse" filename="mesh/Spring_Beam_Centreline.msh" />
<MeshOBJLoader name="meshLoaderCoarseOBJ" filename="mesh/Spring_Beam_Centreline.obj" />

<MeshGmshLoader name="meshLoaderCoarse_crimper" filename="mesh/crimper_for_beam_stent.msh" />
<MeshOBJLoader name="loader_crimper" filename="mesh/crimper_for_beam_stent.obj" />

<Node name="stent_30D">
	<!--<EulerImplicitSolver rayleighStiffness="0" rayleighMass="0" printLog="false" />
	<EigenSparseLU template="CompressedRowSparseMatrixMat3x3d"/>
	
	<EulerImplicitSolver rayleighStiffness="0" printLog="false"  rayleighMass="0.1" />
    <BTDLinearSolver template="BTDMatrix6d" verbose="true" />-->
	
	<EulerImplicitSolver name="cg_odesolver"   rayleighStiffness="0.1" rayleighMass="0.1" />
    <!--<CGLinearSolver name="linear solver" iterations="25" tolerance="1e-09" threshold="1e-09" />-->
	<BTDLinearSolver printLog="false" verbose="false" />
	
    <EdgeSetTopologyContainer name="topo" src="@../meshLoaderCoarse" />
	<EdgeSetTopologyModifier name="colliseEdgeModifier" />		
    <EdgeSetGeometryAlgorithms template="Rigid3d" name="GeomAlgo" drawEdges="true" edgeScale="1.0"/>
	
    <MechanicalObject template="Rigid3" name="MechanicalModel" showObject="1" />
	
	<!--<BeamInterpolation name="BeamInterpolation" radius="0.1"/> 
	<AdaptiveBeamForceFieldAndMass name="BeamForceField"   computeMass="1" massDensity="500"/> 
	<FixedProjectiveConstraint name="FixedProjectiveConstraint" indices="0" />-->
	<BeamFEMForceField name="FEM" radius="0.00025" radiusInner="0" youngModulus="2.03e11" poissonRatio="0.3"/>
	<LinearSolverConstraintCorrection />
	
	<Node name="Visual1" >
        <OglModel name="VisualModel" src="@../../meshLoaderCoarse" />
    </Node>
	
	<Node name="Collision1">
        <MechanicalObject />
		<IdentityMapping  />
    </Node>

</Node>
<BilateralLagrangianConstraint template="Rigid3" object1="@stent_30D/MechanicalModel" object2="@stent_30D/MechanicalModel" first_point="39" second_point="40" />