10 #ifndef ROL_CDFOBJECTIVE_H
11 #define ROL_CDFOBJECTIVE_H
17 #include "ROL_Ptr.hpp"
26 ROL::Ptr<BatchManager<Real> >
bman_;
29 std::vector<ROL::Ptr<Distribution<Real> > >
dist_;
51 wts_[0] = 0.1527533871307258;
pts_[0] = -0.0765265211334973;
52 wts_[1] = 0.1527533871307258;
pts_[1] = 0.0765265211334973;
53 wts_[2] = 0.1491729864726037;
pts_[2] = -0.2277858511416451;
54 wts_[3] = 0.1491729864726037;
pts_[3] = 0.2277858511416451;
55 wts_[4] = 0.1420961093183820;
pts_[4] = -0.3737060887154195;
56 wts_[5] = 0.1420961093183820;
pts_[5] = 0.3737060887154195;
57 wts_[6] = 0.1316886384491766;
pts_[6] = -0.5108670019508271;
58 wts_[7] = 0.1316886384491766;
pts_[7] = 0.5108670019508271;
59 wts_[8] = 0.1181945319615184;
pts_[8] = -0.6360536807265150;
60 wts_[9] = 0.1181945319615184;
pts_[9] = 0.6360536807265150;
61 wts_[10] = 0.1019301198172404;
pts_[10] = -0.7463319064601508;
62 wts_[11] = 0.1019301198172404;
pts_[11] = 0.7463319064601508;
63 wts_[12] = 0.0832767415767048;
pts_[12] = -0.8391169718222188;
64 wts_[13] = 0.0832767415767048;
pts_[13] = 0.8391169718222188;
65 wts_[14] = 0.0626720483341091;
pts_[14] = -0.9122344282513259;
66 wts_[15] = 0.0626720483341091;
pts_[15] = 0.9122344282513259;
67 wts_[16] = 0.0406014298003869;
pts_[16] = -0.9639719272779138;
68 wts_[17] = 0.0406014298003869;
pts_[17] = 0.9639719272779138;
69 wts_[18] = 0.0176140071391521;
pts_[18] = -0.9931285991850949;
70 wts_[19] = 0.0176140071391521;
pts_[19] = 0.9931285991850949;
81 Real val = 0, hs = 0, xpt = 0, xwt = 0, sum = 0, half(0.5), one(1);
82 for (
int k = 0; k < numSamples; k++) {
87 bman_->sumAll(&val,&sum,1);
91 Real
gradientCDF(std::vector<Real> &gradx, std::vector<Real> &gradp,
92 const int dim,
const Real loc,
96 gradx.resize(numSamples,0); gradp.resize(numSamples,0);
97 Real val = 0, hs = 0, xpt = 0, xwt = 0, sum = 0, half(0.5), one(1);
98 for (
int k = 0; k < numSamples; k++) {
103 * std::exp(-std::pow((loc-xpt)/(
sqrt2_*scale_),2));
106 bman_->sumAll(&val,&sum,1);
110 Real
hessVecCDF(std::vector<Real> &hvxx, std::vector<Real> &hvxp, std::vector<Real> &hvpx,
111 std::vector<Real> &gradx, std::vector<Real> &gradp,
112 Real &sumx, Real &sump,
113 const int dim,
const Real loc,
119 hvxx.resize(numSamples,0); hvxp.resize(numSamples,0); hvpx.resize(numSamples,0);
120 gradx.resize(numSamples,0); gradp.resize(numSamples,0);
122 std::vector<Real> psum(3,0), out(3,0);
123 Real val = 0, hs = 0, dval = 0, scale3 = std::pow(
scale_,3);
124 Real xpt = 0, xwt = 0, vpt = 0, vwt = 0, half(0.5), one(1);
125 for (
int k = 0; k < numSamples; k++) {
130 dval = std::exp(-std::pow((loc-xpt)/(
sqrt2_*scale_),2));
133 hvxx[k] = -(xwt/(
sqrt2_*
sqrtpi_*scale3)) * dval * (loc-xpt) * vpt;
136 psum[1] += vpt*gradx[k];
137 psum[2] += vwt*gradp[k];
139 bman_->sumAll(&psum[0],&out[0],3);
140 val = out[0]; sumx = out[1]; sump = out[2];
147 const Real scale = 1.e-2,
148 const bool optProb =
true,
const bool optAtom =
true)
165 Real val(0), diff(0), pt(0), wt(0), meas(0), lb(0), one(1);
169 meas = ((meas > ROL_EPSILON<Real>()) ? meas : one);
171 pt = meas*
pts_[k] + lb;
174 val += wt*std::pow(diff,2);
185 const int numSamples = prob.getNumMyAtoms();
186 std::vector<Real> gradx(numSamples,0.), gradp(numSamples,0);
187 Real diff(0), pt(0), wt(0), meas(0), lb(0), val(0), one(1);
188 std::vector<Real> val_wt(numSamples,0), tmp(
dimension_,0);
189 std::vector<std::vector<Real> > val_pt(numSamples,tmp);
193 meas = ((meas > ROL_EPSILON<Real>()) ? meas : one);
195 pt = meas*
pts_[k] + lb;
198 diff = (val-
dist_[d]->evaluateCDF(pt));
199 for (
int j = 0; j < numSamples; j++) {
200 (val_pt[j])[d] += wt * diff * gradx[j];
201 val_wt[j] += wt * diff * gradp[j];
208 for (
int k = 0; k < numSamples; k++) {
210 gprob.setProbability(k,val_wt[k]);
226 const int numSamples = prob.getNumMyAtoms();
227 std::vector<Real> hvxx(numSamples,0), hvxp(numSamples,0), hvpx(numSamples,0);
228 std::vector<Real> gradx(numSamples,0), gradp(numSamples,0);
229 Real diff(0), pt(0), wt(0), meas(0), lb(0), val(0), sumx(0), sump(0), one(1);
230 std::vector<Real> val_wt(numSamples,0), tmp(
dimension_,0);
231 std::vector<std::vector<Real> > val_pt(numSamples,tmp);
235 meas = ((meas > ROL_EPSILON<Real>()) ? meas : one);
237 pt = meas*
pts_[k] + lb;
239 val =
hessVecCDF(hvxx,hvxp,hvpx,gradx,gradp,sumx,sump,d,pt,prob,atom,vprob,vatom);
240 diff = (val-
dist_[d]->evaluateCDF(pt));
241 for (
int j = 0; j < numSamples; j++) {
242 (val_pt[j])[d] += wt * ( (sump + sumx) * gradx[j] + diff * (hvxx[j] + hvxp[j]) );
243 val_wt[j] += wt * ( (sump + sumx) * gradp[j] + diff * hvpx[j] );
250 for (
int k = 0; k < numSamples; k++) {
252 hprob.setProbability(k,val_wt[k]);
Provides the interface to evaluate objective functions.
std::vector< Real > lowerBound_
Real hessVecCDF(std::vector< Real > &hvxx, std::vector< Real > &hvxp, std::vector< Real > &hvpx, std::vector< Real > &gradx, std::vector< Real > &gradp, Real &sumx, Real &sump, const int dim, const Real loc, const ProbabilityVector< Real > &prob, const AtomVector< Real > &atom, const ProbabilityVector< Real > &vprob, const AtomVector< Real > &vatom) const
const ROL::Ptr< const ProbabilityVector< Real > > getProbabilityVector(void) const
Provides the std::vector implementation of the ROL::Vector interface.
Real value(const Vector< Real > &x, Real &tol)
Compute value.
const Real getProbability(const int i) const
ROL::Ptr< const std::vector< Real > > getAtom(const int i) const
void setAtom(const int i, const std::vector< Real > &pt)
CDFObjective(const std::vector< ROL::Ptr< Distribution< Real > > > &dist, const ROL::Ptr< BatchManager< Real > > &bman, const Real scale=1.e-2, const bool optProb=true, const bool optAtom=true)
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
void initializeQuadrature(void)
Provides the std::vector implementation of the ROL::Vector interface.
Real gradientCDF(std::vector< Real > &gradx, std::vector< Real > &gradp, const int dim, const Real loc, const ProbabilityVector< Real > &prob, const AtomVector< Real > &atom) const
std::vector< Real > upperBound_
Provides the std::vector implementation of the ROL::Vector interface.
std::vector< ROL::Ptr< Distribution< Real > > > dist_
const ROL::Ptr< const AtomVector< Real > > getAtomVector(void) const
ROL::Ptr< BatchManager< Real > > bman_
Real valueCDF(const int dim, const Real loc, const ProbabilityVector< Real > &prob, const AtomVector< Real > &atom) const
int getNumMyAtoms(void) const