Minor bugfix to trust region subproblem solver - upgrade to version 1.0.1

Author: lindonroberts

README.rst

@@ -1,12 +1,14 @@
 ===================================================================================
 Py-BOBYQA: Derivative-Free Solver for Bound-Constrained Minimization |PyPI Version|
 ===================================================================================
-Py-BOBYQA is a flexible package for solving bound-constrained general objective minimization, without requiring derivatives of the objective. It is a Python implementation of the BOBYQA algorithm by Powell.
+Py-BOBYQA is a flexible package for solving bound-constrained general objective minimization, without requiring derivatives of the objective. It is a Python implementation of the BOBYQA algorithm by Powell. Py-BOBYQA is particularly useful when evaluations of the objective function are expensive and/or noisy.
 
 More details about Py-BOBYQA can be found in our paper: C. Cartis, J. Fiala, B. Marteau and L. Roberts, Improving the Flexibility and Robustness of Model-Based Derivative-Free Optimization Solvers, technical report, University of Oxford, (2018).
 
 The original paper by Powell is: M. J. D. Powell, The BOBYQA algorithm for bound constrained optimization without derivatives, technical report DAMTP 2009/NA06, University of Cambridge (2009), and the original Fortran implementation is available `here <http://mat.uc.pt/~zhang/software.html>`_.
 
+If you are interested in solving least-squares minimization problems, you may wish to try `DFO-LS <https://github.com/numericalalgorithmsgroup/dfols>`_, which has the same features as Py-BOBYQA (plus some more), and exploits the least-squares problem structure, so performs better on such problems.
+
 Documentation
 -------------
 See manual.pdf or `here <http://people.maths.ox.ac.uk/robertsl/pybobyqa>`_.
