diff --git a/examples/stlib/PrefabScene_beginner.py b/examples/stlib/PrefabScene_beginner.py new file mode 100644 index 000000000..2c4552c59 --- /dev/null +++ b/examples/stlib/PrefabScene_beginner.py @@ -0,0 +1,30 @@ +from stlib.materials.rigid import Rigid +from stlib.materials.deformable import Deformable +from stlib.geometries.cube import CubeParameters +from stlib.geometries.file import FileParameters +from splib.simulation.headers import setupLagrangianCollision +from splib.simulation.linear_solvers import addLinearSolver +from splib.simulation.ode_solvers import addImplicitODE + +#To be added in splib +def addSolvers(root): + addLinearSolver(root) + addImplicitODE(root) + root.addObject("LinearSolverConstraintCorrection", linearsolver="@LinearSolver") + +def createScene(root): + root.gravity = [0, 0, -9.81] + + setupLagrangianCollision(root) + addSolvers(root) + + rigidParams = Rigid.getParameters() + rigidParams.geometry = CubeParameters([0, 0, 0], 1, 3) + root.add(Rigid,rigidParams) + + deformableParams = Deformable.getParameters() + #Add transformation somewhere here + deformableParams.geometry = FileParameters("SofaScene/Logo.vtk") + root.add(Deformable,deformableParams) + + return root diff --git a/examples/stlib/PrefabScene_expert.py b/examples/stlib/PrefabScene_expert.py new file mode 100644 index 000000000..0bc427d37 --- /dev/null +++ b/examples/stlib/PrefabScene_expert.py @@ -0,0 +1,14 @@ +from splib.helper import exportScene +from stlib.misc.entity import Entity + + +def createScene(root): + Entity.Deformable(Mesh,MechaProperties) + + return root + + +if __name__=="__main__": + root = exportScene() + createScene(root) + pass \ No newline at end of file diff --git a/examples/stlib/PrefabScene_intermediate.py b/examples/stlib/PrefabScene_intermediate.py new file mode 100644 index 000000000..7ce22c745 --- /dev/null +++ b/examples/stlib/PrefabScene_intermediate.py @@ -0,0 +1,8 @@ +from splib.helper import exportScene +from stlib.misc.entity import Entity + + +def createScene(root): + + + return root \ No newline at end of file diff --git a/splib/Testing.py b/splib/Testing.py new file mode 100644 index 000000000..146b4096b --- /dev/null +++ b/splib/Testing.py @@ -0,0 +1,86 @@ +from splib.topology.dynamic import * +from splib.simulation.headers import * +from splib.simulation.ode_solvers import * +from splib.simulation.linear_solvers import * +from splib.mechanics.linear_elasticity import * +from splib.mechanics.mass import * +from splib.mechanics.fix_points import * +from splib.topology.loader import * +from splib.core.node_wrapper import * + +from splib.core.node_wrapper import ChildWrapper, ObjectWrapper +from splib.helper import exportScene + + +def createScene(rootNode): + rootNode.dt = 0.03 + rootNode.gravity = [0,-9.81,0] + + setupLagrangianCollision(rootNode,requiredPlugins={"pluginName":['Sofa.Component.Constraint.Projective', + 'Sofa.Component.Engine.Select', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.Mass', + 'Sofa.Component.ODESolver.Backward', + 'Sofa.Component.SolidMechanics.FEM.Elastic', + 'Sofa.Component.StateContainer', + 'Sofa.Component.Topology.Container.Grid', + 'Sofa.Component.IO.Mesh', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.Topology.Container.Dynamic', + 'Sofa.Component.Visual']}) + # + # + # ## TODO : Being able to call "childNode.addAnything" by using the __getattr__ method + Liver0=rootNode.addChild("Liver0") + Liver0.addObject("EulerImplicitSolver",name="ODESolver",rayleighStiffness="0",rayleighMass="0") + Liver0.addObject("SparseLDLSolver",name="LinearSolver",template="CompressedRowSparseMatrixMat3x3",parallelInverseProduct="False") + Liver0.addObject("MeshGmshLoader",name="meshLoader",filename="mesh/liver.msh") + Liver0.addObject("TetrahedronSetTopologyModifier",name="modifier") + Liver0.addObject("TetrahedronSetTopologyContainer",name="container",src="@meshLoader") + Liver0.addObject("TetrahedronSetGeometryAlgorithms",name="algorithms") + Liver0.addObject("MechanicalObject",name="mstate",template="Vec3d") + Liver0.addObject("LinearSolverConstraintCorrection",name="constraintCorrection") + Liver0.addObject("TetrahedronFEMForceField",name="constitutiveLaw",youngModulus="3000",poissonRatio="0.3",method="large") + Liver0.addObject("MeshMatrixMass",name="mass",massDensity="1") + Liver0.addObject("BoxROI",name="fixedBoxROI",box="0 3 0 2 5 2") + Liver0.addObject("FixedProjectiveConstraint",name="fixedConstraints",indices="@fixedBoxROI.indices") + Visu=Liver0.addChild("Visu") + Visu.addObject("TriangleSetTopologyModifier",name="modifier") + Visu.addObject("TriangleSetTopologyContainer",name="container") + Visu.addObject("TriangleSetGeometryAlgorithms",name="algorithms") + Visu.addObject("Tetra2TriangleTopologicalMapping",name="TopologicalMapping",input="@../container",output="@container") + Visu.addObject("OglModel",name="OglModel",topology="@container",color="[1.0, 0.2, 0.8]") + Visu.addObject("IdentityMapping",name="Mapping",isMechanical="False") + + + SimulatedLiver1 = rootNode.addChild("Liver1") + addImplicitODE(SimulatedLiver1) + addLinearSolver(SimulatedLiver1,iterative=False, template="CompressedRowSparseMatrixMat3x3") + loadMesh(SimulatedLiver1,filename="mesh/liver.msh") + addDynamicTopology(SimulatedLiver1,type=ElementType.TETRAHEDRA,source="@meshLoader") + SimulatedLiver1.addObject("MechanicalObject",name="mstate", template='Vec3d') + SimulatedLiver1.addObject("LinearSolverConstraintCorrection",name="constraintCorrection") + addLinearElasticity(SimulatedLiver1,ElementType.TETRAHEDRA, poissonRatio="0.3", youngModulus="3000", method='large') + addMass(SimulatedLiver1,template='Vec3d',massDensity="2") + addFixation(SimulatedLiver1,ConstraintType.PROJECTIVE,boxROIs=[0, 3, 0, 2, 5, 2]) + + SimulatedLiverVisu = SimulatedLiver1.addChild("Visu") + addDynamicTopology(SimulatedLiverVisu,ElementType.TRIANGLES) + SimulatedLiverVisu.addObject("Tetra2TriangleTopologicalMapping", name="TopologicalMapping", input="@../container", output="@container") + SimulatedLiverVisu.addObject("OglModel", name="OglModel", topology="@container",color=[1.0,0.2,0.8]) + SimulatedLiverVisu.addObject("IdentityMapping",name="Mapping",isMechanical=False) + + + return rootNode + + + +if __name__ == "__main__": + Node = exportScene() + createScene(Node) + + + + + + diff --git a/splib/__init__.py b/splib/__init__.py index 34e9c316c..6b380c11a 100644 --- a/splib/__init__.py +++ b/splib/__init__.py @@ -1,4 +1,2 @@ -"""SPLIB is now relocated at: https://github.com/SofaDefrost/STLIB, please clone and install the plugin to use it""" -raise Exception("SPLIB is now relocated at: https://github.com/SofaDefrost/STLIB, please clone and install the plugin to use it") - +__all__ = ["core","topology", "simulation","modeler","mechanics","Testing","helper"] diff --git a/splib/core/__init__.py b/splib/core/__init__.py new file mode 100644 index 000000000..78deb3fad --- /dev/null +++ b/splib/core/__init__.py @@ -0,0 +1 @@ +__all__ = ["node_wrapper", "utils","enum_types"] diff --git a/splib/core/enum_types.py b/splib/core/enum_types.py new file mode 100644 index 000000000..5d4504311 --- /dev/null +++ b/splib/core/enum_types.py @@ -0,0 +1,51 @@ +from enum import Enum + +class ConstitutiveLaw(Enum): + ELASTIC = 1 + HYPERELASTIC = 2 + +class ODEType(Enum): + EXPLICIT = 1 + IMPLICIT = 2 + +class SolverType(Enum): + DIRECT = 1 + ITERATIVE = 2 + +class MappingType(Enum): + BARYCENTRIC = 1 + IDENTITY = 2 + RIGID = 3 + +class ConstraintType(Enum): + PROJECTIVE = 1 + WEAK = 2 + LAGRANGIAN = 3 + +class CollisionPrimitive(Enum): + POINTS = 1 + LINES = 2 + TRIANGLES = 3 + SPHERES = 4 + +class ElementType(Enum): + POINTS = 1 + EDGES = 2 + TRIANGLES = 3 + QUADS = 4 + TETRAHEDRA = 5 + HEXAHEDRA = 6 + +class StateType(Enum): + VEC3 = 3 + VEC1 = 1 + RIGID = 7 + + def __str__(self): + if self == StateType.VEC3: + return "Vec3" + if self == StateType.VEC1: + return "Vec1" + if self == StateType.RIGID: + return "Rigid3" + return 'Unknown' diff --git a/splib/core/node_wrapper.py b/splib/core/node_wrapper.py new file mode 100644 index 000000000..7646b9e84 --- /dev/null +++ b/splib/core/node_wrapper.py @@ -0,0 +1,73 @@ +from functools import wraps +from splib.core.utils import defaultValueType +# The two classes are not merged because one could want to use a ReusableMethod +# (enabling to pass dictionary fixing params) without wanting to use a full PrefabSimulation + + +class BaseWrapper(object): + + def __init__(self,node): + self.node = node + + def __getattr__(self, item): + return getattr(self.node,item) + + +class ObjectWrapper(BaseWrapper): + def __init__(self,*args,**kwargs): + super().__init__(*args,**kwargs) + + def addObject(self,*args, **kwargs): + parameters = {} + # Expand any parameter that has been set by the user in a custom dictionary + # and having the same key as the component name + if "name" in kwargs: + parameters["name"] = kwargs["name"] + if kwargs["name"] in kwargs: + if isinstance(kwargs[kwargs["name"]], dict): + parameters = {**parameters, **kwargs[kwargs["name"]]} + 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. ") + + # Add all the parameters from kwargs that are not dictionary + for param in kwargs: + if param in parameters: + if not(param == "name"): + print("[Warning] You are redefining the parameter '"+ param + "' of object " + str(args[0])) + elif not(isinstance(kwargs[param], dict)) and not(isinstance(kwargs[param],defaultValueType)): + parameters = {**parameters,param:kwargs[param]} + + return self.node.addObject(*args,**parameters) + + + + +class ChildWrapper(ObjectWrapper): + def __init__(self,*args,**kwargs): + super().__init__(*args,**kwargs) + + # This method enforces that the object that is created by the addChild method keeps the prefab type + def addChild(self,*args, **kwargs): + child = self.node.addChild(*args,**kwargs) + returnObject = self.__new__(type(self)) + returnObject.__init__(child) + return returnObject + + + + +def ReusableMethod(method): + @wraps(method) + def wrapper(*args, **kwargs): + + node = args[0] + # We don't want to wrap an object that is already wrapped + # It wouldn't break anything but nodes might get wrapped and wrapped again multiplying redirections... + if( not isinstance(node,ObjectWrapper) ): + node = ObjectWrapper(node) + + if len(args)>1: + return method(node,*args[1:],**kwargs) + else: + return method(node,**kwargs) + return wrapper diff --git a/splib/core/utils.py b/splib/core/utils.py new file mode 100644 index 000000000..56b91f53c --- /dev/null +++ b/splib/core/utils.py @@ -0,0 +1,36 @@ +from typing import List, Callable, Tuple, Dict +from functools import wraps + +class defaultValueType(): + def __init__(self): + pass + +DEFAULT_VALUE = defaultValueType() + +def isDefault(obj): + return isinstance(obj,defaultValueType) + + +def getParameterSet(name : str,parameterSet : Dict) -> Dict: + if name in parameterSet: + if isinstance(parameterSet[name], dict): + return parameterSet[name] + return {} + + +def MapKeywordArg(objectName,*argumentMaps): + def MapArg(method): + @wraps(method) + def wrapper(*args, **kwargs): + containerParams = getParameterSet(objectName,kwargs) + for argMap in argumentMaps: + if (argMap[0] in kwargs) and not(kwargs[argMap[0]] is None): + containerParams[argMap[1]] = kwargs[argMap[0]] + kwargs[objectName] = containerParams + return method(*args,**kwargs) + return wrapper + return MapArg + + + + diff --git a/splib/helper.py b/splib/helper.py new file mode 100644 index 000000000..b4881af5f --- /dev/null +++ b/splib/helper.py @@ -0,0 +1,35 @@ +class displayNode(): + def __init__(self,_level=0): + self.prefix = "" + for i in range(_level): + self.prefix += "| " + + def addObject(self,type:str,**kwargs): + print(self.prefix + type + " with " + str(kwargs)) + + def addChild(self,name:str): + print(self.prefix + "-> Node : " + name) + return displayNode(len(self.prefix) + 1) + + +class exportScene(): + def __init__(self,name="rootNode"): + self.name = name + + def addObject(self,type:str,**kwargs): + suffix = "" + for i in kwargs: + suffix += "," + str(i) + "=\"" + str(kwargs[i]) + "\"" + print(self.name+".addObject(\"" + type + "\"" + suffix + ")") + + def addChild(self,name:str): + print(name + '=' + self.name+".addChild(\"" + name + "\")") + setattr(self,name,exportScene(name)) + return getattr(self,name) + + def __setattr__(self, key, value): + if(not(key == "name")): + print(self.__dict__["name"] + "." + key + " = " + str(value)) + self.__dict__[key] = value + else: + self.__dict__[key] = value diff --git a/splib/mechanics/__init__.py b/splib/mechanics/__init__.py new file mode 100644 index 000000000..3761021a2 --- /dev/null +++ b/splib/mechanics/__init__.py @@ -0,0 +1 @@ +__all__ = ["linear_elasticity","hyperelasticity","fix_points","collision_model","mass"] \ No newline at end of file diff --git a/splib/mechanics/collision_model.py b/splib/mechanics/collision_model.py new file mode 100644 index 000000000..e1abdf457 --- /dev/null +++ b/splib/mechanics/collision_model.py @@ -0,0 +1,34 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import DEFAULT_VALUE +from splib.core.enum_types import CollisionPrimitive + + +@ReusableMethod +def addCollisionModels(node, primitive : CollisionPrimitive, + topology=DEFAULT_VALUE, + selfCollision=DEFAULT_VALUE, + proximity=DEFAULT_VALUE, + group=DEFAULT_VALUE, + contactStiffness=DEFAULT_VALUE, + contactFriction=DEFAULT_VALUE, + spheresRadius=DEFAULT_VALUE, + **kwargs): + match primitive: + case CollisionPrimitive.POINTS: + node.addObject("PointCollisionModel", name="PointCollision", topology=topology, selfCollision=selfCollision, proximity=proximity, contactStiffness=contactStiffness, contactFriction=contactFriction, group=group, **kwargs) + return + case CollisionPrimitive.LINES: + node.addObject("LineCollisionModel", name="EdgeCollision", topology=topology, selfCollision=selfCollision, proximity=proximity, contactStiffness=contactStiffness, contactFriction=contactFriction, group=group, **kwargs) + return + case CollisionPrimitive.TRIANGLES: + 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) + return + case _: + return + + #Cube and tetra are missing. + + \ No newline at end of file diff --git a/splib/mechanics/fix_points.py b/splib/mechanics/fix_points.py new file mode 100644 index 000000000..ddbbd6308 --- /dev/null +++ b/splib/mechanics/fix_points.py @@ -0,0 +1,30 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import isDefault, DEFAULT_VALUE +from splib.core.enum_types import ConstraintType +from enum import Enum + + +##box +@ReusableMethod +def addFixation(node,type:ConstraintType,boxROIs=DEFAULT_VALUE, sphereROIs=DEFAULT_VALUE, indices=DEFAULT_VALUE, fixAll=DEFAULT_VALUE,**kwargs): + if (isDefault(indices)): + if(not isDefault(boxROIs)): + node.addObject("BoxROI",name='fixedBoxROI',box=boxROIs,**kwargs) + indices="@fixedBoxROI.indices" + if(not isDefault(sphereROIs)): + node.addObject("SphereROI",name='fixedSphereROI',centers=sphereROIs[0],radii=sphereROIs[1],**kwargs) + indices="@fixedSphereROI.indices" + + match type: + case ConstraintType.WEAK: + node.addObject("FixedWeakConstraint",name="fixedConstraints", indices=indices, fixAll=fixAll, **kwargs) + return + case ConstraintType.PROJECTIVE: + node.addObject("FixedProjectiveConstraint",name="fixedConstraints", indices=indices, fixAll=fixAll, **kwargs) + return + case ConstraintType.LAGRANGIAN: + node.addObject("LagrangianFixedConstraint",name="fixedConstraints", indices=indices, fixAll=fixAll, **kwargs) + return + case _: + print('Contraint type is either ConstraintType.PROJECTIVE, ConstraintType.WEAK or ConstraintType.LAGRANGIAN') + return diff --git a/splib/mechanics/hyperelasticity.py b/splib/mechanics/hyperelasticity.py new file mode 100644 index 000000000..fb255df41 --- /dev/null +++ b/splib/mechanics/hyperelasticity.py @@ -0,0 +1,14 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import DEFAULT_VALUE +from splib.core.enum_types import ElementType + + +@ReusableMethod +def addHyperelasticity(node,elem:ElementType,materialName=DEFAULT_VALUE, parameterSet=DEFAULT_VALUE, matrixRegularization=DEFAULT_VALUE,**kwargs): + match elem: + case ElementType.TETRAHEDRA: + node.addObject("TetrahedronHyperelasticityFEMForceField",name="constitutiveLaw", materialName=materialName, parameterSet=parameterSet, matrixRegularization=matrixRegularization, **kwargs) + return + case _: + print('Hyperelasticity model only exist for Tetrahedron elements.') + return \ No newline at end of file diff --git a/splib/mechanics/linear_elasticity.py b/splib/mechanics/linear_elasticity.py new file mode 100644 index 000000000..1e9b30290 --- /dev/null +++ b/splib/mechanics/linear_elasticity.py @@ -0,0 +1,26 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import DEFAULT_VALUE +from splib.core.enum_types import ElementType + + +@ReusableMethod +def addLinearElasticity(node, elementType:ElementType, youngModulus=DEFAULT_VALUE, poissonRatio=DEFAULT_VALUE, method=DEFAULT_VALUE, **kwargs): + match elementType: + case ElementType.EDGES: + node.addObject("BeamFEMForceField", name="constitutiveLaw", youngModulus=youngModulus, poissonRatio=poissonRatio, method=method, **kwargs) + return + case ElementType.TRIANGLES: + node.addObject("TriangleFEMForceField", name="constitutiveLaw", youngModulus=youngModulus, poissonRatio=poissonRatio, method=method, **kwargs) + return + case ElementType.QUADS: + node.addObject("QuadBendingFEMForceField", name="constitutiveLaw", youngModulus=youngModulus, poissonRatio=poissonRatio, method=method, **kwargs) + return + case ElementType.TETRAHEDRA: + node.addObject("TetrahedronFEMForceField", name="constitutiveLaw", youngModulus=youngModulus, poissonRatio=poissonRatio, method=method, **kwargs) + return + case ElementType.HEXAHEDRA: + node.addObject("HexahedronFEMForceField", name="constitutiveLaw", youngModulus=youngModulus, poissonRatio=poissonRatio, method=method, **kwargs) + return + case _: + print('Linear elasticity is only available for topology of type EDGES, TRIANGLES, QUADS, TETRAHEDRON, HEXAHEDRON') + return \ No newline at end of file diff --git a/splib/mechanics/mass.py b/splib/mechanics/mass.py new file mode 100644 index 000000000..c003c3e64 --- /dev/null +++ b/splib/mechanics/mass.py @@ -0,0 +1,23 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import defaultValueType, DEFAULT_VALUE, isDefault +from splib.core.enum_types import ElementType + + +# TODO : use the massDensity ONLY and deduce totalMass if necessary from it + volume + +@ReusableMethod +def addMass(node, elem:ElementType, totalMass=DEFAULT_VALUE, massDensity=DEFAULT_VALUE, lumping=DEFAULT_VALUE, topology=DEFAULT_VALUE, **kwargs): + if (not isDefault(totalMass)) and (not isDefault(massDensity)) : + print("[warning] You defined the totalMass and the massDensity in the same time, only taking massDensity into account") + del kwargs["massDensity"] + + if(elem !=ElementType.POINTS and elem !=ElementType.EDGES): + node.addObject("MeshMatrixMass",name="mass", totalMass=totalMass, massDensity=massDensity, lumping=lumping, topology=topology, **kwargs) + else: + if (not isDefault(massDensity)) : + print("[warning] mass density can only be used on a surface or volumetric topology. Please use totalMass instead") + if (not isDefault(lumping)) : + print("[warning] lumping can only be set for surface or volumetric topology") + + node.addObject("UniformMass",name="mass", totalMass=totalMass, topology=topology,**kwargs) + diff --git a/splib/simulation/__init__.py b/splib/simulation/__init__.py new file mode 100644 index 000000000..d59a665ef --- /dev/null +++ b/splib/simulation/__init__.py @@ -0,0 +1 @@ +__all__ = ["headers","linear_solvers","ode_solvers"] diff --git a/splib/simulation/headers.py b/splib/simulation/headers.py new file mode 100644 index 000000000..cf2f5163d --- /dev/null +++ b/splib/simulation/headers.py @@ -0,0 +1,128 @@ +from splib.core.node_wrapper import ReusableMethod +from enum import Enum + +from splib.core.utils import DEFAULT_VALUE + + +class CollisionType(Enum): + NONE = 1 + PENALITY = 2 + LAGRANGIAN = 3 + + +@ReusableMethod +def setupDefaultHeader(node, displayFlags = "showVisualModels", backgroundColor=[1,1,1,1], parallelComputing=False,**kwargs): + + node.addObject('VisualStyle', displayFlags=displayFlags) + node.addObject('BackgroundSetting', color=backgroundColor) + + node.addObject("RequiredPlugin", name="requiredPlugins", pluginName=['Sofa.Component.Constraint.Projective', + 'Sofa.Component.Engine.Select', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.Mass', + 'Sofa.Component.ODESolver.Backward', + 'Sofa.Component.SolidMechanics.FEM.Elastic', + 'Sofa.Component.StateContainer', + 'Sofa.Component.Topology.Container.Grid', + 'Sofa.Component.IO.Mesh', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.ODESolver.Forward', + 'Sofa.Component.Topology.Container.Dynamic', + 'Sofa.Component.Visual', + ], + **kwargs) + node.addObject('DefaultAnimationLoop',name="animation", parallelODESolving=parallelComputing, **kwargs) + + return node + + +@ReusableMethod +def setupPenalityCollisionHeader(node, displayFlags = "showVisualModels",backgroundColor=[1,1,1,1], stick=False, parallelComputing=False, alarmDistance=DEFAULT_VALUE, contactDistance=DEFAULT_VALUE, **kwargs): + node.addObject('VisualStyle', displayFlags=displayFlags) + node.addObject('BackgroundSetting', color=backgroundColor) + + node.addObject("RequiredPlugin", name="requiredPlugins", pluginName=['Sofa.Component.Constraint.Projective', + 'Sofa.Component.Engine.Select', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.Mass', + 'Sofa.Component.ODESolver.Backward', + 'Sofa.Component.SolidMechanics.FEM.Elastic', + 'Sofa.Component.StateContainer', + 'Sofa.Component.Topology.Container.Grid', + 'Sofa.Component.IO.Mesh', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.ODESolver.Forward', + 'Sofa.Component.Topology.Container.Dynamic', + 'Sofa.Component.Visual', + ], + **kwargs) + + parallelPrefix = "" + if(parallelComputing): + parallelPrefix="Parallel" + + node.addObject('DefaultAnimationLoop',name="animation", **kwargs) + node.addObject('CollisionPipeline', name="collisionPipeline", **kwargs) + node.addObject(parallelPrefix+'BruteForceBroadPhase', name="broadPhase", **kwargs) + node.addObject(parallelPrefix+'BVHNarrowPhase', name="narrowPhase", **kwargs) + + if(stick): + node.addObject('CollisionResponse',name="ContactManager", response="BarycentricStickContact",**kwargs) + else: + node.addObject('CollisionResponse',name="ContactManager", response="BarycentricPenalityContact",**kwargs) + node.addObject('LocalMinDistance' ,name="Distance", alarmDistance=alarmDistance, contactDistance=contactDistance, **kwargs) + + return node + + +@ReusableMethod +def setupLagrangianCollision(node, displayFlags = "showVisualModels",backgroundColor=[1,1,1,1], parallelComputing=False, stick=False, alarmDistance=DEFAULT_VALUE, contactDistance=DEFAULT_VALUE, frictionCoef=0.0, tolerance=0.0, maxIterations=100, **kwargs): + node.addObject('VisualStyle', displayFlags=displayFlags) + node.addObject('BackgroundSetting', color=backgroundColor) + + node.addObject("RequiredPlugin", name="requiredPlugins", pluginName=['Sofa.Component.Constraint.Lagrangian', + 'Sofa.Component.Constraint.Projective', + 'Sofa.Component.Engine.Select', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.Mass', + 'Sofa.Component.ODESolver.Backward', + 'Sofa.Component.SolidMechanics.FEM.Elastic', + 'Sofa.Component.StateContainer', + 'Sofa.Component.Topology.Container.Grid', + 'Sofa.Component.IO.Mesh', + 'Sofa.Component.LinearSolver.Direct', + 'Sofa.Component.ODESolver.Forward', + 'Sofa.Component.Topology.Container.Dynamic', + 'Sofa.Component.Visual', + ], + **kwargs) + + + node.addObject('FreeMotionAnimationLoop',name="animation", + parallelCollisionDetectionAndFreeMotion=parallelComputing, + parallelODESolving=parallelComputing, + **kwargs) + + parallelPrefix = "" + if(parallelComputing): + parallelPrefix="Parallel" + + node.addObject('CollisionPipeline', name="collisionPipeline", + **kwargs) + + node.addObject(parallelPrefix+'BruteForceBroadPhase', name="broadPhase", + **kwargs) + + node.addObject(parallelPrefix+'BVHNarrowPhase', name="narrowPhase", + **kwargs) + + if(stick): + node.addObject('CollisionResponse',name="ContactManager", response="StickContactConstraint", responseParams="tol="+str(tolerance),**kwargs) + else: + node.addObject('CollisionResponse',name="ContactManager", response="FrictionContactConstraint", responseParams="mu="+str(frictionCoef),**kwargs) + + node.addObject('NewProximityIntersection' ,name="Distance", alarmDistance=alarmDistance, contactDistance=contactDistance, **kwargs) + node.addObject('ProjectedGaussSeidelConstraintSolver',name="ConstraintSolver", tolerance=tolerance, maxIterations=maxIterations, multithreading=parallelComputing,**kwargs) + node.addObject("ConstraintAttachButtonSetting") + + return node diff --git a/splib/simulation/linear_solvers.py b/splib/simulation/linear_solvers.py new file mode 100644 index 000000000..b6b90f027 --- /dev/null +++ b/splib/simulation/linear_solvers.py @@ -0,0 +1,31 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import * + +@ReusableMethod +def addLinearSolver(node,iterative=False,iterations=DEFAULT_VALUE,tolerance=DEFAULT_VALUE,threshold=DEFAULT_VALUE,template=DEFAULT_VALUE,constantSparsity=False,parallelInverseProduct=DEFAULT_VALUE,**kwargs): + containerParams = getParameterSet("LinearSolver",kwargs) + if iterative: + if not isDefault(iterations): + containerParams["iterations"] = iterations + if not isDefault(tolerance): + containerParams["tolerance"] = tolerance + if not isDefault(threshold): + containerParams["threshold"] = threshold + else: + if isDefault(template) and not("template" in containerParams): + containerParams["template"] = 'CompressedRowSparseMatrix' + elif not isDefault(template) : + containerParams["template"] = template + kwargs["LinearSolver"] = containerParams + + if(constantSparsity): + node.addObject("ConstantSparsityPatternSystem",name='LinearSystem',**kwargs) + kwargs["LinearSolver"]["template"] = 'CompressedRowSparseMatrix' + kwargs["LinearSolver"]["linearSystem"]="@LinearSystem" + + + + if iterative: + node.addObject('CGLinearSolver', name='LinearSolver', parallelInverseProduct=parallelInverseProduct, **kwargs) + else: + node.addObject('SparseLDLSolver', name='LinearSolver', parallelInverseProduct=parallelInverseProduct, **kwargs) \ No newline at end of file diff --git a/splib/simulation/ode_solvers.py b/splib/simulation/ode_solvers.py new file mode 100644 index 000000000..eed57bab8 --- /dev/null +++ b/splib/simulation/ode_solvers.py @@ -0,0 +1,14 @@ +from splib.core.node_wrapper import ReusableMethod + +@ReusableMethod +def addImplicitODE(node,static=False,**kwargs): + if( not(static) ): + node.addObject("EulerImplicitSolver",name="ODESolver",**kwargs) + else: + node.addObject("StaticSolver",name="ODESolver",**kwargs) + +@ReusableMethod +def addExplicitODE(node,**kwargs): + node.addObject("EulerExplicitSolver",name="ODESolver",**kwargs) + + diff --git a/splib/topology/__init__.py b/splib/topology/__init__.py new file mode 100644 index 000000000..23469f77d --- /dev/null +++ b/splib/topology/__init__.py @@ -0,0 +1 @@ +__all__ = ["static","dynamic","loader"] \ No newline at end of file diff --git a/splib/topology/dynamic.py b/splib/topology/dynamic.py new file mode 100644 index 000000000..1f18be29d --- /dev/null +++ b/splib/topology/dynamic.py @@ -0,0 +1,67 @@ +from splib.core.node_wrapper import ReusableMethod +from enum import Enum +from splib.core.utils import MapKeywordArg, DEFAULT_VALUE +from splib.core.enum_types import ElementType + + + +@ReusableMethod +def addPointTopology(node,position=DEFAULT_VALUE,source=DEFAULT_VALUE, **kwargs): + node.addObject("PointSetTopologyModifier", name="modifier", **kwargs) + node.addObject("PointSetTopologyContainer", name="container", src=source, position=position, **kwargs) + # node.addObject("PointSetGeometryAlgorithms", name="algorithms", **kwargs) + +@ReusableMethod +def addEdgeTopology(node,position=DEFAULT_VALUE,edges=DEFAULT_VALUE,source=DEFAULT_VALUE,**kwargs): + node.addObject("EdgeSetTopologyModifier", name="modifier",**kwargs) + node.addObject("EdgeSetTopologyContainer", name="container", src=source, position=position, edges=edges, **kwargs) + # node.addObject("EdgeSetGeometryAlgorithms", name="algorithms",**kwargs) + +@ReusableMethod +def addTriangleTopology(node,position=DEFAULT_VALUE,edges=DEFAULT_VALUE,triangles=DEFAULT_VALUE,source=DEFAULT_VALUE,**kwargs): + node.addObject("TriangleSetTopologyModifier", name="modifier",**kwargs) + node.addObject("TriangleSetTopologyContainer", name="container", src=source, position=position, edges=edges, triangles=triangles, **kwargs) + # node.addObject("TriangleSetGeometryAlgorithms", name="algorithms",**kwargs) + +@ReusableMethod +def addQuadTopology(node,position=DEFAULT_VALUE,edges=DEFAULT_VALUE,quads=DEFAULT_VALUE,source=DEFAULT_VALUE,**kwargs): + node.addObject("QuadSetTopologyModifier", name="modifier",**kwargs) + node.addObject("QuadSetTopologyContainer", name="container", src=source, position=position, edges=edges, quads=quads, **kwargs) + # node.addObject("QuadSetGeometryAlgorithms", name="algorithms",**kwargs) + +@ReusableMethod +def addTetrahedronTopology(node,position=DEFAULT_VALUE,edges=DEFAULT_VALUE,triangles=DEFAULT_VALUE,tetrahedra=DEFAULT_VALUE,source=DEFAULT_VALUE,**kwargs): + node.addObject("TetrahedronSetTopologyModifier", name="modifier",**kwargs) + node.addObject("TetrahedronSetTopologyContainer", name="container", src=source, position=position, edges=edges, triangles=triangles, tetrahedra=tetrahedra, **kwargs) + # node.addObject("TetrahedronSetGeometryAlgorithms", name="algorithms",**kwargs) + +@ReusableMethod +def addHexahedronTopology(node,position=DEFAULT_VALUE,edges=DEFAULT_VALUE,quads=DEFAULT_VALUE,hexahedra=DEFAULT_VALUE,source=DEFAULT_VALUE,**kwargs): + node.addObject("HexahedronSetTopologyModifier", name="modifier",**kwargs) + node.addObject("HexahedronSetTopologyContainer", name="container", src=source, position=position, edges=edges, quads=quads, hexahedra=hexahedra, **kwargs) + # node.addObject("HexahedronSetGeometryAlgorithms", name="algorithms",**kwargs) + +def addDynamicTopology(node, elementType:ElementType, **kwargs): + match elementType: + case ElementType.POINTS: + addPointTopology(node,**kwargs) + return + case ElementType.EDGES: + addEdgeTopology(node,**kwargs) + return + case ElementType.TRIANGLES: + addTriangleTopology(node,**kwargs) + return + case ElementType.QUADS: + addQuadTopology(node,**kwargs) + return + case ElementType.TETRAHEDRA: + addTetrahedronTopology(node,**kwargs) + return + case ElementType.HEXAHEDRA: + addHexahedronTopology(node,**kwargs) + return + case _: + print('Topology type should be one of the following : "ElementType.POINTS, ElementType.EDGES, ElementType.TRIANGLES, ElementType.QUADS, ElementType.TETRAHEDRA, ElementType.HEXAHEDRA" ') + return + \ No newline at end of file diff --git a/splib/topology/loader.py b/splib/topology/loader.py new file mode 100644 index 000000000..c4f54b01f --- /dev/null +++ b/splib/topology/loader.py @@ -0,0 +1,22 @@ +from splib.core.node_wrapper import ReusableMethod + +@ReusableMethod +def loadMesh(node,filename,**kwargs): + splitedName = filename.split('.') + if len(splitedName) == 1: + print('[Error] : A file name without extension was provided.') + return + + if splitedName[-1] in ['vtk', 'obj', 'stl', 'msh', 'sph']: + if splitedName[-1] == "msh": + return node.addObject("MeshGmshLoader", name="loader",filename=filename, **kwargs) + elif splitedName[-1] == "sph": + return node.addObject("SphereLoader", name="loader",filename=filename, **kwargs) + else: + return node.addObject("Mesh"+splitedName[-1].upper()+"Loader", name="loader", filename=filename, **kwargs) + else: + print('[Error] : File extension ' + splitedName[-1] + ' not recognised.') + + + + diff --git a/splib/topology/static.py b/splib/topology/static.py new file mode 100644 index 000000000..591460425 --- /dev/null +++ b/splib/topology/static.py @@ -0,0 +1,7 @@ +from splib.core.node_wrapper import ReusableMethod +from splib.core.utils import DEFAULT_VALUE + +@ReusableMethod +def addStaticTopology(node, source=DEFAULT_VALUE, **kwargs): + node.addObject("MeshTopology", name="container", src=source, **kwargs) + diff --git a/stlib/README.md b/stlib/README.md new file mode 100644 index 000000000..c4e7dd901 --- /dev/null +++ b/stlib/README.md @@ -0,0 +1,52 @@ +# STLIB + +## Terminology + +| Term | Description | +| -------------- | -------------------------------------------------------------- | +| Component* | Element of the scene hierarchy implementing a given behavior | +| ~~Object~~ | A deprecated synonym of a Component | +| Node* | Element of the scene hierarchy holding other Nodes (often refered as childs) or Components | +| Data* | Attribute of a Component or a Node | +| Prefab | A Node assembling of Components and Nodes (a "fragment" of a scene) | +| Geometry | A Prefab that describes shapes with their topologies (i.e a shapes with their space descritization and their associated connectivity) | +| Entity | A physical Prefab that represents real-world properties and behaviors used in a simulation. An entity should always have a geometry but includes neither a linear solver nor an integration scheme.| +| Parameters | Every Prefab has a set of parameters. These parameters can contain data, links, callable or being composed of other parameters. Some of them can be optional. ~~Must inherit from `stlib.core.baseParameter.BaseParameter` and have a `@dataclasses.dataclass` decorator~~. Must have a `@stlib.parameters` decorator. | + +\*Defined in SOFA documentation [here](https://www.sofa-framework.org/doc/using-sofa/terminology). + + +## Concepts & Structure + +This library is structured to provide a set of _Prefabs_ that can be used to build complex simulations in SOFA. +Prefabs are reusable fragments of a scene that can be specified through Parameters. +We introduce two different concepts, Prefab and Parameters: +- Prefabs defining the logic of instantiation +- Parameters providing the information (e.g data, links, callable) needed by a Prefab for its own instantiation + +We introduce two types of Prefabs: +- __Geometry__: A prefab that describes shapes with their topologies (i.e a shape with its space discretization and its associated connectivity). +- __Entity__: A physical prefab that represents real-world properties and behaviors used in a simulation. An entity should always have a geometry but includes neither a linear solver nor an integration scheme. + + +## Usage + +STLIB has been designed to suit the following levels of use: + +- __Beginners__: + - Create simple simulations using predefined Prefabs. + - Use the provided Prefabs and Parameters without needing to understand the underlying implementation. +- __Intermediate users__: + - Create more complex simulations by combining existing Prefabs. + - Redefine Parameters for their own usage. +- __Advanced users__: + - Create their own Prefabs from scratch or by extending the provided ones. + - Enrich the library with new Prefabs and Parameters. + +## Available Parameters + +Each Prefab comes with a set of Parameters, these have been selected when the following criteria are met: +- Data which corresponds to an intraseque property of the Prefab, required for its instantiation +- Data which does not have a default value +- Data which cannot be inferred +- Data which are commun to all possible usage of the Prefab diff --git a/stlib/__init__.py b/stlib/__init__.py new file mode 100644 index 000000000..0ad659858 --- /dev/null +++ b/stlib/__init__.py @@ -0,0 +1,68 @@ +__all__ = ["core","entities","geometries","materials","collision","visual"] + +import Sofa.Core +from stlib.core.basePrefab import BasePrefab + +def __genericAdd(self : Sofa.Core.Node, typeName, **kwargs): + def findName(cname, names): + """Compute a working unique name in the node""" + rname = cname + for i in range(0, len(names)): + if rname not in names: + return rname + rname = cname + str(i+1) + return rname + + + def checkName(context : Sofa.Core.Node, name): + # Check if the name already exists, if this happens, create a new one. + if name in context.children or name in context.objects: + names = {node.name.value for node in context.children} + names = names.union({object.name.value for object in context.objects}) + name = findName(name, names) + return name + + + # Check if a name is provided, if not, use the one of the class + params = kwargs.copy() + if isinstance(typeName, type) and issubclass(typeName, BasePrefab): #Only for prefabs + if len(params.keys()) > 1 or (len(params.keys()) == 1 and "parameters" not in params): + raise RuntimeError("Invalid argument, a prefab takes only the \"parameters\" kwargs as input") + + elif "name" not in params : #This doesn't apply to prefab + if isinstance(typeName, str): + params["name"] = typeName + elif isinstance(typeName, type) and issubclass(typeName, Sofa.Core.Node): + params["name"] = typeName.__name__ + elif isinstance(typeName, Sofa.Core.Node): + params["name"] = "Node" + elif isinstance(typeName, type) and issubclass(typeName, Sofa.Core.Object): + params["name"] = typeName.name.value + elif isinstance(typeName, type) and issubclass(typeName, Sofa.Core.ObjectDeclaration): + params["name"] = typeName.__name__ + else: + raise RuntimeError("Invalid argument ", typeName) + + if isinstance(typeName, type) and issubclass(typeName, BasePrefab) and len(params.keys()) == 1: + params["parameters"].name = checkName(self, params["parameters"].name) + else: + params["name"] = checkName(self, params["name"]) + + # Dispatch the creation to either addObject or addChild + if isinstance(typeName, type) and issubclass(typeName, BasePrefab): + pref = self.addChild(typeName(**params)) + pref.init() + elif isinstance(typeName, Sofa.Core.Node): + pref = self.addChild(typeName(**params)) + elif isinstance(typeName, type) and issubclass(typeName, Sofa.Core.Object): + pref = self.addObject(typeName(**params)) + elif isinstance(typeName, type) and issubclass(typeName, Sofa.Core.ObjectDeclaration): + pref = self.addObject(typeName.__name__, **params) + elif isinstance(typeName, str): + pref = self.addObject(typeName, **params) + else: + raise RuntimeError("Invalid argument", typeName) + return pref + +# Inject the method so it become available as if it was part of Sofa.Core.Node +Sofa.Core.Node.add = __genericAdd diff --git a/stlib/collision.py b/stlib/collision.py new file mode 100644 index 000000000..b8f059249 --- /dev/null +++ b/stlib/collision.py @@ -0,0 +1,64 @@ +from stlib.core.basePrefab import BasePrefab +from stlib.core.baseParameters import BaseParameters, Optional, dataclasses +from stlib.geometries import Geometry, GeometryParameters +from stlib.geometries.file import FileParameters +from splib.core.enum_types import CollisionPrimitive +from splib.core.utils import DEFAULT_VALUE +from splib.mechanics.collision_model import addCollisionModels +from Sofa.Core import Object + +@dataclasses.dataclass +class CollisionParameters(BaseParameters): + name : str = "Collision" + + primitives : list[CollisionPrimitive] = dataclasses.field(default_factory = lambda :[CollisionPrimitive.TRIANGLES]) + + selfCollision : Optional[bool] = DEFAULT_VALUE + bothSide : Optional[bool] = DEFAULT_VALUE + group : Optional[int] = DEFAULT_VALUE + contactDistance : Optional[float] = DEFAULT_VALUE + + geometry : GeometryParameters = dataclasses.field(default_factory = lambda : GeometryParameters()) + + +class Collision(BasePrefab): + def __init__(self, parameters: CollisionParameters): + BasePrefab.__init__(self, parameters) + + def init(self): + + geom = self.add(Geometry, parameters = self.parameters.geometry) + + self.addObject("MechanicalObject", template="Vec3", position=f"@{self.parameters.geometry.name}/container.position") + for primitive in self.parameters.primitives: + addCollisionModels(self, primitive, + topology=f"@{self.parameters.geometry.name}/container", + selfCollision=self.parameters.selfCollision, + group=self.parameters.group, + **self.parameters.kwargs) + + + @staticmethod + def getParameters(**kwargs) -> CollisionParameters: + return CollisionParameters(**kwargs) + + +def createScene(root): + + root.addObject("VisualStyle", displayFlags="showCollisionModels") + + # Create a visual from a mesh file + parameters = Collision.getParameters() + parameters.group = 1 + parameters.geometry = FileParameters(filename="mesh/cube.obj") + # Expert parameters + # parameters.kwargs = { + # "TriangleCollisionModel":{"contactStiffness": 100.0, "contactFriction": 0.5} + # } + collision = root.add(Collision, parameters) + + # OR set the parameters post creation + # collision.TriangleCollisionModel.contactStiffness = 100.0 + # collision.TriangleCollisionModel.contactFriction = 0.5 + # collision.TriangleCollisionModel.set(contactStiffness=100.0, contactFriction=0.5) # we have information of what is possible + # collision.TriangleCollisionModel.set({"contactStiffness": 100.0, "contactFriction": 0.5}) # we can do n'importe quoi diff --git a/stlib/core/baseEntity.py b/stlib/core/baseEntity.py new file mode 100644 index 000000000..1456436bb --- /dev/null +++ b/stlib/core/baseEntity.py @@ -0,0 +1,11 @@ +import Sofa.Core +from baseParameters import BaseParameters + +class BaseEntity(Sofa.Core.Prefab): + + parameters : BaseParameters + + def __init__(self): + Sofa.Core.Prefab.__init__(self) + + diff --git a/stlib/core/baseGeometry.py b/stlib/core/baseGeometry.py new file mode 100644 index 000000000..e69de29bb diff --git a/stlib/core/baseParameters.py b/stlib/core/baseParameters.py new file mode 100644 index 000000000..be7a633b8 --- /dev/null +++ b/stlib/core/baseParameters.py @@ -0,0 +1,13 @@ +import dataclasses +from splib.core.utils import DEFAULT_VALUE + +import dataclasses +from typing import Callable, Optional, Any + +@dataclasses.dataclass +class BaseParameters(object): + name : str = "Object" + kwargs : dict = dataclasses.field(default_factory=dict) + + def toDict(self): + return dataclasses.asdict(self) diff --git a/stlib/core/basePrefab.py b/stlib/core/basePrefab.py new file mode 100644 index 000000000..83d29f55d --- /dev/null +++ b/stlib/core/basePrefab.py @@ -0,0 +1,23 @@ +import copy +import Sofa +import Sofa.Core +from stlib.core.basePrefabParameters import BasePrefabParameters + +class BasePrefab(Sofa.Core.Node): + """ + A Prefab is a Sofa.Node that assembles a set of components and nodes + """ + + def __init__(self, parameters: BasePrefabParameters): + Sofa.Core.Node.__init__(self, name=parameters.name) + self.parameters = parameters + + def init(self): + raise NotImplemented("To be overridden by child class") + + + def localToGlobalCoordinates(pointCloudInput, pointCloudOutput): + raise NotImplemented("Send an email to Damien, he will help you. Guaranteed :)") + + + diff --git a/stlib/core/basePrefabParameters.py b/stlib/core/basePrefabParameters.py new file mode 100644 index 000000000..d6a69dee9 --- /dev/null +++ b/stlib/core/basePrefabParameters.py @@ -0,0 +1,15 @@ +import dataclasses + +@dataclasses.dataclass +class BasePrefabParameters(object): + name : str = "object" + kwargs : dict = dataclasses.field(default_factory=dict) + + # Transformation information + # TODO: these data are going to be added in Node in SOFA (C++ implementation) + translation : list[float] = dataclasses.field(default_factory = lambda : [0., 0., 0.]) + rotation : list[float] = dataclasses.field(default_factory = lambda : [0., 0., 0.]) + scale : list[float] = dataclasses.field(default_factory = lambda : [1., 1., 1.]) + + def toDict(self): + return dataclasses.asdict(self) diff --git a/stlib/entities/__entity__.py b/stlib/entities/__entity__.py new file mode 100644 index 000000000..9a332dca2 --- /dev/null +++ b/stlib/entities/__entity__.py @@ -0,0 +1,80 @@ +from stlib.core.baseParameters import BaseParameters +from stlib.collision import CollisionParameters, Collision +from stlib.visual import VisualParameters, Visual +from stlib.materials import Material, MaterialParameters +from stlib.geometries import Geometry +import dataclasses +from typing import Callable, Optional +from stlib.geometries import GeometryParameters +from splib.core.enum_types import StateType +from stlib.core.basePrefab import BasePrefab + + +@dataclasses.dataclass +class EntityParameters(BaseParameters): + name : str = "Entity" + + stateType : StateType = StateType.VEC3 + + ### QUID + addCollision : Optional[Callable] = lambda x : Collision(CollisionParameters()) + addVisual : Optional[Callable] = lambda x : Visual(VisualParameters()) + + geometry : GeometryParameters = None + material : MaterialParameters = None + collision : Optional[CollisionParameters] = None + visual : Optional[VisualParameters] = None + + + +class Entity(BasePrefab): + + # A simulated object + material : Material + visual : Visual + collision : Collision + geometry : Geometry + + parameters : EntityParameters + + + def __init__(self, parameters=EntityParameters(), **kwargs): + BasePrefab.__init__(self, parameters) + + + def init(self): + self.geometry = self.add(Geometry, parameters=self.parameters.geometry) + + ### Check compatilibility of Material + if self.parameters.material.stateType != self.parameters.stateType: + print("WARNING: imcompatibility between templates of both the entity and the material") + self.parameters.material.stateType = self.parameters.stateType + + self.material = self.add(Material, parameters=self.parameters.material) + self.material.States.position.parent = self.geometry.container.position.linkpath + + if self.parameters.collision is not None: + self.collision = self.add(Collision, parameters=self.parameters.collision) + self.addMapping(self.collision) + + if self.parameters.visual is not None: + self.visual = self.add(Visual, parameters=self.parameters.visual) + self.addMapping(self.visual) + + + def addMapping(self, destinationPrefab): + + template = f'{self.parameters.stateType},Vec3' # TODO: check that it is always true + + 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, + template=template) + else: + destinationPrefab.addObject("RigidMapping", + output=destinationPrefab.linkpath, + input=self.material.linkpath, + template=template) diff --git a/stlib/entities/__init__.py b/stlib/entities/__init__.py new file mode 100644 index 000000000..bfad7c6c5 --- /dev/null +++ b/stlib/entities/__init__.py @@ -0,0 +1 @@ +from .__entity__ import * \ No newline at end of file diff --git a/stlib/geometries/__geometry__.py b/stlib/geometries/__geometry__.py new file mode 100644 index 000000000..58a249fd8 --- /dev/null +++ b/stlib/geometries/__geometry__.py @@ -0,0 +1,74 @@ +from stlib.core.basePrefab import BasePrefab +from stlib.core.baseParameters import BaseParameters, Optional, dataclasses, Any +from splib.topology.dynamic import addDynamicTopology +from splib.topology.static import addStaticTopology +from splib.core.enum_types import ElementType +from splib.core.utils import DEFAULT_VALUE +from Sofa.Core import Object + +import numpy as np + + +class Geometry(BasePrefab):... + +@dataclasses.dataclass +class InternalDataProvider(object): + position : Any = None + # Topology information + edges : Any = DEFAULT_VALUE + triangles : Any = DEFAULT_VALUE + quads : Any = DEFAULT_VALUE + tetrahedra : Any = DEFAULT_VALUE + hexahedra : Any = DEFAULT_VALUE + + def generateAttribute(self, parent : Geometry): + pass + + +@dataclasses.dataclass +class GeometryParameters(BaseParameters): + name : str = "Geometry" + + # Type of the highest degree element + elementType : Optional[ElementType] = None + data : Optional[InternalDataProvider] = None + + dynamicTopology : bool = False + + def Data(self): + return InternalDataProvider() + + + +class Geometry(BasePrefab): + # container : Object # This should be more specialized into the right SOFA type + # modifier : Optional[Object] + + parameters : GeometryParameters + + def __init__(self, parameters: GeometryParameters): + BasePrefab.__init__(self, parameters) + + + + 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)) + else: + raise ValueError + else: + addStaticTopology(self, + 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 + }) \ No newline at end of file diff --git a/stlib/geometries/__init__.py b/stlib/geometries/__init__.py new file mode 100644 index 000000000..34a466c77 --- /dev/null +++ b/stlib/geometries/__init__.py @@ -0,0 +1 @@ +from .__geometry__ import * \ No newline at end of file diff --git a/stlib/geometries/cube.py b/stlib/geometries/cube.py new file mode 100644 index 000000000..4d5a52ca3 --- /dev/null +++ b/stlib/geometries/cube.py @@ -0,0 +1,14 @@ +from stlib.geometries import GeometryParameters + +class CubeParameters(GeometryParameters): + def __init__(self, center, edgeLength, pointPerEdge, dynamicTopology = False): + + customGeom = CubeParameters.createData(center, edgeLength, pointPerEdge) + GeometryParameters.__init__(data = customGeom, dynamicTopology = dynamicTopology) + + @staticmethod + def createData(center, edgeLength, pointPerEdge) -> GeometryParameters.Data : + data = GeometryParameters.Data() + #Fill data + return data + \ No newline at end of file diff --git a/stlib/geometries/extract.py b/stlib/geometries/extract.py new file mode 100644 index 000000000..7143fb07d --- /dev/null +++ b/stlib/geometries/extract.py @@ -0,0 +1,71 @@ +from stlib.geometries import GeometryParameters, InternalDataProvider, Geometry +from stlib.core.baseParameters import dataclasses +from splib.topology.dynamic import addDynamicTopology +from splib.topology.loader import loadMesh +from splib.core.enum_types import ElementType + +import Sofa +from Sofa.Core import Node + + +class ExtractInternalDataProvider(InternalDataProvider): + destinationType : ElementType + sourceType : ElementType + sourceName : str + + def __init__(self, destinationType : ElementType, sourceType : ElementType, sourceName : str): + self.destinationType = destinationType + self.sourceType = sourceType + self.sourceName = sourceName + + def __post_init__(self): + if(not (self.sourceType == ElementType.TETRAHEDRA and self.destinationType == ElementType.TRIANGLES) + and not (self.sourceType == ElementType.HEXAHEDRA and self.destinationType == ElementType.QUADS) ): + raise ValueError("Only configuration possible are 'Tetrahedra to Triangles' and 'Hexahedra to Quads'") + + InternalDataProvider.__init__(self) + + 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 + # fromLink = parent.parents[0].parents[0].getChild(self.SourceName).container.linkpath + # TODO: the line above cannot work if the nodes and objects are not added to the graph prior the end of __init__() call + # !!! 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) + if self.sourceType == ElementType.TETRAHEDRA: + node.addObject("Tetra2TriangleTopologicalMapping", input=fromLink, output=node.container.linkpath) + elif self.sourceType == ElementType.HEXAHEDRA: + node.addObject("Hexa2QuadTopologicalMapping", input=fromLink, output=node.container.linkpath) + else: + Sofa.msg_error("[stlib/geometry/exctrat.py]", "Element type: " + str(self.sourceType) + " not supported.") + + self.position = node.container.position.linkpath + if node.container.findData("edges") is not None: + self.edges = node.container.edges.linkpath + if node.container.findData("triangles") is not None: + self.triangles = node.container.triangles.linkpath + if node.container.findData("quads") is not None: + self.quads = node.container.quads.linkpath + if node.container.findData("hexahedra") is not None: + self.hexahedra = node.container.hexahedra.linkpath + if node.container.findData("tetras") is not None: + self.tetrahedra = node.container.tetras.linkpath + + + +class ExtractParameters(GeometryParameters): + def __init__(self, + sourceParameters : GeometryParameters, + destinationType : ElementType, + dynamicTopology : bool = False, ): + GeometryParameters.__init__(self, + data = ExtractInternalDataProvider(destinationType = destinationType, + sourceType = sourceParameters.elementType, + sourceName = sourceParameters.name), + dynamicTopology = dynamicTopology, + elementType = destinationType) + diff --git a/stlib/geometries/file.py b/stlib/geometries/file.py new file mode 100644 index 000000000..b057728c1 --- /dev/null +++ b/stlib/geometries/file.py @@ -0,0 +1,34 @@ +from stlib.geometries import GeometryParameters, InternalDataProvider, Geometry +from stlib.core.baseParameters import dataclasses +from splib.topology.loader import loadMesh +from splib.core.enum_types import ElementType + +from Sofa.Core import Node + +@dataclasses.dataclass +class FileInternalDataProvider(InternalDataProvider): + filename : str = "mesh/cube.obj" + + 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) + + + +class FileParameters(GeometryParameters): + + def __init__(self, filename, dynamicTopology = False, elementType : ElementType = None ): + GeometryParameters.__init__(self, + data = FileInternalDataProvider(filename=filename), + dynamicTopology = dynamicTopology, + elementType = elementType) + diff --git a/stlib/geometries/plane.py b/stlib/geometries/plane.py new file mode 100644 index 000000000..2c3c637a4 --- /dev/null +++ b/stlib/geometries/plane.py @@ -0,0 +1,45 @@ +from stlib.geometries import GeometryParameters, InternalDataProvider, Geometry +import dataclasses +import numpy as np + +@dataclasses.dataclass +class PlaneDataProvider(InternalDataProvider): + center : np.ndarray[float] = dataclasses.field(default_factory = lambda : np.array([0,0,0])) + normal : np.ndarray[float] = dataclasses.field(default_factory = lambda : np.array([0,0,1])) + lengthNormal : np.ndarray[float] = dataclasses.field(default_factory = lambda : np.array([1,0,0])) + lengthNbEdge : int = 1 + widthNbEdge : int = 1 + lengthSize : float = 1.0 + widthSize : float = 1.0 + + def __post_init__(self, **kwargs): + InternalDataProvider.__init__(self,**kwargs) + + def generateAttribute(self, parent : Geometry): + + lengthEdgeSize = self.lengthSize / self.lengthNbEdge + widthEdgeSize = self.widthSize / self.widthNbEdge + + self.widthNormal = np.cross(self.normal,self.lengthNormal) + bottomLeftCorner = self.center - self.lengthNormal * self.lengthNbEdge * lengthEdgeSize / 2 - self.widthNormal * self.widthNbEdge * widthEdgeSize / 2 + + self.position = np.array([[ bottomLeftCorner + j * self.widthNormal * widthEdgeSize + i * self.lengthNormal * lengthEdgeSize for j in range(self.widthNbEdge + 1) ] for i in range(self.lengthNbEdge + 1)]) + + self.triangles = np.empty((2 * self.widthNbEdge * self.lengthNbEdge, 3), dtype = int) + for i in range(self.lengthNbEdge): + for j in range(self.widthNbEdge): + self.triangles[i*self.widthNbEdge*2 + j * 2 , 0] = j + i * (self.widthNbEdge + 1) + self.triangles[i*self.widthNbEdge*2 + j * 2 , 1] = j + (i+1) * (self.widthNbEdge + 1) + self.triangles[i*self.widthNbEdge*2 + j * 2 , 2] = j + 1 + i * (self.widthNbEdge + 1) + self.triangles[i*self.widthNbEdge*2 + j * 2 + 1, 0] = j + 1 + i * (self.widthNbEdge + 1) + self.triangles[i*self.widthNbEdge*2 + j * 2 + 1, 1] = j + (i+1) * (self.widthNbEdge + 1) + self.triangles[i*self.widthNbEdge*2 + j * 2 + 1, 2] = j + 1 + (i+1) * (self.widthNbEdge + 1) + + + + +class PlaneParameters(GeometryParameters): + + def __init__(self, center, normal, lengthNormal, lengthNbEdge, widthNbEdge, lengthSize, widthSize, dynamicTopology = False): + GeometryParameters.__init__(self, data = PlaneDataProvider(center=center, normal=normal, lengthNormal=lengthNormal, lengthNbEdge=lengthNbEdge, widthNbEdge=widthNbEdge, lengthSize=lengthSize, widthSize=widthSize), + dynamicTopology = dynamicTopology) diff --git a/stlib/geometries/sphere.py b/stlib/geometries/sphere.py new file mode 100644 index 000000000..2cb64f6ee --- /dev/null +++ b/stlib/geometries/sphere.py @@ -0,0 +1,14 @@ +from stlib.geometries import GeometryParameters + +class SphereParameters(GeometryParameters): + def __init__(self, center, radius, pointPerRad, dynamicTopology = False): + + customGeom = SphereParameters.createData(center, radius, pointPerRad) + GeometryParameters.__init__(data = customGeom, dynamicTopology = dynamicTopology) + + @staticmethod + def createData(center, radius, pointPerRad) -> GeometryParameters.Data : + data = GeometryParameters.Data() + #Fill data + return data + \ No newline at end of file diff --git a/stlib/materials/__init__.py b/stlib/materials/__init__.py new file mode 100644 index 000000000..6a748749a --- /dev/null +++ b/stlib/materials/__init__.py @@ -0,0 +1 @@ +from .__material__ import * \ No newline at end of file diff --git a/stlib/materials/__material__.py b/stlib/materials/__material__.py new file mode 100644 index 000000000..ff88dbd43 --- /dev/null +++ b/stlib/materials/__material__.py @@ -0,0 +1,31 @@ +from stlib.core.baseParameters import BaseParameters, Callable, Optional, dataclasses, Any +from splib.core.utils import defaultValueType, DEFAULT_VALUE, isDefault +from splib.core.enum_types import StateType + +from stlib.core.basePrefab import BasePrefab +from splib.mechanics.mass import addMass + +@dataclasses.dataclass +class MaterialParameters(BaseParameters): + name : str = "Material" + + massDensity : float = DEFAULT_VALUE + massLumping : bool = DEFAULT_VALUE + + stateType : StateType = StateType.VEC3 + + addMaterial : Optional[Callable] = lambda node : addMass(node, node.parameters.stateType, massDensity=node.parameters.massDensity, lumping=node.parameters.massLumping) + + +# TODO : previously called Behavior +class Material(BasePrefab): + + parameters : MaterialParameters + + def __init__(self, parameters: MaterialParameters): + BasePrefab.__init__(self, parameters) + + + def init(self): + self.addObject("MechanicalObject", name="States", template=str(self.parameters.stateType)) + self.parameters.addMaterial(self) diff --git a/stlib/materials/deformable.py b/stlib/materials/deformable.py new file mode 100644 index 000000000..a24763df5 --- /dev/null +++ b/stlib/materials/deformable.py @@ -0,0 +1,52 @@ +from stlib.materials import MaterialParameters +from splib.core.enum_types import ConstitutiveLaw +from stlib.core.baseParameters import Callable, Optional, dataclasses +from splib.mechanics.linear_elasticity import * +from splib.mechanics.hyperelasticity import * +from splib.mechanics.mass import addMass + + +@dataclasses.dataclass +class DeformableBehaviorParameters(MaterialParameters): + + constitutiveLawType : ConstitutiveLaw = ConstitutiveLaw.ELASTIC + elementType : ElementType = ElementType.TETRAHEDRA + parameters : list[float] = dataclasses.field(default_factory=lambda: [1000, 0.45]) # young modulus, poisson ratio + + 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 + massKwargs["geometryState"] = "@States" + + addMass(node, node.parameters.elementType, massDensity=node.parameters.massDensity, lumping=node.parameters.massLumping, topology="@../Geometry/container", mass=massKwargs) + # TODO : change this with inheritance + if(node.parameters.constitutiveLawType == ConstitutiveLaw.HYPERELASTIC): + addHyperelasticity(node, node.parameters.elementType, node.parameters.parameters, topology="@../Geometry/container") + else: + addLinearElasticity(node, node.parameters.elementType, node.parameters.parameters[0], node.parameters.parameters[1], topology="@../Geometry/container") + + addMaterial : Optional[Callable] = __addDeformableMaterial + + +def createScene(root) : + from stlib.entities import Entity, EntityParameters + from stlib.visual import VisualParameters + from stlib.geometries.file import FileParameters + + root.addObject('RequiredPlugin', name='Sofa.Component.Visual') # Needed to use components [VisualStyle] + root.addObject("VisualStyle", displayFlags=["showBehavior"]) + + bunnyParameters = EntityParameters() + bunnyParameters.geometry = FileParameters(filename="mesh/Bunny.vtk") + bunnyParameters.geometry.elementType = ElementType.TETRAHEDRA # TODO: this is required by extract.py. Should it be done automatically in geometry.py ? + bunnyParameters.material = DeformableBehaviorParameters() + bunnyParameters.material.constitutiveLawType = ConstitutiveLaw.ELASTIC + bunnyParameters.material.parameters = [1000, 0.45] + bunnyParameters.visual = VisualParameters() + # bunnyParameters.visual.geometry = ExtractParameters(sourceParameters=bunnyParameters.geometry, + # destinationType=ElementType.TRIANGLES) + bunnyParameters.visual.geometry = FileParameters(filename="mesh/Bunny.stl") + bunnyParameters.visual.color = [1, 1, 1, 0.5] + bunny = root.add(Entity, parameters=bunnyParameters) + # bunny.init() \ No newline at end of file diff --git a/stlib/materials/rigid.py b/stlib/materials/rigid.py new file mode 100644 index 000000000..5668b3562 --- /dev/null +++ b/stlib/materials/rigid.py @@ -0,0 +1,13 @@ +from stlib.core.baseParameters import BaseParameters, Optional, dataclasses +from stlib.geometries import GeometryParameters + + + +@dataclasses.dataclass +class RigidParameters(BaseParameters): + + geometry : GeometryParameters + mass : Optional[float] = None + + def toDict(self): + return dataclasses.asdict(self) diff --git a/stlib/visual.py b/stlib/visual.py new file mode 100644 index 000000000..a0f377f7d --- /dev/null +++ b/stlib/visual.py @@ -0,0 +1,39 @@ +from stlib.core.basePrefab import BasePrefab +from stlib.core.baseParameters import BaseParameters, Optional, dataclasses +from stlib.geometries import Geometry, GeometryParameters +from stlib.geometries.file import FileParameters +from splib.core.utils import DEFAULT_VALUE +from Sofa.Core import Object + +@dataclasses.dataclass +class VisualParameters(BaseParameters): + name : str = "Visual" + + color : Optional[list[float]] = DEFAULT_VALUE + texture : Optional[str] = DEFAULT_VALUE + + geometry : GeometryParameters = dataclasses.field(default_factory = lambda : GeometryParameters()) + + +class Visual(BasePrefab): + + def __init__(self, parameters: VisualParameters): + BasePrefab.__init__(self, parameters) + + def init(self): + self.geometry = self.add(Geometry, parameters=self.parameters.geometry) + self.addObject("OglModel", color=self.parameters.color, src=self.geometry.container.linkpath) + + + @staticmethod + def getParameters(**kwargs) -> VisualParameters: + return VisualParameters(**kwargs) + + +def createScene(root): + + # Create a visual from a mesh file + parameters = Visual.getParameters() + parameters.name = "LiverVisual" + parameters.geometry = FileParameters(filename="mesh/liver.obj") + root.add(Visual, parameters) \ No newline at end of file diff --git a/tests/test_new_add.py b/tests/test_new_add.py new file mode 100644 index 000000000..0c1036b65 --- /dev/null +++ b/tests/test_new_add.py @@ -0,0 +1,84 @@ +import unittest +import Sofa +import SofaRuntime +import Sofa.Core +import stlib + +class ObjectDeclaration(object): + ... + +Sofa.Core.ObjectDeclaration = ObjectDeclaration + +class MechanicalObject(ObjectDeclaration): + pass + +class TestNewAdd(unittest.TestCase): + def test_add_node_with_node_type(self): + root = Sofa.Core.Node("root") + root.add(Sofa.Core.Node, name="aNodeA") + self.assertEqual(len(root.children), 1) + self.assertEqual(root.children[0].name.value, "aNodeA") + + def test_add_node_with_node_instance(self): + root = Sofa.Core.Node("root") + root.add(Sofa.Core.Node("aNodeB")) + self.assertEqual(len(root.children), 1) + self.assertEqual(root.children[0].name.value, "aNodeB") + + def test_add_object_with_string_type(self): + root = Sofa.Core.Node("root") + root.add("MechanicalObject", name="anObject1", position=[[1,2,3]]) + self.assertEqual(len(root.objects), 1) + self.assertEqual(root.objects[0].name.value, "anObject1") + self.assertEqual(root.objects[0].position.value.shape, (1,3)) + + def test_add_object_with_object_type(self): + root = Sofa.Core.Node("root") + root.add(MechanicalObject, name="anObject2", position=[[1,2,3]]) + self.assertEqual(len(root.objects), 1) + self.assertEqual(root.objects[0].name.value, "anObject2") + self.assertEqual(root.objects[0].position.value.shape, (1,3)) + + def test_automatic_name_generation(self): + root = Sofa.Core.Node("root") + root.add(MechanicalObject, position=[[1,2,3]]) + root.add(MechanicalObject, position=[[1,2,3]]) + root.add(MechanicalObject, position=[[1,2,3]]) + self.assertEqual(root.objects[0].name.value, "MechanicalObject") + self.assertEqual(root.objects[1].name.value, "MechanicalObject1") + self.assertEqual(root.objects[2].name.value, "MechanicalObject2") + + root.add(Sofa.Core.Node, name="TestNode") + root.add(Sofa.Core.Node, name="TestNode") + self.assertEqual(root.children[0].name.value, "TestNode") + self.assertEqual(root.children[1].name.value, "TestNode1") + + root.add(Sofa.Core.Node) + root.add(Sofa.Core.Node) + self.assertEqual(root.children[2].name.value, "Node") + self.assertEqual(root.children[3].name.value, "Node1") + + def test_add_node_with_kwargs(self): + root = Sofa.Core.Node("root") + root.add(Sofa.Core.Node, name="aNodeC", gravity=[1,2,3]) + self.assertEqual(root.children[0].gravity.value, [1,2,3]) + + def test_add_instanciated_prefab(self): + root = Sofa.Core.Node("root") + from stlib.entities import Entity, EntityParameters + + bunnyParameters = EntityParameters() + bunny = root.add(Entity(bunnyParameters)) + + def test_add_prefab_with_parameter_object(self): + root = Sofa.Core.Node("root") + from stlib.entities import Entity, EntityParameters + + bunnyParameters = EntityParameters() + bunny = root.add(Entity, bunnyParameters) + + +if __name__ == '__main__': + + SofaRuntime.importPlugin("Sofa.Component.StateContainer") + unittest.main()