diff --git a/CMakeLists.txt b/CMakeLists.txt index 2c0a3bf48..fcfddddee 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -279,6 +279,16 @@ GMDS_ADD_COMPONENT( OFF # must be covered ) +GMDS_ADD_COMPONENT( + MORPHMESH + morphMesh + GMDSmorphMesh + "morphing the mesh" + OFF + OFF # must be covered +) + + #============================================================================== set (GMDS_INCLUDE_DIRS APPEND) diff --git a/math/CMakeLists.txt b/math/CMakeLists.txt index 7614a6751..de84c8500 100644 --- a/math/CMakeLists.txt +++ b/math/CMakeLists.txt @@ -37,6 +37,7 @@ set(GMDS_INC inc/gmds/math/Vector.h inc/gmds/math/VectorDyn.h inc/gmds/math/StreamComputation.h + inc/gmds/math/BezierSurface.h ) set(GMDS_SRC src/AxisAngleRotation.cpp @@ -67,10 +68,11 @@ set(GMDS_SRC src/VectorDyn.cpp src/Orientation.cpp src/Vector.cpp - src/TransfiniteInterpolation.cpp inc/gmds/math/BezierSurface.h src/BezierSurface.cpp) + src/TransfiniteInterpolation.cpp + src/BezierSurface.cpp) #============================================================================== find_package(Eigen3 3.3 REQUIRED) -#============================================================================== +#==============================================================================k add_library(${GMDS_LIB} ${GMDS_INC} ${GMDS_SRC}) #============================================================================== include(GenerateExportHeader) diff --git a/morphMesh/CMakeLists.txt b/morphMesh/CMakeLists.txt new file mode 100644 index 000000000..cde2637d9 --- /dev/null +++ b/morphMesh/CMakeLists.txt @@ -0,0 +1,60 @@ +#============================================================================== +# LIBRARY DEFINTION (SOURCE FILES) +#============================================================================== +# Nommer tout en GMDS_MODULE_NAME, GMDS_SRC, ... dans les composants +set(GMDS_LIB ${LIB_GMDS_MORPHMESH}) +set(GMDS_LIB_PREFIX gmds/morphMesh) + +set(GMDS_INC + ${CMAKE_BINARY_DIR}/exports/${GMDS_LIB}_export.h + inc/gmds/morphMesh/FastLocalize.h + inc/gmds/morphMesh/MorphMesh.h + inc/gmds/morphMesh/EllipticMorph.h + ) +set(GMDS_SRC + src/FastLocalize.cpp + src/MorphMesh.cpp + src/EllipticMorph.cpp + ) +#============================================================================== +add_library(${GMDS_LIB} ${GMDS_INC} ${GMDS_SRC}) +#============================================================================== +include(GenerateExportHeader) +generate_export_header(${GMDS_LIB} + EXPORT_FILE_NAME ${CMAKE_BINARY_DIR}/exports/${GMDS_LIB}_export.h + EXPORT_MACRO_NAME ${GMDS_LIB}_API) +#============================================================================== +# TARGET DEFINITION +#============================================================================== +include(GNUInstallDirs) + +#LIBRARY TO INSTALL +target_link_libraries(${GMDS_LIB} PUBLIC + ${LIB_GMDS_IG} + ${LIB_GMDS_IG_ALGO} + ${LIB_GMDS_IO} ANN) + +target_compile_features(${GMDS_LIB} PUBLIC cxx_std_14) + +# INCLUDE TO INSTALL +target_include_directories(${GMDS_LIB} PUBLIC + $ + $ + ) +set_target_properties(${GMDS_LIB} PROPERTIES PUBLIC_HEADER "${GMDS_INC}") + +install(TARGETS ${GMDS_LIB} + EXPORT GMDS_SUITE + PUBLIC_HEADER DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/${GMDS_LIB_PREFIX} + LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR} + ARCHIVE DESTINATION ${CMAKE_INSTALL_LIBDIR} + RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR}) +#============================================================================== +if(WITH_TEST) + add_subdirectory(tst) +endif(WITH_TEST) +#============================================================================== +add_executable(morphing src/main.cpp) +target_link_libraries(morphing PRIVATE ${GMDS_LIB}) +target_compile_features(morphing PUBLIC cxx_std_14) +install(TARGETS morphing) \ No newline at end of file diff --git a/morphMesh/README.md b/morphMesh/README.md new file mode 100644 index 000000000..c494bd762 --- /dev/null +++ b/morphMesh/README.md @@ -0,0 +1,3 @@ +#morphMesh module + +This module propose algorithm to morph mesh using deformation function. \ No newline at end of file diff --git a/morphMesh/inc/gmds/morphMesh/EllipticMorph.h b/morphMesh/inc/gmds/morphMesh/EllipticMorph.h new file mode 100644 index 000000000..f56251d7c --- /dev/null +++ b/morphMesh/inc/gmds/morphMesh/EllipticMorph.h @@ -0,0 +1,55 @@ + +#ifndef GMDS_ELLIPTICMORPH_H +#define GMDS_ELLIPTICMORPH_H + +/*----------------------------------------------------------------------------*/ +#include "LIB_GMDS_MORPHMESH_export.h" +#include "gmds/ig/Mesh.h" +/*----------------------------------------------------------------------------*/ +namespace gmds { +/*----------------------------------------------------------------------------*/ +namespace morphmesh { +/*----------------------------------------------------------------------------*/ +class LIB_GMDS_MORPHMESH_API EllipticMorph +{ + + public: + /*------------------------------------------------------------------------*/ + /** @brief Constructor. + */ + EllipticMorph(Mesh* AMesh); + + /*------------------------------------------------------------------------*/ + /** @brief Destructor. + */ + ~EllipticMorph(); + + /*------------------------------------------------------------------------*/ + /** @brief Constructor. + */ + void execute(const std::vector>& AEllipses); + + private: + + /*Surement a retirer plus tard*/ + void markLockedCells(); + + + private: + + /* Locked nodes */ + int m_locked; + /* Nodes on the surface of the mesh */ + std::vector m_surfNodes; + std::vector m_locked_faces; + + + /* The mesh to deform */ + Mesh* m_mesh; +}; +/*----------------------------------------------------------------------------*/ +}// namespace morphmesh +/*----------------------------------------------------------------------------*/ +}// namespace gmds +/*----------------------------------------------------------------------------*/ +#endif // GMDS_ELLIPTICMORPH_H diff --git a/morphMesh/inc/gmds/morphMesh/FastLocalize.h b/morphMesh/inc/gmds/morphMesh/FastLocalize.h new file mode 100644 index 000000000..849218c42 --- /dev/null +++ b/morphMesh/inc/gmds/morphMesh/FastLocalize.h @@ -0,0 +1,52 @@ +/*----------------------------------------------------------------------------*/ +#ifndef GMDS_MORPHMESH_FASTLOCALIZE_H +#define GMDS_MORPHMESH_FASTLOCALIZE_H +/*----------------------------------------------------------------------------*/ +#include "LIB_GMDS_MORPHMESH_export.h" +# include +/*----------------------------------------------------------------------------*/ +#include +#include +/*----------------------------------------------------------------------------*/ +namespace gmds { +/*----------------------------------------------------------------------------*/ +class LIB_GMDS_MORPHMESH_API FastLocalize { + public: + + /*-------------------------------------------------------------------*/ + /** @brief Constructor. + * @param AMesh the mesh where we work on + */ + FastLocalize(const std::vector &ANodes); + /*-------------------------------------------------------------------*/ + /** @brief Destructor. + */ + virtual ~FastLocalize(); + + /*-------------------------------------------------------------------*/ + /** @brief Given an input point \p APoint return the cell data that + * indicate what is the closest cell of \APoint (it can be a node or + * a face) + * @param APoint point to localize + * @return the cell data indicating where \p APoint is in this->m_mesh + */ + TCellID find(const math::Point& APoint); + /*-------------------------------------------------------------------*/ + + private: + void buildANNTree(); + private: + /** mesh we work on */ + std::vector m_nodes; + std::map m_ann2gmds_id; + ANNkd_tree* m_kdTree; + ANNidxArray m_nnIdx; // near neighbor indices + ANNdistArray m_dists; // near neighbor distances + ANNpoint m_queryPt; // query point + ANNpointArray m_dataPts; // data points +}; +/*----------------------------------------------------------------------------*/ +} // namespace gmds +/*----------------------------------------------------------------------------*/ +#endif // GMDS_MORPHMESH_FASTLOCALIZE_H +/*----------------------------------------------------------------------------*/ diff --git a/morphMesh/inc/gmds/morphMesh/MorphMesh.h b/morphMesh/inc/gmds/morphMesh/MorphMesh.h new file mode 100644 index 000000000..09d76532f --- /dev/null +++ b/morphMesh/inc/gmds/morphMesh/MorphMesh.h @@ -0,0 +1,66 @@ +/*----------------------------------------------------------------------------*/ +#ifndef GMDS_MORPHMESH_H +#define GMDS_MORPHMESH_H + +/*----------------------------------------------------------------------------*/ +#include "LIB_GMDS_MORPHMESH_export.h" +#include "gmds/ig/Mesh.h" +/*----------------------------------------------------------------------------*/ +namespace gmds { +/*----------------------------------------------------------------------------*/ +namespace morphmesh { +/*----------------------------------------------------------------------------*/ +class LIB_GMDS_MORPHMESH_API MorphMesh +{ + + public: + /*------------------------------------------------------------------------*/ + /** @brief Constructor. + */ + MorphMesh(Mesh* AMesh, const std::vector &APoints, double ARadius); + + /*------------------------------------------------------------------------*/ + /** @brief Destructor. + */ + ~MorphMesh(); + + /*------------------------------------------------------------------------*/ + /** @brief Constructor. + */ + void execute(); + + private: + + /*Surement a retirer plus tard*/ + void markLockedCells(); + + bool computeLocalOrigin(const math::Point& AP, math::Point& AResult); + bool computeLocalOriginZ(const math::Point& AP, math::Point& AResult); + + double findHomothetyRatio(const math::Point& AP, const Node& ANode); + + private: + + /* Locked nodes */ + int m_locked; + Variable* locked_faces; + /* List of "homothetic points" that will define the mesh deformation */ + std::vector m_targets; + /* The radius of an area that will be modified by a homothetic point*/ + double m_radius; + /* Nodes on the surface of the mesh */ + std::vector m_surfNodes; + std::vector m_locked_faces; + //std::map is_locked_node; + //std::map nodes_origins; + + + /* The mesh to deform */ + Mesh* m_mesh; +}; +/*----------------------------------------------------------------------------*/ +}// namespace morphmesh +/*----------------------------------------------------------------------------*/ +}// namespace gmds +/*----------------------------------------------------------------------------*/ +#endif // GMDS_MORPHMESH_H diff --git a/morphMesh/src/EllipticMorph.cpp b/morphMesh/src/EllipticMorph.cpp new file mode 100644 index 000000000..f54472cd5 --- /dev/null +++ b/morphMesh/src/EllipticMorph.cpp @@ -0,0 +1,178 @@ +// +// Created by calderans on 30/08/23. +// + +#include "gmds/morphMesh/EllipticMorph.h" +#include "gmds/io/IGMeshIOService.h" +#include "gmds/morphMesh/FastLocalize.h" +/*----------------------------------------------------------------------------*/ + +using namespace gmds; +using namespace morphmesh; +/*----------------------------------------------------------------------------*/ +EllipticMorph::EllipticMorph(Mesh *AMesh): m_mesh(AMesh) { + + m_locked = m_mesh->newMark(); + markLockedCells(); + +} +/*----------------------------------------------------------------------------*/ +EllipticMorph::~EllipticMorph() {} +/*----------------------------------------------------------------------------*/ +void EllipticMorph::execute(const std::vector>& AEllipses) +{ + + std::vector nodes_int; + std::vector nodes_ext; + + + for(auto n : m_mesh->nodes()){ + Node node = m_mesh->get(n); + if(m_mesh->isMarked(n, m_locked)){ + nodes_int.push_back(node); + } + if((node.getIDs().size() == 4 || node.getIDs().size() == 2) && !m_mesh->isMarked(n,m_locked)){ + nodes_ext.push_back(node); + } + } + std::map vecs; + std::set modified_nodes; + + for(int i = 0; i < AEllipses.size()-1; i++) { + + double Xmin = AEllipses[i][0]; + double Xmax = AEllipses[i+1][0]; + + std::vector execution; + + for(auto n : m_mesh->nodes()){ + + Node node = m_mesh->get(n); + + if (node.X() >= Xmin && node.X() < Xmax && !m_mesh->isMarked(n, m_locked)) { + execution.push_back(node); + } + } + + double coefa_min = AEllipses[i][1]; + double coefb_min = AEllipses[i][2]; + double coefc_min = AEllipses[i][3]; + + double coefa_max = AEllipses[i+1][1]; + double coefb_max = AEllipses[i+1][2]; + double coefc_max = AEllipses[i+1][3]; + + FastLocalize fl_int(nodes_int); + FastLocalize fl_ext(nodes_ext); + + + for (auto n : execution) { + modified_nodes.insert(n.id()); + + math::Point p = n.point(); + + Node n_ext = m_mesh->get(fl_ext.find(p)); + math::Point p_int = m_mesh->get(fl_int.find(p)).point(); + math::Point p_ext = n_ext.point(); + math::Point p_axe(p.X(), 0, 0); + + double distXmin = p.X()-Xmin; + double distXmax = Xmax-Xmin; + double w = distXmin/distXmax; + + // math::Point p_origin = p.distance(p_int) < p.distance(p_axe) ? ; + + double coefy_min = p.Y() >= 0 ? coefa_min : coefc_min; + double coefy_max = p.Y() >= 0 ? coefa_max : coefc_max; + + double coefy_final = coefy_min * (1/(1 + exp(-5*( ( (1 - w ) *2)-1) ) )) + coefy_max * (1/(1 + exp(-5*( ( w *2)-1) ) )); + double coefz_final = coefb_min * (1/(1 + exp(-5*( ( (1 - w ) *2)-1) ) )) + coefb_max * (1/(1 + exp(-5*( ( w *2)-1) ) )); + + //double coefy_final = coefy_min * (1 - w) + coefy_max * (w); + //double coefz_final = coefb_min * (1 - w) + coefb_max * (w); + + + double theta = atan(n.Y() / n.Z()); + if (theta < 0) { + theta *= -1; + } + + double dist1 = p.distance(p_int); + double dist2 = p_ext.distance(p_int); + double omega = dist1 / dist2; + + // double z_new = omega*(coefa * cos(theta)); + // double y_new = omega*(coefb * sin(theta)); + double sinus = sin(theta); + double cosinus = cos(theta); + + if (n.Y() < 0) { + sinus *= -1; + } + if (n.Z() < 0) { + cosinus *= -1; + } + + math::Vector3d vec; + vec.setXYZ(0, (coefy_final * sinus) - p_ext.Y(), (coefz_final * cosinus) - p_ext.Z()); + + vec = omega * vec; + + vecs[n.id()] = vec; + } + } + + for(auto n : modified_nodes){ + + math::Vector3d vec(vecs[n]); + Node node = m_mesh->get(n); + double y = node.Y(); + double z = node.Z(); + + node.setY(y+vec.Y()); + node.setZ(z+vec.Z()); + } + +} +/*----------------------------------------------------------------------------*/ +void EllipticMorph::markLockedCells() +{ + Variable* locked_regions = m_mesh->newVariable("m_locked"); + + + std::set regions; + std::set nodes; + nodes.insert(27006); + m_mesh->mark(27006,m_locked); + for (int i = 0; i < 15; ++i) { + for(auto n : nodes){ + Node node = m_mesh->get(n); + + for(const auto& r : node.get()){ + regions.emplace(r); + } + } + nodes.clear(); + for(const auto& r : regions){ + locked_regions->set(r.id(),1); + for(auto n : r.getIDs()){ + if(!m_mesh->isMarked(n,m_locked)){ + nodes.emplace(n); + m_mesh->mark(n,m_locked); + + } + } + } + regions.clear(); + } + + + + for(auto n : m_mesh->nodes()){ + Node node = m_mesh->get(n); + math::Point p(node.point().X(),0,0); + if(node.point().distance(p) < 0.2){ + m_mesh->mark(n,m_locked); + } + } +} \ No newline at end of file diff --git a/morphMesh/src/FastLocalize.cpp b/morphMesh/src/FastLocalize.cpp new file mode 100644 index 000000000..2a90de9b9 --- /dev/null +++ b/morphMesh/src/FastLocalize.cpp @@ -0,0 +1,70 @@ +/*------------------------------------------------------------------------*/ +#include "gmds/morphMesh/FastLocalize.h" +//#include +/*------------------------------------------------------------------------*/ +using namespace gmds; +/*------------------------------------------------------------------------*/ +FastLocalize::FastLocalize(const std::vector &ANodes): m_nodes(ANodes){ + int nb_pnts = m_nodes.size(); + int k = 20; // max number of nearest neighbors + int dim = 3; // dimension + int maxPts = nb_pnts; // maximum number of data points + + int nPts; // actual number of data points + + m_queryPt = annAllocPt(3); + m_dataPts = annAllocPts(maxPts, dim); // allocate data points + m_nnIdx = new ANNidx[k]; // allocate near neigh indices + m_dists = new ANNdist[k]; // allocate near neighbor dists + + //======================================================== + // (1) Fill in the ANN structure for storing points + // + // Important: Points in APnts and dataPnts are stored with + // same index. + //======================================================== + nPts = 0; + for(auto const &n: m_nodes){ + math::Point p = n.point(); + m_dataPts[nPts][0] = p.X(); + m_dataPts[nPts][1] = p.Y(); + m_dataPts[nPts][2] = p.Z(); + m_ann2gmds_id[nPts]=n.id(); + nPts++; + } + + //======================================================== + // (2) Build the search structure + //======================================================== + m_kdTree = new ANNkd_tree(m_dataPts, // the data points + nPts, // number of points + dim); // dimension of space +} +/*------------------------------------------------------------------------*/ +FastLocalize::~FastLocalize(){ + delete [] m_nnIdx; // clean things up + delete [] m_dists; + delete m_kdTree; + annDeallocPt(m_queryPt); + annDeallocPts(m_dataPts); + annClose(); // done with ANN + +} +/*------------------------------------------------------------------------*/ +TCellID +FastLocalize::find(const math::Point &APoint) +{ + int k = 5; // max number of nearest neighbors + m_queryPt[0]=APoint.X(); + m_queryPt[1]=APoint.Y(); + m_queryPt[2]=APoint.Z(); + + m_kdTree->annkSearch( // search + m_queryPt,// query point + k, + m_nnIdx, + m_dists, + 0.01); + return m_ann2gmds_id[m_nnIdx[0]]; +} +/*------------------------------------------------------------------------*/ \ No newline at end of file diff --git a/morphMesh/src/MorphMesh.cpp b/morphMesh/src/MorphMesh.cpp new file mode 100644 index 000000000..02b1e5ae6 --- /dev/null +++ b/morphMesh/src/MorphMesh.cpp @@ -0,0 +1,329 @@ +/*----------------------------------------------------------------------------*/ +#include "gmds/morphMesh/MorphMesh.h" +#include "gmds/math/Plane.h" +#include "gmds/math/Triangle.h" +#include + +/*----------------------------------------------------------------------------*/ +namespace gmds { +/*----------------------------------------------------------------------------*/ +namespace morphmesh { +/*----------------------------------------------------------------------------*/ +MorphMesh::MorphMesh(Mesh *AMesh, const std::vector &APoints, double ARadius): m_mesh(AMesh) +{ + m_radius = ARadius; + m_targets = APoints; + m_locked = m_mesh->newMark(); + + locked_faces = m_mesh->newVariable("locked_face"); +} +/*----------------------------------------------------------------------------*/ +MorphMesh::~MorphMesh() {} +/*----------------------------------------------------------------------------*/ +void MorphMesh::execute() +{ + std::cout<<"Marking cells"<nodes()){ + Node node = m_mesh->get(n); + if(node.get().size() == 4) + m_surfNodes.push_back(n); + } + + std::vector durations1; + std::vector durations2; + + std::cout<<"Start deformation"<get(n).point(); + double dist = n_point.distance(p); + if(dist < min_dist){ + min_dist = dist; + n0 = n; + } + } + Node node0 = m_mesh->get(n0); + math::Point point0 = node0.point(); + // on calcul le ration d'homothétie pour ce point + //double H = findHomothetyRatio(p,node0); + double Y0 = point0.Y(); + double X0 = point0.X(); + double Z = (Y0 - p.Y())/2; + math::Point tmp; + computeLocalOrigin(point0, tmp); + double omega_height0 = Y0 - tmp.Y(); + math::Vector3d t; + t.setXYZ(0,p.Y()-Y0,0); + math::Vector3d tz; + tz.setXYZ(0,0,Z); + + + for(auto n : m_mesh->nodes()) { + if (!m_mesh->isMarked(n, m_locked)) { + Node node = m_mesh->get(n); + + double dist = node.point().distance(point0); + if (dist <= m_radius) { + time_t start,end; + time(&start); + // Si le point est couvert par la zone à modifier alors on continue le traitement + //double w = 1 - (dist / m_radius); + double sigmoid = 1/(1 + exp(-5*( ( (1 - (dist / m_radius) ) *2)-1) ) ); + + math::Point n_origin; + + bool origin_locked = computeLocalOrigin(node.point(), n_origin); + + double omega_height = Y0 - n_origin.Y(); + double n_height = node.point().Y() - n_origin.Y(); + + math::Vector3d final_t = (n_height / omega_height) * (sigmoid * t); + math::Point result = node.point() + final_t; + + node.setPoint(result); + time(&end); + durations1.emplace_back(double(end - start)); + } + + double distX = fabs(node.X() - X0); + if (distX <= m_radius) { + time_t start,end; + time(&start); + // double w = 0.5+2*cbrt((1 - (distX / m_radius)) - 0.5); + double sigmoid = 1 / (1 + exp(-5 * (((1 - (distX / m_radius)) * 2) - 1))); + + math::Point n_origin; + bool origin_locked = computeLocalOriginZ(node.point(), n_origin); + + double omega_height = fabs(Y0) - fabs(n_origin.Z()); + double n_height = fabs(node.point().Z()) - fabs(n_origin.Z()); + + int sens = node.Z() < 0 ? -1 : 1; + math::Vector3d final_t = (n_height / omega_height) * (sigmoid * (sens * tz)); + + math::Point result = node.point() + final_t; + node.setPoint(result); + time(&end); + durations1.emplace_back(double(end - start)); + } + } + } + /*for(auto n : m_mesh->nodes()) { + if (!m_mesh->isMarked(n, m_locked)) { + Node node = m_mesh->get(n); + double distX = fabs(node.X() - X0); + if (distX <= m_radius) { + // double w = 0.5+2*cbrt((1 - (distX / m_radius)) - 0.5); + double sigmoid = 1 / (1 + exp(-5 * (((1 - (distX / m_radius)) * 2) - 1))); + + math::Point n_origin; + bool origin_locked = computeLocalOriginZ(node.point(), n_origin); + + double omega_height = fabs(Y0) - fabs(n_origin.Z()); + double n_height = fabs(node.point().Z()) - fabs(n_origin.Z()); + + int sens = node.Z() < 0 ? -1 : 1; + math::Vector3d final_t = (n_height / omega_height) * (sigmoid * (sens * tz)); + + math::Point result = node.point() + final_t; + node.setPoint(result); + } + } + }*/ + } + + time(&endglob); + + + double avg = durations1[0]; + for(int i = 1; iunmarkAll(m_locked); + m_mesh->freeMark(m_locked); +} +/*----------------------------------------------------------------------------*/ +void MorphMesh::markLockedCells() +{ + Variable* locked_regions = m_mesh->newVariable("m_locked"); + std::set regions; + std::set nodes; + nodes.insert(27006); + m_mesh->mark(27006,m_locked); + for (int i = 0; i < 15; ++i) { + for(auto n : nodes){ + Node node = m_mesh->get(n); + for(const auto& r : node.get()){ + regions.emplace(r); + } + } + nodes.clear(); + for(const auto& r : regions){ + locked_regions->set(r.id(),1); + for(auto n : r.getIDs()){ + if(!m_mesh->isMarked(n,m_locked)){ + nodes.emplace(n); + m_mesh->mark(n,m_locked); + } + } + } + regions.clear(); + } + + for(auto n : m_mesh->nodes()){ + Node node = m_mesh->get(n); + if(node.Y() <= 0 ){ + m_mesh->mark(n,m_locked); + } + } + + //Test origine locale décallée + /*for(auto r : m_mesh->regions()){ + //if((4440<=r && r<=4639) || (5140<=r && r<=5339) || (5840<=r && r<=6039) || (6540<=r && r<=6739)){ + if((37000<=r && r<=37999) || (39000<=r && r<=39999) || (41000<=r && r<=41999) || (43000<=r && r<=43999)){ + locked_regions->set(r,1); + Region region = m_mesh->get(r); + for(auto n : region.getIDs()){ + m_mesh->mark(n,m_locked); + } + } + }*/ + + + for(auto f : m_mesh->faces()){ + Face face = m_mesh->get(f); + + std::vector f_regions = face.getIDs(); + if(f_regions.size() == 2){ + //Une face entre un hexa lock et un non lock + if((locked_regions->value(f_regions[0]) == 1) && (locked_regions->value(f_regions[1]) == 0) || + (locked_regions->value(f_regions[0]) == 0) && (locked_regions->value(f_regions[1]) == 1)){ + locked_faces->set(f,1); + m_locked_faces.push_back(f); + } + } + } +} +/*----------------------------------------------------------------------------*/ +bool MorphMesh::computeLocalOrigin(const math::Point& AP, math::Point& AResult) +{ + /* + * 1) on crée un segment entre de AP et (X,Y,0) + * 2) On intersecte le segment avec les faces lock + * - si on a un point alors on retourne ce point + * - sinon on retourne (X,Y,0) + */ + std::vector points; + points.emplace_back(AP.X(),0,AP.Z()); + points.emplace_back(AP.X()+0.001,0,AP.Z()+0.001); + points.emplace_back(AP.X()-0.001,0,AP.Z()-0.001); + + for (auto f : m_locked_faces) { + for(const auto& p : points) { + //if (locked_faces->value(f) == 1) { + math::Segment depth(AP, p); + Face face = m_mesh->get(f); + math::Plane plane(face.center(), face.normal()); + // on cherche si la profondeur a une intersection avec un plan d'une face locked + if (depth.intersect(plane, AResult)) { + // si oui on cherche si le point d'intersection est dans la face + // On fait ça en créant deux triangles qui couvrent entièrement la face + // et qui permettent d'utiliser la méthode isIn() + std::vector nodes = face.get(); + math::Triangle tri1(nodes[0].point(), nodes[1].point(), nodes[2].point()); + math::Triangle tri2(nodes[2].point(), nodes[3].point(), nodes[0].point()); + if (tri1.intersect(depth) || tri2.intersect(depth)) { + return true; + } + } + //} + } + } + AResult = {AP.X(),0,AP.Z()}; + + return false; +} +/*----------------------------------------------------------------------------*/ +bool MorphMesh::computeLocalOriginZ(const math::Point& AP, math::Point& AResult) +{ + /* + * 1) on crée un segment entre de AP et (X,Y,0) + * 2) On intersecte le segment avec les faces lock + * - si on a un point alors on retourne ce point + * - sinon on retourne (X,Y,0) + */ + std::vector points; + points.emplace_back(AP.X(),AP.Y(),0); + points.emplace_back(AP.X()+0.001,AP.Y()+0.001,0); + points.emplace_back(AP.X()-0.001,AP.Y()-0.001,0); + + for (auto f : m_locked_faces) { + for(const auto& p : points) { + //if (locked_faces->value(f) == 1) { + math::Segment depth(AP,math::Point(AP.X(),AP.Y(),0)); + Face face = m_mesh->get(f); + math::Plane plane(face.center(), face.normal()); + // on cherche si la profondeur a une intersection avec un plan d'une face locked + if (depth.intersect(plane, AResult)) { + // std::cout<<"test"< nodes = face.get(); + math::Triangle tri1(nodes[0].point(), nodes[1].point(), nodes[2].point()); + math::Triangle tri2(nodes[2].point(), nodes[3].point(), nodes[0].point()); + if (tri1.intersect(depth) || tri2.intersect(depth)) { + return true; + } + } + //} + } + } + + AResult = {AP.X(),AP.Y(),0}; + + return false; +} +/*----------------------------------------------------------------------------*/ +double MorphMesh::findHomothetyRatio(const math::Point &AP, const Node& ANode) +{ + + math::Point p_node = ANode.point(); + // Le vecteur de déformation initial : La translation du noeud du maillage vers le point + math::Vector3d deform(AP - ANode.point()); + + + // On cherche l'origine du noeud initial + math::Point origin; + bool origin_locked = computeLocalOrigin(ANode.point(), origin); + // On fait un vecteur origin->ni pour calculer l'homothétie des deux vecteurs. + // Faire le vecteur dans ce sens là permet de facilement déduire le ration d'homothétie + math::Vector3d height(ANode.point() - origin); + + + double dist1 = origin.distance(ANode.point()); + double dist2 = origin.distance(AP); + + std::cout<<"dist origin -> node "< point "< +#include +/*------------------------------------------------------------------------*/ +using namespace gmds; +/*------------------------------------------------------------------------*/ +int main(int argc, char* argv[]){ + + std::string param_file(argv[1]); + + gmds::Mesh m_mesh(gmds::MeshModel(gmds::DIM3 | gmds::R | gmds::F | gmds::N | gmds::N2R | gmds::R2N | gmds::R2F | gmds::F2N | gmds::F2R)); + gmds::IGMeshIOService ioService(&m_mesh); + gmds::VTKReader vtkReader(&ioService); + vtkReader.setCellOptions(gmds::N|gmds::R); + vtkReader.setDataOptions(gmds::N|gmds::R); + vtkReader.read(param_file); + + gmds::MeshDoctor doc(&m_mesh); + doc.buildN2R(m_mesh.getModel()); + doc.buildFacesAndR2F(); + doc.updateUpwardConnectivity(); + + //gmds::morphmesh::MorphMesh morph(&m_mesh,points,2.5); + //morph.execute(); + + morphmesh::EllipticMorph emorph(&m_mesh); + + std::vector> ellipses; + //ellipses.push_back({2,1,1,1}); + ellipses.push_back({2,1,1,1}); + ellipses.push_back({3,1.2,1.5,1}); + ellipses.push_back({4,1.2,1.5,1}); + ellipses.push_back({5,1.4,1.75,1.5}); + ellipses.push_back({7,1.4,1.75,1.5}); + ellipses.push_back({8,1,1,1}); + //ellipses.push_back({5.6,1.2,1.5,1}); + //ellipses.push_back({8.5,1.2,1.5,1}); + + emorph.execute(ellipses); + + std::cout<<"test"< +#include +#include +#include +#include +#include +#include +#include +/*----------------------------------------------------------------------------*/ +/* TESTS UNITAIRES */ +/*----------------------------------------------------------------------------*/ + +TEST(MorphMesh_FastLocalizeTestSuite, test_FastLocalize_2D) +{ + // WE READ + gmds::Mesh m(gmds::MeshModel(gmds::DIM3 | gmds::N | gmds::F | gmds::F2N)); + + std::string dir(TEST_SAMPLES_DIR); + std::string vtk_file = dir + "/Aero/2D/C1_2D_0.1.vtk"; + + gmds::IGMeshIOService ioService(&m); + gmds::VTKReader vtkReader(&ioService); + vtkReader.setCellOptions(gmds::N | gmds::F); + vtkReader.read(vtk_file); + + // fastlocalize works on a point cloud, no need for the faces + for(auto i: m.faces()) { + m.deleteFace(i); + } + + gmds::FastLocalize fl(&m); + + gmds::Cell::Data data = fl.find(gmds::math::Point({-0.0338,0.74,0})); + ASSERT_EQ(data.dim, 0); + ASSERT_EQ(data.id, 734); + + data = fl.find(gmds::math::Point({-0.2,0.46,0})); + ASSERT_EQ(data.dim, 0); + ASSERT_EQ(data.id, 11); +} + + +TEST(MorphMesh_FastLocalizeTestSuite, test_FastLocalize_3D) +{ + // WE READ + gmds::Mesh m(gmds::MeshModel(gmds::DIM3 | gmds::N | gmds::F | gmds::F2N)); + + std::string dir(TEST_SAMPLES_DIR); + std::string vtk_file = dir + "/Aero/3D/C1_3D_0.1.vtk"; + + gmds::IGMeshIOService ioService(&m); + gmds::VTKReader vtkReader(&ioService); + vtkReader.setCellOptions(gmds::N | gmds::F); + vtkReader.read(vtk_file); + + // fastlocalize works on a point cloud, no need for the faces + for(auto i: m.faces()) { + m.deleteFace(i); + } + + gmds::FastLocalize fl(&m); + + gmds::Cell::Data data = fl.find(gmds::math::Point({-1.0,0.0,0})); + ASSERT_EQ(data.dim, 0); + ASSERT_EQ(data.id, 14322); + + data = fl.find(gmds::math::Point({-1.33,0.678,-1.89})); + ASSERT_EQ(data.dim, 0); + ASSERT_EQ(data.id, 4048); +} + + +TEST(MorphMesh_FastLocalizeTestSuite, test_FastLocalize_nodes) +{ + // WE READ + gmds::Mesh m(gmds::MeshModel(gmds::DIM3 | gmds::N | gmds::F | gmds::F2N)); + + for(int i=0; i<50; i++) { + m.newNode(i,1,1); + } + + gmds::FastLocalize fl(&m); + + gmds::Cell::Data data = fl.find(gmds::math::Point({1.9,0.0,0})); + ASSERT_EQ(data.dim, 0); + ASSERT_EQ(data.id, 2); +} \ No newline at end of file diff --git a/morphMesh/tst/main_test.cpp b/morphMesh/tst/main_test.cpp new file mode 100644 index 000000000..174d6c0cd --- /dev/null +++ b/morphMesh/tst/main_test.cpp @@ -0,0 +1,13 @@ +/*----------------------------------------------------------------------------*/ +#include +/*----------------------------------------------------------------------------*/ +// Files containing the different test suites to launch + +#include "FastLocalizeTestSuite.h" +/*----------------------------------------------------------------------------*/ +int main(int argc, char ** argv) { + ::testing::InitGoogleTest(&argc, argv); + return RUN_ALL_TESTS(); +} +/*----------------------------------------------------------------------------*/ +