mmin_bfgs2.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2017, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_GSL_MMIN_BFGS2_H
24 #define O2SCL_GSL_MMIN_BFGS2_H
25 
26 /** \file mmin_bfgs2.h
27  \brief File defining \ref o2scl::mmin_bfgs2
28 */
29 #include <gsl/gsl_blas.h>
30 #include <gsl/gsl_poly.h>
31 #include <gsl/gsl_multimin.h>
32 #include <o2scl/mmin.h>
33 #include <o2scl/cblas.h>
34 
35 #ifndef DOXYGEN_NO_O2NS
36 namespace o2scl {
37 #endif
38 
39  /** \brief Virtual base for the mmin_bfgs2 wrapper
40 
41  This class is useful so that the gsl_mmin_linmin class doesn't
42  need to depend on any template parameters, even though it will
43  need a wrapping object as an argument for the
44  gsl_mmin_linmin::min() function.
45  */
46  class mmin_wrap_gsl {
47  public:
48  virtual ~mmin_wrap_gsl() {}
49  /// Function
50  virtual double wrap_f(double alpha)=0;
51  /// Derivative
52  virtual double wrap_df(double alpha)=0;
53  /// Function and derivative
54  virtual void wrap_fdf(double alpha, double *f, double *df)=0;
55  };
56 
57  /** \brief Wrapper class for the mmin_bfgs2 minimizer
58 
59  This is a reimplmentation of the internal GSL wrapper for
60  function calls in the BFGS minimizer
61  */
62  template<class func_t=multi_funct11,
64  class dfunc_t=grad_funct11,
65  class auto_grad_t=
68 
69 #ifndef DOXYGEN_INTERNAL
70 
71  protected:
72 
73  /// Function
74  func_t *func;
75 
76  /// Derivative
77  dfunc_t *dfunc;
78 
79  /// The automatic gradient object
80  auto_grad_t *agrad;
81 
82  /// True if the gradient was given by the user
83  bool grad_given;
84 
85  /** \name fixed values
86  */
87  //@{
88  vec_t *x;
89  vec_t *g;
90  vec_t *p;
91  //@}
92 
93  /** \name cached values, for x(alpha)=x + alpha * p
94  */
95  //@{
96  double f_alpha;
97  double df_alpha;
98  //@}
99 
100  /** \name cache keys
101  */
102  //@{
103  double f_cache_key;
104  double df_cache_key;
105  double x_cache_key;
106  double g_cache_key;
107  //@}
108 
109  /// Move to a new point, using the cached value if possible
110  void moveto(double alpha) {
111 
112  if (alpha == x_cache_key) {
113  /* using previously cached position */
114  return;
115  }
116 
117  /* set x_alpha=x + alpha * p */
118 
119  for(size_t i=0;i<dim;i++) {
120  av_x_alpha[i]=(*x)[i];
121  }
122  o2scl_cblas::daxpy(alpha,dim,*p,av_x_alpha);
123 
124  x_cache_key=alpha;
125  }
126 
127  /// Compute the slope
128  double slope() {
129  return o2scl_cblas::ddot(dim,av_g_alpha,*p);
130  }
131 
132  /// Evaluate the function
133  virtual double wrap_f(double alpha) {
134  if (alpha==f_cache_key) {
135  return f_alpha;
136  }
137  moveto(alpha);
138  f_alpha=(*func)(dim,av_x_alpha);
139 
140  f_cache_key=alpha;
141  return f_alpha;
142  }
143 
144  /// Evaluate the derivative
145  virtual double wrap_df(double alpha) {
146 
147  /* using previously cached df(alpha) */
148  if (alpha==df_cache_key) return df_alpha;
149 
150  moveto(alpha);
151  if (alpha!=g_cache_key) {
152  if (grad_given) {
153  (*dfunc)(dim,av_x_alpha,av_g_alpha);
154  } else {
155  (*agrad)(dim,av_x_alpha,av_g_alpha);
156  }
157  g_cache_key=alpha;
158  }
159  df_alpha=slope();
160  df_cache_key=alpha;
161 
162  return df_alpha;
163  }
164 
165  /// Evaluate the function and the derivative
166  virtual void wrap_fdf(double alpha, double *f, double *df) {
167 
168  /* Check for previously cached values */
169  if (alpha == f_cache_key && alpha == df_cache_key) {
170  *f=f_alpha;
171  *df=df_alpha;
172  return;
173  }
174  if (alpha == f_cache_key || alpha == df_cache_key) {
175  *f=wrap_f(alpha);
176  *df=wrap_df(alpha);
177  return;
178  }
179 
180  moveto(alpha);
181  f_alpha=(*func)(dim,av_x_alpha);
182  if (grad_given) {
183  (*dfunc)(dim,av_x_alpha,av_g_alpha);
184  } else {
185  (*agrad)(dim,av_x_alpha,av_g_alpha);
186  }
187  f_cache_key=alpha;
188  g_cache_key=alpha;
189 
190  df_alpha=slope();
191  df_cache_key=alpha;
192 
193  *f=f_alpha;
194  *df=df_alpha;
195 
196  return;
197  }
198 
199 #endif
200 
201  public:
202 
203  /// Temporary storage
204  vec_t av_x_alpha;
205 
206  /// Temporary storage
207  vec_t av_g_alpha;
208 
209  /// Number of minimization dimensions
210  size_t dim;
211 
212  /// Initialize wrapper
213  void prepare_wrapper(func_t &ufunc, dfunc_t *udfunc, vec_t &t_x,
214  double f, vec_t &t_g, vec_t &t_p, auto_grad_t *ag) {
215 
216  func=&ufunc;
217  dfunc=udfunc;
218  if (dfunc!=0) grad_given=true;
219  else grad_given=false;
220  agrad=ag;
221 
222  x=&t_x;
223  g=&t_g;
224  p=&t_p;
225 
226  x_cache_key=0.0;
227  f_cache_key=0.0;
228  g_cache_key=0.0;
229  df_cache_key=0.0;
230 
231  for(size_t i=0;i<dim;i++) {
232  av_x_alpha[i]=(*x)[i];
233  av_g_alpha[i]=(*g)[i];
234  }
235 
236  f_alpha=f;
237  df_alpha=slope();
238 
239  return;
240  }
241 
242  /// Update position
243  void update_position(double alpha, vec_t &t_x, double *t_f, vec_t &t_g) {
244 
245  /* ensure that everything is fully cached */
246  {
247  double t_f_alpha, t_df_alpha;
248  wrap_fdf(alpha,&t_f_alpha,&t_df_alpha);
249  }
250 
251  *t_f=f_alpha;
252  for(size_t i=0;i<dim;i++) {
253  t_x[i]=av_x_alpha[i];
254  t_g[i]=av_g_alpha[i];
255  }
256  }
257 
258  /** \brief Convert cache values to the new minimizer direction
259 
260  Convert the cache values from the end of the current minimisation
261  to those needed for the start of the next minimisation, alpha=0
262  */
264 
265  /* The new x_alpha for alpha=0 is the current position and the
266  new g_alpha for alpha=0 is the current gradient at the
267  endpoint
268  */
269  for(size_t i=0;i<dim;i++) {
270  av_x_alpha[i]=(*x)[i];
271  av_g_alpha[i]=(*g)[i];
272  }
273  x_cache_key=0.0;
274  g_cache_key=0.0;
275 
276  /* The function value does not change */
277  f_cache_key=0.0;
278 
279  /* Calculate the slope along the new direction vector, p */
280  df_alpha=slope();
281  df_cache_key=0.0;
282 
283  return;
284  }
285 
286  };
287 
288  /** \brief The line minimizer for mmin_bfgs2
289  */
291 
292 #ifndef DOXYGEN_INTERNAL
293 
294  protected:
295 
296  /** \brief Minimize the interpolating quadratic
297 
298  Find a minimum in x=[0,1] of the interpolating quadratic through
299  (0,f0) (1,f1) with derivative fp0 at x=0. The interpolating
300  polynomial is q(x)=f0 + fp0 * z + (f1-f0-fp0) * z^2
301  */
302  double interp_quad(double f0, double fp0, double f1, double zl,
303  double zh);
304 
305  /** \brief Minimize the interpolating cubic
306 
307  Find a minimum in x=[0,1] of the interpolating cubic through
308  (0,f0) (1,f1) with derivatives fp0 at x=0 and fp1 at x=1.
309 
310  The interpolating polynomial is:
311 
312  c(x)=f0 + fp0 * z + eta * z^2 + xi * z^3
313 
314  where eta=3*(f1-f0)-2*fp0-fp1, xi=fp0+fp1-2*(f1-f0).
315  */
316  double cubic(double c0, double c1, double c2, double c3, double z);
317 
318  /// Test to see curvature is positive
319  void check_extremum(double c0, double c1, double c2, double c3, double z,
320  double *zmin, double *fmin);
321 
322  /// Interpolate using a cubic
323  double interp_cubic(double f0, double fp0, double f1,
324  double fp1, double zl, double zh);
325 
326  /// Perform the interpolation
327  double interpolate(double a, double fa, double fpa, double b,
328  double fb, double fpb, double xmin, double xmax,
329  int order);
330 #endif
331 
332  public:
333 
334  /** \brief The line minimization
335 
336  Recommended values from \ref Fletcher87 are
337  <tt>rho=0.01, sigma=0.1, tau1=9, tau2=0.05,
338  tau3=0.5 </tt>
339  */
340  int minimize(mmin_wrap_gsl &wrap, double rho, double sigma,
341  double tau1, double tau2, double tau3,
342  int order, double alpha1, double *alpha_new);
343  };
344 
345  /** \brief Multidimensional minimization by the BFGS
346  algorithm (GSL)
347 
348  The functions mmin() and mmin_de() min a given function
349  until the gradient is smaller than the value of mmin::tol_rel
350  (which defaults to \f$ 10^{-4} \f$ ).
351 
352  See an example for the usage of this class in
353  \ref ex_mmin_sect .
354 
355  This class includes the optimizations from the GSL minimizer \c
356  vector_bfgs2.
357 
358  Default template arguments
359  - \c func_t - \ref multi_funct11
360  - \c vec_t - \ref boost::numeric::ublas::vector <double >
361  - \c dfunc_t - \ref mm_funct11
362  - \c auto_grad_t - \ref gradient<func_t,
363  \ref boost::numeric::ublas::vector <double > >
364  - \c def_auto_grad_t - \ref gradient_gsl<func_t,
365  \ref boost::numeric::ublas::vector < double > >
366 
367  \todo While BFGS does well in the \c ex_mmin example with the
368  initial guess of \f$ (1,0,7\pi) \f$ it seems to converge more
369  poorly for the spring function than the other minimizers with
370  other initial guesses, and I think this will happen in the GSL
371  versions too. I need to examine this more closely with some code
372  designed to clearly show this.
373 
374  \future When the bfgs2 line minimizer returns a zero status,
375  the minimization fails. When err_nonconv is false, the
376  minimizer isn't able to update the x vector so the mmin()
377  function doesn't return the best minimum obtained so far.
378  This is a bit confusing, and could be improved.
379  */
380  template<class func_t=multi_funct11,
382  class dfunc_t=grad_funct11,
383  class auto_grad_t=
385  class def_auto_grad_t=
387  class mmin_bfgs2 : public mmin_base<func_t,func_t,vec_t> {
388 
389 #ifndef DOXYGEN_INTERNAL
390 
391  protected:
392 
393  /// \name The original variables from the GSL state structure
394  //@{
395  int iter;
396  double step;
397  double g0norm;
398  double pnorm;
399  double delta_f;
400  /* f'(0) for f(x-alpha*p) */
401  double fp0;
402  vec_t x0;
403  vec_t g0;
404  vec_t p;
405  /* work space */
406  vec_t dx0;
407  vec_t dg0;
408  /* wrapper function */
410  /* minimization parameters */
411  double rho;
412  double sigma;
413  double tau1;
414  double tau2;
415  double tau3;
416  int order;
417  //@}
418 
419  /// The line minimizer
421 
422  /// \name Store the arguments to set() so we can use them for iterate()
423  //@{
424  vec_t *st_x;
425  vec_t st_dx;
426  vec_t st_grad;
427  double st_f;
428  //@}
429 
430  /// Memory size
431  size_t dim;
432 
433  /// Automatic gradient object
434  auto_grad_t *agrad;
435 
436 #endif
437 
438  public:
439 
440  mmin_bfgs2() {
441  // The default values of lmin_tol and step_size from the
442  // example in the GSL documentation
443  lmin_tol=1.0e-4;
444  step_size=0.01;
445  agrad=&def_grad;
446  }
447 
448  virtual ~mmin_bfgs2() {}
449 
450  /// Perform an iteration
451  virtual int iterate() {
452 
453  double alpha=0.0, alpha1;
454 
455  double pg, dir;
456  int status;
457 
458  double f0=st_f;
459 
460  if (pnorm == 0.0 || g0norm == 0.0 || fp0 == 0) {
461  for(size_t i=0;i<this->dim;i++) st_dx[i]=0.0;
462  O2SCL_CONV2("Either pnorm, g0norm, or fp0 vanished in ",
463  "mmin_bfgs2::iterate().",
464  exc_enoprog,this->err_nonconv);
465  // The functions mmin() and mmin_de() use this return value, so
466  // when err_nonconv is false, we still want to return a non-zero
467  // value.
468  return exc_enoprog;
469  }
470 
471  double dbl_eps=std::numeric_limits<double>::epsilon();
472 
473  if (delta_f < 0) {
474  double del;
475  if (-delta_f>10.0*dbl_eps*fabs(f0)) del=-delta_f;
476  else del=10.0*dbl_eps*fabs(f0);
477  if (2.0*del/(-fp0)<1.0) alpha1=2.0*del/(-fp0);
478  else alpha1=1.0;
479  } else {
480  alpha1=fabs(step);
481  }
482 
483  /* line minimisation, with cubic interpolation (order=3) */
484 
485  status=lm.minimize(wrap,rho,sigma,tau1,tau2,tau3,order,
486  alpha1,&alpha);
487 
488  if (status != success) {
489  O2SCL_CONV2("Variable stepb vanished in ",
490  "mmin_bfgs2::iterate().",
491  exc_enoprog,this->err_nonconv);
492 
493  // The functions mmin() and mmin_de() use this return value, so
494  // when err_nonconv is false, we still want to return a non-zero
495  // value.
496  return exc_enoprog;
497  }
498 
499  wrap.update_position(alpha,*st_x,&st_f,st_grad);
500 
501  delta_f=st_f - f0;
502 
503  /* Choose a new direction for the next step */
504 
505  {
506  /* This is the BFGS update: */
507  /* p'=g1 - A dx - B dg */
508  /* A=- (1+ dg.dg/dx.dg) B + dg.g/dx.dg */
509  /* B=dx.g/dx.dg */
510 
511  double dxg, dgg, dxdg, dgnorm, A, B;
512 
513  /* dx0=x - x0 */
514  for(size_t i=0;i<dim;i++) dx0[i]=(*st_x)[i];
515  o2scl_cblas::daxpy(-1.0,dim,x0,dx0);
516 
517  /* keep a copy */
518  for(size_t i=0;i<dim;i++) st_dx[i]=dx0[i];
519 
520  /* dg0=g - g0 */
521  for(size_t i=0;i<dim;i++) dg0[i]=st_grad[i];
522  o2scl_cblas::daxpy(-1.0,dim,g0,dg0);
523 
524  dxg=o2scl_cblas::ddot(dim,dx0,st_grad);
525  dgg=o2scl_cblas::ddot(dim,dg0,st_grad);
526  dxdg=o2scl_cblas::ddot(dim,dx0,dg0);
527 
528  dgnorm=o2scl_cblas::dnrm2(dim,dg0);
529 
530  if (dxdg != 0) {
531  B=dxg / dxdg;
532  A=-(1.0 + dgnorm * dgnorm / dxdg) * B + dgg / dxdg;
533  } else {
534  B=0;
535  A=0;
536  }
537 
538  for(size_t i=0;i<dim;i++) p[i]=st_grad[i];
539  o2scl_cblas::daxpy(-A,dim,dx0,p);
540  o2scl_cblas::daxpy(-B,dim,dg0,p);
541  }
542 
543  for(size_t i=0;i<dim;i++) {
544  g0[i]=st_grad[i];
545  x0[i]=(*st_x)[i];
546  }
547 
548  g0norm=o2scl_cblas::dnrm2(dim,g0);
549  pnorm=o2scl_cblas::dnrm2(dim,p);
550 
551  /* update direction and fp0 */
552 
553  pg=o2scl_cblas::ddot(dim,st_grad,p);
554  dir=(pg >= 0.0) ? -1.0 : +1.0;
555  o2scl_cblas::dscal(dir/pnorm,dim,p);
556  pnorm=o2scl_cblas::dnrm2(dim,p);
557  fp0=o2scl_cblas::ddot(dim,p,g0);
558 
559  wrap.change_direction();
560 
561  return success;
562  }
563 
564  /// Return string denoting type("mmin_bfgs2")
565  virtual const char *type() { return "mmin_bfgs2";}
566 
567  /// Allocate the memory
568  virtual int allocate(size_t n) {
569 
570  p.resize(n);
571  x0.resize(n);
572  g0.resize(n);
573  dx0.resize(n);
574  dg0.resize(n);
575 
576  for(size_t i=0;i<n;i++) {
577  p[i]=0.0;
578  x0[i]=0.0;
579  g0[i]=0.0;
580  dx0[i]=0.0;
581  dg0[i]=0.0;
582  }
583 
584  st_dx.resize(n);
585  st_grad.resize(n);
586 
587  wrap.av_x_alpha.resize(n);
588  wrap.av_g_alpha.resize(n);
589  wrap.dim=n;
590  dim=n;
591 
592  return success;
593  }
594 
595  /// Free the allocated memory
596  virtual int free() {
597  wrap.av_x_alpha.clear();
598  wrap.av_g_alpha.clear();
599  st_grad.clear();
600  dg0.clear();
601  dx0.clear();
602  g0.clear();
603  x0.clear();
604  p.clear();
605  st_dx.clear();
606  wrap.dim=0;
607  dim=0;
608  return 0;
609  }
610 
611  /// Reset the minimizer to use the current point as a new starting point
612  int restart() {
613  iter=0;
614  return success;
615  }
616 
617  /// Set the function and initial guess
618  virtual int set(vec_t &x, double u_step_size, double tol_u,
619  func_t &ufunc) {
620 
621  iter=0;
622  step=u_step_size;
623  delta_f=0;
624 
625  st_x=&x;
626 
627  st_f=ufunc(dim,x);
628  agrad->set_function(ufunc);
629  (*agrad)(dim,x,st_grad);
630 
631  /* Use the gradient as the initial direction */
632 
633  for(size_t i=0;i<dim;i++) {
634  x0[i]=x[i];
635  g0[i]=st_grad[i];
636  }
637  g0norm=o2scl_cblas::dnrm2(dim,g0);
638  for(size_t i=0;i<dim;i++) {
639  p[i]=st_grad[i];
640  }
641  o2scl_cblas::dscal(-1.0/g0norm,dim,p);
642 
643  // pnorm should be 1
644  pnorm=o2scl_cblas::dnrm2(dim,p);
645  fp0=-g0norm;
646 
647  /* Prepare the wrapper */
648 
649  wrap.prepare_wrapper(ufunc,0,x0,st_f,g0,p,agrad);
650 
651  /* Prepare 1d minimisation parameters */
652 
653  rho=0.01;
654  sigma=tol_u;
655  tau1=9;
656  tau2=0.05;
657  tau3=0.5;
658  // use cubic interpolation where possible
659  order=3;
660 
661  return success;
662 
663  }
664 
665  /// Set the function, the gradient, and the initial guess
666  virtual int set_de(vec_t &x, double u_step_size, double tol_u,
667  func_t &ufunc, dfunc_t &udfunc) {
668 
669  iter=0;
670  step=u_step_size;
671  delta_f=0;
672 
673  st_x=&x;
674 
675  st_f=ufunc(dim,x);
676  udfunc(dim,x,st_grad);
677 
678  /* Use the gradient as the initial direction */
679 
680  for(size_t i=0;i<dim;i++) {
681  x0[i]=x[i];
682  g0[i]=st_grad[i];
683  }
684  g0norm=o2scl_cblas::dnrm2(dim,g0);
685  for(size_t i=0;i<dim;i++) {
686  p[i]=st_grad[i];
687  }
688  o2scl_cblas::dscal(-1.0/g0norm,dim,p);
689 
690  // pnorm should be 1
691  pnorm=o2scl_cblas::dnrm2(dim,p);
692  fp0=-g0norm;
693 
694  /* Prepare the wrapper */
695 
696  wrap.prepare_wrapper(ufunc,&udfunc,x0,st_f,g0,p,agrad);
697 
698  /* Prepare 1d minimisation parameters */
699 
700  rho=0.01;
701  sigma=tol_u;
702  tau1=9;
703  tau2=0.05;
704  tau3=0.5;
705  // use cubic interpolation where possible
706  order=3;
707 
708  return success;
709 
710  }
711 
712  /// The size of the first trial step (default 0.01)
713  double step_size;
714 
715  /// The tolerance for the 1-dimensional minimizer
716  double lmin_tol;
717 
718  /// Default automatic gradient object
719  def_auto_grad_t def_grad;
720 
721  /** \brief Calculate the minimum \c min of \c func w.r.t the
722  array \c x of size \c nn.
723  */
724  virtual int mmin(size_t nn, vec_t &xx, double &fmin,
725  func_t &ufunc) {
726 
727  if (nn==0) {
728  O2SCL_ERR2("Tried to min over zero variables ",
729  " in mmin_bfgs2::mmin().",exc_einval);
730  }
731 
732  int xiter=0, status;
733 
734  allocate(nn);
735 
736  set(xx,step_size,lmin_tol,ufunc);
737 
738  do {
739 
740  xiter++;
741 
742  status=iterate();
743 
744  if (status) {
745  break;
746  }
747 
748  // Equivalent to gsl_multimin_test_gradient with
749  // additional code to print out present iteration
750 
751  double norm=o2scl_cblas::dnrm2(nn,st_grad);
752 
753  if(this->verbose>0) {
754  this->print_iter(nn,*st_x,st_f,xiter,
755  norm,this->tol_rel,"mmin_bfgs2");
756  }
757 
758  if (norm<this->tol_rel) {
759  status=success;
760  } else {
761  status=gsl_continue;
762  }
763 
764  } while (status == gsl_continue && xiter < this->ntrial);
765 
766  std::cout << "H: " << (*st_x)[0] << std::endl;
767  for(size_t i=0;i<nn;i++) xx[i]=(*st_x)[i];
768  fmin=st_f;
769 
770  free();
771  this->last_ntrial=xiter;
772 
773  if (status==gsl_continue && xiter==this->ntrial) {
774  O2SCL_CONV_RET("Too many iterations in mmin_bfgs2::mmin().",
775  exc_emaxiter,this->err_nonconv);
776  }
777  return status;
778  }
779 
780  /** \brief Calculate the minimum \c min of \c func w.r.t the
781  array \c x of size \c nn.
782  */
783  virtual int mmin_de(size_t nn, vec_t &xx, double &fmin,
784  func_t &ufunc, dfunc_t &udfunc) {
785 
786  if (nn==0) {
787  O2SCL_ERR2("Tried to min over zero variables ",
788  "in mmin_bfgs2::mmin().",exc_einval);
789  }
790 
791  int xiter=0, status;
792 
793  allocate(nn);
794 
795  set_de(xx,step_size,lmin_tol,ufunc,udfunc);
796 
797  do {
798  xiter++;
799 
800  status=iterate();
801 
802  if (status) {
803  break;
804  }
805 
806  // Equivalent to gsl_multimin_test_gradient with
807  // additional code to print out present iteration
808 
809  double norm=o2scl_cblas::dnrm2(nn,st_grad);
810 
811  if(this->verbose>0) {
812  this->print_iter(nn,*st_x,st_f,xiter,
813  norm,this->tol_rel,"mmin_bfgs2");
814  }
815 
816  if (norm<this->tol_rel) status=success;
817  else status=gsl_continue;
818 
819  } while (status == gsl_continue && xiter < this->ntrial);
820 
821  for(size_t i=0;i<nn;i++) xx[i]=(*st_x)[i];
822  fmin=st_f;
823 
824  free();
825  this->last_ntrial=xiter;
826 
827  if (status==gsl_continue && xiter==this->ntrial) {
828  O2SCL_CONV_RET("Too many iterations in mmin_bfgs2::mmin_de().",
829  exc_emaxiter,this->err_nonconv);
830  }
831 
832  return status;
833  }
834 
835 #ifndef DOXYGEN_INTERNAL
836 
837  private:
838 
842  (const mmin_bfgs2<func_t,vec_t,dfunc_t,auto_grad_t,def_auto_grad_t>&);
843 
844 #endif
845 
846  };
847 
848 #ifndef DOXYGEN_NO_O2NS
849 }
850 #endif
851 
852 #endif
Class for automatically computing gradients [abstract base].
Definition: mmin.h:54
void update_position(double alpha, vec_t &t_x, double *t_f, vec_t &t_g)
Update position.
Definition: mmin_bfgs2.h:243
double step_size
The size of the first trial step (default 0.01)
Definition: mmin_bfgs2.h:713
The main O<span style=&#39;position: relative; top: 0.3em; font-size: 0.8em&#39;>2</span>scl O$_2$scl names...
Definition: anneal.h:42
dfunc_t * dfunc
Derivative.
Definition: mmin_bfgs2.h:77
double lmin_tol
The tolerance for the 1-dimensional minimizer.
Definition: mmin_bfgs2.h:716
#define O2SCL_CONV_RET(d, n, b)
Set a "convergence" error and return the error value.
Definition: err_hnd.h:292
void moveto(double alpha)
Move to a new point, using the cached value if possible.
Definition: mmin_bfgs2.h:110
virtual int mmin_de(size_t nn, vec_t &xx, double &fmin, func_t &ufunc, dfunc_t &udfunc)
Calculate the minimum min of func w.r.t the array x of size nn.
Definition: mmin_bfgs2.h:783
The line minimizer for mmin_bfgs2.
Definition: mmin_bfgs2.h:290
vec_t av_x_alpha
Temporary storage.
Definition: mmin_bfgs2.h:204
invalid argument supplied by user
Definition: err_hnd.h:59
virtual void wrap_fdf(double alpha, double *f, double *df)=0
Function and derivative.
virtual int iterate()
Perform an iteration.
Definition: mmin_bfgs2.h:451
bool grad_given
True if the gradient was given by the user.
Definition: mmin_bfgs2.h:83
exceeded max number of iterations
Definition: err_hnd.h:73
int minimize(mmin_wrap_gsl &wrap, double rho, double sigma, double tau1, double tau2, double tau3, int order, double alpha1, double *alpha_new)
The line minimization.
iteration has not converged
Definition: err_hnd.h:51
void change_direction()
Convert cache values to the new minimizer direction.
Definition: mmin_bfgs2.h:263
virtual int mmin(size_t nn, vec_t &xx, double &fmin, func_t &ufunc)
Calculate the minimum min of func w.r.t the array x of size nn.
Definition: mmin_bfgs2.h:724
virtual int allocate(size_t n)
Allocate the memory.
Definition: mmin_bfgs2.h:568
virtual int set_de(vec_t &x, double u_step_size, double tol_u, func_t &ufunc, dfunc_t &udfunc)
Set the function, the gradient, and the initial guess.
Definition: mmin_bfgs2.h:666
virtual int free()
Free the allocated memory.
Definition: mmin_bfgs2.h:596
double dnrm2(const size_t N, const vec_t &X)
Compute the norm of the vector X.
Definition: cblas_base.h:156
func_t * func
Function.
Definition: mmin_bfgs2.h:74
iteration is not making progress toward solution
Definition: err_hnd.h:105
virtual double wrap_df(double alpha)
Evaluate the derivative.
Definition: mmin_bfgs2.h:145
Multidimensional minimization [abstract base].
Definition: mmin.h:164
vec_t av_g_alpha
Temporary storage.
Definition: mmin_bfgs2.h:207
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
Definition: err_hnd.h:281
def_auto_grad_t def_grad
Default automatic gradient object.
Definition: mmin_bfgs2.h:719
double ddot(const size_t N, const vec_t &X, const vec2_t &Y)
Compute .
Definition: cblas_base.h:123
int restart()
Reset the minimizer to use the current point as a new starting point.
Definition: mmin_bfgs2.h:612
Simple automatic computation of gradient by finite differencing.
Definition: mmin.h:96
virtual double wrap_df(double alpha)=0
Derivative.
virtual double wrap_f(double alpha)
Evaluate the function.
Definition: mmin_bfgs2.h:133
size_t dim
Memory size.
Definition: mmin_bfgs2.h:431
#define O2SCL_CONV2(d, d2, n, b)
Set a "convergence" error, two-string version.
Definition: err_hnd.h:286
virtual void wrap_fdf(double alpha, double *f, double *df)
Evaluate the function and the derivative.
Definition: mmin_bfgs2.h:166
void dscal(const double alpha, const size_t N, vec_t &X)
Compute .
Definition: cblas_base.h:190
auto_grad_t * agrad
Automatic gradient object.
Definition: mmin_bfgs2.h:434
std::function< int(size_t, boost::numeric::ublas::vector< double > &, boost::numeric::ublas::vector< double > &)> grad_funct11
Array of multi-dimensional functions typedef.
Definition: mmin.h:42
double slope()
Compute the slope.
Definition: mmin_bfgs2.h:128
void daxpy(const double alpha, const size_t N, const vec_t &X, vec2_t &Y)
Compute .
Definition: cblas_base.h:98
void prepare_wrapper(func_t &ufunc, dfunc_t *udfunc, vec_t &t_x, double f, vec_t &t_g, vec_t &t_p, auto_grad_t *ag)
Initialize wrapper.
Definition: mmin_bfgs2.h:213
Virtual base for the mmin_bfgs2 wrapper.
Definition: mmin_bfgs2.h:46
Wrapper class for the mmin_bfgs2 minimizer.
Definition: mmin_bfgs2.h:67
Multidimensional minimization by the BFGS algorithm (GSL)
Definition: mmin_bfgs2.h:387
std::function< double(size_t, const boost::numeric::ublas::vector< double > &)> multi_funct11
Multi-dimensional function typedef.
Definition: multi_funct.h:45
Success.
Definition: err_hnd.h:47
virtual const char * type()
Return string denoting type("mmin_bfgs2")
Definition: mmin_bfgs2.h:565
mmin_linmin_gsl lm
The line minimizer.
Definition: mmin_bfgs2.h:420
auto_grad_t * agrad
The automatic gradient object.
Definition: mmin_bfgs2.h:80
virtual double wrap_f(double alpha)=0
Function.
size_t dim
Number of minimization dimensions.
Definition: mmin_bfgs2.h:210

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).