ROL
ROL_LogQuantileQuadrangle.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_LOGQUANTILEQUAD_HPP
45 #define ROL_LOGQUANTILEQUAD_HPP
46 
47 #include "ROL_ExpectationQuad.hpp"
48 #include "ROL_PlusFunction.hpp"
49 
73 namespace ROL {
74 
75 template<class Real>
76 class LogQuantileQuadrangle : public ExpectationQuad<Real> {
77 private:
78 
79  Teuchos::RCP<PlusFunction<Real> > pf_;
80 
81  Real alpha_;
82  Real rate_;
83  Real eps_;
84 
85  void checkInputs(void) const {
86  Real zero(0), one(1);
87  TEUCHOS_TEST_FOR_EXCEPTION((alpha_ < zero) || (alpha_ >= one), std::invalid_argument,
88  ">>> ERROR (ROL::LogQuantileQuadrangle): Linear growth rate must be between 0 and 1!");
89  TEUCHOS_TEST_FOR_EXCEPTION((rate_ <= zero), std::invalid_argument,
90  ">>> ERROR (ROL::LogQuantileQuadrangle): Exponential growth rate must be positive!");
91  TEUCHOS_TEST_FOR_EXCEPTION((eps_ <= zero), std::invalid_argument,
92  ">>> ERROR (ROL::LogQuantileQuadrangle): Smoothing parameter must be positive!");
93  TEUCHOS_TEST_FOR_EXCEPTION(pf_ == Teuchos::null, std::invalid_argument,
94  ">>> ERROR (ROL::LogQuantileQuadrangle): PlusFunction pointer is null!");
95  }
96 
97 public:
105  LogQuantileQuadrangle(Real alpha, Real rate, Real eps,
106  Teuchos::RCP<PlusFunction<Real> > &pf )
107  : ExpectationQuad<Real>(), alpha_(alpha), rate_(rate), eps_(eps), pf_(pf) {
108  checkInputs();
109  }
110 
122  LogQuantileQuadrangle(Teuchos::ParameterList &parlist)
123  : ExpectationQuad<Real>() {
124  Teuchos::ParameterList& list
125  = parlist.sublist("SOL").sublist("Risk Measure").sublist("Log-Quantile Quadrangle");
126  // Check CVaR inputs
127  alpha_ = list.get<Real>("Slope for Linear Growth");
128  rate_ = list.get<Real>("Rate for Exponential Growth");
129  eps_ = list.get<Real>("Smoothing Parameter");
130  // Build plus function
131  pf_ = Teuchos::rcp( new PlusFunction<Real>(list) );
132  checkInputs();
133  }
134 
135  Real error(Real x, int deriv = 0) {
136  Real zero(0), one(1);
137  TEUCHOS_TEST_FOR_EXCEPTION( (deriv > 2), std::invalid_argument,
138  ">>> ERROR (ROL::LogQuantileQuadrangle::error): deriv greater than 2!");
139  TEUCHOS_TEST_FOR_EXCEPTION( (deriv < 0), std::invalid_argument,
140  ">>> ERROR (ROL::LogQuantileQuadrangle::error): deriv less than 0!");
141 
142  Real X = ((deriv == 0) ? x : ((deriv == 1) ? one : zero));
143  return regret(x,deriv) - X;
144  }
145 
146  Real regret(Real x, int deriv = 0) {
147  Real zero(0), one(1);
148  TEUCHOS_TEST_FOR_EXCEPTION( (deriv > 2), std::invalid_argument,
149  ">>> ERROR (ROL::LogQuantileQuadrangle::regret): deriv greater than 2!");
150  TEUCHOS_TEST_FOR_EXCEPTION( (deriv < 0), std::invalid_argument,
151  ">>> ERROR (ROL::LogQuantileQuadrangle::regret): deriv less than 0!");
152 
153  Real arg = std::exp(rate_*x);
154  Real sarg = rate_*arg;
155  Real reg = (pf_->evaluate(arg-one,deriv) *
156  ((deriv == 0) ? one/rate_ : ((deriv == 1) ? arg : sarg*arg))
157  + ((deriv == 2) ? pf_->evaluate(arg-one,deriv-1)*sarg : zero))
158  + ((deriv%2 == 0) ? -one : one) * alpha_ * pf_->evaluate(-x,deriv);
159  return reg;
160  }
161 
162  void checkRegret(void) {
164  // Check v'(eps)
165  Real x = eps_, two(2), p1(0.1), zero(0), one(1);
166  Real vx(0), vy(0);
167  Real dv = regret(x,1);
168  Real t(1), diff(0), err(0);
169  std::cout << std::right << std::setw(20) << "CHECK REGRET: v'(eps) is correct? \n";
170  std::cout << std::right << std::setw(20) << "t"
171  << std::setw(20) << "v'(x)"
172  << std::setw(20) << "(v(x+t)-v(x-t))/2t"
173  << std::setw(20) << "Error"
174  << "\n";
175  for (int i = 0; i < 13; i++) {
176  vy = regret(x+t,0);
177  vx = regret(x-t,0);
178  diff = (vy-vx)/(two*t);
179  err = std::abs(diff-dv);
180  std::cout << std::scientific << std::setprecision(11) << std::right
181  << std::setw(20) << t
182  << std::setw(20) << dv
183  << std::setw(20) << diff
184  << std::setw(20) << err
185  << "\n";
186  t *= p1;
187  }
188  std::cout << "\n";
189  // check v''(eps)
190  vx = zero;
191  vy = zero;
192  dv = regret(x,2);
193  t = one;
194  diff = zero;
195  err = zero;
196  std::cout << std::right << std::setw(20) << "CHECK REGRET: v''(eps) is correct? \n";
197  std::cout << std::right << std::setw(20) << "t"
198  << std::setw(20) << "v''(x)"
199  << std::setw(20) << "(v'(x+t)-v'(x-t))/2t"
200  << std::setw(20) << "Error"
201  << "\n";
202  for (int i = 0; i < 13; i++) {
203  vy = regret(x+t,1);
204  vx = regret(x-t,1);
205  diff = (vy-vx)/(two*t);
206  err = std::abs(diff-dv);
207  std::cout << std::scientific << std::setprecision(11) << std::right
208  << std::setw(20) << t
209  << std::setw(20) << dv
210  << std::setw(20) << diff
211  << std::setw(20) << err
212  << "\n";
213  t *= p1;
214  }
215  std::cout << "\n";
216  // Check v'(0)
217  x = zero;
218  vx = zero;
219  vy = zero;
220  dv = regret(x,1);
221  t = one;
222  diff = zero;
223  err = zero;
224  std::cout << std::right << std::setw(20) << "CHECK REGRET: v'(0) is correct? \n";
225  std::cout << std::right << std::setw(20) << "t"
226  << std::setw(20) << "v'(x)"
227  << std::setw(20) << "(v(x+t)-v(x-t))/2t"
228  << std::setw(20) << "Error"
229  << "\n";
230  for (int i = 0; i < 13; i++) {
231  vy = regret(x+t,0);
232  vx = regret(x-t,0);
233  diff = (vy-vx)/(two*t);
234  err = std::abs(diff-dv);
235  std::cout << std::scientific << std::setprecision(11) << std::right
236  << std::setw(20) << t
237  << std::setw(20) << dv
238  << std::setw(20) << diff
239  << std::setw(20) << err
240  << "\n";
241  t *= p1;
242  }
243  std::cout << "\n";
244  // check v''(eps)
245  vx = zero;
246  vy = zero;
247  dv = regret(x,2);
248  t = one;
249  diff = zero;
250  err = zero;
251  std::cout << std::right << std::setw(20) << "CHECK REGRET: v''(0) is correct? \n";
252  std::cout << std::right << std::setw(20) << "t"
253  << std::setw(20) << "v''(x)"
254  << std::setw(20) << "(v'(x+t)-v'(x-t))/2t"
255  << std::setw(20) << "Error"
256  << "\n";
257  for (int i = 0; i < 13; i++) {
258  vy = regret(x+t,1);
259  vx = regret(x-t,1);
260  diff = (vy-vx)/(two*t);
261  err = std::abs(diff-dv);
262  std::cout << std::scientific << std::setprecision(11) << std::right
263  << std::setw(20) << t
264  << std::setw(20) << dv
265  << std::setw(20) << diff
266  << std::setw(20) << err
267  << "\n";
268  t *= p1;
269  }
270  std::cout << "\n";
271  // Check v'(0)
272  x = -eps_;
273  vx = zero;
274  vy = zero;
275  dv = regret(x,1);
276  t = one;
277  diff = zero;
278  err = zero;
279  std::cout << std::right << std::setw(20) << "CHECK REGRET: v'(-eps) is correct? \n";
280  std::cout << std::right << std::setw(20) << "t"
281  << std::setw(20) << "v'(x)"
282  << std::setw(20) << "(v(x+t)-v(x-t))/2t"
283  << std::setw(20) << "Error"
284  << "\n";
285  for (int i = 0; i < 13; i++) {
286  vy = regret(x+t,0);
287  vx = regret(x-t,0);
288  diff = (vy-vx)/(two*t);
289  err = std::abs(diff-dv);
290  std::cout << std::scientific << std::setprecision(11) << std::right
291  << std::setw(20) << t
292  << std::setw(20) << dv
293  << std::setw(20) << diff
294  << std::setw(20) << err
295  << "\n";
296  t *= p1;
297  }
298  std::cout << "\n";
299  // check v''(eps)
300  vx = zero;
301  vy = zero;
302  dv = regret(x,2);
303  t = one;
304  diff = zero;
305  err = zero;
306  std::cout << std::right << std::setw(20) << "CHECK REGRET: v''(-eps) is correct? \n";
307  std::cout << std::right << std::setw(20) << "t"
308  << std::setw(20) << "v''(x)"
309  << std::setw(20) << "(v'(x+t)-v'(x-t))/2t"
310  << std::setw(20) << "Error"
311  << "\n";
312  for (int i = 0; i < 13; i++) {
313  vy = regret(x+t,1);
314  vx = regret(x-t,1);
315  diff = (vy-vx)/(two*t);
316  err = std::abs(diff-dv);
317  std::cout << std::scientific << std::setprecision(11) << std::right
318  << std::setw(20) << t
319  << std::setw(20) << dv
320  << std::setw(20) << diff
321  << std::setw(20) << err
322  << "\n";
323  t *= p1;
324  }
325  std::cout << "\n";
326  }
327 
328 };
329 
330 }
331 #endif
Provides a general interface for risk measures generated through the expectation risk quadrangle...
virtual void checkRegret(void)
Run default derivative tests for the scalar regret function.
void checkRegret(void)
Run default derivative tests for the scalar regret function.
Teuchos::RCP< PlusFunction< Real > > pf_
Provides an interface for the conditioanl entropic risk using the expectation risk quadrangle...
LogQuantileQuadrangle(Real alpha, Real rate, Real eps, Teuchos::RCP< PlusFunction< Real > > &pf)
Constructor.
LogQuantileQuadrangle(Teuchos::ParameterList &parlist)
Constructor.
Real regret(Real x, int deriv=0)
Evaluate the scalar regret function at x.