diff --git a/src/Extrema/Extrema_GenExtCC.gxx b/src/Extrema/Extrema_GenExtCC.gxx index 969a5b85ae..33826fc5c8 100644 --- a/src/Extrema/Extrema_GenExtCC.gxx +++ b/src/Extrema/Extrema_GenExtCC.gxx @@ -163,7 +163,6 @@ static Standard_Real ProjPOnC(const Pnt& theP, if (aD < aDist) aDist = aD; } - aDist = sqrt(aDist); } return aDist; } @@ -340,10 +339,20 @@ void Extrema_GenExtCC::Perform() aNbInter[1] = anIntervals2->Length() - 1; } } + if (C1.IsClosed() && aNbInter[0] == 1) + { + ChangeIntervals(anIntervals1, 3); + aNbInter[0] = anIntervals1->Length() - 1; + } + if (C2.IsClosed() && aNbInter[1] == 1) + { + ChangeIntervals(anIntervals2, 3); + aNbInter[1] = anIntervals2->Length() - 1; + } // Lipchitz constant computation. const Standard_Real aMaxLC = 10000.; - Standard_Real aLC = 9.0; // Default value. + Standard_Real aLC = 100.0; // Default value. const Standard_Real aMaxDer1 = 1.0 / C1.Resolution(1.0); const Standard_Real aMaxDer2 = 1.0 / C2.Resolution(1.0); Standard_Real aMaxDer = Max(aMaxDer1, aMaxDer2) * Sqrt(2.0); @@ -383,6 +392,43 @@ void Extrema_GenExtCC::Perform() } Extrema_GlobOptFuncCCC2 aFunc (C1, C2); + if (aLC < aMaxLC || aMaxDer > aMaxLC) + { + //Estimation of Lipschitz constant by gradient of optimization function + //using sampling in parameter space. + math_Vector aT(1, 2), aG(1, 2); + Standard_Real aF, aMaxG = 0.; + Standard_Real t1, t2, dt1, dt2; + Standard_Integer n1 = 21, n2 = 21, i1, i2; + dt1 = (C1.LastParameter() - C1.FirstParameter()) / (n1 - 1); + dt2 = (C2.LastParameter() - C2.FirstParameter()) / (n2 - 1); + for (i1 = 1, t1 = C1.FirstParameter(); i1 <= n1; ++i1, t1 += dt1) + { + aT(1) = t1; + for (i2 = 1, t2 = C2.FirstParameter(); i2 <= n2; ++i2, t2 += dt2) + { + aT(2) = t2; + aFunc.Values(aT, aF, aG); + Standard_Real aMod = aG(1)*aG(1) + aG(2)*aG(2); + aMaxG = Max(aMaxG, aMod); + } + } + aMaxG = Sqrt(aMaxG); + if (aMaxG > aMaxDer) + { + aLC = Min(aMaxG, aMaxLC); + isConstLockedFlag = Standard_True; + } + if (aMaxG > 100. * aMaxLC) + { + aLC = 100. * aMaxLC; + isConstLockedFlag = Standard_True; + } + else if (aMaxG < 0.1 * aMaxDer) + { + isConstLockedFlag = Standard_True; + } + } math_GlobOptMin aFinder(&aFunc, myLowBorder, myUppBorder, aLC); aFinder.SetLipConstState(isConstLockedFlag); aFinder.SetContinuity(aContinuity == GeomAbs_C2 ? 2 : 1); @@ -518,7 +564,6 @@ void Extrema_GenExtCC::Perform() aVec(2) = (aCurrent.Y() + aNext.Y()) * 0.5; aFunc.Value(aVec, aVal); - if (Abs(aVal - aF) < Precision::Confusion()) { // It seems the parallel segment is found. diff --git a/src/Extrema/Extrema_GlobOptFuncCC.cxx b/src/Extrema/Extrema_GlobOptFuncCC.cxx index 099f671e87..93b5fd60b0 100644 --- a/src/Extrema/Extrema_GlobOptFuncCC.cxx +++ b/src/Extrema/Extrema_GlobOptFuncCC.cxx @@ -43,7 +43,7 @@ static Standard_Boolean _Value(const Adaptor3d_Curve& C1, return Standard_False; } - F = C2.Value(v).Distance(C1.Value(u)); + F = C2.Value(v).SquareDistance(C1.Value(u)); return Standard_True; } @@ -64,7 +64,7 @@ static Standard_Boolean _Value(const Adaptor2d_Curve2d& C1, return Standard_False; } - F = C2.Value(v).Distance(C1.Value(u)); + F = C2.Value(v).SquareDistance(C1.Value(u)); return Standard_True; } @@ -89,13 +89,14 @@ static Standard_Boolean _Gradient(const Adaptor3d_Curve& C1, C1.D1(X(1), C1D0, C1D1); C2.D1(X(2), C2D0, C2D1); - + G(1) = - (C2D0.X() - C1D0.X()) * C1D1.X() - (C2D0.Y() - C1D0.Y()) * C1D1.Y() - (C2D0.Z() - C1D0.Z()) * C1D1.Z(); G(2) = (C2D0.X() - C1D0.X()) * C2D1.X() + (C2D0.Y() - C1D0.Y()) * C2D1.Y() + (C2D0.Z() - C1D0.Z()) * C2D1.Z(); + G *= 2.; return Standard_True; } @@ -121,8 +122,11 @@ static Standard_Boolean _Gradient(const Adaptor2d_Curve2d& C1, G(1) = - (C2D0.X() - C1D0.X()) * C1D1.X() - (C2D0.Y() - C1D0.Y()) * C1D1.Y(); + G(2) = (C2D0.X() - C1D0.X()) * C2D1.X() + (C2D0.Y() - C1D0.Y()) * C2D1.Y(); + G *= 2.; + return Standard_True; } @@ -166,6 +170,7 @@ static Standard_Boolean _Hessian (const Adaptor3d_Curve& C1, + (C2D0.X() - C1D0.X()) * C2D2.X() + (C2D0.Y() - C1D0.Y()) * C2D2.Y() + (C2D0.Z() - C1D0.Z()) * C2D2.Z(); + H *= 2.; return Standard_True; } @@ -204,10 +209,11 @@ static Standard_Boolean _Hessian (const Adaptor2d_Curve2d& C1, + C2D1.Y() * C2D1.Y() + (C2D0.X() - C1D0.X()) * C2D2.X() + (C2D0.Y() - C1D0.Y()) * C2D2.Y(); + H *= 2.; return Standard_True; } -// C0 +//C0 //======================================================================= //function : Extrema_GlobOptFuncCCC0 @@ -415,6 +421,5 @@ Standard_Boolean Extrema_GlobOptFuncCCC2::Values(const math_Vector& X,Standard_R else isHessianComputed = _Hessian(*myC1_2d, *myC2_2d, X, H); - return (Value(X, F) && Gradient(X, G) && isHessianComputed); } diff --git a/src/IntTools/IntTools_FaceFace.cxx b/src/IntTools/IntTools_FaceFace.cxx index 443783c26e..5e638519ff 100644 --- a/src/IntTools/IntTools_FaceFace.cxx +++ b/src/IntTools/IntTools_FaceFace.cxx @@ -1169,7 +1169,6 @@ void IntTools_FaceFace::MakeCurve(const Standard_Integer Index, tol2d = myTolApprox; } - Standard_Real aReachedTol = Precision::Confusion(); bIsDecomposited = IntTools_WLineTool:: DecompositionOfWLine(WL, myHS1, @@ -1180,7 +1179,6 @@ void IntTools_FaceFace::MakeCurve(const Standard_Integer Index, bAvoidLineConstructor, myTol, aSeqOfL, - aReachedTol, myContext); // aNbSeqOfL=aSeqOfL.Length(); @@ -1188,7 +1186,7 @@ void IntTools_FaceFace::MakeCurve(const Standard_Integer Index, Standard_Real aTolC = 0.; if (bIsDecomposited) { nbiter=aNbSeqOfL; - aTolC = aReachedTol; + aTolC = Precision::Confusion(); } else { nbiter=1; diff --git a/src/IntTools/IntTools_WLineTool.cxx b/src/IntTools/IntTools_WLineTool.cxx index 8e589ffcc8..8d8e1a6b31 100644 --- a/src/IntTools/IntTools_WLineTool.cxx +++ b/src/IntTools/IntTools_WLineTool.cxx @@ -228,203 +228,6 @@ Standard_Boolean IntTools_WLineTool::NotUseSurfacesForApprox(const TopoDS_Face& /////////////////////// DecompositionOfWLine //////////////////////////// -//======================================================================= -//function : CheckTangentZonesExist -//purpose : static subfunction in ComputeTangentZones -//======================================================================= -static -Standard_Boolean CheckTangentZonesExist(const Handle(GeomAdaptor_Surface)& theSurface1, - const Handle(GeomAdaptor_Surface)& theSurface2) -{ - if ( ( theSurface1->GetType() != GeomAbs_Torus ) || - ( theSurface2->GetType() != GeomAbs_Torus ) ) - return Standard_False; - - gp_Torus aTor1 = theSurface1->Torus(); - gp_Torus aTor2 = theSurface2->Torus(); - - if ( aTor1.Location().Distance( aTor2.Location() ) > Precision::Confusion() ) - return Standard_False; - - if ( ( fabs( aTor1.MajorRadius() - aTor2.MajorRadius() ) > Precision::Confusion() ) || - ( fabs( aTor1.MinorRadius() - aTor2.MinorRadius() ) > Precision::Confusion() ) ) - return Standard_False; - - if ( ( aTor1.MajorRadius() < aTor1.MinorRadius() ) || - ( aTor2.MajorRadius() < aTor2.MinorRadius() ) ) - return Standard_False; - - return Standard_True; -} - - -//======================================================================= -//function : ComputeTangentZones -//purpose : static subfunction in DecompositionOfWLine -//======================================================================= -static -Standard_Integer ComputeTangentZones( const Handle(GeomAdaptor_Surface)& theSurface1, - const Handle(GeomAdaptor_Surface)& theSurface2, - const TopoDS_Face& theFace1, - const TopoDS_Face& theFace2, - Handle(TColgp_HArray1OfPnt2d)& theResultOnS1, - Handle(TColgp_HArray1OfPnt2d)& theResultOnS2, - Handle(TColStd_HArray1OfReal)& theResultRadius, - const Handle(IntTools_Context)& aContext) -{ - Standard_Integer aResult = 0; - if ( !CheckTangentZonesExist( theSurface1, theSurface2 ) ) - return aResult; - - - TColgp_SequenceOfPnt2d aSeqResultS1, aSeqResultS2; - TColStd_SequenceOfReal aSeqResultRad; - - gp_Torus aTor1 = theSurface1->Torus(); - gp_Torus aTor2 = theSurface2->Torus(); - - gp_Ax2 anax1( aTor1.Location(), aTor1.Axis().Direction() ); - gp_Ax2 anax2( aTor2.Location(), aTor2.Axis().Direction() ); - Standard_Integer j = 0; - - for ( j = 0; j < 2; j++ ) { - Standard_Real aCoef = ( j == 0 ) ? -1 : 1; - Standard_Real aRadius1 = fabs(aTor1.MajorRadius() + aCoef * aTor1.MinorRadius()); - Standard_Real aRadius2 = fabs(aTor2.MajorRadius() + aCoef * aTor2.MinorRadius()); - - gp_Circ aCircle1( anax1, aRadius1 ); - gp_Circ aCircle2( anax2, aRadius2 ); - - // roughly compute radius of tangent zone for perpendicular case - Standard_Real aCriteria = Precision::Confusion() * 0.5; - - Standard_Real aT1 = aCriteria; - Standard_Real aT2 = aCriteria; - if ( j == 0 ) { - // internal tangency - Standard_Real aR = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius(); - //aT1 = aCriteria * aCriteria + aR * aR - ( aR - aCriteria ) * ( aR - aCriteria ); - aT1 = 2. * aR * aCriteria; - aT2 = aT1; - } - else { - // external tangency - Standard_Real aRb = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius(); - Standard_Real aRm = ( aRadius1 < aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius(); - Standard_Real aDelta = aRb - aCriteria; - aDelta *= aDelta; - aDelta -= aRm * aRm; - aDelta /= 2. * (aRb - aRm); - aDelta -= 0.5 * (aRb - aRm); - - aT1 = 2. * aRm * (aRm - aDelta); - aT2 = aT1; - } - aCriteria = ( aT1 > aT2) ? aT1 : aT2; - if ( aCriteria > 0 ) - aCriteria = sqrt( aCriteria ); - - if ( aCriteria > 0.5 * aTor1.MinorRadius() ) { - // too big zone -> drop to minimum - aCriteria = Precision::Confusion(); - } - - GeomAdaptor_Curve aC1( new Geom_Circle(aCircle1) ); - GeomAdaptor_Curve aC2( new Geom_Circle(aCircle2) ); - Extrema_ExtCC anExtrema(aC1, aC2, 0, 2. * M_PI, 0, 2. * M_PI, - Precision::PConfusion(), Precision::PConfusion()); - - if ( anExtrema.IsDone() ) { - - Standard_Integer i = 0; - for ( i = 1; i <= anExtrema.NbExt(); i++ ) { - if ( anExtrema.SquareDistance(i) > aCriteria * aCriteria ) - continue; - - Extrema_POnCurv P1, P2; - anExtrema.Points( i, P1, P2 ); - - Standard_Boolean bFoundResult = Standard_True; - gp_Pnt2d pr1, pr2; - - Standard_Integer surfit = 0; - for ( surfit = 0; surfit < 2; surfit++ ) { - GeomAPI_ProjectPointOnSurf& aProjector = - (surfit == 0) ? aContext->ProjPS(theFace1) : aContext->ProjPS(theFace2); - - gp_Pnt aP3d = (surfit == 0) ? P1.Value() : P2.Value(); - aProjector.Perform(aP3d); - - if(!aProjector.IsDone()) - bFoundResult = Standard_False; - else { - if(aProjector.LowerDistance() > aCriteria) { - bFoundResult = Standard_False; - } - else { - Standard_Real foundU = 0, foundV = 0; - aProjector.LowerDistanceParameters(foundU, foundV); - if ( surfit == 0 ) - pr1 = gp_Pnt2d( foundU, foundV ); - else - pr2 = gp_Pnt2d( foundU, foundV ); - } - } - } - if ( bFoundResult ) { - aSeqResultS1.Append( pr1 ); - aSeqResultS2.Append( pr2 ); - aSeqResultRad.Append( aCriteria ); - - // torus is u and v periodic - const Standard_Real twoPI = M_PI + M_PI; - Standard_Real arr1tmp[2] = {pr1.X(), pr1.Y()}; - Standard_Real arr2tmp[2] = {pr2.X(), pr2.Y()}; - - // iteration on period bounds - for ( Standard_Integer k1 = 0; k1 < 2; k1++ ) { - Standard_Real aBound = ( k1 == 0 ) ? 0 : twoPI; - Standard_Real aShift = ( k1 == 0 ) ? twoPI : -twoPI; - - // iteration on surfaces - for ( Standard_Integer k2 = 0; k2 < 2; k2++ ) { - Standard_Real* arr1 = ( k2 == 0 ) ? arr1tmp : arr2tmp; - Standard_Real* arr2 = ( k2 != 0 ) ? arr1tmp : arr2tmp; - TColgp_SequenceOfPnt2d& aSeqS1 = ( k2 == 0 ) ? aSeqResultS1 : aSeqResultS2; - TColgp_SequenceOfPnt2d& aSeqS2 = ( k2 != 0 ) ? aSeqResultS1 : aSeqResultS2; - - if (fabs(arr1[0] - aBound) < Precision::PConfusion()) { - aSeqS1.Append( gp_Pnt2d( arr1[0] + aShift, arr1[1] ) ); - aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) ); - aSeqResultRad.Append( aCriteria ); - } - if (fabs(arr1[1] - aBound) < Precision::PConfusion()) { - aSeqS1.Append( gp_Pnt2d( arr1[0], arr1[1] + aShift) ); - aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) ); - aSeqResultRad.Append( aCriteria ); - } - } - } // - } - } - } - } - aResult = aSeqResultRad.Length(); - - if ( aResult > 0 ) { - theResultOnS1 = new TColgp_HArray1OfPnt2d( 1, aResult ); - theResultOnS2 = new TColgp_HArray1OfPnt2d( 1, aResult ); - theResultRadius = new TColStd_HArray1OfReal( 1, aResult ); - - for ( Standard_Integer i = 1 ; i <= aResult; i++ ) { - theResultOnS1->SetValue( i, aSeqResultS1.Value(i) ); - theResultOnS2->SetValue( i, aSeqResultS2.Value(i) ); - theResultRadius->SetValue( i, aSeqResultRad.Value(i) ); - } - } - return aResult; -} - //======================================================================= //function : IsPointOnBoundary //purpose : static subfunction in DecompositionOfWLine @@ -456,27 +259,6 @@ Standard_Boolean IsPointOnBoundary(const Standard_Real theParameter, return bRet; } -//======================================================================= -//function : IsInsideTanZone -//purpose : Check if point is inside a radial tangent zone. -// static subfunction in DecompositionOfWLine and FindPoint -//======================================================================= -static -Standard_Boolean IsInsideTanZone(const gp_Pnt2d& thePoint, - const gp_Pnt2d& theTanZoneCenter, - const Standard_Real theZoneRadius, - Handle(GeomAdaptor_Surface) theGASurface) -{ - Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius ); - Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius ); - Standard_Real aRadiusSQR = ( aUResolution < aVResolution ) ? aUResolution : aVResolution; - aRadiusSQR *= aRadiusSQR; - if ( thePoint.SquareDistance( theTanZoneCenter ) <= aRadiusSQR ) - return Standard_True; - - return Standard_False; -} - //======================================================================= //function : AdjustByNeighbour //purpose : static subfunction in DecompositionOfWLine @@ -651,72 +433,6 @@ Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint, return Standard_False; } -//======================================================================= -//function : FindPoint -//purpose : Find point on the boundary of radial tangent zone -// static subfunction in DecompositionOfWLine -//======================================================================= -static -Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint, - const gp_Pnt2d& theLastPoint, - const Standard_Real theUmin, - const Standard_Real theUmax, - const Standard_Real theVmin, - const Standard_Real theVmax, - const gp_Pnt2d& theTanZoneCenter, - const Standard_Real theZoneRadius, - Handle(GeomAdaptor_Surface) theGASurface, - gp_Pnt2d& theNewPoint) { - theNewPoint = theLastPoint; - - if ( !IsInsideTanZone( theLastPoint, theTanZoneCenter, theZoneRadius, theGASurface) ) - return Standard_False; - - Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius ); - Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius ); - - Standard_Real aRadius = ( aUResolution < aVResolution ) ? aUResolution : aVResolution; - gp_Ax22d anAxis( theTanZoneCenter, gp_Dir2d(1, 0), gp_Dir2d(0, 1) ); - gp_Circ2d aCircle( anAxis, aRadius ); - - // - gp_Vec2d aDir( theLastPoint.XY() - theFirstPoint.XY() ); - Standard_Real aLength = aDir.Magnitude(); - if ( aLength <= gp::Resolution() ) - return Standard_False; - gp_Lin2d aLine( theFirstPoint, aDir ); - - // - Handle(Geom2d_Line) aCLine = new Geom2d_Line( aLine ); - Handle(Geom2d_TrimmedCurve) aC1 = new Geom2d_TrimmedCurve( aCLine, 0, aLength ); - Handle(Geom2d_Circle) aC2 = new Geom2d_Circle( aCircle ); - - Standard_Real aTol = aRadius * 0.001; - aTol = ( aTol < Precision::PConfusion() ) ? Precision::PConfusion() : aTol; - - Geom2dAPI_InterCurveCurve anIntersector; - anIntersector.Init( aC1, aC2, aTol ); - - if ( anIntersector.NbPoints() == 0 ) - return Standard_False; - - Standard_Boolean aFound = Standard_False; - Standard_Real aMinDist = aLength * aLength; - Standard_Integer i = 0; - for ( i = 1; i <= anIntersector.NbPoints(); i++ ) { - gp_Pnt2d aPInt = anIntersector.Point( i ); - if ( aPInt.SquareDistance( theFirstPoint ) < aMinDist ) { - if ( ( aPInt.X() >= theUmin ) && ( aPInt.X() <= theUmax ) && - ( aPInt.Y() >= theVmin ) && ( aPInt.Y() <= theVmax ) ) { - theNewPoint = aPInt; - aFound = Standard_True; - } - } - } - - return aFound; -} - //======================================================================= //function : DecompositionOfWLine //purpose : @@ -731,7 +447,6 @@ Standard_Boolean IntTools_WLineTool:: const Standard_Boolean theAvoidLConstructor, const Standard_Real theTol, IntPatch_SequenceOfLine& theNewLines, - Standard_Real& theReachedTol3d, const Handle(IntTools_Context)& aContext) { Standard_Boolean bRet, bAvoidLineConstructor; @@ -757,13 +472,7 @@ Standard_Boolean IntTools_WLineTool:: TColStd_Array1OfListOfInteger anArrayOfLines(1, aNbPnts); TColStd_Array1OfInteger anArrayOfLineType(1, aNbPnts); TColStd_ListOfInteger aListOfPointIndex; - - Handle(TColgp_HArray1OfPnt2d) aTanZoneS1; - Handle(TColgp_HArray1OfPnt2d) aTanZoneS2; - Handle(TColStd_HArray1OfReal) aTanZoneRadius; - Standard_Integer aNbZone = ComputeTangentZones( theSurface1, theSurface2, theFace1, theFace2, - aTanZoneS1, aTanZoneS2, aTanZoneRadius, aContext); - + // nblines=0; aTol=Precision::Confusion(); @@ -834,24 +543,6 @@ Standard_Boolean IntTools_WLineTool:: bIsCurrentPointOnBoundary = Standard_True; break; } - else { - // check if a point belong to a tangent zone. Begin - Standard_Integer zIt = 0; - for ( zIt = 1; zIt <= aNbZone; zIt++ ) { - gp_Pnt2d aPZone = (i == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt); - Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt); - - if ( IsInsideTanZone(gp_Pnt2d( U, V ), aPZone, aZoneRadius, aGASurface ) ) { - // set boundary flag to split the curve by a tangent zone - bIsPointOnBoundary = Standard_True; - bIsCurrentPointOnBoundary = Standard_True; - if ( theReachedTol3d < aZoneRadius ) { - theReachedTol3d = aZoneRadius; - } - break; - } - } - } }//for(j = 0; j < 2; j++) { if(bIsCurrentPointOnBoundary){ @@ -930,7 +621,7 @@ Standard_Boolean IntTools_WLineTool:: Standard_Integer nbboundaries = 0; Standard_Boolean bIsNearBoundary = Standard_False; - Standard_Integer aZoneIndex = 0; + //Standard_Integer aZoneIndex = 0; Standard_Integer bIsUBoundary = Standard_False; // use if nbboundaries == 1 Standard_Integer bIsFirstBoundary = Standard_False; // use if nbboundaries == 1 @@ -980,32 +671,6 @@ Standard_Boolean IntTools_WLineTool:: } } - // check if a point belong to a tangent zone. Begin - for ( Standard_Integer zIt = 1; zIt <= aNbZone; zIt++ ) { - gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt); - Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt); - - Standard_Integer aneighbourpointindex1 = (j == 0) ? iFirst : iLast; - const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1); - Standard_Real nU1, nV1; - - if(surfit == 0) - aNeighbourPoint.ParametersOnS1(nU1, nV1); - else - aNeighbourPoint.ParametersOnS2(nU1, nV1); - gp_Pnt2d ap1(nU1, nV1); - gp_Pnt2d ap2 = AdjustByNeighbour( ap1, gp_Pnt2d( U, V ), aGASurface ); - - - if ( IsInsideTanZone( ap2, aPZone, aZoneRadius, aGASurface ) ) { - aZoneIndex = zIt; - bIsNearBoundary = Standard_True; - if ( theReachedTol3d < aZoneRadius ) { - theReachedTol3d = aZoneRadius; - } - } - } - // check if a point belong to a tangent zone. End Standard_Boolean bComputeLineEnd = Standard_False; if(nbboundaries == 2) { @@ -1144,20 +809,7 @@ Standard_Boolean IntTools_WLineTool:: gp_Pnt2d ap1(nU1, nV1); gp_Pnt2d ap2; - - if ( aZoneIndex ) { - // exclude point from a tangent zone - anewpoint = AdjustByNeighbour( ap1, anewpoint, aGASurface ); - gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(aZoneIndex) : aTanZoneS2->Value(aZoneIndex); - Standard_Real aZoneRadius = aTanZoneRadius->Value(aZoneIndex); - - if ( FindPoint(ap1, anewpoint, umin, umax, vmin, vmax, - aPZone, aZoneRadius, aGASurface, ap2) ) { - anewpoint = ap2; - found = Standard_True; - } - } - else if ( aGASurface->IsUPeriodic() || aGASurface->IsVPeriodic() ) { + if ( aGASurface->IsUPeriodic() || aGASurface->IsVPeriodic() ) { // re-compute point near boundary if shifted on a period ap2 = AdjustByNeighbour( ap1, anewpoint, aGASurface ); diff --git a/src/IntTools/IntTools_WLineTool.hxx b/src/IntTools/IntTools_WLineTool.hxx index 89bb5a446a..a70fe4d503 100644 --- a/src/IntTools/IntTools_WLineTool.hxx +++ b/src/IntTools/IntTools_WLineTool.hxx @@ -46,7 +46,6 @@ public: const Standard_Boolean theAvoidLConstructor, const Standard_Real theTol, IntPatch_SequenceOfLine& theNewLines, - Standard_Real& theReachedTol3d, const Handle(IntTools_Context)& ); }; diff --git a/src/LocOpe/LocOpe_WiresOnShape.cxx b/src/LocOpe/LocOpe_WiresOnShape.cxx index d3fbe5ef25..3ccfa9f8dc 100644 --- a/src/LocOpe/LocOpe_WiresOnShape.cxx +++ b/src/LocOpe/LocOpe_WiresOnShape.cxx @@ -1295,6 +1295,10 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge, Standard_Real aTolV[2]; aTolV[0] =BRep_Tool::Tolerance(theVertices[0]); aTolV[1] =BRep_Tool::Tolerance(theVertices[1]); + Standard_Real ext = 16.; // = 4 * 4 - to avoid creating microedges, area around vertices is increased + // up to 4 vertex tolerance. Such approach is usual for other topological + // algorithms, for example, Boolean Operations. + Standard_Real aTolVExt[2] = { ext * aTolV[0] * aTolV[0], ext * aTolV[1] * aTolV[1] }; BRepAdaptor_Curve2d thePCurve(theEdge, theFace); Bnd_Box2d theBox; @@ -1304,6 +1308,9 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge, Standard_Real aFpar, aLpar; const Handle(Geom_Curve)& theCurve = BRep_Tool::Curve(theEdge, thePar[0], thePar[1]); GeomAdaptor_Curve theGAcurve(theCurve, thePar[0], thePar[1]); + Standard_Real aTolV2d[2] = { theGAcurve.Resolution(aTolV[0]), theGAcurve.Resolution(aTolV[1]) }; + aTolV2d[0] = Max(aTolV2d[0], Precision::PConfusion()); + aTolV2d[1] = Max(aTolV2d[1], Precision::PConfusion()); Standard_Real aDistMax = Precision::Confusion() * Precision::Confusion(); TopExp_Explorer Explo(theFace, TopAbs_EDGE); for (; Explo.More(); Explo.Next()) @@ -1330,6 +1337,23 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge, isOverlapped = Standard_True; return; } + // Check extremity distances + Standard_Real dists[4]; + gp_Pnt aP11, aP12, aP21, aP22; + anExtrema.TrimmedSquareDistances(dists[0], dists[1], dists[2], dists[3], + aP11, aP12, aP21, aP22); + for (i = 0; i < 4; ++i) + { + if (i < 2) + j = 0; + else + j = 1; + if (dists[i] < aTolVExt[j] / ext) + { + return; + } + } + for (i = 1; i <= aNbExt; i++) { Standard_Real aDist = anExtrema.SquareDistance(i); @@ -1342,7 +1366,7 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge, Standard_Real anIntPar = aPOnC2.Parameter(); for (j = 0; j < 2; j++) //try to find intersection on an extremity of "theEdge" { - if (Abs(theIntPar - thePar[j]) <= Precision::PConfusion()) + if (Abs(theIntPar - thePar[j]) <= aTolV2d[j]) break; } //intersection found in the middle of the edge @@ -1351,10 +1375,10 @@ void FindInternalIntersections(const TopoDS_Edge& theEdge, gp_Pnt aPoint = aCurve->Value(anIntPar); gp_Pnt aPointInt = theCurve->Value(theIntPar); - if (aPointInt.SquareDistance(thePnt[0]) > aTolV[0] * aTolV[0] && - aPointInt.SquareDistance(thePnt[1]) > aTolV[1] * aTolV[1] && - aPoint.SquareDistance(thePnt[0]) > aTolV[0] * aTolV[0] && - aPoint.SquareDistance(thePnt[1]) > aTolV[1] * aTolV[1]) + if (aPointInt.SquareDistance(thePnt[0]) > aTolVExt[0] && + aPointInt.SquareDistance(thePnt[1]) > aTolVExt[1] && + aPoint.SquareDistance(thePnt[0]) > aTolVExt[0] && + aPoint.SquareDistance(thePnt[1]) > aTolVExt[1]) { SplitPars.Append(theIntPar); if( aDist > aDistMax) diff --git a/tests/bugs/modalg_5/bug23706_10 b/tests/bugs/modalg_5/bug23706_10 index d80c11f0e4..6a5d69e19d 100755 --- a/tests/bugs/modalg_5/bug23706_10 +++ b/tests/bugs/modalg_5/bug23706_10 @@ -13,7 +13,7 @@ set info [extrema r3 r4] if {[regexp "ext_1" $info]} { set dist [lindex [length ext_1] end] - if { $dist > 4.0e-13 } { + if { $dist > 5.0e-11 } { puts "Error: Extrema distance is too big" } } else { diff --git a/tests/bugs/modalg_5/bug23706_11 b/tests/bugs/modalg_5/bug23706_11 index b1afdfe5a0..9df8d5c5f4 100755 --- a/tests/bugs/modalg_5/bug23706_11 +++ b/tests/bugs/modalg_5/bug23706_11 @@ -10,8 +10,11 @@ bsplinecurve r1 2 5 1 3 2 1 3 1 4 1 5 3 2 5 3 1 3 7 3 1 4 8 3 1 4 8 3 1 5 9 3 1 bsplinecurve r2 2 5 2 3 2.5 1 3 1 3.5 1 4 3 -1 2 3 1 1 11 3 1 3 9 3 1 3 9 3 1 3 9 3 1 5 7 3 1 7 4 3 1 set info [extrema r1 r2] -if { [llength $info] != 1 } { - puts "Error : Extrema is wrong" +if {[regexp "ext_1" $info]} { + set dist [lindex [length ext_1] end] + if { $dist > 1.0e-8 } { + puts "Error: Extrema distance is too big" + } } else { - puts "OK: Extrema is valid" + puts "Error: Extrema is not found" } diff --git a/tests/bugs/modalg_6/bug27665 b/tests/bugs/modalg_6/bug27665 index 85fe8294f2..7fa5fc9be7 100644 --- a/tests/bugs/modalg_6/bug27665 +++ b/tests/bugs/modalg_6/bug27665 @@ -5,7 +5,7 @@ puts "" ################################################# # BrepExrtrema_DistShapeShape bad performance on OCCT 6.7.0 ################################################# -cpulimit 100 +cpulimit 500 restore [locate_data_file bug27665_wircmpd.brep] w explode w diff --git a/tests/lowalgos/extcc/bug29858_03 b/tests/lowalgos/extcc/bug29858_03 index 7058727e53..51702e229d 100644 --- a/tests/lowalgos/extcc/bug29858_03 +++ b/tests/lowalgos/extcc/bug29858_03 @@ -18,11 +18,12 @@ mkcurve c2 e2 set info [extrema c1 c2] # Check result -regexp {Extrema 1 is point : +([-0-9.+eE]+) +([-0-9.+eE]+) +([-0-9.+eE]+)} $info full x y z -# Point check -set good_x 0.0 -set good_y 0.070710562195021642 -set good_z -0.65305318986891325 -checkreal "Intersection point x:" ${x} ${good_x} 0.01 0.01 -checkreal "Intersection point y:" ${y} ${good_y} 0.01 0.01 -checkreal "Intersection point z:" ${z} ${good_z} 0.01 0.01 \ No newline at end of file +if {[regexp "ext_1" $info]} { + set dist [lindex [length ext_1] end] + if { $dist > 1.0e-10 } { + puts "Error: Extrema distance is too big" + } +} else { + puts "Error: Extrema is not found" +} + diff --git a/tests/lowalgos/extcc/bug32882 b/tests/lowalgos/extcc/bug32882 new file mode 100644 index 0000000000..9e91bbedf6 --- /dev/null +++ b/tests/lowalgos/extcc/bug32882 @@ -0,0 +1,36 @@ +puts "========" +puts "OCC32882: Modeling Data - Extrema curve/curve cannot find all solutions" +puts "========" +puts "" + + +# Read input +restore [locate_data_file bug32882.brep] cc +explode cc + +# Extract geometry from topology +mkcurve c1 cc_1 +mkcurve c2 cc_2 +mkcurve c3 cc_3 + +# Run extrema c1/c3 +set info [extrema c1 c3] + +# Check number of solution +if { [llength $info] != 3 } { + puts "Error: Invalid extrema number in extrema c1-c3 output" +} + +# Check result +checklength ext_1 -l 2.929642751e-14 -eps .01 +checklength ext_2 -l 3.480934286e-14 -eps .01 +checklength ext_3 -l 3.177643716e-14 -eps .01 + +# Run extrema c3/c2 +set info [extrema c3 c2] +# Check number of solutions +if { [llength $info] != 2 } { + puts "Error: Invalid extrema number in extrema c3-c2 output" +} +checklength ext_1 -l 5.684341886e-14 -eps .01 +checklength ext_2 -l 2.929642751e-14 -eps .01