44 #ifndef ROL_BPOEOBJECTIVE_HPP 45 #define ROL_BPOEOBJECTIVE_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);
151 const Real order,
const Real threshold,
155 const bool storage =
true )
156 : ParametrizedObjective_(pObj), order_(order), threshold_(threshold),
157 ValueSampler_(vsampler), GradientSampler_(gsampler), HessianSampler_(hsampler),
158 storage_(storage), initialized_(false) {
159 value_storage_.clear();
160 gradient_storage_.clear();
164 const Real order,
const Real threshold,
167 const bool storage =
true )
168 : ParametrizedObjective_(pObj), order_(order), threshold_(threshold),
169 ValueSampler_(vsampler), GradientSampler_(gsampler), HessianSampler_(gsampler),
170 storage_(storage), initialized_(false) {
171 value_storage_.clear();
172 gradient_storage_.clear();
176 const Real order,
const Real threshold,
178 const bool storage =
true )
179 : ParametrizedObjective_(pObj), order_(order), threshold_(threshold),
180 ValueSampler_(sampler), GradientSampler_(sampler), HessianSampler_(sampler),
181 storage_(storage), initialized_(false) {
182 value_storage_.clear();
183 gradient_storage_.clear();
187 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
189 ParametrizedObjective_->update(*xvec,flag,iter);
190 ValueSampler_->update(*xvec);
192 value_storage_.clear();
195 GradientSampler_->update(*xvec);
196 HessianSampler_->update(*xvec);
198 gradient_storage_.clear();
204 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
207 std::vector<Real> point;
208 Real weight = 0.0, myval = 0.0, pval = 0.0, val = 0.0, bp = 0.0;
209 int start = ValueSampler_->start(), end = ValueSampler_->numMySamples();
210 for (
int i = start; i < end; i++ ) {
211 weight = ValueSampler_->getMyWeight(i);
212 point = ValueSampler_->getMyPoint(i);
218 myval += weight*((order_==1.0) ? bp
219 : std::pow(bp,order_));
223 ValueSampler_->sumAll(&myval,&val,1);
225 return ((order_==1.0) ? val : std::pow(val,1.0/order_));
229 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
234 g.
zero(); sumGrad0_->zero(); pointGrad_->zero();
235 std::vector<Real> point, val(2,0.0), myval(2,0.0);
236 Real weight = 0.0, pval = 0.0, pvalp0 = 0.0, pvalp1 = 0.0, bp = 0.0;
237 int start = GradientSampler_->start(), end = GradientSampler_->numMySamples();
238 for (
int i = start; i < end; i++ ) {
239 weight = GradientSampler_->getMyWeight(i);
240 point = GradientSampler_->getMyPoint(i);
246 pvalp0 = ((order_==1.0) ? bp
247 : std::pow(bp,order_));
248 pvalp1 = ((order_==1.0) ? 1.0
249 : ((order_==2.0) ? bp
250 : std::pow(bp,order_-1.0)));
252 myval[0] += weight*pvalp0;
257 sumGrad0_->axpy(weight*pvalp1,*pointGrad_);
261 GradientSampler_->sumAll(&myval[0],&val[0],2);
263 Real gvar = 0.0; gradient0_->zero();
265 GradientSampler_->sumAll(*sumGrad0_,*gradient0_);
266 Real norm = std::pow(val[0],(order_-1.0)/order_);
267 gradient0_->scale(xvar/norm);
277 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
279 Teuchos::RCP<Vector<Real> > vvec; Real vvar = 0.0;
284 hv.
zero(); sumHess_->zero(); hessvec_->zero();
285 sumGrad0_->zero(); sumGrad1_->zero(); sumGrad2_->zero();
286 gradient0_->zero(); gradient1_->zero(); gradient2_->zero();
287 std::vector<Real> point, val(5,0.0), myval(5,0.0);
288 Real pval = 0.0, pvalp0 = 0.0, pvalp1 = 0.0, pvalp2 = 0.0;
289 Real weight = 0.0, gv = 0.0, bp = 0.0;
290 int start = HessianSampler_->start(), end = HessianSampler_->numMySamples();
291 for (
int i = start; i < end; i++ ) {
293 weight = HessianSampler_->getMyWeight(i);
294 point = HessianSampler_->getMyPoint(i);
300 pvalp0 = ((order_==1.0) ? bp
301 : std::pow(bp,order_));
302 pvalp1 = ((order_==1.0) ? 1.0
303 : ((order_==2.0) ? bp
304 : std::pow(bp,order_-1.0)));
305 pvalp2 = ((order_==1.0) ? 0.0
306 : ((order_==2.0) ? 1.0
307 : ((order_==3.0) ? bp
308 : std::pow(bp,order_-2.0))));
310 myval[0] += weight*pvalp0;
312 myval[2] += weight*pvalp2*(pval-
threshold_)*(pval-threshold_);
315 gv = pointGrad_->dot(vvec->dual());
317 myval[3] += weight*pvalp1*gv;
318 myval[4] += weight*pvalp2*(pval-
threshold_)*gv;
319 sumGrad0_->axpy(weight*pvalp1,*pointGrad_);
320 sumGrad1_->axpy(weight*pvalp2*(pval-threshold_),*pointGrad_);
321 sumGrad2_->axpy(weight*pvalp2*gv,*pointGrad_);
323 getHessVec(*pointHess_,*vvec,*xvec,point,tol);
325 sumHess_->axpy(weight*pvalp1,*pointHess_);
329 HessianSampler_->sumAll(&myval[0],&val[0],5);
330 Real hvar = 0.0; hessvec_->zero();
332 HessianSampler_->sumAll(*sumGrad0_,*gradient0_);
333 HessianSampler_->sumAll(*sumGrad1_,*gradient1_);
334 HessianSampler_->sumAll(*sumGrad2_,*gradient2_);
335 HessianSampler_->sumAll(*sumHess_,*hessvec_);
337 Real norm0 = ((order_==1.0) ? 1.0
338 : ((order_==2.0) ? std::sqrt(val[0])
339 : std::pow(val[0],(order_-1.0)/order_)));
340 Real norm1 = ((order_==1.0) ? val[0]
341 : std::pow(val[0],(2.0*order_-1.0)/order_));
342 hvar = (order_-1.0)*((val[2]/norm0 - val[1]*val[1]/norm1)*vvar
343 +xvar*(val[4]/norm0 - val[3]*val[1]/norm1))
345 hessvec_->scale(xvar/norm0);
346 Real coeff = -(order_-1.0)*xvar*(xvar*val[3]+vvar*val[1])/norm1+vvar/norm0;
347 hessvec_->axpy(coeff,*gradient0_);
348 hessvec_->axpy((order_-1.0)*vvar*xvar/norm0,*gradient1_);
349 hessvec_->axpy((order_-1.0)*xvar*xvar/norm0,*gradient2_);
Provides the interface to evaluate objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Teuchos::RCP< Vector< Real > > gradient1_
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)
Teuchos::RCP< SampleGenerator< Real > > HessianSampler_
Teuchos::RCP< Vector< Real > > pointGrad_
Teuchos::RCP< Vector< Real > > sumGrad1_
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
void initialize(const Vector< Real > &x)
void getGradient(Vector< Real > &g, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
Teuchos::RCP< Vector< Real > > sumGrad0_
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.
void getValue(Real &val, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
Teuchos::RCP< Vector< Real > > pointHess_
std::map< std::vector< Real >, Real > value_storage_
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
std::map< std::vector< Real >, Teuchos::RCP< Vector< Real > > > gradient_storage_
void setVector(const Vector< Real > &vec)
Teuchos::RCP< const Vector< Real > > getVector() const
Teuchos::RCP< SampleGenerator< Real > > ValueSampler_
BPOEObjective(const Teuchos::RCP< ParametrizedObjective< Real > > &pObj, const Real order, const Real threshold, const Teuchos::RCP< SampleGenerator< Real > > &vsampler, const Teuchos::RCP< SampleGenerator< Real > > &gsampler, const Teuchos::RCP< SampleGenerator< Real > > &hsampler, const bool storage=true)
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.
BPOEObjective(const Teuchos::RCP< ParametrizedObjective< Real > > &pObj, const Real order, const Real threshold, const Teuchos::RCP< SampleGenerator< Real > > &vsampler, const Teuchos::RCP< SampleGenerator< Real > > &gsampler, const bool storage=true)
Teuchos::RCP< Vector< Real > > gradient2_
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Teuchos::RCP< Vector< Real > > sumHess_
Teuchos::RCP< Vector< Real > > sumGrad2_
Teuchos::RCP< Vector< Real > > gradient0_
BPOEObjective(const Teuchos::RCP< ParametrizedObjective< Real > > &pObj, const Real order, const Real threshold, const Teuchos::RCP< SampleGenerator< Real > > &sampler, const bool storage=true)
virtual void set(const Vector &x)
Set where .
void getHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
Teuchos::RCP< SampleGenerator< Real > > GradientSampler_
Teuchos::RCP< Vector< Real > > hessvec_
Teuchos::RCP< ParametrizedObjective< Real > > ParametrizedObjective_
static const double ROL_EPSILON
Platform-dependent machine epsilon.