44 #ifndef ROL_MONTECARLOGENERATOR_HPP 45 #define ROL_MONTECARLOGENERATOR_HPP 60 std::vector<std::vector<Real> >
data_;
68 const std::vector<ROL::Ptr<Distribution<Real> > >
dist_;
70 Real
ierf(Real input)
const {
71 std::vector<Real> coeff;
74 Real tmp = c * (std::sqrt(pi)/two * input);
78 while (std::abs(tmp) > tol*std::abs(val)) {
80 for (
unsigned i = 0; i < coeff.size(); i++ ) {
81 c += coeff[i]*coeff[coeff.size()-1-i]/((i+1)*(2*i+1));
83 Real ind =
static_cast<Real
>(cnt);
84 tmp = c/(two*ind+one) * std::pow(std::sqrt(pi)/two*input, two*ind+one);
93 return static_cast<Real
>(rand())/static_cast<Real>(RAND_MAX);
96 std::vector<std::vector<Real> >
sample(
int nSamp,
bool store =
true) {
98 const Real
zero(0), one(1), two(2), tol(0.1);
101 std::vector<Real> pts(nSamp*dim,
zero);
104 for (
int i = 0; i < nSamp; ++i) {
106 for (
int j = 0; j < dim; ++j) {
116 for (
int j = 0; j < dim; ++j) {
118 while (std::abs(pts[j + i*dim]) > tol*ROL_INF<Real>()) {
128 int frac = nSamp / nProc;
129 int rem = nSamp % nProc;
130 int N = frac + ((rank < rem) ? 1 : 0);
132 for (
int i = 0; i < rank; ++i) {
133 offset += frac + ((i < rem) ? 1 : 0);
135 std::vector<std::vector<Real> > mypts;
136 std::vector<Real> pt(dim);
137 for (
int i = 0; i < N; ++i) {
139 for (
int j = 0; j < dim; ++j) {
140 pt[j] = pts[j + I*dim];
145 std::vector<Real> mywts(N, one/static_cast<Real>(nSamp));
160 const bool use_SA =
false,
161 const bool adaptive =
false,
162 const int numNewSamps = 0)
176 ROL_TEST_FOR_EXCEPTION(
nSamp_ < nProc, std::invalid_argument,
177 ">>> ERROR (ROL::MonteCarloGenerator): Total number of samples is less than the number of batches!");
182 std::vector<std::vector<Real> > &bounds,
184 const bool use_SA =
false,
185 const bool adaptive =
false,
186 const int numNewSamps = 0)
199 ROL_TEST_FOR_EXCEPTION(
nSamp_ < nProc, std::invalid_argument,
200 ">>> ERROR (ROL::MonteCarloGenerator): Total number of samples is less than the number of batches!");
201 unsigned dim = bounds.size();
204 for (
unsigned j = 0; j < dim; j++ ) {
205 if ( (bounds[j])[0] > (bounds[j])[1] ) {
206 tmp = (bounds[j])[0];
207 (bounds[j])[0] = (bounds[j])[1];
208 (bounds[j])[1] = tmp;
209 data_.push_back(bounds[j]);
211 data_.push_back(bounds[j]);
217 const std::vector<Real> &mean,
218 const std::vector<Real> &
std,
220 const bool use_SA =
false,
221 const bool adaptive =
false,
222 const int numNewSamps = 0 )
235 ROL_TEST_FOR_EXCEPTION(
nSamp_ < nProc, std::invalid_argument,
236 ">>> ERROR (ROL::MonteCarloGenerator): Total number of samples is less than the number of batches!");
237 unsigned dim = mean.size();
239 std::vector<Real> tmp(2,static_cast<Real>(0));
240 for (
unsigned j = 0; j < dim; j++ ) {
243 data_.push_back(tmp);
262 Real
zero(0), one(1), tol(1e-8);
279 return std::sqrt(var/(
nSamp_))*tol;
283 return static_cast<Real
>(0);
289 Real
zero(0), one(1), tol(1e-4);
294 ng = (vals[cnt])->norm();
308 return std::sqrt(var/(
nSamp_))*tol;
312 return static_cast<Real
>(0);
318 Real
zero(0), one(1);
319 std::vector<std::vector<Real> > pts;
321 for (
int i = 0; i < SampleGenerator<Real>::numMySamples(); i++ ) {
326 pts.insert(pts.end(),pts_new.begin(),pts_new.end());
328 std::vector<Real> wts(pts.size(),one/((Real)
nSamp_));
virtual std::vector< Real > getMyPoint(const int i) const
virtual void update(const Vector< Real > &x)
Real computeError(std::vector< ROL::Ptr< Vector< Real > > > &vals, const Vector< Real > &x)
Real ierf(Real input) const
void sumAll(Real *input, Real *output, int dim) const
Defines the linear algebra or vector space interface.
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
MonteCarloGenerator(const int nSamp, const std::vector< ROL::Ptr< Distribution< Real > > > &dist, const ROL::Ptr< BatchManager< Real > > &bman, const bool use_SA=false, const bool adaptive=false, const int numNewSamps=0)
Real computeError(std::vector< Real > &vals)
int numGlobalSamples(void) const
static constexpr Real pi() noexcept
std::vector< std::vector< Real > > sample(int nSamp, bool store=true)
virtual void refine(void)
int numBatches(void) const
MonteCarloGenerator(const int nSamp, std::vector< std::vector< Real > > &bounds, const ROL::Ptr< BatchManager< Real > > &bman, const bool use_SA=false, const bool adaptive=false, const int numNewSamps=0)
MonteCarloGenerator(const int nSamp, const std::vector< Real > &mean, const std::vector< Real > &std, const ROL::Ptr< BatchManager< Real > > &bman, const bool use_SA=false, const bool adaptive=false, const int numNewSamps=0)
void setPoints(std::vector< std::vector< Real > > &p)
const std::vector< ROL::Ptr< Distribution< Real > > > dist_
void update(const Vector< Real > &x)
void broadcast(Real *input, int cnt, int root) const
void setWeights(std::vector< Real > &w)
std::vector< std::vector< Real > > data_