diff --git a/src/BOPAlgo/BOPAlgo_PaveFiller_5.cxx b/src/BOPAlgo/BOPAlgo_PaveFiller_5.cxx index ae20670cc5..6c4ab28a1a 100644 --- a/src/BOPAlgo/BOPAlgo_PaveFiller_5.cxx +++ b/src/BOPAlgo/BOPAlgo_PaveFiller_5.cxx @@ -424,8 +424,17 @@ void BOPAlgo_PaveFiller::PerformEF() const TopoDS_Vertex& aV = TopoDS::Vertex(myDS->Shape(nV[j])); const gp_Pnt aP = BRep_Tool::Pnt(aV); Standard_Real aDistPP = aP.Distance(aPnew); - UpdateVertex(nV[j], aDistPP); - myVertsToAvoidExtension.Add(nV[j]); + Standard_Real aTol = BRep_Tool::Tolerance(aV); + Standard_Real aMaxDist = 1.e4 * aTol; + if (aTol < .01) + { + aMaxDist = Min(aMaxDist, 0.1); + } + if (aDistPP < aMaxDist) + { + UpdateVertex(nV[j], aDistPP); + myVertsToAvoidExtension.Add(nV[j]); + } } } continue; diff --git a/src/ElSLib/ElSLib.cxx b/src/ElSLib/ElSLib.cxx index e2a815c2cd..86d7f591f5 100644 --- a/src/ElSLib/ElSLib.cxx +++ b/src/ElSLib/ElSLib.cxx @@ -1480,10 +1480,19 @@ void ElSLib::TorusParameters(const gp_Ax3& Pos, Standard_Real cosu = cos(U); Standard_Real sinu = sin(U); gp_Dir dx(cosu,sinu,0.); - gp_Dir dP(x - MajorRadius * cosu, - y - MajorRadius * sinu, - z); - V = dx.AngleWithRef(dP,dx^gp::DZ()); + gp_XYZ dPV(x - MajorRadius * cosu, + y - MajorRadius * sinu, + z); + Standard_Real aMag = dPV.Modulus(); + if (aMag <= gp::Resolution()) + { + V = 0.; + } + else + { + gp_Dir dP(dPV); + V = dx.AngleWithRef(dP, dx^gp::DZ()); + } if (V < -1.e-16) V += PIPI; else if (V < 0) V = 0; } diff --git a/src/Extrema/Extrema_ExtCS.cxx b/src/Extrema/Extrema_ExtCS.cxx index 68cc4138d2..0f980966e2 100644 --- a/src/Extrema/Extrema_ExtCS.cxx +++ b/src/Extrema/Extrema_ExtCS.cxx @@ -39,6 +39,7 @@ #include #include #include +#include Extrema_ExtCS::Extrema_ExtCS() { @@ -177,6 +178,27 @@ void Extrema_ExtCS::Perform(const Adaptor3d_Curve& C, if (myS->IsVPeriodic()) NbV = 13; + if (clast - cfirst <= Precision::Confusion()) + { + Standard_Real aCPar = (cfirst + clast) / 2.; + gp_Pnt aPm = C.Value(aCPar); + Extrema_ExtPS anExtPS(aPm, *myS, ufirst, ulast, + vfirst, vlast, mytolS, mytolS, Extrema_ExtFlag_MIN); + myDone = anExtPS.IsDone(); + if (myDone) { + Standard_Integer NbExt = anExtPS.NbExt(); + Standard_Real T = aCPar, U, V; + Extrema_POnCurv PC; + Extrema_POnSurf PS; + for (i = 1; i <= NbExt; i++) { + PS = anExtPS.Point(i); + PS.Parameter(U, V); + AddSolution(C, T, U, V, PC.Value(), PS.Value(), anExtPS.SquareDistance(i)); + } + } + return; + } + Extrema_GenExtCS Ext(C, *myS, NbT, NbU, NbV, cfirst, clast, ufirst, ulast, vfirst, vlast, mytolC, mytolS); diff --git a/src/Extrema/Extrema_GenExtCS.cxx b/src/Extrema/Extrema_GenExtCS.cxx index e83890b525..bae4a5dcde 100644 --- a/src/Extrema/Extrema_GenExtCS.cxx +++ b/src/Extrema/Extrema_GenExtCS.cxx @@ -15,19 +15,18 @@ // commercial license or contractual agreement. +#include #include #include #include -#include -#include -#include +#include +#include #include +#include +#include #include #include #include -#include -#include -#include #include #include #include @@ -36,14 +35,32 @@ #include #include #include -#include -#include +#include #include +#include const Standard_Real MaxParamVal = 1.0e+10; const Standard_Real aBorderDivisor = 1.0e+4; const Standard_Real HyperbolaLimit = 23.; //ln(MaxParamVal) +static Standard_Boolean IsQuadric(const GeomAbs_SurfaceType theSType) +{ + if (theSType == GeomAbs_Plane) return Standard_True; + if (theSType == GeomAbs_Cylinder) return Standard_True; + if (theSType == GeomAbs_Cone) return Standard_True; + if (theSType == GeomAbs_Sphere) return Standard_True; + if (theSType == GeomAbs_Torus) return Standard_True; + return Standard_False; +} +static Standard_Boolean IsConic(const GeomAbs_CurveType theCType) +{ + if (theCType == GeomAbs_Line) return Standard_True; + if (theCType == GeomAbs_Circle) return Standard_True; + if (theCType == GeomAbs_Ellipse) return Standard_True; + if (theCType == GeomAbs_Hyperbola) return Standard_True; + if (theCType == GeomAbs_Parabola) return Standard_True; + return Standard_False; +} // restrict maximal parameter on hyperbola to avoid FPE static Standard_Real GetCurvMaxParamVal (const Adaptor3d_Curve& theC) { @@ -255,60 +272,106 @@ void Extrema_GenExtCS::Perform (const Adaptor3d_Curve& C, mytmin = -aCMaxVal; } // - math_Vector Tol(1, 3); + Standard_Integer aNbVar = 3; + GeomAbs_SurfaceType aSType = myS->GetType(); + if (IsQuadric(aSType)) + { + aNbVar = 1; + } + else + { + GeomAbs_CurveType aCType = C.GetType(); + if (IsConic(aCType)) + { + aNbVar = 2; + } + } + + math_Vector Tol(1, 3), TUV(1, 3), TUVinf(1, 3), TUVsup(1, 3); + // Tol(1) = mytol1; Tol(2) = mytol2; Tol(3) = mytol2; - math_Vector TUV(1,3), TUVinf(1,3), TUVsup(1,3); + // TUVinf(1) = mytmin; TUVinf(2) = trimumin; TUVinf(3) = trimvmin; - + // TUVsup(1) = mytsup; TUVsup(2) = trimusup; TUVsup(3) = trimvsup; - + // // Number of particles used in PSO algorithm (particle swarm optimization). const Standard_Integer aNbParticles = 48; + // + if (aNbVar == 3) + { + GlobMinGenCS(C, aNbParticles, TUVinf, TUVsup, TUV); + } + else if (aNbVar == 2) + { + GlobMinConicS(C, aNbParticles, TUVinf, TUVsup, TUV); + } + else + { + GlobMinCQuadric(C, aNbParticles, TUVinf, TUVsup, TUV); + } - math_PSOParticlesPool aParticles(aNbParticles, 3); + // Find min approximation + math_FunctionSetRoot anA(myF, Tol); + anA.Perform(myF, TUV, TUVinf, TUVsup); - math_Vector aMinTUV(1,3); - aMinTUV = TUVinf + (TUVsup - TUVinf) / aBorderDivisor; + myDone = Standard_True; +} +//======================================================================= +//function : GlobMinGenCS +//purpose : +//======================================================================= +void Extrema_GenExtCS::GlobMinGenCS(const Adaptor3d_Curve& theC, + const Standard_Integer theNbParticles, + const math_Vector& theTUVinf, + const math_Vector& theTUVsup, + math_Vector& theTUV) +{ + math_PSOParticlesPool aParticles(theNbParticles, 3); - math_Vector aMaxTUV(1,3); - aMaxTUV = TUVsup - (TUVsup - TUVinf) / aBorderDivisor; + math_Vector aMinTUV(1, 3); + aMinTUV = theTUVinf + (theTUVsup - theTUVinf) / aBorderDivisor; + + math_Vector aMaxTUV(1, 3); + aMaxTUV = theTUVsup - (theTUVsup - theTUVinf) / aBorderDivisor; Standard_Real aStepCU = (aMaxTUV(1) - aMinTUV(1)) / mytsample; Standard_Real aStepSU = (aMaxTUV(2) - aMinTUV(2)) / myusample; Standard_Real aStepSV = (aMaxTUV(3) - aMinTUV(3)) / myvsample; // Correct number of curve samples in case of low resolution + Standard_Integer aNewCsample = mytsample; Standard_Real aScaleFactor = 5.0; - Standard_Real aResolutionCU = aStepCU / C.Resolution (1.0); + Standard_Real aResolutionCU = aStepCU / theC.Resolution(1.0); - Standard_Real aMinResolution = aScaleFactor * Min (aResolutionCU, - Min (aStepSU / myS->UResolution (1.0), aStepSV / myS->VResolution (1.0))); + Standard_Real aMinResolution = aScaleFactor * Min(aResolutionCU, + Min(aStepSU / myS->UResolution(1.0), aStepSV / myS->VResolution(1.0))); - if (aMinResolution > Epsilon (1.0)) + if (aMinResolution > Epsilon(1.0)) { if (aResolutionCU > aMinResolution) { const Standard_Integer aMaxNbNodes = 50; - mytsample = Min(aMaxNbNodes, + aNewCsample = Min(aMaxNbNodes, RealToInt(mytsample * aResolutionCU / aMinResolution)); - aStepCU = (aMaxTUV(1) - aMinTUV(1)) / mytsample; + aStepCU = (aMaxTUV(1) - aMinTUV(1)) / aNewCsample; } } // Pre-compute curve sample points. - TColgp_HArray1OfPnt aCurvPnts (0, mytsample); + TColgp_Array1OfPnt aCurvPnts(0, aNewCsample); Standard_Real aCU1 = aMinTUV(1); - for (Standard_Integer aCUI = 0; aCUI <= mytsample; aCUI++, aCU1 += aStepCU) - aCurvPnts.SetValue (aCUI, C.Value (aCU1)); + for (Standard_Integer aCUI = 0; aCUI <= aNewCsample; aCUI++, aCU1 += aStepCU) + aCurvPnts.SetValue(aCUI, theC.Value(aCU1)); PSO_Particle* aParticle = aParticles.GetWorstParticle(); // Select specified number of particles from pre-computed set of samples @@ -319,9 +382,9 @@ void Extrema_GenExtCS::Perform (const Adaptor3d_Curve& C, for (Standard_Integer aSVI = 0; aSVI <= myvsample; aSVI++, aSV += aStepSV) { Standard_Real aCU2 = aMinTUV(1); - for (Standard_Integer aCUI = 0; aCUI <= mytsample; aCUI++, aCU2 += aStepCU) + for (Standard_Integer aCUI = 0; aCUI <= aNewCsample; aCUI++, aCU2 += aStepCU) { - Standard_Real aSqDist = mySurfPnts->Value(aSUI, aSVI).SquareDistance(aCurvPnts.Value(aCUI)); + Standard_Real aSqDist = mySurfPnts->Value(aSUI, aSVI).SquareDistance(aCurvPnts.Value(aCUI)); if (aSqDist < aParticle->Distance) { @@ -333,7 +396,7 @@ void Extrema_GenExtCS::Perform (const Adaptor3d_Curve& C, aParticle->BestPosition[1] = aSU; aParticle->BestPosition[2] = aSV; - aParticle->Distance = aSqDist; + aParticle->Distance = aSqDist; aParticle->BestDistance = aSqDist; aParticle = aParticles.GetWorstParticle(); @@ -342,21 +405,225 @@ void Extrema_GenExtCS::Perform (const Adaptor3d_Curve& C, } } - math_Vector aStep(1,3); + math_Vector aStep(1, 3); aStep(1) = aStepCU; aStep(2) = aStepSU; aStep(3) = aStepSV; // Find min approximation Standard_Real aValue; - Extrema_GlobOptFuncCS aFunc(&C, myS); - math_PSO aPSO(&aFunc, TUVinf, TUVsup, aStep); - aPSO.Perform(aParticles, aNbParticles, aValue, TUV); + Extrema_GlobOptFuncCS aFunc(&theC, myS); + math_PSO aPSO(&aFunc, theTUVinf, theTUVsup, aStep); + aPSO.Perform(aParticles, theNbParticles, aValue, theTUV); - math_FunctionSetRoot anA(myF, Tol); - anA.Perform(myF, TUV, TUVinf, TUVsup); +} +//======================================================================= +//function : GlobMinConicS +//purpose : +//======================================================================= +void Extrema_GenExtCS::GlobMinConicS(const Adaptor3d_Curve& theC, + const Standard_Integer theNbParticles, + const math_Vector& theTUVinf, + const math_Vector& theTUVsup, + math_Vector& theTUV) +{ + Standard_Integer aNbVar = 2; + math_Vector anUVinf(1, aNbVar), anUVsup(1, aNbVar), anUV(1, aNbVar); + Standard_Integer i; + for (i = 1; i <= aNbVar; ++i) + { + anUVinf(i) = theTUVinf(i + 1); + anUVsup(i) = theTUVsup(i + 1); + } + // + math_PSOParticlesPool aParticles(theNbParticles, aNbVar); + + math_Vector aMinUV(1, aNbVar); + aMinUV = anUVinf + (anUVsup - anUVinf) / aBorderDivisor; + + math_Vector aMaxUV(1, aNbVar); + aMaxUV = anUVsup - (anUVsup - anUVinf) / aBorderDivisor; + + //Increase numbers of UV samples to improve searching global minimum + Standard_Integer anAddsample = Max(mytsample / 2, 3); + Standard_Integer anUsample = myusample + anAddsample; + Standard_Integer aVsample = myvsample + anAddsample; + // + Standard_Real aStepSU = (aMaxUV(1) - aMinUV(1)) / anUsample; + Standard_Real aStepSV = (aMaxUV(2) - aMinUV(2)) / aVsample; + // + Extrema_GlobOptFuncConicS aFunc(myS, anUVinf(1), anUVsup(1), anUVinf(2), anUVsup(2)); + aFunc.LoadConic(&theC, theTUVinf(1), theTUVsup(1)); + + + PSO_Particle* aParticle = aParticles.GetWorstParticle(); + // Select specified number of particles from pre-computed set of samples + Standard_Real aSU = aMinUV(1); + + for (Standard_Integer aSUI = 0; aSUI <= anUsample; aSUI++, aSU += aStepSU) + { + anUV(1) = aSU; + Standard_Real aSV = aMinUV(2); + for (Standard_Integer aSVI = 0; aSVI <= aVsample; aSVI++, aSV += aStepSV) + { + anUV(2) = aSV; + Standard_Real aSqDist; + if (!aFunc.Value(anUV, aSqDist)) + { + aSqDist = Precision::Infinite(); + } + + if (aSqDist < aParticle->Distance) + { + aParticle->Position[0] = aSU; + aParticle->Position[1] = aSV; + + aParticle->BestPosition[0] = aSU; + aParticle->BestPosition[1] = aSV; + + aParticle->Distance = aSqDist; + aParticle->BestDistance = aSqDist; + + aParticle = aParticles.GetWorstParticle(); + } + } + } + + math_Vector aStep(1, aNbVar); + aStep(1) = aStepSU; + aStep(2) = aStepSV; + + // Find min approximation + Standard_Real aValue; + math_PSO aPSO(&aFunc, anUVinf, anUVsup, aStep); + aPSO.Perform(aParticles, theNbParticles, aValue, anUV); + // + Standard_Real aCT = aFunc.ConicParameter(anUV); + if (theC.IsPeriodic()) + { + if (aCT < theTUVinf(1) - Precision::PConfusion() || aCT > theTUVsup(1) + Precision::PConfusion()) + { + aCT = ElCLib::InPeriod(aCT, theTUVinf(1), theTUVinf(1) + 2. * M_PI); + } + } + theTUV(1) = aCT; + theTUV(2) = anUV(1); + theTUV(3) = anUV(2); + +} +//======================================================================= +//function : GlobMinCQuadric +//purpose : +//======================================================================= +void Extrema_GenExtCS::GlobMinCQuadric(const Adaptor3d_Curve& theC, + const Standard_Integer theNbParticles, + const math_Vector& theTUVinf, + const math_Vector& theTUVsup, + math_Vector& theTUV) +{ + Standard_Integer aNbVar = 1; + math_Vector aTinf(1, aNbVar), aTsup(1, aNbVar), aT(1, aNbVar); + aTinf(1) = theTUVinf(1); + aTsup(1) = theTUVsup(1); + // + math_PSOParticlesPool aParticles(theNbParticles, aNbVar); + + math_Vector aMinT(1, aNbVar); + aMinT = aTinf + (aTsup - aTinf) / aBorderDivisor; + + math_Vector aMaxT(1, aNbVar); + aMaxT = aTsup - (aTsup - aTinf) / aBorderDivisor; + // + + //Increase numbers of curve samples to improve searching global minimum + //because dimension of optimisation task is redused + const Standard_Integer aMaxNbNodes = 50; + Standard_Integer aNewCsample = mytsample; + Standard_Integer anAddsample = Max(myusample / 2, 3); + aNewCsample += anAddsample; + aNewCsample = Min(aNewCsample, aMaxNbNodes); + // + // Correct number of curve samples in case of low resolution + Standard_Real aStepCT = (aMaxT(1) - aMinT(1)) / aNewCsample; + Standard_Real aStepSU = (theTUVsup(2) - theTUVinf(2)) / myusample; + Standard_Real aStepSV = (theTUVsup(3) - theTUVinf(3)) / myvsample; + Standard_Real aScaleFactor = 5.0; + Standard_Real aResolutionCU = aStepCT / theC.Resolution(1.0); + + Standard_Real aMinResolution = aScaleFactor * Min(aResolutionCU, + Min(aStepSU / myS->UResolution(1.0), aStepSV / myS->VResolution(1.0))); + + if (aMinResolution > Epsilon(1.0)) + { + if (aResolutionCU > aMinResolution) + { + + aNewCsample = Min(aMaxNbNodes, + RealToInt(aNewCsample * aResolutionCU / aMinResolution)); + + aStepCT = (aMaxT(1) - aMinT(1)) / aNewCsample; + } + } + + // + Extrema_GlobOptFuncCQuadric aFunc(&theC, aTinf(1), aTsup(1)); + aFunc.LoadQuad(myS, theTUVinf(2), theTUVsup(2), theTUVinf(3), theTUVsup(3)); + + PSO_Particle* aParticle = aParticles.GetWorstParticle(); + // Select specified number of particles from pre-computed set of samples + Standard_Real aCT = aMinT(1); + for (Standard_Integer aCUI = 0; aCUI <= aNewCsample; aCUI++, aCT += aStepCT) + { + aT(1) = aCT; + Standard_Real aSqDist; + if (!aFunc.Value(aT, aSqDist)) + { + aSqDist = Precision::Infinite(); + } + + if (aSqDist < aParticle->Distance) + { + aParticle->Position[0] = aCT; + + aParticle->BestPosition[0] = aCT; + + aParticle->Distance = aSqDist; + aParticle->BestDistance = aSqDist; + + aParticle = aParticles.GetWorstParticle(); + } + } + // + math_Vector aStep(1, aNbVar); + aStep(1) = aStepCT; + + // Find min approximation + Standard_Real aValue; + math_PSO aPSO(&aFunc, aTinf, aTsup, aStep); + aPSO.Perform(aParticles, theNbParticles, aValue, aT); + // + math_Vector anUV(1, 2); + aFunc.QuadricParameters(aT, anUV); + if (myS->IsUPeriodic()) + { + if (anUV(1) < theTUVinf(2) - Precision::PConfusion() || anUV(1)> theTUVsup(2) + Precision::PConfusion()) + { + anUV(1) = ElCLib::InPeriod(anUV(1), theTUVinf(2), theTUVinf(2) + 2. * M_PI); + } + } + // + if (myS->IsVPeriodic()) + { + if (anUV(2) < theTUVinf(3) - Precision::PConfusion() || anUV(2)> theTUVsup(3) + Precision::PConfusion()) + { + anUV(2) = ElCLib::InPeriod(anUV(2), theTUVinf(3), theTUVinf(3) + 2. * M_PI); + } + } + // + theTUV(1) = aT(1); + theTUV(2) = anUV(1); + theTUV(3) = anUV(2); - myDone = Standard_True; } //======================================================================= diff --git a/src/Extrema/Extrema_GenExtCS.hxx b/src/Extrema/Extrema_GenExtCS.hxx index d2ca3d5432..7f64f7cffd 100644 --- a/src/Extrema/Extrema_GenExtCS.hxx +++ b/src/Extrema/Extrema_GenExtCS.hxx @@ -28,6 +28,7 @@ #include #include #include + class StdFail_NotDone; class Standard_OutOfRange; class Standard_TypeMismatch; @@ -112,6 +113,24 @@ private: Standard_EXPORT Adaptor3d_SurfacePtr BidonSurface() const; + Standard_EXPORT void GlobMinGenCS(const Adaptor3d_Curve& theC, + const Standard_Integer theNbParticles, + const math_Vector& theTUVinf, + const math_Vector& theTUVsup, + math_Vector& theTUV); + + Standard_EXPORT void GlobMinConicS(const Adaptor3d_Curve& theC, + const Standard_Integer theNbParticles, + const math_Vector& theTUVinf, + const math_Vector& theTUVsup, + math_Vector& theTUV); + + Standard_EXPORT void GlobMinCQuadric(const Adaptor3d_Curve& theC, + const Standard_Integer theNbParticles, + const math_Vector& theTUVinf, + const math_Vector& theTUVsup, + math_Vector& theTUV); + Standard_Boolean myDone; Standard_Boolean myInit; diff --git a/src/Extrema/Extrema_GlobOptFuncCQuadric.cxx b/src/Extrema/Extrema_GlobOptFuncCQuadric.cxx new file mode 100644 index 0000000000..26250f8474 --- /dev/null +++ b/src/Extrema/Extrema_GlobOptFuncCQuadric.cxx @@ -0,0 +1,296 @@ +// Copyright (c) 2020 OPEN CASCADE SAS +// +// This file is part of Open CASCADE Technology software library. +// +// This library is free software; you can redistribute it and/or modify it under +// the terms of the GNU Lesser General Public License version 2.1 as published +// by the Free Software Foundation, with special exception defined in the file +// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT +// distribution for complete text of the license and disclaimer of any warranty. +// +// Alternatively, this file may be used under the terms of Open CASCADE +// commercial license or contractual agreement + +#include + +#include +#include +#include + + +//======================================================================= +//function : value +//purpose : +//======================================================================= +void Extrema_GlobOptFuncCQuadric::value(Standard_Real ct, + Standard_Real &F) +{ + Standard_Real u, v; + // + gp_Pnt aCP = myC->Value(ct); + switch (mySType) + { + case GeomAbs_Plane: + ElSLib::Parameters(myPln, aCP, u, v); + break; + case GeomAbs_Cylinder: + ElSLib::Parameters(myCylinder, aCP, u, v); + break; + case GeomAbs_Cone: + ElSLib::Parameters(myCone, aCP, u, v); + break; + case GeomAbs_Sphere: + ElSLib::Parameters(mySphere, aCP, u, v); + break; + case GeomAbs_Torus: + ElSLib::Parameters(myTorus, aCP, u, v); + break; + default: + F = Precision::Infinite(); + return; + } + // + if (mySType != GeomAbs_Plane) + { + if (myUl > 2. * M_PI + Precision::PConfusion()) + { + u += 2. * M_PI; + } + } + if (mySType == GeomAbs_Torus) + { + if (myVl > 2. * M_PI + Precision::PConfusion()) + { + v += 2. * M_PI; + } + } + + F = RealLast(); + if (u >= myUf && u <= myUl && v >= myVf && v <= myVl) + { + gp_Pnt aPS = myS->Value(u, v); + F = Min(F, aCP.SquareDistance(aPS)); + } + Standard_Integer i; + for (i = 0; i < 4; ++i) + { + F = Min(F, aCP.SquareDistance(myPTrim[i])); + } +} + + +//======================================================================= +//function : checkInputData +//purpose : +//======================================================================= +Standard_Boolean Extrema_GlobOptFuncCQuadric::checkInputData(const math_Vector &X, + Standard_Real &ct) +{ + ct = X(X.Lower()); + + if (ct < myTf || ct > myTl ) + { + return Standard_False; + } + return Standard_True; +} + +//======================================================================= +//function : Extrema_GlobOptFuncCQuadric +//purpose : Constructor +//======================================================================= +Extrema_GlobOptFuncCQuadric::Extrema_GlobOptFuncCQuadric(const Adaptor3d_Curve *C, + const Adaptor3d_Surface *S) +: myC(C) +{ + myTf = myC->FirstParameter(); + myTl = myC->LastParameter(); + Standard_Real anUf = S->FirstUParameter(), anUl = S->LastUParameter(); + Standard_Real aVf = S->FirstVParameter(), aVl = S->LastVParameter(); + LoadQuad(S, anUf, anUl, aVf, aVl); +} +//======================================================================= +//function : Extrema_GlobOptFuncCQuadric +//purpose : Constructor +//======================================================================= +Extrema_GlobOptFuncCQuadric::Extrema_GlobOptFuncCQuadric(const Adaptor3d_Curve *C) + : myC(C) +{ + myTf = myC->FirstParameter(); + myTl = myC->LastParameter(); +} +//======================================================================= +//function : Extrema_GlobOptFuncCQuadric +//purpose : Constructor +//======================================================================= +Extrema_GlobOptFuncCQuadric::Extrema_GlobOptFuncCQuadric(const Adaptor3d_Curve *C, + const Standard_Real theTf, const Standard_Real theTl) + : myC(C), myTf(theTf), myTl(theTl) +{ +} + +//======================================================================= +//function : LoadQuad +//purpose : +//======================================================================= +void Extrema_GlobOptFuncCQuadric::LoadQuad( const Adaptor3d_Surface *S, + const Standard_Real theUf, const Standard_Real theUl, + const Standard_Real theVf, const Standard_Real theVl) +{ + myS = S; + myUf = theUf; + myUl = theUl; + myVf = theVf; + myVl = theVl; + // + if (myS->IsUPeriodic()) + { + Standard_Real aTMax = 2. * M_PI + Precision::PConfusion(); + if (myUf > aTMax || myUf < -Precision::PConfusion() || + Abs(myUl - myUf) > aTMax) + { + ElCLib::AdjustPeriodic(0., 2. * M_PI, + Min(Abs(myUl - myUf) / 2, Precision::PConfusion()), + myUf, myUl); + } + } + if (myS->IsVPeriodic()) + { + Standard_Real aTMax = 2. * M_PI + Precision::PConfusion(); + if (myVf > aTMax || myVf < -Precision::PConfusion() || + Abs(myVl - myVf) > aTMax) + { + ElCLib::AdjustPeriodic(0., 2. * M_PI, + Min(Abs(myVl - myVf) / 2, Precision::PConfusion()), + myVf, myVl); + } + } + myPTrim[0] = myS->Value(myUf, myVf); + myPTrim[1] = myS->Value(myUl, myVf); + myPTrim[2] = myS->Value(myUl, myVl); + myPTrim[3] = myS->Value(myUf, myVl); + mySType = S->GetType(); + switch (mySType) + { + case GeomAbs_Plane: + myPln = myS->Plane(); + break; + case GeomAbs_Cylinder: + myCylinder = myS->Cylinder(); + break; + case GeomAbs_Cone: + myCone = myS->Cone(); + break; + case GeomAbs_Sphere: + mySphere = myS->Sphere(); + break; + case GeomAbs_Torus: + myTorus = myS->Torus(); + break; + default: + break; + } + +} + +//======================================================================= +//function : NbVariables +//purpose : +//======================================================================= +Standard_Integer Extrema_GlobOptFuncCQuadric::NbVariables() const +{ + return 1; +} + +//======================================================================= +//function : Value +//purpose : +//======================================================================= +Standard_Boolean Extrema_GlobOptFuncCQuadric::Value(const math_Vector &X, + Standard_Real &F) +{ + Standard_Real ct; + if (!checkInputData(X, ct)) + return Standard_False; + + value(ct, F); + if (Precision::IsInfinite(F)) + { + return Standard_False; + } + return Standard_True; +} + +//======================================================================= +//function : QuadricParameters +//purpose : +//======================================================================= +void Extrema_GlobOptFuncCQuadric::QuadricParameters(const math_Vector& theCT, + math_Vector& theUV ) const +{ + Standard_Real u, v; + // + //Arrays of extremity points parameters correspond to array of corner + //points myPTrim[] + Standard_Real uext[4] = { myUf, myUl, myUl, myUf }; + Standard_Real vext[4] = { myVf, myVf, myVl, myVl }; + gp_Pnt aCP = myC->Value(theCT(1)); + switch (mySType) + { + case GeomAbs_Plane: + ElSLib::Parameters(myPln, aCP, u, v); + break; + case GeomAbs_Cylinder: + ElSLib::Parameters(myCylinder, aCP, u, v); + break; + case GeomAbs_Cone: + ElSLib::Parameters(myCone, aCP, u, v); + break; + case GeomAbs_Sphere: + ElSLib::Parameters(mySphere, aCP, u, v); + break; + case GeomAbs_Torus: + ElSLib::Parameters(myTorus, aCP, u, v); + break; + default: + theUV(1) = myUf; + theUV(2) = myUl; + return; + } + // + if (mySType != GeomAbs_Plane) + { + if (myUl > 2. * M_PI + Precision::PConfusion()) + { + u += 2. * M_PI; + } + } + if (mySType == GeomAbs_Torus) + { + if (myVl > 2. * M_PI + Precision::PConfusion()) + { + v += 2. * M_PI; + } + } + + Standard_Real F = RealLast(); + if (u >= myUf && u <= myUl && v >= myVf && v <= myVl) + { + gp_Pnt aPS = myS->Value(u, v); + F = aCP.SquareDistance(aPS); + } + Standard_Integer i; + for (i = 0; i < 4; ++i) + { + Standard_Real Fi = aCP.SquareDistance(myPTrim[i]); + if (Fi < F) + { + F = Fi; + u = uext[i]; + v = vext[i]; + } + } + theUV(1) = u; + theUV(2) = v; +} + diff --git a/src/Extrema/Extrema_GlobOptFuncCQuadric.hxx b/src/Extrema/Extrema_GlobOptFuncCQuadric.hxx new file mode 100644 index 0000000000..e4cdb8b5f7 --- /dev/null +++ b/src/Extrema/Extrema_GlobOptFuncCQuadric.hxx @@ -0,0 +1,88 @@ +// Copyright (c) 2020 OPEN CASCADE SAS +// +// This file is part of Open CASCADE Technology software library. +// +// This library is free software; you can redistribute it and/or modify it under +// the terms of the GNU Lesser General Public License version 2.1 as published +// by the Free Software Foundation, with special exception defined in the file +// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT +// distribution for complete text of the license and disclaimer of any warranty. +// +// Alternatively, this file may be used under the terms of Open CASCADE +// commercial license or contractual agreement + +#ifndef _Extrema_GlobOptFuncCQuadric_HeaderFile +#define _Extrema_GlobOptFuncCQuadric_HeaderFile + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +//! This class implements function which calculate square Eucluidean distance +//! between point on surface and nearest point on Conic. +class Extrema_GlobOptFuncCQuadric : public math_MultipleVarFunction +{ +public: + + //! Curve and surface should exist during all the lifetime of Extrema_GlobOptFuncCQuadric. + Standard_EXPORT Extrema_GlobOptFuncCQuadric(const Adaptor3d_Curve *C); + + Standard_EXPORT Extrema_GlobOptFuncCQuadric(const Adaptor3d_Curve *C, + const Standard_Real theTf, + const Standard_Real theTl); + + Standard_EXPORT Extrema_GlobOptFuncCQuadric(const Adaptor3d_Curve *C, + const Adaptor3d_Surface *S); + + Standard_EXPORT void LoadQuad(const Adaptor3d_Surface *S, + const Standard_Real theUf, + const Standard_Real theUl, + const Standard_Real theVf, + const Standard_Real theVl); + + Standard_EXPORT virtual Standard_Integer NbVariables() const; + + Standard_EXPORT virtual Standard_Boolean Value(const math_Vector &theX, + Standard_Real &theF); + //! Parameters of quadric for point on curve defined by theCT + Standard_EXPORT void QuadricParameters(const math_Vector& theCT, + math_Vector& theUV) const; + +private: + + Standard_Boolean checkInputData(const math_Vector &X, + Standard_Real &ct); + + void value(Standard_Real ct, + Standard_Real &F); + + + const Adaptor3d_Curve *myC; + const Adaptor3d_Surface *myS; + GeomAbs_SurfaceType mySType; + gp_Pln myPln; + gp_Cone myCone; + gp_Cylinder myCylinder; + gp_Sphere mySphere; + gp_Torus myTorus; + gp_Pnt myPTrim[4]; + // Boundaries + Standard_Real myTf; + Standard_Real myTl; + Standard_Real myUf; + Standard_Real myUl; + Standard_Real myVf; + Standard_Real myVl; + +}; + +#endif diff --git a/src/Extrema/Extrema_GlobOptFuncConicS.cxx b/src/Extrema/Extrema_GlobOptFuncConicS.cxx new file mode 100644 index 0000000000..8ea81b15be --- /dev/null +++ b/src/Extrema/Extrema_GlobOptFuncConicS.cxx @@ -0,0 +1,259 @@ +// Copyright (c) 2020 OPEN CASCADE SAS +// +// This file is part of Open CASCADE Technology software library. +// +// This library is free software; you can redistribute it and/or modify it under +// the terms of the GNU Lesser General Public License version 2.1 as published +// by the Free Software Foundation, with special exception defined in the file +// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT +// distribution for complete text of the license and disclaimer of any warranty. +// +// Alternatively, this file may be used under the terms of Open CASCADE +// commercial license or contractual agreement + +#include + +#include +#include + +//F(u, v) = Conic.SquareDistance(myS(u, v)) + +//======================================================================= +//function : value +//purpose : +//======================================================================= +void Extrema_GlobOptFuncConicS::value(Standard_Real su, + Standard_Real sv, + Standard_Real &F) +{ + Standard_Real ct; + gp_Pnt aPS = myS->Value(su, sv); + switch (myCType) + { + case GeomAbs_Line: + ct = ElCLib::Parameter(myLin, aPS); + break; + case GeomAbs_Circle: + ct = ElCLib::Parameter(myCirc, aPS); + break; + case GeomAbs_Ellipse: + ct = ElCLib::Parameter(myElips, aPS); + break; + case GeomAbs_Hyperbola: + ct = ElCLib::Parameter(myHypr, aPS); + break; + case GeomAbs_Parabola: + ct = ElCLib::Parameter(myParab, aPS); + break; + default: + F = Precision::Infinite(); + return; + } + // + if (myCType == GeomAbs_Circle || myCType == GeomAbs_Ellipse) + { + if (myTl > 2. * M_PI + Precision::PConfusion()) + { + ct += 2. * M_PI; + } + } + F = RealLast(); + if (ct >= myTf && ct <= myTl) + { + gp_Pnt aPC = myC->Value(ct); + F = Min(F, aPS.SquareDistance(aPC)); + } + F = Min(F, aPS.SquareDistance(myCPf)); + F = Min(F, aPS.SquareDistance(myCPl)); + +} + + +//======================================================================= +//function : checkInputData +//purpose : +//======================================================================= +Standard_Boolean Extrema_GlobOptFuncConicS::checkInputData(const math_Vector &X, + Standard_Real &su, + Standard_Real &sv) +{ + Standard_Integer aStartIndex = X.Lower(); + su = X(aStartIndex); + sv = X(aStartIndex + 1); + + if (su < myUf || su > myUl || + sv < myVf || sv > myVl) + { + return Standard_False; + } + return Standard_True; +} + +//======================================================================= +//function : Extrema_GlobOptFuncConicS +//purpose : Constructor +//======================================================================= +Extrema_GlobOptFuncConicS::Extrema_GlobOptFuncConicS(const Adaptor3d_Surface *S, + const Standard_Real theUf, const Standard_Real theUl, + const Standard_Real theVf, const Standard_Real theVl) +: myS(S), myUf(theUf), myUl(theUl), + myVf(theVf), myVl(theVl) +{ + +} + +//======================================================================= +//function : Extrema_GlobOptFuncConicS +//purpose : Constructor +//======================================================================= +Extrema_GlobOptFuncConicS::Extrema_GlobOptFuncConicS(const Adaptor3d_Surface *S) + : myS(S), myUf(S->FirstUParameter()), myUl(S->LastUParameter()), + myVf(S->FirstVParameter()), myVl(S->LastVParameter()) +{ + +} +//======================================================================= +//function : Extrema_GlobOptFuncConicS +//purpose : Constructor +//======================================================================= +Extrema_GlobOptFuncConicS::Extrema_GlobOptFuncConicS(const Adaptor3d_Curve *C, + const Adaptor3d_Surface *S) + : myS(S), myUf(S->FirstUParameter()), myUl(S->LastUParameter()), + myVf(S->FirstVParameter()), myVl(S->LastVParameter()) +{ + Standard_Real aCTf = C->FirstParameter(); + Standard_Real aCTl = C->LastParameter(); + LoadConic(C, aCTf, aCTl); +} + +//======================================================================= +//function : LoadConic +//purpose : +//======================================================================= +void Extrema_GlobOptFuncConicS::LoadConic(const Adaptor3d_Curve *C, + const Standard_Real theTf, const Standard_Real theTl) +{ + myC = C; + myTf = theTf; + myTl = theTl; + if (myC->IsPeriodic()) + { + Standard_Real aTMax = 2. * M_PI + Precision::PConfusion(); + if (myTf > aTMax || myTf < -Precision::PConfusion() || + Abs(myTl - myTf) > aTMax) + { + ElCLib::AdjustPeriodic(0., 2. * M_PI, + Min(Abs(myTl - myTf) / 2, Precision::PConfusion()), + myTf, myTl); + } + } + myCPf = myC->Value(myTf); + myCPl = myC->Value(myTl); + myCType = myC->GetType(); + switch (myCType) + { + case GeomAbs_Line: + myLin = myC->Line(); + break; + case GeomAbs_Circle: + myCirc = myC->Circle(); + break; + case GeomAbs_Ellipse: + myElips = myC->Ellipse(); + break; + case GeomAbs_Hyperbola: + myHypr = myC->Hyperbola(); + break; + case GeomAbs_Parabola: + myParab = myC->Parabola(); + break; + default: + break; + } +} + +//======================================================================= +//function : NbVariables +//purpose : +//======================================================================= +Standard_Integer Extrema_GlobOptFuncConicS::NbVariables() const +{ + return 2; +} + +//======================================================================= +//function : Value +//purpose : +//======================================================================= +Standard_Boolean Extrema_GlobOptFuncConicS::Value(const math_Vector &X, + Standard_Real &F) +{ + Standard_Real su, sv; + if (!checkInputData(X, su, sv)) + return Standard_False; + + value(su, sv, F); + if (Precision::IsInfinite(F)) + { + return Standard_False; + } + return Standard_True; +} + +//======================================================================= +//function : ConicParameter +//purpose : +//======================================================================= +Standard_Real Extrema_GlobOptFuncConicS::ConicParameter(const math_Vector& theUV) const +{ + Standard_Real ct; + gp_Pnt aPS = myS->Value(theUV(1), theUV(2)); + switch (myCType) + { + case GeomAbs_Line: + ct = ElCLib::Parameter(myLin, aPS); + break; + case GeomAbs_Circle: + ct = ElCLib::Parameter(myCirc, aPS); + break; + case GeomAbs_Ellipse: + ct = ElCLib::Parameter(myElips, aPS); + break; + case GeomAbs_Hyperbola: + ct = ElCLib::Parameter(myHypr, aPS); + break; + case GeomAbs_Parabola: + ct = ElCLib::Parameter(myParab, aPS); + break; + default: + ct = myTf; + return ct; + } + // + if (myCType == GeomAbs_Circle || myCType == GeomAbs_Ellipse) + { + if (myTl > 2. * M_PI + Precision::PConfusion()) + { + ct += 2. * M_PI; + } + } + Standard_Real F = RealLast(); + if (ct >= myTf && ct <= myTl) + { + gp_Pnt aPC = myC->Value(ct); + F = Min(F, aPS.SquareDistance(aPC)); + } + Standard_Real Fext = aPS.SquareDistance(myCPf); + if (Fext < F) + { + F = Fext; + ct = myTf; + } + Fext = aPS.SquareDistance(myCPl); + if (Fext < F) + { + F = Fext; + ct = myTl; + } + return ct; +} \ No newline at end of file diff --git a/src/Extrema/Extrema_GlobOptFuncConicS.hxx b/src/Extrema/Extrema_GlobOptFuncConicS.hxx new file mode 100644 index 0000000000..296cb67029 --- /dev/null +++ b/src/Extrema/Extrema_GlobOptFuncConicS.hxx @@ -0,0 +1,90 @@ + +// Copyright (c) 2020 OPEN CASCADE SAS +// +// This file is part of Open CASCADE Technology software library. +// +// This library is free software; you can redistribute it and/or modify it under +// the terms of the GNU Lesser General Public License version 2.1 as published +// by the Free Software Foundation, with special exception defined in the file +// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT +// distribution for complete text of the license and disclaimer of any warranty. +// +// Alternatively, this file may be used under the terms of Open CASCADE +// commercial license or contractual agreement + +#ifndef _Extrema_GlobOptFuncConicS_HeaderFile +#define _Extrema_GlobOptFuncConicS_HeaderFile + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +//! This class implements function which calculate square Eucluidean distance +//! between point on surface and nearest point on Conic. +class Extrema_GlobOptFuncConicS : public math_MultipleVarFunction +{ +public: + + //! Curve and surface should exist during all the lifetime of Extrema_GlobOptFuncConicS. + Standard_EXPORT Extrema_GlobOptFuncConicS(const Adaptor3d_Curve *C, + const Adaptor3d_Surface *S); + + Standard_EXPORT Extrema_GlobOptFuncConicS(const Adaptor3d_Surface *S); + + Standard_EXPORT Extrema_GlobOptFuncConicS(const Adaptor3d_Surface *S, + const Standard_Real theUf, + const Standard_Real theUl, + const Standard_Real theVf, + const Standard_Real theVl); + + Standard_EXPORT void LoadConic(const Adaptor3d_Curve *S, const Standard_Real theTf, const Standard_Real theTl); + + Standard_EXPORT virtual Standard_Integer NbVariables() const; + + Standard_EXPORT virtual Standard_Boolean Value(const math_Vector &theX, + Standard_Real &theF); + + //! Parameter of conic for point on surface defined by theUV + Standard_EXPORT Standard_Real ConicParameter(const math_Vector& theUV) const; + +private: + + Standard_Boolean checkInputData(const math_Vector &X, + Standard_Real &su, + Standard_Real &sv); + + void value(Standard_Real su, + Standard_Real sv, + Standard_Real &F); + + + const Adaptor3d_Curve *myC; + const Adaptor3d_Surface *myS; + GeomAbs_CurveType myCType; + gp_Lin myLin; + gp_Circ myCirc; + gp_Elips myElips; + gp_Hypr myHypr; + gp_Parab myParab; + gp_Pnt myCPf; + gp_Pnt myCPl; + //Boundaries + Standard_Real myTf; + Standard_Real myTl; + Standard_Real myUf; + Standard_Real myUl; + Standard_Real myVf; + Standard_Real myVl; + +}; + +#endif diff --git a/src/Extrema/FILES b/src/Extrema/FILES index 34ec6a72e0..1e3d13ee43 100644 --- a/src/Extrema/FILES +++ b/src/Extrema/FILES @@ -99,6 +99,10 @@ Extrema_GlobOptFuncCC.cxx Extrema_GlobOptFuncCC.hxx Extrema_GlobOptFuncCS.cxx Extrema_GlobOptFuncCS.hxx +Extrema_GlobOptFuncConicS.cxx +Extrema_GlobOptFuncConicS.hxx +Extrema_GlobOptFuncCQuadric.cxx +Extrema_GlobOptFuncCQuadric.hxx Extrema_GLocateExtPC.gxx Extrema_HArray1OfPOnCurv.hxx Extrema_HArray1OfPOnCurv2d.hxx diff --git a/src/math/math_PSOParticlesPool.cxx b/src/math/math_PSOParticlesPool.cxx index f3b7bf6ca6..8fd5dc28cb 100644 --- a/src/math/math_PSOParticlesPool.cxx +++ b/src/math/math_PSOParticlesPool.cxx @@ -31,7 +31,7 @@ math_PSOParticlesPool::math_PSOParticlesPool(const Standard_Integer theParticles { myParticlesCount = theParticlesCount; myDimensionCount = theDimensionCount; - + myMemory.Init(0.); // Pointers adjusting. Standard_Integer aParIdx, aShiftIdx; for(aParIdx = 1; aParIdx <= myParticlesCount; ++aParIdx) diff --git a/tests/bugs/modalg_5/bug25232_8 b/tests/bugs/modalg_5/bug25232_8 index 49d0e7f63a..61892b8be4 100644 --- a/tests/bugs/modalg_5/bug25232_8 +++ b/tests/bugs/modalg_5/bug25232_8 @@ -27,5 +27,5 @@ mkvolume result fcon3 fp checkprops result -s 1706.51 checkshape result -checknbshapes result -vertex 4 -edge 5 -wire 2 -face 2 -shell 1 -solid 1 -compsolid 0 -compound 0 -shape 15 +checknbshapes result -vertex 5 -edge 6 -wire 2 -face 2 -shell 1 -solid 1 -compsolid 0 -compound 0 -shape 17 checkview -display result -2d -path ${imagedir}/${test_image}.png diff --git a/tests/bugs/modalg_7/bug27928 b/tests/bugs/modalg_7/bug27928 index d48df7980e..24073c9b6c 100644 --- a/tests/bugs/modalg_7/bug27928 +++ b/tests/bugs/modalg_7/bug27928 @@ -36,7 +36,7 @@ checkprops r2 -s 2075.44 -v 1489.33 checknbshapes r3 -wire 5 -face 5 -shell 1 -solid 1 checkprops r3 -s 2521.83 -v 1824.69 -checknbshapes r4 -vertex 12 -edge 18 -t +checknbshapes r4 -vertex 13 -edge 19 -t checkprops r4 -l 825.645 checkview -display r0 -2d -path ${imagedir}/${test_image}.png diff --git a/tests/bugs/modalg_7/bug29580_1 b/tests/bugs/modalg_7/bug29580_1 index 260b75ebc2..010433c134 100644 --- a/tests/bugs/modalg_7/bug29580_1 +++ b/tests/bugs/modalg_7/bug29580_1 @@ -9,6 +9,7 @@ puts "" restore [locate_data_file bug29580_Cylinder.brep] b1 restore [locate_data_file bug29580_Solid.brep] b2 +settolerance b2 e 1.e-5 bfuse result b1 b2 foreach f [explode result f] { @@ -21,7 +22,7 @@ foreach f [explode result f] { } checkshape result -checknbshapes result -wire 11 -face 10 -shell 1 -solid 1 -checkprops result -s 865.851 -v 1622.17 +checknbshapes result -wire 14 -face 13 -shell 1 -solid 1 +checkprops result -s 866.155 -v 1622.85 checkview -display result -2d -path ${imagedir}/${test_image}.png diff --git a/tests/bugs/modalg_7/bug29839 b/tests/bugs/modalg_7/bug29839 new file mode 100644 index 0000000000..3b0d828bcb --- /dev/null +++ b/tests/bugs/modalg_7/bug29839 @@ -0,0 +1,22 @@ +puts "=================================================================" +puts "OCC29839: Unexpected extrema behavior" +puts "=================================================================" +puts "" + +restore [locate_data_file bug29839.brep] s + +explode s +distmini d s_1 s_2 +checknbshapes d -edge 0 -vertex 1 + +set dist [dval d_val] + +if { ${dist} > 2.e-10} { + puts "Error: bad distance" +} + +smallview +X+Y +donly s_1 s_2 d +fit + +checkview -screenshot -2d -path ${imagedir}/${test_image}.png