if (ext.NbExt()!=0){
Extrema_POnCurv POnC, POnL;
ext.Points(1, POnC, POnL);
- param.ChangeValue(nb) =POnC.Parameter();
+ if (POnC.Value().Distance(POnL.Value()) < Precision::Confusion())
+ param.ChangeValue(nb) =POnC.Parameter();
+ else
+ {
+ if (!cproj.Value(nb).IsNull()) {
+ cproj.Value(nb)->D0(cproj.Value(nb)->LastParameter(),p01);
+ }
+ else if (!cproj.Value(nb+1).IsNull()) {
+ cproj.Value(nb+1)->D0(cproj.Value(nb+1)->FirstParameter(),p01);
+ }
+ }
}
}
if (!ext.IsDone()||ext.NbExt()==0) {
-- generic classes for CURVE-CURVE extremas:
----------------------------------------------
private generic class FuncExtCC, SeqPOnC;
- generic class GenExtCC, CCF;
+ generic class GenExtCC;
generic class GenLocateExtCC, CCLocF;
- generic class CurveCache;
----------------------------------------------
-- Curve-Curve:
-- 3d:
class ExtCC;
-
- class CCache instantiates CurveCache from Extrema
- (Curve from Adaptor3d,
- Pnt from gp,
- HArray1OfPnt from TColgp);
class ECC instantiates GenExtCC from Extrema
(Curve from Adaptor3d,
CurveTool from Extrema,
Curve from Adaptor3d,
CurveTool from Extrema,
- CCache from Extrema,
HArray1OfPnt from TColgp,
POnCurv from Extrema,
Pnt from gp,
Vec from gp);
class LocateExtCC;
-
- class LCCache instantiates CurveCache from Extrema
- (Curve from Adaptor3d,
- Pnt from gp,
- HArray1OfPnt from TColgp);
class ELCC instantiates GenExtCC from Extrema
(Curve from Adaptor3d,
CurveTool from Extrema,
Curve from Adaptor3d,
CurveTool from Extrema,
- LCCache from Extrema,
HArray1OfPnt from TColgp,
POnCurv from Extrema,
Pnt from gp,
-- 2d:
class ExtCC2d;
- class CCache2d instantiates CurveCache from Extrema
- (Curve2d from Adaptor2d,
- Pnt2d from gp,
- HArray1OfPnt2d from TColgp);
-
class ECC2d instantiates GenExtCC from Extrema
(Curve2d from Adaptor2d,
Curve2dTool from Extrema,
Curve2d from Adaptor2d,
Curve2dTool from Extrema,
- CCache2d from Extrema,
HArray1OfPnt2d from TColgp,
POnCurv2d from Extrema,
Pnt2d from gp,
class LocateExtCC2d;
- class LCCache2d instantiates CurveCache from Extrema
- (Curve2d from Adaptor2d,
- Pnt2d from gp,
- HArray1OfPnt2d from TColgp);
class ELCC2d instantiates GenExtCC from Extrema
(Curve2d from Adaptor2d,
Curve2dTool from Extrema,
Curve2d from Adaptor2d,
Curve2dTool from Extrema,
- LCCache2d from Extrema,
HArray1OfPnt2d from TColgp,
POnCurv2d from Extrema,
Pnt2d from gp,
+++ /dev/null
--- Created on: 2008-12-28
--- Created by: Roman Lygin
--- Copyright (c) 2008-2014 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.
-
--- roman.lygin@gmail.com
-
-generic class CurveCache from Extrema
- (Curve as any; -- Adaptor3d_Curve or Adaptor2d_Curve2d
- Pnt as any; -- gp_Pnt or gp_Pnt2d
- ArrayOfPnt as Transient from Standard) -- TColgp_HArray1OfPnt or TColgp_HArray1OfPnt2d
-inherits Transient from Standard
-
- ---Purpose:
-
-uses
- Array1OfReal from TColStd,
- HArray1OfReal from TColStd
-
-raises NotDone from StdFail
-
-is
-
- Create returns mutable CurveCache from Extrema;
-
- Create (theC: Curve; theUFirst, theULast: Real; theNbSamples: Integer;
- theToCalculate: Boolean)returns mutable CurveCache from Extrema;
-
- Create (theC: Curve; theUFirst, theULast: Real;
- IntervalsCN: Array1OfReal from TColStd;
- StartIndex, EndIndex: Integer;
- Coeff: Real)
- returns mutable CurveCache from Extrema;
-
- SetCurve (me: mutable; theC: Curve; theNbSamples: Integer;
- theToCalculate: Boolean);
- ---Purpose: Sets the \a theRank-th curve (1 or 2) and its parameters.
- -- Caches sample points for further reuse in Perform().
-
- SetCurve (me: mutable; theC: Curve;
- theUFirst, theULast: Real; theNbSamples: Integer; theToCalculate: Boolean);
- ---Purpose:
-
- SetRange (me: mutable; Uinf, Usup: Real; theToCalculate: Boolean);
- ---Purpose:
-
- CalculatePoints (me: mutable);
- ---Purpose: Calculates points along the curve and stores
- -- them in internal array for further reuse in Perform().
-
- IsValid (me) returns Boolean;
- ---C++: inline
- ---Purpose: Returns True if the points array has already been calculated
- -- for specified curve and range
-
- Parameters (me) returns HArray1OfReal from TColStd
- raises NotDone from StdFail;
- ---C++: inline
- ---C++: return const &
- ---Purpose: Raises StdFail_NotDone if points have not yet been calculated.
-
- Points (me) returns ArrayOfPnt
- raises NotDone from StdFail;
- ---C++: inline
- ---C++: return const &
- ---Purpose: Raises StdFail_NotDone if points have not yet been calculated.
-
- CurvePtr (me) returns Address from Standard;
- ---C++: inline
- ---Purpose:
-
- NbSamples (me) returns Integer;
- ---C++: inline
- ---Purpose:
-
- FirstParameter (me) returns Real;
- ---C++: inline
- ---Purpose:
-
- LastParameter (me) returns Real;
- ---C++: inline
- ---Purpose:
-
- TrimFirstParameter (me) returns Real;
- ---C++: inline
- ---Purpose: For bounded curves returns FirstParameter(), for unbounded - -1e-10.
-
- TrimLastParameter (me) returns Real;
- ---C++: inline
- ---Purpose: For bounded curves returns LastParameter(), for unbounded - 1e-10.
-
-fields
-
- myC : Address from Standard;
- myFirst : Real;
- myLast : Real;
- myTrimFirst : Real;
- myTrimLast : Real;
- myNbSamples : Integer;
- myParamArray : HArray1OfReal from TColStd;
- myPntArray : ArrayOfPnt;
- myIsArrayValid: Boolean;
-
-end;
+++ /dev/null
-// Created on: 2008-12-28
-// Created by: Roman Lygin
-// Copyright (c) 2008-2014 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.
-
-// roman.lygin@gmail.com
-
-#include <Precision.hxx>
-
-//=======================================================================
-//function : Extrema_CurveCache
-//purpose :
-//=======================================================================
-
-Extrema_CurveCache::Extrema_CurveCache(const Curve& theC,
- const Standard_Real theUFirst,
- const Standard_Real theULast,
- const Standard_Integer theNbSamples,
- const Standard_Boolean theToCalculate) :
- myC (0), myNbSamples (-1), myIsArrayValid (Standard_False)
-{
- SetCurve (theC, theUFirst, theULast, theNbSamples, theToCalculate);
-}
-
-//=======================================================================
-//function : Extrema_CurveCache
-//purpose : with sampling by knots and between them
-//=======================================================================
-
-Extrema_CurveCache::Extrema_CurveCache(const Curve& theC,
- const Standard_Real theUFirst,
- const Standard_Real theULast,
- const TColStd_Array1OfReal& IntervalsCN,
- const Standard_Integer StartIndex,
- const Standard_Integer EndIndex,
- const Standard_Real Coeff)
-{
- myC = (Standard_Address)&theC;
- myIsArrayValid = Standard_False;
- myParamArray.Nullify();
- myPntArray.Nullify();
-
- myTrimFirst = myFirst = theUFirst;
- myTrimLast = myLast = theULast;
-
- Standard_Integer Nbp = 3;
- if (2 * Coeff < 10000.0)
- Nbp = Max((Standard_Integer) (2 * Coeff), Nbp);
- myNbSamples = (EndIndex - StartIndex)*Nbp + 1;
-
- const Standard_Integer aNbSTresh = 10000;
- if (myNbSamples > aNbSTresh)
- {
- Nbp = aNbSTresh / (EndIndex - StartIndex);
- myNbSamples = (EndIndex - StartIndex)*Nbp + 1;
- }
-
- //Cache points
- myParamArray = new TColStd_HArray1OfReal(1, myNbSamples);
- myPntArray = new ArrayOfPnt (1, myNbSamples);
-
- const Curve& aCurve = *((Curve*)myC);
-
- Standard_Integer i, j, k = 1;
- Standard_Real aPar;
- for (i = StartIndex; i < EndIndex; i++)
- {
- Standard_Real delta = (IntervalsCN(i+1) - IntervalsCN(i)) / Nbp;
- for (j = 0; j < Nbp; j++)
- {
- aPar = IntervalsCN(i) + j*delta;
- myParamArray->SetValue(k, aPar);
- myPntArray->SetValue(k++, aCurve.Value(aPar));
- }
- }
- Standard_Real aDelta = (myTrimLast - myTrimFirst) / myNbSamples / 200.;
- myParamArray->SetValue(1, myTrimFirst + aDelta);
- myPntArray->SetValue(1, aCurve.Value(myTrimFirst + aDelta));
- myParamArray->SetValue(myNbSamples, myTrimLast - aDelta);
- myPntArray->SetValue(myNbSamples, aCurve.Value(myTrimLast - aDelta));
-
- myIsArrayValid = Standard_True; //cache is available now
-}
-
-//=======================================================================
-//function : Extrema_CurveCache
-//purpose :
-//=======================================================================
-
-Extrema_CurveCache::Extrema_CurveCache() : myC (0), myNbSamples (-1),
- myIsArrayValid (Standard_False)
-{
-}
-
-//=======================================================================
-//function : SetCurve
-//purpose :
-//=======================================================================
-
-void Extrema_CurveCache::SetCurve (const Curve& theC,
- const Standard_Integer theNbSamples,
- const Standard_Boolean theToCalculate)
-{
- myC = (Standard_Address)&theC;
- myNbSamples = theNbSamples;
- myIsArrayValid = Standard_False;
- myParamArray.Nullify();
- myPntArray.Nullify();
- if (theToCalculate) {
- CalculatePoints();
- }
-}
-
-//=======================================================================
-//function : SetCurve
-//purpose :
-//=======================================================================
-
-void Extrema_CurveCache::SetCurve (const Curve& theC,
- const Standard_Real theUFirst,
- const Standard_Real theULast,
- const Standard_Integer theNbSamples,
- const Standard_Boolean theToCalculate)
-{
- SetCurve (theC, theNbSamples, Standard_False); //no calculation
- SetRange (theUFirst, theULast, theToCalculate);
-}
-
-//=======================================================================
-//function : SetRange
-//purpose :
-//=======================================================================
-
-void Extrema_CurveCache::SetRange (const Standard_Real theUFirst,
- const Standard_Real theULast,
- const Standard_Boolean theToCalculate)
-{
- //myTrimFirst and myTrimLast are used to compute values on unlimited curves
- myTrimFirst = myFirst = theUFirst;
- if (Precision::IsInfinite(myTrimFirst)){
- myTrimFirst = -1.0e+10;
- }
- myTrimLast = myLast = theULast;
- if (Precision::IsInfinite(myTrimLast)){
- myTrimLast = 1.0e+10;
- }
-
- myIsArrayValid = Standard_False;
- myParamArray.Nullify();
- myPntArray.Nullify();
- if (theToCalculate) {
- CalculatePoints();
- }
-}
-
-//=======================================================================
-//function : CalculatePoints
-//purpose :
-//=======================================================================
-
-void Extrema_CurveCache::CalculatePoints()
-{
- if (myIsArrayValid) return; //no need to recalculate if nothing has changed
- const Curve& aCurve = *((Curve*)myC);
-
- // compute myNbSamples point along the [myTrimFirst, myTrimLast] range
-
- Standard_Real aDelta = myTrimLast - myTrimFirst;
- Standard_Real aPar0 = aDelta / myNbSamples / 100.;
- aDelta = (aDelta - aPar0) / (myNbSamples - 1);
- aPar0 = myTrimFirst + (aPar0/2.);
-
- //Cache points
-
- myParamArray = new TColStd_HArray1OfReal(1, myNbSamples);
- myPntArray = new ArrayOfPnt (1, myNbSamples);
-
- Standard_Integer i;
- Standard_Real aPar;
- for (i = 1, aPar = aPar0; i <= myNbSamples; i++, aPar += aDelta) {
- myParamArray->SetValue(i, aPar);
- myPntArray->SetValue (i, aCurve.Value (aPar));
- }
-
- myIsArrayValid = Standard_True; //cache is available now
-}
+++ /dev/null
-// Created on: 2008-12-28
-// Created by: Roman Lygin
-// Copyright (c) 2008-2014 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.
-
-// roman.lygin@gmail.com
-
-#include <StdFail_NotDone.hxx>
-#include <TColStd_HArray1OfReal.hxx>
-
-//=======================================================================
-//function : IsValid
-//purpose :
-//=======================================================================
-
-inline Standard_Boolean Extrema_CurveCache::IsValid() const
-{
- return myIsArrayValid;
-}
-
-//=======================================================================
-//function : Parameters
-//purpose :
-//=======================================================================
-
-inline const Handle(TColStd_HArray1OfReal)& Extrema_CurveCache::Parameters() const
-{
- StdFail_NotDone_Raise_if (!myIsArrayValid, "Extrema_CurveCache::Parameters()")
- return myParamArray;
-}
-
-//=======================================================================
-//function : Points
-//purpose :
-//=======================================================================
-
-inline const Handle(ArrayOfPnt)& Extrema_CurveCache::Points() const
-{
- StdFail_NotDone_Raise_if (!myIsArrayValid, "Extrema_CurveCache::Points()")
- return myPntArray;
-}
-
-//=======================================================================
-//function : CurvePtr
-//purpose :
-//=======================================================================
-
-inline Standard_Address Extrema_CurveCache::CurvePtr() const
-{
- return myC;
-}
-
-//=======================================================================
-//function : NbSamples
-//purpose :
-//=======================================================================
-
-inline Standard_Integer Extrema_CurveCache::NbSamples() const
-{
- return myNbSamples;
-}
-
-//=======================================================================
-//function : FirstParameter
-//purpose :
-//=======================================================================
-
-inline Standard_Real Extrema_CurveCache::FirstParameter() const
-{
- return myFirst;
-}
-
-//=======================================================================
-//function : LastParameter
-//purpose :
-//=======================================================================
-
-inline Standard_Real Extrema_CurveCache::LastParameter() const
-{
- return myLast;
-}
-
-//=======================================================================
-//function : TrimFirstParameter
-//purpose :
-//=======================================================================
-
-inline Standard_Real Extrema_CurveCache::TrimFirstParameter() const
-{
- return myTrimFirst;
-}
-
-//=======================================================================
-//function : TrimLastParameter
-//purpose :
-//=======================================================================
-
-inline Standard_Real Extrema_CurveCache::TrimLastParameter() const
-{
- return myTrimLast;
-}
myInf: Real [2];
mySup: Real [2];
myTol: Real [2];
- myCacheLists: ListOfTransient from TColStd [2]; -- lists of Handle(Extrema_CCache)
P1f: Pnt;
P1l: Pnt;
P2f: Pnt;
#include <Adaptor3d_Curve.hxx>
#include <Extrema_CurveTool.hxx>
-#include <Extrema_CCache.hxx>
//=======================================================================
//function : Extrema_ExtCC
const Standard_Real V1,
const Standard_Real V2,
const Standard_Real TolC1,
- const Standard_Real TolC2) :
- myDone (Standard_False)
+ const Standard_Real TolC2)
+: myECC(C1, C2, U1, U2, V1, V2),
+ myDone (Standard_False)
{
SetCurve (1, C1, U1, U2);
SetCurve (2, C2, V1, V2);
Extrema_ExtCC::Extrema_ExtCC(const Adaptor3d_Curve& C1,
const Adaptor3d_Curve& C2,
const Standard_Real TolC1,
- const Standard_Real TolC2) :
- myDone (Standard_False)
+ const Standard_Real TolC2)
+: myECC(C1, C2),
+ myDone (Standard_False)
{
SetCurve (1, C1, C1.FirstParameter(), C1.LastParameter());
SetCurve (2, C2, C2.FirstParameter(), C2.LastParameter());
Standard_OutOfRange_Raise_if (theRank < 1 || theRank > 2, "Extrema_ExtCC::SetCurve()")
Standard_Integer anInd = theRank - 1;
myC[anInd] = (Standard_Address)&C;
-
- //clear the previous cache to rebuild it later in Perform()
- myCacheLists[anInd].Clear();
}
//=======================================================================
void Extrema_ExtCC::Perform()
{
Standard_NullObject_Raise_if (!myC[0] || !myC[1], "Extrema_ExtCC::Perform()")
+ myECC.SetParams(*((Adaptor3d_Curve*)myC[0]),
+ *((Adaptor3d_Curve*)myC[1]), myInf[0], mySup[0], myInf[1], mySup[1]);
myDone = Standard_False;
mypoints.Clear();
mySqDist.Clear();
if (Precision::IsInfinite(U12) || Precision::IsInfinite(U22)) mydist22 = RealLast();
else mydist22 = P1l.SquareDistance(P2l);
- myECC.SetTolerance (Tol);
-
//Depending on the types of curves, the algorithm is chosen:
//- _ExtElC, when one of the curves is a line and the other is elementary,
// or there are two circles;
Results(CCXtrem, U11, U12, U21, U22);
}
else {
- Standard_Integer i;
- Standard_Integer aNbS = 32; //default number of sample points per interval (why 32?)
- for (i = 0; i < 2; i++) {
- TColStd_ListOfTransient& aCacheList = myCacheLists[i];
- if (aCacheList.IsEmpty()) {
- //no caches, build them
- Adaptor3d_Curve& aC = *(Adaptor3d_Curve*)myC[i];
- //single interval from myInf[i] to mySup[i]
- Handle(Extrema_CCache) aCache = new Extrema_CCache(aC, myInf[i], mySup[i], aNbS, Standard_True);
- aCacheList.Append (aCache);
- }
- }
- Handle(Extrema_CCache) aCache1 = Handle(Extrema_CCache)::DownCast (myCacheLists[0].First());
- myECC.SetCurveCache (1, aCache1);
- Handle(Extrema_CCache) aCache2 = Handle(Extrema_CCache)::DownCast (myCacheLists[1].First());
- myECC.SetCurveCache (2, aCache2);
myECC.Perform();
- Results (myECC, U11, U12, U21, U22);
+ Results(myECC, U11, U12, U21, U22);
}
} else {
- //common case - use _GenExtCC
- //1. check and prepare caches
-
- Standard_Integer i;
- Standard_Integer aNbS[2] = {32, 32}, aNbInter[2] = {1, 1};
- Standard_Real Coeff[2] = {1., 1.};
- Standard_Real rf = 0., rl = 0., LL[2] = {-1., -1.};
- Standard_Boolean KnotSampling[2] = {Standard_False, Standard_False};
- for (i = 0; i < 2; i++) {
- TColStd_ListOfTransient& aCacheList = myCacheLists[i];
- if (aCacheList.IsEmpty()) {
- //no caches, build them
- Adaptor3d_Curve& aC = *(Adaptor3d_Curve*)myC[i];
-
- Standard_Real du1 = 0., t = 0.;
- gp_Pnt P1, P2;
- GeomAbs_CurveType aType = aC.GetType();
- switch (aType) {
- case GeomAbs_BezierCurve:
- aNbS[i] = aC.NbPoles() * 2;
- break;
- case GeomAbs_BSplineCurve:
- if (aC.Degree() <= 3 &&
- aC.NbKnots() > 2)
- KnotSampling[i] = Standard_True;
-
- aNbS[i] = aC.NbPoles() * 2;
- rf = (Extrema_CurveTool::BSpline(*((Adaptor3d_Curve*)myC[i])))->FirstParameter();
- rl = (Extrema_CurveTool::BSpline(*((Adaptor3d_Curve*)myC[i])))->LastParameter();
- aNbS[i] = (Standard_Integer) ( aNbS[i] * ((mySup[i] - myInf[i]) / (rl - rf)) + 1 );
- case GeomAbs_OtherCurve:
- case GeomAbs_Line:
- //adjust number of sample points for Lines, B-Splines and Other curves
- aNbInter[i] = aC.NbIntervals (GeomAbs_C2);
- aNbS[i] = Max(aNbS[i] / aNbInter[i], 3);
- LL[i] = 0.;
- du1 = (mySup[i] - myInf[i]) / ((aNbS[i]-1)*aNbInter[i]);
- P1 = Extrema_CurveTool::Value(*((Adaptor3d_Curve*)myC[i]), myInf[i]);
- for(t = myInf[i] + du1; t <= mySup[i]; t += du1) {
- P2 = Extrema_CurveTool::Value(*((Adaptor3d_Curve*)myC[i]), t);
- LL[i] += P1.Distance(P2);
- P1 = P2;
- }
- LL[i] /= ((aNbS[i]-1)*aNbInter[i]);
- break;
- default:
- break;
- }
- }
- }
-
- if(LL[0] > 0. && LL[1] > 0.) {
- if(LL[0] > 4.*LL[1]) {
- Coeff[0] = LL[0]/LL[1]/2.;
- if (aNbS[0] * Coeff[0] <= 10000.0)
- aNbS[0] = (Standard_Integer) ( aNbS[0] * Coeff[0] );
- }
- else if(LL[1] > 4.*LL[0]) {
- Coeff[1] = LL[1]/LL[0]/2.;
- if (aNbS[1] * Coeff[1] <= 10000.0)
- aNbS[1] = (Standard_Integer) (aNbS[1] * Coeff[1] );
- }
- }
- //modified by NIZNHY-PKV Tue Apr 17 10:01:32 2012f
- Standard_Integer aNbSTresh;
- //
- aNbSTresh=10000;
- //
- for (i = 0; i < 2; ++i) {
- if (aNbS[i]>aNbSTresh) {
- aNbS[i]=aNbSTresh;
- }
- }
- //modified by NIZNHY-PKV Tue Apr 17 10:01:34 2012t
- for (i = 0; i < 2; i++) {
- TColStd_ListOfTransient& aCacheList = myCacheLists[i];
- if (aCacheList.IsEmpty()) {
- //no caches, build them
- Adaptor3d_Curve& aC = *(Adaptor3d_Curve*)myC[i];
-
- if (aC.GetType() >= GeomAbs_BSplineCurve)
- {
- //can be multiple intervals, one cache per one C2 interval
- TColStd_Array1OfReal anArr (1, aNbInter[i] + 1);
- aC.Intervals (anArr, GeomAbs_C2);
-
- if (KnotSampling[i])
- {
- Standard_Integer NbIntervCN = aC.NbIntervals(GeomAbs_CN);
- TColStd_Array1OfReal IntervalsCN(1, NbIntervCN + 1);
- aC.Intervals(IntervalsCN, GeomAbs_CN);
-
- Standard_Integer j = 1, start_j = 1, k;
- Standard_Real NextKnot;
- for (k = 1; k <= aNbInter[i]; k++)
- {
- do
- {
- NextKnot = IntervalsCN(j+1);
- j++;
- }
- while (NextKnot != anArr(k+1));
-
- Handle(Extrema_CCache) aCache =
- new Extrema_CCache(aC, anArr(k), anArr(k+1),
- IntervalsCN, start_j, j, Coeff[i]);
- aCacheList.Append (aCache);
-
- start_j = j;
- }
- }
- else
- {
- for (Standard_Integer k = 1; k <= aNbInter[i]; k++) {
- Handle(Extrema_CCache) aCache =
- new Extrema_CCache(aC, anArr(k), anArr(k+1), aNbS[i], Standard_True);
- aCacheList.Append (aCache);
- }
- }
- }
- else
- {
- //single interval from myInf[i] to mySup[i]
- Handle(Extrema_CCache) aCache = new Extrema_CCache (aC,
- myInf[i], mySup[i], aNbS[i], Standard_True);
- aCacheList.Append (aCache);
- }
- }
- }
-
- //2. process each cache from one list with each cache from the other
- TColStd_ListIteratorOfListOfTransient anIt1 (myCacheLists[0]);
- for (; anIt1.More(); anIt1.Next()) {
- Handle(Extrema_CCache) aCache1 = Handle(Extrema_CCache)::DownCast (anIt1.Value());
- myECC.SetCurveCache (1, aCache1);
- TColStd_ListIteratorOfListOfTransient anIt2 (myCacheLists[1]);
- for (; anIt2.More(); anIt2.Next()) {
- Handle(Extrema_CCache) aCache2 = Handle(Extrema_CCache)::DownCast (anIt2.Value());
- myECC.SetCurveCache (2, aCache2);
- myECC.Perform();
- Results(myECC, U11, U12, U21, U22);
- }
- }
+ myECC.Perform();
+ Results(myECC, U11, U12, U21, U22);
}
}
SequenceOfReal from TColStd,
Curve2d from Adaptor2d,
Curve2dTool from Extrema,
- CCache2d from Extrema,
ECC2d from Extrema
is static protected;
--- modified by NIZHNY-EAP Thu Jan 27 16:53:25 2000 ___BEGIN___
- Results(me: in out;AlgExt: ECC2d from Extrema; C : Curve2d from Adaptor2d;
--- modified by NIZHNY-EAP Thu Jan 27 16:53:26 2000 ___END___
+ Results(me: in out;AlgExt: ECC2d from Extrema;
Ut11, Ut12, Ut21, Ut22: Real;
Period1 : Real from Standard = 0.0;
Period2 : Real from Standard = 0.0)
#include <Extrema_Curve2dTool.hxx>
-//=======================================================================
-//function : IsParallelDot
-//purpose : This function returns True if angle between <theV1> and
-//<theV2> vectors or between <theV1> and <-theV2> vectors is less than
-//AngTol.
-//=======================================================================
-static Standard_Boolean IsParallelDot( gp_Vec2d theV1,
- gp_Vec2d theV2,
- Standard_Real AngTol)
- {
- return Abs(theV1.Dot(theV2)) >=
- theV1.Magnitude()*theV2.Magnitude()*cos(AngTol);
- }
-
-Extrema_ExtCC2d::Extrema_ExtCC2d() {}
+Extrema_ExtCC2d::Extrema_ExtCC2d()
+{
+}
Extrema_ExtCC2d::Extrema_ExtCC2d(const Adaptor2d_Curve2d& C1,
{
mypoints.Clear();
mySqDist.Clear();
- Standard_Integer NbU = 32, NbV = 32;
GeomAbs_CurveType type1 = Extrema_Curve2dTool::GetType(C1), type2 = Extrema_Curve2dTool::GetType(*((Adaptor2d_Curve2d*)myC));
Standard_Real U11, U12, U21, U22, Tol = Min(mytolc1, mytolc2);
// Extrema_POnCurv2d P1, P2;
case GeomAbs_BezierCurve:
case GeomAbs_OtherCurve:
case GeomAbs_BSplineCurve: {
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
Standard_Real Period2 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(*((Adaptor2d_Curve2d*)myC))) Period2 = Extrema_Curve2dTool::Period(*((Adaptor2d_Curve2d*)myC));
- Results(Xtrem, C1, U11, U12, U21, U22, 2*M_PI,Period2);
+ Results(Xtrem, U11, U12, U21, U22, 2*M_PI,Period2);
}
break;
case GeomAbs_Line: {
break;
case GeomAbs_Ellipse: {
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Ellipse(C1), Extrema_Curve2dTool::Ellipse(*((Adaptor2d_Curve2d*)myC)));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22,2*M_PI, 2*M_PI);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22,2*M_PI, 2*M_PI);
}
break;
case GeomAbs_Parabola: {
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Ellipse(C1), Extrema_Curve2dTool::Parabola(*((Adaptor2d_Curve2d*)myC)));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 2*M_PI, 0.);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 2*M_PI, 0.);
}
break;
case GeomAbs_Hyperbola: {
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Ellipse(C1), Extrema_Curve2dTool::Hyperbola(*((Adaptor2d_Curve2d*)myC)));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 2*M_PI, 0.);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 2*M_PI, 0.);
}
break;
case GeomAbs_BezierCurve:
case GeomAbs_OtherCurve:
case GeomAbs_BSplineCurve: {
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
Standard_Real Period2 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(*((Adaptor2d_Curve2d*)myC))) Period2 = Extrema_Curve2dTool::Period(*((Adaptor2d_Curve2d*)myC));
- Results(Xtrem, C1, U11, U12, U21, U22, 2*M_PI,Period2);
+ Results(Xtrem, U11, U12, U21, U22, 2*M_PI,Period2);
}
break;
case GeomAbs_Line: {
case GeomAbs_Ellipse: {
//inverse = Standard_True;
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Ellipse(*((Adaptor2d_Curve2d*)myC)), Extrema_Curve2dTool::Parabola(C1));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 0., 2*M_PI);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 0., 2*M_PI);
}
break;
case GeomAbs_Parabola: {
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Parabola(C1), Extrema_Curve2dTool::Parabola(*((Adaptor2d_Curve2d*)myC)));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 0., 0.);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 0., 0.);
}
break;
case GeomAbs_Hyperbola: {
//inverse = Standard_True;
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Hyperbola(*((Adaptor2d_Curve2d*)myC)), Extrema_Curve2dTool::Parabola(C1));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 0., 0.);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 0., 0.);
}
break;
case GeomAbs_BezierCurve:
case GeomAbs_OtherCurve:
case GeomAbs_BSplineCurve: {
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
Standard_Real Period2 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(*((Adaptor2d_Curve2d*)myC))) Period2 = Extrema_Curve2dTool::Period(*((Adaptor2d_Curve2d*)myC));
- Results(Xtrem, C1, U11, U12, U21, U22, 0., Period2);
+ Results(Xtrem, U11, U12, U21, U22, 0., Period2);
}
break;
case GeomAbs_Line: {
case GeomAbs_Ellipse: {
//inverse = Standard_True;
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Ellipse(*((Adaptor2d_Curve2d*)myC)), Extrema_Curve2dTool::Hyperbola(C1));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 0., 2*M_PI );
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 0., 2*M_PI );
}
break;
case GeomAbs_Parabola: {
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Hyperbola(C1), Extrema_Curve2dTool::Parabola(*((Adaptor2d_Curve2d*)myC)));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 0., 0.);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 0., 0.);
}
break;
case GeomAbs_Hyperbola: {
//Extrema_ExtElC2d Xtrem(Extrema_Curve2dTool::Hyperbola(C1), Extrema_Curve2dTool::Hyperbola(*((Adaptor2d_Curve2d*)myC)));
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
- Results(Xtrem, C1, U11, U12, U21, U22, 0., 0.);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
+ Results(Xtrem, U11, U12, U21, U22, 0., 0.);
}
break;
case GeomAbs_OtherCurve:
case GeomAbs_BezierCurve:
case GeomAbs_BSplineCurve: {
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC)
- , NbU, NbV, mytolc1, mytolc2);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
Standard_Real Period2 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(*((Adaptor2d_Curve2d*)myC))) Period2 = Extrema_Curve2dTool::Period(*((Adaptor2d_Curve2d*)myC));
- Results(Xtrem, C1, U11, U12, U21, U22, 0., Period2);
+ Results(Xtrem, U11, U12, U21, U22, 0., Period2);
}
break;
case GeomAbs_Line: {
case GeomAbs_BezierCurve:
case GeomAbs_OtherCurve:
case GeomAbs_BSplineCurve: {
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
Standard_Real Period1 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(C1)) Period1 = Extrema_Curve2dTool::Period(C1);
Standard_Real Period2 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(*((Adaptor2d_Curve2d*)myC))) Period2 = Extrema_Curve2dTool::Period(*((Adaptor2d_Curve2d*)myC));
- Results(Xtrem, C1, U11, U12, U21, U22, Period1, Period2);
+ Results(Xtrem, U11, U12, U21, U22, Period1, Period2);
}
break;
case GeomAbs_BezierCurve:
case GeomAbs_OtherCurve:
case GeomAbs_BSplineCurve: {
- Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC),
- NbU, NbV, mytolc1, mytolc2);
+ Extrema_ECC2d Xtrem(C1, *((Adaptor2d_Curve2d*)myC));
+ Xtrem.Perform();
Standard_Real Period2 = 0.;
if (Extrema_Curve2dTool::IsPeriodic(*((Adaptor2d_Curve2d*)myC))) Period2 = Extrema_Curve2dTool::Period(*((Adaptor2d_Curve2d*)myC));
- Results(Xtrem, C1, U11, U12, U21, U22, 0., Period2);
+ Results(Xtrem, U11, U12, U21, U22, 0., Period2);
}
break;
case GeomAbs_Line: {
void Extrema_ExtCC2d::Results(const Extrema_ECC2d& AlgExt,
-// modified by NIZHNY-EAP Wed Feb 23 14:51:24 2000 ___BEGIN___
- const Adaptor2d_Curve2d& C1,
-// modified by NIZHNY-EAP Wed Feb 23 14:51:26 2000 ___END___
- const Standard_Real Ut11,
- const Standard_Real Ut12,
- const Standard_Real Ut21,
- const Standard_Real Ut22,
- const Standard_Real Period1,
- const Standard_Real Period2)
+ const Standard_Real Ut11,
+ const Standard_Real Ut12,
+ const Standard_Real Ut21,
+ const Standard_Real Ut22,
+ const Standard_Real Period1,
+ const Standard_Real Period2)
{
Standard_Integer i, NbExt;
Standard_Real Val, U, U2;
Extrema_POnCurv2d P1, P2;
-
+
myDone = AlgExt.IsDone();
- if (myDone) {
- if (!myIsPar) {
+ if (myDone)
+ {
+ if (!myIsPar)
+ {
NbExt = AlgExt.NbExt();
- for (i = 1; i <= NbExt; i++) {
- // Verification de la validite des parametres pour le cas trimme:
- AlgExt.Points(i, P1, P2);
- U = P1.Parameter();
- if (Period1 != 0.0) U = ElCLib::InPeriod(U,Ut11,Ut11+Period1);
- U2 = P2.Parameter();
- if (Period2 != 0.0) U2 = ElCLib::InPeriod(U2,Ut21,Ut21+Period2);
-
- if ((U >= Ut11 - Precision::PConfusion()) &&
- (U <= Ut12 + Precision::PConfusion()) &&
- (U2 >= Ut21 - Precision::PConfusion()) &&
- (U2 <= Ut22 + Precision::PConfusion())) {
-// modified by NIZHNY-EAP Thu Jan 27 16:40:55 2000 ___BEGIN___
- // to be sure that it's a real extrema
- gp_Pnt2d p;
- gp_Vec2d v1, v2;
- Extrema_Curve2dTool::D1(C1,U,p, v1);
- Extrema_Curve2dTool::D1(*((Adaptor2d_Curve2d*)myC),U2,p, v2);
- if (IsParallelDot(v1, v2, Precision::Angular()))
- {
- mynbext++;
- Val = AlgExt.SquareDistance(i);
- P1.SetValues(U, P1.Value());
- P2.SetValues(U2, P2.Value());
- mySqDist.Append(Val);
- mypoints.Append(P1);
- mypoints.Append(P2);
- }
-// modified by NIZHNY-EAP Thu Jan 27 16:41:00 2000 ___END___
- }
+ for (i = 1; i <= NbExt; i++)
+ {
+ // Verification de la validite des parametres pour le cas trimme:
+ AlgExt.Points(i, P1, P2);
+ U = P1.Parameter();
+ if (Period1 != 0.0)
+ U = ElCLib::InPeriod(U,Ut11,Ut11+Period1);
+ U2 = P2.Parameter();
+ if (Period2 != 0.0)
+ U2 = ElCLib::InPeriod(U2,Ut21,Ut21+Period2);
+
+ if ((U >= Ut11 - Precision::PConfusion()) &&
+ (U <= Ut12 + Precision::PConfusion()) &&
+ (U2 >= Ut21 - Precision::PConfusion()) &&
+ (U2 <= Ut22 + Precision::PConfusion()))
+ {
+ mynbext++;
+ Val = AlgExt.SquareDistance(i);
+ P1.SetValues(U, P1.Value());
+ P2.SetValues(U2, P2.Value());
+ mySqDist.Append(Val);
+ mypoints.Append(P1);
+ mypoints.Append(P2);
+ }
}
}
Tool1 as any; -- as ToolCurve(Curve1)
Curve2 as any;
Tool2 as any; -- as ToolCurve(Curve2)
- Cache as Transient from Standard; -- CurveCache from Extrema
ArrayOfPnt as Transient from Standard; -- as returned by Extrema_CurveCache::Points()
POnC as any;
Pnt as any;
---Purpose: It calculates all the distance between two curves.
-- These distances can be maximum or minimum.
-
+uses SequenceOfReal from TColStd,
+ Vector from math
+
raises InfiniteSolutions from StdFail,
NotDone from StdFail,
OutOfRange from Standard
-
-private class CCF instantiates FuncExtCC from Extrema (Curve1, Tool1,
- Curve2, Tool2,
- POnC, Pnt, Vec);
is
-- between Uinf and Usup for C1 and between Vinf and Vsup
-- for C2.
- Create (C1: Curve1; C2: Curve2; NbU,NbV: Integer; TolU,TolV: Real) returns GenExtCC;
+ Create (C1: Curve1; C2: Curve2) returns GenExtCC;
---Purpose: It calculates all the distances.
-- The function F(u,v)=distance(C1(u),C2(v)) has an
- -- extremum when gradient(f)=0. The algorithm searchs
- -- all the zeros.
- -- NbU,NbV are used to locate the close points to
- -- find the zeros. They must be great enough such that
- -- if there is N extrema, there will be N extrema
- -- between the segment and the grid.
- -- TolU and TolV are used to determinethe conditions
- -- to stop the iterations; at the iteration number n:
- -- (Un - Un-1) < TolU and (Vn - Vn-1) < TolV .
+ -- extremum when gradient(f)=0. The algorithm uses
+ -- Evtushenko's global optimization solver..
- Create (C1: Curve1; C2: Curve2; Uinf, Usup, Vinf, Vsup: Real;
- NbU,NbV: Integer; TolU,TolV: Real) returns GenExtCC;
+ Create (C1: Curve1; C2: Curve2; Uinf, Usup, Vinf, Vsup: Real) returns GenExtCC;
---Purpose: Calculates all the distances as above
-- between Uinf and Usup for C1 and between Vinf and Vsup
-- for C2.
- SetCurveCache (me: in out; theRank: Integer; theCache: Cache);
- ---Purpose:
-
- SetTolerance (me: in out; Tol: Real);
- ---Purpose:
+ SetParams(me: in out; C1: Curve1; C2: Curve2; Uinf, Usup, Vinf, Vsup: Real)
+ ---Purpose: Set params in case of empty constructor is usage.
+ is static;
Perform(me: in out) is static;
---Purpose: Performs calculations.
is static;
fields
- myF : CCF from Extrema;
+ myLowBorder : Vector from math;
+ myUppBorder : Vector from math;
+ mySolCount : Integer from Standard;
+ myPoints1 : SequenceOfReal from TColStd;
+ myPoints2 : SequenceOfReal from TColStd;
+ myC : Address from Standard [2];
myDone : Boolean;
- myCache: Cache [2];
end GenExtCC;
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
-#include <StdFail_NotDone.hxx>
-#include <math_FunctionSetRoot.hxx>
-#include <math_NewtonFunctionSetRoot.hxx>
-#include <TColStd_Array2OfInteger.hxx>
-#include <TColStd_Array2OfReal.hxx>
+#include <Extrema_GlobOptFuncCC.hxx>
+#include <math_GlobOptMin.hxx>
#include <Standard_NullObject.hxx>
#include <Standard_OutOfRange.hxx>
+#include <StdFail_NotDone.hxx>
+#include <TColStd_Array1OfReal.hxx>
#include <Precision.hxx>
-Extrema_GenExtCC::Extrema_GenExtCC () : myDone (Standard_False)
+Extrema_GenExtCC::Extrema_GenExtCC()
+: myLowBorder(1,2),
+ myUppBorder(1,2),
+ myDone(Standard_False)
{
}
Extrema_GenExtCC::Extrema_GenExtCC (const Curve1& C1,
- const Curve2& C2,
- const Standard_Integer NbU,
- const Standard_Integer NbV,
- const Standard_Real TolU,
- const Standard_Real TolV) : myF (C1,C2, Min(TolU, TolV)), myDone (Standard_False)
+ const Curve2& C2)
+: myLowBorder(1,2),
+ myUppBorder(1,2),
+ myDone(Standard_False)
{
- SetCurveCache (1, new Cache (C1, C1.FirstParameter(), C1.LastParameter(), NbU, Standard_True));
- SetCurveCache (2, new Cache (C2, C2.FirstParameter(), C2.LastParameter(), NbV, Standard_True));
- Perform();
+ myC[0] = (Standard_Address)&C1;
+ myC[1] = (Standard_Address)&C2;
+ myLowBorder(1) = C1.FirstParameter();
+ myLowBorder(2) = C2.FirstParameter();
+ myUppBorder(1) = C1.LastParameter();
+ myUppBorder(2) = C2.LastParameter();
}
Extrema_GenExtCC::Extrema_GenExtCC (const Curve1& C1,
- const Curve2& C2,
- const Standard_Real Uinf,
- const Standard_Real Usup,
- const Standard_Real Vinf,
- const Standard_Real Vsup,
- const Standard_Integer NbU,
- const Standard_Integer NbV,
- const Standard_Real /*TolU*/,
- const Standard_Real /*TolV*/) : myF (C1,C2), myDone (Standard_False)
-{
- SetCurveCache (1, new Cache (C1, Uinf, Usup, NbU, Standard_True));
- SetCurveCache (2, new Cache (C2, Vinf, Vsup, NbV, Standard_True));
- Perform();
-}
-
-void Extrema_GenExtCC::SetCurveCache (const Standard_Integer theRank,
- const Handle(Cache)& theCache)
+ const Curve2& C2,
+ const Standard_Real Uinf,
+ const Standard_Real Usup,
+ const Standard_Real Vinf,
+ const Standard_Real Vsup)
+: myLowBorder(1,2),
+ myUppBorder(1,2),
+ myDone(Standard_False)
{
- Standard_OutOfRange_Raise_if (theRank < 1 || theRank > 2, "Extrema_GenExtCC::SetCurveCache()")
- myF.SetCurve (theRank, *(Curve1*)theCache->CurvePtr());
- Standard_Integer anInd = theRank - 1;
- myCache[anInd] = theCache;
+ myC[0] = (Standard_Address)&C1;
+ myC[1] = (Standard_Address)&C2;
+ myLowBorder(1) = Uinf;
+ myLowBorder(2) = Vinf;
+ myUppBorder(1) = Usup;
+ myUppBorder(2) = Vsup;
}
-void Extrema_GenExtCC::SetTolerance (const Standard_Real Tol)
+void Extrema_GenExtCC::SetParams (const Curve1& C1,
+ const Curve2& C2,
+ const Standard_Real Uinf,
+ const Standard_Real Usup,
+ const Standard_Real Vinf,
+ const Standard_Real Vsup)
{
- myF.SetTolerance (Tol);
+ myC[0] = (Standard_Address)&C1;
+ myC[1] = (Standard_Address)&C2;
+ myLowBorder(1) = Uinf;
+ myLowBorder(2) = Vinf;
+ myUppBorder(1) = Usup;
+ myUppBorder(2) = Vsup;
}
-
//=============================================================================
void Extrema_GenExtCC::Perform ()
-/*-----------------------------------------------------------------------------
-Fonction:
- Recherche de toutes les distances extremales entre les courbes C1 et C2.
- a partir de 2 echantillonnages (NbU,NbV).
-
-Methode:
- L'algorithme part de l'hypothese que les echantillonnages sont suffisamment
- fins pour que, s'il existe N distances extremales entre les 2 courbes,
- alors il existe aussi N extrema entre les 2 ensembles de points.
- Ainsi, l'algorithme consiste a partir des extrema des echantillons
- pour trouver les extrema des courbes.
- Les extrema sont calcules par l'algorithme math_FunctionSetRoot avec les
- arguments suivants:
- - myF: Extrema_FuncExtCC cree a partir de C1 et C2,
- - UV: math_Vector dont les composantes sont les parametres des points de
- l'extremum sur les ensembles de points,
- - Tol: Min(TolU,TolV), (Prov.:math_FunctionSetRoot n'autorise pas un vecteur)
- - UVinf: math_Vector dont les composantes sont les bornes inferieures de u et
- v,
- - UVsup: math_Vector dont les composantes sont les bornes superieures de u et
- v.
-
-Traitement:
- a- Constitution du tableau des square distances (TbDist2(0,NbU+1,0,NbV+1)):
- Le tableau est volontairement etendu; les lignes 0 et NbU+1 et les
- colonnes 0 et NbV+1 seront initialisees a RealFirst() ou RealLast()
- pour simplifier les tests effectues dans l'etape b
- (on n'a pas besoin de tester si le point est sur une extremite).
- b- Calcul des extrema:
- On recherche d'abord les minima et ensuite les maxima. Ces 2 traitements
- se passent de facon similaire:
- b.a- Initialisations:
- - des 'bords' du tableau TbDist2 (a RealLast() dans le cas des minima
- et a RealLast() dans le cas des maxima),
- - du tableau TbSel(0,NbU+1,0,NbV+1) de selection des points pour un
- calcul d'extremum local (a 0). Lorsqu'un couple de points sera
- selectionne, il ne sera plus selectionnable, ainsi que les couples
- adjacents (8 au maximum).
- Les adresses correspondantes seront mises a 1.
- b.b- Calcul des minima (ou maxima):
- On boucle sur toutes les square distances du tableau TbDist2:
- - recherche d'un minimum (ou maximum) sur les echantillons,
- - calcul de l'extremum sur les courbes,
- - mise a jour du tableau TbSel.
------------------------------------------------------------------------------*/
{
myDone = Standard_False;
- const Handle(Cache)& aCache1 = myCache[0];
- const Handle(Cache)& aCache2 = myCache[1];
- Standard_NullObject_Raise_if ((aCache1.IsNull() || aCache2.IsNull()),
- "Extrema_GenExtCC::Perform()")
-
- Standard_Integer aNbU = aCache1->NbSamples(), aNbV = aCache2->NbSamples();
- Standard_OutOfRange_Raise_if ((aNbU < 2 ||aNbV < 2), "Extrema_GenExtCC::Perform()")
-
-/*
-a- Constitution du tableau des distances (TbDist2(0,NbU+1,0,NbV+1)):
- ---------------------------------------------------------------
-*/
-
- //ensure that caches have been calculated
- if (!aCache1->IsValid())
- aCache1->CalculatePoints();
- if (!aCache2->IsValid())
- aCache2->CalculatePoints();
-
-// Calcul des distances
- const Handle(ArrayOfPnt)& aPntArray1 = aCache1->Points();
- const Handle(ArrayOfPnt)& aPntArray2 = aCache2->Points();
- Standard_Integer NoU, NoV;
- TColStd_Array2OfReal TbDist2(0, aNbU + 1, 0, aNbV + 1);
- for (NoU = 1; NoU <= aNbU; NoU++) {
- const Pnt& P1 = aPntArray1->Value (NoU);
- for (NoV = 1; NoV <= aNbV; NoV++) {
- const Pnt& P2 = aPntArray2->Value (NoV);
- TbDist2(NoU,NoV) = P1.SquareDistance(P2);
+ Curve1 &C1 = *(Curve1*)myC[0];
+ Curve2 &C2 = *(Curve2*)myC[1];
+
+ Standard_Integer aNbInter[2];
+ aNbInter[0] = C1.NbIntervals(GeomAbs_C2);
+ aNbInter[1] = C2.NbIntervals(GeomAbs_C2);
+ TColStd_Array1OfReal anIntervals1(1, aNbInter[0] + 1);
+ TColStd_Array1OfReal anIntervals2(1, aNbInter[1] + 1);
+ C1.Intervals(anIntervals1, GeomAbs_C2);
+ C2.Intervals(anIntervals2, GeomAbs_C2);
+
+ math_MultipleVarFunction *aFunc = new Extrema_GlobOptFuncCCC2(C1, C2);
+ math_GlobOptMin aFinder(aFunc, myLowBorder, myUppBorder);
+
+ Standard_Integer i,j,k;
+ math_Vector aFirstBorderInterval(1,2);
+ math_Vector aSecondBorderInterval(1,2);
+ Standard_Real aF = RealLast(); // best functional value
+ Standard_Real aCurrF = RealLast(); // current functional value computed on current interval
+ for(i = 1; i <= aNbInter[0]; i++)
+ {
+ for(j = 1; j <= aNbInter[1]; j++)
+ {
+ aFirstBorderInterval(1) = anIntervals1(i);
+ aFirstBorderInterval(2) = anIntervals2(j);
+ aSecondBorderInterval(1) = anIntervals1(i + 1);
+ aSecondBorderInterval(2) = anIntervals2(j + 1);
+
+ aFinder.SetLocalParams(aFirstBorderInterval, aSecondBorderInterval);
+ aFinder.Perform();
+
+ aCurrF = aFinder.GetF();
+ if (aCurrF < aF + Precision::Confusion())
+ {
+ if (aCurrF > Abs(aF - Precision::Confusion()) || (aCurrF < 1.0e-15 && aF < 1.0e-15))
+ {
+ Standard_Integer myTmpSolCount = aFinder.NbExtrema();
+ math_Vector sol(1,2);
+ for(k = 1; k <= myTmpSolCount; k++)
+ {
+ aFinder.Points(k, sol);
+ myPoints1.Append(sol(1));
+ myPoints2.Append(sol(2));
+ }
+ mySolCount += myTmpSolCount;
+ } // if (aCurrF > aF - Precision::Confusion())
+ else
+ {
+ aF = aCurrF;
+ mySolCount = aFinder.NbExtrema();
+ myPoints1.Clear();
+ myPoints2.Clear();
+ math_Vector sol(1,2);
+ for(k = 1; k <= mySolCount; k++)
+ {
+ aFinder.Points(k, sol);
+ myPoints1.Append(sol(1));
+ myPoints2.Append(sol(2));
+ }
+ } // else
+ } //if (aCurrF < aF + Precision::Confusion())
}
}
-/*
-b- Calcul des minima:
- -----------------
- b.a) Initialisations:
-*/
-// - generales
- math_Vector Tol(1, 2);
- Tol(1) = myF.Tolerance();
- Tol(2) = myF.Tolerance();
- math_Vector UV(1,2), UVinf(1,2), UVsup(1,2);
- UVinf(1) = aCache1->TrimFirstParameter();
- UVinf(2) = aCache2->TrimFirstParameter();
- UVsup(1) = aCache1->TrimLastParameter();
- UVsup(2) = aCache2->TrimLastParameter();
-
- myF.SubIntervalInitialize(UVinf,UVsup);
-
-// - des 'bords' du tableau TbDist2
- for (NoV = 0; NoV <= aNbV+1; NoV++) {
- TbDist2(0,NoV) = RealLast();
- TbDist2(aNbU+1,NoV) = RealLast();
- }
- for (NoU = 1; NoU <= aNbU; NoU++) {
- TbDist2(NoU,0) = RealLast();
- TbDist2(NoU,aNbV+1) = RealLast();
- }
-
-// - du tableau TbSel(0,aNbU+1,0,aNbV+1) de selection des points
- TColStd_Array2OfInteger TbSel(0,aNbU+1,0,aNbV+1);
- TbSel.Init(0);
-
-/*
- b.b) Calcul des minima:
-*/
-// - recherche d un minimum sur la grille
- Standard_Integer Nu, Nv;
- Standard_Real Dist2;
- const Handle(TColStd_HArray1OfReal)& aParamArray1 = aCache1->Parameters();
- const Handle(TColStd_HArray1OfReal)& aParamArray2 = aCache2->Parameters();
- for (NoU = 1; NoU <= aNbU; NoU++) {
- for (NoV = 1; NoV <= aNbV; NoV++) {
- if (TbSel(NoU,NoV) == 0) {
- Dist2 = TbDist2(NoU,NoV);
- if ((TbDist2(NoU-1,NoV-1) >= Dist2) &&
- (TbDist2(NoU-1,NoV ) >= Dist2) &&
- (TbDist2(NoU-1,NoV+1) >= Dist2) &&
- (TbDist2(NoU ,NoV-1) >= Dist2) &&
- (TbDist2(NoU ,NoV+1) >= Dist2) &&
- (TbDist2(NoU+1,NoV-1) >= Dist2) &&
- (TbDist2(NoU+1,NoV ) >= Dist2) &&
- (TbDist2(NoU+1,NoV+1) >= Dist2)) {
-
-// - calcul de l extremum sur la surface:
-
- UV(1) = aParamArray1->Value(NoU);
- UV(2) = aParamArray2->Value(NoV);
-
- math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup);
-
-
-// - mise a jour du tableau TbSel
- for (Nu = NoU-1; Nu <= NoU+1; Nu++) {
- for (Nv = NoV-1; Nv <= NoV+1; Nv++) {
- TbSel(Nu,Nv) = 1;
- }
- }
- }
- } // if (TbSel(NoU,NoV)
- } // for (NoV = 1; ...
- } // for (NoU = 1; ...
-/*
-c- Calcul des maxima:
- -----------------
- c.a) Initialisations:
-*/
-// - des 'bords' du tableau TbDist2
- for (NoV = 0; NoV <= aNbV+1; NoV++) {
- TbDist2(0,NoV) = RealFirst();
- TbDist2(aNbU+1,NoV) = RealFirst();
- }
- for (NoU = 1; NoU <= aNbU; NoU++) {
- TbDist2(NoU,0) = RealFirst();
- TbDist2(NoU,aNbV+1) = RealFirst();
- }
-
-// - du tableau TbSel(0,aNbU+1,0,aNbV+1) de selection des points
- TbSel.Init(0);
- /*for (NoU = 0; NoU <= aNbU+1; NoU++) {
- for (NoV = 0; NoV <= aNbV+1; NoV++) {
- TbSel(NoU,NoV) = 0;
+ // Clear solutions clusters if it is necessary.
+ for(i = 1; i <= mySolCount - 1; i++)
+ {
+ for(j = i + 1; j <= mySolCount; j++)
+ {
+ if ((myPoints1(i) - myPoints1(j)) < (myUppBorder(1) - myLowBorder(1)) * Precision::Confusion() &&
+ (myPoints2(i) - myPoints2(j)) < (myUppBorder(2) - myLowBorder(2)) * Precision::Confusion())
+ {
+ // Points with indexes i and j is in same cluster, delete j point from it.
+ myPoints1.Remove(j);
+ myPoints2.Remove(j);
+ j--;
+ mySolCount--;
+ }
}
- }*/
-/*
- c.b) Calcul des maxima:
-*/
-// - recherche d un maximum sur la grille
- for (NoU = 1; NoU <= aNbU; NoU++) {
- for (NoV = 1; NoV <= aNbV; NoV++) {
- if (TbSel(NoU,NoV) == 0) {
- Dist2 = TbDist2(NoU,NoV);
- if ((TbDist2(NoU-1,NoV-1) <= Dist2) &&
- (TbDist2(NoU-1,NoV ) <= Dist2) &&
- (TbDist2(NoU-1,NoV+1) <= Dist2) &&
- (TbDist2(NoU ,NoV-1) <= Dist2) &&
- (TbDist2(NoU ,NoV+1) <= Dist2) &&
- (TbDist2(NoU+1,NoV-1) <= Dist2) &&
- (TbDist2(NoU+1,NoV ) <= Dist2) &&
- (TbDist2(NoU+1,NoV+1) <= Dist2)) {
-
-// - calcul de l extremum sur la surface:
-
- UV(1) = aParamArray1->Value(NoU);
- UV(2) = aParamArray2->Value(NoV);
+ }
- math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup);
-
-// - mise a jour du tableau TbSel
- for (Nu = NoU-1; Nu <= NoU+1; Nu++) {
- for (Nv = NoV-1; Nv <= NoV+1; Nv++) {
- TbSel(Nu,Nv) = 1;
- }
- }
- }
- } // if (TbSel(NoU,NoV))
- } // for (NoV = 1; ...)
- } // for (NoU = 1; ...)
+ delete aFunc;
myDone = Standard_True;
}
//=============================================================================
-Standard_Boolean Extrema_GenExtCC::IsDone () const { return myDone; }
+Standard_Boolean Extrema_GenExtCC::IsDone () const
+{
+ return myDone;
+}
//=============================================================================
Standard_Integer Extrema_GenExtCC::NbExt () const
{
StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::NbExt()")
- return myF.NbExt();
+
+ return mySolCount;
}
//=============================================================================
{
StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
- return myF.SquareDistance(N);
+
+ return Tool1::Value(*((Curve1*)myC[0]), myPoints1(N)).SquareDistance(Tool2::Value(*((Curve2*)myC[1]), myPoints2(N)));
}
//=============================================================================
void Extrema_GenExtCC::Points (const Standard_Integer N,
- POnC& P1, POnC& P2) const
+ POnC& P1,
+ POnC& P2) const
{
- StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
- Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
- myF.Points(N,P1,P2);
-}
+ StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::Points()")
+ Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::Points()")
+
+ P1.SetValues(myPoints1(N), Tool1::Value(*((Curve1*)myC[0]), myPoints1(N)));
+ P2.SetValues(myPoints2(N), Tool2::Value(*((Curve2*)myC[1]), myPoints2(N)));
+}
\ No newline at end of file
--- /dev/null
+// Created on: 2014-01-20
+// Created by: Alexaner Malyshev
+// Copyright (c) 2014-2014 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 <Extrema_GlobOptFuncCC.hxx>
+
+#include <gp_Pnt.hxx>
+#include <gp_Pnt2d.hxx>
+#include <gp_Vec.hxx>
+#include <gp_Vec2d.hxx>
+#include <math_Vector.hxx>
+#include <Standard_Integer.hxx>
+#include <Standard_OutOfRange.hxx>
+
+static Standard_Integer _NbVariables()
+{
+ return 2;
+}
+
+// 3d _Value
+static Standard_Boolean _Value(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2,
+ const math_Vector& X,
+ Standard_Real& F)
+{
+ Standard_Real u = X(1);
+ Standard_Real v = X(2);
+
+ if (u < C1.FirstParameter() ||
+ u > C1.LastParameter() ||
+ v < C2.FirstParameter() ||
+ v > C2.LastParameter())
+ {
+ return Standard_False;
+ }
+
+ F = C2.Value(v).Distance(C1.Value(u));
+ return Standard_True;
+}
+
+// 2d _Value
+static Standard_Boolean _Value(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2,
+ const math_Vector& X,
+ Standard_Real& F)
+{
+ Standard_Real u = X(1);
+ Standard_Real v = X(2);
+
+ if (u < C1.FirstParameter() ||
+ u > C1.LastParameter() ||
+ v < C2.FirstParameter() ||
+ v > C2.LastParameter())
+ {
+ return Standard_False;
+ }
+
+ F = C2.Value(v).Distance(C1.Value(u));
+ return Standard_True;
+}
+
+//! F = (x2(v) - x1(u))^2 + (y2(v) - y1(u))^2 + (z2(v) - z1(u))^2
+
+// 3d _Gradient
+static Standard_Boolean _Gradient(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2,
+ const math_Vector& X,
+ math_Vector& G)
+{
+ gp_Pnt C1D0, C2D0;
+ gp_Vec C1D1, C2D1;
+
+ if(X(1) < C1.FirstParameter() ||
+ X(1) > C1.LastParameter() ||
+ X(2) < C2.FirstParameter() ||
+ X(2) > C2.LastParameter())
+ {
+ return Standard_False;
+ }
+
+ 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();
+ return Standard_True;
+}
+
+// 2d _Graient
+static Standard_Boolean _Gradient(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2,
+ const math_Vector& X,
+ math_Vector& G)
+{
+ gp_Pnt2d C1D0, C2D0;
+ gp_Vec2d C1D1, C2D1;
+
+ if(X(1) < C1.FirstParameter() ||
+ X(1) > C1.LastParameter() ||
+ X(2) < C2.FirstParameter() ||
+ X(2) > C2.LastParameter())
+ {
+ return Standard_False;
+ }
+
+ 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();
+ G(2) = (C2D0.X() - C1D0.X()) * C2D1.X()
+ + (C2D0.Y() - C1D0.Y()) * C2D1.Y();
+ return Standard_True;
+}
+
+// 3d _Hessian
+static Standard_Boolean _Hessian (const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2,
+ const math_Vector& X,
+ math_Matrix & H)
+{
+ gp_Pnt C1D0, C2D0;
+ gp_Vec C1D1, C2D1;
+ gp_Vec C1D2, C2D2;
+
+ if(X(1) < C1.FirstParameter() ||
+ X(1) > C1.LastParameter() ||
+ X(2) < C2.FirstParameter() ||
+ X(2) > C2.LastParameter())
+ {
+ return Standard_False;
+ }
+
+ C1.D2(X(1), C1D0, C1D1, C1D2);
+ C2.D2(X(2), C2D0, C2D1, C2D2);
+
+ H(1, 1) = C1D1.X() * C1D1.X()
+ + C1D1.Y() * C1D1.Y()
+ + C1D1.Z() * C1D1.Z()
+ - (C2D0.X() - C1D0.X()) * C1D2.X()
+ - (C2D0.Y() - C1D0.Y()) * C1D2.Y()
+ - (C2D0.Z() - C1D0.Z()) * C1D2.Z();
+
+ H(1, 2) = - C2D1.X() * C1D1.X()
+ - C2D1.Y() * C1D1.Y()
+ - C2D1.Z() * C1D1.Z();
+
+ H(2,1) = H(1,2);
+
+ H(2,2) = C2D1.X() * C2D1.X()
+ + C2D1.Y() * C2D1.Y()
+ + C2D1.Z() * C2D1.Z()
+ + (C2D0.X() - C1D0.X()) * C2D2.X()
+ + (C2D0.Y() - C1D0.Y()) * C2D2.Y()
+ + (C2D0.Z() - C1D0.Z()) * C2D2.Z();
+ return Standard_True;
+}
+
+// 2d _Hessian
+static Standard_Boolean _Hessian (const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2,
+ const math_Vector& X,
+ math_Matrix & H)
+{
+ gp_Pnt2d C1D0, C2D0;
+ gp_Vec2d C1D1, C2D1;
+ gp_Vec2d C1D2, C2D2;
+
+ if(X(1) < C1.FirstParameter() ||
+ X(1) > C1.LastParameter() ||
+ X(2) < C2.FirstParameter() ||
+ X(2) > C2.LastParameter())
+ {
+ return Standard_False;
+ }
+
+ C1.D2(X(1), C1D0, C1D1, C1D2);
+ C2.D2(X(2), C2D0, C2D1, C2D2);
+
+ H(1, 1) = C1D1.X() * C1D1.X()
+ + C1D1.Y() * C1D1.Y()
+ - (C2D0.X() - C1D0.X()) * C1D2.X()
+ - (C2D0.Y() - C1D0.Y()) * C1D2.Y();
+
+ H(1, 2) = - C2D1.X() * C1D1.X()
+ - C2D1.Y() * C1D1.Y();
+
+ H(2,1) = H(1,2);
+
+ H(2,2) = C2D1.X() * C2D1.X()
+ + C2D1.Y() * C2D1.Y()
+ + (C2D0.X() - C1D0.X()) * C2D2.X()
+ + (C2D0.Y() - C1D0.Y()) * C2D2.Y();
+ return Standard_True;
+}
+
+// C0
+
+//=======================================================================
+//function : Extrema_GlobOptFuncCCC0
+//purpose : Constructor
+//=======================================================================
+Extrema_GlobOptFuncCCC0::Extrema_GlobOptFuncCCC0(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2)
+: myC1_3d(&C1),
+ myC2_3d(&C2)
+{
+ myType = 1;
+}
+
+//=======================================================================
+//function : Extrema_GlobOptFuncCCC0
+//purpose : Constructor
+//=======================================================================
+Extrema_GlobOptFuncCCC0::Extrema_GlobOptFuncCCC0(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2)
+: myC1_2d(&C1),
+ myC2_2d(&C2)
+{
+ myType = 2;
+}
+
+
+//=======================================================================
+//function : NbVariables
+//purpose :
+//=======================================================================
+Standard_Integer Extrema_GlobOptFuncCCC0::NbVariables() const
+{
+ return _NbVariables();
+}
+
+//=======================================================================
+//function : Value
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC0::Value(const math_Vector& X,Standard_Real& F)
+{
+ if (myType == 1)
+ return _Value(*myC1_3d, *myC2_3d, X, F);
+ else
+ return _Value(*myC1_2d, *myC2_2d, X, F);
+}
+
+// C1
+
+//=======================================================================
+//function : Extrema_GlobOptFuncCCC1
+//purpose : Constructor
+//=======================================================================
+Extrema_GlobOptFuncCCC1::Extrema_GlobOptFuncCCC1(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2)
+: myC1_3d(&C1),
+ myC2_3d(&C2)
+{
+ myType = 1;
+}
+
+//=======================================================================
+//function : Extrema_GlobOptFuncCCC1
+//purpose : Constructor
+//=======================================================================
+Extrema_GlobOptFuncCCC1::Extrema_GlobOptFuncCCC1(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2)
+: myC1_2d(&C1),
+ myC2_2d(&C2)
+{
+ myType = 2;
+}
+
+//=======================================================================
+//function : NbVariables
+//purpose :
+//=======================================================================
+Standard_Integer Extrema_GlobOptFuncCCC1::NbVariables() const
+{
+ return _NbVariables();
+}
+
+//=======================================================================
+//function : Value
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC1::Value(const math_Vector& X,Standard_Real& F)
+{
+ if (myType == 1)
+ return _Value(*myC1_3d, *myC2_3d, X, F);
+ else
+ return _Value(*myC1_2d, *myC2_2d, X, F);
+}
+
+//=======================================================================
+//function : Gradient
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC1::Gradient(const math_Vector& X,math_Vector& G)
+{
+ if (myType == 1)
+ return _Gradient(*myC1_3d, *myC2_3d, X, G);
+ else
+ return _Gradient(*myC1_2d, *myC2_2d, X, G);
+}
+
+//=======================================================================
+//function : Values
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC1::Values(const math_Vector& X,Standard_Real& F,math_Vector& G)
+{
+ return (Value(X, F) && Gradient(X, G));
+}
+
+// C2
+
+//=======================================================================
+//function : Extrema_GlobOptFuncCCC2
+//purpose : Constructor
+//=======================================================================
+Extrema_GlobOptFuncCCC2::Extrema_GlobOptFuncCCC2(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2)
+: myC1_3d(&C1),
+ myC2_3d(&C2)
+{
+ myType = 1;
+}
+
+//=======================================================================
+//function : Extrema_GlobOptFuncCCC2
+//purpose : Constructor
+//=======================================================================
+Extrema_GlobOptFuncCCC2::Extrema_GlobOptFuncCCC2(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2)
+: myC1_2d(&C1),
+ myC2_2d(&C2)
+{
+ myType = 2;
+}
+
+//=======================================================================
+//function : NbVariables
+//purpose :
+//=======================================================================
+Standard_Integer Extrema_GlobOptFuncCCC2::NbVariables() const
+{
+ return _NbVariables();
+}
+
+//=======================================================================
+//function : Value
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC2::Value(const math_Vector& X,Standard_Real& F)
+{
+ if (myType == 1)
+ return _Value(*myC1_3d, *myC2_3d, X, F);
+ else
+ return _Value(*myC1_2d, *myC2_2d, X, F);
+}
+
+//=======================================================================
+//function : Gradient
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC2::Gradient(const math_Vector& X,math_Vector& G)
+{
+ if (myType == 1)
+ return _Gradient(*myC1_3d, *myC2_3d, X, G);
+ else
+ return _Gradient(*myC1_2d, *myC2_2d, X, G);
+}
+
+//=======================================================================
+//function : Values
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC2::Values(const math_Vector& X,Standard_Real& F,math_Vector& G)
+{
+ return (Value(X, F) && Gradient(X, G));
+}
+
+//=======================================================================
+//function : Values
+//purpose :
+//=======================================================================
+Standard_Boolean Extrema_GlobOptFuncCCC2::Values(const math_Vector& X,Standard_Real& F,math_Vector& G,math_Matrix& H)
+{
+ Standard_Boolean isHessianComputed = Standard_False;
+ if (myType == 1)
+ isHessianComputed = _Hessian(*myC1_3d, *myC2_3d, X, H);
+ else
+ isHessianComputed = _Hessian(*myC1_2d, *myC2_2d, X, H);
+
+
+ return (Value(X, F) && Gradient(X, G) && isHessianComputed);
+}
--- /dev/null
+// Created on: 2014-01-20
+// Created by: Alexaner Malyshev
+// Copyright (c) 2014-2014 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_GlobOptFuncCC_HeaderFile
+#define _Extrema_GlobOptFuncCC_HeaderFile
+
+#include <Adaptor2d_Curve2d.hxx>
+#include <Adaptor3d_Curve.hxx>
+#include <math_Matrix.hxx>
+#include <math_Vector.hxx>
+#include <math_MultipleVarFunction.hxx>
+#include <math_MultipleVarFunctionWithGradient.hxx>
+#include <math_MultipleVarFunctionWithHessian.hxx>
+
+//! This class implements function which operates with Eucluidean distance <br>
+//! between 2 curves in case of C1 and C2 continuity is C0.
+class Extrema_GlobOptFuncCCC0 : public math_MultipleVarFunction
+{
+public:
+
+ Standard_EXPORT Extrema_GlobOptFuncCCC0(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2);
+
+ Standard_EXPORT Extrema_GlobOptFuncCCC0(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2);
+
+
+
+ Standard_EXPORT virtual Standard_Integer NbVariables() const;
+
+ Standard_EXPORT virtual Standard_Boolean Value(const math_Vector& X,Standard_Real& F);
+
+
+private:
+
+ Extrema_GlobOptFuncCCC0 & operator = (const Extrema_GlobOptFuncCCC0 & theOther);
+
+ const Adaptor3d_Curve *myC1_3d, *myC2_3d;
+ const Adaptor2d_Curve2d *myC1_2d, *myC2_2d;
+ Standard_Integer myType;
+};
+
+
+//! This class implements function which operates with Eucluidean distance <br>
+//! between 2 curves in case of C1 and C2 continuity is C1.
+class Extrema_GlobOptFuncCCC1 : public math_MultipleVarFunctionWithGradient
+{
+public:
+
+ Standard_EXPORT Extrema_GlobOptFuncCCC1(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2);
+
+ Standard_EXPORT Extrema_GlobOptFuncCCC1(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2);
+
+ Standard_EXPORT virtual Standard_Integer NbVariables() const;
+
+ Standard_EXPORT virtual Standard_Boolean Value(const math_Vector& X,Standard_Real& F);
+
+ Standard_EXPORT virtual Standard_Boolean Gradient(const math_Vector& X,math_Vector& G);
+
+ Standard_EXPORT virtual Standard_Boolean Values(const math_Vector& X,Standard_Real& F,math_Vector& G);
+
+
+private:
+
+ Extrema_GlobOptFuncCCC1 & operator = (const Extrema_GlobOptFuncCCC1 & theOther);
+
+ const Adaptor3d_Curve *myC1_3d, *myC2_3d;
+ const Adaptor2d_Curve2d *myC1_2d, *myC2_2d;
+ Standard_Integer myType;
+};
+
+
+//! This class implements function which operates with Eucluidean distance <br>
+//! between 2 curves in case of C1 and C2 continuity is C2 or more.
+class Extrema_GlobOptFuncCCC2 : public math_MultipleVarFunctionWithHessian
+{
+public:
+
+ Standard_EXPORT Extrema_GlobOptFuncCCC2(const Adaptor3d_Curve& C1,
+ const Adaptor3d_Curve& C2);
+
+ Standard_EXPORT Extrema_GlobOptFuncCCC2(const Adaptor2d_Curve2d& C1,
+ const Adaptor2d_Curve2d& C2);
+
+ Standard_EXPORT virtual Standard_Integer NbVariables() const;
+
+ Standard_EXPORT virtual Standard_Boolean Value(const math_Vector& X,Standard_Real& F);
+
+ Standard_EXPORT virtual Standard_Boolean Gradient(const math_Vector& X,math_Vector& G);
+
+ Standard_EXPORT virtual Standard_Boolean Values(const math_Vector& X,Standard_Real& F,math_Vector& G);
+
+ Standard_EXPORT virtual Standard_Boolean Values(const math_Vector& X,Standard_Real& F,math_Vector& G,math_Matrix& H);
+
+
+private:
+
+ Extrema_GlobOptFuncCCC2 & operator = (const Extrema_GlobOptFuncCCC2 & theOther);
+
+ const Adaptor3d_Curve *myC1_3d, *myC2_3d;
+ const Adaptor2d_Curve2d *myC1_2d, *myC2_2d;
+ Standard_Integer myType;
+};
+
+#endif
Vec2d from gp,
HArray1OfPnt2d from TColgp,
Curve2d from Adaptor2d,
- Curve2dTool from Extrema,
- LCCache2d from Extrema
+ Curve2dTool from Extrema
raises DomainError from Standard,
NotDone from StdFail
Extrema_HUBTreeOfSphere.hxx
+Extrema_GlobOptFuncCC.hxx
+Extrema_GlobOptFuncCC.cxx
}
}
}
- if (!IntersFound) //intersection is inside "theEdge" => split
+ if (!IntersFound && aDist <= Precision::Confusion()) //intersection is inside "theEdge" => split
{
gp_Pnt aPoint = aCurve->Value(anIntPar);
if (aPoint.Distance(thePnt[0]) > BRep_Tool::Tolerance(theVertices[0]) &&
math_Vector.cxx
math_IntegerVector.hxx
math_IntegerVector.cxx
-math_SingleTab.hxx
\ No newline at end of file
+math_SingleTab.hxx
+math_GlobOptMin.hxx
+math_GlobOptMin.cxx
\ No newline at end of file
deferred class MultipleVarFunctionWithHessian;
deferred class FunctionSet;
deferred class FunctionSetWithDerivatives;
+ imported GlobOptMin;
class IntegerRandom;
class Gauss;
class GaussLeastSquare;
--- /dev/null
+// Created on: 2014-01-20
+// Created by: Alexaner Malyshev
+// Copyright (c) 2014-2014 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 <math_GlobOptMin.hxx>
+
+#include <math_BFGS.hxx>
+#include <math_Matrix.hxx>
+#include <math_MultipleVarFunctionWithGradient.hxx>
+#include <math_MultipleVarFunctionWithHessian.hxx>
+#include <math_NewtonMinimum.hxx>
+#include <math_Powell.hxx>
+#include <Precision.hxx>
+#include <Standard_Integer.hxx>
+#include <Standard_Real.hxx>
+
+const Handle(Standard_Type)& STANDARD_TYPE(math_GlobOptMin)
+{
+ static Handle(Standard_Type) _atype = new Standard_Type ("math_GlobOptMin", sizeof (math_GlobOptMin));
+ return _atype;
+}
+
+//=======================================================================
+//function : math_GlobOptMin
+//purpose : Constructor
+//=======================================================================
+math_GlobOptMin::math_GlobOptMin(math_MultipleVarFunction* theFunc,
+ const math_Vector& theA,
+ const math_Vector& theB,
+ Standard_Real theC)
+: myN(theFunc->NbVariables()),
+ myA(1, myN),
+ myB(1, myN),
+ myGlobA(1, myN),
+ myGlobB(1, myN),
+ myX(1, myN),
+ myTmp(1, myN),
+ myV(1, myN)
+{
+ Standard_Integer i;
+
+ myFunc = theFunc;
+ myC = theC;
+ myZ = -1;
+ mySolCount = 0;
+
+ for(i = 1; i <= myN; i++)
+ {
+ myGlobA(i) = theA(i);
+ myGlobB(i) = theB(i);
+
+ myA(i) = theA(i);
+ myB(i) = theB(i);
+ }
+
+ myDone = Standard_False;
+}
+
+//=======================================================================
+//function : SetGlobalParams
+//purpose : Set params without memory allocation.
+//=======================================================================
+void math_GlobOptMin::SetGlobalParams(math_MultipleVarFunction* theFunc,
+ const math_Vector& theA,
+ const math_Vector& theB,
+ Standard_Real theC)
+{
+ Standard_Integer i;
+
+ myFunc = theFunc;
+ myC = theC;
+ myZ = -1;
+ mySolCount = 0;
+
+ for(i = 1; i <= myN; i++)
+ {
+ myGlobA(i) = theA(i);
+ myGlobB(i) = theB(i);
+
+ myA(i) = theA(i);
+ myB(i) = theB(i);
+ }
+
+ myDone = Standard_False;
+}
+
+//=======================================================================
+//function : SetLocalParams
+//purpose : Set params without memory allocation.
+//=======================================================================
+void math_GlobOptMin::SetLocalParams(const math_Vector& theLocalA,
+ const math_Vector& theLocalB)
+{
+ Standard_Integer i;
+
+ myZ = -1;
+ mySolCount = 0;
+
+ for(i = 1; i <= myN; i++)
+ {
+ myA(i) = theLocalA(i);
+ myB(i) = theLocalB(i);
+ }
+
+ myDone = Standard_False;
+}
+
+//=======================================================================
+//function : ~math_GlobOptMin
+//purpose :
+//=======================================================================
+math_GlobOptMin::~math_GlobOptMin()
+{
+}
+
+//=======================================================================
+//function : Perform
+//purpose : Compute Global extremum point
+//=======================================================================
+// In this algo indexes started from 1, not from 0.
+void math_GlobOptMin::Perform()
+{
+ Standard_Integer i;
+
+ // Compute parameters range
+ Standard_Real minLength = RealLast();
+ Standard_Real maxLength = RealFirst();
+ for(i = 1; i <= myN; i++)
+ {
+ Standard_Real currentLength = myB(i) - myA(i);
+ if (currentLength < minLength)
+ minLength = currentLength;
+ if (currentLength > maxLength)
+ maxLength = currentLength;
+ }
+
+ Standard_Real myTol = 1e-2;
+
+ myE1 = minLength * myTol / myC;
+ myE2 = 2.0 * myTol / myC * maxLength;
+ myE3 = - maxLength * myTol / 4;
+
+ // Compure start point.
+ math_Vector aPnt(1,2);
+ for(i = 1; i <= myN; i++)
+ {
+ Standard_Real currCentral = (myA(i) + myB(i)) / 2.0;
+ aPnt(i) = currCentral;
+ myY.Append(currCentral);
+ }
+
+ myFunc->Value(aPnt, myF);
+ mySolCount++;
+
+ computeGlobalExtremum(myN);
+
+ myDone = Standard_True;
+
+}
+
+//=======================================================================
+//function : computeLocalExtremum
+//purpose :
+//=======================================================================
+Standard_Boolean math_GlobOptMin::computeLocalExtremum(const math_Vector& thePnt,
+ Standard_Real& theVal,
+ math_Vector& theOutPnt)
+{
+ Standard_Integer i;
+
+ //Newton method
+ if (dynamic_cast<math_MultipleVarFunctionWithHessian*>(myFunc))
+ {
+ math_MultipleVarFunctionWithHessian* myTmp =
+ dynamic_cast<math_MultipleVarFunctionWithHessian*> (myFunc);
+
+ math_NewtonMinimum newtonMinimum(*myTmp, thePnt);
+ if (newtonMinimum.IsDone())
+ {
+ newtonMinimum.Location(theOutPnt);
+ theVal = newtonMinimum.Minimum();
+ }
+ else return Standard_False;
+ } else
+
+ // BFGS method used.
+ if (dynamic_cast<math_MultipleVarFunctionWithGradient*>(myFunc))
+ {
+ math_MultipleVarFunctionWithGradient* myTmp =
+ dynamic_cast<math_MultipleVarFunctionWithGradient*> (myFunc);
+ math_BFGS bfgs(*myTmp, thePnt);
+ if (bfgs.IsDone())
+ {
+ bfgs.Location(theOutPnt);
+ theVal = bfgs.Minimum();
+ }
+ else return Standard_False;
+ } else
+
+ // Powell method used.
+ if (dynamic_cast<math_MultipleVarFunction*>(myFunc))
+ {
+ math_Matrix m(1, myN, 1, myN, 0.0);
+ for(i = 1; i <= myN; i++)
+ m(1, 1) = 1.0;
+
+ math_Powell powell(*myFunc, thePnt, m, 1e-10);
+
+ if (powell.IsDone())
+ {
+ powell.Location(theOutPnt);
+ theVal = powell.Minimum();
+ }
+ else return Standard_False;
+ }
+
+ if (isInside(theOutPnt))
+ return Standard_True;
+ else
+ return Standard_False;
+}
+
+//=======================================================================
+//function : ComputeGlobalExtremum
+//purpose :
+//=======================================================================
+void math_GlobOptMin::computeGlobalExtremum(Standard_Integer j)
+{
+ Standard_Integer i;
+ Standard_Real d; // Functional in moved point.
+ Standard_Real val = RealLast(); // Local extrema computed in moved point.
+ Standard_Real stepBestValue = RealLast();
+ math_Vector stepBestPoint(1,2);
+ Standard_Boolean isInside = Standard_False;
+ Standard_Real r;
+
+ for(myX(j) = myA(j) + myE1; myX(j) < myB(j) + myE1; myX(j) += myV(j))
+ {
+ if (myX(j) > myB(j))
+ myX(j) = myB(j);
+
+ if (j == 1)
+ {
+ isInside = Standard_False;
+ myFunc->Value(myX, d);
+ r = (d - myF) * myZ;
+ if(r > myE3)
+ {
+ isInside = computeLocalExtremum(myX, val, myTmp);
+ }
+ stepBestValue = (isInside && (val < d))? val : d;
+ stepBestPoint = (isInside && (val < d))? myTmp : myX;
+
+ // Solutions are close to each other.
+ if (Abs(stepBestValue - myF) < Precision::Confusion() * 0.01)
+ {
+ if (!isStored(stepBestPoint))
+ {
+ if ((stepBestValue - myF) * myZ > 0.0)
+ myF = stepBestValue;
+ for(i = 1; i <= myN; i++)
+ myY.Append(stepBestPoint(i));
+ mySolCount++;
+ }
+ }
+
+ // New best solution.
+ if ((stepBestValue - myF) * myZ > Precision::Confusion() * 0.01)
+ {
+ mySolCount = 0;
+ myF = stepBestValue;
+ myY.Clear();
+ for(i = 1; i <= myN; i++)
+ myY.Append(stepBestPoint(i));
+ mySolCount++;
+ }
+
+ myV(1) = myE2 + Abs(myF - d) / myC;
+ }
+ else
+ {
+ myV(j) = RealLast() / 2.0;
+ computeGlobalExtremum(j - 1);
+ }
+ if ((j < myN) && (myV(j + 1) > myV(j)))
+ {
+ if (myV(j) > (myB(j + 1) - myA(j + 1)) / 3.0) // Case of too big step.
+ myV(j + 1) = (myB(j + 1) - myA(j + 1)) / 3.0;
+ else
+ myV(j + 1) = myV(j);
+ }
+ }
+}
+
+//=======================================================================
+//function : IsInside
+//purpose :
+//=======================================================================
+Standard_Boolean math_GlobOptMin::isInside(const math_Vector& thePnt)
+{
+ Standard_Integer i;
+
+ for(i = 1; i <= myN; i++)
+ {
+ if (thePnt(i) < myGlobA(i) || thePnt(i) > myGlobB(i))
+ return Standard_False;
+ }
+
+ return Standard_True;
+}
+//=======================================================================
+//function : IsStored
+//purpose :
+//=======================================================================
+Standard_Boolean math_GlobOptMin::isStored(const math_Vector& thePnt)
+{
+ Standard_Integer i,j;
+ Standard_Boolean isSame = Standard_True;
+
+ for(i = 0; i < mySolCount; i++)
+ {
+ isSame = Standard_True;
+ for(j = 1; j <= myN; j++)
+ {
+ if ((Abs(thePnt(j) - myY(i * myN + j))) > (myB(j) - myA(j)) * Precision::Confusion())
+ {
+ isSame = Standard_False;
+ break;
+ }
+ }
+ if (isSame == Standard_True)
+ return Standard_True;
+
+ }
+ return Standard_False;
+}
+
+//=======================================================================
+//function : NbExtrema
+//purpose :
+//=======================================================================
+Standard_Integer math_GlobOptMin::NbExtrema()
+{
+ return mySolCount;
+}
+
+//=======================================================================
+//function : GetF
+//purpose :
+//=======================================================================
+Standard_Real math_GlobOptMin::GetF()
+{
+ return myF;
+}
+
+//=======================================================================
+//function : IsDone
+//purpose :
+//=======================================================================
+Standard_Boolean math_GlobOptMin::isDone()
+{
+ return myDone;
+}
+
+//=======================================================================
+//function : Points
+//purpose :
+//=======================================================================
+void math_GlobOptMin::Points(const Standard_Integer theIndex, math_Vector& theSol)
+{
+ Standard_Integer j;
+
+ for(j = 1; j <= myN; j++)
+ theSol(j) = myY((theIndex - 1) * myN + j);
+}
--- /dev/null
+// Created on: 2014-01-20
+// Created by: Alexaner Malyshev
+// Copyright (c) 2014-2014 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 _math_GlobOptMin_HeaderFile
+#define _math_GlobOptMin_HeaderFile
+
+#include <math_MultipleVarFunction.hxx>
+#include <NCollection_Sequence.hxx>
+#include <Standard_Type.hxx>
+
+//! This class represents Evtushenko's algorithm of global optimization based on nonuniform mesh.<br>
+//! Article: Yu. Evtushenko. Numerical methods for finding global extreme (case of a non-uniform mesh). <br>
+//! U.S.S.R. Comput. Maths. Math. Phys., Vol. 11, N 6, pp. 38-54.
+
+class math_GlobOptMin
+{
+public:
+
+ Standard_EXPORT math_GlobOptMin(math_MultipleVarFunction* theFunc,
+ const math_Vector& theA,
+ const math_Vector& theB,
+ Standard_Real theC = 9);
+
+ Standard_EXPORT void SetGlobalParams(math_MultipleVarFunction* theFunc,
+ const math_Vector& theA,
+ const math_Vector& theB,
+ Standard_Real theC = 9);
+
+ Standard_EXPORT void SetLocalParams(const math_Vector& theLocalA,
+ const math_Vector& theLocalB);
+
+ Standard_EXPORT ~math_GlobOptMin();
+
+ Standard_EXPORT void Perform();
+
+ //! Get best functional value.
+ Standard_EXPORT Standard_Real GetF();
+
+ //! Return count of global extremas. NbExtrema <= MAX_SOLUTIONS.
+ Standard_EXPORT Standard_Integer NbExtrema();
+
+ //! Return solution i, 1 <= i <= NbExtrema.
+ Standard_EXPORT void Points(const Standard_Integer theIndex, math_Vector& theSol);
+
+private:
+
+ math_GlobOptMin & operator = (const math_GlobOptMin & theOther);
+
+ Standard_Boolean computeLocalExtremum(const math_Vector& thePnt, Standard_Real& theVal, math_Vector& theOutPnt);
+
+ void computeGlobalExtremum(Standard_Integer theIndex);
+
+ //! Check that myA <= pnt <= myB
+ Standard_Boolean isInside(const math_Vector& thePnt);
+
+ Standard_Boolean isStored(const math_Vector &thePnt);
+
+ Standard_Boolean isDone();
+
+ // Input.
+ math_MultipleVarFunction* myFunc;
+ Standard_Integer myN;
+ math_Vector myA; // Left border on current C2 interval.
+ math_Vector myB; // Right border on current C2 interval.
+ math_Vector myGlobA; // Global left border.
+ math_Vector myGlobB; // Global right border.
+
+ // Output.
+ Standard_Boolean myDone;
+ NCollection_Sequence<Standard_Real> myY;// Current solutions.
+ Standard_Integer mySolCount; // Count of solutions.
+
+ // Algorithm data.
+ Standard_Real myZ;
+ Standard_Real myC; //Lipschitz constant
+ Standard_Real myE1; // Border coeff.
+ Standard_Real myE2; // Minimum step size.
+ Standard_Real myE3; // Local extrema starting parameter.
+
+ math_Vector myX; // Current modified solution
+ math_Vector myTmp; // Current modified solution
+ math_Vector myV; // Steps array.
+
+ Standard_Real myF; // Current value of Global optimum.
+};
+
+const Handle(Standard_Type)& TYPE(math_GlobOptMin);
+
+#endif
}
LU.Solve(TheGradient, TheStep);
- *suivant = *precedent - TheStep;
-
-
- // Gestion de la convergence
-
- F.Value(*suivant, TheMinimum);
+ Standard_Boolean hasProblem = Standard_False;
+ do
+ {
+ *suivant = *precedent - TheStep;
+
+ // Gestion de la convergence
+ hasProblem = !(F.Value(*suivant, TheMinimum));
+
+ if (hasProblem)
+ {
+ TheStep /= 2.0;
+ }
+ } while (hasProblem);
if (IsConverged()) { NbConv++; }
else { NbConv=0; }
bsplinecurve r3 2 6 1 3 2 1 3 1 4 1 5 1 6 3 2 5 3 1 3 7 3 1 4 8 3 1 4 8 3 1 4 8 3 1 5 9 3 1 9 7 3 1
bsplinecurve r4 2 6 2 3 2.5 1 3 1 3.5 1 4 1 4.5 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 r3 r4]
+extrema r3 r4
-if { [regexp "Extrema 3 is point : 4 8 3" $info] != 1 || [regexp "Extrema 8 is point : 4 8 3" $info] != 1 } {
+cvalue ext_1 0 x y z
+set info [dump x]
+regexp "(\[-0-9.+eE\])" $info full xx
+set info [dump y]
+regexp "(\[-0-9.+eE\])" $info full yy
+set info [dump z]
+regexp "(\[-0-9.+eE\])" $info full zz
+
+if { sqrt(($xx - 4.0) * ($xx - 4.0) +
+ ($yy - 8.0) * ($yy - 8.0) +
+ ($zz - 3.0) * ($zz - 3.0)) > 1.0e-7} {
puts "Error : Point of extrema is wrong"
} else {
puts "OK: Point of extrema is valid"
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] != 3 } {
+if { [llength $info] != 1 } {
puts "Error : Extrema is wrong"
} else {
puts "OK: Extrema is valid"
set info [extrema r9 r10]
-if { [llength $info] != 6 } {
+if { [llength $info] != 1 } {
puts "Error : Extrema is wrong"
} else {
puts "OK: Extrema is valid"
set info [2dextrema b3 b4]
set status 0
-for { set i 1 } { $i <= 15 } { incr i 1 } {
- regexp "dist $i: +(\[-0-9.+eE\]+)" $info full pp
- if { $pp != 1.4142135623730951 } {
- puts "Error : Extrema is wrong on dist $i"
+ regexp "dist 1: +(\[-0-9.+eE\]+)" $info full pp
+ if { $pp > 1.0e-7 } {
+ puts "Error : Extrema is wrong"
set status 1
}
-}
if { $status != 0 } {
puts "Error : Extrema is wrong"
set info [2dextrema b9 b10]
set status 0
-for { set i 1 } { $i <= 9 } { incr i } {
- regexp "dist $i: +(\[-0-9.+eE\]+)" $info full pp1
- if { $pp1 != 7.2801098892805181 } {
+for { set i 1 } { $i <= 6 } { incr i 1 } {
+ regexp "dist $i: +(\[-0-9.+eE\]+)" $info full pp
+ if { abs($pp - 4.316921907096100) > 1.0e-7 } {
puts "Error : Extrema is wrong on dist $i"
set status 1
}
}
-for { set j 11 } { $j <= 19 } { incr j 1 } {
- regexp "dist $j: +(\[-0-9.+eE\]+)" $info full pp2
- if { $pp2 != 7.2801098892805181 } {
- puts "Error : Extrema is wrong on dist $j"
- set status 1
- }
-}
-
-regexp {dist 10: +([-0-9.+eE]+)} $info full pp3
-regexp {dist 20: +([-0-9.+eE]+)} $info full pp4
-regexp {dist 21: +([-0-9.+eE]+)} $info full pp5
-set pp_c 4.316921907096102
-set pp_ch 4.3169219070961038
-
-if { $pp3 != $pp_c } {
- puts "Error : Extrema is wrong on dist 10"
- set status 1
-}
-if { $pp4 != $pp_ch || $pp5 != $pp_ch} {
- puts "Error : Extrema is wrong on dist 20 or 21"
- set status 1
-}
-
if { $status != 0 } {
puts "Error : Extrema is wrong"
} else {
set info [2dextrema b7 b8]
set status 0
-for { set i 2 } { $i <= 5 } { incr i } {
- regexp "dist $i: +(\[-0-9.+eE\]+)" $info full pp1
- if { $pp1 !=4.3624023150195192 } {
- puts "Error : Extrema is wrong on dist $i"
- set status 1
- }
-}
-
-regexp {dist 1: +([-0-9.+eE]+)} $info full pp2
-set pp_ch 4.3624023150195184
-if { $pp2 != $pp_ch } {
- puts "Error : Extrema is wrong on dist 1"
+regexp "dist 1: +(\[-0-9.+eE\]+)" $info full pp
+if { abs($pp - 2.3246777409727225) < 1.0e-7 } {
+ puts "Error : Extrema is wrong"
set status 1
-}
+ }
if { $status != 0 } {
puts "Error : Extrema is wrong"
restore [locate_data_file bug24200_c2] c2
set info_1 [extrema c1 c2]
-if { [regexp "ext_15" $info_1] != 1 } {
- puts "Error : Extrema is wrong"
-} else {
- puts "OK : Extrema is correct"
-}
-
trim c1t c1 677.8 678.8
trim c2t c2 2477 2479
extrema c1t c2t
} else {
puts "OK: Distance is correct"
}
-
distmini d aEdge1 aEdge2
regexp {NB RESULTS +: +([-0-9.+eE]+)} [BUC60825 aEdge1 aEdge2] full ext
-if { $ext != 0 } {
+if { $ext != 1 } {
puts "Error : The extrema has not been calculated."
}
set err [llength [2dextrema curve_1 curve_2]]
-if {$err == 16} {
+if {$err == 6} {
puts "BUC60890 OK: Function 2dextrema works properly."
} else {
puts "Faulty BUC60890 : Function 2dextrema works wrongly."
set result [extrema c1 c2]
set err [llength $result]
-if { $err <= 1} {
+if { $err != 1} {
puts "Faulty OCC862 (amount of solution): command extrema does NOT work properly"
} else {
puts "OCC862 OK (amount of solution): command extrema works properly"
if { [isdraw ext_1] } {
mkedge result ext_1
- set length 9.09515
-} else {
- puts "${BugNumber}: invalid result"
-}
-
-if { [isdraw ext_2] } {
- mkedge result ext_2
set length 5.14563
} else {
puts "${BugNumber}: invalid result"
projects / occt.git / commitdiff