44 #ifndef ROL_STOCHASTICOBJECTIVE_HPP
45 #define ROL_STOCHASTICOBJECTIVE_HPP
60 Ptr<RandVarFunctional<Real>>
rvf_;
83 if (xstat == nullPtr) {
84 xstat = makePtr<const std::vector<Real>>(0);
92 if (xstat == nullPtr) {
93 xstat = makePtr<std::vector<Real>>(0);
106 const bool storage =
true,
107 const int comp = 0,
const int index = 0 )
111 rvf->useStorage(storage);
118 const bool storage =
true,
119 const int comp = 0,
const int index = 0 )
125 const bool storage =
true,
126 const int comp = 0,
const int index = 0 )
130 ParameterList &parlist,
134 const int comp = 0,
const int index = 0 )
138 std::string name, type = parlist.sublist(
"SOL").get(
"Type",
"Risk Averse");
139 if (type ==
"Risk Averse")
140 name = parlist.sublist(
"SOL").sublist(
"Risk Measure").get(
"Name",
"CVaR");
142 if (type ==
"Risk Averse" && name ==
"Convex Combination Risk Measure")
143 rvf_ = makePtr<ConvexCombinationRiskMeasure<Real>>(parlist);
145 rvf_ = RandVarFunctionalFactory<Real>(parlist);
147 bool storage = parlist.sublist(
"SOL").get(
"Store Sampled Value and Gradient",
true);
148 rvf_->useStorage(storage);
152 ROL::ParameterList &parlist,
155 const int comp = 0,
const int index = 0 )
159 ROL::ParameterList &parlist,
161 const int comp = 0,
const int index = 0 )
166 return rvf_->computeStatistic(xstat);
176 rvf_->resetStorage(type);
178 obj_->update(*x0,type,iter);
190 rvf_->resetStorage(flag);
192 obj_->update(*x0,flag,iter);
204 rvf_->initialize(*x0);
206 for (
int i = 0; i <
vsampler_->numMySamples(); i++ ) {
208 rvf_->updateValue(*
obj_,*x0,*xstat,tol);
219 Ptr<std::vector<Real>> gstat =
getStat(g);
220 rvf_->initialize(*x0);
221 for (
int i = 0; i <
gsampler_->numMySamples(); i++ ) {
223 rvf_->updateGradient(*
obj_,*x0,*xstat,tol);
236 Ptr<std::vector<Real>> hvstat =
getStat(hv);
237 rvf_->initialize(*x0);
238 for (
int i = 0; i <
hsampler_->numMySamples(); i++ ) {
240 rvf_->updateHessVec(*
obj_,*v0,*vstat,*x0,*xstat,tol);
242 rvf_->getHessVec(*hv0,*hvstat,*v0,*vstat,*x0,*xstat,*
hsampler_);
Provides the interface to evaluate objective functions.
Real computeStatistic(const Vector< Real > &x) const
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
StochasticObjective(const Ptr< Objective< Real >> &obj, const Ptr< RandVarFunctional< Real >> &rvf, const Ptr< SampleGenerator< Real >> &vsampler, const Ptr< SampleGenerator< Real >> &gsampler, const bool storage=true, const int comp=0, const int index=0)
Ptr< Objective< Real > > obj_
virtual void zero()
Set to zero vector.
virtual ~StochasticObjective()
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.
Ptr< std::vector< Real > > getStat(Vector< Real > &x) const
Ptr< const Vector< Real > > getConstVector(const Vector< Real > &x) const
StochasticObjective(const Ptr< Objective< Real >> &obj, ROL::ParameterList &parlist, const Ptr< SampleGenerator< Real >> &sampler, const int comp=0, const int index=0)
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
StochasticObjective(const Ptr< Objective< Real >> &obj, const Ptr< RandVarFunctional< Real >> &rvf, const Ptr< SampleGenerator< Real >> &sampler, const bool storage=true, const int comp=0, const int index=0)
Ptr< std::vector< Real > > getStatistic(const int comp=0, const int index=0)
Ptr< SampleGenerator< Real > > hsampler_
Ptr< SampleGenerator< Real > > vsampler_
Ptr< RandVarFunctional< Real > > rvf_
StochasticObjective(const Ptr< Objective< Real >> &obj, const Ptr< RandVarFunctional< Real >> &rvf, const Ptr< SampleGenerator< Real >> &vsampler, const Ptr< SampleGenerator< Real >> &gsampler, const Ptr< SampleGenerator< Real >> &hsampler, const bool storage=true, const int comp=0, const int index=0)
Ptr< const std::vector< Real > > getConstStat(const Vector< Real > &x) const
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Ptr< const Vector< Real > > getVector(void) const
Provides the interface to implement any functional that maps a random variable to a (extended) real n...
Ptr< Vector< Real > > getVector(Vector< Real > &x) const
virtual void set(const Vector &x)
Set where .
StochasticObjective(const Ptr< Objective< Real >> &obj, ROL::ParameterList &parlist, const Ptr< SampleGenerator< Real >> &vsampler, const Ptr< SampleGenerator< Real >> &gsampler, const int comp=0, const int index=0)
void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
StochasticObjective(const Ptr< Objective< Real >> &obj, ParameterList &parlist, const Ptr< SampleGenerator< Real >> &vsampler, const Ptr< SampleGenerator< Real >> &gsampler, const Ptr< SampleGenerator< Real >> &hsampler, const int comp=0, const int index=0)
Ptr< SampleGenerator< Real > > gsampler_
Real value(const Vector< Real > &x, Real &tol)
Compute value.