ROL
ROL_lDFP.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_LDFP_H
11 #define ROL_LDFP_H
12 
17 #include "ROL_Secant.hpp"
18 
19 namespace ROL {
20 
21 template<class Real>
22 class lDFP : public Secant<Real> {
23 private:
27 
28 public:
29  lDFP(int M, bool useDefaultScaling = true, Real Bscaling = Real(1))
30  : Secant<Real>(M,useDefaultScaling,Bscaling) {}
31 
32  // Apply lBFGS Approximate Inverse Hessian
33  void applyH( Vector<Real> &Hv, const Vector<Real> &v ) const {
34  const Real one(1);
35 
36  // Apply initial Hessian approximation to v
37  applyH0(Hv,v);
38 
39  std::vector<Ptr<Vector<Real>>> a(state_->current+1);
40  std::vector<Ptr<Vector<Real>>> b(state_->current+1);
41  Real bv(0), av(0), bs(0), as(0);
42  for (int i = 0; i <= state_->current; i++) {
43  b[i] = Hv.clone();
44  b[i]->set(*(state_->iterDiff[i]));
45  b[i]->scale(1.0/sqrt(state_->product[i]));
46  //bv = b[i]->dot(v.dual());
47  bv = b[i]->apply(v);
48  Hv.axpy(bv,*b[i]);
49 
50  a[i] = Hv.clone();
51  applyH0(*a[i],*(state_->gradDiff[i]));
52 
53  for (int j = 0; j < i; j++) {
54  //bs = b[j]->dot((state_->gradDiff[i])->dual());
55  bs = b[j]->apply(*(state_->gradDiff[i]));
56  a[i]->axpy(bs,*b[j]);
57  //as = a[j]->dot((state_->gradDiff[i])->dual());
58  as = a[j]->apply(*(state_->gradDiff[i]));
59  a[i]->axpy(-as,*a[j]);
60  }
61  //as = a[i]->dot((state_->gradDiff[i])->dual());
62  as = a[i]->apply(*(state_->gradDiff[i]));
63  a[i]->scale(one/sqrt(as));
64  //av = a[i]->dot(v.dual());
65  av = a[i]->apply(v);
66  Hv.axpy(-av,*a[i]);
67  }
68  }
69 
70  // Apply Initial Secant Approximate Hessian
71  virtual void applyH0( Vector<Real> &Hv, const Vector<Real> &v ) const {
72  Hv.set(v.dual());
73  if (useDefaultScaling_) {
74  if (state_->iter != 0 && state_->current != -1) {
75  Real ss = state_->iterDiff[state_->current]->dot(*(state_->iterDiff[state_->current]));
76  Hv.scale(state_->product[state_->current]/ss);
77  }
78  }
79  else {
80  Hv.scale(static_cast<Real>(1)/Bscaling_);
81  }
82  }
83 
84  // Apply lBFGS Approximate Hessian
85  void applyB( Vector<Real> &Bv, const Vector<Real> &v ) const {
86  const Real zero(0);
87 
88  Bv.set(v.dual());
89  std::vector<Real> alpha(state_->current+1,zero);
90  for (int i = state_->current; i>=0; i--) {
91  alpha[i] = state_->gradDiff[i]->dot(Bv);
92  alpha[i] /= state_->product[i];
93  Bv.axpy(-alpha[i],(state_->iterDiff[i])->dual());
94  }
95 
96  // Apply initial inverse Hessian approximation to v
97  Ptr<Vector<Real>> tmp = Bv.clone();
98  applyB0(*tmp,Bv.dual());
99  Bv.set(*tmp);
100 
101  Real beta(0);
102  for (int i = 0; i <= state_->current; i++) {
103  //beta = state_->iterDiff[i]->dot(Bv.dual());
104  beta = state_->iterDiff[i]->apply(Bv);
105  beta /= state_->product[i];
106  Bv.axpy((alpha[i]-beta),*(state_->gradDiff[i]));
107  }
108  }
109 
110  // Apply Initial Secant Approximate Hessian
111  virtual void applyB0( Vector<Real> &Bv, const Vector<Real> &v ) const {
112  Bv.set(v.dual());
113  if (useDefaultScaling_) {
114  if (state_->iter != 0 && state_->current != -1) {
115  Real ss = state_->iterDiff[state_->current]->dot(*(state_->iterDiff[state_->current]));
116  Bv.scale(ss/state_->product[state_->current]);
117  }
118  }
119  else {
120  Bv.scale(Bscaling_);
121  }
122  }
123 };
124 
125 }
126 
127 #endif
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:119
Provides definitions for limited-memory DFP operators.
Definition: ROL_lDFP.hpp:22
virtual void applyH(Vector< Real > &Hv, const Vector< Real > &v) const =0
bool useDefaultScaling_
Definition: ROL_Secant.hpp:50
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
Real Bscaling_
Definition: ROL_Secant.hpp:51
useDefaultScaling
Definition: ROL_lDFP.hpp:30