0026022: Extrema_ExtCC gives not precise solution

Conditional optimization added to Newton optimization algorithm.
Test case for issue CR26022

diff --git a/src/math/math_GlobOptMin.cxx b/src/math/math_GlobOptMin.cxx
index 2213926e34107aa5a24c0c6f146b86c9e3eb3269..059992bb2125d8ea368494b54709abe7db68c777 100644 (file)
--- a/src/math/math_GlobOptMin.cxx
+++ b/src/math/math_GlobOptMin.cxx
@@ -258,8 +258,8 @@ Standard_Boolean math_GlobOptMin::computeLocalExtremum(const math_Vector& thePnt
   {
     math_MultipleVarFunctionWithHessian* myTmp = 
       dynamic_cast<math_MultipleVarFunctionWithHessian*> (myFunc);
-    
     math_NewtonMinimum newtonMinimum(*myTmp);
+    newtonMinimum.SetBoundary(myGlobA, myGlobB);
     newtonMinimum.Perform(*myTmp, thePnt);
 
     if (newtonMinimum.IsDone())

diff --git a/src/math/math_NewtonMinimum.cdl b/src/math/math_NewtonMinimum.cdl
index d4a0b886aefea886f5e5f8c96820ad92fe2aadce..3649d79f3b0964e0d8f3b8466ac8d4324869ed24 100644 (file)
--- a/src/math/math_NewtonMinimum.cdl
+++ b/src/math/math_NewtonMinimum.cdl
@@ -97,8 +97,13 @@ is
     raises DimensionError,
           NotDone
     is static;
-    
-    
+
+    SetBoundary(me: in out; 
+                theLeftBorder : in Vector;
+                theRightBorder: in Vector)
+    ---Purpose: Set boundaries.
+    is static;
+
     Minimum(me)
        ---Purpose: returns the value of the minimum.
         -- Exception NotDone is raised if the minimum was not found.
@@ -166,6 +171,9 @@ CTol:            Real    is protected;
 nbiter:          Integer is protected;
 NoConvexTreatement: Boolean is protected;
 Convex           : Boolean is protected;
+myIsBoundsDefined  : Boolean is protected;
+myLeft             : Vector is protected;
+myRight            : Vector is protected;
 Itermax:         Integer;
 
 

diff --git a/src/math/math_NewtonMinimum.cxx b/src/math/math_NewtonMinimum.cxx
index 7d3ebc4ac4ee05e231776c24768c77e45af087c3..32f9714e72d379dbc3d518742ac3006db40771af 100644 (file)
--- a/src/math/math_NewtonMinimum.cxx
+++ b/src/math/math_NewtonMinimum.cxx
@@ -21,6 +21,7 @@
 //#endif
 
 #include <math_NewtonMinimum.ixx>
+#include <Precision.hxx>
 
 #include <math_Gauss.hxx>
 #include <math_Jacobi.hxx>
@@ -49,6 +50,9 @@ math_NewtonMinimum::math_NewtonMinimum(
   nbiter            (0),
   NoConvexTreatement(theWithSingularity),
   Convex            (Standard_True),
+  myIsBoundsDefined(Standard_False),
+  myLeft(1, theFunction.NbVariables(), 0.0),
+  myRight(1, theFunction.NbVariables(), 0.0),
   Done              (Standard_False),
   Itermax           (theNbIterations)
 {
@@ -62,6 +66,18 @@ math_NewtonMinimum::~math_NewtonMinimum()
 {
 }
 
+//=======================================================================
+//function : SetBoundary
+//purpose  : Set boundaries for conditional optimization
+//=======================================================================
+void math_NewtonMinimum::SetBoundary(const math_Vector& theLeftBorder,
+                                     const math_Vector& theRightBorder)
+{
+  myLeft = theLeftBorder;
+  myRight = theRightBorder;
+  myIsBoundsDefined = Standard_True;
+}
+
 //=======================================================================
 //function : Perform
 //purpose  : 
@@ -134,8 +150,52 @@ void math_NewtonMinimum::Perform(math_MultipleVarFunctionWithHessian& F,
       TheStatus = math_DirectionSearchError;
       return;
     }
-   
     LU.Solve(TheGradient, TheStep);
+
+    if (myIsBoundsDefined)
+    {
+      // Project point on bounds or nullify TheStep coords if point lies on boudary.
+
+      *suivant = *precedent - TheStep;
+      Standard_Real aMult = RealLast();
+      for(Standard_Integer anIdx = 1; anIdx <= myLeft.Upper(); anIdx++)
+      {
+        if (suivant->Value(anIdx) < myLeft(anIdx))
+        {
+          Standard_Real aValue = Abs(precedent->Value(anIdx) - myLeft(anIdx)) / Abs(TheStep(anIdx));
+          aMult = Min (aValue, aMult);
+        }
+
+        if (suivant->Value(anIdx) > myRight(anIdx))
+        {
+          Standard_Real aValue = Abs(precedent->Value(anIdx) - myRight(anIdx)) / Abs(TheStep(anIdx));
+          aMult = Min (aValue, aMult);
+        }
+      }
+
+      if (aMult != RealLast())
+      {
+        if (aMult > Precision::PConfusion())
+        {
+          // Project point into param space.
+          TheStep *= aMult;
+        }
+        else
+        {
+          // Old point on border and new point out of border:
+          // Nullify corresponding TheStep indexes.
+          for(Standard_Integer anIdx = 1; anIdx <= myLeft.Upper(); anIdx++)
+          {
+            if (Abs(precedent->Value(anIdx) - myRight(anIdx)) < Precision::PConfusion() ||
+                Abs(precedent->Value(anIdx) - myLeft(anIdx) ) < Precision::PConfusion())
+            {
+              TheStep(anIdx) = 0.0;
+            }
+          }
+        }
+      }
+    }
+
     Standard_Boolean hasProblem = Standard_False;
     do
     {

diff --git a/tests/bugs/fclasses/bug25635_1 b/tests/bugs/fclasses/bug25635_1
index db42b001ba7849371bd52e59919bbdfc9f2cc535..2ae860f420caafbb4c28df5f5d7e3c6c5d08b753 100755 (executable)
--- a/tests/bugs/fclasses/bug25635_1
+++ b/tests/bugs/fclasses/bug25635_1
@@ -36,7 +36,7 @@ checkreal "Parameter2" ${Parameter2} ${expected_Parameter2} ${tol_abs} ${tol_rel
 set expected_OriginX 2.
 checkreal "OriginX" ${OriginX} ${expected_OriginX} ${tol_abs} ${tol_rel}
 
-set expected_OriginY -6.e-05 
+set expected_OriginY 7.e-05 
 checkreal "OriginY" ${OriginY} ${expected_OriginY} ${tol_abs} ${tol_rel}
 
 set expected_AxisX 1.

diff --git a/tests/bugs/fclasses/bug26022 b/tests/bugs/fclasses/bug26022
new file mode 100644 (file)
index 0000000..c6d15e9
--- /dev/null
+++ b/tests/bugs/fclasses/bug26022
@@ -0,0 +1,40 @@
+puts "========"
+puts "OCC26022"
+puts "========"
+puts ""
+##############################################
+# Extrema_ExtCC gives not precise solution
+##############################################
+
+restore [locate_data_file bug26022_splitnoproblem671_notria.brep] a
+
+explode a
+explode a_2 e
+subshape a_2 e 2
+mkcurve c1 a_1
+mkcurve c2 a_2_2
+
+extrema c1 c2 1
+
+cvalue c1 prm_1_1 x1 y1 z1
+bounds c1 u1 u2
+cvalue c1 u2 x2 y2 z2
+
+regexp {is ([\d.\-e]+)} [length ext_1] str dist
+set dist_to_end [dval sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)+(z1-z2)*(z1-z2))]
+set parametric_dist_to_end [dval u2-prm_1_1]
+
+puts "dist of solution $dist"
+puts "parametric dist to end of curve $parametric_dist_to_end"
+puts "dist to end of curve $dist_to_end"
+
+set tol_abs 4.0e-12
+set tol_rel 1.0e-2
+
+set expected_dist 0.0
+set expected_parametric_dist_to_end 0.0
+set expected_dist_to_end 0.0
+
+checkreal "dist of solution" ${dist} ${expected_dist} ${tol_abs} ${tol_rel}
+checkreal "parametric dist to end of curve" ${parametric_dist_to_end} ${expected_parametric_dist_to_end} ${tol_abs} ${tol_rel}
+checkreal "dist to end of curve" ${dist_to_end} ${expected_dist_to_end} ${tol_abs} ${tol_rel}

diff --git a/tests/bugs/modalg_5/bug23706_14 b/tests/bugs/modalg_5/bug23706_14
index 72c0e5d4e61d437e489c3775b96a961b55c4a178..516eaee789da328c126b6f1321b4a7dbe6518d20 100755 (executable)
--- a/tests/bugs/modalg_5/bug23706_14
+++ b/tests/bugs/modalg_5/bug23706_14
@@ -13,7 +13,7 @@ set info [2dextrema b9 b10]
 set status 0
 for { set i 1 } { $i <= 1 } { incr i 1 } {
     regexp "dist $i: +(\[-0-9.+eE\]+)" $info full pp
-    if { abs($pp - 3.6710601078037173) > 1.0e-7 } {
+    if { abs($pp - 3.6712618987696386) > 1.0e-7 } {
        puts "Error : Extrema is wrong on dist $i"
        set status 1
     }