Added S, but still need to make the free motion not crash

Author: bakpaul

examples/stlib/SofaScene.py

@@ -1,19 +1,27 @@
 from stlib.geometries.plane import PlaneParameters
-from stlib.collision import CollisionParameters
-from stlib.collision import Collision
-from stlib.visual import Visual
-from splib.core.enum_types import CollisionPrimitive
-from splib.simulation.headers import setupLagrangianCollision
+from stlib.geometries.file import FileParameters
+from stlib.geometries.extract import ExtractParameters
+from stlib.materials.deformable import DeformableBehaviorParameters
+from stlib.collision import Collision, CollisionParameters
+from stlib.entities import Entity, EntityParameters
+from stlib.visual import Visual, VisualParameters
+from splib.core.enum_types import CollisionPrimitive, ElementType, ConstitutiveLaw
+from splib.simulation.headers import setupLagrangianCollision, setupDefaultHeader
+from splib.simulation.ode_solvers import addImplicitODE
+from splib.simulation.linear_solvers import addLinearSolver
 import dataclasses
 import numpy as np
 
 
 
 def createScene(root):
+    root.gravity=[0,0,9.81]
     ##Solvers
-    setupLagrangianCollision(root, displayFlags = "showVisualModels",backgroundColor=[0.8, 0.8, 0.8, 1],
-                             parallelComputing = True,alarmDistance=0.3, contactDistance=0.02,
-                             frictionCoef=0.5, tolerance=1.0e-4, maxIterations=20)
+    setupDefaultHeader(root, displayFlags = "showVisualModels",backgroundColor=[0.8, 0.8, 0.8, 1],
+                             parallelComputing = True)
+    # setupLagrangianCollision(root, displayFlags = "showVisualModels",backgroundColor=[0.8, 0.8, 0.8, 1],
+    #                          parallelComputing = True,alarmDistance=0.3, contactDistance=0.02,
+    #                          frictionCoef=0.5, tolerance=1.0e-4, maxIterations=20)
 
     ##Environement
     planes_lengthNormal = np.array([0, 1, 0])
@@ -29,7 +37,8 @@ def createScene(root):
     plane1_collisionParams.geometry = PlaneParameters(np.array([15,0,1]), np.array([0,0,-1]),
                                                       planes_lengthNormal, planes_lengthNbEdge, planes_widthNbEdge, planes_lengthSize, planes_widthSize)
     plane1 = root.add(Collision, parameters = plane1_collisionParams)
-    #TODO being able to reuse already loaded geometry of current prefab to add any new sub prefab
+    # TODO being able to reuse already loaded geometry of current prefab to add any new sub prefab
+    # We need to enable to directly pass a link to an already existing prefab in place of a prefab parameter object
     plane1_visu = plane1.addChild("Visu")
     plane1_visu.addObject("OglModel", name="VisualModel", src="@../Geometry/container")
 
@@ -46,3 +55,68 @@ def createScene(root):
 
 
     ## Real models
+    Beam = root.addChild("Beam")
+
+    VolumetricObjects = root.addChild("VolumetricObjects")
+    addImplicitODE(VolumetricObjects)
+    addLinearSolver(VolumetricObjects, constantSparsity=True)
+
+    ### Logo
+    LogoParams = EntityParameters()
+    LogoParams.name = "Logo"
+    LogoParams.geometry = FileParameters(filename="mesh/SofaScene/Logo.vtk")
+    LogoParams.geometry.elementType = ElementType.TETRAHEDRA
+    LogoParams.material = DeformableBehaviorParameters()
+    LogoParams.material.constitutiveLawType = ConstitutiveLaw.ELASTIC
+    LogoParams.material.parameters = [200, 0.4]
+
+    def logoAddMaterial(node):
+        DeformableBehaviorParameters.addDeformableMaterial(node)
+        node.addObject("ConstantForceField", name="ConstantForceUpwards", totalForce=[0, 0, -5.0])
+        #TODO deal with that is a more smooth way in the material directly
+        node.addObject("LinearSolverConstraintCorrection", name="ConstraintCorrection", linearSolver=VolumetricObjects.LinearSolver.linkpath, ODESolver=VolumetricObjects.ODESolver.linkpath)
+
+
+    LogoParams.material.addMaterial = logoAddMaterial
+    LogoParams.material.massDensity = 0.003261
+    LogoParams.collision = CollisionParameters()
+    LogoParams.collision.primitives = [CollisionPrimitive.SPHERES]
+    LogoParams.collision.geometry = FileParameters(filename="mesh/SofaScene/LogoColli.sph")
+    #TODO make this flawless with spheres. Here collisions elements are not in the topology and a link is to be made between the loader and the collision object
+    LogoParams.collision.kwargs = {"SphereCollision" : {"radius" : "@Geometry/loader.listRadius"}}
+    LogoParams.visual = VisualParameters()
+    LogoParams.visual.geometry = FileParameters(filename="mesh/SofaScene/LogoVisu.obj")
+    LogoParams.visual.color = [0.7, .35, 0, 0.8]
+
+    Logo = VolumetricObjects.add(Entity, parameters = LogoParams)
+
+    ### S
+    SParams = EntityParameters()
+    SParams.name = "S"
+    SParams.geometry = FileParameters(filename="mesh/SofaScene/S.vtk")
+    SParams.geometry.elementType = ElementType.TETRAHEDRA
+    SParams.material = DeformableBehaviorParameters()
+    SParams.material.constitutiveLawType = ConstitutiveLaw.ELASTIC
+    SParams.material.parameters = [200, 0.45]
+
+    def SAddMaterial(node):
+        DeformableBehaviorParameters.addDeformableMaterial(node)
+        #TODO deal with that is a more smooth way in the material directly
+        node.addObject("LinearSolverConstraintCorrection", name="ConstraintCorrection", linearSolver=VolumetricObjects.LinearSolver.linkpath, ODESolver=VolumetricObjects.ODESolver.linkpath)
+
+    SParams.material.addMaterial = SAddMaterial
+    SParams.material.massDensity = 0.011021
+    SParams.collision = CollisionParameters()
+    SParams.collision.primitives = [CollisionPrimitive.TRIANGLES]
+    # #TODO: to fix link issues for extracted geometry, it might be better to give source geometry relative link + parameters
+    SParams.collision.geometry = ExtractParameters(destinationType=ElementType.TRIANGLES, sourceParameters=SParams.geometry )
+    SParams.visual = VisualParameters()
+    SParams.visual.geometry = FileParameters(filename="mesh/SofaScene/SVisu.obj")
+    SParams.visual.color = [0.7, .7, 0.7, 0.8]
+
+    S = VolumetricObjects.add(Entity, parameters = SParams)
+
+
+    SDensity = 0.011021
+    ODensity = SDensity
+    ADensity = 0.00693695

