Salome HOME
0022523: [CEA 1096] Add a colomn "biquadratic" in "Mesh Information"
[modules/smesh.git] / src / StdMeshers / StdMeshers_Quadrangle_2D.hxx
1 // Copyright (C) 2007-2014  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   : StdMeshers_Quadrangle_2D.hxx
23 //           Moved here from SMESH_Quadrangle_2D.hxx
24 //  Author : Paul RASCLE, EDF
25 //  Module : SMESH
26
27 #ifndef _SMESH_QUADRANGLE_2D_HXX_
28 #define _SMESH_QUADRANGLE_2D_HXX_
29
30 #include "SMESH_Algo.hxx"
31 #include "SMESH_ProxyMesh.hxx"
32 #include "SMESH_StdMeshers.hxx"
33 #include "StdMeshers_FaceSide.hxx"
34 #include "StdMeshers_QuadrangleParams.hxx"
35
36 #include <TopoDS_Face.hxx>
37 #include <Bnd_B2d.hxx>
38
39 class SMDS_MeshNode;
40 class SMESH_Mesh;
41 class SMESH_MesherHelper;
42 class SMESH_ProxyMesh;
43 struct uvPtStruct;
44
45
46 enum TSideID { QUAD_BOTTOM_SIDE=0, QUAD_RIGHT_SIDE, QUAD_TOP_SIDE, QUAD_LEFT_SIDE, NB_QUAD_SIDES };
47
48 typedef uvPtStruct UVPtStruct;
49 struct FaceQuadStruct
50 {
51   struct Side // a side of FaceQuadStruct
52   {
53     struct Contact // contact of two sides
54     {
55       int   point; // index of a grid point of this side where two sides meat
56       Side* other_side;
57       int   other_point;
58     };
59     StdMeshers_FaceSidePtr grid;
60     int                    from, to;     // indices of grid points used by the quad
61     int                    di;           // +1 or -1 depending on IsReversed()
62     std::set<int>          forced_nodes; // indices of forced grid points
63     std::vector<Contact>   contacts;     // contacts with sides of other quads
64     int                    nbNodeOut;    // nb of missing nodes on an opposite shorter side
65
66     Side(StdMeshers_FaceSidePtr theGrid = StdMeshers_FaceSidePtr());
67     Side& operator=(const Side& otherSide);
68     operator StdMeshers_FaceSidePtr() { return grid; }
69     operator const StdMeshers_FaceSidePtr() const { return grid; }
70     void AddContact( int ip, Side* side, int iop );
71     int  ToSideIndex( int quadNodeIndex ) const;
72     int  ToQuadIndex( int sideNodeIndex ) const;
73     bool IsForced( int nodeIndex ) const;
74     bool IsReversed() const { return nbNodeOut ? false : to < from; }
75     bool Reverse(bool keepGrid);
76     int  NbPoints() const { return Abs( to - from ); }
77     double Param( int nodeIndex ) const;
78     double Length( int from=-1, int to=-1) const;
79     gp_XY Value2d( double x ) const;
80     const UVPtStruct& First() const { return GetUVPtStruct()[ from ]; }
81     const UVPtStruct& Last()  const {
82       return GetUVPtStruct()[ to-nbNodeOut-(IsReversed() ? -1 : +1)];
83     }
84     // some sortcuts
85     const vector<UVPtStruct>& GetUVPtStruct(bool isXConst=0, double constValue=0) const
86     { return nbNodeOut ?
87         grid->SimulateUVPtStruct( NbPoints()-nbNodeOut-1, isXConst, constValue ) :
88         grid->GetUVPtStruct( isXConst, constValue );
89     }
90   };
91   struct SideIterator // iterator on UVPtStruct of a Side
92   {
93     const UVPtStruct *uvPtr, *uvEnd;
94     int               dPtr, counter;
95     SideIterator(): uvPtr(0), uvEnd(0), dPtr(0), counter(0) {}
96     void Init( const Side& side ) {
97       dPtr  = counter = 0;
98       uvPtr = uvEnd = 0;
99       if ( side.NbPoints() > 0 ) {
100         uvPtr = & side.First();
101         uvEnd = & side.Last();
102         dPtr  = ( uvEnd > uvPtr ) ? +1 : -1;
103         uvEnd += dPtr;
104       }
105     }
106     bool More() const { return uvPtr != uvEnd; }
107     void Next() { uvPtr += dPtr; ++counter; }
108     UVPtStruct& UVPt() const { return (UVPtStruct&) *uvPtr; }
109     UVPtStruct& operator[](int i) { return (UVPtStruct&) uvPtr[ i*dPtr]; }
110     int  Count() const { return counter; }
111   };
112
113   std::vector< Side >      side;
114   std::vector< UVPtStruct> uv_grid;
115   int                      iSize, jSize;
116   TopoDS_Face              face;
117   Bnd_B2d                  uv_box;
118   std::string              name; // to ease debugging
119
120   FaceQuadStruct ( const TopoDS_Face& F = TopoDS_Face(), const std::string& nm="main" );
121   UVPtStruct& UVPt( int i, int j ) { return uv_grid[ i + j * iSize ]; }
122   void  shift    ( size_t nb, bool keepUnitOri, bool keepGrid=false );
123   int & nbNodeOut( int iSide ) { return side[ iSide ].nbNodeOut; }
124   bool  findCell ( const gp_XY& uv, int & i, int & j );
125   bool  isNear   ( const gp_XY& uv, int & i, int & j, int nbLoops=1 );
126   bool  isEqual  ( const gp_XY& uv, int   i, int   j );
127   void  normPa2IJ( double x, double y, int & i, int & j );
128   void  updateUV ( const gp_XY& uv, int   i, int   j, bool isVertical );
129
130   typedef boost::shared_ptr<FaceQuadStruct> Ptr;
131 };
132
133 class STDMESHERS_EXPORT StdMeshers_Quadrangle_2D: public SMESH_2D_Algo
134 {
135  public:
136   StdMeshers_Quadrangle_2D(int hypId, int studyId, SMESH_Gen* gen);
137   virtual ~StdMeshers_Quadrangle_2D();
138
139   virtual bool CheckHypothesis(SMESH_Mesh&         aMesh,
140                                const TopoDS_Shape& aShape,
141                                Hypothesis_Status&  aStatus);
142
143   virtual bool Compute(SMESH_Mesh&         aMesh,
144                        const TopoDS_Shape& aShape);
145
146   virtual bool Evaluate(SMESH_Mesh &         aMesh,
147                         const TopoDS_Shape & aShape,
148                         MapShapeNbElems&     aResMap);
149
150   FaceQuadStruct::Ptr CheckAnd2Dcompute(SMESH_Mesh&         aMesh,
151                                         const TopoDS_Shape& aShape,
152                                         const bool          CreateQuadratic);
153
154   FaceQuadStruct::Ptr CheckNbEdges(SMESH_Mesh&         aMesh,
155                                    const TopoDS_Shape& aShape,
156                                    const bool          considerMesh=false);
157
158   static bool IsApplicable(const TopoDS_Shape & aShape, bool toCheckAll);
159
160  protected:
161
162   bool checkNbEdgesForEvaluate(SMESH_Mesh& aMesh,
163                                const TopoDS_Shape & aShape,
164                                MapShapeNbElems& aResMap,
165                                std::vector<int>& aNbNodes,
166                                bool& IsQuadratic);
167
168   bool setNormalizedGrid(FaceQuadStruct::Ptr quad);
169
170   void splitQuadFace(SMESHDS_Mesh *       theMeshDS,
171                      const int            theFaceID,
172                      const SMDS_MeshNode* theNode1,
173                      const SMDS_MeshNode* theNode2,
174                      const SMDS_MeshNode* theNode3,
175                      const SMDS_MeshNode* theNode4);
176
177   bool computeQuadDominant(SMESH_Mesh&         aMesh,
178                            const TopoDS_Face&  aFace);
179
180   bool computeQuadDominant(SMESH_Mesh&         aMesh,
181                            const TopoDS_Face&  aFace,
182                            FaceQuadStruct::Ptr quad);
183
184   bool computeQuadPref(SMESH_Mesh&         aMesh,
185                        const TopoDS_Face&  aFace,
186                        FaceQuadStruct::Ptr quad);
187
188   bool computeTriangles(SMESH_Mesh&         aMesh,
189                         const TopoDS_Face&  aFace,
190                         FaceQuadStruct::Ptr quad);
191
192   bool evaluateQuadPref(SMESH_Mesh&         aMesh,
193                         const TopoDS_Shape& aShape,
194                         std::vector<int>&   aNbNodes,
195                         MapShapeNbElems&    aResMap,
196                         bool                isQuadratic);
197
198   bool computeReduced (SMESH_Mesh&         aMesh,
199                        const TopoDS_Face&  aFace,
200                        FaceQuadStruct::Ptr quad);
201
202   void updateDegenUV(FaceQuadStruct::Ptr quad);
203
204   void smooth (FaceQuadStruct::Ptr quad);
205
206   bool check();
207
208   int getCorners(const TopoDS_Face&          theFace,
209                  SMESH_Mesh &                theMesh,
210                  std::list<TopoDS_Edge>&     theWire,
211                  std::vector<TopoDS_Vertex>& theVertices,
212                  int &                       theNbDegenEdges,
213                  const bool                  considerMesh);
214
215   bool getEnforcedUV();
216
217   bool addEnforcedNodes();
218
219   int splitQuad(FaceQuadStruct::Ptr quad, int i, int j);
220
221   void shiftQuad(FaceQuadStruct::Ptr& quad, const int num );
222
223   typedef std::map< StdMeshers_FaceSidePtr, std::vector< FaceQuadStruct::Ptr > > TQuadsBySide;
224   void updateSideUV( FaceQuadStruct::Side&  side,
225                      int                    iForced,
226                      const TQuadsBySide&    quads,
227                      int *                  iNext=NULL);
228
229
230  protected: // Fields
231
232   bool myQuadranglePreference;
233   bool myTrianglePreference;
234   int  myTriaVertexID;
235   bool myNeedSmooth, myCheckOri;
236   const StdMeshers_QuadrangleParams* myParams;
237   StdMeshers_QuadType                myQuadType;
238
239   SMESH_MesherHelper*                myHelper;
240   SMESH_ProxyMesh::Ptr               myProxyMesh;
241   std::list< FaceQuadStruct::Ptr >   myQuadList;
242
243   struct ForcedPoint
244   {
245     gp_XY         uv;
246     gp_XYZ        xyz;
247     TopoDS_Vertex vertex;
248
249     double U() const { return uv.X(); }
250     double V() const { return uv.Y(); }
251     operator const gp_XY& () { return uv; }
252   };
253   std::vector< ForcedPoint >         myForcedPnts;
254 };
255
256 #endif