+ const bool allHex = ( theHexas.size() == myMesh->NbHexas() );
+
+ while ( startHex )
+ {
+ // move in two directions from startHex via facetID
+ for ( int is2nd = 0; is2nd < 2; ++is2nd )
+ {
+ curHex = startHex;
+ int curFacet = facetID;
+ if ( is2nd ) // do not treat startHex twice
+ {
+ vTool.Set( curHex );
+ if ( vTool.IsFreeFace( curFacet, &curHex ))
+ {
+ curHex = 0;
+ }
+ else
+ {
+ vTool.GetFaceNodes( curFacet, facetNodes );
+ vTool.Set( curHex );
+ curFacet = vTool.GetFaceIndex( facetNodes );
+ }
+ }
+ while ( curHex )
+ {
+ // store a facet to split
+ if ( curHex->GetGeomType() != SMDSGeom_HEXA )
+ {
+ theFacets.insert( make_pair( curHex, -1 ));
+ break;
+ }
+ if ( !allHex && !theHexas.count( curHex ))
+ break;
+
+ pair< TFacetOfElem::iterator, bool > facetIt2isNew =
+ theFacets.insert( make_pair( curHex, curFacet ));
+ if ( !facetIt2isNew.second )
+ break;
+
+ // remember not-to-split facets in facetsToCheck
+ int oppFacet = vTool.GetOppFaceIndexOfHex( curFacet );
+ for ( int iF = 0; iF < vTool.NbFaces(); ++iF )
+ {
+ if ( iF == curFacet && iF == oppFacet )
+ continue;
+ TVolumeFaceKey facetKey ( vTool, iF );
+ TElemFacets elemFacet( facetIt2isNew.first, iF );
+ pair< TFacetMap::iterator, bool > it2isnew =
+ facetsToCheck.insert( make_pair( facetKey, elemFacet ));
+ if ( !it2isnew.second )
+ facetsToCheck.erase( it2isnew.first ); // adjacent hex already checked
+ }
+ // pass to a volume adjacent via oppFacet
+ if ( vTool.IsFreeFace( oppFacet, &curHex ))
+ {
+ curHex = 0;
+ }
+ else
+ {
+ // get a new curFacet
+ vTool.GetFaceNodes( oppFacet, facetNodes );
+ vTool.Set( curHex );
+ curFacet = vTool.GetFaceIndex( facetNodes, /*hint=*/curFacet );
+ }
+ }
+ } // move in two directions from startHex via facetID
+
+ // Find a new startHex by facetsToCheck
+
+ startHex = 0;
+ facetID = -1;
+ TFacetMap::iterator fIt = facetsToCheck.begin();
+ while ( !startHex && fIt != facetsToCheck.end() )
+ {
+ const TElemFacets& elemFacets = fIt->second;
+ const SMDS_MeshElement* hex = elemFacets.first->first;
+ int splitFacet = elemFacets.first->second;
+ int lateralFacet = elemFacets.second;
+ facetsToCheck.erase( fIt );
+ fIt = facetsToCheck.begin();
+
+ vTool.Set( hex );
+ if ( vTool.IsFreeFace( lateralFacet, &curHex ) ||
+ curHex->GetGeomType() != SMDSGeom_HEXA )
+ continue;
+ if ( !allHex && !theHexas.count( curHex ))
+ continue;
+
+ startHex = curHex;
+
+ // find a facet of startHex to split
+
+ set<const SMDS_MeshNode*> lateralNodes;
+ vTool.GetFaceNodes( lateralFacet, lateralNodes );
+ vTool.GetFaceNodes( splitFacet, facetNodes );
+ int oppLateralFacet = vTool.GetOppFaceIndexOfHex( lateralFacet );
+ vTool.Set( startHex );
+ lateralFacet = vTool.GetFaceIndex( lateralNodes, oppLateralFacet );
+
+ // look for a facet of startHex having common nodes with facetNodes
+ // but not lateralFacet
+ for ( int iF = 0; iF < vTool.NbFaces(); ++iF )
+ {
+ if ( iF == lateralFacet )
+ continue;
+ int nbCommonNodes = 0;
+ const SMDS_MeshNode** nn = vTool.GetFaceNodes( iF );
+ for ( int iN = 0, nbN = vTool.NbFaceNodes( iF ); iN < nbN; ++iN )
+ nbCommonNodes += facetNodes.count( nn[ iN ]);
+
+ if ( nbCommonNodes >= 2 )
+ {
+ facetID = iF;
+ break;
+ }
+ }
+ if ( facetID < 0 )
+ throw SALOME_Exception( THIS_METHOD "facet of a new startHex not found");
+ }
+ } // while ( startHex )
+
+ return;
+}
+
+namespace
+{
+ //================================================================================
+ /*!
+ * \brief Selects nodes of several elements according to a given interlace
+ * \param [in] srcNodes - nodes to select from
+ * \param [out] tgtNodesVec - array of nodes of several elements to fill in
+ * \param [in] interlace - indices of nodes for all elements
+ * \param [in] nbElems - nb of elements
+ * \param [in] nbNodes - nb of nodes in each element
+ * \param [in] mesh - the mesh
+ * \param [out] elemQueue - a list to push elements found by the selected nodes
+ * \param [in] type - type of elements to look for
+ */
+ //================================================================================
+
+ void selectNodes( const vector< const SMDS_MeshNode* >& srcNodes,
+ vector< const SMDS_MeshNode* >* tgtNodesVec,
+ const int* interlace,
+ const int nbElems,
+ const int nbNodes,
+ SMESHDS_Mesh* mesh = 0,
+ list< const SMDS_MeshElement* >* elemQueue=0,
+ SMDSAbs_ElementType type=SMDSAbs_All)
+ {
+ for ( int iE = 0; iE < nbElems; ++iE )
+ {
+ vector< const SMDS_MeshNode* >& elemNodes = tgtNodesVec[iE];
+ const int* select = & interlace[iE*nbNodes];
+ elemNodes.resize( nbNodes );
+ for ( int iN = 0; iN < nbNodes; ++iN )
+ elemNodes[iN] = srcNodes[ select[ iN ]];
+ }
+ const SMDS_MeshElement* e;
+ if ( elemQueue )
+ for ( int iE = 0; iE < nbElems; ++iE )
+ if (( e = mesh->FindElement( tgtNodesVec[iE], type, /*noMedium=*/false)))
+ elemQueue->push_back( e );
+ }
+}
+
+//=======================================================================
+/*
+ * Split bi-quadratic elements into linear ones without creation of additional nodes
+ * - bi-quadratic triangle will be split into 3 linear quadrangles;
+ * - bi-quadratic quadrangle will be split into 4 linear quadrangles;
+ * - tri-quadratic hexahedron will be split into 8 linear hexahedra;
+ * Quadratic elements of lower dimension adjacent to the split bi-quadratic element
+ * will be split in order to keep the mesh conformal.
+ * \param elems - elements to split
+ */
+//=======================================================================
+
+void SMESH_MeshEditor::SplitBiQuadraticIntoLinear(TIDSortedElemSet& theElems)
+{
+ vector< const SMDS_MeshNode* > elemNodes(27), subNodes[12], splitNodes[8];
+ vector<const SMDS_MeshElement* > splitElems;
+ list< const SMDS_MeshElement* > elemQueue;
+ list< const SMDS_MeshElement* >::iterator elemIt;
+
+ SMESHDS_Mesh * mesh = GetMeshDS();
+ ElemFeatures *elemType, hexaType(SMDSAbs_Volume), quadType(SMDSAbs_Face), segType(SMDSAbs_Edge);
+ int nbElems, nbNodes;
+
+ TIDSortedElemSet::iterator elemSetIt = theElems.begin();
+ for ( ; elemSetIt != theElems.end(); ++elemSetIt )
+ {
+ elemQueue.clear();
+ elemQueue.push_back( *elemSetIt );
+ for ( elemIt = elemQueue.begin(); elemIt != elemQueue.end(); ++elemIt )
+ {
+ const SMDS_MeshElement* elem = *elemIt;
+ switch( elem->GetEntityType() )
+ {
+ case SMDSEntity_TriQuad_Hexa: // HEX27
+ {
+ elemNodes.assign( elem->begin_nodes(), elem->end_nodes() );
+ nbElems = nbNodes = 8;
+ elemType = & hexaType;
+
+ // get nodes for new elements
+ static int vInd[8][8] = {{ 0,8,20,11, 16,21,26,24 },
+ { 1,9,20,8, 17,22,26,21 },
+ { 2,10,20,9, 18,23,26,22 },
+ { 3,11,20,10, 19,24,26,23 },
+ { 16,21,26,24, 4,12,25,15 },
+ { 17,22,26,21, 5,13,25,12 },
+ { 18,23,26,22, 6,14,25,13 },
+ { 19,24,26,23, 7,15,25,14 }};
+ selectNodes( elemNodes, & splitNodes[0], &vInd[0][0], nbElems, nbNodes );
+
+ // add boundary faces to elemQueue
+ static int fInd[6][9] = {{ 0,1,2,3, 8,9,10,11, 20 },
+ { 4,5,6,7, 12,13,14,15, 25 },
+ { 0,1,5,4, 8,17,12,16, 21 },
+ { 1,2,6,5, 9,18,13,17, 22 },
+ { 2,3,7,6, 10,19,14,18, 23 },
+ { 3,0,4,7, 11,16,15,19, 24 }};
+ selectNodes( elemNodes, & subNodes[0], &fInd[0][0], 6,9, mesh, &elemQueue, SMDSAbs_Face );
+
+ // add boundary segments to elemQueue
+ static int eInd[12][3] = {{ 0,1,8 }, { 1,2,9 }, { 2,3,10 }, { 3,0,11 },
+ { 4,5,12}, { 5,6,13}, { 6,7,14 }, { 7,4,15 },
+ { 0,4,16}, { 1,5,17}, { 2,6,18 }, { 3,7,19 }};
+ selectNodes( elemNodes, & subNodes[0], &eInd[0][0], 12,3, mesh, &elemQueue, SMDSAbs_Edge );
+ break;
+ }
+ case SMDSEntity_BiQuad_Triangle: // TRIA7
+ {
+ elemNodes.assign( elem->begin_nodes(), elem->end_nodes() );
+ nbElems = 3;
+ nbNodes = 4;
+ elemType = & quadType;
+
+ // get nodes for new elements
+ static int fInd[3][4] = {{ 0,3,6,5 }, { 1,4,6,3 }, { 2,5,6,4 }};
+ selectNodes( elemNodes, & splitNodes[0], &fInd[0][0], nbElems, nbNodes );
+
+ // add boundary segments to elemQueue
+ static int eInd[3][3] = {{ 0,1,3 }, { 1,2,4 }, { 2,0,5 }};
+ selectNodes( elemNodes, & subNodes[0], &eInd[0][0], 3,3, mesh, &elemQueue, SMDSAbs_Edge );
+ break;
+ }
+ case SMDSEntity_BiQuad_Quadrangle: // QUAD9
+ {
+ elemNodes.assign( elem->begin_nodes(), elem->end_nodes() );
+ nbElems = 4;
+ nbNodes = 4;
+ elemType = & quadType;
+
+ // get nodes for new elements
+ static int fInd[4][4] = {{ 0,4,8,7 }, { 1,5,8,4 }, { 2,6,8,5 }, { 3,7,8,6 }};
+ selectNodes( elemNodes, & splitNodes[0], &fInd[0][0], nbElems, nbNodes );
+
+ // add boundary segments to elemQueue
+ static int eInd[4][3] = {{ 0,1,4 }, { 1,2,5 }, { 2,3,6 }, { 3,0,7 }};
+ selectNodes( elemNodes, & subNodes[0], &eInd[0][0], 4,3, mesh, &elemQueue, SMDSAbs_Edge );
+ break;
+ }
+ case SMDSEntity_Quad_Edge:
+ {
+ if ( elemIt == elemQueue.begin() )
+ continue; // an elem is in theElems
+ elemNodes.assign( elem->begin_nodes(), elem->end_nodes() );
+ nbElems = 2;
+ nbNodes = 2;
+ elemType = & segType;
+
+ // get nodes for new elements
+ static int eInd[2][2] = {{ 0,2 }, { 2,1 }};
+ selectNodes( elemNodes, & splitNodes[0], &eInd[0][0], nbElems, nbNodes );
+ break;
+ }
+ default: continue;
+ } // switch( elem->GetEntityType() )
+
+ // Create new elements
+
+ SMESHDS_SubMesh* subMesh = mesh->MeshElements( elem->getshapeId() );
+
+ splitElems.clear();
+
+ //elemType->SetID( elem->GetID() ); // create an elem with the same ID as a removed one
+ mesh->RemoveFreeElement( elem, subMesh, /*fromGroups=*/false );
+ //splitElems.push_back( AddElement( splitNodes[ 0 ], *elemType ));
+ //elemType->SetID( -1 );
+
+ for ( int iE = 0; iE < nbElems; ++iE )
+ splitElems.push_back( AddElement( splitNodes[ iE ], *elemType ));
+
+
+ ReplaceElemInGroups( elem, splitElems, mesh );
+
+ if ( subMesh )
+ for ( size_t i = 0; i < splitElems.size(); ++i )
+ subMesh->AddElement( splitElems[i] );
+ }
+ }
+}
+
+//=======================================================================
+//function : AddToSameGroups
+//purpose : add elemToAdd to the groups the elemInGroups belongs to
+//=======================================================================
+
+void SMESH_MeshEditor::AddToSameGroups (const SMDS_MeshElement* elemToAdd,
+ const SMDS_MeshElement* elemInGroups,
+ SMESHDS_Mesh * aMesh)
+{
+ const set<SMESHDS_GroupBase*>& groups = aMesh->GetGroups();
+ if (!groups.empty()) {
+ set<SMESHDS_GroupBase*>::const_iterator grIt = groups.begin();
+ for ( ; grIt != groups.end(); grIt++ ) {
+ SMESHDS_Group* group = dynamic_cast<SMESHDS_Group*>( *grIt );
+ if ( group && group->Contains( elemInGroups ))
+ group->SMDSGroup().Add( elemToAdd );
+ }
+ }
+}
+
+
+//=======================================================================
+//function : RemoveElemFromGroups
+//purpose : Remove removeelem to the groups the elemInGroups belongs to
+//=======================================================================
+void SMESH_MeshEditor::RemoveElemFromGroups (const SMDS_MeshElement* removeelem,
+ SMESHDS_Mesh * aMesh)
+{
+ const set<SMESHDS_GroupBase*>& groups = aMesh->GetGroups();
+ if (!groups.empty())
+ {
+ set<SMESHDS_GroupBase*>::const_iterator GrIt = groups.begin();
+ for (; GrIt != groups.end(); GrIt++)
+ {
+ SMESHDS_Group* grp = dynamic_cast<SMESHDS_Group*>(*GrIt);
+ if (!grp || grp->IsEmpty()) continue;
+ grp->SMDSGroup().Remove(removeelem);
+ }
+ }
+}
+
+//================================================================================
+/*!
+ * \brief Replace elemToRm by elemToAdd in the all groups
+ */
+//================================================================================
+
+void SMESH_MeshEditor::ReplaceElemInGroups (const SMDS_MeshElement* elemToRm,
+ const SMDS_MeshElement* elemToAdd,
+ SMESHDS_Mesh * aMesh)
+{
+ const set<SMESHDS_GroupBase*>& groups = aMesh->GetGroups();
+ if (!groups.empty()) {
+ set<SMESHDS_GroupBase*>::const_iterator grIt = groups.begin();
+ for ( ; grIt != groups.end(); grIt++ ) {
+ SMESHDS_Group* group = dynamic_cast<SMESHDS_Group*>( *grIt );
+ if ( group && group->SMDSGroup().Remove( elemToRm ) && elemToAdd )
+ group->SMDSGroup().Add( elemToAdd );
+ }
+ }
+}
+
+//================================================================================
+/*!
+ * \brief Replace elemToRm by elemToAdd in the all groups
+ */
+//================================================================================
+
+void SMESH_MeshEditor::ReplaceElemInGroups (const SMDS_MeshElement* elemToRm,
+ const vector<const SMDS_MeshElement*>& elemToAdd,
+ SMESHDS_Mesh * aMesh)
+{
+ const set<SMESHDS_GroupBase*>& groups = aMesh->GetGroups();
+ if (!groups.empty())
+ {
+ set<SMESHDS_GroupBase*>::const_iterator grIt = groups.begin();
+ for ( ; grIt != groups.end(); grIt++ ) {
+ SMESHDS_Group* group = dynamic_cast<SMESHDS_Group*>( *grIt );
+ if ( group && group->SMDSGroup().Remove( elemToRm ) )
+ for ( int i = 0; i < elemToAdd.size(); ++i )
+ group->SMDSGroup().Add( elemToAdd[ i ] );
+ }
+ }
+}
+
+//=======================================================================
+//function : QuadToTri
+//purpose : Cut quadrangles into triangles.
+// theCrit is used to select a diagonal to cut
+//=======================================================================
+
+bool SMESH_MeshEditor::QuadToTri (TIDSortedElemSet & theElems,
+ const bool the13Diag)
+{
+ myLastCreatedElems.Clear();
+ myLastCreatedNodes.Clear();
+
+ MESSAGE( "::QuadToTri()" );
+
+ SMESHDS_Mesh * aMesh = GetMeshDS();
+
+ Handle(Geom_Surface) surface;
+ SMESH_MesherHelper helper( *GetMesh() );
+
+ TIDSortedElemSet::iterator itElem;
+ for ( itElem = theElems.begin(); itElem != theElems.end(); itElem++ ) {
+ const SMDS_MeshElement* elem = *itElem;
+ if ( !elem || elem->GetType() != SMDSAbs_Face )
+ continue;
+ bool isquad = elem->NbNodes()==4 || elem->NbNodes()==8;
+ if(!isquad) continue;
+
+ if(elem->NbNodes()==4) {
+ // retrieve element nodes
+ const SMDS_MeshNode* aNodes [4];
+ SMDS_ElemIteratorPtr itN = elem->nodesIterator();
+ int i = 0;
+ while ( itN->more() )
+ aNodes[ i++ ] = static_cast<const SMDS_MeshNode*>( itN->next() );