splib/core/node_wrapper.py

@@ -25,7 +25,9 @@ def addObject(self,*args, **kwargs):
             parameters["name"] = kwargs["name"]
             if kwargs["name"] in kwargs:
                 if isinstance(kwargs[kwargs["name"]], dict):
-                    parameters = {**parameters, **kwargs[kwargs["name"]]}
+                    for param in kwargs[kwargs["name"]]:
+                        if not(isinstance(kwargs[kwargs["name"]][param],defaultValueType)):
+                            parameters = {**parameters, param : kwargs[kwargs["name"]][param]}
                 else:
                     print("[Warning] You are passing a keyword arg with the same name as one obj without it being a Dict, it will not be used. ")
 

splib/mechanics/collision_model.py

@@ -24,7 +24,7 @@ def addCollisionModels(node, primitive : CollisionPrimitive,
             node.addObject("TriangleCollisionModel", name="TriangleCollision", topology=topology, selfCollision=selfCollision, proximity=proximity, contactStiffness=contactStiffness, contactFriction=contactFriction, group=group,**kwargs)
             return
         case CollisionPrimitive.SPHERES:
-            node.addObject("SphereCollisionModel", name="SphereCollision", topology=topology, selfCollision=selfCollision, proximity=proximity, contactStiffness=contactStiffness, contactFriction=contactFriction, group=group, radius=spheresRadius, **kwargs)
+            node.addObject("SphereCollisionModel", name="SphereCollision", selfCollision=selfCollision, proximity=proximity, contactStiffness=contactStiffness, contactFriction=contactFriction, group=group, radius=spheresRadius, **kwargs)
             return
         case _:
             return

stlib/entities/__entity__.py

@@ -65,16 +65,17 @@ def init(self):
     def addMapping(self, destinationPrefab):
 
         template = f'{self.parameters.stateType},Vec3' # TODO: check that it is always true
-        
+
+        #TODO: all paths are expecting Geometry to be called Geomtry and so on. We need to robustify this by using the name parameter somehow
         if( self.parameters.stateType == StateType.VEC3):
             destinationPrefab.addObject("BarycentricMapping", 
                                         output=destinationPrefab.linkpath, 
-                                        output_topology=destinationPrefab.geometry.container.linkpath,
-                                        input=self.material.linkpath, 
-                                        input_topology=self.geometry.container.linkpath, 
+                                        output_topology=destinationPrefab.Geometry.container.linkpath,
+                                        input=self.Material.linkpath,
+                                        input_topology=self.Geometry.container.linkpath,
                                         template=template)
         else:
             destinationPrefab.addObject("RigidMapping", 
                                          output=destinationPrefab.linkpath, 
-                                        input=self.material.linkpath, 
+                                        input=self.Material.linkpath,
                                         template=template)

stlib/geometries/__geometry__.py

@@ -56,9 +56,17 @@ def init(self):
         # Generate attribute (positions, edges, triangles, quads, tetrahedra, hexahedra) from the internal data provider
         if isinstance(self.parameters.data, InternalDataProvider) :
             self.parameters.data.generateAttribute(self)
+
         if self.parameters.dynamicTopology :
             if self.parameters.elementType is not None :
-                addDynamicTopology(self, container = dataclasses.asdict(self.parameters.data))
+                addDynamicTopology(self, elementType=self.parameters.elementType, container = {
+                                                                                                "position": self.parameters.data.position,
+                                                                                                "edges": self.parameters.data.edges,
+                                                                                                "triangles": self.parameters.data.triangles,
+                                                                                                "quads": self.parameters.data.quads,
+                                                                                                "tetrahedra": self.parameters.data.tetrahedra,
+                                                                                                "hexahedra": self.parameters.data.hexahedra
+                                                                                            })
             else:
                 raise ValueError
         else:

stlib/geometries/extract.py

@@ -25,7 +25,7 @@ def __post_init__(self):
 
         InternalDataProvider.__init__(self)
 
-    def generateAttribute(self, parent : Geometry):                
+    def generateAttribute(self, parent : Geometry):
         node = parent.addChild("ExtractedGeometry")
 
         #TODO: Specify somewhere in the doc that this should only be used for mapped topologies that extract parent topology surface
@@ -34,8 +34,8 @@ def generateAttribute(self, parent : Geometry):
         # !!! also, on a fail, nothing is added to the graph, which makes things harder to debug
         # !!! also, does not work because of the function canCreate(), which checks the input (not yet created?)
         # this is all related
-        fromLink = "@../../Geometry.container" # TODO: can we do better than this?
-        addDynamicTopology(node, elementType=self.sourceType)
+        fromLink = "@../../../Geometry/container" # TODO: can we do better than this?
+        addDynamicTopology(node, elementType=self.destinationType)
         if self.sourceType == ElementType.TETRAHEDRA:
             node.addObject("Tetra2TriangleTopologicalMapping", input=fromLink, output=node.container.linkpath)
         elif self.sourceType == ElementType.HEXAHEDRA:
@@ -60,12 +60,11 @@ def generateAttribute(self, parent : Geometry):
 class ExtractParameters(GeometryParameters):
     def __init__(self, 
                  sourceParameters : GeometryParameters, 
-                 destinationType : ElementType,  
-                 dynamicTopology : bool = False, ):
+                 destinationType : ElementType,):
         GeometryParameters.__init__(self,
                                     data = ExtractInternalDataProvider(destinationType = destinationType, 
                                                                        sourceType = sourceParameters.elementType,
                                                                        sourceName = sourceParameters.name), 
-                                    dynamicTopology = dynamicTopology, 
+                                    dynamicTopology = True,
                                     elementType = destinationType)
 

