44 #ifndef ROL_HMCROBJECTIVE_HPP 45 #define ROL_HMCROBJECTIVE_HPP 47 #include "Teuchos_RCP.hpp" 86 pointGrad_ = x.
dual().clone();
87 pointHess_ = x.
dual().clone();
88 gradient0_ = x.
dual().clone();
89 sumGrad0_ = x.
dual().clone();
90 gradient1_ = x.
dual().clone();
91 sumGrad1_ = x.
dual().clone();
92 gradient2_ = x.
dual().clone();
93 sumGrad2_ = x.
dual().clone();
94 hessvec_ = x.
dual().clone();
95 sumHess_ = x.
dual().clone();
106 if ( !initialized_ ) {
112 const std::vector<Real> ¶m, Real &tol) {
113 if ( storage_ && value_storage_.count(param) ) {
114 val = value_storage_[param];
117 ParametrizedObjective_->setParameter(param);
118 val = ParametrizedObjective_->value(x,tol);
120 value_storage_.insert(std::pair<std::vector<Real>,Real>(param,val));
126 const std::vector<Real> ¶m, Real &tol) {
127 if ( storage_ && gradient_storage_.count(param) ) {
128 g.
set(*(gradient_storage_[param]));
131 ParametrizedObjective_->setParameter(param);
132 ParametrizedObjective_->gradient(g,x,tol);
134 Teuchos::RCP<Vector<Real> > tmp = g.
clone();
135 gradient_storage_.insert(std::pair<std::vector<Real>,Teuchos::RCP<
Vector<Real> > >(param,tmp));
136 gradient_storage_[param]->set(g);
142 const std::vector<Real> ¶m, Real &tol) {
143 ParametrizedObjective_->setParameter(param);
144 ParametrizedObjective_->hessVec(hv,v,x,tol);
152 Real order, Real prob,
156 bool storage =
true )
157 : ParametrizedObjective_(pObj),
158 ValueSampler_(vsampler), GradientSampler_(gsampler), HessianSampler_(hsampler),
159 initialized_(false), storage_(storage) {
160 order_ = ((order < 1.0) ? 1.0 : order);
161 prob_ = ((prob > 1.0) ? 1.0 : ((prob < 0.0) ? 0.0 : prob));
162 value_storage_.clear();
163 gradient_storage_.clear();
167 Real order, Real prob,
170 bool storage =
true )
171 : ParametrizedObjective_(pObj),
172 ValueSampler_(vsampler), GradientSampler_(gsampler), HessianSampler_(gsampler),
173 initialized_(false), storage_(storage) {
174 order_ = ((order < 1.0) ? 1.0 : order);
175 prob_ = ((prob > 1.0) ? 1.0 : ((prob < 0.0) ? 0.0 : prob));
176 value_storage_.clear();
177 gradient_storage_.clear();
181 Real order, Real prob,
183 bool storage =
true )
184 : ParametrizedObjective_(pObj), order_(order), prob_(prob),
185 ValueSampler_(sampler), GradientSampler_(sampler), HessianSampler_(sampler),
186 initialized_(false), storage_(storage) {
187 order_ = ((order < 1.0) ? 1.0 : order);
188 prob_ = ((prob > 1.0) ? 1.0 : ((prob < 0.0) ? 0.0 : prob));
189 value_storage_.clear();
190 gradient_storage_.clear();
194 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
196 ParametrizedObjective_->update(*xvec,flag,iter);
197 ValueSampler_->update(*xvec);
199 value_storage_.clear();
202 GradientSampler_->update(*xvec);
203 HessianSampler_->update(*xvec);
205 gradient_storage_.clear();
211 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
214 std::vector<Real> point;
215 Real weight = 0.0, myval = 0.0, pval = 0.0, val = 0.0;
216 int start = ValueSampler_->start(), end = ValueSampler_->numMySamples();
217 for (
int i = start; i < end; i++ ) {
218 weight = ValueSampler_->getMyWeight(i);
219 point = ValueSampler_->getMyPoint(i);
224 myval += weight*((order_ == 1.0) ? pval-xvar
225 : std::pow(pval-xvar,order_));
229 ValueSampler_->sumAll(&myval,&val,1);
234 return xvar + ((order_ == 1.0) ? val
235 : ((order_ == 2.0) ? std::sqrt(val)
236 : std::pow(val,1.0/order_)))/(1.0-
prob_);
240 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
245 g.
zero(); sumGrad0_->zero();
246 std::vector<Real> point, val(2,0.0), myval(2,0.0);
247 Real weight = 0.0, pval = 0.0, pvalp0 = 0.0, pvalp1 = 0.0;
248 int start = GradientSampler_->start(), end = GradientSampler_->numMySamples();
249 for (
int i = start; i < end; i++ ) {
250 weight = GradientSampler_->getMyWeight(i);
251 point = GradientSampler_->getMyPoint(i);
256 pvalp0 = ((order_ == 1.0) ? pval-xvar
257 : std::pow(pval-xvar,order_));
258 pvalp1 = ((order_ == 1.0) ? 1.0
259 : ((order_ == 2.0) ? pval-xvar
260 : std::pow(pval-xvar,order_-1.0)));
262 myval[0] += weight*pvalp0;
263 myval[1] += weight*pvalp1;
267 sumGrad0_->axpy(weight*pvalp1,*pointGrad_);
270 Real gvar = 1.0; gradient0_->zero();
272 GradientSampler_->sumAll(&myval[0],&val[0],2);
274 GradientSampler_->sumAll(*sumGrad0_,*gradient0_);
276 Real norm = ((order_ == 1.0) ? 1.0
277 : ((order_ == 2.0) ? std::sqrt(val[0])
278 : std::pow(val[0],(order_-1.0)/order_)));
279 gvar -= val[1]/((1.0-
prob_)*norm);
280 gradient0_->scale(1.0/((1.0-prob_)*norm));
289 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
291 Teuchos::RCP<Vector<Real> > vvec; Real vvar = 0.0;
297 sumGrad0_->zero(); sumGrad1_->zero(); sumGrad2_->zero(); sumHess_->zero();
298 gradient0_->zero(); gradient1_->zero(); gradient2_->zero();
300 std::vector<Real> point;
301 Real pval = 0.0, pvalp0 = 0.0, pvalp1 = 0.0, pvalp2 = 0.0, gv = 0.0;
302 std::vector<Real> val(5,0.0), myval(5,0.0);
303 int start = HessianSampler_->start(), end = HessianSampler_->numMySamples();
304 for (
int i = start; i < end; i++ ) {
306 weight = HessianSampler_->getMyWeight(i);
307 point = HessianSampler_->getMyPoint(i);
312 pvalp0 = ((order_ == 1.0) ? pval-xvar
313 : std::pow(pval-xvar,order_));
314 pvalp1 = ((order_ == 1.0) ? 1.0
315 : ((order_ == 2.0) ? pval-xvar
316 : std::pow(pval-xvar,order_-1.0)));
317 pvalp2 = ((order_ == 1.0) ? 0.0
318 : ((order_ == 2.0) ? 1.0
319 : ((order_ == 3.0) ? pval-xvar
320 : std::pow(pval-xvar,order_-2.0))));
322 myval[0] += weight*pvalp0;
323 myval[1] += weight*pvalp1;
324 myval[2] += weight*pvalp2;
327 gv = pointGrad_->dot(vvec->dual());
329 myval[3] += weight*pvalp1*gv;
330 myval[4] += weight*pvalp2*gv;
331 sumGrad0_->axpy(weight*pvalp1,*pointGrad_);
332 sumGrad1_->axpy(weight*pvalp2,*pointGrad_);
333 sumGrad2_->axpy(weight*pvalp2*gv,*pointGrad_);
335 getHessVec(*pointHess_,*vvec,*xvec,point,tol);
337 sumHess_->axpy(weight*pvalp1,*pointHess_);
340 Real hvar = 0.0; hessvec_->zero();
341 HessianSampler_->sumAll(&myval[0],&val[0],5);
344 HessianSampler_->sumAll(*sumGrad0_,*gradient0_);
345 HessianSampler_->sumAll(*sumGrad1_,*gradient1_);
346 HessianSampler_->sumAll(*sumGrad2_,*gradient2_);
347 HessianSampler_->sumAll(*sumHess_,*hessvec_);
349 Real norm0 = (1.0-
prob_)*((order_ == 1.0) ? 1.0
350 : ((order_ == 2.0) ? std::sqrt(val[0])
351 : std::pow(val[0],(order_-1.0)/order_)));
352 Real norm1 = (1.0-
prob_)*((order_ == 1.0) ? 1.0
353 : std::pow(val[0],(2.0*order_-1.0)/order_));
354 hvar = (order_-1.0)*((val[2]/norm0 - val[1]*val[1]/norm1)*vvar
355 -(val[4]/norm0 - val[3]*val[1]/norm1));
356 hessvec_->scale(1.0/norm0);
357 hessvec_->axpy(-(order_-1.0)*vvar/norm0,*gradient1_);
358 hessvec_->axpy((order_-1.0)*(vvar*val[1]-val[3])/norm1,*gradient0_);
359 hessvec_->axpy((order_-1.0)/norm0,*gradient2_);
void getValue(Real &val, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
Provides the interface to evaluate objective functions.
Teuchos::RCP< SampleGenerator< Real > > GradientSampler_
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Teuchos::RCP< ParametrizedObjective< Real > > ParametrizedObjective_
Teuchos::RCP< SampleGenerator< Real > > ValueSampler_
void getHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
Teuchos::RCP< Vector< Real > > sumGrad1_
Teuchos::RCP< Vector< Real > > gradient1_
Teuchos::RCP< Vector< Real > > hessvec_
Teuchos::RCP< Vector< Real > > sumGrad2_
Teuchos::RCP< Vector< Real > > gradient0_
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply preconditioner to vector.
const Real getStatistic(const int i=0) const
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
Teuchos::RCP< Vector< Real > > pointGrad_
std::map< std::vector< Real >, Real > value_storage_
Teuchos::RCP< SampleGenerator< Real > > HessianSampler_
std::map< std::vector< Real >, Teuchos::RCP< Vector< Real > > > gradient_storage_
void setVector(const Vector< Real > &vec)
Teuchos::RCP< const Vector< Real > > getVector() const
void initialize(const Vector< Real > &x)
void setStatistic(const Real stat)
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Teuchos::RCP< Vector< Real > > sumHess_
HMCRObjective(Teuchos::RCP< ParametrizedObjective< Real > > &pObj, Real order, Real prob, Teuchos::RCP< SampleGenerator< Real > > &sampler, bool storage=true)
Teuchos::RCP< Vector< Real > > pointHess_
virtual void set(const Vector &x)
Set where .
HMCRObjective(Teuchos::RCP< ParametrizedObjective< Real > > &pObj, Real order, Real prob, Teuchos::RCP< SampleGenerator< Real > > &vsampler, Teuchos::RCP< SampleGenerator< Real > > &gsampler, bool storage=true)
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Teuchos::RCP< Vector< Real > > gradient2_
void getGradient(Vector< Real > &g, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
void unwrap_const_CVaR_vector(Teuchos::RCP< Vector< Real > > &xvec, Real &xvar, const Vector< Real > &x)
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
Teuchos::RCP< Vector< Real > > sumGrad0_
static const double ROL_EPSILON
Platform-dependent machine epsilon.
HMCRObjective(Teuchos::RCP< ParametrizedObjective< Real > > &pObj, Real order, Real prob, Teuchos::RCP< SampleGenerator< Real > > &vsampler, Teuchos::RCP< SampleGenerator< Real > > &gsampler, Teuchos::RCP< SampleGenerator< Real > > &hsampler, bool storage=true)