X-Git-Url: http://git.salome-platform.org/gitweb/?p=modules%2Fsmesh.git;a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_Quadrangle_2D.cxx;h=da85567a93dbb69534f8036d213ddeb61303c06f;hp=016a5098db2c27e8f35e096d8077e519502b4483;hb=d9faba6c847c1c1a4d4f501ca5ac5725a25a8236;hpb=209c9bd9c33255f46bfa60acc76c4c3878f815aa diff --git a/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx b/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx index 016a5098d..da85567a9 100644 --- a/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx +++ b/src/StdMeshers/StdMeshers_Quadrangle_2D.cxx @@ -67,9 +67,7 @@ #ifndef StdMeshers_Array2OfNode_HeaderFile #define StdMeshers_Array2OfNode_HeaderFile typedef const SMDS_MeshNode* SMDS_MeshNodePtr; -DEFINE_BASECOLLECTION (StdMeshers_BaseCollectionNodePtr, SMDS_MeshNodePtr) -DEFINE_ARRAY2(StdMeshers_Array2OfNode, - StdMeshers_BaseCollectionNodePtr, SMDS_MeshNodePtr) +typedef NCollection_Array2 StdMeshers_Array2OfNode; #endif using namespace std; @@ -339,7 +337,7 @@ bool StdMeshers_Quadrangle_2D::computeTriangles(SMESH_Mesh& aMesh, FaceQuadStruct::Ptr newQuad = myQuadList.back(); if ( quad != newQuad ) // split done { - { + { // update left side limit till where to make triangles FaceQuadStruct::Ptr botQuad = // a bottom part ( quad->side[ QUAD_LEFT_SIDE ].from == 0 ) ? quad : newQuad; if ( botQuad->nbNodeOut( QUAD_LEFT_SIDE ) > 0 ) @@ -359,12 +357,33 @@ bool StdMeshers_Quadrangle_2D::computeTriangles(SMESH_Mesh& aMesh, { splitQuad( quad, quad->iSize-2, 0 ); } - if ( quad->nbNodeOut( QUAD_LEFT_SIDE )) + if ( quad->nbNodeOut( QUAD_LEFT_SIDE )) { splitQuad( quad, 1, 0 ); + + if ( quad->nbNodeOut( QUAD_TOP_SIDE )) + { + newQuad = myQuadList.back(); + if ( newQuad == quad ) // too narrow to split + { + // update left side limit till where to make triangles + quad->side[ QUAD_LEFT_SIDE ].to--; + } + else + { + FaceQuadStruct::Ptr leftQuad = + ( quad->side[ QUAD_BOTTOM_SIDE ].from == 0 ) ? quad : newQuad; + leftQuad->nbNodeOut( QUAD_TOP_SIDE ) = 0; + } + } } - return computeQuadDominant( aMesh, aFace ); + if ( ! computeQuadDominant( aMesh, aFace )) + return false; + + // try to fix zero-area triangles near straight-angle corners + + return true; } //================================================================================ @@ -591,17 +610,59 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh& aMesh, iup = nbhoriz - 1; int stop = 0; - // if left edge is out, we will stop at a second node - //if (quad->nbNodeOut(3)) stop++; - if ( quad->nbNodeOut( QUAD_RIGHT_SIDE )) - quad->UVPt( nbhoriz-1, 0 ).node = uv_e1[ nbright-2 ].node; - if ( quad->nbNodeOut( QUAD_LEFT_SIDE )) - quad->UVPt( 0, 0 ).node = uv_e3[ nbleft-2 ].node; + if ( quad->side[3].grid->Edge(0).IsNull() ) // left side is simulated one + { + // quad divided at I but not at J, as nbvertic==nbright==2 + stop++; // we stop at a second node + } + else + { + if ( quad->nbNodeOut( QUAD_RIGHT_SIDE )) + quad->UVPt( nbhoriz-1, 0 ).node = uv_e1[ nbright-2 ].node; + if ( quad->nbNodeOut( QUAD_LEFT_SIDE )) + quad->UVPt( 0, 0 ).node = uv_e3[ nbleft-2 ].node; + if ( nbright > 2 ) // there was a split at J + quad->nbNodeOut( QUAD_LEFT_SIDE ) = 0; + } + const SMDS_MeshNode *a, *b, *c, *d; + i = nbup - 1; + // avoid creating zero-area triangles near a straight-angle corner + { + a = uv_e2[i].node; + b = uv_e2[i-1].node; + c = uv_e1[nbright-2].node; + SMESH_TNodeXYZ pa( a ), pb( b ), pc( c ); + double area = 0.5 * (( pb - pa ) ^ ( pc - pa )).Modulus(); + if ( Abs( area ) < 1e-20 ) + { + --g; + d = quad->UVPt( g, nbvertic-2 ).node; + if ( myTrianglePreference ) + { + if ( SMDS_MeshFace* face = myHelper->AddFace(a, d, c)) + meshDS->SetMeshElementOnShape(face, geomFaceID); + } + else + { + if ( SMDS_MeshFace* face = myHelper->AddFace(a, b, d, c)) + { + meshDS->SetMeshElementOnShape(face, geomFaceID); + SMESH_ComputeErrorPtr& err = aMesh.GetSubMesh( aFace )->GetComputeError(); + if ( !err || err->IsOK() || err->myName < COMPERR_WARNING ) + { + err.reset( new SMESH_ComputeError( COMPERR_WARNING, + "Bad quality quad created")); + err->myBadElements.push_back( face ); + } + } + --i; + } + } + } // for each node of the up edge find nearest node // in the first row of the regular grid and link them - for (i = nbup - 1; i > stop; i--) { - const SMDS_MeshNode *a, *b, *c, *d; + for ( ; i > stop; i--) { a = uv_e2[i].node; b = uv_e2[i - 1].node; gp_Pnt pb (b->X(), b->Y(), b->Z()); @@ -745,8 +806,43 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh& aMesh, int g = nbvertic - 1; // last processed node in the grid int stop = 0; i = quad->side[ QUAD_LEFT_SIDE ].to-1; // nbleft - 1; - for (; i > stop; i--) { - const SMDS_MeshNode *a, *b, *c, *d; + + const SMDS_MeshNode *a, *b, *c, *d; + // avoid creating zero-area triangles near a straight-angle corner + { + a = uv_e3[i].node; + b = uv_e3[i-1].node; + c = quad->UVPt( 1, g ).node; + SMESH_TNodeXYZ pa( a ), pb( b ), pc( c ); + double area = 0.5 * (( pb - pa ) ^ ( pc - pa )).Modulus(); + if ( Abs( area ) < 1e-20 ) + { + --g; + d = quad->UVPt( 1, g ).node; + if ( myTrianglePreference ) + { + if ( SMDS_MeshFace* face = myHelper->AddFace(a, d, c)) + meshDS->SetMeshElementOnShape(face, geomFaceID); + } + else + { + if ( SMDS_MeshFace* face = myHelper->AddFace(a, b, d, c)) + { + meshDS->SetMeshElementOnShape(face, geomFaceID); + SMESH_ComputeErrorPtr& err = aMesh.GetSubMesh( aFace )->GetComputeError(); + if ( !err || err->IsOK() || err->myName < COMPERR_WARNING ) + { + err.reset( new SMESH_ComputeError( COMPERR_WARNING, + "Bad quality quad created")); + err->myBadElements.push_back( face ); + } + } + --i; + } + } + } + for (; i > stop; i--) // loop on nodes on the left side + { a = uv_e3[i].node; b = uv_e3[i - 1].node; gp_Pnt pb (b->X(), b->Y(), b->Z()); @@ -756,12 +852,13 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh& aMesh, if (i == stop + 1) { // down bondary reached c = quad->uv_grid[nbhoriz*jlow + 1].node; near = jlow; - } else { + } + else { double mind = RealLast(); for (int k = g; k >= jlow; k--) { const SMDS_MeshNode *nk; if (k > jup) - nk = uv_e2[1].node; + nk = quad->uv_grid[nbhoriz*jup + 1].node; //uv_e2[1].node; else nk = quad->uv_grid[nbhoriz*k + 1].node; gp_Pnt pnk (nk->X(), nk->Y(), nk->Z()); @@ -782,11 +879,10 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh& aMesh, } else { // make quadrangle if (near + 1 > jup) - d = uv_e2[1].node; + d = quad->uv_grid[nbhoriz*jup + 1].node; //uv_e2[1].node; else d = quad->uv_grid[nbhoriz*(near + 1) + 1].node; - //SMDS_MeshFace* face = meshDS->AddFace(a, b, c, d); - if (!myTrianglePreference){ + if (!myTrianglePreference) { SMDS_MeshFace* face = myHelper->AddFace(a, b, c, d); if (face) meshDS->SetMeshElementOnShape(face, geomFaceID); } @@ -798,7 +894,7 @@ bool StdMeshers_Quadrangle_2D::computeQuadDominant(SMESH_Mesh& aMesh, for (int k = near + 1; k < g; k++) { c = quad->uv_grid[nbhoriz*k + 1].node; if (k + 1 > jup) - d = uv_e2[1].node; + d = quad->uv_grid[nbhoriz*jup + 1].node; //uv_e2[1].node; else d = quad->uv_grid[nbhoriz*(k + 1) + 1].node; SMDS_MeshFace* face = myHelper->AddFace(a, c, d); @@ -1036,45 +1132,40 @@ FaceQuadStruct::Ptr StdMeshers_Quadrangle_2D::CheckNbEdges(SMESH_Mesh & { list< TopoDS_Edge > sideEdges; TopoDS_Vertex nextSideV = corners[( iSide + 1 - nbUsedDegen ) % corners.size() ]; - while ( edgeIt != edges.end() && - !nextSideV.IsSame( myHelper->IthVertex( 0, *edgeIt ))) + bool nextSideVReached = false; + do { - if ( SMESH_Algo::isDegenerated( *edgeIt ) ) + const TopoDS_Edge& edge = *edgeIt; + nextSideVReached = nextSideV.IsSame( myHelper->IthVertex( 1, edge )); + if ( SMESH_Algo::isDegenerated( edge )) { - if ( myNeedSmooth ) - { - ++edgeIt; // no side on the degenerated EDGE - } - else + if ( !myNeedSmooth ) // need to make a side on a degen edge { if ( sideEdges.empty() ) { + sideEdges.push_back( edge ); ++nbUsedDegen; - sideEdges.push_back( *edgeIt++ ); // a degenerated side - break; + nextSideVReached = true; } else { - break; // do not append a degenerated EDGE to a regular side + break; } } } else { - sideEdges.push_back( *edgeIt++ ); + sideEdges.push_back( edge ); } + ++edgeIt; } + while ( edgeIt != edges.end() && !nextSideVReached ); + if ( !sideEdges.empty() ) { - quad->side.push_back( StdMeshers_FaceSide::New( F, sideEdges, &aMesh, iSide < QUAD_TOP_SIDE, - ignoreMediumNodes, myProxyMesh )); - ++iSide; - } - else if ( !SMESH_Algo::isDegenerated( *edgeIt ) && // closed EDGE - myHelper->IthVertex( 0, *edgeIt ).IsSame( myHelper->IthVertex( 1, *edgeIt ))) - { - quad->side.push_back( StdMeshers_FaceSide::New( F, *edgeIt++, &aMesh, iSide < QUAD_TOP_SIDE, - ignoreMediumNodes, myProxyMesh)); + quad->side.push_back + ( StdMeshers_FaceSide::New( F, sideEdges, &aMesh, iSide < QUAD_TOP_SIDE, + ignoreMediumNodes, myProxyMesh )); ++iSide; } if ( quad->side.size() == 4 ) @@ -3750,7 +3841,7 @@ void StdMeshers_Quadrangle_2D::updateDegenUV(FaceQuadStruct::Ptr quad) uv1.v = uv2.v = 0.5 * ( uv1.v + uv2.v ); } - else if ( quad->side.size() == 4 && myQuadType == QUAD_STANDARD) + else if ( quad->side.size() == 4 /*&& myQuadType == QUAD_STANDARD*/) // Set number of nodes on a degenerated side to be same as on an opposite side // ---------------------------------------------------------------------------- @@ -3977,8 +4068,9 @@ bool StdMeshers_Quadrangle_2D::check() int iPrev = myHelper->WrapIndex( i-1, wire->NbEdges() ); const TopoDS_Edge& e1 = wire->Edge( iPrev ); const TopoDS_Edge& e2 = wire->Edge( i ); - double angle = myHelper->GetAngle( e1, e2, geomFace ); - if ( maxAngle < angle && angle < 0.9 * M_PI ) + double angle = myHelper->GetAngle( e1, e2, geomFace, wire->FirstVertex( i )); + if (( maxAngle < angle ) && + ( 5.* M_PI/180 < angle && angle < 175.* M_PI/180 )) { maxAngle = angle; iVertex = i; @@ -4130,7 +4222,7 @@ int StdMeshers_Quadrangle_2D::getCorners(const TopoDS_Face& theFace, TopoDS_Vertex v = helper.IthVertex( 0, *edge ); if ( !theConsiderMesh || SMESH_Algo::VertexNode( v, helper.GetMeshDS() )) { - double angle = SMESH_MesherHelper::GetAngle( prevE, *edge, theFace ); + double angle = SMESH_MesherHelper::GetAngle( prevE, *edge, theFace, v ); vertexByAngle.insert( make_pair( angle, v )); angleByVertex.Bind( v, angle ); } @@ -4182,7 +4274,7 @@ int StdMeshers_Quadrangle_2D::getCorners(const TopoDS_Face& theFace, if ( nbCorners == 3 ) vMap.Add( triaVertex ); multimap::reverse_iterator a2v = vertexByAngle.rbegin(); - for ( ; a2v != vertexByAngle.rend() && vMap.Extent() < nbCorners; ++a2v ) + for ( int iC = 0; a2v != vertexByAngle.rend() && iC < nbCorners; ++a2v, ++iC ) vMap.Add( (*a2v).second ); // check if there are possible variations in choosing corners @@ -4707,7 +4799,7 @@ int StdMeshers_Quadrangle_2D::splitQuad(FaceQuadStruct::Ptr quad, int I, int J) myQuadList.push_back( FaceQuadStruct::Ptr( newQuad )); vector points; - if ( I > 0 ) + if ( I > 0 && I <= quad->iSize-2 ) { points.reserve( quad->jSize ); for ( int jP = 0; jP < quad->jSize; ++jP ) @@ -4746,7 +4838,7 @@ int StdMeshers_Quadrangle_2D::splitQuad(FaceQuadStruct::Ptr quad, int I, int J) return QUAD_LEFT_SIDE; } - else if ( J > 0 ) //// split horizontally, a new quad is below an old one + else if ( J > 0 && J <= quad->jSize-2 ) //// split horizontally, a new quad is below an old one { points.reserve( quad->iSize ); for ( int iP = 0; iP < quad->iSize; ++iP )