#include <Standard_Real.hxx>
#include <Precision.hxx>
+//=======================================================================
+//function : DistanceToBorder
+//purpose :
+//=======================================================================
+static Standard_Real DistanceToBorder(const math_Vector & theX,
+ const math_Vector & theMin,
+ const math_Vector & theMax)
+{
+ Standard_Real aDist = RealLast();
+
+ for (Standard_Integer anIdx = theMin.Lower(); anIdx <= theMin.Upper(); ++anIdx)
+ {
+ const Standard_Real aDist1 = Abs (theX(anIdx) - theMin(anIdx));
+ const Standard_Real aDist2 = Abs (theX(anIdx) - theMax(anIdx));
+
+ aDist = Min (aDist, Min (aDist1, aDist2));
+ }
+
+ return aDist;
+}
+
//=======================================================================
//function : math_GlobOptMin
myTmp(1, myN),
myV(1, myN),
myMaxV(1, myN),
- myExpandCoeff(1, myN),
myCellSize(0, myN - 1),
myFilter(theFunc->NbVariables()),
myCont(2)
myMaxV(i) = (myB(i) - myA(i)) / 3.0;
}
- myExpandCoeff(1) = 1.0;
- for(i = 2; i <= myN; i++)
- {
- myExpandCoeff(i) = (myB(i) - myA(i)) / (myB(i - 1) - myA(i - 1));
- }
-
myTol = theDiscretizationTol;
mySameTol = theSameTol;
myMaxV(i) = (myB(i) - myA(i)) / 3.0;
}
- myExpandCoeff(1) = 1.0;
- for(i = 2; i <= myN; i++)
- {
- myExpandCoeff(i) = (myB(i) - myA(i)) / (myB(i - 1) - myA(i - 1));
- }
-
myTol = theDiscretizationTol;
mySameTol = theSameTol;
myMaxV(i) = (myB(i) - myA(i)) / 3.0;
}
- myExpandCoeff(1) = 1.0;
- for(i = 2; i <= myN; i++)
- {
- myExpandCoeff(i) = (myB(i) - myA(i)) / (myB(i - 1) - myA(i - 1));
- }
-
myDone = Standard_False;
}
theSameTol = mySameTol;
}
-//=======================================================================
-//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(const Standard_Boolean isFindSingleSolution)
{
- Standard_Integer i;
+ myDone = Standard_False;
// Compute parameters range
Standard_Real minLength = RealLast();
Standard_Real maxLength = RealFirst();
- for(i = 1; i <= myN; i++)
+ for(Standard_Integer i = 1; i <= myN; i++)
{
Standard_Real currentLength = myB(i) - myA(i);
if (currentLength < minLength)
minLength = currentLength;
if (currentLength > maxLength)
maxLength = currentLength;
+
+ myV(i) = 0.0;
}
if (minLength < Precision::PConfusion())
void math_GlobOptMin::computeGlobalExtremum(Standard_Integer j)
{
Standard_Integer i;
- Standard_Real d; // Functional in moved point.
+ Standard_Real d = RealLast(), aPrevVal; // Functional in original and moved points.
Standard_Real val = RealLast(); // Local extrema computed in moved point.
Standard_Real aStepBestValue = RealLast();
math_Vector aStepBestPoint(1, myN);
- Standard_Boolean isInside = Standard_False;
- Standard_Real r;
- Standard_Boolean isReached = Standard_False;
+ Standard_Boolean isInside = Standard_False,
+ isReached = Standard_False;
+ Standard_Real r1, r2, r;
- for(myX(j) = myA(j) + myE1;
- (myX(j) < myB(j) + myE1) && (!isReached);
- myX(j) += myV(j))
+ for(myX(j) = myA(j) + myE1; !isReached; myX(j) += myV(j))
{
if (myX(j) > myB(j))
{
if (j == 1)
{
isInside = Standard_False;
+ aPrevVal = d;
myFunc->Value(myX, d);
- r = (d + myZ * myC * myLastStep - myF) * myZ;
+ r1 = (d + myZ * myC * myLastStep - myF) * myZ; // Evtushenko estimation.
+ r2 = ((d + aPrevVal - myC * myLastStep) * 0.5 - myF) * myZ; // Shubert / Piyavsky estimation.
+ r = Min(r1, r2);
if(r > myE3)
{
- isInside = computeLocalExtremum(myX, val, myTmp);
+ Standard_Real aSaveParam = myX(1);
+
+ // Piyavsky midpoint estimation.
+ Standard_Real aParam = (2 * myX(1) - myV(1) ) * 0.5 + (aPrevVal - d) * 0.5 / myC;
+ if (Precision::IsInfinite(aPrevVal))
+ aParam = myX(1) - myV(1) * 0.5; // Protection from upper dimension step.
+
+ myX(1) = aParam;
+ Standard_Real aVal = 0;
+ myFunc->Value(myX, aVal);
+ myX(1) = aSaveParam;
+
+ if ( (aVal < d && aVal < aPrevVal) ||
+ DistanceToBorder(myX, myA, myB) < myE1 ) // Condition optimization case near the border.
+ {
+ isInside = computeLocalExtremum(myX, val, myTmp);
+ }
}
aStepBestValue = (isInside && (val < d))? val : d;
aStepBestPoint = (isInside && (val < d))? myTmp : myX;
for(i = 1; i < j; i++)
myV(i) = 0.0;
}
- // Compute step in (j + 1) dimension according to scale.
if (j < myN)
{
- Standard_Real aUpperDimStep = myV(j) * myExpandCoeff(j + 1);
+ Standard_Real aUpperDimStep = Max(myV(j), myE2);
if (myV(j + 1) > aUpperDimStep)
{
if (aUpperDimStep > myMaxV(j + 1)) // Case of too big step.
return Standard_False;
}
-//=======================================================================
-//function : NbExtrema
-//purpose :
-//=======================================================================
-Standard_Integer math_GlobOptMin::NbExtrema()
-{
- return mySolCount;
-}
-
-//=======================================================================
-//function : GetF
-//purpose :
-//=======================================================================
-Standard_Real math_GlobOptMin::GetF()
-{
- return myF;
-}
-
-//=======================================================================
-//function : SetFunctionalMinimalValue
-//purpose :
-//=======================================================================
-void math_GlobOptMin::SetFunctionalMinimalValue(const Standard_Real theMinimalValue)
-{
- myFunctionalMinimalValue = theMinimalValue;
-}
-
-//=======================================================================
-//function : GetFunctionalMinimalValue
-//purpose :
-//=======================================================================
-Standard_Real math_GlobOptMin::GetFunctionalMinimalValue()
-{
- return myFunctionalMinimalValue;
-}
-
-//=======================================================================
-//function : IsDone
-//purpose :
-//=======================================================================
-Standard_Boolean math_GlobOptMin::isDone()
-{
- return myDone;
-}
-
//=======================================================================
//function : Points
//purpose :
for(i = 1; i <= myN; i++)
myY.Append(aBestPnt(i));
mySolCount = 1;
-
- myDone = Standard_False;
}
-
-//=======================================================================
-//function : SetLipConstState
-//purpose :
-//=======================================================================
-void math_GlobOptMin::SetLipConstState(const Standard_Boolean theFlag)
-{
- myIsConstLocked = theFlag;
-}
\ No newline at end of file
#include <NCollection_Sequence.hxx>
#include <Standard_Type.hxx>
-class NCollection_CellFilter_Inspector
-{
-public:
-
- //! Points and target type
- typedef math_Vector Point;
- typedef math_Vector Target;
-
- NCollection_CellFilter_Inspector(const Standard_Integer theDim,
- const Standard_Real theTol)
- : myCurrent(1, theDim)
- {
- myTol = theTol * theTol;
- myIsFind = Standard_False;
- Dimension = theDim;
- }
-
- //! Access to co-ordinate
- static Standard_Real Coord (int i, const Point &thePnt)
- {
- return thePnt(i + 1);
- }
-
- //! Auxiliary method to shift point by each coordinate on given value;
- //! useful for preparing a points range for Inspect with tolerance
- void Shift (const Point& thePnt,
- const NCollection_Array1<Standard_Real> &theTol,
- Point& theLowPnt,
- Point& theUppPnt) const
- {
- for(Standard_Integer anIdx = 1; anIdx <= Dimension; anIdx++)
- {
- theLowPnt(anIdx) = thePnt(anIdx) - theTol(anIdx - 1);
- theUppPnt(anIdx) = thePnt(anIdx) + theTol(anIdx - 1);
- }
- }
-
- void ClearFind()
- {
- myIsFind = Standard_False;
- }
-
- Standard_Boolean isFind()
- {
- return myIsFind;
- }
-
- //! Set current point to search for coincidence
- void SetCurrent (const math_Vector& theCurPnt)
- {
- myCurrent = theCurPnt;
- }
-
- //! Implementation of inspection method
- NCollection_CellFilter_Action Inspect (const Target& theObject)
- {
- Standard_Real aSqDist = (myCurrent - theObject).Norm2();
-
- if(aSqDist < myTol)
- {
- myIsFind = Standard_True;
- }
-
- return CellFilter_Keep;
- }
-
-private:
- Standard_Real myTol;
- math_Vector myCurrent;
- Standard_Boolean myIsFind;
- Standard_Integer Dimension;
-};
-
-//! 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>
+//! This class represents Evtushenko's algorithm of global optimization based on non-uniform mesh.
+//! Article: Yu. Evtushenko. Numerical methods for finding global extreme (case of a non-uniform mesh).
//! U.S.S.R. Comput. Maths. Math. Phys., Vol. 11, N 6, pp. 38-54.
-
+//!
+//! This method performs search on non-uniform mesh. The search space is a box in R^n space.
+//! The default behavior is to find all minimums in that box. Computation of maximums is not supported.
+//!
+//! The search box can be split into smaller boxes by discontinuity criteria.
+//! This functionality is covered by SetGlobalParams and SetLocalParams API.
+//!
+//! It is possible to set continuity of the local boxes.
+//! Such option can forcibly change local extrema search.
+//! In other words if theFunc can be casted to the function with Hessian but, continuity is set to 1
+//! Gradient based local optimization method will be used, not Hessian based method.
+//! This functionality is covered by SetContinuity and GetContinuity API.
+//!
+//! It is possible to freeze Lipschitz const to avoid internal modifications on it.
+//! This functionality is covered by SetLipConstState and GetLipConstState API.
+//!
+//! It is possible to perform single solution search.
+//! This functionality is covered by first parameter in Perform method.
+//!
+//! It is possible to set / get minimal value of the functional.
+//! It works well together with single solution search.
+//! This functionality is covered by SetFunctionalMinimalValue and GetFunctionalMinimalValue API.
class math_GlobOptMin
{
public:
+ //! Constructor. Perform method is not called from it.
+ //! @param theFunc - objective functional.
+ //! @param theLowerBorder - lower corner of the search box.
+ //! @param theUpperBorder - upper corner of the search box.
+ //! @param theC - Lipschitz constant.
+ //! @param theDiscretizationTol - parameter space discretization tolerance.
+ //! @param theSameTol - functional value space indifference tolerance.
Standard_EXPORT math_GlobOptMin(math_MultipleVarFunction* theFunc,
const math_Vector& theLowerBorder,
const math_Vector& theUpperBorder,
const Standard_Real theDiscretizationTol = 1.0e-2,
const Standard_Real theSameTol = 1.0e-7);
+ //! @param theFunc - objective functional.
+ //! @param theLowerBorder - lower corner of the search box.
+ //! @param theUpperBorder - upper corner of the search box.
+ //! @param theC - Lipschitz constant.
+ //! @param theDiscretizationTol - parameter space discretization tolerance.
+ //! @param theSameTol - functional value space indifference tolerance.
Standard_EXPORT void SetGlobalParams(math_MultipleVarFunction* theFunc,
const math_Vector& theLowerBorder,
const math_Vector& theUpperBorder,
const Standard_Real theDiscretizationTol = 1.0e-2,
const Standard_Real theSameTol = 1.0e-7);
+ //! Method to reduce bounding box. Perform will use this box.
+ //! @param theLocalA - lower corner of the local box.
+ //! @param theLocalB - upper corner of the local box.
Standard_EXPORT void SetLocalParams(const math_Vector& theLocalA,
const math_Vector& theLocalB);
+ //! Method to set tolerances.
+ //! @param theDiscretizationTol - parameter space discretization tolerance.
+ //! @param theSameTol - functional value space indifference tolerance.
Standard_EXPORT void SetTol(const Standard_Real theDiscretizationTol,
const Standard_Real theSameTol);
+ //! Method to get tolerances.
+ //! @param theDiscretizationTol - parameter space discretization tolerance.
+ //! @param theSameTol - functional value space indifference tolerance.
Standard_EXPORT void GetTol(Standard_Real& theDiscretizationTol,
Standard_Real& theSameTol);
- Standard_EXPORT ~math_GlobOptMin();
-
//! @param isFindSingleSolution - defines whether to find single solution or all solutions.
Standard_EXPORT void Perform(const Standard_Boolean isFindSingleSolution = Standard_False);
+ //! Return solution theIndex, 1 <= theIndex <= NbExtrema.
+ Standard_EXPORT void Points(const Standard_Integer theIndex, math_Vector& theSol);
+
+ //! Set / Get continuity of local borders splits (0 ~ C0, 1 ~ C1, 2 ~ C2).
+ inline void SetContinuity(const Standard_Integer theCont) { myCont = theCont; }
+ inline Standard_Integer GetContinuity() const {return myCont; }
+
+ //! Set / Get functional minimal value.
+ inline void SetFunctionalMinimalValue(const Standard_Real theMinimalValue)
+ { myFunctionalMinimalValue = theMinimalValue; }
+ inline Standard_Real GetFunctionalMinimalValue() const {return myFunctionalMinimalValue;}
+
+ //! Set / Get Lipchitz constant modification state.
+ //! True means that the constant is locked and unlocked otherwise.
+ inline void SetLipConstState(const Standard_Boolean theFlag) {myIsConstLocked = theFlag; }
+ inline Standard_Boolean GetLipConstState() const { return myIsConstLocked; }
+
+ //! Return computation state of the algorithm.
+ inline Standard_Boolean isDone() const {return myDone; }
+
//! Get best functional value.
- Standard_EXPORT Standard_Real GetF();
+ inline Standard_Real GetF() const {return myF;}
//! Return count of global extremas.
- Standard_EXPORT Standard_Integer NbExtrema();
+ inline Standard_Integer NbExtrema() const {return mySolCount;}
- //! Return solution theIndex, 1 <= theIndex <= NbExtrema.
- Standard_EXPORT void Points(const Standard_Integer theIndex, math_Vector& theSol);
+private:
- Standard_EXPORT Standard_Boolean isDone();
+ //! Class for duplicate fast search. For internal usage only.
+ class NCollection_CellFilter_Inspector
+ {
+ public:
- //! Set functional minimal value.
- Standard_EXPORT void SetFunctionalMinimalValue(const Standard_Real theMinimalValue);
+ //! Points and target type
+ typedef math_Vector Point;
+ typedef math_Vector Target;
- //! Lock/Unlock Lipchitz constant for internal modifications.
- Standard_EXPORT void SetLipConstState(const Standard_Boolean theFlag);
+ NCollection_CellFilter_Inspector(const Standard_Integer theDim,
+ const Standard_Real theTol)
+ : myCurrent(1, theDim)
+ {
+ myTol = theTol * theTol;
+ myIsFind = Standard_False;
+ Dimension = theDim;
+ }
- //! Get functional minimal value.
- Standard_EXPORT Standard_Real GetFunctionalMinimalValue();
+ //! Access to coordinate.
+ static Standard_Real Coord (int i, const Point &thePnt)
+ {
+ return thePnt(i + 1);
+ }
- //! Get continuity of local borders splits.
- inline Standard_Integer GetContinuity() const {return myCont; }
+ //! Auxiliary method to shift point by each coordinate on given value;
+ //! useful for preparing a points range for Inspect with tolerance.
+ void Shift (const Point& thePnt,
+ const NCollection_Array1<Standard_Real> &theTol,
+ Point& theLowPnt,
+ Point& theUppPnt) const
+ {
+ for(Standard_Integer anIdx = 1; anIdx <= Dimension; anIdx++)
+ {
+ theLowPnt(anIdx) = thePnt(anIdx) - theTol(anIdx - 1);
+ theUppPnt(anIdx) = thePnt(anIdx) + theTol(anIdx - 1);
+ }
+ }
- //! Set continuity of local borders splits.
- inline void SetContinuity(const Standard_Integer theCont) { myCont = theCont; }
+ void ClearFind()
+ {
+ myIsFind = Standard_False;
+ }
+
+ Standard_Boolean isFind()
+ {
+ return myIsFind;
+ }
+
+ //! Set current point to search for coincidence
+ void SetCurrent (const math_Vector& theCurPnt)
+ {
+ myCurrent = theCurPnt;
+ }
+
+ //! Implementation of inspection method
+ NCollection_CellFilter_Action Inspect (const Target& theObject)
+ {
+ Standard_Real aSqDist = (myCurrent - theObject).Norm2();
+
+ if(aSqDist < myTol)
+ {
+ myIsFind = Standard_True;
+ }
+
+ return CellFilter_Keep;
+ }
+
+ private:
+ Standard_Real myTol;
+ math_Vector myCurrent;
+ Standard_Boolean myIsFind;
+ Standard_Integer Dimension;
+ };
-private:
// Compute cell size.
void initCellSize();
void computeGlobalExtremum(Standard_Integer theIndex);
//! Check possibility to stop computations.
- //! Find single solution + in neighbourhood of best possible solution.
+ //! Find single solution + in neighborhood of best possible solution.
Standard_Boolean CheckFunctionalStopCriteria();
//! Computes starting value / approximation:
//! myF - initial best value.
//! myY - initial best point.
//! myC - approximation of Lipschitz constant.
- //! to imporve convergence speed.
+ //! to improve convergence speed.
void computeInitialValues();
//! Check that myA <= thePnt <= myB
// Algorithm data.
Standard_Real myZ;
- Standard_Real myE1; // Border coeff.
+ Standard_Real myE1; // Border coefficient.
Standard_Real myE2; // Minimum step size.
Standard_Real myE3; // Local extrema starting parameter.
math_Vector myV; // Steps array.
math_Vector myMaxV; // Max Steps array.
Standard_Real myLastStep; // Last step.
- math_Vector myExpandCoeff; // Define expand coefficient between neighboring indiced dimensions.
NCollection_Array1<Standard_Real> myCellSize;
Standard_Integer myMinCellFilterSol;