SALOME - SMESH
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
StdMeshers_Prism_3D.hxx
Go to the documentation of this file.
00001 //  Copyright (C) 2007-2008  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 //  Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
00004 //  CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
00005 //
00006 //  This library is free software; you can redistribute it and/or
00007 //  modify it under the terms of the GNU Lesser General Public
00008 //  License as published by the Free Software Foundation; either
00009 //  version 2.1 of the License.
00010 //
00011 //  This library is distributed in the hope that it will be useful,
00012 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014 //  Lesser General Public License for more details.
00015 //
00016 //  You should have received a copy of the GNU Lesser General Public
00017 //  License along with this library; if not, write to the Free Software
00018 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00019 //
00020 //  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00021 //
00022 //  SMESH SMESH : implementaion of SMESH idl descriptions
00023 //  File   : StdMeshers_Prism_3D.hxx
00024 //  Module : SMESH
00025 //
00026 #ifndef _SMESH_Prism_3D_HXX_
00027 #define _SMESH_Prism_3D_HXX_
00028 
00029 #include "SMESH_StdMeshers.hxx"
00030 
00031 #include "SMESH_3D_Algo.hxx"
00032 #include "SMDS_TypeOfPosition.hxx"
00033 #include "SMDS_MeshNode.hxx"
00034 #include "SMESH_Block.hxx"
00035 #include "SMESH_Mesh.hxx"
00036 #include "SMESHDS_Mesh.hxx"
00037 #include "SMESH_subMesh.hxx"
00038 #include "SMESH_MesherHelper.hxx"
00039 #include "SMESH_Comment.hxx"
00040 
00041 #include <vector>
00042 
00043 #include <Adaptor3d_Curve.hxx>
00044 #include <Adaptor3d_Surface.hxx>
00045 #include <Adaptor2d_Curve2d.hxx>
00046 #include <BRepAdaptor_Surface.hxx>
00047 #include <TopTools_IndexedMapOfOrientedShape.hxx>
00048 #include <gp_XYZ.hxx>
00049 
00050 
00051 class SMESHDS_SubMesh;
00052 class TopoDS_Edge;
00053 class TopoDS_Faces;
00054 struct TNode;
00055 
00056 //typedef std::map<const SMDS_MeshNode*, const SMDS_MeshNode*> TNodeNodeMap;
00057 typedef std::vector<const SMDS_MeshNode* > TNodeColumn;
00058 
00059 // map of bottom nodes to the column of nodes above them
00060 // (the column includes the bottom nodes)
00061 typedef std::map< TNode, TNodeColumn > TNode2ColumnMap;
00062 typedef std::map< double, TNodeColumn > TParam2ColumnMap;
00063 typedef std::map< double, TNodeColumn >::const_iterator TParam2ColumnIt;
00064 
00065 typedef TopTools_IndexedMapOfOrientedShape TBlockShapes;
00066 
00067 // ===============================================
00071 // ===============================================
00072 
00073 struct TNode
00074 {
00075   const SMDS_MeshNode* myNode;
00076   gp_XYZ               myParams;
00077 
00078   gp_XYZ GetCoords() const { return gp_XYZ( myNode->X(), myNode->Y(), myNode->Z() ); }
00079   gp_XYZ GetParams() const { return myParams; }
00080   gp_XYZ& ChangeParams() { return myParams; }
00081   bool HasParams() const { return myParams.X() >= 0.0; }
00082   SMDS_TypeOfPosition GetPositionType() const
00083   { return myNode ? myNode->GetPosition()->GetTypeOfPosition() : SMDS_TOP_UNSPEC; }
00084   bool IsNeighbor( const TNode& other ) const;
00085 
00086   TNode(const SMDS_MeshNode* node = 0): myNode(node), myParams(-1,-1,-1) {}
00087   bool operator < (const TNode& other) const { return myNode->GetID() < other.myNode->GetID(); }
00088 };
00089 
00090 // ===============================================================
00097 // ===============================================================
00098 
00099 class STDMESHERS_EXPORT StdMeshers_PrismAsBlock: public SMESH_Block
00100 {
00101 public:
00105   StdMeshers_PrismAsBlock();
00106 
00107   ~StdMeshers_PrismAsBlock();
00108 
00118   bool Init(SMESH_MesherHelper* helper, const TopoDS_Shape& shape3D);
00119 
00123   SMESH_ComputeErrorPtr GetError() const { return myError; }
00124 
00129   int VerticalSize() const { return myParam2ColumnMaps[0].begin()->second.size(); }
00130 
00131   bool HasNotQuadElemOnTop() const { return myNotQuadOnTop; }
00132 
00138   const TNodeColumn* GetNodeColumn(const SMDS_MeshNode* node) const;
00139 
00146   const TParam2ColumnMap& GetParam2ColumnMap(const int baseEdgeID,
00147                                              bool &    isReverse)
00148   {
00149     std::pair< TParam2ColumnMap*, bool > & col_frw =
00150       myShapeIndex2ColumnMap[ baseEdgeID ];
00151     isReverse = !col_frw.second;
00152     return * col_frw.first;
00153   }
00154   
00159   SMESH_Mesh* Mesh() const { return myHelper->GetMesh(); }
00160 
00165   SMESHDS_Mesh* MeshDS() const { return Mesh()->GetMeshDS(); }
00166 
00172   SMESH_subMesh* SubMesh(const int shapeID) const
00173   { return Mesh()->GetSubMesh( Shape( shapeID )); }
00174 
00180   SMESHDS_SubMesh* SubMeshDS(const int shapeID) const
00181   { return SubMesh(shapeID)->GetSubMeshDS(); }
00182 
00188   const TopoDS_Shape& Shape(const int shapeID) const
00189   { return myShapeIDMap( shapeID ); }
00190 
00196   int ShapeID(const TopoDS_Shape& shape) const
00197   { return myShapeIDMap.FindIndex( shape ); }
00198 
00207   static bool IsForwardEdge(SMESHDS_Mesh*           meshDS,
00208                             const TParam2ColumnMap& columnsMap,
00209                             const TopoDS_Edge &     bottomEdge,
00210                             const int               sideFaceID);
00219   static bool GetWallFaces( SMESH_Mesh*                     mesh,
00220                             const TopoDS_Shape &            mainShape,
00221                             const TopoDS_Shape &            bottomFace,
00222                             const std::list< TopoDS_Edge >& bottomEdges,
00223                             std::list< TopoDS_Face >&       wallFaces);
00224 
00225 private:
00226 
00227   // --------------------------------------------------------------------
00235   // --------------------------------------------------------------------
00236   class TSideFace: public Adaptor3d_Surface
00237   {
00238     int                             myID; 
00239     // map used to find out real UV by it's normalized UV
00240     TParam2ColumnMap*               myParamToColumnMap;
00241     BRepAdaptor_Surface             mySurface;
00242     TopoDS_Edge                     myBaseEdge;
00243     // first and last normalized params and orientaion for each component or it-self
00244     std::vector< std::pair< double, double> > myParams;
00245     bool                            myIsForward;
00246     std::vector< TSideFace* >       myComponents;
00247     SMESH_MesherHelper *            myHelper;
00248   public:
00249     TSideFace( SMESH_MesherHelper* helper,
00250                const int           faceID,
00251                const TopoDS_Face&  face,
00252                const TopoDS_Edge&  baseEdge,
00253                TParam2ColumnMap*   columnsMap,
00254                const double        first = 0.0,
00255                const double        last = 1.0);
00256     TSideFace( const std::vector< TSideFace* >&             components,
00257                const std::vector< std::pair< double, double> > & params);
00258     TSideFace( const TSideFace& other );
00259     ~TSideFace();
00260     bool IsComplex() const
00261     { return ( NbComponents() > 0 || myParams[0].first != 0. || myParams[0].second != 1. ); }
00262     int FaceID() const { return myID; }
00263     TParam2ColumnMap* GetColumns() const { return myParamToColumnMap; }
00264     gp_XY GetNodeUV(const TopoDS_Face& F, const SMDS_MeshNode* n) const
00265     { return myHelper->GetNodeUV( F, n ); }
00266     const TopoDS_Edge & BaseEdge() const { return myBaseEdge; }
00267     int ColumnHeight() const {
00268       if ( NbComponents() ) return GetComponent(0)->GetColumns()->begin()->second.size();
00269       else                  return GetColumns()->begin()->second.size(); }
00270     double GetColumns(const double U, TParam2ColumnIt & col1, TParam2ColumnIt& col2 ) const;
00271     int NbComponents() const { return myComponents.size(); }
00272     TSideFace* GetComponent(const int i) const { return myComponents.at( i ); }
00273     void SetComponent(const int i, TSideFace* c)
00274     { if ( myComponents[i] ) delete myComponents[i]; myComponents[i]=c; }
00275     TSideFace* GetComponent(const double U, double& localU) const;
00276     bool IsForward() const { return myIsForward; }
00277     // boundary geometry for a face
00278     Adaptor3d_Surface* Surface() const { return new TSideFace( *this ); }
00279     bool GetPCurves(Adaptor2d_Curve2d* pcurv[4]) const;
00280     Adaptor2d_Curve2d* HorizPCurve(const bool isTop, const TopoDS_Face& horFace) const;
00281     Adaptor3d_Curve* HorizCurve(const bool isTop) const;
00282     Adaptor3d_Curve* VertiCurve(const bool isMax) const;
00283     TopoDS_Edge GetEdge( const int edge ) const;
00284     int InsertSubShapes( TBlockShapes& shapeMap ) const;
00285     // redefine Adaptor methods
00286     gp_Pnt Value(const Standard_Real U,const Standard_Real V) const;
00287   };
00288 
00289   // --------------------------------------------------------------------
00293   // --------------------------------------------------------------------
00294   class STDMESHERS_EXPORT TVerticalEdgeAdaptor: public Adaptor3d_Curve
00295   {
00296     const TNodeColumn* myNodeColumn;
00297   public:
00298     TVerticalEdgeAdaptor( const TParam2ColumnMap* columnsMap, const double parameter );
00299     gp_Pnt Value(const Standard_Real U) const;
00300     Standard_Real FirstParameter() const { return 0; }
00301     Standard_Real LastParameter() const { return 1; }
00302   };
00303 
00304   // --------------------------------------------------------------------
00308   // --------------------------------------------------------------------
00309   class STDMESHERS_EXPORT THorizontalEdgeAdaptor: public Adaptor3d_Curve
00310   {
00311     const TSideFace* mySide;
00312     double           myV;
00313   public:
00314     THorizontalEdgeAdaptor( const TSideFace* sideFace, const bool isTop)
00315       :mySide(sideFace), myV( isTop ? 1.0 : 0.0 ) {}
00316     gp_Pnt Value(const Standard_Real U) const;
00317     Standard_Real FirstParameter() const { return 0; }
00318     Standard_Real LastParameter() const { return 1; }
00319   };
00320 
00321   // --------------------------------------------------------------------
00325   // --------------------------------------------------------------------
00326   class STDMESHERS_EXPORT TPCurveOnHorFaceAdaptor: public Adaptor2d_Curve2d
00327   {
00328     const TSideFace*  mySide;
00329     int               myZ;
00330     TopoDS_Face       myFace;
00331   public:
00332     TPCurveOnHorFaceAdaptor( const TSideFace*   sideFace,
00333                              const bool         isTop,
00334                              const TopoDS_Face& horFace)
00335       : mySide(sideFace), myFace(horFace), myZ(isTop ? mySide->ColumnHeight() - 1 : 0 ) {}
00336     gp_Pnt2d Value(const Standard_Real U) const;
00337     Standard_Real FirstParameter() const { return 0; }
00338     Standard_Real LastParameter() const { return 1; }
00339   };
00340   // --------------------------------------------------------------------
00341 
00342   bool myNotQuadOnTop;
00343   SMESH_MesherHelper* myHelper;
00344   TBlockShapes myShapeIDMap;
00345 
00346   // container of 4 side faces
00347   TSideFace*                 mySide; 
00348   // node columns for each base edge
00349   std::vector< TParam2ColumnMap > myParam2ColumnMaps;
00350   // to find a column for a node by edge SMESHDS Index
00351   std::map< int, std::pair< TParam2ColumnMap*, bool > > myShapeIndex2ColumnMap;
00352 
00353   SMESH_ComputeErrorPtr myError;
00357   bool error(int error, const SMESH_Comment& comment = "") {
00358     myError = SMESH_ComputeError::New(error,comment);
00359     return myError->IsOK();
00360   }
00361   //std::vector< SMESH_subMesh* >           mySubMeshesVec; // submesh by in-block id
00362 };
00363 
00364 // =============================================
00368 // =============================================
00369 
00370 class STDMESHERS_EXPORT StdMeshers_Prism_3D: public SMESH_3D_Algo
00371 {
00372 public:
00373   StdMeshers_Prism_3D(int hypId, int studyId, SMESH_Gen* gen);
00374   virtual ~StdMeshers_Prism_3D();
00375 
00376   virtual bool CheckHypothesis(SMESH_Mesh&                          aMesh,
00377                                const TopoDS_Shape&                  aShape,
00378                                SMESH_Hypothesis::Hypothesis_Status& aStatus);
00379 
00380   virtual bool Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape);
00381 
00382   virtual bool Evaluate(SMESH_Mesh & aMesh, const TopoDS_Shape & aShape,
00383                         MapShapeNbElems& aResMap);
00384 
00393   void ProjectTriangles() { myProjectTriangles = true; }
00394 
00400   static void AddPrisms( std::vector<const TNodeColumn*> & nodeColumns,
00401                          SMESH_MesherHelper*          helper);
00402 
00403 private:
00404 
00411   bool assocOrProjBottom2Top();
00412 
00418   bool projectBottomToTop();
00419 
00426   bool setFaceAndEdgesXYZ( const int faceID, const gp_XYZ& params, int z );
00427 
00428 private:
00429 
00430   bool myProjectTriangles;
00431 
00432   StdMeshers_PrismAsBlock myBlock;
00433   SMESH_MesherHelper*     myHelper;
00434 
00435   std::vector<gp_XYZ>     myShapeXYZ; // point on each sub-shape
00436 
00437   // map of bottom nodes to the column of nodes above them
00438   // (the column includes the bottom node)
00439   TNode2ColumnMap         myBotToColumnMap;
00440 };
00441 
00442 #endif