+ _LayerEdge* leOnV[2] = { getLEdgeOnV(0), getLEdgeOnV(1) };
+
+ // get cosin to use in findEdgesToSmooth()
+ _edgeDir[0] = getEdgeDir( E, leOnV[0]->_nodes[0], data.GetHelper() );
+ _edgeDir[1] = getEdgeDir( E, leOnV[1]->_nodes[0], data.GetHelper() );
+ _leOnV[0]._cosin = Abs( leOnV[0]->_cosin );
+ _leOnV[1]._cosin = Abs( leOnV[1]->_cosin );
+ if ( _eos._sWOL.IsNull() ) // 3D
+ for ( int iEnd = 0; iEnd < 2; ++iEnd )
+ _leOnV[iEnd]._cosin = Abs( _edgeDir[iEnd].Normalized() * leOnV[iEnd]->_normal );
+
+ if ( isAnalytic() )
+ return;
+
+ // divide E to have offset segments with low deflection
+ BRepAdaptor_Curve c3dAdaptor( E );
+ const double curDeflect = 0.1; //0.01; // Curvature deflection == |p1p2]*sin(p1p2,p1pM)
+ const double angDeflect = 0.1; //0.09; // Angular deflection == sin(p1pM,pMp2)
+ GCPnts_TangentialDeflection discret(c3dAdaptor, angDeflect, curDeflect);
+ if ( discret.NbPoints() <= 2 )
+ {
+ _anaCurve = new Geom_Line( gp::OX() ); // only type does matter
+ return;
+ }
+
+ const double u0 = c3dAdaptor.FirstParameter();
+ gp_Pnt p; gp_Vec tangent;
+ if ( discret.NbPoints() >= (int) _eos.size() + 2 )
+ {
+ _offPoints.resize( discret.NbPoints() );
+ for ( size_t i = 0; i < _offPoints.size(); i++ )
+ {
+ double u = discret.Parameter( i+1 );
+ c3dAdaptor.D1( u, p, tangent );
+ _offPoints[i]._xyz = p.XYZ();
+ _offPoints[i]._edgeDir = tangent.XYZ();
+ _offPoints[i]._param = GCPnts_AbscissaPoint::Length( c3dAdaptor, u0, u ) / _curveLen;
+ }
+ }
+ else
+ {
+ std::vector< double > params( _eos.size() + 2 );
+
+ params[0] = data.GetHelper().GetNodeU( E, leOnV[0]->_nodes[0] );
+ params.back() = data.GetHelper().GetNodeU( E, leOnV[1]->_nodes[0] );
+ for ( size_t i = 0; i < _eos.size(); i++ )
+ params[i+1] = data.GetHelper().GetNodeU( E, _eos[i]->_nodes[0] );
+
+ if ( params[1] > params[ _eos.size() ] )
+ std::reverse( params.begin() + 1, params.end() - 1 );
+
+ _offPoints.resize( _eos.size() + 2 );
+ for ( size_t i = 0; i < _offPoints.size(); i++ )
+ {
+ const double u = params[i];
+ c3dAdaptor.D1( u, p, tangent );
+ _offPoints[i]._xyz = p.XYZ();
+ _offPoints[i]._edgeDir = tangent.XYZ();
+ _offPoints[i]._param = GCPnts_AbscissaPoint::Length( c3dAdaptor, u0, u ) / _curveLen;
+ }
+ }
+
+ // set _2edges
+ _offPoints [0]._2edges.set( &_leOnV[0], &_leOnV[0], 0.5, 0.5 );
+ _offPoints.back()._2edges.set( &_leOnV[1], &_leOnV[1], 0.5, 0.5 );
+ _2NearEdges tmp2edges;
+ tmp2edges._edges[1] = _eos._edges[0];
+ _leOnV[0]._2neibors = & tmp2edges;
+ _leOnV[0]._nodes = leOnV[0]->_nodes;
+ _leOnV[1]._nodes = leOnV[1]->_nodes;
+ _LayerEdge* eNext, *ePrev = & _leOnV[0];
+ for ( size_t iLE = 0, i = 1; i < _offPoints.size()-1; i++ )
+ {
+ // find _LayerEdge's located before and after an offset point
+ // (_eos._edges[ iLE ] is next after ePrev)
+ while ( iLE < _eos._edges.size() && _offPoints[i]._param > _leParams[ iLE ] )
+ ePrev = _eos._edges[ iLE++ ];
+ eNext = ePrev->_2neibors->_edges[1];
+
+ gp_Pnt p0 = SMESH_TNodeXYZ( ePrev->_nodes[0] );
+ gp_Pnt p1 = SMESH_TNodeXYZ( eNext->_nodes[0] );
+ double r = p0.Distance( _offPoints[i]._xyz ) / p0.Distance( p1 );
+ _offPoints[i]._2edges.set( ePrev, eNext, 1-r, r );
+ }
+
+ // replace _LayerEdge's on VERTEX by _leOnV in _offPoints._2edges
+ for ( size_t i = 0; i < _offPoints.size(); i++ )
+ if ( _offPoints[i]._2edges._edges[0] == leOnV[0] )
+ _offPoints[i]._2edges._edges[0] = & _leOnV[0];
+ else break;
+ for ( size_t i = _offPoints.size()-1; i > 0; i-- )
+ if ( _offPoints[i]._2edges._edges[1] == leOnV[1] )
+ _offPoints[i]._2edges._edges[1] = & _leOnV[1];
+ else break;
+
+ // set _normal of _leOnV[0] and _leOnV[1] to be normal to the EDGE
+
+ int iLBO = _offPoints.size() - 2; // last but one
+
+ if ( leOnV[ 0 ]->Is( _LayerEdge::MULTI_NORMAL ))
+ _leOnV[ 0 ]._normal = getNormalNormal( _eos._edges[1]->_normal, _edgeDir[0] );
+ else
+ _leOnV[ 0 ]._normal = getNormalNormal( leOnV[0]->_normal, _edgeDir[0] );
+ if ( leOnV[ 1 ]->Is( _LayerEdge::MULTI_NORMAL ))
+ _leOnV[ 1 ]._normal = getNormalNormal( _eos._edges.back()->_normal, _edgeDir[1] );
+ else
+ _leOnV[ 1 ]._normal = getNormalNormal( leOnV[1]->_normal, _edgeDir[1] );
+ _leOnV[ 0 ]._len = 0;
+ _leOnV[ 1 ]._len = 0;
+ _leOnV[ 0 ]._lenFactor = _offPoints[1 ]._2edges._edges[1]->_lenFactor;
+ _leOnV[ 1 ]._lenFactor = _offPoints[iLBO]._2edges._edges[0]->_lenFactor;
+
+ _iSeg[0] = 0;
+ _iSeg[1] = _offPoints.size()-2;
+
+ // initialize OffPnt::_len
+ for ( size_t i = 0; i < _offPoints.size(); ++i )
+ _offPoints[i]._len = 0;
+
+ if ( _eos._edges[0]->NbSteps() > 1 ) // already inflated several times, init _xyz
+ {
+ _leOnV[0]._len = leOnV[0]->_len;
+ _leOnV[1]._len = leOnV[1]->_len;
+ for ( size_t i = 0; i < _offPoints.size(); i++ )
+ {
+ _LayerEdge* e0 = _offPoints[i]._2edges._edges[0];
+ _LayerEdge* e1 = _offPoints[i]._2edges._edges[1];
+ const double w0 = _offPoints[i]._2edges._wgt[0];
+ const double w1 = _offPoints[i]._2edges._wgt[1];
+ double avgLen = ( e0->_len * w0 + e1->_len * w1 );
+ gp_XYZ avgXYZ = ( SMESH_TNodeXYZ( e0->_nodes.back() ) * w0 +
+ SMESH_TNodeXYZ( e1->_nodes.back() ) * w1 );
+ _offPoints[i]._xyz = avgXYZ;
+ _offPoints[i]._len = avgLen;
+ }
+ }
+}
+
+//================================================================================
+/*!
+ * \brief return _normal of _leOnV[is2nd] normal to the EDGE
+ */
+//================================================================================
+
+gp_XYZ _Smoother1D::getNormalNormal( const gp_XYZ & normal,
+ const gp_XYZ& edgeDir)
+{
+ gp_XYZ cross = normal ^ edgeDir;
+ gp_XYZ norm = edgeDir ^ cross;
+ double size = norm.Modulus();
+
+ // if ( size == 0 ) // MULTI_NORMAL _LayerEdge
+ // return gp_XYZ( 1e-100, 1e-100, 1e-100 );
+
+ return norm / size;
+}
+
+//================================================================================
+/*!
+ * \brief Sort _LayerEdge's by a parameter on a given EDGE
+ */
+//================================================================================
+
+void _SolidData::SortOnEdge( const TopoDS_Edge& E,
+ vector< _LayerEdge* >& edges)
+{
+ map< double, _LayerEdge* > u2edge;
+ for ( size_t i = 0; i < edges.size(); ++i )
+ u2edge.insert( u2edge.end(),
+ make_pair( _helper->GetNodeU( E, edges[i]->_nodes[0] ), edges[i] ));
+
+ ASSERT( u2edge.size() == edges.size() );
+ map< double, _LayerEdge* >::iterator u2e = u2edge.begin();
+ for ( size_t i = 0; i < edges.size(); ++i, ++u2e )
+ edges[i] = u2e->second;
+
+ Sort2NeiborsOnEdge( edges );
+}
+
+//================================================================================
+/*!
+ * \brief Set _2neibors according to the order of _LayerEdge on EDGE
+ */
+//================================================================================
+
+void _SolidData::Sort2NeiborsOnEdge( vector< _LayerEdge* >& edges )
+{
+ if ( edges.size() < 2 || !edges[0]->_2neibors ) return;
+
+ for ( size_t i = 0; i < edges.size()-1; ++i )
+ if ( edges[i]->_2neibors->tgtNode(1) != edges[i+1]->_nodes.back() )
+ edges[i]->_2neibors->reverse();
+
+ const size_t iLast = edges.size() - 1;
+ if ( edges.size() > 1 &&
+ edges[iLast]->_2neibors->tgtNode(0) != edges[iLast-1]->_nodes.back() )
+ edges[iLast]->_2neibors->reverse();
+}
+
+//================================================================================
+/*!
+ * \brief Return _EdgesOnShape* corresponding to the shape
+ */
+//================================================================================
+
+_EdgesOnShape* _SolidData::GetShapeEdges(const TGeomID shapeID )
+{
+ if ( shapeID < (int)_edgesOnShape.size() &&
+ _edgesOnShape[ shapeID ]._shapeID == shapeID )
+ return _edgesOnShape[ shapeID ]._subMesh ? & _edgesOnShape[ shapeID ] : 0;
+
+ for ( size_t i = 0; i < _edgesOnShape.size(); ++i )
+ if ( _edgesOnShape[i]._shapeID == shapeID )
+ return _edgesOnShape[i]._subMesh ? & _edgesOnShape[i] : 0;
+
+ return 0;
+}
+
+//================================================================================
+/*!
+ * \brief Return _EdgesOnShape* corresponding to the shape
+ */
+//================================================================================
+
+_EdgesOnShape* _SolidData::GetShapeEdges(const TopoDS_Shape& shape )
+{
+ SMESHDS_Mesh* meshDS = _proxyMesh->GetMesh()->GetMeshDS();
+ return GetShapeEdges( meshDS->ShapeToIndex( shape ));
+}
+
+//================================================================================
+/*!
+ * \brief Prepare data of the _LayerEdge for smoothing on FACE
+ */
+//================================================================================
+
+void _SolidData::PrepareEdgesToSmoothOnFace( _EdgesOnShape* eos, bool substituteSrcNodes )
+{
+ SMESH_MesherHelper helper( *_proxyMesh->GetMesh() );
+
+ set< TGeomID > vertices;
+ TopoDS_Face F;
+ if ( eos->ShapeType() == TopAbs_FACE )
+ {
+ // check FACE concavity and get concave VERTEXes
+ F = TopoDS::Face( eos->_shape );
+ if ( isConcave( F, helper, &vertices ))
+ _concaveFaces.insert( eos->_shapeID );
+
+ // set eos._eosConcaVer
+ eos->_eosConcaVer.clear();
+ eos->_eosConcaVer.reserve( vertices.size() );
+ for ( set< TGeomID >::iterator v = vertices.begin(); v != vertices.end(); ++v )
+ {
+ _EdgesOnShape* eov = GetShapeEdges( *v );
+ if ( eov && eov->_edges.size() == 1 )
+ {
+ eos->_eosConcaVer.push_back( eov );
+ for ( size_t i = 0; i < eov->_edges[0]->_neibors.size(); ++i )
+ eov->_edges[0]->_neibors[i]->Set( _LayerEdge::DIFFICULT );
+ }
+ }
+
+ // SetSmooLen() to _LayerEdge's on FACE
+ for ( size_t i = 0; i < eos->_edges.size(); ++i )
+ {
+ eos->_edges[i]->SetSmooLen( Precision::Infinite() );
+ }
+ SMESH_subMeshIteratorPtr smIt = eos->_subMesh->getDependsOnIterator(/*includeSelf=*/false);
+ while ( smIt->more() ) // loop on sub-shapes of the FACE
+ {
+ _EdgesOnShape* eoe = GetShapeEdges( smIt->next()->GetId() );
+ if ( !eoe ) continue;
+
+ vector<_LayerEdge*>& eE = eoe->_edges;
+ for ( size_t iE = 0; iE < eE.size(); ++iE ) // loop on _LayerEdge's on EDGE or VERTEX
+ {
+ if ( eE[iE]->_cosin <= theMinSmoothCosin )
+ continue;
+
+ SMDS_ElemIteratorPtr segIt = eE[iE]->_nodes[0]->GetInverseElementIterator(SMDSAbs_Edge);
+ while ( segIt->more() )
+ {
+ const SMDS_MeshElement* seg = segIt->next();
+ if ( !eos->_subMesh->DependsOn( seg->getshapeId() ))
+ continue;
+ if ( seg->GetNode(0) != eE[iE]->_nodes[0] )
+ continue; // not to check a seg twice
+ for ( size_t iN = 0; iN < eE[iE]->_neibors.size(); ++iN )
+ {
+ _LayerEdge* eN = eE[iE]->_neibors[iN];
+ if ( eN->_nodes[0]->getshapeId() != eos->_shapeID )
+ continue;
+ double dist = SMESH_MeshAlgos::GetDistance( seg, SMESH_TNodeXYZ( eN->_nodes[0] ));
+ double smooLen = getSmoothingThickness( eE[iE]->_cosin, dist );
+ eN->SetSmooLen( Min( smooLen, eN->GetSmooLen() ));
+ eN->Set( _LayerEdge::NEAR_BOUNDARY );
+ }
+ }
+ }
+ }
+ } // if ( eos->ShapeType() == TopAbs_FACE )
+
+ for ( size_t i = 0; i < eos->_edges.size(); ++i )
+ {
+ eos->_edges[i]->_smooFunction = 0;
+ eos->_edges[i]->Set( _LayerEdge::TO_SMOOTH );
+ }
+ bool isCurved = false;
+ for ( size_t i = 0; i < eos->_edges.size(); ++i )
+ {
+ _LayerEdge* edge = eos->_edges[i];
+
+ // get simplices sorted
+ _Simplex::SortSimplices( edge->_simplices );
+
+ // smoothing function
+ edge->ChooseSmooFunction( vertices, _n2eMap );
+
+ // set _curvature
+ double avgNormProj = 0, avgLen = 0;
+ for ( size_t iS = 0; iS < edge->_simplices.size(); ++iS )
+ {
+ _Simplex& s = edge->_simplices[iS];
+
+ gp_XYZ vec = edge->_pos.back() - SMESH_TNodeXYZ( s._nPrev );
+ avgNormProj += edge->_normal * vec;
+ avgLen += vec.Modulus();
+ if ( substituteSrcNodes )
+ {
+ s._nNext = _n2eMap[ s._nNext ]->_nodes.back();
+ s._nPrev = _n2eMap[ s._nPrev ]->_nodes.back();
+ }
+ }
+ avgNormProj /= edge->_simplices.size();
+ avgLen /= edge->_simplices.size();
+ if (( edge->_curvature = _Curvature::New( avgNormProj, avgLen )))
+ {
+ isCurved = true;
+ SMDS_FacePosition* fPos = dynamic_cast<SMDS_FacePosition*>( edge->_nodes[0]->GetPosition() );
+ if ( !fPos )
+ for ( size_t iS = 0; iS < edge->_simplices.size() && !fPos; ++iS )
+ fPos = dynamic_cast<SMDS_FacePosition*>( edge->_simplices[iS]._nPrev->GetPosition() );
+ if ( fPos )
+ edge->_curvature->_uv.SetCoord( fPos->GetUParameter(), fPos->GetVParameter() );
+ }
+ }
+
+ // prepare for putOnOffsetSurface()
+ if (( eos->ShapeType() == TopAbs_FACE ) &&
+ ( isCurved || !eos->_eosConcaVer.empty() ))
+ {
+ eos->_offsetSurf = helper.GetSurface( TopoDS::Face( eos->_shape ));
+ eos->_edgeForOffset = 0;
+
+ double maxCosin = -1;
+ for ( TopExp_Explorer eExp( eos->_shape, TopAbs_EDGE ); eExp.More(); eExp.Next() )
+ {
+ _EdgesOnShape* eoe = GetShapeEdges( eExp.Current() );
+ if ( !eoe || eoe->_edges.empty() ) continue;
+
+ vector<_LayerEdge*>& eE = eoe->_edges;
+ _LayerEdge* e = eE[ eE.size() / 2 ];
+ if ( e->_cosin > maxCosin )
+ {
+ eos->_edgeForOffset = e;
+ maxCosin = e->_cosin;
+ }
+ }
+ }
+}
+
+//================================================================================
+/*!
+ * \brief Add faces for smoothing
+ */
+//================================================================================
+
+void _SolidData::AddShapesToSmooth( const set< _EdgesOnShape* >& eosToSmooth,
+ const set< _EdgesOnShape* >* edgesNoAnaSmooth )
+{
+ set< _EdgesOnShape * >::const_iterator eos = eosToSmooth.begin();
+ for ( ; eos != eosToSmooth.end(); ++eos )
+ {
+ if ( !*eos || (*eos)->_toSmooth ) continue;
+
+ (*eos)->_toSmooth = true;
+
+ if ( (*eos)->ShapeType() == TopAbs_FACE )
+ {
+ PrepareEdgesToSmoothOnFace( *eos, /*substituteSrcNodes=*/false );
+ (*eos)->_toSmooth = true;
+ }
+ }
+
+ // avoid _Smoother1D::smoothAnalyticEdge() of edgesNoAnaSmooth
+ if ( edgesNoAnaSmooth )
+ for ( eos = edgesNoAnaSmooth->begin(); eos != edgesNoAnaSmooth->end(); ++eos )
+ {
+ if ( (*eos)->_edgeSmoother )
+ (*eos)->_edgeSmoother->_anaCurve.Nullify();
+ }
+}
+
+//================================================================================
+/*!
+ * \brief Limit _LayerEdge::_maxLen according to local curvature
+ */
+//================================================================================
+
+void _ViscousBuilder::limitMaxLenByCurvature( _SolidData& data, SMESH_MesherHelper& helper )
+{
+ // find intersection of neighbor _LayerEdge's to limit _maxLen
+ // according to local curvature (IPAL52648)
+
+ // This method must be called after findCollisionEdges() where _LayerEdge's
+ // get _lenFactor initialized in the case of eos._hyp.IsOffsetMethod()
+
+ for ( size_t iS = 0; iS < data._edgesOnShape.size(); ++iS )
+ {
+ _EdgesOnShape& eosI = data._edgesOnShape[iS];
+ if ( eosI._edges.empty() ) continue;
+ if ( !eosI._hyp.ToSmooth() )
+ {
+ for ( size_t i = 0; i < eosI._edges.size(); ++i )
+ {
+ _LayerEdge* eI = eosI._edges[i];
+ for ( size_t iN = 0; iN < eI->_neibors.size(); ++iN )
+ {
+ _LayerEdge* eN = eI->_neibors[iN];
+ if ( eI->_nodes[0]->GetID() < eN->_nodes[0]->GetID() ) // treat this pair once
+ {
+ _EdgesOnShape* eosN = data.GetShapeEdges( eN );
+ limitMaxLenByCurvature( eI, eN, eosI, *eosN, helper );
+ }
+ }
+ }
+ }
+ else if ( eosI.ShapeType() == TopAbs_EDGE )
+ {
+ const TopoDS_Edge& E = TopoDS::Edge( eosI._shape );
+ if ( SMESH_Algo::IsStraight( E, /*degenResult=*/true )) continue;
+
+ _LayerEdge* e0 = eosI._edges[0];
+ for ( size_t i = 1; i < eosI._edges.size(); ++i )
+ {
+ _LayerEdge* eI = eosI._edges[i];
+ limitMaxLenByCurvature( eI, e0, eosI, eosI, helper );
+ e0 = eI;
+ }
+ }
+ }
+}
+
+//================================================================================
+/*!
+ * \brief Limit _LayerEdge::_maxLen according to local curvature
+ */
+//================================================================================
+
+void _ViscousBuilder::limitMaxLenByCurvature( _LayerEdge* e1,
+ _LayerEdge* e2,
+ _EdgesOnShape& eos1,
+ _EdgesOnShape& eos2,
+ SMESH_MesherHelper& helper )
+{
+ gp_XYZ plnNorm = e1->_normal ^ e2->_normal;
+ double norSize = plnNorm.SquareModulus();
+ if ( norSize < std::numeric_limits<double>::min() )
+ return; // parallel normals
+
+ // find closest points of skew _LayerEdge's
+ SMESH_TNodeXYZ src1( e1->_nodes[0] ), src2( e2->_nodes[0] );
+ gp_XYZ dir12 = src2 - src1;
+ gp_XYZ perp1 = e1->_normal ^ plnNorm;
+ gp_XYZ perp2 = e2->_normal ^ plnNorm;
+ double dot1 = perp2 * e1->_normal;
+ double dot2 = perp1 * e2->_normal;
+ double u1 = ( perp2 * dir12 ) / dot1;
+ double u2 = - ( perp1 * dir12 ) / dot2;
+ if ( u1 > 0 && u2 > 0 )
+ {
+ double ovl = ( u1 * e1->_normal * dir12 -
+ u2 * e2->_normal * dir12 ) / dir12.SquareModulus();
+ if ( ovl > theSmoothThickToElemSizeRatio )
+ {
+ e1->_maxLen = Min( e1->_maxLen, 0.75 * u1 / e1->_lenFactor );
+ e2->_maxLen = Min( e2->_maxLen, 0.75 * u2 / e2->_lenFactor );
+ }
+ }
+}
+
+//================================================================================
+/*!
+ * \brief Fill data._collisionEdges
+ */
+//================================================================================
+
+void _ViscousBuilder::findCollisionEdges( _SolidData& data, SMESH_MesherHelper& helper )
+{
+ data._collisionEdges.clear();
+
+ // set the full thickness of the layers to LEs
+ for ( size_t iS = 0; iS < data._edgesOnShape.size(); ++iS )
+ {
+ _EdgesOnShape& eos = data._edgesOnShape[iS];
+ if ( eos._edges.empty() ) continue;
+ if ( eos.ShapeType() != TopAbs_EDGE && eos.ShapeType() != TopAbs_VERTEX ) continue;
+
+ for ( size_t i = 0; i < eos._edges.size(); ++i )
+ {
+ if ( eos._edges[i]->Is( _LayerEdge::BLOCKED )) continue;
+ double maxLen = eos._edges[i]->_maxLen;
+ eos._edges[i]->_maxLen = Precision::Infinite(); // avoid blocking
+ eos._edges[i]->SetNewLength( 1.5 * maxLen, eos, helper );
+ eos._edges[i]->_maxLen = maxLen;
+ }
+ }
+
+ // make temporary quadrangles got by extrusion of
+ // mesh edges along _LayerEdge._normal's
+
+ vector< const SMDS_MeshElement* > tmpFaces;
+
+ for ( size_t iS = 0; iS < data._edgesOnShape.size(); ++iS )
+ {
+ _EdgesOnShape& eos = data._edgesOnShape[ iS ];
+ if ( eos.ShapeType() != TopAbs_EDGE )
+ continue;
+ if ( eos._edges.empty() )
+ {
+ _LayerEdge* edge[2] = { 0, 0 }; // LE of 2 VERTEX'es
+ SMESH_subMeshIteratorPtr smIt = eos._subMesh->getDependsOnIterator(/*includeSelf=*/false);
+ while ( smIt->more() )
+ if ( _EdgesOnShape* eov = data.GetShapeEdges( smIt->next()->GetId() ))
+ if ( eov->_edges.size() == 1 )
+ edge[ bool( edge[0]) ] = eov->_edges[0];
+
+ if ( edge[1] )
+ {
+ _TmpMeshFaceOnEdge* f = new _TmpMeshFaceOnEdge( edge[0], edge[1], --_tmpFaceID );
+ tmpFaces.push_back( f );
+ }
+ }
+ for ( size_t i = 0; i < eos._edges.size(); ++i )
+ {
+ _LayerEdge* edge = eos._edges[i];
+ for ( int j = 0; j < 2; ++j ) // loop on _2NearEdges
+ {
+ const SMDS_MeshNode* src2 = edge->_2neibors->srcNode(j);
+ if ( src2->GetPosition()->GetDim() > 0 &&
+ src2->GetID() < edge->_nodes[0]->GetID() )
+ continue; // avoid using same segment twice
+
+ // a _LayerEdge containg tgt2
+ _LayerEdge* neiborEdge = edge->_2neibors->_edges[j];
+
+ _TmpMeshFaceOnEdge* f = new _TmpMeshFaceOnEdge( edge, neiborEdge, --_tmpFaceID );
+ tmpFaces.push_back( f );
+ }
+ }
+ }
+
+ // Find _LayerEdge's intersecting tmpFaces.
+
+ SMDS_ElemIteratorPtr fIt( new SMDS_ElementVectorIterator( tmpFaces.begin(),
+ tmpFaces.end()));
+ SMESHUtils::Deleter<SMESH_ElementSearcher> searcher
+ ( SMESH_MeshAlgos::GetElementSearcher( *getMeshDS(), fIt ));
+
+ double dist1, dist2, segLen, eps = 0.5;
+ _CollisionEdges collEdges;
+ vector< const SMDS_MeshElement* > suspectFaces;
+ const double angle45 = Cos( 45. * M_PI / 180. );
+
+ for ( size_t iS = 0; iS < data._edgesOnShape.size(); ++iS )
+ {
+ _EdgesOnShape& eos = data._edgesOnShape[ iS ];
+ if ( eos.ShapeType() == TopAbs_FACE || !eos._sWOL.IsNull() )
+ continue;
+ // find sub-shapes whose VL can influence VL on eos
+ set< TGeomID > neighborShapes;
+ PShapeIteratorPtr fIt = helper.GetAncestors( eos._shape, *_mesh, TopAbs_FACE );
+ while ( const TopoDS_Shape* face = fIt->next() )
+ {
+ TGeomID faceID = getMeshDS()->ShapeToIndex( *face );
+ if ( _EdgesOnShape* eof = data.GetShapeEdges( faceID ))
+ {
+ SMESH_subMeshIteratorPtr subIt = eof->_subMesh->getDependsOnIterator(/*includeSelf=*/false);
+ while ( subIt->more() )
+ neighborShapes.insert( subIt->next()->GetId() );
+ }
+ }
+ if ( eos.ShapeType() == TopAbs_VERTEX )
+ {
+ PShapeIteratorPtr eIt = helper.GetAncestors( eos._shape, *_mesh, TopAbs_EDGE );
+ while ( const TopoDS_Shape* edge = eIt->next() )
+ neighborShapes.erase( getMeshDS()->ShapeToIndex( *edge ));
+ }
+ // find intersecting _LayerEdge's
+ for ( size_t i = 0; i < eos._edges.size(); ++i )
+ {
+ if ( eos._edges[i]->Is( _LayerEdge::MULTI_NORMAL )) continue;
+ _LayerEdge* edge = eos._edges[i];
+ gp_Ax1 lastSegment = edge->LastSegment( segLen, eos );
+ segLen *= 1.2;
+
+ gp_Vec eSegDir0, eSegDir1;
+ if ( edge->IsOnEdge() )
+ {
+ SMESH_TNodeXYZ eP( edge->_nodes[0] );
+ eSegDir0 = SMESH_TNodeXYZ( edge->_2neibors->srcNode(0) ) - eP;
+ eSegDir1 = SMESH_TNodeXYZ( edge->_2neibors->srcNode(1) ) - eP;
+ }
+ suspectFaces.clear();
+ searcher->GetElementsInSphere( SMESH_TNodeXYZ( edge->_nodes.back()), edge->_len * 2,
+ SMDSAbs_Face, suspectFaces );
+ collEdges._intEdges.clear();
+ for ( size_t j = 0 ; j < suspectFaces.size(); ++j )
+ {
+ const _TmpMeshFaceOnEdge* f = (const _TmpMeshFaceOnEdge*) suspectFaces[j];
+ if ( f->_le1 == edge || f->_le2 == edge ) continue;
+ if ( !neighborShapes.count( f->_le1->_nodes[0]->getshapeId() )) continue;
+ if ( !neighborShapes.count( f->_le2->_nodes[0]->getshapeId() )) continue;
+ if ( edge->IsOnEdge() ) {
+ if ( edge->_2neibors->include( f->_le1 ) ||
+ edge->_2neibors->include( f->_le2 )) continue;
+ }
+ else {
+ if (( f->_le1->IsOnEdge() && f->_le1->_2neibors->include( edge )) ||
+ ( f->_le2->IsOnEdge() && f->_le2->_2neibors->include( edge ))) continue;
+ }
+ dist1 = dist2 = Precision::Infinite();
+ if ( !edge->SegTriaInter( lastSegment, f->_nn[0], f->_nn[1], f->_nn[2], dist1, eps ))
+ dist1 = Precision::Infinite();
+ if ( !edge->SegTriaInter( lastSegment, f->_nn[3], f->_nn[2], f->_nn[0], dist2, eps ))
+ dist2 = Precision::Infinite();
+ if (( dist1 > segLen ) && ( dist2 > segLen ))
+ continue;
+
+ if ( edge->IsOnEdge() )
+ {
+ // skip perpendicular EDGEs
+ gp_Vec fSegDir = SMESH_TNodeXYZ( f->_nn[0] ) - SMESH_TNodeXYZ( f->_nn[3] );
+ bool isParallel = ( isLessAngle( eSegDir0, fSegDir, angle45 ) ||
+ isLessAngle( eSegDir1, fSegDir, angle45 ) ||
+ isLessAngle( eSegDir0, fSegDir.Reversed(), angle45 ) ||
+ isLessAngle( eSegDir1, fSegDir.Reversed(), angle45 ));
+ if ( !isParallel )
+ continue;
+ }
+
+ // either limit inflation of edges or remember them for updating _normal
+ // double dot = edge->_normal * f->GetDir();
+ // if ( dot > 0.1 )
+ {
+ collEdges._intEdges.push_back( f->_le1 );
+ collEdges._intEdges.push_back( f->_le2 );
+ }