Mercurial > octave

changeset 26332:d4211810202a

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
qp and sqp: Non-fixed tolerance for qp (bug #53506). * __qp__.cc (F__qp__): Add the tolerance as an input parameter. * qp.m: Add tolerance "TolFun" to the options struct. * sqp.m: The qp subproblem now accepts the tolerance passed to sqp as a parameter.
author Maor Shutman <maorus12@gmail.com>
date Fri, 04 May 2018 19:18:46 +0300
parents 0a78614a2484
children 4bad0d5b97b3
files libinterp/corefcn/__qp__.cc scripts/optimization/qp.m scripts/optimization/sqp.m
diffstat 3 files changed, 12 insertions(+), 12 deletions(-) [+]
line wrap: on
line diff
--- a/libinterp/corefcn/__qp__.cc	Sun Dec 30 20:53:52 2018 -0800
+++ b/libinterp/corefcn/__qp__.cc	Fri May 04 19:18:46 2018 +0300
@@ -89,15 +89,13 @@
 qp (const Matrix& H, const ColumnVector& q,
     const Matrix& Aeq, const ColumnVector& beq,
     const Matrix& Ain, const ColumnVector& bin,
-    int maxit,
+    int maxit, double rtol,
     ColumnVector& x, ColumnVector& lambda, int& iter)
 {
   int info = 0;
 
   iter = 0;
 
-  double rtol = sqrt (std::numeric_limits<double>::epsilon ());
-
   // Problem dimension.
   octave_idx_type n = x.numel ();
 
@@ -493,7 +491,7 @@
 Undocumented internal function.
 @end deftypefn */)
 {
-  if (args.length () != 8)
+  if (args.length () != 9)
     print_usage ();
 
   const ColumnVector x0  (args(0).vector_value ());
@@ -504,7 +502,8 @@
   const Matrix Ain       (args(5).matrix_value ());
   const ColumnVector bin (args(6).vector_value ());
   const int maxit        (args(7).int_value ());
-
+  const double rtol      (args(8).double_value());
+    
   int iter = 0;
 
   // Copy the initial guess into the working variable
@@ -513,7 +512,7 @@
   // Reordering the Lagrange multipliers
   ColumnVector lambda;
 
-  int info = qp (H, q, Aeq, beq, Ain, bin, maxit, x, lambda, iter);
+  int info = qp (H, q, Aeq, beq, Ain, bin, maxit, rtol, x, lambda, iter);
 
   return ovl (x, lambda, info, iter);
 }
--- a/scripts/optimization/qp.m	Sun Dec 30 20:53:52 2018 -0800
+++ b/scripts/optimization/qp.m	Fri May 04 19:18:46 2018 +0300
@@ -169,6 +169,7 @@
   endif
 
   maxit = optimget (options, "MaxIter", 200);
+  tol = optimget (options, "TolX", sqrt (eps));
 
   ## Validate the quadratic penalty.
   if (! issquare (H))
@@ -243,7 +244,7 @@
     endif
 
     if (! isempty (lb) && ! isempty (ub))
-      rtol = sqrt (eps);
+      rtol = tol;
       for i = 1:n
         if (abs (lb (i) - ub(i)) < rtol*(1 + max (abs (lb(i) + ub(i)))))
           ## These are actually an equality constraint
@@ -291,7 +292,7 @@
       endif
 
       if (! isempty (A_lb) && ! isempty (A_ub))
-        rtol = sqrt (eps);
+        rtol = tol;
         for i = 1:dimA_in
           if (abs (A_lb(i) - A_ub(i))
               < rtol*(1 + max (abs (A_lb(i) + A_ub(i)))))
@@ -331,7 +332,7 @@
       || isa (A, "single") || isa (b, "single"))
     rtol = sqrt (eps ("single"));
   else
-    rtol = sqrt (eps);
+    rtol = tol;
   endif
 
   eq_infeasible = (n_eq > 0 && norm (A*x0-b) > rtol*(1+abs (b)));
@@ -407,10 +408,10 @@
       x0 = xbar;
     endif
   endif
-
+  
   if (info == 0)
     ## The initial (or computed) guess is feasible.  Call the solver.
-    [x, lambda, info, iter] = __qp__ (x0, H, q, A, b, Ain, bin, maxit);
+    [x, lambda, info, iter] = __qp__ (x0, H, q, A, b, Ain, bin, maxit, rtol);
   else
     iter = 0;
     x = x0;
--- a/scripts/optimization/sqp.m	Sun Dec 30 20:53:52 2018 -0800
+++ b/scripts/optimization/sqp.m	Fri May 04 19:18:46 2018 +0300
@@ -417,7 +417,7 @@
 
     old_lambda = lambda;
     [p, obj_qp, INFO, lambda] = qp (x, B, c, F, g, [], [], d, C,
-                                    Inf (size (d)));
+                                    Inf (size (d)), struct ("TolX", tol));
 
     info = INFO.info;