]> SALOME platform Git repositories - modules/smesh.git/blob - src/SMESHUtils/SMESH_MeshAlgos.hxx
Salome HOME
Merge eap/23491 branch.
[modules/smesh.git] / src / SMESHUtils / SMESH_MeshAlgos.hxx
1 // Copyright (C) 2007-2016  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
5 //
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
10 //
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14 // Lesser General Public License for more details.
15 //
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
19 //
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 //
22 // File      : SMESH_MeshAlgos.hxx
23 // Created   : Tue Apr 30 18:00:36 2013
24 // Author    : Edward AGAPOV (eap)
25
26 // This file holds some low level algorithms extracted from SMESH_MeshEditor
27 // to make them accessible from Controls package
28
29
30 #ifndef __SMESH_MeshAlgos_HXX__
31 #define __SMESH_MeshAlgos_HXX__
32
33 #include "SMESH_Utils.hxx"
34
35 #include "SMDSAbs_ElementType.hxx"
36 #include "SMDS_ElemIterator.hxx"
37 #include "SMESH_TypeDefs.hxx"
38
39 #include <TopAbs_State.hxx>
40 #include <vector>
41
42 class gp_Pnt;
43 class gp_Ax1;
44 class SMDS_MeshNode;
45 class SMDS_MeshElement;
46 class SMDS_Mesh;
47
48 //=======================================================================
49 /*!
50  * \brief Searcher for the node closest to a point
51  */
52 //=======================================================================
53
54 struct SMESHUtils_EXPORT SMESH_NodeSearcher
55 {
56   virtual const SMDS_MeshNode* FindClosestTo( const gp_Pnt& pnt ) = 0;
57   virtual void MoveNode( const SMDS_MeshNode* node, const gp_Pnt& toPnt ) = 0;
58   virtual int  FindNearPoint(const gp_Pnt&                        point,
59                              const double                         tolerance,
60                              std::vector< const SMDS_MeshNode* >& foundNodes) = 0;
61   virtual ~SMESH_NodeSearcher() {}
62 };
63
64 //=======================================================================
65 /*!
66  * \brief Searcher for elements
67  */
68 //=======================================================================
69
70 struct SMESHUtils_EXPORT SMESH_ElementSearcher
71 {
72   /*!
73    * \brief Find elements of given type where the given point is IN or ON.
74    *        Returns nb of found elements and elements them-selves.
75    *
76    * 'ALL' type means elements of any type excluding nodes and 0D elements
77    */
78   virtual int FindElementsByPoint(const gp_Pnt&                           point,
79                                   SMDSAbs_ElementType                     type,
80                                   std::vector< const SMDS_MeshElement* >& foundElems) = 0;
81   /*!
82    * \brief Return an element most close to the given point
83    */
84   virtual const SMDS_MeshElement* FindClosestTo( const gp_Pnt&       point,
85                                                  SMDSAbs_ElementType type) = 0;
86   /*!
87    * \brief Return elements possibly intersecting the line
88    */
89   virtual void GetElementsNearLine( const gp_Ax1&                           line,
90                                     SMDSAbs_ElementType                     type,
91                                     std::vector< const SMDS_MeshElement* >& foundElems) = 0;
92   /*!
93    * \brief Return elements whose bounding box intersects a sphere
94    */
95   virtual void GetElementsInSphere( const gp_XYZ&                           center,
96                                     const double                            radius,
97                                     SMDSAbs_ElementType                     type,
98                                     std::vector< const SMDS_MeshElement* >& foundElems) = 0;
99   /*!
100    * \brief Find out if the given point is out of closed 2D mesh.
101    */
102   virtual TopAbs_State GetPointState(const gp_Pnt& point) = 0;
103
104   /*!
105    * \brief Return a projection of a given point to a 2D mesh.
106    *        Optionally return the closest face
107    */
108   virtual gp_XYZ Project(const gp_Pnt&            point,
109                          SMDSAbs_ElementType      type,
110                          const SMDS_MeshElement** closestFace= 0) = 0;
111
112   virtual ~SMESH_ElementSearcher();
113 };
114
115 namespace SMESH_MeshAlgos
116 {
117   /*!
118    * \brief Return true if the point is IN or ON of the element
119    */
120   SMESHUtils_EXPORT
121   bool IsOut( const SMDS_MeshElement* element, const gp_Pnt& point, double tol );
122
123   SMESHUtils_EXPORT
124   double GetDistance( const SMDS_MeshElement* elem, const gp_Pnt& point, gp_XYZ* closestPnt = 0 );
125
126   SMESHUtils_EXPORT
127   double GetDistance( const SMDS_MeshEdge* edge, const gp_Pnt& point, gp_XYZ* closestPnt = 0 );
128
129   SMESHUtils_EXPORT
130   double GetDistance( const SMDS_MeshFace* face, const gp_Pnt& point, gp_XYZ* closestPnt = 0 );
131
132   SMESHUtils_EXPORT
133   double GetDistance( const SMDS_MeshVolume* volume, const gp_Pnt& point, gp_XYZ* closestPnt = 0 );
134
135   SMESHUtils_EXPORT
136   void GetBarycentricCoords( const gp_XY& point,
137                              const gp_XY& t0, const gp_XY& t1, const gp_XY& t2,
138                              double &    bc0, double &    bc1);
139
140   /*!
141    * Return a face having linked nodes n1 and n2 and which is
142    * - not in avoidSet,
143    * - in elemSet provided that !elemSet.empty()
144    * i1 and i2 optionally returns indices of n1 and n2
145    */
146   SMESHUtils_EXPORT
147   const SMDS_MeshElement* FindFaceInSet(const SMDS_MeshNode*    n1,
148                                         const SMDS_MeshNode*    n2,
149                                         const TIDSortedElemSet& elemSet,
150                                         const TIDSortedElemSet& avoidSet,
151                                         int*                    i1=0,
152                                         int*                    i2=0);
153   /*!
154    * \brief Calculate normal of a mesh face
155    */
156   SMESHUtils_EXPORT
157   bool FaceNormal(const SMDS_MeshElement* F, gp_XYZ& normal, bool normalized=true);
158
159   /*!
160    * \brief Return nodes common to two elements
161    */
162   SMESHUtils_EXPORT
163   std::vector< const SMDS_MeshNode*> GetCommonNodes(const SMDS_MeshElement* e1,
164                                                     const SMDS_MeshElement* e2);
165   /*!
166    * \brief Return true if node1 encounters first in the face and node2, after.
167    *        The nodes are supposed to be neighbor nodes in the face.
168    */
169   SMESHUtils_EXPORT
170   bool IsRightOrder( const SMDS_MeshElement* face,
171                      const SMDS_MeshNode*    node0,
172                      const SMDS_MeshNode*    node1 );
173
174   /*!
175    * \brief Mark elements given by SMDS_Iterator
176    */
177   template< class ElemIter >
178   void MarkElems( ElemIter it, const bool isMarked )
179   {
180     while ( it->more() ) it->next()->setIsMarked( isMarked );
181   }
182   /*!
183    * \brief Mark elements given by std iterators
184    */
185   template< class ElemIter >
186   void MarkElems( ElemIter it, ElemIter end, const bool isMarked )
187   {
188     for ( ; it != end; ++it ) (*it)->setIsMarked( isMarked );
189   }
190   /*!
191    * \brief Mark nodes of elements given by SMDS_Iterator
192    */
193   template< class ElemIter >
194   void MarkElemNodes( ElemIter it, const bool isMarked, const bool markElem = false )
195   {
196     if ( markElem )
197       while ( it->more() ) {
198         const SMDS_MeshElement* e = it->next();
199         e->setIsMarked( isMarked );
200         MarkElems( e->nodesIterator(), isMarked );
201       }
202     else
203       while ( it->more() )
204         MarkElems( it->next()->nodesIterator(), isMarked );
205   }
206   /*!
207    * \brief Mark elements given by std iterators
208    */
209   template< class ElemIter >
210   void MarkElemNodes( ElemIter it, ElemIter end, const bool isMarked, const bool markElem = false )
211   {
212     if ( markElem )
213       for ( ; it != end; ++it ) {
214         (*it)->setIsMarked( isMarked );
215         MarkElems( (*it)->nodesIterator(), isMarked );
216       }
217     else
218       for ( ; it != end; ++it )
219         MarkElems( (*it)->nodesIterator(), isMarked );
220   }
221
222   /*!
223    * \brief Return SMESH_NodeSearcher. The caller is responsible for deleting it
224    */
225   SMESHUtils_EXPORT
226   SMESH_NodeSearcher* GetNodeSearcher( SMDS_Mesh& mesh );
227
228   SMESHUtils_EXPORT
229   SMESH_NodeSearcher* GetNodeSearcher( SMDS_ElemIteratorPtr elemIt );
230
231   /*!
232    * \brief Return SMESH_ElementSearcher. The caller is responsible for deleting it
233    */
234   SMESHUtils_EXPORT
235   SMESH_ElementSearcher* GetElementSearcher( SMDS_Mesh& mesh,
236                                              double     tolerance=-1.);
237   SMESHUtils_EXPORT
238   SMESH_ElementSearcher* GetElementSearcher( SMDS_Mesh& mesh,
239                                              SMDS_ElemIteratorPtr elemIt,
240                                              double     tolerance=-1. );
241
242
243
244   typedef std::vector<const SMDS_MeshNode*> TFreeBorder;
245   typedef std::vector<TFreeBorder>          TFreeBorderVec;
246   struct TFreeBorderPart
247   {
248     int _border; // border index within a TFreeBorderVec
249     int _node1;  // node index within the border-th TFreeBorder
250     int _node2;
251     int _nodeLast;
252   };
253   typedef std::vector<TFreeBorderPart>  TCoincidentGroup;
254   typedef std::vector<TCoincidentGroup> TCoincidentGroupVec;
255   struct CoincidentFreeBorders
256   {
257     TFreeBorderVec      _borders;          // nodes of all free borders
258     TCoincidentGroupVec _coincidentGroups; // groups of coincident parts of borders
259   };
260
261   /*!
262    * Returns TFreeBorder's coincident within the given tolerance.
263    * If the tolerance <= 0.0 then one tenth of an average size of elements adjacent
264    * to free borders being compared is used.
265    *
266    * (Implemented in ./SMESH_FreeBorders.cxx)
267    */
268   SMESHUtils_EXPORT
269   void FindCoincidentFreeBorders(SMDS_Mesh&              mesh,
270                                  double                  tolerance,
271                                  CoincidentFreeBorders & foundFreeBordes);
272   /*!
273    * Returns all or only closed TFreeBorder's.
274    * Optionally check if the mesh is manifold and if faces are correctly oriented.
275    *
276    * (Implemented in ./SMESH_FreeBorders.cxx)
277    */
278   SMESHUtils_EXPORT
279   void FindFreeBorders(SMDS_Mesh&       mesh,
280                        TFreeBorderVec & foundFreeBordes,
281                        const bool       closedOnly,
282                        bool*            isManifold = 0,
283                        bool*            isGoodOri = 0);
284   /*!
285    * Fill a hole defined by a TFreeBorder with 2D elements.
286    *
287    * (Implemented in ./SMESH_FillHole.cxx)
288    */
289   SMESHUtils_EXPORT
290   void FillHole(const TFreeBorder &                   freeBorder,
291                 SMDS_Mesh&                            mesh,
292                 std::vector<const SMDS_MeshElement*>& newFaces);
293
294
295   /*!
296    * \brief Find nodes whose merge makes the element invalid
297    *
298    * (Implemented in SMESH_DeMerge.cxx)
299    */
300   SMESHUtils_EXPORT
301   void DeMerge(const SMDS_MeshElement*              elem,
302                std::vector< const SMDS_MeshNode* >& newNodes,
303                std::vector< const SMDS_MeshNode* >& noMergeNodes);
304
305 } // namespace SMESH_MeshAlgos
306
307 #endif