44 #ifndef ROL_RISKAVERSEOBJECTIVE_HPP
45 #define ROL_RISKAVERSEOBJECTIVE_HPP
47 #include "Teuchos_RefCountPtr.hpp"
60 Teuchos::RCP<ParametrizedObjective<Real> >
pObj_;
63 Teuchos::RCP<RiskMeasure<Real> >
rm_;
75 bool storage =
true ) :
storage_(storage) {
76 pObj_ = Teuchos::rcp(&pObj,
false);
77 vsampler_ = Teuchos::rcp(&vsampler,
false);
78 gsampler_ = Teuchos::rcp(&gsampler,
false);
79 rm_ = Teuchos::rcp(&rm,
false);
93 pObj_ = Teuchos::rcp(&pObj,
false);
96 rm_ = Teuchos::rcp(&rm,
false);
103 pObj_->update(x,flag,iter);
118 Teuchos::RCP<Vector<Real> > x0;
120 for (
int i = 0; i <
vsampler_->numMySamples(); i++ ) {
126 val =
pObj_->value(*x0,tol);
139 Teuchos::RCP<Vector<Real> > x0;
141 Teuchos::RCP<Vector<Real> > g0 = x0->clone();
143 for (
int i = 0; i <
gsampler_->numMySamples(); i++ ) {
149 val =
pObj_->value(*x0,tol);
158 pObj_->gradient(*g0,*x0,tol);
160 Teuchos::RCP<Vector<Real> > tmp = g0->clone();
173 Teuchos::RCP<Vector<Real> > x0;
174 Teuchos::RCP<Vector<Real> > v0;
175 rm_->reset(x0,x,v0,v);
178 Teuchos::RCP<Vector<Real> > g0 = x0->clone();
179 Teuchos::RCP<Vector<Real> > h0 = x0->clone();
180 for (
int i = 0; i <
gsampler_->numMySamples(); i++ ) {
186 val =
pObj_->value(*x0,tol);
195 pObj_->gradient(*g0,*x0,tol);
197 Teuchos::RCP<Vector<Real> > tmp = g0->clone();
203 pObj_->hessVec(*h0,*v0,*x0,tol);
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
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< SampleGenerator< Real > > vsampler_
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
Teuchos::RCP< ParametrizedObjective< Real > > pObj_
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply preconditioner to vector.
RiskAverseObjective(ParametrizedObjective< Real > &pObj, RiskMeasure< Real > &rm, SampleGenerator< Real > &sampler, bool storage=true)
RiskAverseObjective(ParametrizedObjective< Real > &pObj, RiskMeasure< Real > &rm, SampleGenerator< Real > &vsampler, SampleGenerator< Real > &gsampler, bool storage=true)
std::map< std::vector< Real >, Teuchos::RCP< Vector< Real > > > gradient_storage_
std::map< std::vector< Real >, Real > value_storage_
Teuchos::RCP< SampleGenerator< Real > > gsampler_
virtual void set(const Vector &x)
Set where .
Teuchos::RCP< RiskMeasure< Real > > rm_
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
virtual ~RiskAverseObjective()
virtual Real value(const Vector< Real > &x, Real &tol)
Compute value.