ROL
ROL_Parabolic.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef ROL_PARABOLIC_HPP
11 #define ROL_PARABOLIC_HPP
12 
13 #include "ROL_Distribution.hpp"
14 #include "ROL_ParameterList.hpp"
15 
16 namespace ROL {
17 
18 template<class Real>
19 class Parabolic : public Distribution<Real> {
20 private:
21  Real a_;
22  Real b_;
23 
24 public:
25  Parabolic(const Real a = 0., const Real b = 1.)
26  : a_(std::min(a,b)), b_(std::max(a,b)) {}
27 
28  Parabolic(ROL::ParameterList &parlist) {
29  a_ = parlist.sublist("SOL").sublist("Distribution").sublist("Parabolic").get("Lower Bound",0.);
30  b_ = parlist.sublist("SOL").sublist("Distribution").sublist("Parabolic").get("Upper Bound",1.);
31  Real tmp = a_;
32  a_ = std::min(a_,b_);
33  b_ = std::max(b_,tmp);
34  }
35 
36  Real evaluatePDF(const Real input) const {
37  Real scale = 6.0/std::pow(b_-a_,3.0);
38  return ((input >= a_ && input <= b_) ? scale*(input-a_)*(b_-input) : 0.);
39  }
40 
41  Real evaluateCDF(const Real input) const {
42  Real d1 = b_-a_, d2 = d1*d1, d3 = d2*d1;
43  Real v1 = input-a_, v2 = v1*v1, v3 = v1*v2;
44  return ((input < a_) ? 0. : ((input > b_) ? 1. :
45  3.0*v2/d2 - 2.0*v3/d3));
46  }
47 
48  Real integrateCDF(const Real input) const {
49  Real d0 = b_+a_, d1 = b_-a_, d2 = d1*d1, d3 = d2*d1;
50  Real v1 = input-a_, v2 = v1*v1, v3 = v1*v2, v4 = v1*v3;
51  return ((input < a_) ? 0. :
52  ((input > b_) ? input - 0.5*d0 :
53  v3/d2 - 0.5*v4/d3));
54  }
55 
56  Real invertCDF(const Real input) const {
57  Real a = a_-b_, b = a_+b_, c = 0.;
58  Real fa = evaluateCDF(a) - input;
59  Real fc = 0.;
60  Real sa = ((fa < 0.) ? -1. : ((fa > 0.) ? 1. : 0.));
61  Real sc = 0.;
62  for (size_t i = 0; i < 100; i++) {
63  c = (a+b)*0.5;
64  fc = evaluateCDF(c) - input;
65  sc = ((fc < 0.) ? -1. : ((fc > 0.) ? 1. : 0.));
66  if ( fc == 0. || (b-a)*0.5 < ROL_EPSILON<Real>() ) {
67  break;
68  }
69  if ( sc == sa ) { a = c; fa = fc; sa = sc; }
70  else { b = c; }
71  }
72  return c;
73  }
74 
75  Real moment(const size_t m) const {
76  Real p = (Real)m;
77  Real a1 = std::pow(a_,p+1), b1 = std::pow(b_,p+1);
78  Real a2 = a1*a_, b2 = b1*b_;
79  Real a3 = a2*a_, b3 = b2*b_;
80  return 6./std::pow(b_-a_,3)
81  * (-(b3-a3)/(p+3) + (a_+b_)*(b2-a2)/(p+2) - a_*b_*(b1-a1)/(p+1));
82  }
83 
84  Real lowerBound(void) const {
85  return a_;
86  }
87 
88  Real upperBound(void) const {
89  return b_;
90  }
91 
92  void test(std::ostream &outStream = std::cout ) const {
93  size_t size = 5;
94  std::vector<Real> X(size,0.);
95  std::vector<int> T(size,0);
96  X[0] = a_-4.0*(Real)rand()/(Real)RAND_MAX;
97  T[0] = 0;
98  X[1] = a_;
99  T[1] = 1;
100  X[2] = (b_-a_)*(Real)rand()/(Real)RAND_MAX + a_;
101  T[2] = 0;
102  X[3] = b_;
103  T[3] = 1;
104  X[4] = b_+4.0*(Real)rand()/(Real)RAND_MAX;
105  T[4] = 0;
106  Distribution<Real>::test(X,T,outStream);
107  }
108 };
109 
110 }
111 
112 #endif
Parabolic(const Real a=0., const Real b=1.)
Real moment(const size_t m) const
Real lowerBound(void) const
Real integrateCDF(const Real input) const
Real evaluatePDF(const Real input) const
virtual void test(std::ostream &outStream=std::cout) const
Real invertCDF(const Real input) const
Real evaluateCDF(const Real input) const
void test(std::ostream &outStream=std::cout) const
Parabolic(ROL::ParameterList &parlist)
Real upperBound(void) const