1 // SMESH SMESH : implementaion of SMESH idl descriptions
3 // Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
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.
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.
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
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
24 // File : StdMeshers_RadialPrism_3D.cxx
26 // Created : Fri Oct 20 11:37:07 2006
27 // Author : Edward AGAPOV (eap)
30 #include "StdMeshers_RadialPrism_3D.hxx"
32 #include "StdMeshers_ProjectionUtils.hxx"
33 #include "StdMeshers_NumberOfLayers.hxx"
34 #include "StdMeshers_LayerDistribution.hxx"
35 #include "StdMeshers_Prism_3D.hxx"
36 #include "StdMeshers_Regular_1D.hxx"
38 #include "SMDS_MeshNode.hxx"
39 #include "SMESHDS_SubMesh.hxx"
40 #include "SMESH_Gen.hxx"
41 #include "SMESH_Mesh.hxx"
42 #include "SMESH_MeshEditor.hxx"
43 #include "SMESH_MesherHelper.hxx"
44 #include "SMESH_subMesh.hxx"
46 #include "utilities.h"
48 #include <TopoDS_Solid.hxx>
49 #include <TopoDS_Shell.hxx>
50 #include <BRepTools.hxx>
51 #include <BRepAdaptor_Curve.hxx>
52 #include <TopTools_ListIteratorOfListOfShape.hxx>
53 #include <BRepBuilderAPI_MakeEdge.hxx>
60 #define RETURN_BAD_RESULT(msg) { MESSAGE(msg); return false; }
61 #define gpXYZ(n) gp_XYZ(n->X(),n->Y(),n->Z())
63 typedef StdMeshers_ProjectionUtils TAssocTool;
65 //=======================================================================
66 //function : StdMeshers_RadialPrism_3D
68 //=======================================================================
70 StdMeshers_RadialPrism_3D::StdMeshers_RadialPrism_3D(int hypId, int studyId, SMESH_Gen* gen)
71 :SMESH_3D_Algo(hypId, studyId, gen)
73 _name = "RadialPrism_3D";
74 _shapeType = (1 << TopAbs_SOLID); // 1 bit per shape type
76 _compatibleHypothesis.push_back("LayerDistribution");
77 _compatibleHypothesis.push_back("NumberOfLayers");
79 myDistributionHypo = 0;
82 //================================================================================
86 //================================================================================
88 StdMeshers_RadialPrism_3D::~StdMeshers_RadialPrism_3D()
91 //=======================================================================
92 //function : CheckHypothesis
94 //=======================================================================
96 bool StdMeshers_RadialPrism_3D::CheckHypothesis(SMESH_Mesh& aMesh,
97 const TopoDS_Shape& aShape,
98 SMESH_Hypothesis::Hypothesis_Status& aStatus)
100 // check aShape that must have 2 shells
101 if ( TAssocTool::Count( aShape, TopAbs_SOLID, 0 ) != 1 ||
102 TAssocTool::Count( aShape, TopAbs_SHELL, 0 ) != 2 )
104 aStatus = HYP_BAD_GEOMETRY;
109 myDistributionHypo = 0;
111 list <const SMESHDS_Hypothesis * >::const_iterator itl;
113 const list <const SMESHDS_Hypothesis * >&hyps = GetUsedHypothesis(aMesh, aShape);
114 if ( hyps.size() == 0 )
116 aStatus = SMESH_Hypothesis::HYP_MISSING;
117 return false; // can't work with no hypothesis
120 if ( hyps.size() > 1 )
122 aStatus = SMESH_Hypothesis::HYP_ALREADY_EXIST;
126 const SMESHDS_Hypothesis *theHyp = hyps.front();
128 string hypName = theHyp->GetName();
130 if (hypName == "NumberOfLayers")
132 myNbLayerHypo = static_cast<const StdMeshers_NumberOfLayers *>(theHyp);
133 aStatus = SMESH_Hypothesis::HYP_OK;
136 if (hypName == "LayerDistribution")
138 myDistributionHypo = static_cast<const StdMeshers_LayerDistribution *>(theHyp);
139 aStatus = SMESH_Hypothesis::HYP_OK;
142 aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
146 //=======================================================================
149 //=======================================================================
151 bool StdMeshers_RadialPrism_3D::Compute(SMESH_Mesh& aMesh, const TopoDS_Shape& aShape)
154 SMESHDS_Mesh * meshDS = aMesh.GetMeshDS();
156 myHelper = new SMESH_MesherHelper( aMesh );
157 myHelper->IsQuadraticSubMesh( aShape );
158 // to delete helper at exit from Compute()
159 std::auto_ptr<SMESH_MesherHelper> helperDeleter( myHelper );
162 TopoDS_Solid solid = TopoDS::Solid( aShape );
163 TopoDS_Shell outerShell = BRepTools::OuterShell( solid );
164 TopoDS_Shape innerShell;
166 for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells )
167 if ( !outerShell.IsSame( It.Value() ))
168 innerShell = It.Value();
170 return error(COMPERR_BAD_SHAPE, SMESH_Comment("Must be 2 shells but not")<<nbShells);
172 // ----------------------------------
173 // Associate subshapes of the shells
174 // ----------------------------------
176 TAssocTool::TShapeShapeMap shape2ShapeMap;
177 if ( !TAssocTool::FindSubShapeAssociation( outerShell, &aMesh,
180 return error(COMPERR_BAD_SHAPE,"Topology of inner and outer shells seems different" );
182 // ------------------
184 // ------------------
186 TNode2ColumnMap node2columnMap;
187 myLayerPositions.clear();
189 for ( exp.Init( outerShell, TopAbs_FACE ); exp.More(); exp.Next() )
191 // Corresponding subshapes
192 TopoDS_Face outFace = TopoDS::Face( exp.Current() );
194 if ( !shape2ShapeMap.IsBound( outFace )) {
195 return error(dfltErr(),SMESH_Comment("Corresponding inner face not found for face #" )
196 << meshDS->ShapeToIndex( outFace ));
198 inFace = TopoDS::Face( shape2ShapeMap( outFace ));
201 // Find matching nodes of in and out faces
202 TNodeNodeMap nodeIn2OutMap;
203 if ( ! TAssocTool::FindMatchingNodesOnFaces( inFace, &aMesh, outFace, &aMesh,
204 shape2ShapeMap, nodeIn2OutMap ))
205 return error(COMPERR_BAD_INPUT_MESH,SMESH_Comment("Mesh on faces #")
206 << meshDS->ShapeToIndex( outFace ) << " and "
207 << meshDS->ShapeToIndex( inFace ) << " seems different" );
211 SMDS_ElemIteratorPtr faceIt = meshDS->MeshElements( inFace )->GetElements();
212 while ( faceIt->more() ) // loop on faces on inFace
214 const SMDS_MeshElement* face = faceIt->next();
215 if ( !face || face->GetType() != SMDSAbs_Face )
217 int nbNodes = face->NbNodes();
218 if ( face->IsQuadratic() )
221 // find node columns for each node
222 vector< const TNodeColumn* > columns( nbNodes );
223 for ( int i = 0; i < nbNodes; ++i )
225 const SMDS_MeshNode* n = face->GetNode( i );
226 TNode2ColumnMap::iterator n_col = node2columnMap.find( n );
227 if ( n_col != node2columnMap.end() )
228 columns[ i ] = & n_col->second;
230 columns[ i ] = makeNodeColumn( node2columnMap, n, nodeIn2OutMap[ n ] );
233 StdMeshers_Prism_3D::AddPrisms( columns, myHelper );
235 } // loop on faces of out shell
240 //================================================================================
242 * \brief Create a column of nodes from outNode to inNode
243 * \param n2ColMap - map of node columns to add a created column
244 * \param outNode - botton node of a column
245 * \param inNode - top node of a column
246 * \retval const TNodeColumn* - a new column pointer
248 //================================================================================
250 TNodeColumn* StdMeshers_RadialPrism_3D::makeNodeColumn( TNode2ColumnMap& n2ColMap,
251 const SMDS_MeshNode* outNode,
252 const SMDS_MeshNode* inNode)
254 SMESHDS_Mesh * meshDS = myHelper->GetMeshDS();
255 int shapeID = myHelper->GetSubShapeID();
257 if ( myLayerPositions.empty() ) {
258 gp_Pnt pIn = gpXYZ( inNode ), pOut = gpXYZ( outNode );
259 computeLayerPositions( pIn, pOut );
261 int nbSegments = myLayerPositions.size() + 1;
263 TNode2ColumnMap::iterator n_col =
264 n2ColMap.insert( make_pair( outNode, TNodeColumn() )).first;
265 TNodeColumn & column = n_col->second;
266 column.resize( nbSegments + 1 );
267 column.front() = outNode;
268 column.back() = inNode;
270 gp_XYZ p1 = gpXYZ( outNode );
271 gp_XYZ p2 = gpXYZ( inNode );
272 for ( int z = 1; z < nbSegments; ++z )
274 double r = myLayerPositions[ z - 1 ];
275 gp_XYZ p = ( 1 - r ) * p1 + r * p2;
276 SMDS_MeshNode* n = meshDS->AddNode( p.X(), p.Y(), p.Z() );
277 meshDS->SetNodeInVolume( n, shapeID );
284 //================================================================================
285 //================================================================================
287 * \brief Class computing layers distribution using data of
288 * StdMeshers_LayerDistribution hypothesis
290 //================================================================================
291 //================================================================================
293 class TNodeDistributor: public StdMeshers_Regular_1D
295 list <const SMESHDS_Hypothesis *> myUsedHyps;
297 // -----------------------------------------------------------------------------
298 static TNodeDistributor* GetDistributor(SMESH_Mesh& aMesh)
300 const int myID = -1000;
301 map < int, SMESH_1D_Algo * > & algoMap = aMesh.GetGen()->_map1D_Algo;
302 map < int, SMESH_1D_Algo * >::iterator id_algo = algoMap.find( myID );
303 if ( id_algo == algoMap.end() )
304 return new TNodeDistributor( myID, 0, aMesh.GetGen() );
305 return static_cast< TNodeDistributor* >( id_algo->second );
307 // -----------------------------------------------------------------------------
308 bool Compute( vector< double > & positions,
312 const StdMeshers_LayerDistribution* hyp)
314 double len = pIn.Distance( pOut );
315 if ( len <= DBL_MIN ) return error(dfltErr(),"Too close points of inner and outer shells");
317 if ( !hyp || !hyp->GetLayerDistribution() )
318 return error(dfltErr(), "Invalid LayerDistribution hypothesis");
320 myUsedHyps.push_back( hyp->GetLayerDistribution() );
322 TopoDS_Edge edge = BRepBuilderAPI_MakeEdge( pIn, pOut );
323 SMESH_Hypothesis::Hypothesis_Status aStatus;
324 if ( !StdMeshers_Regular_1D::CheckHypothesis( aMesh, edge, aStatus ))
325 return error(dfltErr(), "StdMeshers_Regular_1D::CheckHypothesis() failed"
326 "with LayerDistribution hypothesis");
328 BRepAdaptor_Curve C3D(edge);
329 double f = C3D.FirstParameter(), l = C3D.LastParameter();
330 list< double > params;
331 if ( !StdMeshers_Regular_1D::computeInternalParameters( C3D, len, f, l, params, false ))
332 return error(dfltErr(),"StdMeshers_Regular_1D failed to compute layers distribution");
335 positions.reserve( params.size() );
336 for (list<double>::iterator itU = params.begin(); itU != params.end(); itU++)
337 positions.push_back( *itU / len );
341 // -----------------------------------------------------------------------------
342 TNodeDistributor( int hypId, int studyId, SMESH_Gen* gen)
343 : StdMeshers_Regular_1D( hypId, studyId, gen)
346 // -----------------------------------------------------------------------------
347 virtual const list <const SMESHDS_Hypothesis *> &
348 GetUsedHypothesis(SMESH_Mesh &, const TopoDS_Shape &, const bool)
352 // -----------------------------------------------------------------------------
355 //================================================================================
357 * \brief Compute positions of nodes between the internal and the external surfaces
358 * \retval bool - is a success
360 //================================================================================
362 bool StdMeshers_RadialPrism_3D::computeLayerPositions(const gp_Pnt& pIn,
367 int nbSegments = myNbLayerHypo->GetNumberOfLayers();
368 myLayerPositions.resize( nbSegments - 1 );
369 for ( int z = 1; z < nbSegments; ++z )
370 myLayerPositions[ z - 1 ] = double( z )/ double( nbSegments );
373 if ( myDistributionHypo ) {
374 SMESH_Mesh * mesh = myHelper->GetMesh();
375 if ( !TNodeDistributor::GetDistributor(*mesh)->Compute( myLayerPositions, pIn, pOut,
376 *mesh, myDistributionHypo ))
378 error( TNodeDistributor::GetDistributor(*mesh)->GetComputeError() );
382 RETURN_BAD_RESULT("Bad hypothesis");