44 #ifndef ROL_HMCROBJECTIVE_HPP 45 #define ROL_HMCROBJECTIVE_HPP 47 #include "Teuchos_RCP.hpp" 108 const std::vector<Real> ¶m, Real &tol) {
116 value_storage_.insert(std::pair<std::vector<Real>,Real>(param,val));
122 const std::vector<Real> ¶m, Real &tol) {
130 Teuchos::RCP<Vector<Real> > tmp = g.
clone();
138 const std::vector<Real> ¶m, Real &tol) {
148 Real order, Real prob,
152 bool storage =
true )
156 order_ = ((order < 1.0) ? 1.0 : order);
157 prob_ = ((prob > 1.0) ? 1.0 : ((prob < 0.0) ? 0.0 : prob));
163 Real order, Real prob,
166 bool storage =
true )
170 order_ = ((order < 1.0) ? 1.0 : order);
171 prob_ = ((prob > 1.0) ? 1.0 : ((prob < 0.0) ? 0.0 : prob));
177 Real order, Real prob,
179 bool storage =
true )
183 order_ = ((order < 1.0) ? 1.0 : order);
184 prob_ = ((prob > 1.0) ? 1.0 : ((prob < 0.0) ? 0.0 : prob));
190 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
207 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
210 std::vector<Real> point;
211 Real weight = 0.0, myval = 0.0, pval = 0.0, val = 0.0;
213 for (
int i = start; i < end; i++ ) {
220 myval += weight*((
order_ == 1.0) ? pval-xvar
221 : std::pow(pval-xvar,
order_));
227 if (std::abs(val) < ROL_EPSILON<Real>()) {
230 return xvar + ((
order_ == 1.0) ? val
231 : ((
order_ == 2.0) ? std::sqrt(val)
236 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
241 std::vector<Real> point, val(2,0.0), myval(2,0.0);
242 Real weight = 0.0, pval = 0.0, pvalp0 = 0.0, pvalp1 = 0.0;
244 for (
int i = start; i < end; i++ ) {
251 pvalp0 = ((
order_ == 1.0) ? pval-xvar
252 : std::pow(pval-xvar,
order_));
253 pvalp1 = ((
order_ == 1.0) ? 1.0
254 : ((
order_ == 2.0) ? pval-xvar
255 : std::pow(pval-xvar,
order_-1.0)));
257 myval[0] += weight*pvalp0;
258 myval[1] += weight*pvalp1;
268 if (std::abs(val[0]) >= ROL_EPSILON<Real>()) {
271 Real norm = ((
order_ == 1.0) ? 1.0
272 : ((
order_ == 2.0) ? std::sqrt(val[0])
274 gvar -= val[1]/((1.0-
prob_)*norm);
284 Teuchos::RCP<Vector<Real> > xvec; Real xvar = 0.0;
286 Teuchos::RCP<Vector<Real> > vvec; Real vvar = 0.0;
294 std::vector<Real> point;
295 Real pval = 0.0, pvalp0 = 0.0, pvalp1 = 0.0, pvalp2 = 0.0, gv = 0.0;
296 std::vector<Real> val(5,0.0), myval(5,0.0);
298 for (
int i = start; i < end; i++ ) {
306 pvalp0 = ((
order_ == 1.0) ? pval-xvar
307 : std::pow(pval-xvar,
order_));
308 pvalp1 = ((
order_ == 1.0) ? 1.0
309 : ((
order_ == 2.0) ? pval-xvar
310 : std::pow(pval-xvar,
order_-1.0)));
311 pvalp2 = ((
order_ == 1.0) ? 0.0
313 : ((
order_ == 3.0) ? pval-xvar
314 : std::pow(pval-xvar,
order_-2.0))));
316 myval[0] += weight*pvalp0;
317 myval[1] += weight*pvalp1;
318 myval[2] += weight*pvalp2;
323 myval[3] += weight*pvalp1*gv;
324 myval[4] += weight*pvalp2*gv;
336 if (std::abs(val[0]) >= ROL_EPSILON<Real>()) {
344 : ((
order_ == 2.0) ? std::sqrt(val[0])
348 hvar = (
order_-1.0)*((val[2]/norm0 - val[1]*val[1]/norm1)*vvar
349 -(val[4]/norm0 - val[3]*val[1]/norm1));
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_
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.
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)
void initialize(const Vector< Real > &x)
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
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.
HMCRObjective(Teuchos::RCP< Objective< Real > > &pObj, Real order, Real prob, Teuchos::RCP< SampleGenerator< Real > > &sampler, bool storage=true)
Teuchos::RCP< Vector< Real > > sumHess_
HMCRObjective(Teuchos::RCP< Objective< Real > > &pObj, Real order, Real prob, Teuchos::RCP< SampleGenerator< Real > > &vsampler, Teuchos::RCP< SampleGenerator< Real > > &gsampler, Teuchos::RCP< SampleGenerator< Real > > &hsampler, bool storage=true)
Teuchos::RCP< Objective< Real > > ParametrizedObjective_
Teuchos::RCP< Vector< Real > > pointHess_
virtual void set(const Vector &x)
Set where .
HMCRObjective(Teuchos::RCP< Objective< 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_