86 :
RiskMeasure<Real>(), plusFunction_(pf), xvar_(0.0), vvar_(0.0),
87 pnorm_(0.0), dpnorm_(0.0), dpnorm1_(0.0), pgv_(0.0), pgv1_(0.0),
89 prob_ = ((prob >= 0.0) ? ((prob <= 1.0) ? prob : 0.5) : 0.5);
90 coeff_ = ((coeff >= 0.0) ? ((coeff <= 1.0) ? coeff : 1.0) : 1.0);
91 order_ = ((order < 2) ? 2 : order);
94 HMCR( Teuchos::ParameterList &parlist )
96 pnorm_(0.0), dpnorm_(0.0), dpnorm1_(0.0), pgv_(0.0), pgv1_(0.0),
98 Teuchos::ParameterList &list
99 = parlist.sublist(
"SOL").sublist(
"Risk Measure").sublist(
"HMCR");
101 Real prob = list.get(
"Confidence Level",0.5);
102 prob_ = ((prob >= 0.0) ? ((prob <= 1.0) ? prob : 0.5) : 0.5);
103 Real coeff = list.get(
"Convex Combination Parameter",0.5);
104 coeff_ = ((coeff >= 0.0) ? ((coeff <= 1.0) ? coeff : 1.0) : 1.0);
105 unsigned order = list.get(
"Order",2);
106 order_ = ((order < 2) ? 2 : order);
117 dualVector_ = (x0->dual()).clone();
118 pHMdualVec0_ = (x0->dual()).clone();
119 HMdualVec0_ = (x0->dual()).clone();
120 pHMdualVec1_ = (x0->dual()).clone();
121 HMdualVec1_ = (x0->dual()).clone();
122 pHMdualVec2_ = (x0->dual()).clone();
123 HMdualVec2_ = (x0->dual()).clone();
124 pHMdualVec3_ = (x0->dual()).clone();
125 HMdualVec3_ = (x0->dual()).clone();
130 pHMdualVec0_->zero(); HMdualVec0_->zero();
144 pHMdualVec1_->zero(); HMdualVec1_->zero();
145 pHMdualVec2_->zero(); HMdualVec2_->zero();
146 pHMdualVec3_->zero(); HMdualVec3_->zero();
153 void update(
const Real val,
const Real weight) {
157 Real pf = plusFunction_->evaluate(val-xvar_,0);
158 pnorm_ += weight*std::pow(pf,(Real)order_);
165 Real pf0 = plusFunction_->evaluate(val-xvar_,0);
166 Real pf1 = plusFunction_->evaluate(val-xvar_,1);
168 Real rorder0 = (Real)order_;
169 Real rorder1 = (Real)order_-1.0;
171 Real pf0p0 = std::pow(pf0,rorder0);
172 Real pf0p1 = std::pow(pf0,rorder1);
174 pnorm_ += weight*pf0p0;
175 dpnorm_ += weight*pf0p1*pf1;
177 pHMdualVec0_->axpy(weight*pf0p1*pf1,g);
185 Real pf0 = plusFunction_->evaluate(val-xvar_,0);
186 Real pf1 = plusFunction_->evaluate(val-xvar_,1);
187 Real pf2 = plusFunction_->evaluate(val-xvar_,2);
189 Real rorder0 = (Real)order_;
190 Real rorder1 = (Real)order_-1.0;
191 Real rorder2 = (Real)order_-2.0;
193 Real pf0p0 = std::pow(pf0,rorder0);
194 Real pf0p1 = std::pow(pf0,rorder1);
195 Real pf0p2 = std::pow(pf0,rorder2);
197 Real coeff1 = pf0p1*pf1;
198 Real coeff2 = rorder1*pf0p2*pf1*pf1 + pf0p1*pf2;
200 pnorm_ += weight*pf0p0;
201 dpnorm_ += weight*coeff1;
202 dpnorm1_ += weight*coeff2;
203 pgv_ += weight*coeff1*gv;
204 pgv1_ += weight*coeff2*gv;
206 pHMdualVec0_->axpy(weight*coeff1,g);
207 pHMdualVec1_->axpy(weight*coeff2,g);
208 pHMdualVec2_->axpy(weight*coeff1,hv);
209 pHMdualVec3_->axpy(weight*coeff2*gv,g);
213 std::vector<Real> val_in(2,0.0), val_out(2,0.0);
216 sampler.
sumAll(&val_in[0],&val_out[0],2);
217 return (1.0-coeff_)*val_out[0]
218 + coeff_*(xvar_ + std::pow(val_out[1],1.0/(Real)order_)/(1.0-
prob_));
222 std::vector<Real> val_in(3,0.0), val_out(3,0.0);
226 sampler.
sumAll(&val_in[0],&val_out[0],3);
228 dualVector_->scale(1.0-coeff_);
231 if ( val_in[1] > 0. ) {
232 sampler.
sumAll(*pHMdualVec0_,*HMdualVec0_);
234 Real denom = (1.0-
prob_)*std::pow(val_out[1],((Real)order_-1.0)/(Real)order_);
236 var -= coeff_*((denom > 0.) ? val_out[2]/denom : 0.);
238 dualVector_->axpy(coeff_/denom,*HMdualVec0_);
246 std::vector<Real> val_in(6,0.0), val_out(6,0.0);
251 sampler.
sumAll(&val_in[0],&val_out[0],6);
255 dualVector_->scale(1.0-coeff_);
257 if ( val_out[1] > 0. ) {
258 sampler.
sumAll(*pHMdualVec0_,*HMdualVec0_);
259 sampler.
sumAll(*pHMdualVec1_,*HMdualVec1_);
260 sampler.
sumAll(*pHMdualVec2_,*HMdualVec2_);
261 sampler.
sumAll(*pHMdualVec3_,*HMdualVec3_);
263 Real rorder0 = (Real)order_;
264 Real rorder1 = (Real)order_-1.0;
265 Real rorder2 = (Real)(2*order_)-1.0;
267 Real denom1 = (1.0-
prob_)*std::pow(val_out[1],rorder1/rorder0);
268 Real denom2 = (1.0-
prob_)*std::pow(val_out[1],rorder2/rorder0);
270 var = coeff_*((val_out[3]/denom1 - rorder1*val_out[2]*val_out[2]/denom2)*vvar_
271 -(val_out[5]/denom1 - rorder1*val_out[4]*val_out[2]/denom2));
273 dualVector_->axpy(coeff_*(-vvar_/denom1),*HMdualVec1_);
274 dualVector_->axpy(coeff_*(vvar_*rorder1*val_out[2]/denom2),*HMdualVec0_);
275 dualVector_->axpy(coeff_/denom1,*HMdualVec3_);
276 dualVector_->axpy(coeff_/denom1,*HMdualVec2_);
277 dualVector_->axpy(coeff_*(-rorder1*val_out[4]/denom2),*HMdualVec0_);
HMCR(Real prob, Real coeff, unsigned order, Teuchos::RCP< PlusFunction< Real > > &pf)
Teuchos::RCP< PlusFunction< Real > > plusFunction_
Teuchos::RCP< Vector< Real > > dualVector_
void reset(Teuchos::RCP< Vector< Real > > &x0, const Vector< Real > &x)
Teuchos::RCP< Vector< Real > > HMdualVec1_
Teuchos::RCP< Vector< Real > > pHMdualVec2_
const Real getStatistic(const int i=0) const
HMCR(Teuchos::ParameterList &parlist)
void getGradient(Vector< Real > &g, SampleGenerator< Real > &sampler)
Defines the linear algebra or vector space interface.
void sumAll(Real *input, Real *output, int dim) const
Teuchos::RCP< const Vector< Real > > getVector() const
Teuchos::RCP< Vector< Real > > HMdualVec2_
Teuchos::RCP< Vector< Real > > HMdualVec0_
Teuchos::RCP< Vector< Real > > HMdualVec3_
Teuchos::RCP< Vector< Real > > pHMdualVec1_
void reset(Teuchos::RCP< Vector< Real > > &x0, const Vector< Real > &x, Teuchos::RCP< Vector< Real > > &v0, const Vector< Real > &v)
void update(const Real val, const Real weight)
void update(const Real val, const Vector< Real > &g, const Real weight)
void update(const Real val, const Vector< Real > &g, const Real gv, const Vector< Real > &hv, const Real weight)
Teuchos::RCP< Vector< Real > > pHMdualVec0_
virtual void update(const Real val, const Real weight)
Real getValue(SampleGenerator< Real > &sampler)
virtual void reset(Teuchos::RCP< Vector< Real > > &x0, const Vector< Real > &x)
void getHessVec(Vector< Real > &hv, SampleGenerator< Real > &sampler)
Teuchos::RCP< Vector< Real > > pHMdualVec3_