44 #ifndef ROL_LOGQUANTILEQUAD_HPP 45 #define ROL_LOGQUANTILEQUAD_HPP 79 Teuchos::RCP<PlusFunction<Real> >
pf_;
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!");
124 Teuchos::ParameterList& list
125 = parlist.sublist(
"SOL").sublist(
"Risk Measure").sublist(
"Log-Quantile Quadrangle");
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");
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!");
142 Real X = ((deriv == 0) ? x : ((deriv == 1) ? one : zero));
143 return regret(x,deriv) - X;
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!");
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);
165 Real x =
eps_, two(2), p1(0.1), zero(0), one(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" 175 for (
int i = 0; i < 13; i++) {
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
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" 202 for (
int i = 0; i < 13; i++) {
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
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" 230 for (
int i = 0; i < 13; i++) {
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
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" 257 for (
int i = 0; i < 13; i++) {
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
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" 285 for (
int i = 0; i < 13; i++) {
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
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" 312 for (
int i = 0; i < 13; i++) {
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
Real error(Real x, int deriv=0)
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.
void checkInputs(void) const
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.