10 #ifndef ROL_SPECTRALRISK_HPP
11 #define ROL_SPECTRALRISK_HPP
57 ROL::Ptr<MixedCVaR<Real> >
mqq_;
64 ROL_TEST_FOR_EXCEPTION(
plusFunction_ == ROL::nullPtr, std::invalid_argument,
65 ">>> ERROR (ROL::SpectralRisk): PlusFunction pointer is null!");
66 if ( dist != ROL::nullPtr) {
67 Real lb = dist->lowerBound();
68 Real ub = dist->upperBound();
69 ROL_TEST_FOR_EXCEPTION(lb < static_cast<Real>(0), std::invalid_argument,
70 ">>> ERROR (ROL::SpectralRisk): Distribution lower bound less than zero!");
71 ROL_TEST_FOR_EXCEPTION(ub > static_cast<Real>(1), std::invalid_argument,
72 ">>> ERROR (ROL::SpectralRisk): Distribution upper bound greater than one!");
79 pts_.clear();
pts_.assign(pts.begin(),pts.end());
80 wts_.clear();
wts_.assign(wts.begin(),wts.end());
82 mqq_ = ROL::makePtr<MixedCVaR<Real>>(pts,wts,pf);
87 const Real lo = dist->lowerBound(), hi = dist->upperBound();
88 const Real half(0.5), one(1), N(nQuad);
89 wts.clear(); wts.resize(nQuad);
90 pts.clear(); pts.resize(nQuad);
92 wts[0] = half/(N-half);
94 for (
int i = 1; i < nQuad; ++i) {
95 wts[i] = one/(N-half);
96 pts[i] = dist->invertCDF(static_cast<Real>(i)/N);
100 wts[0] = half/(N-one);
102 for (
int i = 1; i < nQuad-1; ++i) {
103 wts[i] = one/(N-one);
104 pts[i] = dist->invertCDF(static_cast<Real>(i)/N);
106 wts[nQuad-1] = half/(N-one);
112 const std::vector<Real> &wts,
113 const bool print =
false)
const {
115 const int nQuad = wts.size();
116 std::cout << std::endl;
117 std::cout << std::scientific << std::setprecision(15);
118 std::cout << std::setw(25) << std::left <<
"Points"
119 << std::setw(25) << std::left <<
"Weights"
121 for (
int i = 0; i < nQuad; ++i) {
122 std::cout << std::setw(25) << std::left << pts[i]
123 << std::setw(25) << std::left << wts[i]
126 std::cout << std::endl;
139 std::vector<Real> wts(nQuad), pts(nQuad);
150 ROL::ParameterList &list
151 = parlist.sublist(
"SOL").sublist(
"Risk Measure").sublist(
"Spectral Risk");
152 int nQuad = list.get(
"Number of Quadrature Points",5);
153 bool print = list.get(
"Print Quadrature to Screen",
false);
155 ROL::Ptr<Distribution<Real> > dist = DistributionFactory<Real>(list);
157 ROL::Ptr<PlusFunction<Real> > pf = ROL::makePtr<PlusFunction<Real>>(list);
159 std::vector<Real> wts(nQuad), pts(nQuad);
168 SpectralRisk(
const std::vector<Real> &pts,
const std::vector<Real> &wts,
179 mqq_->setStorage(value_storage,gradient_storage);
185 mqq_->setHessVecStorage(gradvec_storage,hessvec_storage);
188 void setSample(
const std::vector<Real> &point,
const Real weight) {
190 mqq_->setSample(point,weight);
195 mqq_->resetStorage(flag);
200 mqq_->resetStorage(type);
209 return mqq_->computeStatistic(xstat);
214 const std::vector<Real> &xstat,
216 mqq_->updateValue(obj,x,xstat,tol);
221 const std::vector<Real> &xstat,
223 mqq_->updateGradient(obj,x,xstat,tol);
228 const std::vector<Real> &vstat,
230 const std::vector<Real> &xstat,
232 mqq_->updateHessVec(obj,v,vstat,x,xstat,tol);
236 const std::vector<Real> &xstat,
238 return mqq_->getValue(x,xstat,sampler);
242 std::vector<Real> &gstat,
244 const std::vector<Real> &xstat,
246 mqq_->getGradient(g,gstat,x,xstat,sampler);
250 std::vector<Real> &hvstat,
252 const std::vector<Real> &vstat,
254 const std::vector<Real> &xstat,
256 mqq_->getHessVec(hv,hvstat,v,vstat,x,xstat,sampler);
Provides the interface to evaluate objective functions.
virtual void setSample(const std::vector< Real > &point, const Real weight)
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 buildMixedQuantile(const std::vector< Real > &pts, const std::vector< Real > &wts, const ROL::Ptr< PlusFunction< Real > > &pf)
void checkInputs(const ROL::Ptr< Distribution< Real > > &dist=ROL::nullPtr) const
virtual void setStorage(const Ptr< ScalarController< Real >> &value_storage, const Ptr< VectorController< Real >> &gradient_storage)
Defines the linear algebra or vector space interface.
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.
virtual void resetStorage(bool flag=true)
Reset internal storage.
Real computeStatistic(const Ptr< const std::vector< Real >> &xstat) const override
Compute statistic.
void resetStorage(bool flag=true)
Reset internal storage.
void setHessVecStorage(const Ptr< ScalarController< Real >> &gradvec_storage, const Ptr< VectorController< Real >> &hessvec_storage)
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 setSample(const std::vector< Real > &point, const Real weight)
void setStorage(const Ptr< ScalarController< Real >> &value_storage, const Ptr< VectorController< Real >> &gradient_storage)
ROL::Ptr< MixedCVaR< Real > > mqq_
Provides an interface for spectral risk measures.
void buildQuadFromDist(std::vector< Real > &pts, std::vector< Real > &wts, const int nQuad, const ROL::Ptr< Distribution< Real > > &dist) const
ROL::Ptr< PlusFunction< Real > > plusFunction_
SpectralRisk(ROL::ParameterList &parlist)
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.
SpectralRisk(const ROL::Ptr< Distribution< Real > > &dist, const int nQuad, const ROL::Ptr< PlusFunction< Real > > &pf)
void updateValue(Objective< Real > &obj, const Vector< Real > &x, const std::vector< Real > &xstat, Real &tol)
Update internal storage for value computation.
Real getValue(const Vector< Real > &x, const std::vector< Real > &xstat, SampleGenerator< Real > &sampler)
Return risk measure value.
void printQuad(const std::vector< Real > &pts, const std::vector< Real > &wts, const bool print=false) const
Provides the interface to implement any functional that maps a random variable to a (extended) real n...
void initialize(const Vector< Real > &x)
Initialize temporary variables.
virtual void setHessVecStorage(const Ptr< ScalarController< Real >> &gradvec_storage, const Ptr< VectorController< Real >> &hessvec_storage)
SpectralRisk(const std::vector< Real > &pts, const std::vector< Real > &wts, const ROL::Ptr< PlusFunction< Real > > &pf)
virtual void initialize(const Vector< Real > &x)
Initialize temporary variables.
void resetStorage(UpdateType type)