ROL
ROL_PD_HMCR2.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_PD_HMCR2_HPP
11 #define ROL_PD_HMCR2_HPP
12 
14 #include "ROL_Types.hpp"
15 
16 namespace ROL {
17 
18 template<class Real>
19 class PD_HMCR2 : public PD_RandVarFunctional<Real> {
20 private:
21  Real beta_;
23 
24  Ptr<ScalarController<Real>> values_;
25  Ptr<ScalarController<Real>> gradvecs_;
26  Ptr<VectorController<Real>> gradients_;
27  Ptr<VectorController<Real>> hessvecs_;
28 
34 
37 
42 
47 
48  void initializeStorage(void) {
49  values_ = makePtr<ScalarController<Real>>();
50  gradvecs_ = makePtr<ScalarController<Real>>();
51  gradients_ = makePtr<VectorController<Real>>();
52  hessvecs_ = makePtr<VectorController<Real>>();
53 
55  RandVarFunctional<Real>::setHessVecStorage(gradvecs_,hessvecs_);
56  }
57 
58  void clear(void) {
59  gradvecs_->reset();
60  hessvecs_->reset();
61  }
62 
63  void checkInputs(void) {
64  Real zero(0), one(1);
65  ROL_TEST_FOR_EXCEPTION((beta_ < zero) || (beta_ >= one), std::invalid_argument,
66  ">>> ERROR (ROL::PD_HMCR2): Confidence parameter beta is out of range!");
68  }
69 
70 public:
71  PD_HMCR2(const Real beta)
72  : PD_RandVarFunctional<Real>(), beta_(beta) {
73  checkInputs();
74  }
75 
76  void setStorage(const Ptr<ScalarController<Real>> &value_storage,
77  const Ptr<VectorController<Real>> &gradient_storage) {
78  values_ = value_storage;
79  gradients_ = gradient_storage;
81  }
82 
83  void setHessVecStorage(const Ptr<ScalarController<Real>> &gradvec_storage,
84  const Ptr<VectorController<Real>> &hessvec_storage) {
85  gradvecs_ = gradvec_storage;
86  hessvecs_ = hessvec_storage;
88  }
89 
91  const Real zero(0), two(2);
92  Real val(0), lold(0), lnew(0), mdiff(0), gdiff(0), sum(0), gsum(0);
93  for (int i = sampler.start(); i < sampler.numMySamples(); ++i) {
94  values_->get(val, sampler.getMyPoint(i));
95  getMultiplier(lold, sampler.getMyPoint(i));
96  lnew = std::max(zero, getPenaltyParameter()*val+lold);
97  sum += sampler.getMyWeight(i) * std::pow(lnew,two);
98  }
99  sampler.sumAll(&sum,&gsum,1);
100  gsum = std::sqrt(gsum);
101  for (int i = sampler.start(); i < sampler.numMySamples(); ++i) {
102  values_->get(val, sampler.getMyPoint(i));
103  getMultiplier(lold, sampler.getMyPoint(i));
104  lnew = std::max(zero, getPenaltyParameter()*val+lold)/gsum;
105  mdiff += sampler.getMyWeight(i) * std::pow(lnew-lold,2);
106  setMultiplier(lnew, sampler.getMyPoint(i));
107  }
108  sampler.sumAll(&mdiff,&gdiff,1);
109  gdiff = std::sqrt(gdiff);
110  return gdiff;
111  }
112 
113  void initialize(const Vector<Real> &x) {
115  mScalar1_ = static_cast<Real>(0);
116  mScalar2_ = static_cast<Real>(0);
117  clear();
118  }
119 
121  const Vector<Real> &x,
122  const std::vector<Real> &xstat,
123  Real &tol) {
124  const Real zero(0), two(2);
125  Real lam(0);
126  getMultiplier(lam, point_);
127  Real val = computeValue(obj, x, tol);
128  Real arg = val - xstat[0];
129  Real pf = std::max(zero, arg + lam/getPenaltyParameter());
130  val_ += weight_ * std::pow(pf,two);
131  setValue(arg, point_);
132  }
133 
134  Real getValue(const Vector<Real> &x,
135  const std::vector<Real> &xstat,
136  SampleGenerator<Real> &sampler) {
137  const Real half(0.5), one(1);
138  Real ev(0);
139  sampler.sumAll(&val_, &ev, 1);
140  Real norm = std::sqrt(ev);
141  Real sig = one/(one-beta_);
142  Real val = (norm <= sig/getPenaltyParameter()
143  ? half * getPenaltyParameter() * ev
144  : sig * (norm - sig*half/getPenaltyParameter()));
145  return xstat[0] + val;
146  }
147 
149  const Vector<Real> &x,
150  const std::vector<Real> &xstat,
151  Real &tol) {
152  const Real zero(0), two(2);
153  Real lam(0);
154  getMultiplier(lam, point_);
155  Real val = computeValue(obj, x, tol);
156  Real arg = val - xstat[0];
157  Real pf = std::max(zero, arg + lam/getPenaltyParameter());
158  if ( pf > zero ) {
159  val_ += weight_ * pf;
160  gv_ += weight_ * std::pow(pf,two);
161  computeGradient(*dualVector_, obj, x, tol);
162  g_->axpy(weight_ * pf, *dualVector_);
163  }
164  }
165 
167  std::vector<Real> &gstat,
168  const Vector<Real> &x,
169  const std::vector<Real> &xstat,
170  SampleGenerator<Real> &sampler) {
171  const Real one(1);
172  std::vector<Real> mv = {val_, gv_};
173  std::vector<Real> ev(2,0);
174  sampler.sumAll(&mv[0], &ev[0], 2);
175  Real norm = std::sqrt(ev[1]);
176  Real sig = one/(one-beta_);
177  Real scal = (norm <= sig/getPenaltyParameter()
179  : sig/norm);
180  gstat[0] = one - scal * ev[0];
181  sampler.sumAll(*g_, g);
182  g.scale(scal);
183  }
184 
186  const Vector<Real> &v,
187  const std::vector<Real> &vstat,
188  const Vector<Real> &x,
189  const std::vector<Real> &xstat,
190  Real &tol) {
191  const Real zero(0), two(2);
192  Real lam(0);
193  getMultiplier(lam, point_);
194  Real val = computeValue(obj, x, tol);
195  Real arg = val - xstat[0];
196  Real pf = std::max(zero, arg + lam/getPenaltyParameter());
197  if ( pf > zero ) {
198  val_ += weight_ * std::pow(pf,two);
199  mScalar1_ += weight_ * pf;
200 
201  Real gv = computeGradVec(*dualVector_, obj, v, x, tol);
202  mScalar2_ += weight_ * pf * gv;
203  gv_ += weight_ * (vstat[0] - gv);
204  g_->axpy(weight_ * pf, *dualVector_);
205  hv_->axpy(weight_ * (gv - vstat[0]), *dualVector_);
206  computeHessVec(*dualVector_, obj, v, x, tol);
207  hv_->axpy(weight_ * pf, *dualVector_);
208  }
209  }
210 
212  std::vector<Real> &hvstat,
213  const Vector<Real> &v,
214  const std::vector<Real> &vstat,
215  const Vector<Real> &x,
216  const std::vector<Real> &xstat,
217  SampleGenerator<Real> &sampler) {
218  const Real one(1);
219  std::vector<Real> mv = {val_, gv_, mScalar1_, mScalar2_};
220  std::vector<Real> ev(4,0);
221  sampler.sumAll(&mv[0],&ev[0],4);
222  Real norm = std::sqrt(ev[0]);
223  Real sig = one/(one-beta_);
224  Real scal = (norm <= sig/getPenaltyParameter()
226  : sig/norm);
227  hvstat[0] = scal * ev[1];
228  sampler.sumAll(*hv_,hv);
229  hv.scale(scal);
230  if (norm > sig/getPenaltyParameter()) {
231  Real norm3 = ev[0]*norm;
232  hvstat[0] += sig/norm3 * (ev[3] - ev[2]*vstat[0]) * ev[2];
233  dualVector_->zero();
234  sampler.sumAll(*g_,*dualVector_);
235  hv.axpy(sig/norm3 * (ev[2]*vstat[0] - ev[3]),*dualVector_);
236  }
237  }
238 };
239 
240 }
241 
242 #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)
virtual void scale(const Real alpha)=0
Compute where .
Ptr< Vector< Real > > g_
virtual void setHessVecStorage(const Ptr< ScalarController< Real >> &gradvec_storage, const Ptr< VectorController< Real >> &hessvec_storage)
void setStorage(const Ptr< ScalarController< Real >> &value_storage, const Ptr< VectorController< Real >> &gradient_storage)
Real computeValue(Objective< Real > &obj, const Vector< Real > &x, Real &tol)
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:119
Ptr< Vector< Real > > hv_
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.
Contains definitions of custom data types in ROL.
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 setMultiplier(Real &lam, const std::vector< Real > &pt)
virtual std::vector< Real > getMyPoint(const int i) const
virtual Real getMyWeight(const int i) const
Ptr< Vector< Real > > dualVector_
virtual void setStorage(const Ptr< ScalarController< Real >> &value_storage, const Ptr< VectorController< Real >> &gradient_storage)
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
virtual int numMySamples(void) const
void sumAll(Real *input, Real *output, int dim) const
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Real computeDual(SampleGenerator< Real > &sampler)
PD_HMCR2(const Real beta)
void checkInputs(void)
void initializeStorage(void)
Ptr< ScalarController< Real > > values_
virtual void setStorage(const Ptr< ScalarController< Real >> &value_storage, const Ptr< VectorController< Real >> &gradient_storage)
void initialize(const Vector< Real > &x)
Initialize temporary variables.
Ptr< VectorController< Real > > gradients_
Real getValue(const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure value.
void getMultiplier(Real &lam, const std::vector< Real > &pt) const
void computeGradient(Vector< Real > &g, Objective< Real > &obj, const Vector< Real > &x, Real &tol)
Ptr< VectorController< Real > > hessvecs_
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.
void setHessVecStorage(const Ptr< ScalarController< Real >> &gradvec_storage, const Ptr< VectorController< Real >> &hessvec_storage)
void updateGradient(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal risk measure storage for gradient computation.
virtual void initialize(const Vector< Real > &x)
Initialize temporary variables.
void clear(void)
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 setValue(const Real val, const std::vector< Real > &pt)
void updateValue(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal storage for value computation.
Ptr< ScalarController< Real > > gradvecs_
virtual void setHessVecStorage(const Ptr< ScalarController< Real >> &gradvec_storage, const Ptr< VectorController< Real >> &hessvec_storage)