<HexahedronElasticForce />¶
Doxygen:
SofaCaribou::forcefield::HexahedronElasticForce
Implementation of a corotational linear elasticity forcefield for hexahedral topologies.
Requires a mechanical object. Requires a topology container.
Attribute |
Format |
Default |
Description |
---|---|---|---|
printLog |
bool |
false |
Output informative messages at the initialization and during the simulation. |
youngModulus |
float |
1000 |
Young’s modulus of the material |
poissonRatio |
float |
0.3 |
Poisson’s ratio of the material |
corotated |
bool |
true |
Whether or not to use corotated elements for the strain computation. The rotation is viewed as constant on the element and is extracted at its center point. |
topology_container |
path |
Path to a topology container (or path to a mesh) that contains the hexahedral elements. |
|
integration_method |
option |
Regular |
Integration method used to integrate the stiffness matrix.
|
Quick example¶
XML
<Node>
<RegularGridTopology name="grid" min="-7.5 -7.5 0" max="7.5 7.5 80" n="9 9 21" />
<MechanicalObject src="@grid" />
<HexahedronSetTopologyContainer name="topology" src="@grid" />
<HexahedronElasticForce topology_container="@topology" youngModulus="3000" poissonRatio="0.49" corotated="1" printLog="1" />
</Node>
Python
node.addObject("RegularGridTopology", name="grid", min=[-7.5, -7.5, 0], max=[7.5, 7.5, 80], n=[9, 9, 21])
node.addObject("MechanicalObject", src="@grid")
node.addObject("HexahedronSetTopologyContainer", name="topology", src="@grid")
node.addObject("HexahedronElasticForce", topology_container="@topology", youngModulus=3000, poissonRatio=0.49, corotated=True, printLog=True)
Available python bindings¶
- class HexahedronElasticForce¶
- class GaussNode¶
- Members
weight : Gauss node’s weight.
numpy.double
jacobian_determinant : Gauss node’s Jacobian determinant.
numpy.double
dN_dx : Gauss node’s shape derivatives w.r.t. the current position vector.
numpy.ndarray
F : Gauss node’s strain tensor.
numpy.ndarray
- gauss_nodes_of(hexahedron_id)¶
- stiffness_matrix_of(hexahedron_id)¶
- Parameters
hexahedron_id (int) – Index of the hexahedron in the topology container.
- Returns
Reference to the elemental 24x24 tangent stiffness matrix
- Return type
- Note
No copy involved.
Get the elemental 24x24 tangent stiffness matrix of the given hexahedron.
- K()¶
- Returns
Reference to the forcefield tangent stiffness matrix
- Return type
- Note
No copy involved.
Get the tangent stiffness matrix of the force field as a compressed sparse column major matrix.
- cond()¶
- Returns
Condition number of the forcefield’s tangent stiffness matrix
- Return type
- eigenvalues()¶
- Returns
Reference to the eigen values of the forcefield’s tangent stiffness matrix.
- Return type
list[numpy.double]