ROL
ROL_lBFGS.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef ROL_LBFGS_H
11 #define ROL_LBFGS_H
12 
17 #include "ROL_Secant.hpp"
18 
19 namespace ROL {
20 
21 template<class Real>
22 class lBFGS : public Secant<Real> {
23 private:
25 
26 public:
27  lBFGS(int M, bool useDefaultScaling = true, Real Bscaling = Real(1))
28  : Secant<Real>(M,useDefaultScaling,Bscaling) {}
29 
30  // Apply lBFGS Approximate Inverse Hessian
31  void applyH( Vector<Real> &Hv, const Vector<Real> &v ) const {
32  const Real zero(0);
33 
34  auto tmp = v.clone();
35  tmp->set(v);
36  std::vector<Real> alpha(state_->current+1,zero);
37  for (int i = state_->current; i>=0; i--) {
38  alpha[i] = state_->iterDiff[i]->apply(*tmp);
39  alpha[i] /= state_->product[i];
40  tmp->axpy(-alpha[i],*state_->gradDiff[i]);
41  }
42 
43  // Apply initial inverse Hessian approximation to v
44  Secant<Real>::applyH0(Hv,*tmp);
45 
46  Real beta(0);
47  for (int i = 0; i <= state_->current; i++) {
48  //beta = Hv.dot((state_->gradDiff[i])->dual());
49  beta = Hv.apply(*state_->gradDiff[i]);
50  beta /= state_->product[i];
51  Hv.axpy((alpha[i]-beta),*(state_->iterDiff[i]));
52  }
53  }
54 
55  // Apply lBFGS Approximate Hessian
56  void applyB( Vector<Real> &Bv, const Vector<Real> &v ) const {
57  const Real one(1);
58 
59  // Apply initial Hessian approximation to v
61 
62  std::vector<Ptr<Vector<Real>>> a(state_->current+1);
63  std::vector<Ptr<Vector<Real>>> b(state_->current+1);
64  Real bv(0), av(0), bs(0), as(0);
65  for (int i = 0; i <= state_->current; i++) {
66  b[i] = Bv.clone();
67  b[i]->set(*(state_->gradDiff[i]));
68  b[i]->scale(one/sqrt(state_->product[i]));
69  //bv = v.dot(b[i]->dual());
70  bv = v.apply(*b[i]);
71  Bv.axpy(bv,*b[i]);
72 
73  a[i] = Bv.clone();
74  Secant<Real>::applyB0(*a[i],*(state_->iterDiff[i]));
75 
76  for (int j = 0; j < i; j++) {
77  //bs = (state_->iterDiff[i])->dot(b[j]->dual());
78  bs = (state_->iterDiff[i])->apply(*b[j]);
79  a[i]->axpy(bs,*b[j]);
80  //as = (state_->iterDiff[i])->dot(a[j]->dual());
81  as = (state_->iterDiff[i])->apply(*a[j]);
82  a[i]->axpy(-as,*a[j]);
83  }
84  //as = (state_->iterDiff[i])->dot(a[i]->dual());
85  as = (state_->iterDiff[i])->apply(*a[i]);
86  a[i]->scale(one/sqrt(as));
87  //av = v.dot(a[i]->dual());
88  av = v.apply(*a[i]);
89  Bv.axpy(-av,*a[i]);
90  }
91  }
92 };
93 
94 }
95 
96 
97 #endif
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual Real apply(const Vector< Real > &x) const
Apply to a dual vector. This is equivalent to the call .
Definition: ROL_Vector.hpp:204
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:119
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
Definition: ROL_Secant.hpp:164
virtual void applyH(Vector< Real > &Hv, const Vector< Real > &v) const =0
Provides definitions for limited-memory BFGS operators.
Definition: ROL_lBFGS.hpp:22
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
virtual void applyB0(Vector< Real > &Bv, const Vector< Real > &v) const
Definition: ROL_Secant.hpp:127
Provides interface for and implements limited-memory secant operators.
Definition: ROL_Secant.hpp:45
virtual void applyB(Vector< Real > &Bv, const Vector< Real > &v) const =0
const Ptr< SecantState< Real > > state_
Definition: ROL_Secant.hpp:48
virtual void applyH0(Vector< Real > &Hv, const Vector< Real > &v) const
Definition: ROL_Secant.hpp:110