X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FStdMeshers%2FStdMeshers_Projection_2D.cxx;h=199a266ed1af7be2abad2a1f12eb7c8ad4e550b2;hb=fb01f8afbf165cba78e3c941d453c1aa29107cab;hp=a31487eefe6b40827458d56e6beb185294fc7b83;hpb=25b744c11fd2621e0f0b2c7a064655dde46e9623;p=modules%2Fsmesh.git diff --git a/src/StdMeshers/StdMeshers_Projection_2D.cxx b/src/StdMeshers/StdMeshers_Projection_2D.cxx index a31487eef..199a266ed 100644 --- a/src/StdMeshers/StdMeshers_Projection_2D.cxx +++ b/src/StdMeshers/StdMeshers_Projection_2D.cxx @@ -1,32 +1,30 @@ -// SMESH SMESH : implementaion of SMESH idl descriptions +// Copyright (C) 2007-2008 CEA/DEN, EDF R&D, OPEN CASCADE +// +// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License. // -// Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, -// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS -// -// This library is free software; you can redistribute it and/or -// modify it under the terms of the GNU Lesser General Public -// License as published by the Free Software Foundation; either -// version 2.1 of the License. -// -// This library is distributed in the hope that it will be useful, -// but WITHOUT ANY WARRANTY; without even the implied warranty of -// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -// Lesser General Public License for more details. -// -// You should have received a copy of the GNU Lesser General Public -// License along with this library; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA -// -// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. // +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA // +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com // +// SMESH SMESH : implementaion of SMESH idl descriptions // File : StdMeshers_Projection_2D.cxx // Module : SMESH // Created : Fri Oct 20 11:37:07 2006 // Author : Edward AGAPOV (eap) - - +// #include "StdMeshers_Projection_2D.hxx" #include "StdMeshers_ProjectionSource2D.hxx" @@ -37,7 +35,7 @@ #include "SMESH_Block.hxx" #include "SMESH_Gen.hxx" #include "SMESH_Mesh.hxx" -#include "SMESH_MeshEditor.hxx" +#include "SMESH_MesherHelper.hxx" #include "SMESH_Pattern.hxx" #include "SMESH_subMesh.hxx" #include "SMESH_subMeshEventListener.hxx" @@ -295,16 +293,14 @@ namespace { // Find a new node connected to nV1 and belonging to edge submesh; const SMDS_MeshNode* nE = 0; SMESHDS_SubMesh* smDS = sm->GetSubMeshDS(); - SMDS_ElemIteratorPtr vElems = nV1->GetInverseElementIterator(); + SMDS_ElemIteratorPtr vElems = nV1->GetInverseElementIterator(SMDSAbs_Face); while ( vElems->more() && !nE ) { const SMDS_MeshElement* elem = vElems->next(); - if ( elem->GetType() != SMDSAbs_Face ) - continue; // new nodes are shared by faces int nbNodes = elem->NbNodes(); if ( elem->IsQuadratic() ) nbNodes /= 2; int iV1 = elem->GetNodeIndex( nV1 ); - // try next aftre nV1 + // try next after nV1 int iE = SMESH_MesherHelper::WrapIndex( iV1 + 1, nbNodes ); if ( smDS->Contains( elem->GetNode( iE ) )) nE = elem->GetNode( iE ); @@ -358,6 +354,100 @@ namespace { } // bool getBoundaryNodes() + //================================================================================ + /*! + * \brief Preform projection in case if tgtFace.IsPartner( srcFace ) + * \param tgtFace - target face + * \param srcFace - source face + * \param tgtMesh - target mesh + * \param srcMesh - source mesh + * \retval bool - true if succeeded + */ + //================================================================================ + + bool projectPartner(const TopoDS_Face& tgtFace, + const TopoDS_Face& srcFace, + SMESH_Mesh * tgtMesh, + SMESH_Mesh * srcMesh, + const TAssocTool::TShapeShapeMap& shape2ShapeMap) + { + if ( !tgtFace.IsPartner( srcFace )) + return false; + + // Fill map of src to tgt nodes with nodes on edges + + map src2tgtNodes; + map::iterator srcN_tgtN; + + for ( TopExp_Explorer srcEdge( srcFace, TopAbs_EDGE); srcEdge.More(); srcEdge.Next() ) + { + const TopoDS_Shape& tgtEdge = shape2ShapeMap( srcEdge.Current() ); + if ( !tgtEdge.IsPartner( srcEdge.Current() )) + return false; + + map< double, const SMDS_MeshNode* > srcNodes, tgtNodes; + if ( !SMESH_Algo::GetSortedNodesOnEdge( srcMesh->GetMeshDS(), + TopoDS::Edge( srcEdge.Current() ), + /*ignoreMediumNodes = */true, + srcNodes ) + || + !SMESH_Algo::GetSortedNodesOnEdge( tgtMesh->GetMeshDS(), + TopoDS::Edge( tgtEdge ), + /*ignoreMediumNodes = */true, + tgtNodes ) + || + srcNodes.size() != tgtNodes.size()) + return false; + + map< double, const SMDS_MeshNode* >::iterator u_tn = tgtNodes.begin(); + map< double, const SMDS_MeshNode* >::iterator u_sn = srcNodes.begin(); + for ( ; u_tn != tgtNodes.end(); ++u_tn, ++u_sn) + src2tgtNodes.insert( make_pair( u_sn->second, u_tn->second )); + } + + // Make new faces + + // transformation to get location of target nodes from source ones + gp_Trsf srcTrsf = srcFace.Location(); + gp_Trsf tgtTrsf = tgtFace.Location(); + gp_Trsf trsf = srcTrsf.Inverted() * tgtTrsf; + + // prepare the helper adding quadratic elements if necessary + SMESH_MesherHelper helper( *tgtMesh ); + helper.IsQuadraticSubMesh( tgtFace ); + helper.SetElementsOnShape( true ); + + const SMDS_MeshNode* nullNode = 0; + + SMESHDS_SubMesh* srcSubDS = srcMesh->GetMeshDS()->MeshElements( srcFace ); + SMDS_ElemIteratorPtr elemIt = srcSubDS->GetElements(); + while ( elemIt->more() ) // loop on all mesh faces on srcFace + { + const SMDS_MeshElement* elem = elemIt->next(); + vector< const SMDS_MeshNode* > tgtFaceNodes; + tgtFaceNodes.reserve( elem->NbNodes() ); + SMDS_ElemIteratorPtr nodeIt = elem->nodesIterator(); + while ( nodeIt->more() ) // loop on nodes of the source element + { + const SMDS_MeshNode* srcNode = (const SMDS_MeshNode*) nodeIt->next(); + srcN_tgtN = src2tgtNodes.insert( make_pair( srcNode, nullNode )).first; + if ( srcN_tgtN->second == nullNode ) + { + // create a new node + gp_Pnt tgtP = gp_Pnt(srcNode->X(),srcNode->Y(),srcNode->Z()).Transformed( trsf ); + srcN_tgtN->second = helper.AddNode( tgtP.X(), tgtP.Y(), tgtP.Z() ); + } + tgtFaceNodes.push_back( srcN_tgtN->second ); + } + // create a new face (with reversed orientation) + if ( tgtFaceNodes.size() == 3 ) + helper.AddFace( tgtFaceNodes[0],tgtFaceNodes[2],tgtFaceNodes[1]); + else + helper.AddFace( tgtFaceNodes[0],tgtFaceNodes[3],tgtFaceNodes[2],tgtFaceNodes[1]); + } + return true; + } + } // namespace //======================================================================= @@ -409,6 +499,10 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& return error(COMPERR_BAD_INPUT_MESH,"Source mesh not computed"); } + // try to project from same face with different location + if ( projectPartner( tgtFace, srcFace, tgtMesh, srcMesh, shape2ShapeMap )) + return true; + // -------------------- // Prepare to mapping // -------------------- @@ -419,18 +513,20 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& // Check if node projection to a face is needed Bnd_B2d uvBox; SMDS_ElemIteratorPtr faceIt = srcSubMesh->GetSubMeshDS()->GetElements(); - for ( int nbN = 0; nbN < 3 && faceIt->more(); ) { + int nbFaceNodes = 0; + for ( ; nbFaceNodes < 3 && faceIt->more(); ) { const SMDS_MeshElement* face = faceIt->next(); SMDS_ElemIteratorPtr nodeIt = face->nodesIterator(); while ( nodeIt->more() ) { const SMDS_MeshNode* node = static_cast( nodeIt->next() ); if ( node->GetPosition()->GetTypeOfPosition() == SMDS_TOP_FACE ) { - nbN++; + nbFaceNodes++; uvBox.Add( helper.GetNodeUV( srcFace, node )); } } } - const bool toProjectNodes = ( uvBox.IsVoid() || uvBox.SquareExtent() < DBL_MIN ); + const bool toProjectNodes = + ( nbFaceNodes > 0 && ( uvBox.IsVoid() || uvBox.SquareExtent() < DBL_MIN )); // Load pattern from the source face SMESH_Pattern mapper; @@ -515,7 +611,7 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& // Sort new and old nodes of a submesh separately - bool isSeam = helper.IsSeamShape( sm->GetId() ); + bool isSeam = helper.IsRealSeam( sm->GetId() ); enum { NEW_NODES = 0, OLD_NODES }; map< double, const SMDS_MeshNode* > u2nodesMaps[2], u2nodesOnSeam; @@ -537,7 +633,7 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& continue; // node is already in the map } - // sort nodes on edges by its position + // sort nodes on edges by their position map< double, const SMDS_MeshNode* > & pos2nodes = u2nodesMaps[isOld ? OLD_NODES : NEW_NODES]; switch ( node->GetPosition()->GetTypeOfPosition() ) { @@ -558,12 +654,13 @@ bool StdMeshers_Projection_2D::Compute(SMESH_Mesh& theMesh, const TopoDS_Shape& } if ( u2nodesMaps[ NEW_NODES ].size() != u2nodesMaps[ OLD_NODES ].size() ) { - if ( u2nodesMaps[ NEW_NODES ].size() == 0 && - sm->GetSubShape().ShapeType() == TopAbs_EDGE && - BRep_Tool::Degenerated( TopoDS::Edge( sm->GetSubShape() ))) + if ( u2nodesMaps[ NEW_NODES ].size() == 0 && + sm->GetSubShape().ShapeType() == TopAbs_EDGE && + helper.IsDegenShape( sm->GetId() ) ) // NPAL15894 (tt88bis.py) - project mesh built by NETGEN_1d_2D that - // does not make segments/nodes on degenerated edges + // does not make segments/nodes on degenerated edges continue; + RETURN_BAD_RESULT("Different nb of old and new nodes on shape #"<< sm->GetId() <<" "<< u2nodesMaps[ OLD_NODES ].size() << " != " << u2nodesMaps[ NEW_NODES ].size());