ROL
ROL_QuantileQuadrangle.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Rapid Optimization Library (ROL) Package
5 // Copyright (2014) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact lead developers:
38 // Drew Kouri (dpkouri@sandia.gov) and
39 // Denis Ridzal (dridzal@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
44 #ifndef ROL_SMOOTHCVARQUAD_HPP
45 #define ROL_SMOOTHCVARQUAD_HPP
46 
47 #include "ROL_ExpectationQuad.hpp"
48 #include "ROL_PlusFunction.hpp"
49 
94 namespace ROL {
95 
96 template<class Real>
97 class QuantileQuadrangle : public ExpectationQuad<Real> {
98 private:
99 
100  Teuchos::RCP<PlusFunction<Real> > pf_;
101 
102  Real prob_;
103  Real lam_;
104  Real eps_;
105 
106  Real alpha_;
107  Real beta_;
108 
109  void checkInputs(void) const {
110  Real zero(0), one(1);
111  TEUCHOS_TEST_FOR_EXCEPTION((prob_ <= zero) || (prob_ >= one), std::invalid_argument,
112  ">>> ERROR (ROL::QuantileQuadrangle): Confidence level must be between 0 and 1!");
113  TEUCHOS_TEST_FOR_EXCEPTION((lam_ < zero) || (lam_ > one), std::invalid_argument,
114  ">>> ERROR (ROL::QuantileQuadrangle): Convex combination parameter must be positive!");
115  TEUCHOS_TEST_FOR_EXCEPTION((eps_ <= zero), std::invalid_argument,
116  ">>> ERROR (ROL::QuantileQuadrangle): Smoothing parameter must be positive!");
117  TEUCHOS_TEST_FOR_EXCEPTION(pf_ == Teuchos::null, std::invalid_argument,
118  ">>> ERROR (ROL::QuantileQuadrangle): PlusFunction pointer is null!");
119  }
120 
121  void setParameters(void) {
122  Real one(1);
123  alpha_ = lam_;
124  beta_ = (one-alpha_*prob_)/(one-prob_);
125  }
126 
127 public:
134  QuantileQuadrangle(Real prob, Real eps, Teuchos::RCP<PlusFunction<Real> > &pf )
135  : ExpectationQuad<Real>(), prob_(prob), lam_(0), eps_(eps), pf_(pf) {
136  checkInputs();
137  setParameters();
138  }
139 
149  QuantileQuadrangle(Real prob, Real lam, Real eps,
150  Teuchos::RCP<PlusFunction<Real> > &pf )
151  : ExpectationQuad<Real>(), prob_(prob), lam_(lam), eps_(eps), pf_(pf) {
152  checkInputs();
153  setParameters();
154  }
155 
167  QuantileQuadrangle(Teuchos::ParameterList &parlist) : ExpectationQuad<Real>() {
168  Teuchos::ParameterList& list
169  = parlist.sublist("SOL").sublist("Risk Measure").sublist("Quantile-Based Quadrangle");
170  // Get CVaR parameters
171  prob_ = list.get<Real>("Confidence Level");
172  lam_ = list.get<Real>("Convex Combination Parameter");
173  eps_ = list.get<Real>("Smoothing Parameter");
174  // Build plus function
175  pf_ = Teuchos::rcp( new PlusFunction<Real>(list) );
176  // Check inputs
177  checkInputs();
178  setParameters();
179  }
180 
181  Real error(Real x, int deriv = 0) {
182  Real one(1);
183  Real err = (beta_-one)*pf_->evaluate(x,deriv)
184  + ((deriv%2) ? -one : one)*(one-alpha_)*pf_->evaluate(-x,deriv);
185  return err;
186  }
187 
188  Real regret(Real x, int deriv = 0) {
189  Real zero(0), one(1);
190  Real X = ((deriv==0) ? x : ((deriv==1) ? one : zero));
191  Real reg = error(x,deriv) + X;
192  return reg;
193  }
194 
195  void checkRegret(void) {
197  // Check v'(eps)
198  Real x = eps_, two(2), p1(0.1), one(1), zero(0);
199  Real vx(0), vy(0);
200  Real dv = regret(x,1);
201  Real t(1), diff(0), err(0);
202  std::cout << std::right << std::setw(20) << "CHECK REGRET: v'(eps) is correct? \n";
203  std::cout << std::right << std::setw(20) << "t"
204  << std::setw(20) << "v'(x)"
205  << std::setw(20) << "(v(x+t)-v(x-t))/2t"
206  << std::setw(20) << "Error"
207  << "\n";
208  for (int i = 0; i < 13; i++) {
209  vy = regret(x+t,0);
210  vx = regret(x-t,0);
211  diff = (vy-vx)/(two*t);
212  err = std::abs(diff-dv);
213  std::cout << std::scientific << std::setprecision(11) << std::right
214  << std::setw(20) << t
215  << std::setw(20) << dv
216  << std::setw(20) << diff
217  << std::setw(20) << err
218  << "\n";
219  t *= p1;
220  }
221  std::cout << "\n";
222  // check v''(eps)
223  vx = zero;
224  vy = zero;
225  dv = regret(x,2);
226  t = one;
227  diff = zero;
228  err = zero;
229  std::cout << std::right << std::setw(20) << "CHECK REGRET: v''(eps) is correct? \n";
230  std::cout << std::right << std::setw(20) << "t"
231  << std::setw(20) << "v''(x)"
232  << std::setw(20) << "(v'(x+t)-v'(x-t))/2t"
233  << std::setw(20) << "Error"
234  << "\n";
235  for (int i = 0; i < 13; i++) {
236  vy = regret(x+t,1);
237  vx = regret(x-t,1);
238  diff = (vy-vx)/(two*t);
239  err = std::abs(diff-dv);
240  std::cout << std::scientific << std::setprecision(11) << std::right
241  << std::setw(20) << t
242  << std::setw(20) << dv
243  << std::setw(20) << diff
244  << std::setw(20) << err
245  << "\n";
246  t *= p1;
247  }
248  std::cout << "\n";
249  // Check v'(0)
250  x = zero;
251  vx = zero;
252  vy = zero;
253  dv = regret(x,1);
254  t = one;
255  diff = zero;
256  err = zero;
257  std::cout << std::right << std::setw(20) << "CHECK REGRET: v'(0) is correct? \n";
258  std::cout << std::right << std::setw(20) << "t"
259  << std::setw(20) << "v'(x)"
260  << std::setw(20) << "(v(x+t)-v(x-t))/2t"
261  << std::setw(20) << "Error"
262  << "\n";
263  for (int i = 0; i < 13; i++) {
264  vy = regret(x+t,0);
265  vx = regret(x-t,0);
266  diff = (vy-vx)/(two*t);
267  err = std::abs(diff-dv);
268  std::cout << std::scientific << std::setprecision(11) << std::right
269  << std::setw(20) << t
270  << std::setw(20) << dv
271  << std::setw(20) << diff
272  << std::setw(20) << err
273  << "\n";
274  t *= p1;
275  }
276  std::cout << "\n";
277  // check v''(eps)
278  vx = zero;
279  vy = zero;
280  dv = regret(x,2);
281  t = one;
282  diff = zero;
283  err = zero;
284  std::cout << std::right << std::setw(20) << "CHECK REGRET: v''(0) is correct? \n";
285  std::cout << std::right << std::setw(20) << "t"
286  << std::setw(20) << "v''(x)"
287  << std::setw(20) << "(v'(x+t)-v'(x-t))/2t"
288  << std::setw(20) << "Error"
289  << "\n";
290  for (int i = 0; i < 13; i++) {
291  vy = regret(x+t,1);
292  vx = regret(x-t,1);
293  diff = (vy-vx)/(two*t);
294  err = std::abs(diff-dv);
295  std::cout << std::scientific << std::setprecision(11) << std::right
296  << std::setw(20) << t
297  << std::setw(20) << dv
298  << std::setw(20) << diff
299  << std::setw(20) << err
300  << "\n";
301  t *= p1;
302  }
303  std::cout << "\n";
304  // Check v'(0)
305  x = -eps_;
306  vx = zero;
307  vy = zero;
308  dv = regret(x,1);
309  t = one;
310  diff = zero;
311  err = zero;
312  std::cout << std::right << std::setw(20) << "CHECK REGRET: v'(-eps) is correct? \n";
313  std::cout << std::right << std::setw(20) << "t"
314  << std::setw(20) << "v'(x)"
315  << std::setw(20) << "(v(x+t)-v(x-t))/2t"
316  << std::setw(20) << "Error"
317  << "\n";
318  for (int i = 0; i < 13; i++) {
319  vy = regret(x+t,0);
320  vx = regret(x-t,0);
321  diff = (vy-vx)/(two*t);
322  err = std::abs(diff-dv);
323  std::cout << std::scientific << std::setprecision(11) << std::right
324  << std::setw(20) << t
325  << std::setw(20) << dv
326  << std::setw(20) << diff
327  << std::setw(20) << err
328  << "\n";
329  t *= p1;
330  }
331  std::cout << "\n";
332  // check v''(eps)
333  vx = zero;
334  vy = zero;
335  dv = regret(x,2);
336  t = one;
337  diff = zero;
338  err = zero;
339  std::cout << std::right << std::setw(20) << "CHECK REGRET: v''(-eps) is correct? \n";
340  std::cout << std::right << std::setw(20) << "t"
341  << std::setw(20) << "v''(x)"
342  << std::setw(20) << "(v'(x+t)-v'(x-t))/2t"
343  << std::setw(20) << "Error"
344  << "\n";
345  for (int i = 0; i < 13; i++) {
346  vy = regret(x+t,1);
347  vx = regret(x-t,1);
348  diff = (vy-vx)/(two*t);
349  err = std::abs(diff-dv);
350  std::cout << std::scientific << std::setprecision(11) << std::right
351  << std::setw(20) << t
352  << std::setw(20) << dv
353  << std::setw(20) << diff
354  << std::setw(20) << err
355  << "\n";
356  t *= p1;
357  }
358  std::cout << "\n";
359  }
360 
361 };
362 
363 }
364 #endif
QuantileQuadrangle(Teuchos::ParameterList &parlist)
Constructor.
Provides a general interface for risk measures generated through the expectation risk quadrangle...
Provides an interface for a convex combination of the expected value and the conditional value-at-ris...
virtual void checkRegret(void)
Run default derivative tests for the scalar regret function.
Real regret(Real x, int deriv=0)
Evaluate the scalar regret function at x.
Real error(Real x, int deriv=0)
Teuchos::RCP< PlusFunction< Real > > pf_
QuantileQuadrangle(Real prob, Real eps, Teuchos::RCP< PlusFunction< Real > > &pf)
Constructor.
void checkRegret(void)
Run default derivative tests for the scalar regret function.
QuantileQuadrangle(Real prob, Real lam, Real eps, Teuchos::RCP< PlusFunction< Real > > &pf)
Constructor.