stlib/geometries/file.py

@@ -13,14 +13,20 @@ def __post_init__(self, **kwargs):
         InternalDataProvider.__init__(self,**kwargs)
 
     def generateAttribute(self, parent : Geometry):    
-        loadMesh(parent, self.filename)      
-
-        self.position = str(parent.loader.position.linkpath)
-        self.edges = str(parent.loader.edges.linkpath)
-        self.triangles = str(parent.loader.triangles.linkpath)
-        self.quads = str(parent.loader.quads.linkpath)
-        self.hexahedra = str(parent.loader.hexahedra.linkpath)
-        self.tetrahedra = str(parent.loader.tetras.linkpath)
+        loadMesh(parent, self.filename)
+
+        if hasattr(parent.loader, 'position'):
+            self.position = str(parent.loader.position.linkpath)
+        if hasattr(parent.loader, 'edges'):
+            self.edges = str(parent.loader.edges.linkpath)
+        if hasattr(parent.loader, 'triangles'):
+            self.triangles = str(parent.loader.triangles.linkpath)
+        if hasattr(parent.loader, 'quads'):
+            self.quads = str(parent.loader.quads.linkpath)
+        if hasattr(parent.loader, 'hexahedra'):
+            self.hexahedra = str(parent.loader.hexahedra.linkpath)
+        if hasattr(parent.loader, 'tetrahedra'):
+            self.tetrahedra = str(parent.loader.tetrahedra.linkpath)
 
 
 

stlib/materials/deformable.py

@@ -13,7 +13,7 @@ class DeformableBehaviorParameters(MaterialParameters):
     elementType : ElementType = ElementType.TETRAHEDRA
     parameters : list[float] = dataclasses.field(default_factory=lambda: [1000, 0.45])  # young modulus, poisson ratio
 
-    def __addDeformableMaterial(node):
+    def addDeformableMaterial(node):
 
         massKwargs = {}
         if node.parameters.elementType != ElementType.EDGES: #If we use the MeshMatrixMass, then the mass will need us to specify the mstate to use
@@ -26,7 +26,7 @@ def __addDeformableMaterial(node):
         else:
             addLinearElasticity(node, node.parameters.elementType, node.parameters.parameters[0], node.parameters.parameters[1], topology="@../Geometry/container")
 
-    addMaterial : Optional[Callable] = __addDeformableMaterial
+    addMaterial : Optional[Callable] = addDeformableMaterial
 
 
 def createScene(root) :