diff --git a/CMakeLists.txt b/CMakeLists.txt index c6569a363..f1455dab3 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,12 @@ 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/BaseRodSectionMaterial.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/src/BeamAdapter/component/engine/WireRestShape.h b/src/BeamAdapter/component/engine/WireRestShape.h index 8c77de028..34598fcc1 100644 --- a/src/BeamAdapter/component/engine/WireRestShape.h +++ b/src/BeamAdapter/component/engine/WireRestShape.h @@ -97,10 +97,10 @@ class WireRestShape : public core::objectmodel::BaseObject 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) diff --git a/src/BeamAdapter/component/engine/WireRestShape.inl b/src/BeamAdapter/component/engine/WireRestShape.inl index dda4ad4b5..c2c89726c 100644 --- a/src/BeamAdapter/component/engine/WireRestShape.inl +++ b/src/BeamAdapter/component/engine/WireRestShape.inl @@ -431,7 +431,7 @@ void WireRestShape::getRestTransformOnX(Transform &global_H_local, co template -void WireRestShape::getYoungModulusAtX(const Real& x_curv, Real& youngModulus, Real& cPoisson) +void WireRestShape::getYoungModulusAtX(const Real& x_curv, Real& youngModulus, Real& cPoisson) const { //Initialization Real _E1, _E2; @@ -462,7 +462,7 @@ void WireRestShape::getYoungModulusAtX(const Real& x_curv, Real& youn template -void WireRestShape::getInterpolationParam(const Real& x_curv, Real &_rho, Real &_A, Real &_Iy , Real &_Iz, Real &_Asy, Real &_Asz, Real &_J) +void WireRestShape::getInterpolationParam(const Real& x_curv, Real &_rho, Real &_A, Real &_Iy , Real &_Iz, Real &_Asy, Real &_Asz, Real &_J) const { if(x_curv <= this->d_straightLength.getValue()) { diff --git a/src/BeamAdapter/component/model/BaseRodSectionMaterial.cpp b/src/BeamAdapter/component/model/BaseRodSectionMaterial.cpp new file mode 100644 index 000000000..eec87a3bc --- /dev/null +++ b/src/BeamAdapter/component/model/BaseRodSectionMaterial.cpp @@ -0,0 +1,34 @@ +/****************************************************************************** +* 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_BASERODSECTIONMATERIAL_CPP + +#include +#include + +namespace sofa::beamadapter +{ + +using namespace sofa::defaulttype; + +template class SOFA_BEAMADAPTER_API BaseRodSectionMaterial; + +}// namespace sofa::beamadapter diff --git a/src/BeamAdapter/component/model/BaseRodSectionMaterial.h b/src/BeamAdapter/component/model/BaseRodSectionMaterial.h new file mode 100644 index 000000000..357a4a501 --- /dev/null +++ b/src/BeamAdapter/component/model/BaseRodSectionMaterial.h @@ -0,0 +1,122 @@ +/****************************************************************************** +* 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; +}; + +#if !defined(SOFA_PLUGIN_BEAMADAPTER_BASERODSECTIONMATERIAL_CPP) +extern template class SOFA_BEAMADAPTER_API BaseRodSectionMaterial; +#endif + +} // 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..1e23ec6f4 --- /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..85ea072d5 --- /dev/null +++ b/src/BeamAdapter/component/model/RodMeshSection.h @@ -0,0 +1,91 @@ +/****************************************************************************** +* 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)); + + using Real = typename DataTypes::Real; + using Transform = typename sofa::defaulttype::SolidTypes::Transform; + using Coord = typename DataTypes::Coord; + using Quat = sofa::type::Quat; + + /// 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 type::Quat& input, type::Vec3& x, type::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..b4de00487 --- /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 /(this->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 type::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; + } + + type::Vec3 PosEndCurve = p.getCenter(); + Quat ExtremityQuat = p.getOrientation(); + type::Vec3 ExtremityPos = PosEndCurve + type::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(type::Vec3(0, 0, 0)); + m_localRestTransforms[0].setOrientation(input); + + for (unsigned int i = 0; i < m_localRestPositions.size() - 1; i++) + { + type::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); + + type::Vec3 localPos = m_localRestPositions[i + 1] - m_localRestPositions[0]; + + m_localRestTransforms[i + 1].setOrigin(localPos); + + m_curvAbs.push_back(tot); + } + m_absOfGeometry = tot; + + this->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, type::Vec3& x, Quat& output) +{ + x.normalize(); + type::Vec3 x0 = input.inverseRotate(x); + + Real cTheta = x0[0]; + Real theta; + if (cTheta > (1 - EPSILON)) + { + output = input; + } + else + { + theta = acos(cTheta); + // axis of rotation + type::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..e0dd278db --- /dev/null +++ b/src/BeamAdapter/component/model/RodSpireSection.cpp @@ -0,0 +1,40 @@ +/****************************************************************************** +* 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 +#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..20b057afb --- /dev/null +++ b/src/BeamAdapter/component/model/RodSpireSection.h @@ -0,0 +1,69 @@ +/****************************************************************************** +* 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)); + + using Real = typename DataTypes::Real; + using Transform = typename sofa::defaulttype::SolidTypes::Transform; + using Quat = sofa::type::Quat; + + /// 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..79e3be5a2 --- /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 = this->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"; + this->d_nbEdgesVisu.setValue(10); + } + + + if (int nbEdgesCollis = this->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"; + this->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(type::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(type::Vec3(0, 1, 0), theta); + Quat newSpireQuat = Qtheta * Qphi; + + // computation of the position + Real radius = d_spireDiameter.getValue() / 2.0; + type::Vec3 PosEndCurve(radius * sin(theta), numSpire * d_spireHeight.getValue(), radius * (cos(theta) - 1)); + type::Vec3 SpirePos = PosEndCurve + type::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..4a433a145 --- /dev/null +++ b/src/BeamAdapter/component/model/RodStraightSection.cpp @@ -0,0 +1,41 @@ +/****************************************************************************** +* 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 + +#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..0e2935d15 --- /dev/null +++ b/src/BeamAdapter/component/model/RodStraightSection.h @@ -0,0 +1,62 @@ +/****************************************************************************** +* 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)); + + using Real = typename DataTypes::Real; + using Transform = typename sofa::defaulttype::SolidTypes::Transform; + using Quat = sofa::type::Quat; + + /// 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..2cc4b445b --- /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 = this->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"; + this->d_nbEdgesVisu.setValue(10); + } + + + if (int nbEdgesCollis = this->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"; + this->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(type::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 >() +#endif + ; + +const int CudaRodSpireSectionClass = core::RegisterObject("Class defining a rod spire section, defining material and geometry parameters using CUDA.") + .add< sofa::beamadapter::RodSpireSection >() +#ifdef SOFA_GPU_CUDA_DOUBLE + .add< sofa::beamadapter::RodSpireSection >() +#endif +; + +const int CudaRodStraightSectionClass = core::RegisterObject("Class defining a rod straight section Material, defining material and geometry parameters using CUDA.") + .add< sofa::beamadapter::RodStraightSection >() +#ifdef SOFA_GPU_CUDA_DOUBLE + .add< sofa::beamadapter::RodStraightSection >() +#endif +; + + + } // namespace sofa::gpu::cuda