From f054d7e85af85b67f78e42c156f1af60d10abb9d Mon Sep 17 00:00:00 2001
From: stevenj <stevenj@alum.mit.edu>
Date: Sun, 26 Aug 2007 10:49:57 -0400
Subject: [PATCH] added option to StoGO (currently disabled) to use NLopt LBFGS
 for local minimization instead of StoGO's version

darcs-hash:20070826144957-c8de0-e4a51850c52e7b716ab954be83869a5e714de814.gz
---
 stogo/Makefile.am |  2 +-
 stogo/linalg.h    |  5 ++--
 stogo/local.cc    | 71 ++++++++++++++++++++++++++++++++++++++++++++---
 3 files changed, 71 insertions(+), 7 deletions(-)

diff --git a/stogo/Makefile.am b/stogo/Makefile.am
index dd715b9..3defc18 100644
--- a/stogo/Makefile.am
+++ b/stogo/Makefile.am
@@ -1,4 +1,4 @@
-AM_CPPFLAGS = -I$(top_srcdir)/util
+AM_CPPFLAGS = -I$(top_srcdir)/util -I$(top_srcdir)/api
 
 noinst_LTLIBRARIES = libstogo.la
 libstogo_la_SOURCES = global.cc linalg.cc local.cc stogo.cc tools.cc	\
diff --git a/stogo/linalg.h b/stogo/linalg.h
index 8e889a7..d5b8d7b 100644
--- a/stogo/linalg.h
+++ b/stogo/linalg.h
@@ -25,10 +25,11 @@ double eps() ;
 
 class RVector{
 protected:
-  double*  elements;  // array of values
-  int      len;       // size of array
 
  public:
+  int      len;       // size of array
+  double*  elements;  // array of values
+
   RVector() ;
   RVector(int);       // Constructor
   RVector(RCRVector); // copy constructor
diff --git a/stogo/local.cc b/stogo/local.cc
index c32d010..6f506e2 100644
--- a/stogo/local.cc
+++ b/stogo/local.cc
@@ -19,6 +19,37 @@
 #  define IF_NLOPT_CHECK_EVALS 
 #endif
 
+////////////////////////////////////////////////////////////////////////
+// SGJ, 2007: allow local to use local optimizers in NLopt, to compare
+// to the BFGS code below
+#if 0
+#include "nlopt.h"
+
+typedef struct {
+  Global *glob;
+  double maxgrad;
+  nlopt_stopping *stop;
+} f_local_data;
+
+static double f_local(int n, const double *x, double *grad, void *data_)
+{
+  f_local_data *data = (f_local_data *) data_;
+  double f;
+  RVector xv, gv;
+  // hack to avoid pointless copy of x and grad
+  xv.len = gv.len = n;
+  xv.elements = const_cast<double *>(x);
+  gv.elements = grad;
+  f=data->glob->ObjectiveGradient(xv, gv,
+				   grad?OBJECTIVE_AND_GRADIENT:OBJECTIVE_ONLY);
+  if (grad) data->maxgrad = max(data->maxgrad, normInf(gv));
+  xv.elements = gv.elements = 0; // prevent deallocation
+  data->stop->nevals++;
+  return f;
+}
+#endif
+////////////////////////////////////////////////////////////////////////
+
 int local(Trial &T, TBox &box, TBox &domain, double eps_cl, double *mgr,
           Global &glob, int axis, RCRVector x_av
 #ifdef NLOPT_UTIL_H
@@ -28,7 +59,7 @@ int local(Trial &T, TBox &box, TBox &domain, double eps_cl, double *mgr,
 
   int n=box.GetDim();
   RVector x(n);
-  double tmp;
+  double tmp, f;
 
   x=T.xvals ;
 
@@ -51,10 +82,39 @@ int local(Trial &T, TBox &box, TBox &domain, double eps_cl, double *mgr,
      return LS_Old ;
    } 
 
-  int k_max, info, outside ;
+#ifdef NLOPT_H
+
+  if (axis != -1) {
+    cout << "NLopt code only works with axis == -1, exiting...\n" ;
+    exit(EXIT_FAILURE);
+  }
+  f_local_data data;
+  data.glob = &glob;
+  data.maxgrad = *mgr;
+  data.stop = stop;
+  nlopt_result ret = nlopt_minimize(NLOPT_LOCAL_LBFGS, n, f_local, &data,
+				    box.lb.raw_data(), box.ub.raw_data(),
+				    x.raw_data(), &f, 
+				    stop->fmin_max,
+				    stop->ftol_rel, stop->ftol_abs,
+				    stop->xtol_rel, stop->xtol_abs,
+				    stop->maxeval - stop->nevals,
+				    stop->maxtime - stop->start);
+  *mgr = data.maxgrad;
+  T.xvals=x ; T.objval=f ;
+  if (ret == NLOPT_MAXEVAL_REACHED || ret == NLOPT_MAXTIME_REACHED)
+    return LS_MaxEvalTime;
+  else if (ret > 0)
+    return LS_New;
+  else
+    return LS_Out; // failure
+  
+#else /* not using NLopt local optimizer ... use original STOgo BFGS code */
+
+  int k_max, info, outside = 0;
   int k, i, good_enough, iTmp ;
 
-  double maxgrad, delta, f, f_new;
+  double maxgrad, delta, f_new;
   double alpha, gamma, beta, d2, s2, nom, den, ro ;
   double nrm_sd, nrm_hn, snrm_hn, nrm_dl ;
   RVector g(n), h_sd(n), h_dl(n), h_n(n), x_new(n), g_new(n) ;
@@ -337,10 +397,13 @@ int local(Trial &T, TBox &box, TBox &domain, double eps_cl, double *mgr,
     exit(1);
   }
 
-  T.xvals=x ; T.objval=f ;
   *mgr=maxgrad ;
+
+  T.xvals=x ; T.objval=f ;
   if (outside>0)
     return LS_Out ;
   else
     return info ;
+
+#endif
 }
-- 
2.30.2