ROL
ROL_QuantileRadius.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_QUANTILERADIUSQUADRANGLE_HPP
11 #define ROL_QUANTILERADIUSQUADRANGLE_HPP
12 
14 #include "ROL_PlusFunction.hpp"
15 
16 #include "ROL_ParameterList.hpp"
17 
18 namespace ROL {
19 
20 template<class Real>
21 class QuantileRadius : public RandVarFunctional<Real> {
22 private:
23  Ptr<PlusFunction<Real> > plusFunction_;
24  Real prob_;
25  Real coeff_;
26  std::vector<Real> vec_;
27 
33 
36 
41 
42  void initializeQR(void) {
43  Real zero(0);
44  // Initialize temporary storage
45  vec_.clear(); vec_.resize(2,zero);
46  }
47 
48  void checkInputs(void) {
49  Real zero(0), one(1);
50  // Check inputs
51  ROL_TEST_FOR_EXCEPTION((prob_>one || prob_<zero), std::invalid_argument,
52  ">>> ERROR (ROL::QuantileRadius): Confidence level out of range!");
53  ROL_TEST_FOR_EXCEPTION((coeff_<zero), std::invalid_argument,
54  ">>> ERROR (ROL::QuantileRadius): Coefficient is negative!");
55  initializeQR();
56  }
57 
58 public:
59 
60  QuantileRadius( ROL::ParameterList &parlist )
61  : RandVarFunctional<Real>() {
62  ROL::ParameterList &list
63  = parlist.sublist("SOL").sublist("Risk Measure").sublist("Quantile Radius");
64  // Grab probability and coefficient arrays
65  prob_ = list.get<Real>("Confidence Level");
66  coeff_ = list.get<Real>("Coefficient");
67  // Build (approximate) plus function
68  plusFunction_ = makePtr<PlusFunction<Real>>(list);
69  checkInputs();
70  }
71 
72  QuantileRadius(const Real prob, const Real coeff,
73  const Ptr<PlusFunction<Real> > &pf)
74  : RandVarFunctional<Real>(), plusFunction_(pf), prob_(prob), coeff_(coeff) {
75  checkInputs();
76  }
77 
78  void initialize(const Vector<Real> &x) {
80  vec_.assign(2,static_cast<Real>(0));
81  }
82 
83  Real computeStatistic(const Ptr<const std::vector<Real>> &xstat) const override {
84  Real stat(0), half(0.5);
85  if (xstat != nullPtr) {
86  stat = half*((*xstat)[0] + (*xstat)[1]);
87  }
88  return stat;
89  }
90 
92  const Vector<Real> &x,
93  const std::vector<Real> &xstat,
94  Real &tol) {
95  const Real half(0.5), one(1);
96  Real val = computeValue(obj,x,tol);
97  Real pf1 = plusFunction_->evaluate(val-xstat[0],0);
98  Real pf2 = plusFunction_->evaluate(-val-xstat[1],0);
99  RandVarFunctional<Real>::val_ += weight_*(val + half*coeff_/(one-prob_)*(pf1 + pf2));
100  }
101 
102  Real getValue(const Vector<Real> &x,
103  const std::vector<Real> &xstat,
104  SampleGenerator<Real> &sampler) {
105  const Real half(0.5);
106  Real cvar(0);
107  sampler.sumAll(&val_,&cvar,1);
108  cvar += half*coeff_*(xstat[0] + xstat[1]);
109  return cvar;
110  }
111 
113  const Vector<Real> &x,
114  const std::vector<Real> &xstat,
115  Real &tol) {
116  const Real half(0.5), one(1);
117  Real val = computeValue(obj,x,tol);
118  Real pf1 = plusFunction_->evaluate(val-xstat[0],1);
119  Real pf2 = plusFunction_->evaluate(-val-xstat[1],1);
120  Real c = half*weight_*coeff_/(one-prob_);
121  vec_[0] -= c*pf1;
122  vec_[1] -= c*pf2;
123  computeGradient(*dualVector_,obj,x,tol);
124  g_->axpy(weight_ + c * (pf1 - pf2),*dualVector_);
125  }
126 
128  std::vector<Real> &gstat,
129  const Vector<Real> &x,
130  const std::vector<Real> &xstat,
131  SampleGenerator<Real> &sampler) {
132  const Real half(0.5);
133  sampler.sumAll(&vec_[0],&gstat[0],2);
134  sampler.sumAll(*g_,g);
135  gstat[0] += half*coeff_;
136  gstat[1] += half*coeff_;
137  }
138 
140  const Vector<Real> &v,
141  const std::vector<Real> &vstat,
142  const Vector<Real> &x,
143  const std::vector<Real> &xstat,
144  Real &tol) {
145  const Real half(0.5), one(1);
146  Real val = computeValue(obj,x,tol);
147  Real pf11 = plusFunction_->evaluate(val-xstat[0],1);
148  Real pf12 = plusFunction_->evaluate(val-xstat[0],2);
149  Real pf21 = plusFunction_->evaluate(-val-xstat[1],1);
150  Real pf22 = plusFunction_->evaluate(-val-xstat[1],2);
151  Real c = half*weight_*coeff_/(one-prob_);
152  Real gv = computeGradVec(*dualVector_,obj,v,x,tol);
153  vec_[0] -= c*pf12*(gv-vstat[0]);
154  vec_[1] += c*pf22*(gv+vstat[1]);
155  hv_->axpy(c*(pf12*(gv-vstat[0]) + pf22*(gv+vstat[1])),*dualVector_);
156  computeHessVec(*dualVector_,obj,v,x,tol);
157  hv_->axpy(weight_ + c * (pf11 - pf21),*dualVector_);
158  }
159 
161  std::vector<Real> &hvstat,
162  const Vector<Real> &v,
163  const std::vector<Real> &vstat,
164  const Vector<Real> &x,
165  const std::vector<Real> &xstat,
166  SampleGenerator<Real> &sampler) {
167  sampler.sumAll(&vec_[0],&hvstat[0],2);
168  sampler.sumAll(*hv_,hv);
169  }
170 };
171 
172 }
173 
174 #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)
Ptr< Vector< Real > > g_
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.
Real computeValue(Objective< Real > &obj, const Vector< Real > &x, Real &tol)
Ptr< Vector< Real > > hv_
QuantileRadius(const Real prob, const Real coeff, const Ptr< PlusFunction< Real > > &pf)
Ptr< Vector< Real > > dualVector_
Real computeStatistic(const Ptr< const std::vector< Real >> &xstat) const override
Compute statistic.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
void sumAll(Real *input, Real *output, int dim) const
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Real getValue(const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure value.
std::vector< Real > vec_
Ptr< PlusFunction< Real > > plusFunction_
void computeGradient(Vector< Real > &g, Objective< Real > &obj, const Vector< Real > &x, Real &tol)
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.
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...
void updateGradient(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal risk measure storage for gradient computation.
QuantileRadius(ROL::ParameterList &parlist)
void updateValue(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal storage for value computation.
void initialize(const Vector< Real > &x)
Initialize temporary variables.
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.
virtual void initialize(const Vector< Real > &x)
Initialize temporary variables.