44 #ifndef ROL_SPECTRALRISK_HPP 45 #define ROL_SPECTRALRISK_HPP 91 Teuchos::RCP<MixedQuantileQuadrangle<Real> >
mqq_;
98 TEUCHOS_TEST_FOR_EXCEPTION(
plusFunction_ == Teuchos::null, std::invalid_argument,
99 ">>> ERROR (ROL::SpectralRisk): PlusFunction pointer is null!");
100 if ( dist != Teuchos::null) {
101 Real lb = dist->lowerBound();
102 Real ub = dist->upperBound();
103 TEUCHOS_TEST_FOR_EXCEPTION(lb < static_cast<Real>(0), std::invalid_argument,
104 ">>> ERROR (ROL::SpectralRisk): Distribution lower bound less than zero!");
105 TEUCHOS_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());
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;
172 std::vector<Real> wts(nQuad), pts(nQuad);
183 Teuchos::ParameterList &list
184 = parlist.sublist(
"SOL").sublist(
"Risk Measure").sublist(
"Spectral Risk");
185 int nQuad = list.get(
"Number of Quadrature Points",5);
186 bool print = list.get(
"Print Quadrature to Screen",
false);
188 Teuchos::RCP<Distribution<Real> > dist = DistributionFactory<Real>(list);
192 std::vector<Real> wts(nQuad), pts(nQuad);
201 SpectralRisk(
const std::vector<Real> &pts,
const std::vector<Real> &wts,
210 std::vector<Real> xstat;
213 int nQuad =
static_cast<int>(
wts_.size());
214 for (
int i = 0; i < nQuad; ++i) {
215 stat +=
wts_[i] * xstat[i];
226 mqq_->reset(x0,x,v0,v);
229 void update(
const Real val,
const Real weight) {
230 mqq_->update(val,weight);
234 mqq_->update(val,g,weight);
239 mqq_->update(val,g,gv,hv,weight);
243 return mqq_->getValue(sampler);
247 mqq_->getGradient(g,sampler);
251 mqq_->getHessVec(hv,sampler);
void printQuad(const std::vector< Real > &pts, const std::vector< Real > &wts, const bool print=false) const
Provides an interface for a convex combination of conditional value-at-risks.
void reset(Teuchos::RCP< Vector< Real > > &x0, const Vector< Real > &x)
Reset internal risk measure storage. Called for value and gradient computation.
SpectralRisk(Teuchos::ParameterList &parlist)
Real computeStatistic(const Vector< Real > &x) const
void getGradient(Vector< Real > &g, SampleGenerator< Real > &sampler)
Return risk measure (sub)gradient.
Defines the linear algebra or vector space interface.
void buildQuadFromDist(std::vector< Real > &pts, std::vector< Real > &wts, const int nQuad, const Teuchos::RCP< Distribution< Real > > &dist) const
SpectralRisk(const Teuchos::RCP< Distribution< Real > > &dist, const int nQuad, const Teuchos::RCP< PlusFunction< Real > > &pf)
Provides an interface for spectral risk measures.
void update(const Real val, const Real weight)
Update internal risk measure storage for value computation.
void reset(Teuchos::RCP< Vector< Real > > &x0, const Vector< Real > &x, Teuchos::RCP< Vector< Real > > &v0, const Vector< Real > &v)
Reset internal risk measure storage. Called for Hessian-times-a-vector computation.
Teuchos::RCP< PlusFunction< Real > > plusFunction_
void checkInputs(Teuchos::RCP< Distribution< Real > > &dist=Teuchos::null) const
void buildMixedQuantile(const std::vector< Real > &pts, const std::vector< Real > &wts, const Teuchos::RCP< PlusFunction< Real > > &pf)
void update(const Real val, const Vector< Real > &g, const Real gv, const Vector< Real > &hv, const Real weight)
Update internal risk measure storage for Hessian-time-a-vector computation.
Teuchos::RCP< MixedQuantileQuadrangle< Real > > mqq_
void update(const Real val, const Vector< Real > &g, const Real weight)
Update internal risk measure storage for gradient computation.
SpectralRisk(const std::vector< Real > &pts, const std::vector< Real > &wts, const Teuchos::RCP< PlusFunction< Real > > &pf)
Provides the interface to implement risk measures.
Real getValue(SampleGenerator< Real > &sampler)
Return risk measure value.
void getHessVec(Vector< Real > &hv, SampleGenerator< Real > &sampler)
Return risk measure Hessian-times-a-vector.