@@ -25,7 +27,7 @@ Additionally, the following python packages should be installed (these will be i
 
 Installation using pip
 ----------------------
-For easy installation, use *pip* (http://www.pip-installer.org/) as root::
+For easy installation, use `pip <http://www.pip-installer.org/>`_ as root::
 
     $ [sudo] pip install Py-BOBYQA
 

docs/history.rst

@@ -0,0 +1,12 @@
+Version History
+===============
+This section lists the different versions of Py-BOBYQA and the updates between them.
+
+Version 1.0 (6 Feb 2018)
+------------------------
+* Initial release of Py-BOBYQA
+
+Version 1.0.1 (20 Feb 2018)
+---------------------------
+* Minor bug fix to trust region subproblem solver (the output :code:`crvmin` is calculated correctly) - this has minimal impact on the performance of Py-BOBYQA.
+

docs/index.rst

@@ -12,7 +12,7 @@ Py-BOBYQA: Derivative-Free Optimizer for Bound-Constrained Minimization
 
 **Author:** `Lindon Roberts <[email protected]>`_
 
-Py-BOBYQA is a flexible package for finding local solutions to nonlinear, nonconvex minimization problems (with optional bound constraints), without requiring any derivatives of the objective. Py-BOBYQA is a Python implementation of the `BOBYQA <http://mat.uc.pt/~zhang/software.html#powell_software>`_ solver by Powell (documentation `here <http://www.damtp.cam.ac.uk/user/na/NA_papers/NA2009_06.pdf>`_).
+Py-BOBYQA is a flexible package for finding local solutions to nonlinear, nonconvex minimization problems (with optional bound constraints), without requiring any derivatives of the objective. Py-BOBYQA is a Python implementation of the `BOBYQA <http://mat.uc.pt/~zhang/software.html#powell_software>`_ solver by Powell (documentation `here <http://www.damtp.cam.ac.uk/user/na/NA_papers/NA2009_06.pdf>`_). It is particularly useful when evaluations of the objective function are expensive and/or noisy.
 
 That is, Py-BOBYQA solves
 
@@ -23,6 +23,8 @@ That is, Py-BOBYQA solves
 
 Full details of the Py-BOBYQA algorithm are given in our paper: C. Cartis, J. Fiala, B. Marteau and L. Roberts, Improving the Flexibility and Robustness of Model-Based Derivative-Free Optimization Solvers, technical report, University of Oxford, (2018).
 
+If you are interested in solving least-squares minimization problems, you may wish to try `DFO-LS <https://github.com/numericalalgorithmsgroup/dfols>`_, which has the same features as Py-BOBYQA (plus some more), and exploits the least-squares problem structure, so performs better on such problems.
+
 Py-BOBYQA is released under the GNU General Public License. Please `contact NAG <http://www.nag.com/content/worldwide-contact-information>`_ for alternative licensing.
 
 .. toctree::
@@ -34,6 +36,7 @@ Py-BOBYQA is released under the GNU General Public License. Please `contact NAG
    userguide
    advanced
    diagnostic
+   history
 
 Acknowledgements
 ----------------

docs/userguide.rst

@@ -124,11 +124,11 @@ Note that Py-BOBYQA is a randomized algorithm: in its first phase, it builds an
   
       ****** Py-BOBYQA Results ******
       Solution xmin = [ 1.  1.]
-      Objective value f(xmin) = 2.790048765e-19
-      Needed 211 objective evaluations (at 211 points)
-      Approximate gradient = [ -1.57869863e-08   7.96417081e-09]
-      Approximate Hessian = [[ 803.2582815  -399.50964113]
-       [-399.50964113  199.12578708]]
+      Objective value f(xmin) = 2.964036794e-19
+      Needed 213 objective evaluations (at 213 points)
+      Approximate gradient = [ -2.57280154e-08   1.26855793e-08]
+      Approximate Hessian = [[ 802.90904563 -400.46022134]
+       [-400.46022134  200.23335154]]
       Exit flag = 0
       Success: rho has reached rhoend
       ******************************
@@ -155,10 +155,10 @@ Py-BOBYQA correctly finds the solution to the constrained problem:
       ****** Py-BOBYQA Results ******
       Solution xmin = [ 0.9   0.81]
       Objective value f(xmin) = 0.01
-      Needed 129 objective evaluations (at 129 points)
-      Approximate gradient = [ -2.00000002e-01  -4.01842008e-09]
-      Approximate Hessian = [[ 649.34361541 -364.85712365]
-       [-364.85712365  200.53013142]]
+      Needed 134 objective evaluations (at 134 points)
+      Approximate gradient = [ -1.99999226e-01  -4.31078784e-07]
+      Approximate Hessian = [[ 649.6790222  -360.18361979]
+       [-360.18361979  200.00205196]]
       Exit flag = 0
       Success: rho has reached rhoend
       ******************************
@@ -189,8 +189,8 @@ And we can now see each evaluation of :code:`objfun`:
       Function eval 2 at point 2 has f = 14.337296 at x = [-1.08  0.85]
       Function eval 3 at point 3 has f = 55.25 at x = [-1.2   0.73]
       ...
-      Function eval 128 at point 128 has f = 0.0100000000000225 at x = [ 0.9         0.81000002]
-      Function eval 129 at point 129 has f = 0.00999999999999997 at x = [ 0.9   0.81]
+      Function eval 133 at point 133 has f = 0.0100000000000165 at x = [ 0.9         0.81000001]
+      Function eval 134 at point 134 has f = 0.00999999999999997 at x = [ 0.9   0.81]
       Did a total of 1 run(s)
 
 If we wanted to save this output to a file, we could replace the above call to :code:`logging.basicConfig()` with
@@ -287,18 +287,18 @@ As noted above, Py-BOBYQA has an input parameter :code:`objfun_has_noise` to ind
   
       soln = pybobyqa.solve(rosenbrock_noisy, x0, objfun_has_noise=True)
 
-Using this setting, we get a more accurate solution, and better estimates of the gradient and Hessian:
+This time, we find the true solution, and better estimates of the gradient and Hessian:
 
   .. code-block:: none
   
       ****** Py-BOBYQA Results ******
       Solution xmin = [ 1.  1.]
-      Objective value f(xmin) = 3.559647071e-18
+      Objective value f(xmin) = 3.418770987e-18
       Needed 300 objective evaluations (at 300 points)
       Did a total of 4 runs
-      Approximate gradient = [ -3.70447710e-07   1.81205404e-07]
-      Approximate Hessian = [[ 804.44834579 -394.71053944]
-       [-394.71053944  194.52019795]]
+      Approximate gradient = [ -1.36175005e-08   2.12249758e-09]
+      Approximate Hessian = [[ 805.93202374 -394.16671315]
+       [-394.16671315  192.99451721]]
       Exit flag = 1
       Warning (max evals): Objective has been called MAXFUN times
       ******************************

examples/rosenbrock_noisy.py

@@ -23,8 +23,8 @@ def rosenbrock_noisy(x):
 print("")
 
 # Call Py-BOBYQA
-soln = pybobyqa.solve(rosenbrock_noisy, x0)
-#soln = pybobyqa.solve(rosenbrock_noisy, x0, objfun_has_noise=True)
+#soln = pybobyqa.solve(rosenbrock_noisy, x0)
+soln = pybobyqa.solve(rosenbrock_noisy, x0, objfun_has_noise=True)
 
 # Display output
 print(soln)

manual.pdf

@@ -121,7 +121,7 @@ def runTest(self):
         soln = pybobyqa.solve(objfun, x0)
         self.assertTrue(array_compare(soln.x, xmin, thresh=1e-2), "Wrong xmin")
         self.assertTrue(array_compare(soln.gradient, gradfun(soln.x), thresh=1e-2), "Wrong gradient")
-        self.assertTrue(array_compare(soln.hessian, hessfun, thresh=1e-1), "Wrong Hessian")
+        self.assertTrue(array_compare(soln.hessian, hessfun, thresh=0.5), "Wrong Hessian")
         self.assertTrue(abs(soln.f - fmin) < 1e-4, "Wrong fmin")
 
 

pybobyqa/tests/test_solver.py

@@ -92,7 +92,8 @@ def runTest(self):
         s_cauchy, red_cauchy, crvmin_cauchy = cauchy_pt(g, hess, Delta)
         self.assertTrue(est_min <= red_cauchy, 'Cauchy reduction not achieved')
         self.assertTrue(np.all(gnew == g + hess.vec_mul(d)), 'Wrong gnew')
-        self.assertAlmostEqual(crvmin, -1.0, 'Wrong crvmin')
+        print(crvmin)
+        self.assertAlmostEqual(crvmin, 1.2, 'Wrong crvmin')
 
 
 class TestUncBdry(unittest.TestCase):
@@ -210,7 +211,8 @@ def runTest(self):
         # print(d)
         self.assertTrue(est_min <= red_cauchy, 'Cauchy reduction not achieved')
         self.assertTrue(np.all(gnew == g + hess.vec_mul(d)), 'Wrong gnew')
-        self.assertAlmostEqual(crvmin, 1.5, 'Wrong crvmin')
+        print(crvmin)
+        self.assertAlmostEqual(crvmin, -1.0, 'Wrong crvmin')
 
 
 class TestConBdry(unittest.TestCase):
@@ -233,7 +235,8 @@ def runTest(self):
         s_cauchy, red_cauchy, crvmin_cauchy = cauchy_pt_box(g, hess, Delta, sl - xopt, su - xopt)
         self.assertTrue(est_min <= red_cauchy, 'Cauchy reduction not achieved')
         self.assertTrue(np.max(np.abs(gnew - g - hess.vec_mul(d))) < 1e-10, 'Wrong gnew')
-        self.assertAlmostEqual(crvmin, 1.0, 'Wrong crvmin')
+        print(crvmin)
+        self.assertAlmostEqual(crvmin, -1.0, 'Wrong crvmin')
         # self.assertAlmostEqual(crvmin, crvmin_cauchy, 'Wrong crvmin')
 
 

pybobyqa/tests/test_trust_region.py

@@ -140,7 +140,7 @@ def trsbox(xopt, g, hess, sl, su, delta):
 
         # Reduce STPLEN if necessary in order to preserve the simple bounds,
         # letting IACT be the index of the new constrained variable.
-        iact = -1
+        iact = None
         for i in range(n):
             if s[i] != 0.0:
                 temp = (su[i] - xopt[i] - d[i] if s[i] > 0.0 else sl[i] - xopt[i] - d[i]) / s[i]
@@ -153,8 +153,8 @@ def trsbox(xopt, g, hess, sl, su, delta):
         if stplen > 0.0:
             iterc += 1
             temp = shs / stepsq
-            if iact == 0 and temp > 0.0:
-                crvmin = (min(crvmin, temp) if crvmin != -1.0 else temp)
+            if iact is None and temp > 0.0:
+                crvmin = min(crvmin, temp) if crvmin != -1.0 else temp
             ggsav = gredsq
             gnew += stplen * hs
             d += stplen * s
@@ -163,7 +163,7 @@ def trsbox(xopt, g, hess, sl, su, delta):
             qred += sdec
 
         # Restart the conjugate gradient method if it has hit a new bound.
-        if iact > -1:
+        if iact is not None:
             nact += 1
             xbdi[iact] = (1 if s[iact] >= 0.0 else -1)
             delsq = delsq - d[iact] ** 2
@@ -240,7 +240,7 @@ def alt_trust_step(n, xopt, hess, sl, su, d, xbdi, nact, gnew, qred):
             # bound, there is a branch back to label 100 after fixing that variable.
             free_variable_reached_bound = False
             angbd = 1.0
-            iact = -1
+            iact = None
             for i in range(n):
                 if xbdi[i] == 0:
                     tempa = xopt[i] + d[i] - sl[i]
@@ -334,7 +334,7 @@ def alt_trust_step(n, xopt, hess, sl, su, d, xbdi, nact, gnew, qred):
             hred = cth * hred + sth * hs
 
             qred += sdec
-            if iact > -1 and isav == iu - 1:
+            if iact is not None and isav == iu - 1:
                 nact += 1
                 xbdi[iact] = xsav
                 restart_alt_loop = True

pybobyqa/trust_region.py

@@ -22,4 +22,4 @@
 
 """
 
-__version__ = '1.0'
+__version__ = '1.0.1'