diff --git a/CMakeLists.txt b/CMakeLists.txt
index c6569a363..2ecb75afd 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -58,10 +58,19 @@ set(HEADER_FILES
${BEAMADAPTER_SRC}/component/mapping/MultiAdaptiveBeamMapping.h
${BEAMADAPTER_SRC}/component/mapping/MultiAdaptiveBeamMapping.inl
+ ${BEAMADAPTER_SRC}/component/model/BaseRodSectionMaterial.h
+ ${BEAMADAPTER_SRC}/component/model/BaseRodSectionMaterial.inl
+ ${BEAMADAPTER_SRC}/component/model/RodMeshSection.h
+ ${BEAMADAPTER_SRC}/component/model/RodMeshSection.inl
+ ${BEAMADAPTER_SRC}/component/model/RodSpireSection.h
+ ${BEAMADAPTER_SRC}/component/model/RodSpireSection.inl
+ ${BEAMADAPTER_SRC}/component/model/RodStraightSection.h
+ ${BEAMADAPTER_SRC}/component/model/RodStraightSection.inl
+
${BEAMADAPTER_SRC}/utils/BeamSection.h
${BEAMADAPTER_SRC}/utils/BeamActions.h
${BEAMADAPTER_SRC}/utils/deprecatedcomponent.h
- )
+)
set(SOURCE_FILES
${BEAMADAPTER_SRC}/initBeamAdapter.cpp
@@ -86,7 +95,11 @@ set(SOURCE_FILES
${BEAMADAPTER_SRC}/component/mapping/AdaptiveBeamMapping.cpp
${BEAMADAPTER_SRC}/component/mapping/BeamLengthMapping.cpp
${BEAMADAPTER_SRC}/component/mapping/MultiAdaptiveBeamMapping.cpp
- )
+
+ ${BEAMADAPTER_SRC}/component/model/RodMeshSection.cpp
+ ${BEAMADAPTER_SRC}/component/model/RodSpireSection.cpp
+ ${BEAMADAPTER_SRC}/component/model/RodStraightSection.cpp
+)
if(SofaImplicitField_FOUND)
set(HEADER_FILES ${HEADER_FILES}
diff --git a/doc/AdaptiveBeamMapping_00000001.png b/doc/AdaptiveBeamMapping_00000001.png
new file mode 100644
index 000000000..3137425a2
Binary files /dev/null and b/doc/AdaptiveBeamMapping_00000001.png differ
diff --git a/doc/BeamAdapter.pdf b/doc/BeamAdapter.pdf
index 3ebc90859..43ccf3cb9 100644
Binary files a/doc/BeamAdapter.pdf and b/doc/BeamAdapter.pdf differ
diff --git a/examples/3instruments.scn b/examples/3instruments.scn
index 2e1843149..8d876cdad 100644
--- a/examples/3instruments.scn
+++ b/examples/3instruments.scn
@@ -24,25 +24,34 @@
-
+
+
+
+
+
-
+
+
+
+
+
-
-
+
+
+
+
+
+
diff --git a/examples/3instruments_collis.scn b/examples/3instruments_collis.scn
index 4a5511faa..7c30267b6 100644
--- a/examples/3instruments_collis.scn
+++ b/examples/3instruments_collis.scn
@@ -38,24 +38,33 @@
-
+
+
+
+
+
-
+
+
+
+
+
-
+
+
+
+
+
@@ -83,7 +92,7 @@
-
diff --git a/examples/SingleBeam.scn b/examples/SingleBeam.scn
index 8899cf388..4f9d73e1c 100644
--- a/examples/SingleBeam.scn
+++ b/examples/SingleBeam.scn
@@ -6,7 +6,9 @@
+
+
@@ -15,8 +17,8 @@
-
-
+
+
diff --git a/examples/SingleBeamDeployment.scn b/examples/SingleBeamDeployment.scn
index ca9739617..6b6287acc 100644
--- a/examples/SingleBeamDeployment.scn
+++ b/examples/SingleBeamDeployment.scn
@@ -15,11 +15,10 @@
-
+
+
+
+
diff --git a/examples/SingleBeamDeploymentCollision.scn b/examples/SingleBeamDeploymentCollision.scn
index b455d0f44..4797b11e3 100644
--- a/examples/SingleBeamDeploymentCollision.scn
+++ b/examples/SingleBeamDeploymentCollision.scn
@@ -32,11 +32,11 @@
-
+
+
+
+
+
diff --git a/examples/Tool_from_loader.scn b/examples/Tool_from_loader.scn
index 544bcba57..e76e8fd76 100644
--- a/examples/Tool_from_loader.scn
+++ b/examples/Tool_from_loader.scn
@@ -23,11 +23,12 @@
-
-
+
+
+
+
+
+
diff --git a/src/BeamAdapter/component/controller/InterventionalRadiologyController.inl b/src/BeamAdapter/component/controller/InterventionalRadiologyController.inl
index 1a6026511..8ccbb5ab7 100644
--- a/src/BeamAdapter/component/controller/InterventionalRadiologyController.inl
+++ b/src/BeamAdapter/component/controller/InterventionalRadiologyController.inl
@@ -248,6 +248,9 @@ void InterventionalRadiologyController::onMouseEvent(MouseEvent * mev
template
void InterventionalRadiologyController::onKeyPressedEvent(KeypressedEvent *kev)
{
+
+ return;
+
/// Control keys for interventonal Radiology simulations:
switch(kev->getKey())
{
@@ -308,56 +311,56 @@ template
void InterventionalRadiologyController::onBeginAnimationStep(const double dt)
{
SOFA_UNUSED(dt);
-
- BaseContext* context = getContext();
- auto xInstrTip = sofa::helper::getWriteOnlyAccessor(d_xTip);
- if(m_FF || m_RW)
- {
- int id = d_controlledInstrument.getValue();
- if (id >= (int)xInstrTip.size())
- {
- msg_warning()<<"Controlled Instument num "<getDt();
- else
- {
- unsigned int newSensorData = m_currentSensorData + 1;
-
- while( m_sensorMotionData[newSensorData][0] < context->getTime() )
- {
- m_currentSensorData = newSensorData;
- newSensorData++;
- }
- if(newSensorData >= m_sensorMotionData.size())
- {
- xInstrTip[id] = 0;
- }
- else
- {
- xInstrTip[id] += m_sensorMotionData[m_currentSensorData][1];
- }
- }
- }
- if (m_RW)
- {
- xInstrTip[id] -= d_speed.getValue()* context->getDt();
- // verif min x :
- if ( xInstrTip[id] < 0.0)
- {
- xInstrTip[id] = 0.0;
- m_RW = false;
- }
- }
- }
-
- /// The tip of the instrument can not be further than its total length
- for (unsigned int i=0; i m_instrumentsList[i]->getRestTotalLength() )
- xInstrTip[i] = m_instrumentsList[i]->getRestTotalLength();
+
+ //BaseContext* context = getContext();
+ //auto xInstrTip = sofa::helper::getWriteOnlyAccessor(d_xTip);
+ //if(m_FF || m_RW)
+ //{
+ // int id = d_controlledInstrument.getValue();
+ // if (id >= (int)xInstrTip.size())
+ // {
+ // msg_warning()<<"Controlled Instument num "<getDt();
+ // else
+ // {
+ // unsigned int newSensorData = m_currentSensorData + 1;
+
+ // while( m_sensorMotionData[newSensorData][0] < context->getTime() )
+ // {
+ // m_currentSensorData = newSensorData;
+ // newSensorData++;
+ // }
+ // if(newSensorData >= m_sensorMotionData.size())
+ // {
+ // xInstrTip[id] = 0;
+ // }
+ // else
+ // {
+ // xInstrTip[id] += m_sensorMotionData[m_currentSensorData][1];
+ // }
+ // }
+ // }
+ // if (m_RW)
+ // {
+ // xInstrTip[id] -= d_speed.getValue()* context->getDt();
+ // // verif min x :
+ // if ( xInstrTip[id] < 0.0)
+ // {
+ // xInstrTip[id] = 0.0;
+ // m_RW = false;
+ // }
+ // }
+ //}
+
+ ///// The tip of the instrument can not be further than its total length
+ //for (unsigned int i=0; i m_instrumentsList[i]->getRestTotalLength() )
+ // xInstrTip[i] = m_instrumentsList[i]->getRestTotalLength();
applyInterventionalRadiologyController();
}
diff --git a/src/BeamAdapter/component/engine/SteerableCatheter.h b/src/BeamAdapter/component/engine/SteerableCatheter.h
index 17ceebb98..27819a66a 100644
--- a/src/BeamAdapter/component/engine/SteerableCatheter.h
+++ b/src/BeamAdapter/component/engine/SteerableCatheter.h
@@ -104,9 +104,10 @@ class SteerableCatheter : public WireRestShape
/// Bring inherited attributes and function in the current lookup context.
/// otherwise any access to the base::attribute would require
/// the "this->" approach.
- using Inherit1::d_spireDiameter;
- using Inherit1::d_length;
- using Inherit1::d_straightLength;
+ Data d_length;
+ Data d_straightLength;
+ Data d_spireDiameter;
+ Data d_spireHeight;
///////////////////////////////////////////////////////////////////////////
};
diff --git a/src/BeamAdapter/component/engine/WireRestShape.cpp b/src/BeamAdapter/component/engine/WireRestShape.cpp
index dacba394e..34572378a 100644
--- a/src/BeamAdapter/component/engine/WireRestShape.cpp
+++ b/src/BeamAdapter/component/engine/WireRestShape.cpp
@@ -52,7 +52,7 @@ using namespace sofa::defaulttype;
///
////////////////////////////////////////////////////////////////////////////////////////////////////
-static int WireRestShapeClass = core::RegisterObject("Describe the shape functions on multiple segments using curvilinear abscissa")
+const int WireRestShapeClass = core::RegisterObject("Describe the shape functions on multiple segments using curvilinear abscissa")
.add< WireRestShape >(true)
;
diff --git a/src/BeamAdapter/component/engine/WireRestShape.h b/src/BeamAdapter/component/engine/WireRestShape.h
index 8c77de028..2ab349ce9 100644
--- a/src/BeamAdapter/component/engine/WireRestShape.h
+++ b/src/BeamAdapter/component/engine/WireRestShape.h
@@ -33,11 +33,11 @@
#include
-#include
+#include
+
#include
#include
#include
-#include
#include
namespace sofa::component::engine
@@ -48,14 +48,17 @@ namespace _wirerestshape_
using sofa::core::topology::TopologyContainer;
using sofa::component::topology::container::dynamic::EdgeSetTopologyModifier;
-using sofa::component::topology::mapping::Edge2QuadTopologicalMapping;
using sofa::core::loader::MeshLoader;
+using namespace sofa::beamadapter;
+
/**
* \class WireRestShape
* \brief Describe the shape functions on multiple segments
- *
- * Describe the shape functions on multiple segments using curvilinear abscissa
+ *
+ * Describe the full shape of a Wire with a given set of @sa BaseRodSectionMaterial. The wire is discretized by a set of beams (given by the keyPoints and the relatives Beam density)
+ * @sa d_keyPoints and @d_density are computed by method @sa initLengths using the set of rod sections description.
+ * This component compute the beam discretization and the shape functions on multiple segments using curvilinear abscissa.
*/
template
class WireRestShape : public core::objectmodel::BaseObject
@@ -69,8 +72,6 @@ class WireRestShape : public core::objectmodel::BaseObject
using Vec3 = sofa::type::Vec<3, Real>;
using Quat = sofa::type::Quat;
- using BeamSection = sofa::beamadapter::BeamSection;
-
/**
* @brief Default Constructor.
*/
@@ -82,25 +83,19 @@ class WireRestShape : public core::objectmodel::BaseObject
virtual ~WireRestShape() = default;
/////////////////////////// Inherited from BaseObject //////////////////////////////////////////
- void parse(core::objectmodel::BaseObjectDescription* arg) override;
void init() override ;
-
- void draw(const core::visual::VisualParams * vparams) override ;
-
- /////////////////////////// Methods of WireRestShape //////////////////////////////////////////
- /// For coils: a part of the coil instrument can be brokenIn2 (by default the point of release is the end of the straight length)
- Real getReleaseCurvAbs() const {return d_straightLength.getValue();}
+ /////////////////////////// Methods of WireRestShape //////////////////////////////////////////
/// This function is called by the force field to evaluate the rest position of each beam
void getRestTransformOnX(Transform &global_H_local, const Real &x);
/// This function gives the Young modulus and Poisson's coefficient of the beam depending on the beam position
- void getYoungModulusAtX(const Real& x_curv, Real& youngModulus, Real& cPoisson);
+ void getYoungModulusAtX(const Real& x_curv, Real& youngModulus, Real& cPoisson) const;
/// This function gives the mass density and the BeamSection data depending on the beam position
- void getInterpolationParam(const Real& x_curv, Real &_rho, Real &_A, Real &_Iy , Real &_Iz, Real &_Asy, Real &_Asz, Real &_J);
+ void getInterpolationParam(const Real& x_curv, Real &_rho, Real &_A, Real &_Iy , Real &_Iz, Real &_Asy, Real &_Asz, Real &_J) const;
/**
* This function provides a type::vector with the curviliar abscissa of the noticeable point(s)
@@ -110,77 +105,48 @@ class WireRestShape : public core::objectmodel::BaseObject
/// Functions enabling to load and use a geometry given from OBJ external file
- void initRestConfig();
- void getRestPosNonProcedural(Real& abs, Coord &p);
void computeOrientation(const Vec3& AB, const Quat& Q, Quat &result);
- void initFromLoader();
- bool checkTopology();
-
- [[nodiscard]] bool fillTopology();
+
+
Real getLength() ;
- void getCollisionSampling(Real &dx, const Real &x_curv) ;
+ void getCollisionSampling(Real &dx, const Real &x_curv);
void getNumberOfCollisionSegment(Real &dx, unsigned int &numLines) ;
- //TODO(dmarchal 2017-05-17) Please specify who and when it will be done either a time after wich
- //we can remove the todo.
- // todo => topological change !
- void releaseWirePart();
- void rotateFrameForAlignX(const Quat &input, Vec3 &x, Quat &output);
+
+ /////////////////////////// Deprecated Methods //////////////////////////////////////////
+
+ /// For coils: a part of the coil instrument can be brokenIn2 (by default the point of release is the end of the straight length)
+ Real getReleaseCurvAbs() const {
+ msg_warning() << "Releasing catheter or brokenIn2 mode is not anymore supported. Feature has been removed after release v23.06";
+ return 0.0;
+ }
+
+ void releaseWirePart() {
+ msg_warning() << "Releasing catheter or brokenIn2 mode is not anymore supported. Feature has been removed after release v23.06";
+ }
+
+protected:
+ /// Internal method to init Lengths vector @sa d_keyPoints using the length of each materials @sa l_sectionMaterials.
+ void initLengths();
+ /// Internal method to init Edge Topology @sa _topology using the list of materials @sa l_sectionMaterials. Returns false if init can't be performed.
+ bool initTopology();
public:
- /// Analitical creation of wire shape...
- Data d_isAProceduralShape;
- Data d_nonProceduralScale;
- Data d_length;
- Data d_straightLength;
- Data d_spireDiameter;
- Data d_spireHeight;
Data > d_density;
Data > d_keyPoints;
- Data< int > d_numEdges;
- Data > d_numEdgesCollis;
-
- /// User Data about the Young modulus
- Data d_poissonRatio;
- Data d_youngModulus1;
- Data d_youngModulus2;
-
- /// Radius
- Data d_radius1;
- Data d_radius2;
- Data d_innerRadius1;
- Data d_innerRadius2;
-
- Data d_massDensity1;
- Data d_massDensity2;
-
- /// broken in 2 case
- Data d_brokenIn2;
- Data d_drawRestShape;
-
-private:
- /// Data required for the File loading
- type::vector m_localRestPositions;
- type::vector m_localRestTransforms;
- type::vector m_curvAbs ;
- double m_absOfGeometry {0};
- BeamSection beamSection1;
- BeamSection beamSection2;
+ /// Vector or links to the Wire section material. The order of the linked material will define the WireShape structure.
+ MultiLink, BaseRodSectionMaterial, BaseLink::FLAG_STOREPATH | BaseLink::FLAG_STRONGLINK> l_sectionMaterials;
+private:
/// Link to be set to the topology container in the component graph.
SingleLink, TopologyContainer, BaseLink::FLAG_STOREPATH | BaseLink::FLAG_STRONGLINK> l_topology;
/// Pointer to the topology container, should be set using @sa l_topology, otherwise will search for one in current Node.
TopologyContainer* _topology{ nullptr };
/// Pointer to the topology modifier. Will be set at init by searching one in @sa _topology context.
EdgeSetTopologyModifier* edgeMod{ nullptr };
-
- /// Link to be set to the topology container in the component graph.
- SingleLink, MeshLoader, BaseLink::FLAG_STOREPATH | BaseLink::FLAG_STRONGLINK> l_loader;
- /// Pointer to the MeshLoader, should be set using @sa l_loader, otherwise will search for one in current Node.
- MeshLoader* loader{ nullptr };
};
diff --git a/src/BeamAdapter/component/engine/WireRestShape.inl b/src/BeamAdapter/component/engine/WireRestShape.inl
index dda4ad4b5..1ae3484e6 100644
--- a/src/BeamAdapter/component/engine/WireRestShape.inl
+++ b/src/BeamAdapter/component/engine/WireRestShape.inl
@@ -35,12 +35,11 @@
#include
#include
-#include
#include
#include
-#define EPSILON 0.0000000001
+#define EPSILON 0.0001
#define VERIF 1
namespace sofa::component::engine
@@ -57,83 +56,15 @@ using sofa::core::objectmodel::BaseContext ;
* @brief Default Constructor.
*/
template
-WireRestShape::WireRestShape() :
- //TODO(dmarchal 2017-05-17) not sure that procedural & nonProceduralScale are very understandable name...are they exclusives ?
- //if so have look in my comment in the init section.
- d_isAProceduralShape( initData(&d_isAProceduralShape,(bool)true,"isAProceduralShape","is the guidewire shape mathemetically defined ?") )
- , d_nonProceduralScale( initData ( &d_nonProceduralScale, (Real)1.0, "nonProceduralScale", "scale of the model defined by file" ) )
- , d_length(initData(&d_length, (Real)1.0, "length", "total length of the wire instrument"))
- , d_straightLength(initData(&d_straightLength, (Real)0.0, "straightLength", "length of the initial straight shape"))
- , d_spireDiameter(initData(&d_spireDiameter, (Real)0.1, "spireDiameter", "diameter of the spire"))
- , d_spireHeight(initData(&d_spireHeight, (Real)0.01, "spireHeight", "height between each spire"))
- , d_density(initData(&d_density, "densityOfBeams", "density of beams between key points"))
- , d_keyPoints(initData(&d_keyPoints,"keyPoints","key points of the shape (curv absc)"))
- , d_numEdges(initData(&d_numEdges, 10, "numEdges","number of Edges for the visual model"))
- , d_numEdgesCollis(initData(&d_numEdgesCollis,"numEdgesCollis", "number of Edges for the collision model" ))
- , d_poissonRatio(initData(&d_poissonRatio,(Real)0.49,"poissonRatio","Poisson Ratio"))
- , d_youngModulus1(initData(&d_youngModulus1,(Real)5000,"youngModulus","Young Modulus"))
- , d_youngModulus2(initData(&d_youngModulus2,(Real)3000,"youngModulusExtremity","youngModulus for beams at the extremity\nonly if not straight"))
- , d_radius1(initData(&d_radius1,(Real)1.0f,"radius","radius"))
- , d_radius2(initData(&d_radius2,(Real)1.0f,"radiusExtremity","radius for beams at the extremity\nonly if not straight"))
- , d_innerRadius1(initData(&d_innerRadius1,(Real)0.0f,"innerRadius","inner radius if it applies"))
- , d_innerRadius2(initData(&d_innerRadius2,(Real)0.0f,"innerRadiusExtremity","inner radius for beams at the extremity\nonly if not straight"))
- , d_massDensity1(initData(&d_massDensity1,(Real)1.0,"massDensity", "Density of the mass (usually in kg/m^3)" ))
- , d_massDensity2(initData(&d_massDensity2,(Real)1.0,"massDensityExtremity", "Density of the mass at the extremity\nonly if not straight" ))
- , d_brokenIn2(initData(&d_brokenIn2, (bool)false, "brokenIn2", ""))
- , d_drawRestShape(initData(&d_drawRestShape, (bool)false, "draw", "draw rest shape"))
- , l_topology(initLink("topology", "link to the topology container"))
- , l_loader(initLink("loader", "link to the MeshLoader"))
+WireRestShape::WireRestShape()
+ : d_density(initData(&d_density, "densityOfBeams", "density of beams between key points"))
+ , d_keyPoints(initData(&d_keyPoints,"keyPoints","key points of the shape (curv absc)"))
+ , l_sectionMaterials(initLink("wireMaterials", "link to Wire Section Materials (to be ordered according to the instrument, from handle to tip)"))
+ , l_topology(initLink("topology", "link to the topology container"))
{
- d_spireDiameter.setGroup("Procedural");
- d_spireHeight.setGroup("Procedural");
-}
-
-template
-void WireRestShape::rotateFrameForAlignX(const Quat &input, Vec3 &x, Quat &output)
-{
- x.normalize();
- Vec3 x0=input.inverseRotate(x);
-
- Real cTheta=x0[0];
- Real theta;
- if (cTheta>(1-EPSILON))
- {
- output = input;
- }
- else
- {
- theta=acos(cTheta);
- // axis of rotation
- Vec3 dw(0,-x0[2],x0[1]);
- dw.normalize();
-
- // computation of the rotation
- Quat inputRoutput;
- inputRoutput.axisToQuat(dw, theta);
- output=input*inputRoutput;
- }
}
-template
-void WireRestShape::parse(core::objectmodel::BaseObjectDescription* args)
-{
- const char* arg = args->getAttribute("procedural") ;
- if(arg)
- {
- msg_warning() << "The attribute 'procedural' has been renamed into 'isAProceduralShape'. " << msgendl
- << "To remove this warning you need to update your scene and replace 'procedural' with 'isAProceduralShape'" ;
-
- /// As arg is owned by the "procedural" attribute it cannot be removed before
- /// being copied in the "isAProceduralShape". So please keep the ordering of the
- /// two following functions.
- args->setAttribute("isAProceduralShape", arg) ;
- args->removeAttribute("procedural") ;
-
- }
-
- Inherit1::parse(args) ;
-}
template
void WireRestShape::init()
@@ -161,146 +92,95 @@ void WireRestShape::init()
return;
}
- if (!d_isAProceduralShape.getValue())
- {
- // Get meshLoader, check first if loader has been set using link. Otherwise will search in current context.
- loader = l_loader.get();
-
- if (!loader)
- this->getContext()->get(loader);
- if (!loader) {
- msg_error() << "Cannot find a mesh loader. Please insert a MeshObjLoader in the same node or use l_loader to specify the path in the scene graph.";
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
- return;
- }
- else
- {
- msg_info() << "Found a mesh with " << loader->d_edges.getValue().size() << " edges";
- initFromLoader();
- }
- }
- else
+ if (l_sectionMaterials.empty())
{
- if (!fillTopology())
- {
- msg_error() << "Error while trying to fill the associated topology, setting the state to Invalid";
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
- return;
- }
+ msg_error() << "No BaseRodSectionMaterial set. At least one material should be set and link using wireMaterials.";
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
+ return;
}
+
+ ////////////////////////////////////////////////////////
+ ////////// keyPoint list and Density Assignement ///////
+ ////////////////////////////////////////////////////////
+
+ initLengths();
+
// Get pointer to the topology Modifier (for topological changes)
_topology->getContext()->get(edgeMod);
-
if (edgeMod == nullptr)
{
msg_warning() << "No EdgeSetTopologyModifier found in the same node as the topology container: " << _topology->getName() << ". This wire won't support topological changes.";
}
+ initTopology();
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Valid);
+ msg_info() << "WireRestShape end init";
+}
- ////////////////////////////////////////////////////////
- ////////// keyPoint list and Density Assignement ///////
- ////////////////////////////////////////////////////////
- auto keyPointList = sofa::helper::getWriteOnlyAccessor(d_keyPoints);
- if(!keyPointList.size())
- {
- keyPointList.push_back(0.0);
- if(d_straightLength.getValue()>= 0.001*this->d_length.getValue() && d_straightLength.getValue() <= 0.999*d_length.getValue())
- keyPointList.push_back(d_straightLength.getValue());
- keyPointList.push_back(d_length.getValue());
- }
- if( d_density.getValue().size() != keyPointList.size()-1)
- {
- auto densityList = sofa::helper::getWriteOnlyAccessor(d_density);
+template
+void WireRestShape::initLengths()
+{
+ auto keyPointList = sofa::helper::getWriteOnlyAccessor(d_keyPoints);
+ auto densityList = sofa::helper::getWriteOnlyAccessor(d_density);
- if(densityList.size() > keyPointList.size()-1 )
- densityList.resize(keyPointList.size()-1);
- else
- {
- densityList.clear();
-
- if(d_straightLength.getValue()>= 0.001*this->d_length.getValue() )
- {
- int numNodes = (int) floor(5.0*d_straightLength.getValue() / d_length.getValue() );
- densityList.push_back(numNodes);
- }
- if( d_straightLength.getValue() <= 0.999*d_length.getValue())
- {
- int numNodes = (int) floor(20.0*(1.0 - d_straightLength.getValue() / d_length.getValue()) );
- densityList.push_back(numNodes);
- }
- }
- }
+ keyPointList.resize(l_sectionMaterials.size() + 1);
+ keyPointList[0] = Real(0.0);
- if(!d_numEdgesCollis.getValue().size())
+ densityList.resize(l_sectionMaterials.size());
+
+ for (unsigned int i = 0; i < l_sectionMaterials.size(); ++i)
{
- auto densityCol = sofa::helper::getWriteOnlyAccessor(d_numEdgesCollis);
- densityCol.resize(keyPointList.size()-1);
- for (unsigned int i=0; igetLength();
+ densityList[i] = rodSection->getNbCollisionEdges();
}
-
- msg_info() <<"WireRestShape end init" ;
-
- // Prepare beam sections
- double r = this->d_radius1.getValue();
- double rInner = this->d_innerRadius1.getValue();
- this->beamSection1._r = r;
- this->beamSection1._rInner = rInner;
- this->beamSection1._Iz = M_PI*(r*r*r*r - rInner*rInner*rInner*rInner)/4.0;
- this->beamSection1._Iy = this->beamSection1._Iz ;
- this->beamSection1._J = this->beamSection1._Iz + this->beamSection1._Iy;
- this->beamSection1._A = M_PI*(r*r - rInner*rInner);
- this->beamSection1._Asy = 0.0;
- this->beamSection1._Asz = 0.0;
-
- r = this->d_radius2.getValue();
- rInner = this->d_innerRadius2.getValue();
- this->beamSection2._r = r;
- this->beamSection2._rInner = rInner;
- this->beamSection2._Iz = M_PI*(r*r*r*r - rInner*rInner*rInner*rInner)/4.0;
- this->beamSection2._Iy = this->beamSection2._Iz ;
- this->beamSection2._J = this->beamSection2._Iz + this->beamSection2._Iy;
- this->beamSection2._A = M_PI*(r*r - rInner*rInner);
- this->beamSection2._Asy = 0.0;
- this->beamSection2._Asz = 0.0;
-
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Valid);
}
template
-void WireRestShape::releaseWirePart(){
-
- d_brokenIn2.setValue(true);
+bool WireRestShape::initTopology()
+{
+ /// fill topology :
+ _topology->clear();
+ _topology->cleanup();
- if ( edgeMod == nullptr )
+ const type::vector& keyPts = d_keyPoints.getValue();
+ if (l_sectionMaterials.size() != keyPts.size() - 1)
{
- msg_error() << "no edgeSetModifier in the node -> cannot do the topological change";
- return;
+ msg_error() << "Wrong number of inputs. Component can't be init. Number of input materials: " << l_sectionMaterials.size() << ", should be equal to keyPointList.size()-1. keyPointList.size() is equal to: " << keyPts.size();
+ return false;
}
- ///////// remove the edge that is cut //////
- for ( sofa::Size i=0; i<_topology->getNbPoints(); i++)
- {
- if( _topology->getPX(i) > this->getReleaseCurvAbs() + EPSILON )
- {
- type::vector edge_remove;
- edge_remove.push_back( i-1 );
-
- msg_info() << "releaseWirePart() -> remove edge number "<< i ;
-
- edgeMod->removeEdges(edge_remove,false); // remove the single edge and do not remove any point...
+
+ Real prev_length = 0.0;
+ int prev_edges = 0;
+ int startPtId = 0;
+ for (auto i = 0; i < l_sectionMaterials.size(); ++i)
+ {
+ // Add topology of the material
+ int nbrVisuEdges = l_sectionMaterials.get(i)->getNbVisualEdges();
+ Real length = fabs(keyPts[i + 1] - keyPts[i]);
+ Real dx = length / nbrVisuEdges;
+
+ // add points from the material
+ for (int i = startPtId; i < nbrVisuEdges + 1; i++) {
+ _topology->addPoint(prev_length + i * dx, 0, 0);
+ }
- msg_info() << "WireRestShape _topology name="<<_topology->getName()<<" - numEdges ="<<_topology->getNbEdges() ;
-
- return;
+ // add segments from the material
+ for (int i = prev_edges; i < prev_edges + nbrVisuEdges; i++) {
+ _topology->addEdge(i, i + 1);
}
+
+ prev_length = length;
+ prev_edges = nbrVisuEdges;
+ startPtId = 1; // Assume the last point of mat[n] == first point of mat[n+1]
}
- dmsg_info() <<" Wire Part is brokenIn2... should implement a topo change !" ;
+ return true;
}
@@ -308,64 +188,55 @@ template
void WireRestShape::getSamplingParameters(type::vector& xP_noticeable,
type::vector& nbP_density) const
{
-
- xP_noticeable.clear();
- nbP_density.clear();
-
- if (d_brokenIn2.getValue())
- {
- for (unsigned int i=0; i getReleaseCurvAbs() )
- break;
- xP_noticeable.push_back(x);
- nbP_density.push_back(d_density.getValue()[i]);
- }
- xP_noticeable.push_back( getReleaseCurvAbs());
-
- dmsg_info() <<"getSamplingParameters brokenIn2 detected - return xP_noticeable ="<& keyPts = d_keyPoints.getValue();
+
+ // verify that size of number of materials == size of keyPoints-1
+ if (l_sectionMaterials.size() != keyPts.size() - 1)
+ {
+ msg_error() << "Problem size of number of materials: " << l_sectionMaterials.size()
+ << " != size of keyPoints-1 " << keyPts.size()-1
+ << ". Returning default values.";
+ numLines = 20;
+ dx = totalLength / numLines;
return;
}
-
-
- for (unsigned int i=1; id_keyPoints.getValue().size(); i++)
+
+ // Check in which section x_used belongs to and get access to this section material
+ for (auto i = 1; i< keyPts.size(); ++i)
{
- if( x_used < this->d_keyPoints.getValue()[i] )
+ if (x_used <= keyPts[i])
{
- numLines = (unsigned int)d_numEdgesCollis.getValue()[i-1];
- dx=(this->d_keyPoints.getValue()[i] - this->d_keyPoints.getValue()[i-1])/numLines;
+ numLines = l_sectionMaterials.get(i-1)->getNbCollisionEdges();
+
+ Real length = fabs(keyPts[i] - keyPts[i-1]);
+ dx = length / numLines;
return;
}
}
- dx=d_length.getValue()/20;
- msg_error() << " problem is getCollisionSampling : x_curv "<::getRestTransformOnX(Transform &global_H_local, co
{
Real x_used = x - EPSILON;
- if(x_used>d_length.getValue())
- x_used=d_length.getValue();
-
- if(x_used<0.0)
- x_used=0.0;
-
- if( x_used < d_straightLength.getValue())
- {
- global_H_local.set(Vec3(x_used, 0.0, 0.0 ), Quat());
- return;
- }
-
- if(d_isAProceduralShape.getValue())
- {
- Real projetedLength = d_spireDiameter.getValue()*M_PI;
- Real lengthSpire=sqrt(d_spireHeight.getValue()*d_spireHeight.getValue() + projetedLength*projetedLength );
- // angle in the z direction
- Real phi= atan(d_spireHeight.getValue()/projetedLength);
-
- Quat Qphi;
- Qphi.axisToQuat(Vec3(0,0,1),phi);
-
- // spire angle (if theta=2*PI, there is a complete spire between startx and x_used)
- Real lengthCurve= x_used-d_straightLength.getValue();
- Real numSpire=lengthCurve/lengthSpire;
- Real theta= 2*M_PI*numSpire;
-
- // computation of the Quat
- Quat Qtheta;
- Qtheta.axisToQuat(Vec3(0,1,0),theta);
- Quat newSpireQuat = Qtheta*Qphi;
-
-
- // computation of the position
- Real radius=d_spireDiameter.getValue()/2.0;
- Vec3 PosEndCurve(radius*sin(theta), numSpire*d_spireHeight.getValue(), radius*(cos(theta)-1) );
- Vec3 SpirePos=PosEndCurve + Vec3(d_straightLength.getValue(),0,0);
-
- global_H_local.set(SpirePos,newSpireQuat);
- }
- else
- {
- x_used = x_used - d_straightLength.getValue();
- x_used = x_used/(d_length.getValue()-d_straightLength.getValue()) * m_absOfGeometry;
-
- Coord p;
- this->getRestPosNonProcedural(x_used,p);
- Vec3 PosEndCurve = p.getCenter();
- Quat ExtremityQuat = p.getOrientation();
- Vec3 ExtremityPos = PosEndCurve + Vec3(d_straightLength.getValue(),0,0);
-
- global_H_local.set(ExtremityPos,ExtremityQuat);
- }
-}
-
-
-template
-void WireRestShape::getYoungModulusAtX(const Real& x_curv, Real& youngModulus, Real& cPoisson)
-{
- //Initialization
- Real _E1, _E2;
- youngModulus = 0.0;
- cPoisson = 0.0;
-
- //Get the two possible values of the Young modulus
- _E1 = this->d_youngModulus1.getValue();
- _E2 = this->d_youngModulus2.getValue();
-
- //Get User data
- cPoisson = this->d_poissonRatio.getValue();
-
- //Depending on the position of the beam, determine the Young modulus
- if(x_curv <= this->d_straightLength.getValue())
- {
- youngModulus = _E1;
- }
- else
- {
- if(_E2 == 0.0)
- youngModulus = _E1;
- else
- youngModulus = _E2;
+ const Real totalLength = this->getLength();
+ if (x_used > totalLength) {
+ x_used = totalLength;
}
- return;
-}
-
-
-template
-void WireRestShape::getInterpolationParam(const Real& x_curv, Real &_rho, Real &_A, Real &_Iy , Real &_Iz, Real &_Asy, Real &_Asz, Real &_J)
-{
- if(x_curv <= this->d_straightLength.getValue())
- {
- if(d_massDensity1.isSet())
- _rho = d_massDensity1.getValue();
-
- if(d_radius1.isSet())
- {
- _A =beamSection1._A;
- _Iy =beamSection1._Iy;
- _Iz =beamSection1._Iz;
- _Asy =beamSection1._Asy;
- _Asz =beamSection1._Asz;
- _J =beamSection1._J;
- }
+ else if (x_used < 0.0) {
+ x_used = 0.0;
}
- else
+
+ const type::vector& keyPts = d_keyPoints.getValue();
+ for (auto i = 1; i < keyPts.size(); ++i)
{
- if(d_massDensity2.isSet())
- _rho = d_massDensity2.getValue();
- else if(d_massDensity1.isSet())
- _rho = d_massDensity1.getValue();
-
- if(d_radius2.isSet())
- {
- _A =beamSection2._A;
- _Iy =beamSection2._Iy;
- _Iz =beamSection2._Iz;
- _Asy =beamSection2._Asy;
- _Asz =beamSection2._Asz;
- _J =beamSection2._J;
- }
- else if(d_radius1.isSet())
+ if (x_used <= keyPts[i])
{
- _A =beamSection1._A;
- _Iy =beamSection1._Iy;
- _Iz =beamSection1._Iz;
- _Asy =beamSection1._Asy;
- _Asz =beamSection1._Asz;
- _J =beamSection1._J;
+ return l_sectionMaterials.get(i - 1)->getRestTransformOnX(global_H_local, x_used, keyPts[i - 1]);
}
}
-}
-
-template
-bool WireRestShape::checkTopology()
-{
- if (!loader->d_edges.getValue().size())
- {
- msg_error() << "There is no edges in the topology loaded by " << loader->getName() ;
- return false;
- }
- if (loader->d_triangles.getValue().size())
- {
- msg_error() << "There are triangles in the topology loaded by " << loader->getName() ;
- return false;
- }
-
- if (loader->d_quads.getValue().size())
- {
- msg_error() << "There are quads in the topology loaded by " << loader->getName() ;
- return false;
- }
-
- if (loader->d_polygons.getValue().size())
- {
- msg_error() << "There are polygons in the topology loaded by " << loader->getName() ;
- return false;
- }
-
- //TODO(dmarchal 2017-05-17) when writing a TODO please specify:
- // who will do that
- // when it will be done
- /// \todo check if the topology is like a wire
-
-
- return true;
+ msg_error() << " problem in getRestTransformOnX : x_curv " << x_used << " is not between keyPoints" << d_keyPoints.getValue();
}
template
-bool WireRestShape::fillTopology()
+void WireRestShape::getYoungModulusAtX(const Real& x_curv, Real& youngModulus, Real& cPoisson) const
{
- if (!_topology)
- {
- msg_error() << "Topology is null";
- return false;
- }
-
- const auto length = this->d_length.getValue();
- if (length <= Real(0.0))
- {
- msg_error() << "Length is 0 (or negative), check if d_length has been given or computed.";
- return false;
- }
+ const Real x_used = x_curv - Real(EPSILON);
+ const type::vector& keyPts = d_keyPoints.getValue();
- int nbrEdges = d_numEdges.getValue();
- if (nbrEdges <= 0)
+ // Depending on the position of the beam, determine the corresponding section material and returning its Young modulus
+ for (auto i = 1; i < keyPts.size(); ++i)
{
- msg_warning() << "Number of edges has been set to an invalid value: " << nbrEdges << ". Value should be a positive integer. Setting to default value: 10";
- nbrEdges = 10;
+ if (x_used <= keyPts[i])
+ {
+ return l_sectionMaterials.get(i - 1)->getYoungModulusAtX(youngModulus, cPoisson);
+ }
}
- /// fill topology :
- _topology->clear();
- _topology->cleanup();
-
- Real dx = this->d_length.getValue() / nbrEdges;
-
- /// add points
- for (int i = 0; i < d_numEdges.getValue() + 1; i++)
- _topology->addPoint(i * dx, 0, 0);
-
- /// add segments
- for (int i = 0; i < d_numEdges.getValue(); i++)
- _topology->addEdge(i, i + 1);
-
- return true;
+ msg_error() << " problem in getYoungModulusAtX : x_curv " << x_curv << " is not between keyPoints" << keyPts;
}
-
template
-void WireRestShape::initFromLoader()
+void WireRestShape::getInterpolationParam(const Real& x_curv, Real &_rho, Real &_A, Real &_Iy , Real &_Iz, Real &_Asy, Real &_Asz, Real &_J) const
{
- if (!checkTopology())
- {
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
- return;
- }
-
- // this is... dirty but d_straightLength itself is not relevant for loader-created wire
- // but the code uses it actively and expects it to not be zero.
- if (d_straightLength.getValue() <= 0.0)
- {
- msg_warning() << "straightLength cannot be 0 (or negative...). Setting a minimum value.";
- d_straightLength.setValue(0.0001);
- }
-
- type::vector vertices;
- sofa::core::topology::BaseMeshTopology::SeqEdges edges;
-
- //get the topology position
- auto topoVertices = sofa::helper::getReadAccessor(loader->d_positions);
-
- //copy the topology edges in a local vector
- auto topoEdges = sofa::helper::getReadAccessor(loader->d_edges);
- edges = topoEdges.ref();
-
- /** renumber the vertices **/
- type::vector verticesConnexion; //gives the number of edges connected to a vertex
- for(unsigned int i =0; i < topoVertices.size(); i++)
- verticesConnexion.push_back(2);
-
- for(const auto& ed : edges)
- {
- verticesConnexion[ed[0]]--;
- verticesConnexion[ed[1]]--;
- }
+ const Real x_used = x_curv - Real(EPSILON);
+ const type::vector& keyPts = d_keyPoints.getValue();
- msg_info() << "Successfully compute the vertex connexion" ;
-
- // check for the first corner of the edge
- unsigned int firstIndex = 0;
- bool found = false;
- while((firstIndex < verticesConnexion.size()) && !found)
- {
- if(verticesConnexion[firstIndex] == 1)
- found = true;
- else
- firstIndex++;
- }
-
- if(firstIndex == verticesConnexion.size())
- {
- msg_error() << "The first vertex of the beam structure is not found, probably because of a closed structure" ;
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
- return;
- }
-
- vertices.push_back(topoVertices[firstIndex]);
-
- while(edges.size() > 0)
+ // Check in which section x_used belongs to and get access to this section material
+ for (auto i = 1; i < keyPts.size(); ++i)
{
- auto it = edges.begin();
- auto end = edges.end();
-
- bool notFound = true;
- while (notFound && (it != end))
+ if (x_used <= keyPts[i])
{
- const auto& ed = (*it);
- auto toDel = it;
- it++;
- if(ed[0] == firstIndex)
- {
- vertices.push_back(topoVertices[ed[1]]);
- firstIndex = ed[1];
- edges.erase(toDel);
- notFound = false;
-
- }
- else if(ed[1] == firstIndex)
- {
- vertices.push_back(topoVertices[ed[0]]);
- firstIndex = ed[0];
- edges.erase(toDel);
- notFound = false;
- }
+ return l_sectionMaterials.get(i - 1)->getInterpolationParam(_rho, _A, _Iy, _Iz, _Asy, _Asz, _J);
}
}
- msg_info() << "Successfully computed the topology" ;
-
- m_localRestPositions = vertices;
-
- for(unsigned int i = 0; i < m_localRestPositions.size() - 1; i++)
- m_localRestPositions[i] *= d_nonProceduralScale.getValue();
-
- initRestConfig();
- // TODO epernod 2022-08-05: Init from loader seems quite buggy, need to check if this is still needed and working
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Valid);
-}
-
-
-template
-void WireRestShape::initRestConfig()
-{
- m_curvAbs.clear();
- double tot = 0;
- m_curvAbs.push_back(0);
- Quat input, output;
- input.identity();
- m_localRestTransforms.resize(m_localRestPositions.size());
- m_localRestTransforms[0].setOrigin(Vec3(0,0,0));
- m_localRestTransforms[0].setOrientation(input);
-
- for(unsigned int i = 0; i < m_localRestPositions.size() - 1; i++)
- {
- Vec3 vec = m_localRestPositions[i+1] - m_localRestPositions[i];
- double norm = vec.norm();
- tot += norm;
-
- this->rotateFrameForAlignX(input, vec, output);
-
- input = output;
-
- m_localRestTransforms[i+1].setOrientation(output);
-
- Vec3 localPos = m_localRestPositions[i+1] - m_localRestPositions[0];
-
- m_localRestTransforms[i+1].setOrigin(localPos);
-
- m_curvAbs.push_back(tot);
- }
- m_absOfGeometry = tot;
-
- Real newLength = d_straightLength.getValue() + m_absOfGeometry;
- d_length.setValue(newLength);
-
- msg_info() <<"Length of the loaded shape = "<< m_absOfGeometry << ", total length with straight length = " << newLength ;
-
- if (!fillTopology())
- {
- msg_error() << "Error while trying to fill the associated topology, setting the state to Invalid";
- this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
- return;
- }
+ msg_error() << " problem in getInterpolationParam : x_curv " << x_curv << " is not between keyPoints" << keyPts;
}
-template
-void WireRestShape::getRestPosNonProcedural(Real& abs, Coord &p)
-{
- /*** find the range which includes the "requested" abs ***/
- double startingAbs = 0; unsigned int index = 0;
-
- while ((startingAbs < abs) && (index < m_localRestPositions.size()))
- {
- index++;
- startingAbs = m_curvAbs[index];
- }
-
- /*** OOB ***/
- if(abs > startingAbs)
- {
- msg_error() << "abs = "<
typename WireRestShape::Real WireRestShape::getLength()
{
- if(d_brokenIn2.getValue())
- return d_straightLength.getValue();
- else
- return d_length.getValue();
+ return d_keyPoints.getValue().back();
}
+
template
void WireRestShape::getNumberOfCollisionSegment(Real &dx, unsigned int &numLines)
{
numLines = 0;
- for (unsigned i=0; igetNbCollisionEdges();
}
- dx=d_length.getValue()/numLines;
+ dx = getLength() / numLines;
}
+
template
void WireRestShape::computeOrientation(const Vec3& AB, const Quat& Q, Quat &result)
{
@@ -824,26 +355,6 @@ void WireRestShape::computeOrientation(const Vec3& AB, const Quat& Q,
}
-template
-void WireRestShape::draw(const core::visual::VisualParams* vparams)
-{
- if (!d_drawRestShape.getValue())
- return;
-
- vparams->drawTool()->saveLastState();
- vparams->drawTool()->setLightingEnabled(false);
-
- std::vector< sofa::type::Vec3 > points;
- points.reserve(m_localRestPositions.size());
-
- for (unsigned int i = 0; i < m_localRestPositions.size(); i++)
- {
- points.emplace_back(m_localRestPositions[i][0], m_localRestPositions[i][1], m_localRestPositions[i][2]);
- }
-
- vparams->drawTool()->drawPoints(points, 10, sofa::type::RGBAColor(1, 0.5, 0.5, 1));
- vparams->drawTool()->restoreLastState();
-}
} // namespace _wirerestshape_
using _wirerestshape_::WireRestShape;
diff --git a/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.h b/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.h
index ac681e6b8..d2a0b910f 100644
--- a/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.h
+++ b/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.h
@@ -154,6 +154,8 @@ class AdaptiveBeamMapping : public Mapping
Data> d_segmentsCurvAbs; /*!< (output) the abscissa of each created point on the collision model */
Data d_parallelMapping; /*!< flag to enable parallel internal computation of apply/applyJ */
+ Data d_onlyVisual;
+
SingleLink,
BInterpolation, BaseLink::FLAG_STOREPATH|BaseLink::FLAG_STRONGLINK> l_adaptativebeamInterpolation;
diff --git a/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.inl b/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.inl
index 8a9029c12..f4df8bdd9 100644
--- a/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.inl
+++ b/src/BeamAdapter/component/mapping/AdaptiveBeamMapping.inl
@@ -77,6 +77,7 @@ AdaptiveBeamMapping::AdaptiveBeamMapping(State< In >* from, State< Out
, d_nbPointsPerBeam(initData(&d_nbPointsPerBeam, 0.0, "nbPointsPerBeam", "if non zero, we will adapt the points depending on the discretization, with this num of points per beam (compatible with useCurvAbs)"))
, d_segmentsCurvAbs(initData(&d_segmentsCurvAbs, "segmentsCurvAbs", "the abscissa of each point on the collision model", true, true))
, d_parallelMapping(initData(&d_parallelMapping, false, "parallelMapping", "flag to enable parallel internal computation"))
+ , d_onlyVisual(initData(&d_onlyVisual, (bool)false, "onlyVisual", "Really not mechanical mapping"))
, l_adaptativebeamInterpolation(initLink("interpolation", "Path to the Interpolation component on scene"), interpolation)
, m_inputMapping(nullptr)
, m_isSubMapping(isSubMapping)
@@ -292,6 +293,9 @@ void AdaptiveBeamMapping< TIn, TOut>::apply(const MechanicalParams* mparams, Dat
template
void AdaptiveBeamMapping< TIn, TOut>::applyJ(const core::MechanicalParams* mparams, Data& dOut, const Data& dIn)
{
+ if (d_onlyVisual.getValue())
+ return;
+
SOFA_UNUSED(mparams);
ScopedAdvancedTimer timer("AdaptiveBeamMapping_applyJ");
@@ -377,6 +381,9 @@ void AdaptiveBeamMapping< TIn, TOut>::applyJ(const core::MechanicalParams* mpara
template
void AdaptiveBeamMapping< TIn, TOut>::applyJT(const core::MechanicalParams* mparams, Data& dOut, const Data& dIn)
{
+ if (d_onlyVisual.getValue())
+ return;
+
SOFA_UNUSED(mparams);
ScopedAdvancedTimer timer("AdaptiveBeamMapping_applyJT");
diff --git a/src/BeamAdapter/component/model/BaseRodSectionMaterial.h b/src/BeamAdapter/component/model/BaseRodSectionMaterial.h
new file mode 100644
index 000000000..d023b6a69
--- /dev/null
+++ b/src/BeamAdapter/component/model/BaseRodSectionMaterial.h
@@ -0,0 +1,120 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+#include
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+using sofa::core::loader::MeshLoader;
+
+/**
+ * \class BaseRodSectionMaterial
+ * \brief Base class describing a Rod section which will define a set of beam elements.
+ *
+ * This class provide an api to define a rod/wire section using physical and geometry parameters.
+ * The section will then be modellized by a set of beam elements. Inheriting class should provide the geometry structure:
+ * @sa RodMeshSection to define a rod using a mesh file, @sa RodSpireSection or @sa RodStraightSection to define procedural shapes.
+ * Method @sa initSection and @sa getRestTransformOnX should be overriden to provide the correct creation and interpolation.
+ *
+ * The rod section is described by:
+ * - Topology parameters: vertices and edges @sa d_nbEdgesVisu and @sa d_nbEdgesCollis
+ * - Geometry parameters: radius @sa d_radius, @sa d_innerRadius and length @sa d_length
+ * - Mechanical parameters: @sa d_poissonRatio and @sa d_youngModulus
+ */
+template
+class BaseRodSectionMaterial : public core::objectmodel::BaseObject
+{
+public:
+ SOFA_CLASS(SOFA_TEMPLATE(BaseRodSectionMaterial, DataTypes), core::objectmodel::BaseObject);
+
+ using Coord = typename DataTypes::Coord;
+ using Real = typename Coord::value_type;
+ using Transform = typename sofa::defaulttype::SolidTypes::Transform;
+ using Vec3 = sofa::type::Vec<3, Real>;
+ using Quat = sofa::type::Quat;
+ using Size = sofa::Size;
+
+ /////////////////////////// Inherited from BaseObject //////////////////////////////////////////
+
+ /// Default Constructor
+ BaseRodSectionMaterial();
+
+ /// init method from BaseObject API. Will call internal @see initSection to be overriden by children
+ void init() override;
+
+
+ /////////////////////////// Geometry and physics Getter //////////////////////////////////////////
+
+ /// Returns the number of visual edges of this section. To be set or computed by child.
+ [[nodiscard]] int getNbVisualEdges() const { return d_nbEdgesVisu.getValue(); }
+
+ /// Returns the number of collision edges of this section. To be set or computed by child.
+ [[nodiscard]] int getNbCollisionEdges() const { return d_nbEdgesCollis.getValue(); }
+
+ /// Returns the total length of this section. To be set or computed by child.
+ [[nodiscard]] Real getLength() const { return d_length.getValue(); }
+
+ /// Returns the Young modulus and Poisson's coefficient of this section
+ void getYoungModulusAtX(Real& youngModulus, Real& cPoisson) const;
+
+ /// Returns the mass density and the BeamSection of this section
+ void getInterpolationParam(Real& _rho, Real& _A, Real& _Iy, Real& _Iz, Real& _Asy, Real& _Asz, Real& _J) const;
+
+
+
+ /// This function is called to get the rest position of the beam depending on the current curved abscisse given in parameter
+ virtual void getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start)
+ {
+ SOFA_UNUSED(global_H_local);
+ SOFA_UNUSED(x_used);
+ SOFA_UNUSED(x_start);
+ }
+
+protected:
+ /// Internal method to init the section. to be overidden by child.
+ virtual bool initSection() { return false; }
+
+public:
+ Data d_poissonRatio; ///< Data defining the mehcanical Poisson ratio of this section
+ Data d_youngModulus; ///< Data defining the mehcanical Young Modulus of this section
+ Data d_massDensity; ///< Data defining the mehcanical mass density of this section
+
+ Data d_radius; ///< Data defining the geometry radius of this section
+ Data d_innerRadius; ///< Data defining the geometry internal radius of this section is hollow
+ Data d_length; ///< Data defining the geometry length of this section
+
+ Data d_nbEdgesVisu; ///< Data defining the number of visual edges composing this section
+ Data d_nbEdgesCollis; ///< Data defining the number of collision edges composing this section
+
+private:
+ /// Internal structure to store physical parameter of the a beam section
+ BeamSection beamSection;
+};
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/BaseRodSectionMaterial.inl b/src/BeamAdapter/component/model/BaseRodSectionMaterial.inl
new file mode 100644
index 000000000..ffef90e69
--- /dev/null
+++ b/src/BeamAdapter/component/model/BaseRodSectionMaterial.inl
@@ -0,0 +1,96 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+
+namespace sofa::beamadapter
+{
+
+template
+BaseRodSectionMaterial::BaseRodSectionMaterial()
+ : d_poissonRatio(initData(&d_poissonRatio, (Real)0.49, "poissonRatio", "Poisson Ratio of this section"))
+ , d_youngModulus(initData(&d_youngModulus, (Real)5000, "youngModulus", "Young Modulus of this section"))
+ , d_massDensity(initData(&d_massDensity, (Real)1.0, "massDensity", "Density of the mass (usually in kg/m^3)"))
+ , d_radius(initData(&d_radius, (Real)1.0, "radius", "Full radius of this section"))
+ , d_innerRadius(initData(&d_innerRadius, (Real)0.0, "innerRadius", "Inner radius of this section if hollow"))
+ , d_length(initData(&d_length, (Real)1.0, "length", "Total length of this section"))
+ , d_nbEdgesVisu(initData(&d_nbEdgesVisu, (Size)10, "nbEdgesVisu", "number of Edges for the visual model"))
+ , d_nbEdgesCollis(initData(&d_nbEdgesCollis, (Size)20, "nbEdgesCollis", "number of Edges for the collision model"))
+{
+
+}
+
+
+template
+void BaseRodSectionMaterial::init()
+{
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Loading);
+
+ // Prepare beam sections
+ double r = this->d_radius.getValue();
+ double rInner = this->d_innerRadius.getValue();
+ this->beamSection._r = r;
+ this->beamSection._rInner = rInner;
+ this->beamSection._Iz = M_PI * (r * r * r * r - rInner * rInner * rInner * rInner) / 4.0;
+ this->beamSection._Iy = this->beamSection._Iz;
+ this->beamSection._J = this->beamSection._Iz + this->beamSection._Iy;
+ this->beamSection._A = M_PI * (r * r - rInner * rInner);
+ this->beamSection._Asy = 0.0;
+ this->beamSection._Asz = 0.0;
+
+ // call delegate method to init the section
+ bool res = initSection();
+
+ if (res)
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Valid);
+ else
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
+}
+
+
+template
+void BaseRodSectionMaterial::getYoungModulusAtX(Real& youngModulus, Real& cPoisson) const
+{
+ youngModulus = this->d_youngModulus.getValue();
+ cPoisson = this->d_poissonRatio.getValue();
+}
+
+
+template
+void BaseRodSectionMaterial::getInterpolationParam(Real& _rho, Real& _A, Real& _Iy, Real& _Iz, Real& _Asy, Real& _Asz, Real& _J) const
+{
+ if (d_massDensity.isSet())
+ _rho = d_massDensity.getValue();
+
+ if (d_radius.isSet())
+ {
+ _A = beamSection._A;
+ _Iy = beamSection._Iy;
+ _Iz = beamSection._Iz;
+ _Asy = beamSection._Asy;
+ _Asz = beamSection._Asz;
+ _J = beamSection._J;
+ }
+}
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodMeshSection.cpp b/src/BeamAdapter/component/model/RodMeshSection.cpp
new file mode 100644
index 000000000..58b506884
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodMeshSection.cpp
@@ -0,0 +1,39 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#define SOFA_PLUGIN_BEAMADAPTER_RODMESHSECTION_CPP
+
+#include
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+using namespace sofa::defaulttype;
+
+const int RodMeshSectionClass = core::RegisterObject("Class defining a Rod Section using a MeshLoader and material parameters.")
+ .add< RodMeshSection >(true);
+
+template class SOFA_BEAMADAPTER_API RodMeshSection;
+
+}// namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodMeshSection.h b/src/BeamAdapter/component/model/RodMeshSection.h
new file mode 100644
index 000000000..0c75b9c8e
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodMeshSection.h
@@ -0,0 +1,86 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+using sofa::core::loader::MeshLoader;
+
+/**
+ * \class RodMeshSection
+ * \brief Specialization class of @sa BaseRodSectionMaterial describing a rod section created using a Mesh file
+ *
+ * This class will describe a rod section defined by a mesh file structure using the link @sa l_loader
+ * Method @sa initFromLoader and @sa initRestConfig will define the beam structure using the geometry of the given mesh
+ * as well as the Length. Mechanical parameters are set using the @sa BaseRodSectionMaterial Data
+ * Method @sa getRestTransformOnX will return the current position of the curviline abscisse along the mesh structure.
+ */
+template
+class RodMeshSection : public sofa::beamadapter::BaseRodSectionMaterial
+{
+public:
+ SOFA_CLASS(SOFA_TEMPLATE(RodMeshSection, DataTypes), SOFA_TEMPLATE(BaseRodSectionMaterial, DataTypes));
+
+ /// Default Constructor
+ RodMeshSection();
+
+ /// Override method to get the rest position of the beam. In this implementation, it will interpolate along the loaded mesh geometry
+ void getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start) override;
+
+protected:
+ /// Internal method to init the section. Called by @sa BaseRodSectionMaterial::init() method
+ bool initSection() override;
+
+ /// Internal method called by initSection to init from a linked MeshLoader @sa l_loader
+ bool initFromLoader();
+ /// Internal method called by initFromLoader to compute @sa m_localRestPositions and @sa m_localRestTransforms given the mesh structure
+ void initRestConfig();
+ /// Method to check if the given loader has a edge set structure
+ bool checkLoaderTopology();
+
+ /// Tool method to rotate the input frame @param input given an axis @param x. Result is set in @param output
+ void rotateFrameForAlignX(const Quat& input, Vec3& x, Quat& output);
+
+public:
+ /// Link to be set to the topology container in the component graph.
+ SingleLink, MeshLoader, BaseLink::FLAG_STOREPATH | BaseLink::FLAG_STRONGLINK> l_loader;
+
+private:
+ /// Pointer to the MeshLoader, should be set using @sa l_loader, otherwise will search for one in current Node.
+ MeshLoader* p_loader{ nullptr };
+
+ type::vector m_localRestPositions; ///< rest position of the key points interpolated on the mesh geometry
+ type::vector m_localRestTransforms; ///< rest transform of the key points interpolated on the mesh geometry
+ type::vector m_curvAbs; ///< set of absciss curviline points
+ Real m_absOfGeometry{ 0 }; ///< max curv absciss of this mesh structure
+};
+
+#if !defined(SOFA_PLUGIN_BEAMADAPTER_RODMESHSECTION_CPP)
+extern template class SOFA_BEAMADAPTER_API RodMeshSection;
+#endif
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodMeshSection.inl b/src/BeamAdapter/component/model/RodMeshSection.inl
new file mode 100644
index 000000000..1867c195b
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodMeshSection.inl
@@ -0,0 +1,298 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+#include
+
+#define EPSILON 0.0001
+
+namespace sofa::beamadapter
+{
+
+template
+RodMeshSection::RodMeshSection()
+ : BaseRodSectionMaterial()
+ , l_loader(initLink("loader", "link to the MeshLoader"))
+{
+
+}
+
+
+template
+bool RodMeshSection::initSection()
+{
+ // Get meshLoader, check first if loader has been set using link. Otherwise will search in current context.
+ p_loader = l_loader.get();
+
+ if (!p_loader)
+ this->getContext()->get(p_loader);
+
+ if (!p_loader) {
+ msg_error() << "Cannot find a mesh loader. Please insert a MeshObjLoader in the same node or use l_loader to specify the path in the scene graph.";
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
+ return false;
+ }
+ else
+ {
+ msg_info() << "Found a mesh with " << p_loader->d_edges.getValue().size() << " edges";
+ return initFromLoader();
+ }
+}
+
+
+template
+void RodMeshSection::getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start)
+{
+ Real abs_curr = x_used - x_start;
+ abs_curr = abs_curr /(d_length.getValue()) * m_absOfGeometry;
+
+ Coord p;
+
+ /*** find the range which includes the "requested" abs ***/
+ Real startingAbs = 0;
+ unsigned int index = 0;
+ while ((startingAbs < abs_curr) && (index < m_localRestPositions.size()))
+ {
+ index++;
+ startingAbs = m_curvAbs[index];
+ }
+
+ /*** OOB ***/
+ if (abs_curr > startingAbs)
+ {
+ msg_error() << "Out of bound position requested= " << abs_curr << " with startingAbs = " << startingAbs;
+ return;
+ }
+ else /*** Expected case ***/
+ {
+ const Real alpha = (abs_curr - m_curvAbs[index - 1]) / (m_curvAbs[index] - m_curvAbs[index - 1]);
+ const Real one_minus_alpha = 1 - alpha;
+ const Vec3 result = m_localRestTransforms[index - 1].getOrigin() * one_minus_alpha + m_localRestTransforms[index].getOrigin() * alpha;
+ Quat slerp;
+ slerp.slerp(m_localRestTransforms[index - 1].getOrientation(), m_localRestTransforms[index].getOrientation(), alpha, true);
+ slerp.normalize();
+
+ p.getCenter() = result;
+ p.getOrientation() = slerp;
+ }
+
+ Vec3 PosEndCurve = p.getCenter();
+ Quat ExtremityQuat = p.getOrientation();
+ Vec3 ExtremityPos = PosEndCurve + Vec3(x_start, 0, 0);
+
+ global_H_local.set(ExtremityPos, ExtremityQuat);
+}
+
+
+template
+bool RodMeshSection::initFromLoader()
+{
+ if (!checkLoaderTopology())
+ {
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
+ return false;
+ }
+
+ type::vector vertices;
+ sofa::core::topology::BaseMeshTopology::SeqEdges edges;
+
+ //get the topology position
+ auto topoVertices = sofa::helper::getReadAccessor(p_loader->d_positions);
+
+ //copy the topology edges in a local vector
+ auto topoEdges = sofa::helper::getReadAccessor(p_loader->d_edges);
+ edges = topoEdges.ref();
+
+ /** renumber the vertices **/
+ type::vector verticesConnexion; //gives the number of edges connected to a vertex
+ for (unsigned int i = 0; i < topoVertices.size(); i++)
+ verticesConnexion.push_back(2);
+
+ for (const auto& ed : edges)
+ {
+ verticesConnexion[ed[0]]--;
+ verticesConnexion[ed[1]]--;
+ }
+
+ msg_info() << "Successfully compute the vertex connexion";
+
+ // check for the first corner of the edge
+ unsigned int firstIndex = 0;
+ bool found = false;
+ while ((firstIndex < verticesConnexion.size()) && !found)
+ {
+ if (verticesConnexion[firstIndex] == 1)
+ found = true;
+ else
+ firstIndex++;
+ }
+
+ if (firstIndex == verticesConnexion.size())
+ {
+ msg_error() << "The first vertex of the beam structure is not found, probably because of a closed structure";
+ this->d_componentState.setValue(sofa::core::objectmodel::ComponentState::Invalid);
+ return false;
+ }
+
+ vertices.push_back(topoVertices[firstIndex]);
+
+ while (edges.size() > 0)
+ {
+ auto it = edges.begin();
+ auto end = edges.end();
+
+ bool notFound = true;
+ while (notFound && (it != end))
+ {
+ const auto& ed = (*it);
+ auto toDel = it;
+ it++;
+ if (ed[0] == firstIndex)
+ {
+ vertices.push_back(topoVertices[ed[1]]);
+ firstIndex = ed[1];
+ edges.erase(toDel);
+ notFound = false;
+
+ }
+ else if (ed[1] == firstIndex)
+ {
+ vertices.push_back(topoVertices[ed[0]]);
+ firstIndex = ed[0];
+ edges.erase(toDel);
+ notFound = false;
+ }
+ }
+ }
+
+ msg_info() << "Successfully computed the topology";
+
+ m_localRestPositions = vertices;
+
+ initRestConfig();
+
+ return true;
+}
+
+
+template
+void RodMeshSection::initRestConfig()
+{
+ m_curvAbs.clear();
+ Real tot = 0;
+ m_curvAbs.push_back(0);
+ Quat input, output;
+ input.identity();
+ m_localRestTransforms.resize(m_localRestPositions.size());
+ m_localRestTransforms[0].setOrigin(Vec3(0, 0, 0));
+ m_localRestTransforms[0].setOrientation(input);
+
+ for (unsigned int i = 0; i < m_localRestPositions.size() - 1; i++)
+ {
+ Vec3 vec = m_localRestPositions[i + 1] - m_localRestPositions[i];
+ Real norm = vec.norm();
+ tot += norm;
+
+ this->rotateFrameForAlignX(input, vec, output);
+
+ input = output;
+
+ m_localRestTransforms[i + 1].setOrientation(output);
+
+ Vec3 localPos = m_localRestPositions[i + 1] - m_localRestPositions[0];
+
+ m_localRestTransforms[i + 1].setOrigin(localPos);
+
+ m_curvAbs.push_back(tot);
+ }
+ m_absOfGeometry = tot;
+
+ d_length.setValue(m_absOfGeometry);
+
+ msg_info() << "Length of the loaded shape = " << m_absOfGeometry;
+}
+
+
+template
+bool RodMeshSection::checkLoaderTopology()
+{
+ if (!p_loader->d_edges.getValue().size())
+ {
+ msg_error() << "There is no edges in the topology loaded by " << p_loader->getName();
+ return false;
+ }
+
+ if (p_loader->d_triangles.getValue().size())
+ {
+ msg_error() << "There are triangles in the topology loaded by " << p_loader->getName();
+ return false;
+ }
+
+ if (p_loader->d_quads.getValue().size())
+ {
+ msg_error() << "There are quads in the topology loaded by " << p_loader->getName();
+ return false;
+ }
+
+ if (p_loader->d_polygons.getValue().size())
+ {
+ msg_error() << "There are polygons in the topology loaded by " << p_loader->getName();
+ return false;
+ }
+
+ return true;
+}
+
+
+template
+void RodMeshSection::rotateFrameForAlignX(const Quat& input, Vec3& x, Quat& output)
+{
+ x.normalize();
+ Vec3 x0 = input.inverseRotate(x);
+
+ Real cTheta = x0[0];
+ Real theta;
+ if (cTheta > (1 - EPSILON))
+ {
+ output = input;
+ }
+ else
+ {
+ theta = acos(cTheta);
+ // axis of rotation
+ Vec3 dw(0, -x0[2], x0[1]);
+ dw.normalize();
+
+ // computation of the rotation
+ Quat inputRoutput;
+ inputRoutput.axisToQuat(dw, theta);
+
+ output = input * inputRoutput;
+ }
+}
+
+
+
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodSpireSection.cpp b/src/BeamAdapter/component/model/RodSpireSection.cpp
new file mode 100644
index 000000000..04b64d83a
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodSpireSection.cpp
@@ -0,0 +1,39 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#define SOFA_PLUGIN_BEAMADAPTER_RODSPIRESECTION_CPP
+
+#include
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+using namespace sofa::defaulttype;
+
+const int RodSpireSectionClass = core::RegisterObject("Class defining a rod spire section, defining material and geometry parameters.")
+ .add< RodSpireSection >(true);
+
+template class SOFA_BEAMADAPTER_API RodSpireSection;
+
+}// namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodSpireSection.h b/src/BeamAdapter/component/model/RodSpireSection.h
new file mode 100644
index 000000000..91a793d88
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodSpireSection.h
@@ -0,0 +1,65 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+using sofa::core::loader::MeshLoader;
+
+/**
+ * \class RodSpireSection
+ * \brief Specialization class of @sa BaseRodSectionMaterial describing a rod spire section.
+ *
+ * This class will describe a rod spire section using spire diameter and height between each spire. Length and mechanical
+ * parameters are the same as @sa BaseRodSectionMaterial Data
+ * Method @sa getRestTransformOnX will return the current position of the curviline abscisse along the spire.
+ */
+template
+class RodSpireSection : public sofa::beamadapter::BaseRodSectionMaterial
+{
+public:
+ SOFA_CLASS(SOFA_TEMPLATE(RodSpireSection, DataTypes), SOFA_TEMPLATE(BaseRodSectionMaterial, DataTypes));
+
+ /// Default Constructor
+ RodSpireSection();
+
+ /// Override method to get the rest position of the beam. In this implementation, it will compute the current position given the spire parameters
+ void getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start) override;
+
+protected:
+ /// Internal method to init the section. Called by @sa BaseRodSectionMaterial::init() method
+ bool initSection() override;
+
+public:
+ Data d_spireDiameter; ///< Data defining the diameter of the spire
+ Data d_spireHeight; ///< Data defining the height between each spire
+};
+
+#if !defined(SOFA_PLUGIN_BEAMADAPTER_RODSPIRESECTION_CPP)
+extern template class SOFA_BEAMADAPTER_API RodSpireSection;
+#endif
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodSpireSection.inl b/src/BeamAdapter/component/model/RodSpireSection.inl
new file mode 100644
index 000000000..cb9ff3230
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodSpireSection.inl
@@ -0,0 +1,97 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+template
+RodSpireSection::RodSpireSection()
+ : BaseRodSectionMaterial()
+ , d_spireDiameter(initData(&d_spireDiameter, (Real)0.1, "spireDiameter", "diameter of the spire"))
+ , d_spireHeight(initData(&d_spireHeight, (Real)0.01, "spireHeight", "height between each spire"))
+{
+
+}
+
+
+template
+bool RodSpireSection::initSection()
+{
+ const auto length = this->d_length.getValue();
+ if (length <= Real(0.0))
+ {
+ msg_error() << "Length is 0 (or negative), check if d_length has been given or well computed.";
+ return false;
+ }
+
+ if (int nbrEdgesVisu = d_nbEdgesVisu.getValue() <= 0)
+ {
+ msg_warning() << "Number of visual edges has been set to an invalid value: " << nbrEdgesVisu << ". Value should be a positive integer. Setting to default value: 10";
+ d_nbEdgesVisu.setValue(10);
+ }
+
+
+ if (int nbEdgesCollis = d_nbEdgesCollis.getValue() <= 0)
+ {
+ msg_warning() << "Number of collision edges has been set to an invalid value: " << nbEdgesCollis << ". Value should be a positive integer. Setting to default value: 20";
+ d_nbEdgesCollis.setValue(10);
+ }
+
+ return true;
+}
+
+
+template
+void RodSpireSection::getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start)
+{
+ Real projetedLength = d_spireDiameter.getValue() * M_PI;
+ Real lengthSpire = sqrt(d_spireHeight.getValue() * d_spireHeight.getValue() + projetedLength * projetedLength);
+ // angle in the z direction
+ Real phi = atan(d_spireHeight.getValue() / projetedLength);
+
+ Quat Qphi;
+ Qphi.axisToQuat(Vec3(0, 0, 1), phi);
+
+ // spire angle (if theta=2*PI, there is a complete spire between startx and x)
+ Real lengthCurve = x_used - x_start;
+ Real numSpire = lengthCurve / lengthSpire;
+ Real theta = 2 * M_PI * numSpire;
+
+ // computation of the Quat
+ Quat Qtheta;
+ Qtheta.axisToQuat(Vec3(0, 1, 0), theta);
+ Quat newSpireQuat = Qtheta * Qphi;
+
+ // computation of the position
+ Real radius = d_spireDiameter.getValue() / 2.0;
+ Vec3 PosEndCurve(radius * sin(theta), numSpire * d_spireHeight.getValue(), radius * (cos(theta) - 1));
+ Vec3 SpirePos = PosEndCurve + Vec3(x_start, 0, 0);
+
+ global_H_local.set(SpirePos, newSpireQuat);
+}
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodStraightSection.cpp b/src/BeamAdapter/component/model/RodStraightSection.cpp
new file mode 100644
index 000000000..777466bf0
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodStraightSection.cpp
@@ -0,0 +1,39 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#define SOFA_PLUGIN_BEAMADAPTER_RODSTRAIGHTSECTION_CPP
+
+#include
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+using namespace sofa::defaulttype;
+
+const int RodStraightSectionClass = core::RegisterObject("Class defining a rod straight section Material, defining material and geometry parameters.")
+ .add< RodStraightSection >(true);
+
+template class SOFA_BEAMADAPTER_API RodStraightSection;
+
+}// namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodStraightSection.h b/src/BeamAdapter/component/model/RodStraightSection.h
new file mode 100644
index 000000000..629c946fb
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodStraightSection.h
@@ -0,0 +1,58 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+/**
+ * \class RodStraightSection
+ * \brief Specialization class of @sa BaseRodSectionMaterial describing a rod straight section.
+ *
+ * This class will describe a rod straight section which will parametrized only using the @sa BaseRodSectionMaterial Data
+ * Method @sa getRestTransformOnX will return: Vec3(current_x, 0 0)
+ */
+template
+class RodStraightSection : public sofa::beamadapter::BaseRodSectionMaterial
+{
+public:
+ SOFA_CLASS(SOFA_TEMPLATE(RodStraightSection, DataTypes), SOFA_TEMPLATE(BaseRodSectionMaterial, DataTypes));
+
+ /// Default Constructor
+ RodStraightSection();
+
+ /// Override method to get the rest position of the beam. In this implementation, it will basically returns Vec3(x_start + x_used, 0 0)
+ void getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start) override;
+
+protected:
+ /// Internal method to init the section. Called by @sa BaseRodSectionMaterial::init() method
+ bool initSection() override;
+};
+
+#if !defined(SOFA_PLUGIN_BEAMADAPTER_RODSTRAIGHTSECTION_CPP)
+extern template class SOFA_BEAMADAPTER_API RodStraightSection;
+#endif
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/component/model/RodStraightSection.inl b/src/BeamAdapter/component/model/RodStraightSection.inl
new file mode 100644
index 000000000..cc499ca9b
--- /dev/null
+++ b/src/BeamAdapter/component/model/RodStraightSection.inl
@@ -0,0 +1,73 @@
+/******************************************************************************
+* BeamAdapter plugin *
+* (c) 2006 Inria, University of Lille, CNRS *
+* *
+* This program is free software; you can redistribute it and/or modify it *
+* under the terms of the GNU Lesser General Public License as published by *
+* the Free Software Foundation; either version 2.1 of the License, or (at *
+* your option) any later version. *
+* *
+* This program is distributed in the hope that it will be useful, but WITHOUT *
+* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or *
+* FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License *
+* for more details. *
+* *
+* You should have received a copy of the GNU Lesser General Public License *
+* along with this program. If not, see . *
+*******************************************************************************
+* Authors: see Authors.md *
+* *
+* Contact information: contact@sofa-framework.org *
+******************************************************************************/
+#pragma once
+
+#include
+#include
+#include
+
+namespace sofa::beamadapter
+{
+
+template
+RodStraightSection::RodStraightSection()
+ : BaseRodSectionMaterial()
+{
+
+}
+
+
+template
+bool RodStraightSection::initSection()
+{
+ const auto length = this->d_length.getValue();
+ if (length <= Real(0.0))
+ {
+ msg_error() << "Length is 0 (or negative), check if d_length has been given or well computed.";
+ return false;
+ }
+
+ if (int nbrEdgesVisu = d_nbEdgesVisu.getValue() <= 0)
+ {
+ msg_warning() << "Number of visual edges has been set to an invalid value: " << nbrEdgesVisu << ". Value should be a positive integer. Setting to default value: 10";
+ d_nbEdgesVisu.setValue(10);
+ }
+
+
+ if (int nbEdgesCollis = d_nbEdgesCollis.getValue() <= 0)
+ {
+ msg_warning() << "Number of collision edges has been set to an invalid value: " << nbEdgesCollis << ". Value should be a positive integer. Setting to default value: 20";
+ d_nbEdgesCollis.setValue(10);
+ }
+
+ return true;
+}
+
+
+template
+void RodStraightSection::getRestTransformOnX(Transform& global_H_local, const Real& x_used, const Real& x_start)
+{
+ global_H_local.set(Vec3(x_start + x_used, 0.0, 0.0), Quat());
+}
+
+
+} // namespace sofa::beamadapter
diff --git a/src/BeamAdapter/gpu/cuda/CudaInstantiations.cpp b/src/BeamAdapter/gpu/cuda/CudaInstantiations.cpp
index 5441f7f7d..94079918e 100644
--- a/src/BeamAdapter/gpu/cuda/CudaInstantiations.cpp
+++ b/src/BeamAdapter/gpu/cuda/CudaInstantiations.cpp
@@ -28,6 +28,9 @@
#include
#include
#include
+#include
+#include
+#include
#include
#include
@@ -86,6 +89,19 @@ namespace sofa::component::mapping
#endif
} // namespace sofa::component::mapping
+namespace sofa::beamadapter
+{
+ template class SOFA_BEAMADAPTER_API RodMeshSection;
+ template class SOFA_BEAMADAPTER_API RodSpireSection;
+ template class SOFA_BEAMADAPTER_API RodStraightSection;
+
+#ifdef SOFA_GPU_CUDA_DOUBLE
+ template class SOFA_BEAMADAPTER_API RodMeshSection;
+ template class SOFA_BEAMADAPTER_API RodSpireSection;
+ template class SOFA_BEAMADAPTER_API RodStraightSection;
+#endif
+} // namespace sofa::beamadapter
+
namespace sofa::gpu::cuda
{
@@ -139,6 +155,29 @@ int CudaMultiAdaptiveBeamMappingClass = core::RegisterObject("Set the positions
#endif
;
+const int CudaRodMeshSectionClass = core::RegisterObject("Class defining a Rod Section using a MeshLoader and material parameters using CUDA.")
+ .add< sofa::beamadapter::RodMeshSection >()
+#ifdef SOFA_GPU_CUDA_DOUBLE
+ .add< sofa::beamadapter::RodMeshSection