44#ifndef ROL_SPECTRALRISK_HPP
45#define ROL_SPECTRALRISK_HPP
91 ROL::Ptr<MixedCVaR<Real> >
mqq_;
98 ROL_TEST_FOR_EXCEPTION(
plusFunction_ == ROL::nullPtr, std::invalid_argument,
99 ">>> ERROR (ROL::SpectralRisk): PlusFunction pointer is null!");
100 if ( dist != ROL::nullPtr) {
101 Real lb = dist->lowerBound();
102 Real ub = dist->upperBound();
103 ROL_TEST_FOR_EXCEPTION(lb <
static_cast<Real
>(0), std::invalid_argument,
104 ">>> ERROR (ROL::SpectralRisk): Distribution lower bound less than zero!");
105 ROL_TEST_FOR_EXCEPTION(ub >
static_cast<Real
>(1), std::invalid_argument,
106 ">>> ERROR (ROL::SpectralRisk): Distribution upper bound greater than one!");
113 pts_.clear();
pts_.assign(pts.begin(),pts.end());
114 wts_.clear();
wts_.assign(wts.begin(),wts.end());
116 mqq_ = ROL::makePtr<MixedCVaR<Real>>(pts,wts,pf);
121 const Real lo = dist->lowerBound(), hi = dist->upperBound();
122 const Real half(0.5), one(1), N(nQuad);
123 wts.clear(); wts.resize(nQuad);
124 pts.clear(); pts.resize(nQuad);
126 wts[0] = half/(N-half);
128 for (
int i = 1; i < nQuad; ++i) {
129 wts[i] = one/(N-half);
130 pts[i] = dist->invertCDF(
static_cast<Real
>(i)/N);
134 wts[0] = half/(N-one);
136 for (
int i = 1; i < nQuad-1; ++i) {
137 wts[i] = one/(N-one);
138 pts[i] = dist->invertCDF(
static_cast<Real
>(i)/N);
140 wts[nQuad-1] = half/(N-one);
146 const std::vector<Real> &wts,
147 const bool print =
false)
const {
149 const int nQuad = wts.size();
150 std::cout << std::endl;
151 std::cout << std::scientific << std::setprecision(15);
152 std::cout << std::setw(25) << std::left <<
"Points"
153 << std::setw(25) << std::left <<
"Weights"
155 for (
int i = 0; i < nQuad; ++i) {
156 std::cout << std::setw(25) << std::left << pts[i]
157 << std::setw(25) << std::left << wts[i]
160 std::cout << std::endl;
173 std::vector<Real> wts(nQuad), pts(nQuad);
184 ROL::ParameterList &list
185 = parlist.sublist(
"SOL").sublist(
"Risk Measure").sublist(
"Spectral Risk");
186 int nQuad = list.get(
"Number of Quadrature Points",5);
187 bool print = list.get(
"Print Quadrature to Screen",
false);
189 ROL::Ptr<Distribution<Real> > dist = DistributionFactory<Real>(list);
191 ROL::Ptr<PlusFunction<Real> > pf = ROL::makePtr<PlusFunction<Real>>(list);
193 std::vector<Real> wts(nQuad), pts(nQuad);
202 SpectralRisk(
const std::vector<Real> &pts,
const std::vector<Real> &wts,
212 int nQuad =
static_cast<int>(
wts_.size());
213 for (
int i = 0; i < nQuad; ++i) {
214 stat +=
wts_[i] * xstat[i];
222 mqq_->setStorage(value_storage,gradient_storage);
228 mqq_->setHessVecStorage(gradvec_storage,hessvec_storage);
231 void setSample(
const std::vector<Real> &point,
const Real weight) {
233 mqq_->setSample(point,weight);
238 mqq_->resetStorage(flag);
243 mqq_->resetStorage(type);
252 return mqq_->computeStatistic(xstat);
257 const std::vector<Real> &xstat,
259 mqq_->updateValue(obj,x,xstat,tol);
264 const std::vector<Real> &xstat,
266 mqq_->updateGradient(obj,x,xstat,tol);
271 const std::vector<Real> &vstat,
273 const std::vector<Real> &xstat,
275 mqq_->updateHessVec(obj,v,vstat,x,xstat,tol);
279 const std::vector<Real> &xstat,
281 return mqq_->getValue(x,xstat,sampler);
285 std::vector<Real> &gstat,
287 const std::vector<Real> &xstat,
289 mqq_->getGradient(g,gstat,x,xstat,sampler);
293 std::vector<Real> &hvstat,
295 const std::vector<Real> &vstat,
297 const std::vector<Real> &xstat,
299 mqq_->getHessVec(hv,hvstat,v,vstat,x,xstat,sampler);
Provides the interface to evaluate objective functions.
Provides the interface to implement any functional that maps a random variable to a (extended) real n...
virtual void setSample(const std::vector< Real > &point, const Real weight)
virtual void initialize(const Vector< Real > &x)
Initialize temporary variables.
virtual void setStorage(const Ptr< ScalarController< Real > > &value_storage, const Ptr< VectorController< Real > > &gradient_storage)
virtual void setHessVecStorage(const Ptr< ScalarController< Real > > &gradvec_storage, const Ptr< VectorController< Real > > &hessvec_storage)
virtual void resetStorage(bool flag=true)
Reset internal storage.
Provides an interface for spectral risk measures.
Real computeStatistic(const Ptr< std::vector< Real > > &xstat) const
SpectralRisk(const ROL::Ptr< Distribution< Real > > &dist, const int nQuad, const ROL::Ptr< PlusFunction< Real > > &pf)
void updateHessVec(Objective< Real > &obj, const Vector< Real > &v, const std::vector< Real > &vstat, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal risk measure storage for Hessian-time-a-vector computation.
void setStorage(const Ptr< ScalarController< Real > > &value_storage, const Ptr< VectorController< Real > > &gradient_storage)
void resetStorage(bool flag=true)
Reset internal storage.
void printQuad(const std::vector< Real > &pts, const std::vector< Real > &wts, const bool print=false) const
void getGradient(Vector< Real > &g, std::vector< Real > &gstat, const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure (sub)gradient.
void buildQuadFromDist(std::vector< Real > &pts, std::vector< Real > &wts, const int nQuad, const ROL::Ptr< Distribution< Real > > &dist) const
void resetStorage(UpdateType type)
Real computeStatistic(const std::vector< Real > &xstat)
void updateGradient(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal risk measure storage for gradient computation.
void setHessVecStorage(const Ptr< ScalarController< Real > > &gradvec_storage, const Ptr< VectorController< Real > > &hessvec_storage)
ROL::Ptr< PlusFunction< Real > > plusFunction_
void buildMixedQuantile(const std::vector< Real > &pts, const std::vector< Real > &wts, const ROL::Ptr< PlusFunction< Real > > &pf)
void getHessVec(Vector< Real > &hv, std::vector< Real > &hvstat, const Vector< Real > &v, const std::vector< Real > &vstat, const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure Hessian-times-a-vector.
void updateValue(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal storage for value computation.
ROL::Ptr< MixedCVaR< Real > > mqq_
void checkInputs(ROL::Ptr< Distribution< Real > > &dist=ROL::nullPtr) const
SpectralRisk(ROL::ParameterList &parlist)
void initialize(const Vector< Real > &x)
Initialize temporary variables.
void setSample(const std::vector< Real > &point, const Real weight)
SpectralRisk(const std::vector< Real > &pts, const std::vector< Real > &wts, const ROL::Ptr< PlusFunction< Real > > &pf)
Real getValue(const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure value.
Defines the linear algebra or vector space interface.