ROL
ROL_PrimalDualInteriorPointOperator.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_PRIMALDUALINTERIORPOINTOPERATOR_H
11 #define ROL_PRIMALDUALINTERIORPOINTOPERATOR_H
12 
13 #include "ROL_LinearOperator.hpp"
14 
15 
16 namespace ROL {
17 
18 
19 
20 template<class Real>
22 
23  typedef Vector<Real> V;
27 
28  static const size_type OPT = 0;
29  static const size_type EQUAL = 1;
30  static const size_type LOWER = 0;
31  static const size_type UPPER = 1;
32 
33  ROL::Ptr<const V> x_; // Optimization vector
34  ROL::Ptr<const V> l_; // constraint multiplier
35 
36  ROL::Ptr<V> scratch_;
37 
38  Real delta_; // Initial correction
39 
40 
41 public:
42 
43  PrimalDualInteriorPointBlock11( ROL::Ptr<OBJ> &obj, ROL::Ptr<CON> &con,
44  const V &x, ROL::Ptr<V> & scratch,
45  Real delta=0 ) :
46  obj_(obj), con_(con), scratch_(scratch), delta_(delta) {
47 
48  const PV &x_pv = dynamic_cast<const PV&>(x);
49 
50  x_ = x_pv.get(OPT);
51  l_ = x_pv.get(EQUAL);
52  }
53 
54  void update( const Vector<Real> &x, bool flag = true, int iter = -1 ) {
55 
56  const PV &x_pv = dynamic_cast<const PV&>(x);
57 
58  x_ = x_pv.get(OPT);
59  l_ = x_pv.get(EQUAL);
60 
61  obj_->update(*x_,flag,true);
62  con_->update(*x_,flag,true);
63  }
64 
65  void apply( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
66 
67 
68 
69  PV &Hv_pv = dynamic_cast<PV&>(Hv);
70  const PV &v_pv = dynamic_cast<const PV&>(v);
71 
72  // output vector components
73  ROL::Ptr<V> Hvx = Hv_pv.get(OPT);
74  ROL::Ptr<V> Hvl = Hv_pv.get(EQUAL);
75 
76  // input vector components
77  ROL::Ptr<const V> vx = v_pv.get(OPT);
78  ROL::Ptr<const V> vl = v_pv.get(EQUAL);
79 
80  obj_->hessVec(*jvx,*vx,*x_,tol);
81  con_->applyAdjointHessian(*scratch_,*l_,*vx,*x_,tol);
82  jvx->plus(*scratch_);
83  con_->applyAdjointJacobian(*scratch_,*vl,*x_,tol);
84  jvx->plus(*scratch_);
85 
86  // Inertial correction
87  if( delta_ != 0 ) {
88  jvx->axpy(delta_,*vx);
89  }
90 
91  con_->applyJacobian(*jvl,*vx,*x_,tol);
92 
93  }
94 
95  void applyInverse( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
96  ROL_TEST_FOR_EXCEPTION( true , std::logic_error,
97  ">>> ERROR (ROL_PrimalDualInteriorPointBlock11, applyInverse): "
98  "Not implemented.");
99  }
100 
101  void setInertia( Real delta ) {
102  delta_ = delta;
103  }
104 
105 
106 }; // class PrimalDualInteriorPointBlock11
107 
108 
109 template<class Real>
111 
112  typedef Vector<Real> V;
114 
115  static const size_type OPT = 0;
116  static const size_type EQUAL = 1;
117  static const size_type LOWER = 0;
118  static const size_type UPPER = 1;
119 
120 public:
121 
122  void apply( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
123 
124 
125  PV &Hv_pv = dynamic_cast<PV&>(Hv);
126  const PV &v_pv = dynamic_cast<const PV&>(v);
127 
128  // output vector components
129  ROL::Ptr<V> Hvx = Hv_pv.get(OPT);
130  ROL::Ptr<V> Hvl = Hv_pv.get(EQUAL);
131 
132  // input vector components
133  ROL::Ptr<const V> vzl = v_pv.get(LOWER);
134  ROL::Ptr<const V> vzu = v_pv.get(UPPER);
135 
136  Hvx->set(*vzu);
137  Hvx->axpy(-1.0,*vzl);
138  Hvl->zero();
139 
140  }
141 
142  void applyInverse( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
143  ROL_TEST_FOR_EXCEPTION( true , std::logic_error,
144  ">>> ERROR (ROL_PrimalDualInteriorPointBlock12, applyInverse): "
145  "Not implemented.");
146  }
147 
148 
149 }; // class PrimalDualInteriorPointBlock12
150 
151 
152 template<class Real>
154 
155  typedef Vector<Real> V;
157 
158  static const size_type OPT = 0;
159  static const size_type EQUAL = 1;
160  static const size_type LOWER = 0;
161  static const size_type UPPER = 1;
162 
163  ROL::Ptr<const V> zl_;
164  ROL::Ptr<const V> zu_;
165 
166 public:
167 
169  const PV &z_pv = dynamic_cast<const PV&>(z);
170  zl_ = z_pv.get(LOWER);
171  zu_ = z_pv.get(UPPER);
172  }
173 
174  void update( const Vector<Real> &z, bool flag = true, int iter = -1 ) {
175  const PV &z_pv = dynamic_cast<const PV&>(z);
176  zl_ = z_pv.get(LOWER);
177  zu_ = z_pv.get(UPPER);
178  }
179 
180  virtual void apply( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
181 
182 
183  PV &Hv_pv = dynamic_cast<PV&>(Hv);
184  const PV &v_pv = dynamic_cast<const PV&>(v);
185 
186  // output vector components
187  ROL::Ptr<V> Hvzl = Hv_pv.get(LOWER);
188  ROL::Ptr<V> Hvzu = Hv_pv.get(UPPER);
189 
190  // input vector components
191  ROL::Ptr<const V> vx = v_pv.get(OPT);
192  ROL::Ptr<const V> vl = v_pv.get(EQUAL);
193 
194  Hvzl->set(*vx);
195  Hvzl->applyBinary(mult_,*zl_);
196 
197  Hvzu->set(*vx);
198  Hvzu->applyBinary(mult_,*zu_);
199  Hvzu->scale(-1.0);
200  }
201 
202  virtual void applyInverse( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
203  ROL_TEST_FOR_EXCEPTION( true , std::logic_error,
204  ">>> ERROR (ROL_PrimalDualInteriorPointBlock21, applyInverse): "
205  "Not implemented.");
206  }
207 
208 }; // class PrimalDualInteriorPointBlock21
209 
210 
211 template<class Real>
213 
214  typedef Vector<Real> V;
217 
218  static const size_type OPT = 0;
219  static const size_type LOWER = 0;
220  static const size_type UPPER = 1;
221 
222  ROL::Ptr<const V> x_;
223  ROL::Ptr<const V> xl_;
224  ROL::Ptr<const V> xu_;
225 
226  Elementwise::Multiply<Real> mult_;
227  Elementwise::Multiply<Real> divinv_;
228 
229 
230 public:
231 
232  PrimalDualInteriorPointBlock22( const ROL::Ptr<BND> &bnd, const Vector<Real> &x ) {
233 
234  const PV &x_pv = dynamic_cast<const PV&>(x);
235 
236  x_ = x_pv.get(OPT);
237  xl_ = bnd.getLowerBound();
238  xu_ = bnd.getUpperBound();
239 
240  }
241 
242  virtual void update( const Vector<Real> &x, bool flag = true, int iter = -1 ) {
243 
244  const PV &x_pv = dynamic_cast<const PV&>(x);
245  x_ = x_pv.get(OPT);
246  }
247 
248  virtual void apply( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
249 
250 
251  PV &Hv_pv = dynamic_cast<PV&>(Hv);
252  const PV &v_pv = dynamic_cast<const PV&>(v);
253 
254  // output vector components
255  ROL::Ptr<V> Hvzl = Hv_pv.get(LOWER);
256  ROL::Ptr<V> Hvzu = Hv_pv.get(UPPER);
257 
258  // input vector components
259  ROL::Ptr<const V> vzl = v_pv.get(LOWER);
260  ROL::Ptr<const V> vzu = v_pv.get(UPPER);
261 
262  Hvzl->set(*x_);
263  Hvzl->axpy(-1.0,*xl_);
264  Hvzl->applyBinary(mult_,*vzl);
265 
266  Hvzu->set(*xu_);
267  Hvzu->axpy(-1.0,*x_);
268  Hvzu->applyBinary(mult_,*vzu);
269 
270  }
271 
272  virtual void applyInverse( Vector<Real> &Hv, const Vector<Real> &v, Real &tol ) const {
273 
274 
275 
276  PV &Hv_pv = dynamic_cast<PV&>(Hv);
277  const PV &v_pv = dynamic_cast<const PV&>(v);
278 
279  // output vector components
280  ROL::Ptr<V> Hvzl = Hv_pv.get(LOWER);
281  ROL::Ptr<V> Hvzu = Hv_pv.get(UPPER);
282 
283  // input vector components
284  ROL::Ptr<const V> vzl = v_pv.get(LOWER);
285  ROL::Ptr<const V> vzu = v_pv.get(UPPER);
286 
287  Hvzl->set(*x_);
288  Hvzl->axpy(-1.0,*xl_);
289  Hvzl->applyBinary(divinv_,*vzl); // Hvzl = vzl/(x-xl)
290 
291  Hvzu->set(*xu_);
292  Hvzu->axpy(-1.0,*x_);
293  Hvzu->applyBinary(divinv_,*vzu); // Hvzu = vzu/(xu-x)
294  }
295 
296 }; // class PrimalDualInteriorPointBlock22
297 
298 
299 } // namespace ROL
300 
301 
302 
303 #endif // ROL_PRIMALDUALINTERIORPOINTOPERATOR_H
304 
Provides the interface to evaluate objective functions.
typename PV< Real >::size_type size_type
PrimalDualInteriorPointBlock22(const ROL::Ptr< BND > &bnd, const Vector< Real > &x)
virtual void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
ROL::Ptr< const Vector< Real > > get(size_type i) const
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update linear operator.
Defines the linear algebra of vector space on a generic partitioned vector.
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update linear operator.
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
void update(const Vector< Real > &z, bool flag=true, int iter=-1)
Update linear operator.
virtual void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
virtual void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
Provides the interface to apply a linear operator.
Provides the interface to apply upper and lower bound constraints.
void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
virtual void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
PrimalDualInteriorPointBlock11(ROL::Ptr< OBJ > &obj, ROL::Ptr< CON > &con, const V &x, ROL::Ptr< V > &scratch, Real delta=0)
const Ptr< Obj > obj_
Defines the general constraint operator interface.