1 // File: StdMeshers_Helper.cxx
2 // Created: 15.02.06 15:22:41
4 // Copyright: Open CASCADE 2006
7 #include "StdMeshers_Helper.hxx"
9 #include "SMDS_FacePosition.hxx"
10 #include "SMDS_EdgePosition.hxx"
11 #include "SMESH_subMesh.hxx"
12 #include "SMESH_MeshEditor.hxx"
14 #include <BRepAdaptor_Surface.hxx>
15 #include <BRepTools.hxx>
16 #include <BRep_Tool.hxx>
17 #include <Geom2d_Curve.hxx>
18 #include <Geom_Curve.hxx>
19 #include <Geom_Surface.hxx>
20 #include <TopExp_Explorer.hxx>
21 #include <TopTools_MapOfShape.hxx>
22 #include <gp_Pnt2d.hxx>
25 //=======================================================================
26 //function : CheckShape
28 //=======================================================================
30 bool StdMeshers_Helper::IsQuadraticSubMesh(const TopoDS_Shape& aSh)
32 SMESHDS_Mesh* meshDS = GetMesh()->GetMeshDS();
33 myShapeID = meshDS->ShapeToIndex(aSh);
34 // we can create quadratic elements only if all elements
35 // created on subshapes of given shape are quadratic
36 // also we have to fill myNLinkNodeMap
37 myCreateQuadratic = true;
38 mySeamShapeIds.clear();
39 TopAbs_ShapeEnum subType( aSh.ShapeType()==TopAbs_FACE ? TopAbs_EDGE : TopAbs_FACE );
40 SMDSAbs_ElementType elemType( subType==TopAbs_FACE ? SMDSAbs_Face : SMDSAbs_Edge );
42 TopExp_Explorer exp( aSh, subType );
43 for (; exp.More() && myCreateQuadratic; exp.Next()) {
44 if ( SMESHDS_SubMesh * subMesh = meshDS->MeshElements( exp.Current() )) {
45 if ( SMDS_ElemIteratorPtr it = subMesh->GetElements() ) {
47 const SMDS_MeshElement* e = it->next();
48 if ( e->GetType() != elemType || !e->IsQuadratic() ) {
49 myCreateQuadratic = false;
54 switch ( e->NbNodes() ) {
56 AddNLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(2)); break;
58 AddNLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(3));
59 AddNLinkNode(e->GetNode(1),e->GetNode(2),e->GetNode(4));
60 AddNLinkNode(e->GetNode(2),e->GetNode(0),e->GetNode(5)); break;
62 AddNLinkNode(e->GetNode(0),e->GetNode(1),e->GetNode(4));
63 AddNLinkNode(e->GetNode(1),e->GetNode(2),e->GetNode(5));
64 AddNLinkNode(e->GetNode(2),e->GetNode(3),e->GetNode(6));
65 AddNLinkNode(e->GetNode(3),e->GetNode(0),e->GetNode(7));
68 myCreateQuadratic = false;
77 if(!myCreateQuadratic) {
78 myNLinkNodeMap.clear();
81 // treatment of periodic faces
82 if ( aSh.ShapeType() == TopAbs_FACE ) {
83 const TopoDS_Face& face = TopoDS::Face( aSh );
84 BRepAdaptor_Surface surface( face );
85 if ( surface.IsUPeriodic() || surface.IsVPeriodic() ) {
86 // look for a seam edge
87 for ( exp.Init( face, TopAbs_EDGE ); exp.More(); exp.Next()) {
88 const TopoDS_Edge& edge = TopoDS::Edge( exp.Current() );
89 if ( BRep_Tool::IsClosed( edge, face )) {
90 // initialize myPar1, myPar2 and myParIndex
91 if ( mySeamShapeIds.empty() ) {
93 BRep_Tool::UVPoints( edge, face, uv1, uv2 );
94 if ( Abs( uv1.Coord(1) - uv2.Coord(1) ) < Abs( uv1.Coord(2) - uv2.Coord(2) ))
96 myParIndex = 1; // U periodic
97 myPar1 = surface.FirstUParameter();
98 myPar2 = surface.LastUParameter();
101 myParIndex = 2; // V periodic
102 myPar1 = surface.FirstVParameter();
103 myPar2 = surface.LastVParameter();
106 // store shapes indices
107 mySeamShapeIds.insert( meshDS->ShapeToIndex( exp.Current() ));
108 for ( TopExp_Explorer v( exp.Current(), TopAbs_VERTEX ); v.More(); v.Next() )
109 mySeamShapeIds.insert( meshDS->ShapeToIndex( v.Current() ));
115 return myCreateQuadratic;
119 //=======================================================================
120 //function : IsMedium
122 //=======================================================================
124 bool StdMeshers_Helper::IsMedium(const SMDS_MeshNode* node,
125 const SMDSAbs_ElementType typeToCheck)
127 return SMESH_MeshEditor::IsMedium( node, typeToCheck );
130 //=======================================================================
131 //function : AddNLinkNode
133 //=======================================================================
135 * Auxilary function for filling myNLinkNodeMap
137 void StdMeshers_Helper::AddNLinkNode(const SMDS_MeshNode* n1,
138 const SMDS_MeshNode* n2,
139 const SMDS_MeshNode* n12)
141 NLink link( n1, n2 );
142 if ( n1 > n2 ) link = NLink( n2, n1 );
143 // add new record to map
144 myNLinkNodeMap.insert( make_pair(link,n12));
147 //=======================================================================
149 * \brief Select UV on either of 2 pcurves of a seam edge, closest to the given UV
150 * \param uv1 - UV on the seam
151 * \param uv2 - UV within a face
152 * \retval gp_Pnt2d - selected UV
154 //=======================================================================
156 gp_Pnt2d StdMeshers_Helper::GetUVOnSeam( const gp_Pnt2d& uv1, const gp_Pnt2d& uv2 ) const
158 double p1 = uv1.Coord( myParIndex );
159 double p2 = uv2.Coord( myParIndex );
160 double p3 = ( Abs( p1 - myPar1 ) < Abs( p1 - myPar2 )) ? myPar2 : myPar1;
161 if ( Abs( p2 - p1 ) > Abs( p2 - p3 ))
163 gp_Pnt2d result = uv1;
164 result.SetCoord( myParIndex, p1 );
168 //=======================================================================
170 * \brief Return node UV on face
171 * \param F - the face
172 * \param n - the node
173 * \param n2 - a medium node will be placed between n and n2
174 * \retval gp_XY - resulting UV
176 * Auxilary function called form GetMediumNode()
178 //=======================================================================
180 gp_XY StdMeshers_Helper::GetNodeUV(const TopoDS_Face& F,
181 const SMDS_MeshNode* n,
182 const SMDS_MeshNode* n2)
185 const SMDS_PositionPtr Pos = n->GetPosition();
186 if(Pos->GetTypeOfPosition()==SMDS_TOP_FACE) {
187 // node has position on face
188 const SMDS_FacePosition* fpos =
189 static_cast<const SMDS_FacePosition*>(n->GetPosition().get());
190 uv = gp_Pnt2d(fpos->GetUParameter(),fpos->GetVParameter());
192 else if(Pos->GetTypeOfPosition()==SMDS_TOP_EDGE) {
193 // node has position on edge => it is needed to find
194 // corresponding edge from face, get pcurve for this
195 // edge and recieve value from this pcurve
196 const SMDS_EdgePosition* epos =
197 static_cast<const SMDS_EdgePosition*>(n->GetPosition().get());
198 SMESHDS_Mesh* meshDS = GetMesh()->GetMeshDS();
199 int edgeID = Pos->GetShapeId();
200 TopoDS_Edge E = TopoDS::Edge(meshDS->IndexToShape(edgeID));
203 Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface(E, F, f, l);
204 uv = C2d->Value( epos->GetUParameter() );
205 // for a node on a seam edge select one of UVs on 2 pcurves
206 if ( n2 && mySeamShapeIds.find( edgeID ) != mySeamShapeIds.end() )
207 uv = GetUVOnSeam( uv, GetNodeUV( F, n2, 0 ));
209 else if(Pos->GetTypeOfPosition()==SMDS_TOP_VERTEX) {
210 SMESHDS_Mesh * meshDS = GetMesh()->GetMeshDS();
211 int vertexID = n->GetPosition()->GetShapeId();
212 const TopoDS_Vertex& V = TopoDS::Vertex(meshDS->IndexToShape(vertexID));
213 uv = BRep_Tool::Parameters( V, F );
214 if ( n2 && mySeamShapeIds.find( vertexID ) != mySeamShapeIds.end() )
215 uv = GetUVOnSeam( uv, GetNodeUV( F, n2, 0 ));
221 //=======================================================================
222 //function : GetMediumNode
224 //=======================================================================
226 * Special function for search or creation medium node
228 const SMDS_MeshNode* StdMeshers_Helper::GetMediumNode(const SMDS_MeshNode* n1,
229 const SMDS_MeshNode* n2,
234 NLink link(( n1 < n2 ? n1 : n2 ), ( n1 < n2 ? n2 : n1 ));
235 ItNLinkNode itLN = myNLinkNodeMap.find( link );
236 if ( itLN != myNLinkNodeMap.end() ) {
237 return (*itLN).second;
240 // create medium node
242 SMESHDS_Mesh* meshDS = GetMesh()->GetMeshDS();
244 // we try to create medium node using UV parameters of
245 // nodes, else - medium between corresponding 3d points
246 const SMDS_PositionPtr Pos1 = n1->GetPosition();
247 const SMDS_PositionPtr Pos2 = n2->GetPosition();
249 if( Pos1->GetTypeOfPosition()==SMDS_TOP_FACE ) {
250 faceID = Pos1->GetShapeId();
252 else if( Pos2->GetTypeOfPosition()==SMDS_TOP_FACE ) {
253 faceID = Pos2->GetShapeId();
256 TopoDS_Face F = TopoDS::Face(meshDS->IndexToShape(faceID));
257 gp_XY p1 = GetNodeUV(F,n1,n2);
258 gp_XY p2 = GetNodeUV(F,n2,n1);
259 double u = (p1.X()+p2.X())/2.;
260 double v = (p1.Y()+p2.Y())/2.;
261 Handle(Geom_Surface) S = BRep_Tool::Surface(F);
262 gp_Pnt P = S->Value(u, v);
263 n12 = meshDS->AddNode(P.X(), P.Y(), P.Z());
264 meshDS->SetNodeOnFace(n12, faceID, u, v);
265 myNLinkNodeMap.insert(NLinkNodeMap::value_type(link,n12));
270 double x = ( n1->X() + n2->X() )/2.;
271 double y = ( n1->Y() + n2->Y() )/2.;
272 double z = ( n1->Z() + n2->Z() )/2.;
273 n12 = meshDS->AddNode(x,y,z);
274 meshDS->SetNodeInVolume(n12, myShapeID);
275 myNLinkNodeMap.insert(NLinkNodeMap::value_type(link,n12));
281 //=======================================================================
284 //=======================================================================
286 * Special function for creation quadratic triangle
288 SMDS_MeshFace* StdMeshers_Helper::AddFace(const SMDS_MeshNode* n1,
289 const SMDS_MeshNode* n2,
290 const SMDS_MeshNode* n3)
292 SMESHDS_Mesh * meshDS = GetMesh()->GetMeshDS();
293 if(!myCreateQuadratic) {
294 return meshDS->AddFace(n1, n2, n3);
297 const SMDS_MeshNode* n12 = GetMediumNode(n1,n2,false);
298 const SMDS_MeshNode* n23 = GetMediumNode(n2,n3,false);
299 const SMDS_MeshNode* n31 = GetMediumNode(n3,n1,false);
301 return meshDS->AddFace(n1, n2, n3, n12, n23, n31);
305 //=======================================================================
308 //=======================================================================
310 * Special function for creation quadratic quadrangle
312 SMDS_MeshFace* StdMeshers_Helper::AddFace(const SMDS_MeshNode* n1,
313 const SMDS_MeshNode* n2,
314 const SMDS_MeshNode* n3,
315 const SMDS_MeshNode* n4)
317 SMESHDS_Mesh * meshDS = GetMesh()->GetMeshDS();
318 if(!myCreateQuadratic) {
319 return meshDS->AddFace(n1, n2, n3, n4);
322 const SMDS_MeshNode* n12 = GetMediumNode(n1,n2,false);
323 const SMDS_MeshNode* n23 = GetMediumNode(n2,n3,false);
324 const SMDS_MeshNode* n34 = GetMediumNode(n3,n4,false);
325 const SMDS_MeshNode* n41 = GetMediumNode(n4,n1,false);
327 return meshDS->AddFace(n1, n2, n3, n4, n12, n23, n34, n41);
331 //=======================================================================
332 //function : AddVolume
334 //=======================================================================
336 * Special function for creation quadratic volume
338 SMDS_MeshVolume* StdMeshers_Helper::AddVolume(const SMDS_MeshNode* n1,
339 const SMDS_MeshNode* n2,
340 const SMDS_MeshNode* n3,
341 const SMDS_MeshNode* n4,
342 const SMDS_MeshNode* n5,
343 const SMDS_MeshNode* n6)
345 SMESHDS_Mesh * meshDS = GetMesh()->GetMeshDS();
346 if(!myCreateQuadratic) {
347 return meshDS->AddVolume(n1, n2, n3, n4, n5, n6);
350 const SMDS_MeshNode* n12 = GetMediumNode(n1,n2,true);
351 const SMDS_MeshNode* n23 = GetMediumNode(n2,n3,true);
352 const SMDS_MeshNode* n31 = GetMediumNode(n3,n1,true);
354 const SMDS_MeshNode* n45 = GetMediumNode(n4,n5,true);
355 const SMDS_MeshNode* n56 = GetMediumNode(n5,n6,true);
356 const SMDS_MeshNode* n64 = GetMediumNode(n6,n4,true);
358 const SMDS_MeshNode* n14 = GetMediumNode(n1,n4,true);
359 const SMDS_MeshNode* n25 = GetMediumNode(n2,n5,true);
360 const SMDS_MeshNode* n36 = GetMediumNode(n3,n6,true);
362 return meshDS->AddVolume(n1, n2, n3, n4, n5, n6,
363 n12, n23, n31, n45, n56, n64, n14, n25, n36);
367 //=======================================================================
368 //function : AddVolume
370 //=======================================================================
372 * Special function for creation quadratic volume
374 SMDS_MeshVolume* StdMeshers_Helper::AddVolume(const SMDS_MeshNode* n1,
375 const SMDS_MeshNode* n2,
376 const SMDS_MeshNode* n3,
377 const SMDS_MeshNode* n4)
379 SMESHDS_Mesh * meshDS = GetMesh()->GetMeshDS();
380 if(!myCreateQuadratic) {
381 return meshDS->AddVolume(n1, n2, n3, n4);
384 const SMDS_MeshNode* n12 = GetMediumNode(n1,n2,true);
385 const SMDS_MeshNode* n23 = GetMediumNode(n2,n3,true);
386 const SMDS_MeshNode* n31 = GetMediumNode(n3,n1,true);
388 const SMDS_MeshNode* n14 = GetMediumNode(n1,n4,true);
389 const SMDS_MeshNode* n24 = GetMediumNode(n2,n4,true);
390 const SMDS_MeshNode* n34 = GetMediumNode(n3,n4,true);
392 return meshDS->AddVolume(n1, n2, n3, n4, n12, n23, n31, n14, n24, n34);
396 //=======================================================================
397 //function : AddVolume
399 //=======================================================================
401 * Special function for creation quadratic volume
403 SMDS_MeshVolume* StdMeshers_Helper::AddVolume(const SMDS_MeshNode* n1,
404 const SMDS_MeshNode* n2,
405 const SMDS_MeshNode* n3,
406 const SMDS_MeshNode* n4,
407 const SMDS_MeshNode* n5,
408 const SMDS_MeshNode* n6,
409 const SMDS_MeshNode* n7,
410 const SMDS_MeshNode* n8)
412 SMESHDS_Mesh * meshDS = GetMesh()->GetMeshDS();
413 if(!myCreateQuadratic) {
414 return meshDS->AddVolume(n1, n2, n3, n4, n5, n6, n7, n8);
417 const SMDS_MeshNode* n12 = GetMediumNode(n1,n2,true);
418 const SMDS_MeshNode* n23 = GetMediumNode(n2,n3,true);
419 const SMDS_MeshNode* n34 = GetMediumNode(n3,n4,true);
420 const SMDS_MeshNode* n41 = GetMediumNode(n4,n1,true);
422 const SMDS_MeshNode* n56 = GetMediumNode(n5,n6,true);
423 const SMDS_MeshNode* n67 = GetMediumNode(n6,n7,true);
424 const SMDS_MeshNode* n78 = GetMediumNode(n7,n8,true);
425 const SMDS_MeshNode* n85 = GetMediumNode(n8,n5,true);
427 const SMDS_MeshNode* n15 = GetMediumNode(n1,n5,true);
428 const SMDS_MeshNode* n26 = GetMediumNode(n2,n6,true);
429 const SMDS_MeshNode* n37 = GetMediumNode(n3,n7,true);
430 const SMDS_MeshNode* n48 = GetMediumNode(n4,n8,true);
432 return meshDS->AddVolume(n1, n2, n3, n4, n5, n6, n7, n8,
433 n12, n23, n34, n41, n56, n67,
434 n78, n85, n15, n26, n37, n48);