diff --git a/examples/Components/animationloop/MechanicalMatrixMapper.pyscn b/examples/Components/animationloop/MechanicalMatrixMapper.pyscn index 5da8c0744a9..5b3cb6dce64 100644 --- a/examples/Components/animationloop/MechanicalMatrixMapper.pyscn +++ b/examples/Components/animationloop/MechanicalMatrixMapper.pyscn @@ -1,55 +1,93 @@ +# This Python scene requires the SofaPython3 plugin +# +# This scene shows the use of MechanicalMatrixMapper. A beam made of hexahedrons is simulated. One part is +# deformable, while the other is rigid. +# +# To achieve this behavior, the two parts are defined independently. The deformable part (respectively the rigid part) +# is defined with vec3 (respectively rigid) degrees of freedom. The two parts are joined as a single object using a +# mapping, which has both parts as input. The FEM is really computed in the mapped node. The deformable part is just a +# subset of the FEM degrees of freedom, while the rigid part is mapped on the rigid degree of freedom. +# +# SOFA does not transfer the mapped stiffness matrix naturally (on an assembled linear system). The FEM in the mapped +# node must be transferred to the global matrix system using MechanicalMatrixMapper. Otherwise, the derivatives of the +# FEM forces are considered null. +# + import Sofa def createScene(rootNode): - rootNode.createObject('RequiredPlugin', name='SofaSparseSolver') - rootNode.createObject('RequiredPlugin', name='SofaBoundaryCondition') - rootNode.createObject('RequiredPlugin', name='SofaDeformable') - rootNode.createObject('RequiredPlugin', name='SofaEngine') - rootNode.createObject('RequiredPlugin', name='SofaGeneralAnimationLoop') - rootNode.createObject('RequiredPlugin', name='SofaImplicitOdeSolver') - rootNode.createObject('RequiredPlugin', name='SofaRigid') - rootNode.createObject('RequiredPlugin', name='SofaMiscMapping') - rootNode.createObject('RequiredPlugin', name='SofaSimpleFem') - - rootNode.createObject('VisualStyle', displayFlags='showCollisionModels hideMappings showForceFields') - - meshDivision = rootNode.createChild('meshDivision') - meshDivision.createObject('RegularGridTopology', name='topology', n=[2, 2, 4] , min=[-1, -1, -2], max=[1, 1, 2]) - meshOfStructure = meshDivision.createObject('MeshTopology',name='beamMesh', src="@topology") - meshDivision.createObject('MechanicalObject', template='Vec3d') - boxFixed= meshDivision.createObject('BoxROI',template="Vec3d", name='box_roi_fix',box=[-1, -1, -2.05, 1, 1, -1.5] , position="@beamMesh.position" ) - boxDeformable= meshDivision.createObject('BoxROI',template="Vec3d", name='box_roi_deformable',box=[-1, -1, -2.05, 1, 1, 1.5], position="@beamMesh.position" ) - boxRigid= meshDivision.createObject('BoxROI',template="Vec3d", name='box_roi_rigid',box=[-1, -1, 1.500001, 1, 1, 2.05], position="@beamMesh.position") - - SolverNode= rootNode.createChild('SolverNode') - SolverNode.createObject('EulerImplicitSolver') - SolverNode.createObject('SparseLDLSolver', name="ldlsolveur") - SolverNode.createObject('MechanicalMatrixMapper', template='Vec3d,Rigid3d', object1='@./deformablePartNode/beamPart1Mech', object2='@./RigidNode/rigid1', nodeToParse='@./deformablePartNode/FEMNode') + rootNode.addObject('RequiredPlugin', name='SofaSparseSolver') + rootNode.addObject('RequiredPlugin', name='SofaBoundaryCondition') + rootNode.addObject('RequiredPlugin', name='SofaDeformable') + rootNode.addObject('RequiredPlugin', name='SofaEngine') + rootNode.addObject('RequiredPlugin', name='SofaGeneralAnimationLoop') + rootNode.addObject('RequiredPlugin', name='SofaImplicitOdeSolver') + rootNode.addObject('RequiredPlugin', name='SofaRigid') + rootNode.addObject('RequiredPlugin', name='SofaMiscMapping') + rootNode.addObject('RequiredPlugin', name='SofaSimpleFem') + + rootNode.addObject('VisualStyle', displayFlags='showCollisionModels hideMappings showForceFields') + + meshDivision = rootNode.addChild('meshDivision') + meshDivision.addObject('RegularGridTopology', name='topology', n=[4, 4, 8] , min=[-1, -1, -2], max=[1, 1, 2]) + meshOfStructure = meshDivision.addObject('MeshTopology',name='beamMesh', src="@topology") + meshDivision.addObject('MechanicalObject', template='Vec3d') + boxDeformable= meshDivision.addObject('BoxROI',template="Vec3d", name='box_roi_deformable', box=[-1, -1, -2.05, 1, 1, 1.5], position="@beamMesh.position" ) + boxRigid= meshDivision.addObject('BoxROI',template="Vec3d", name='box_roi_rigid', box=[-1, -1, 1.500001, 1, 1, 2.05], position="@beamMesh.position") + + SolverNode= rootNode.addChild('SolverNode') + SolverNode.addObject('EulerImplicitSolver') + SolverNode.addObject('SparseLDLSolver', name="ldlsolveur") + SolverNode.addObject('MechanicalMatrixMapper', template='Vec3d,Rigid3d', + object1='@./deformablePartNode/beamPart1Mech', + object2='@./RigidNode/rigid1', + nodeToParse='@./deformablePartNode/FEMNode') ##### Deformable Part of the BEAM (Main Body) - deformablePartNode= SolverNode.createChild('deformablePartNode') - deformablePartNode.createObject('MechanicalObject', template='Vec3d',name='beamPart1Mech', position="@"+boxDeformable.getPathName()+".pointsInROI") + deformablePartNode= SolverNode.addChild('deformablePartNode') + deformablePartNode.addObject('MechanicalObject', template='Vec3d',name='beamPart1Mech', position="@"+boxDeformable.getPathName()+".pointsInROI") + deformablePartNode.addObject('BoxROI',template="Vec3d", name='box_roi_fix', box=[-1, -1, -2.05, 1, 1, -1.5], position="@beamPart1Mech.position" ) + deformablePartNode.addObject('FixedConstraint', indices='@box_roi_fix.indices') ##### Rigid Part of the BEAM (Top) - RigidNode= SolverNode.createChild('RigidNode') - RigidNode.createObject("MechanicalObject",template="Rigid3d",name="rigid1", position=[0, 0, 2, 0, 0, 0, 1], showObject=True, showObjectScale=0.5) - RigidifiedNode= RigidNode.createChild('RigidifiedNode') + RigidNode= SolverNode.addChild('RigidNode') + RigidNode.addObject("MechanicalObject",template="Rigid3d",name="rigid1", position=[0, 0, 2, 0, 0, 0, 1], showObject=True, showObjectScale=0.5) + RigidifiedNode= RigidNode.addChild('RigidifiedNode') - RigidifiedNode.createObject("MechanicalObject",name="rigidMecha",template="Vec3d", position="@"+boxRigid.getPathName()+".pointsInROI") - RigidifiedNode.createObject("RigidMapping",globalToLocalCoords="true") + RigidifiedNode.addObject("MechanicalObject",name="rigidMecha",template="Vec3d", position="@"+boxRigid.getPathName()+".pointsInROI") + RigidifiedNode.addObject("RigidMapping",globalToLocalCoords="true") -##### COMBINED BEAM - FEMNode= deformablePartNode.createChild('FEMNode') +##### COMBINED BEAM + FEMNode= deformablePartNode.addChild('FEMNode') RigidifiedNode.addChild(FEMNode) - FEMNode.createObject('MeshTopology',name='meshInput',src="@"+meshOfStructure.getPathName()) - FEMNode.createObject('MechanicalObject', template='Vec3d',name='beamMecha') - FEMNode.createObject('HexahedronFEMForceField', name='HexaFF', src="@meshInput", poissonRatio=0.49, youngModulus=2000) - FEMNode.createObject('RestShapeSpringsForceField', name='restShapeFF', points="@"+boxFixed.getPathName()+'.indices', stiffness=10000, angularStiffness=10000) - FEMNode.createObject('UniformMass', totalMass=0.1) - FEMNode.createObject('ConstantForceField', name='xMoins', indices=15, forces=[1000, 0, 0]) - FEMNode.createObject("SubsetMultiMapping",name="subsetMapping",template="Vec3d,Vec3d", input="@../beamPart1Mech @../../RigidNode/RigidifiedNode/rigidMecha",output="@./beamMecha", indexPairs="0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 0 10 0 11 1 0 1 1 1 2 1 3") + FEMNode.addObject('MeshTopology',name='meshInput',src="@"+meshOfStructure.getPathName()) + FEMNode.addObject('MechanicalObject', template='Vec3d',name='beamMecha') + FEMNode.addObject('HexahedronFEMForceField', name='HexaFF', src="@meshInput", poissonRatio=0.49, youngModulus=2000) + FEMNode.addObject('UniformMass', totalMass=0.1) + + FEMNode.addObject('BoxROI', name='corner', box=[1, 1, 2, 1, 1, 2]) + FEMNode.addObject('ConstantForceField', name='xMoins', indices='@corner.indices', forces=[1000, 0, 0]) + + + meshDivision.init() + boxDeformable.init() + boxRigid.init() + index_pairs = [ [0, 0] ] * len(meshOfStructure.position) + + i = 0 + for index in boxDeformable.indices.value: + index_pairs[index] = [0, i] + i += 1 + i = 0 + for index in boxRigid.indices.value: + index_pairs[index] = [1, i] + i += 1 + + FEMNode.addObject("SubsetMultiMapping",name="subsetMapping",template="Vec3d,Vec3d", + input="@../beamPart1Mech @../../RigidNode/RigidifiedNode/rigidMecha", + output="@./beamMecha", + indexPairs=index_pairs) - return rootNode - + return rootNode \ No newline at end of file