updates following recent writeup, using node with best LLB as best node

Author: thartland

src/hiopbbpy/opt/bnbalgorithm.py

@@ -724,6 +724,7 @@ def initialize(self, l0 = None, u0 = None, queue = None):
     # Least upper bound (LUB)
     self.best_node = root
     self.LUB = self.best_node.aq_U
+    self.LLB = self.best_node.aq_L
     if queue is not None:
       print("root node acqf_L, acqf_U = ", root.aq_L, root.aq_U)
       for _, _, node in queue:
@@ -827,19 +828,23 @@ def bnboptimize(self, l_init, u_init):
             print("upper-bound = ", child.aq_U)
           assert child.aq_U >= child.aq_L, "ERROR: child upper bound < child lower bound for child"
           if child.aq_U <= self.LUB:
+            self.LUB = child.aq_U
+          if child.aq_L <= self.LLB:
             self.best_node = child
-            self.LUB = self.best_node.aq_U
+            self.LLB = self.best_node.aq_L
             updated_best_node = True
-        gap_history.append(self.best_node.aq_U - self.best_node.aq_L)
         if not updated_best_node:
           print("best node not updated")
-          print("min |child.aq_U - LUB| = ", min([abs(child.aq_U - self.LUB) for child in children]))
+          print("min |child.aq_L - LLB| = ", min([abs(child.aq_L - self.LLB) for child in children]))
           if self.pure_BBS and self.sync_mode:
-            idx = np.argmin([child.aq_U for child in children])
+            print("forcing best node update")
+            idx = np.argmin([child.aq_L for child in children])
             self.best_node = children[idx]
-            self.LUB = self.best_node.aq_U
+            self.LLB = self.best_node.aq_L
+            updated_best_node = True
         else:
           print("best node updated")
+        gap_history.append(self.best_node.aq_U - self.best_node.aq_L)
 
         # pre-prune
         children_lower_bounds = [child.aq_L for child in children]
@@ -1387,7 +1392,6 @@ def compute_acqf_bounds(self, l, u, skip_LB=False):
     acqf_solve_success = False 
     if not self.acqf_UB_solver == "MINEVAL": # local gradient-based optimization method
       constraints = []
-      #box_bounds = [[l[i], u[i]] for i in range(len(l))]
       box_bounds = np.array([l, u]).T
       acqf_callback = {'obj' : self.acqf.scalar_evaluate}
       if self.acqf.has_gradient:

src/hiopbbpy/opt/boalgorithm.py

@@ -8,13 +8,16 @@
 import numpy as np
 from numpy.random import uniform
 from scipy.stats import qmc
+from scipy.optimize import minimize
 from ..surrogate_modeling.gp import GaussianProcess
 from .acquisition import LCBacquisition, EIacquisition
 from ..problems.problem import Problem
 from ..utils.util import Evaluator, Logger
 from .bnbalgorithm import BnBAlgorithm
 from .opt_utils import minimizer_wrapper
 
+
+
 # A base class defining a general framework for Bayesian Optimization
 class BOAlgorithmBase:
   def __init__(self):

src/hiopbbpy/problems/BraninProblem.py

@@ -29,7 +29,7 @@ def _evaluate(self, x: np.ndarray) -> np.ndarray:
     ne, nx = x.shape
     assert nx == self.ndim
 
-    y = np.zeros((ne, 1), complex)
+    y = np.zeros((ne, 1))
     b = 5.1 / (4.0 * (np.pi) ** 2)
     c = 5.0 / np.pi
     r = 6.0

src/hiopbbpy/surrogate_modeling/krg.py

@@ -12,23 +12,26 @@
 
 
 class smtKRG(GaussianProcess):
-  def __init__(self, theta, xlimits, ndim, corr="pow_exp", pow_exp_power=1.0, noise0=0.0, nugget = 100. * np.finfo(np.double).eps, random_state=None, hyper_opt="TNC", eval_noise=False):
+  def __init__(self, theta0, xlimits, ndim, corr="pow_exp", pow_exp_power=1.0, noise0=0.0, nugget = 100. * np.finfo(np.double).eps, random_state=None, hyper_opt="TNC", eval_noise=False, fix_theta = False, theta_bounds=[0.1,10.]):
     super().__init__(ndim, xlimits)
     if random_state is None:
       random_state = 42
     design_space = DesignSpace(xlimits, random_state=random_state)
 
-    theta_bounds = [(10.0e0)**(-1.), 10.0e0]
+    if fix_theta:
+      theta_bounds_ = [theta0-1.e-8, theta0+1.e-8]
+    else:
+      theta_bounds_ = theta_bounds
     self.surrogatesmt = KRG(design_space=design_space,
                             print_global=False,
                             noise0=[noise0],
                             eval_noise=False,
                             corr=corr,
                             pow_exp_power=pow_exp_power,
                             hyper_opt=hyper_opt,
-                            theta0 = [theta],
+                            theta0 = [theta0] * ndim,
                             nugget=nugget,
-                            theta_bounds=theta_bounds,
+                            theta_bounds=theta_bounds_,
                             )
     self.trained = False