ROL
ROL_ExpectationQuadDeviation.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_EXPECTATIONQUADDEVIATION_HPP
11 #define ROL_EXPECTATIONQUADDEVIATION_HPP
12 
13 #include "ROL_ExpectationQuad.hpp"
15 #include "ROL_Types.hpp"
16 
52 namespace ROL {
53 
54 template<class Real>
56 private:
57  Ptr<ExpectationQuad<Real>> eq_;
58 
64 
67 
72 
73 public:
75  : RandVarFunctional<Real>(), eq_(eq) {}
76 
79  void checkRegret(void) {
80  eq_->check();
81  }
82 
84  const Vector<Real> &x,
85  const std::vector<Real> &xstat,
86  Real &tol) {
87  Real val = computeValue(obj,x,tol);
88  Real r = eq_->error(val-xstat[0],0);
89  val_ += weight_ * r;
90  }
91 
93  const Vector<Real> &x,
94  const std::vector<Real> &xstat,
95  Real &tol) {
96  Real val = computeValue(obj,x,tol);
97  Real r = eq_->error(val-xstat[0],1);
98  if (std::abs(r) >= ROL_EPSILON<Real>()) {
99  val_ -= weight_ * r;
100  computeGradient(*dualVector_,obj,x,tol);
101  g_->axpy(weight_*r,*dualVector_);
102  }
103  }
104 
106  const Vector<Real> &v,
107  const std::vector<Real> &vstat,
108  const Vector<Real> &x,
109  const std::vector<Real> &xstat,
110  Real &tol) {
111  Real val = computeValue(obj,x,tol);
112  Real r1 = eq_->error(val-xstat[0],1);
113  Real r2 = eq_->error(val-xstat[0],2);
114  if (std::abs(r2) >= ROL_EPSILON<Real>()) {
115  Real gv = computeGradVec(*dualVector_,obj,v,x,tol);
116  val_ += weight_ * r2 * (vstat[0] - gv);
117  hv_->axpy(weight_*r2*(gv-vstat[0]),*dualVector_);
118  }
119  if (std::abs(r1) >= ROL_EPSILON<Real>()) {
120  computeHessVec(*dualVector_,obj,v,x,tol);
121  hv_->axpy(weight_*r1,*dualVector_);
122  }
123  }
124 
125  Real getValue(const Vector<Real> &x,
126  const std::vector<Real> &xstat,
127  SampleGenerator<Real> &sampler) {
128  Real val(0);
129  sampler.sumAll(&val_,&val,1);
130  return val;
131  }
132 
134  std::vector<Real> &gstat,
135  const Vector<Real> &x,
136  const std::vector<Real> &xstat,
137  SampleGenerator<Real> &sampler) {
138  Real stat(0);
139  sampler.sumAll(&val_,&stat,1);
140  gstat[0] = stat;
141  sampler.sumAll(*g_,g);
142  }
143 
145  std::vector<Real> &hvstat,
146  const Vector<Real> &v,
147  const std::vector<Real> &vstat,
148  const Vector<Real> &x,
149  const std::vector<Real> &xstat,
150  SampleGenerator<Real> &sampler) {
151  Real stat(0);
152  sampler.sumAll(&val_,&stat,1);
153  hvstat[0] = stat;
154  sampler.sumAll(*hv_,hv);
155  }
156 };
157 
158 }
159 
160 #endif
Provides the interface to evaluate objective functions.
void computeHessVec(Vector< Real > &hv, Objective< Real > &obj, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
void checkRegret(void)
Run derivative tests for the scalar regret function.
Ptr< Vector< Real > > g_
Provides a general interface for risk and error measures generated through the expectation risk quadr...
Real computeValue(Objective< Real > &obj, const Vector< Real > &x, Real &tol)
Ptr< Vector< Real > > hv_
void updateHessVec(Objective< Real > &obj, const Vector< Real > &v, const std::vector< Real > &vstat, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal risk measure storage for Hessian-time-a-vector computation.
void updateValue(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal storage for value computation.
Contains definitions of custom data types in ROL.
Ptr< Vector< Real > > dualVector_
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
void sumAll(Real *input, Real *output, int dim) const
void updateGradient(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal risk measure storage for gradient computation.
void computeGradient(Vector< Real > &g, Objective< Real > &obj, const Vector< Real > &x, Real &tol)
void getGradient(Vector< Real > &g, std::vector< Real > &gstat, const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure (sub)gradient.
Real computeGradVec(Vector< Real > &g, Objective< Real > &obj, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Provides the interface to implement any functional that maps a random variable to a (extended) real n...
Ptr< ExpectationQuad< Real > > eq_
void getHessVec(Vector< Real > &hv, std::vector< Real > &hvstat, const Vector< Real > &v, const std::vector< Real > &vstat, const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure Hessian-times-a-vector.
ExpectationQuadDeviation(const Ptr< ExpectationQuad< Real >> &eq)
Real getValue(const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure value.