44 #ifndef ROL_TRUNCATEDGAUSSIAN_HPP 45 #define ROL_TRUNCATEDGAUSSIAN_HPP 49 #include "Teuchos_RCP.hpp" 50 #include "Teuchos_ParameterList.hpp" 62 Teuchos::RCP<Gaussian<Real> >
gauss_;
72 const Real mean = 0.,
const Real sdev = 1.)
73 :
a_((lo < up) ? lo : up),
b_((up > lo) ? up : lo),
83 Teuchos::ParameterList TGlist
84 = parlist.sublist(
"SOL").sublist(
"Distribution").sublist(
"Truncated Gaussian");
85 a_ = TGlist.get(
"Lower Bound",-1.);
86 b_ = TGlist.get(
"Upper Bound",1.);
89 b_ = std::max(
b_,tmp);
91 mean_ = TGlist.get(
"Mean",0.);
92 sdev_ = TGlist.get(
"Standard Deviation",1.);
104 return ((input <=
a_) ? 0.0 : ((input >=
b_) ? 0.0 :
110 return ((input <=
a_) ? 0.0 : ((input >=
b_) ? 1.0 :
115 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
116 ">>> ERROR (ROL::TruncatedGaussian): Truncated Gaussian integrateCDF not implemented!");
117 return ((input < 0.5*(
a_+
b_)) ? 0.0 : input - 0.5*(
a_+
b_));
132 case 1: val = mean;
break;
133 case 2: val = var*(1.+(
alpha_*phiA-
beta_*phiB)/
Z_-std::pow((phiA-phiB)/
Z_,2))+mean*mean;
break;
135 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
136 ">>> ERROR (ROL::TruncatedGaussian): Truncated Gaussian moment not implemented for m > 2!");
149 void test(std::ostream &outStream = std::cout )
const {
151 std::vector<Real> X(size,0.);
152 std::vector<int> T(size,0);
153 X[0] =
a_-4.0*(Real)rand()/(Real)RAND_MAX;
157 X[2] = (
b_-
a_)*(Real)rand()/(Real)RAND_MAX +
a_;
161 X[4] =
b_+4.0*(Real)rand()/(Real)RAND_MAX;
Real moment(const size_t m) const
Real invertCDF(const Real input) const
Teuchos::RCP< Gaussian< Real > > gauss_
void test(std::ostream &outStream=std::cout) const
virtual void test(std::ostream &outStream=std::cout) const
Real evaluateCDF(const Real input) const
Real integrateCDF(const Real input) const
Real lowerBound(void) const
TruncatedGaussian(Teuchos::ParameterList &parlist)
Real upperBound(void) const
TruncatedGaussian(const Real lo=-1., const Real up=1., const Real mean=0., const Real sdev=1.)
Real evaluatePDF(const Real input) const