Salome HOME
PAL16774,PAL16631(SALOME crash after a mesh computation that failed because of lack...
[modules/smesh.git] / src / SMESH / SMESH_Block.hxx
1 //  Copyright (C) 2003  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
2 //  CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS 
3 // 
4 //  This library is free software; you can redistribute it and/or 
5 //  modify it under the terms of the GNU Lesser General Public 
6 //  License as published by the Free Software Foundation; either 
7 //  version 2.1 of the License. 
8 // 
9 //  This library is distributed in the hope that it will be useful, 
10 //  but WITHOUT ANY WARRANTY; without even the implied warranty of 
11 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU 
12 //  Lesser General Public License for more details. 
13 // 
14 //  You should have received a copy of the GNU Lesser General Public 
15 //  License along with this library; if not, write to the Free Software 
16 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA 
17 // 
18 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19
20 // File      : SMESH_Block.hxx
21 // Created   : Tue Nov 30 12:42:18 2004
22 // Author    : Edward AGAPOV (eap)
23
24
25 #ifndef SMESH_Block_HeaderFile
26 #define SMESH_Block_HeaderFile
27
28 #include <TopExp.hxx>
29 #include <TopTools_IndexedMapOfOrientedShape.hxx>
30 #include <TopoDS_Edge.hxx>
31 #include <TopoDS_Face.hxx>
32 #include <TopoDS_Shell.hxx>
33 #include <TopoDS_Vertex.hxx>
34 #include <gp_XY.hxx>
35 #include <gp_XYZ.hxx>
36 #include <math_FunctionSetWithDerivatives.hxx>
37
38 #include <ostream>
39 #include <vector>
40 #include <list>
41
42 class SMDS_MeshVolume;
43 class SMDS_MeshNode;
44 class Adaptor3d_Surface;
45 class Adaptor2d_Curve2d;
46 class Adaptor3d_Curve;
47 class gp_Pnt;
48
49 // =========================================================
50 // class calculating coordinates of 3D points by normalized
51 // parameters inside the block and vice versa
52 // =========================================================
53
54 class SMESH_Block: public math_FunctionSetWithDerivatives
55 {
56  public:
57   enum TShapeID {
58     // ----------------------------
59     // Ids of the block sub-shapes
60     // ----------------------------
61     ID_NONE = 0,
62
63     ID_V000 = 1, ID_V100, ID_V010, ID_V110, ID_V001, ID_V101, ID_V011, ID_V111,
64
65     ID_Ex00, ID_Ex10, ID_Ex01, ID_Ex11,
66     ID_E0y0, ID_E1y0, ID_E0y1, ID_E1y1,
67     ID_E00z, ID_E10z, ID_E01z, ID_E11z,
68
69     ID_Fxy0, ID_Fxy1, ID_Fx0z, ID_Fx1z, ID_F0yz, ID_F1yz,
70
71     ID_Shell
72     };
73   enum { // to use TShapeID for indexing certain type subshapes
74
75     ID_FirstV = ID_V000, ID_FirstE = ID_Ex00, ID_FirstF = ID_Fxy0
76
77   };
78
79
80  public:
81   // -------------------------------------------------
82   // Block topology in terms of block sub-shapes' ids
83   // -------------------------------------------------
84
85   static int NbVertices()  { return  8; }
86   static int NbEdges()     { return 12; }
87   static int NbFaces()     { return  6; }
88   static int NbSubShapes() { return ID_Shell; }
89   // to avoid magic numbers when allocating memory for subshapes
90
91   static inline bool IsVertexID( int theShapeID )
92   { return ( theShapeID >= ID_V000 && theShapeID <= ID_V111 ); }
93
94   static inline bool IsEdgeID( int theShapeID )
95   { return ( theShapeID >= ID_Ex00 && theShapeID <= ID_E11z ); }
96
97   static inline bool IsFaceID( int theShapeID )
98   { return ( theShapeID >= ID_Fxy0 && theShapeID <= ID_F1yz ); }
99
100   static int ShapeIndex( int theShapeID )
101   {
102     if ( IsVertexID( theShapeID )) return theShapeID - ID_V000;
103     if ( IsEdgeID( theShapeID ))   return theShapeID - ID_Ex00;
104     if ( IsFaceID( theShapeID ))   return theShapeID - ID_Fxy0;
105     return 0;
106   }
107   // return index [0-...] for each type of sub-shapes,
108   // for example :
109   // ShapeIndex( ID_Ex00 ) == 0
110   // ShapeIndex( ID_Ex10 ) == 1
111
112   static void GetFaceEdgesIDs (const int faceID, vector< int >& edgeVec );
113   // return edges IDs of a face in the order u0, u1, 0v, 1v
114
115   static void GetEdgeVertexIDs (const int edgeID, vector< int >& vertexVec );
116   // return vertex IDs of an edge
117
118   static int GetCoordIndOnEdge (const int theEdgeID)
119   { return (theEdgeID < ID_E0y0) ? 1 : (theEdgeID < ID_E00z) ? 2 : 3; }
120   // return an index of a coordinate which varies along the edge
121
122   static double* GetShapeCoef (const int theShapeID);
123   // for theShapeID( TShapeID ), returns 3 coefficients used
124   // to compute an addition of an on-theShape point to coordinates
125   // of an in-shell point. If an in-shell point has parameters (Px,Py,Pz),
126   // then the addition of a point P is computed as P*kx*ky*kz and ki is
127   // defined by the returned coef like this:
128   // ki = (coef[i] == 0) ? 1 : (coef[i] < 0) ? 1 - Pi : Pi
129
130   static int GetShapeIDByParams ( const gp_XYZ& theParams );
131   // define an id of the block sub-shape by point parameters
132
133   static ostream& DumpShapeID (const int theBlockShapeID, ostream& stream);
134   // DEBUG: dump an id of a block sub-shape
135
136
137  public:
138   // ---------------
139   // Initialization
140   // ---------------
141
142   SMESH_Block();
143
144   bool LoadBlockShapes(const TopoDS_Shell&         theShell,
145                        const TopoDS_Vertex&        theVertex000,
146                        const TopoDS_Vertex&        theVertex001,
147                        TopTools_IndexedMapOfOrientedShape& theShapeIDMap );
148   // Initialize block geometry with theShell,
149   // add sub-shapes of theBlock to theShapeIDMap so that they get
150   // IDs acoording to enum TShapeID
151
152   bool LoadBlockShapes(const TopTools_IndexedMapOfOrientedShape& theShapeIDMap);
153   // Initialize block geometry with shapes from theShapeIDMap
154
155   bool LoadMeshBlock(const SMDS_MeshVolume*        theVolume,
156                      const int                     theNode000Index,
157                      const int                     theNode001Index,
158                      vector<const SMDS_MeshNode*>& theOrderedNodes);
159   // prepare to work with theVolume and
160   // return nodes in theVolume corners in the order of TShapeID enum
161
162   bool LoadFace(const TopoDS_Face& theFace,
163                 const int          theFaceID,
164                 const TopTools_IndexedMapOfOrientedShape& theShapeIDMap);
165   // Load face geometry.
166   // It is enough to compute params or coordinates on the face.
167   // Face subshapes must be loaded into theShapeIDMap before
168
169   static bool Insert(const TopoDS_Shape& theShape,
170                      const int           theShapeID,
171                      TopTools_IndexedMapOfOrientedShape& theShapeIDMap);
172   // Insert theShape into theShapeIDMap with theShapeID,
173   // Not yet set shapes preceding theShapeID are filled with compounds
174   // Return true if theShape was successfully bound to theShapeID
175
176   static bool FindBlockShapes(const TopoDS_Shell&         theShell,
177                               const TopoDS_Vertex&        theVertex000,
178                               const TopoDS_Vertex&        theVertex001,
179                               TopTools_IndexedMapOfOrientedShape& theShapeIDMap );
180   // add sub-shapes of theBlock to theShapeIDMap so that they get
181   // IDs acoording to enum TShapeID
182
183 public:
184   // ---------------------------------
185   // Define coordinates by parameters
186   // ---------------------------------
187
188   bool VertexPoint( const int theVertexID, gp_XYZ& thePoint ) const {
189     if ( !IsVertexID( theVertexID ))           return false;
190     thePoint = myPnt[ theVertexID - ID_FirstV ]; return true;
191   }
192   // return vertex coordinates, parameters are defined by theVertexID
193
194   bool EdgePoint( const int theEdgeID, const gp_XYZ& theParams, gp_XYZ& thePoint ) const {
195     if ( !IsEdgeID( theEdgeID ))                                 return false;
196     thePoint = myEdge[ theEdgeID - ID_FirstE ].Point( theParams ); return true;
197   }
198   // return coordinates of a point on edge
199
200   bool EdgeU( const int theEdgeID, const gp_XYZ& theParams, double& theU ) const {
201     if ( !IsEdgeID( theEdgeID ))                              return false;
202     theU = myEdge[ theEdgeID - ID_FirstE ].GetU( theParams ); return true;
203   }
204   // return parameter on edge by in-block parameters
205
206   bool FacePoint( const int theFaceID, const gp_XYZ& theParams, gp_XYZ& thePoint ) const {
207     if ( !IsFaceID ( theFaceID ))                                return false;
208     thePoint = myFace[ theFaceID - ID_FirstF ].Point( theParams ); return true;
209   }
210   // return coordinates of a point on face
211
212   bool FaceUV( const int theFaceID, const gp_XYZ& theParams, gp_XY& theUV ) const {
213     if ( !IsFaceID ( theFaceID ))                               return false;
214     theUV = myFace[ theFaceID - ID_FirstF ].GetUV( theParams ); return true;
215   }
216   // return UV coordinates on a face by in-block parameters
217
218   bool ShellPoint( const gp_XYZ& theParams, gp_XYZ& thePoint ) const;
219   // return coordinates of a point in shell
220
221   static bool ShellPoint(const gp_XYZ&         theParams,
222                          const vector<gp_XYZ>& thePointOnShape,
223                          gp_XYZ&               thePoint );
224   // computes coordinates of a point in shell by points on sub-shapes
225   // and point parameters.
226   // thePointOnShape[ subShapeID ] must be a point on a subShape;
227   // thePointOnShape.size() == ID_Shell, thePointOnShape[0] not used
228
229
230  public:
231   // ---------------------------------
232   // Define parameters by coordinates
233   // ---------------------------------
234
235   bool ComputeParameters (const gp_Pnt& thePoint,
236                           gp_XYZ&       theParams,
237                           const int     theShapeID    = ID_Shell,
238                           const gp_XYZ& theParamsHint = gp_XYZ(-1,-1,-1));
239   // compute point parameters in the block.
240   // Note: for edges, it is better to use EdgeParameters()
241
242   bool VertexParameters(const int theVertexID, gp_XYZ& theParams);
243   // return parameters of a vertex given by TShapeID
244
245   bool EdgeParameters(const int theEdgeID, const double theU, gp_XYZ& theParams);
246   // return parameters of a point given by theU on edge
247
248
249  public:
250   // ---------------
251   // Block geomerty
252   // ---------------
253
254   
255   
256  public:
257   // ---------
258   // Services
259   // ---------
260
261   static bool IsForwardEdge (const TopoDS_Edge &                       theEdge,
262                              const TopTools_IndexedMapOfOrientedShape& theShapeIDMap) {
263     int v1ID = theShapeIDMap.FindIndex( TopExp::FirstVertex( theEdge ).Oriented( TopAbs_FORWARD ));
264     int v2ID = theShapeIDMap.FindIndex( TopExp::LastVertex( theEdge ).Oriented( TopAbs_FORWARD ));
265     return ( v1ID < v2ID );
266   }
267   // Return true if an in-block parameter increases along theEdge curve
268
269   static int GetOrderedEdges (const TopoDS_Face&        theFace,
270                               TopoDS_Vertex             theFirstVertex,
271                               std::list< TopoDS_Edge >& theEdges,
272                               std::list< int >  &       theNbVertexInWires);
273   // Return nb wires and a list of oredered edges.
274   // It is used to assign indices to subshapes.
275   // theFirstVertex may be NULL.
276   // Always try to set a seam edge first
277
278  public:
279   // -----------------------------------------------------------
280   // Methods of math_FunctionSetWithDerivatives used internally
281   // to define parameters by coordinates
282   // -----------------------------------------------------------
283   Standard_Integer NbVariables() const;
284   Standard_Integer NbEquations() const;
285   Standard_Boolean Value(const math_Vector& X,math_Vector& F) ;
286   Standard_Boolean Derivatives(const math_Vector& X,math_Matrix& D) ;
287   Standard_Boolean Values(const math_Vector& X,math_Vector& F,math_Matrix& D) ;
288   Standard_Integer GetStateNumber ();
289
290  protected:
291
292   /*!
293    * \brief Call it after geometry initialisation
294    */
295   void init();
296
297   // Note: to compute params of a point on a face, it is enough to set
298   // TFace, TEdge's and points for that face only
299
300   // Note 2: curve adaptors need to have only Value(double), FirstParameter() and
301   // LastParameter() defined to be used by Block algoritms
302
303   class TEdge {
304     int                myCoordInd;
305     double             myFirst;
306     double             myLast;
307     Adaptor3d_Curve*   myC3d;
308     // if mesh volume
309     gp_XYZ             myNodes[2];
310   public:
311     void Set( const int edgeID, Adaptor3d_Curve* curve, const bool isForward );
312     void Set( const int edgeID, const gp_XYZ& node1, const gp_XYZ& node2 );
313     Adaptor3d_Curve* GetCurve() const { return myC3d; }
314     double EndParam(int i) const { return i ? myLast : myFirst; }
315     int CoordInd() const { return myCoordInd; }
316     const gp_XYZ& NodeXYZ(int i) const { return i ? myNodes[1] : myNodes[0]; }
317     gp_XYZ Point( const gp_XYZ& theParams ) const; // Return coord by params
318     double GetU( const gp_XYZ& theParams ) const;  // Return U by params
319     TEdge(): myC3d(0) {}
320     ~TEdge();
321   };
322
323   class TFace {
324     // 4 edges in the order u0, u1, 0v, 1v
325     int                  myCoordInd[ 4 ];
326     double               myFirst   [ 4 ];
327     double               myLast    [ 4 ];
328     Adaptor2d_Curve2d*   myC2d     [ 4 ];
329     // 4 corner points in the order 00, 10, 11, 01
330     gp_XY                myCorner  [ 4 ];
331     // surface
332     Adaptor3d_Surface*   myS;
333     // if mesh volume
334     gp_XYZ               myNodes[4];
335   public:
336     void Set( const int faceID, Adaptor3d_Surface* S, // must be in GetFaceEdgesIDs() order:
337               Adaptor2d_Curve2d* c2d[4], const bool isForward[4] );
338     void Set( const int faceID, const TEdge& edgeU0, const TEdge& edgeU1 );
339     gp_XY  GetUV( const gp_XYZ& theParams ) const;
340     gp_XYZ Point( const gp_XYZ& theParams ) const;
341     int GetUInd() const { return myCoordInd[ 0 ]; }
342     int GetVInd() const { return myCoordInd[ 2 ]; }
343     void GetCoefs( int i, const gp_XYZ& theParams, double& eCoef, double& vCoef ) const;
344     TFace(): myS(0) { myC2d[0]=myC2d[1]=myC2d[2]=myC2d[3]=0; }
345     ~TFace();
346   };
347
348   // geometry in the order as in TShapeID:
349   // 8 vertices
350   gp_XYZ myPnt[ 8 ];
351   // 12 edges
352   TEdge  myEdge[ 12 ];
353   // 6 faces
354   TFace  myFace[ 6 ];
355
356   // for param computation
357
358   enum { SQUARE_DIST = 0, DRV_1, DRV_2, DRV_3 };
359   double distance () const { return sqrt( myValues[ SQUARE_DIST ]); }
360   double funcValue(double sqDist) const { return mySquareFunc ? sqDist : sqrt(sqDist); }
361   bool computeParameters(const gp_Pnt& thePoint, gp_XYZ& theParams, const gp_XYZ& theParamsHint);
362
363   int      myFaceIndex;
364   double   myFaceParam;
365   int      myNbIterations;
366   double   mySumDist;
367   double   myTolerance;
368   bool     mySquareFunc;
369
370   gp_XYZ   myPoint; // the given point
371   gp_XYZ   myParam; // the best parameters guess
372   double   myValues[ 4 ]; // values computed at myParam: square distance and 3 derivatives
373
374   typedef pair<gp_XYZ,gp_XYZ> TxyzPair;
375   TxyzPair my3x3x3GridNodes[ 27 ]; // to compute the first param guess
376   bool     myGridComputed;
377 };
378
379
380 #endif