44 #ifndef ROL_RISKNEUTRALOBJECTIVE_HPP 45 #define ROL_RISKNEUTRALOBJECTIVE_HPP 47 #include "Teuchos_RefCountPtr.hpp" 75 const std::vector<Real> ¶m, Real &tol) {
76 if ( storage_ && value_storage_.count(param) ) {
77 val = value_storage_[param];
80 ParametrizedObjective_->setParameter(param);
81 val = ParametrizedObjective_->value(x,tol);
83 value_storage_.insert(std::pair<std::vector<Real>,Real>(param,val));
89 const std::vector<Real> ¶m, Real &tol) {
90 if ( storage_ && gradient_storage_.count(param) ) {
91 g.
set(*(gradient_storage_[param]));
94 ParametrizedObjective_->setParameter(param);
95 ParametrizedObjective_->gradient(g,x,tol);
97 Teuchos::RCP<Vector<Real> > tmp = g.
clone();
98 gradient_storage_.insert(std::pair<std::vector<Real>,Teuchos::RCP<
Vector<Real> > >(param,tmp));
99 gradient_storage_[param]->set(g);
105 const std::vector<Real> ¶m, Real &tol) {
106 ParametrizedObjective_->setParameter(param);
107 ParametrizedObjective_->hessVec(hv,v,x,tol);
118 const bool storage =
true )
119 : ParametrizedObjective_(pObj),
120 ValueSampler_(vsampler), GradientSampler_(gsampler), HessianSampler_(hsampler),
121 firstUpdate_(true), storage_(true) {
122 value_storage_.clear();
123 gradient_storage_.clear();
129 const bool storage =
true )
130 : ParametrizedObjective_(pObj),
131 ValueSampler_(vsampler), GradientSampler_(gsampler), HessianSampler_(gsampler),
132 firstUpdate_(true), storage_(true) {
133 value_storage_.clear();
134 gradient_storage_.clear();
139 const bool storage =
true )
140 : ParametrizedObjective_(pObj),
141 ValueSampler_(sampler), GradientSampler_(sampler), HessianSampler_(sampler),
142 firstUpdate_(true), storage_(true) {
143 value_storage_.clear();
144 gradient_storage_.clear();
148 if ( firstUpdate_ ) {
149 gradient_ = (x.
dual()).clone();
150 pointDual_ = (x.
dual()).clone();
151 sumDual_ = (x.
dual()).clone();
152 firstUpdate_ =
false;
154 ParametrizedObjective_->update(x,flag,iter);
155 ValueSampler_->update(x);
158 value_storage_.clear();
161 GradientSampler_->update(x);
162 HessianSampler_->update(x);
165 gradient_storage_.clear();
171 Real myval = 0.0, ptval = 0.0, val = 0.0, error = 2.0*tol + 1.0;
172 std::vector<Real> ptvals;
173 while ( error > tol ) {
174 ValueSampler_->refine();
175 for (
int i = ValueSampler_->start(); i < ValueSampler_->numMySamples(); i++ ) {
176 getValue(ptval,x,ValueSampler_->getMyPoint(i),tol);
177 myval += ValueSampler_->getMyWeight(i)*ptval;
178 ptvals.push_back(ptval);
180 error = ValueSampler_->computeError(ptvals);
183 ValueSampler_->sumAll(&myval,&val,1);
185 ValueSampler_->setSamples();
190 g.
zero(); pointDual_->zero(); sumDual_->zero();
191 std::vector<Teuchos::RCP<Vector<Real> > > ptgs;
192 Real error = 2.0*tol + 1.0;
193 while ( error > tol ) {
194 GradientSampler_->refine();
195 for (
int i = GradientSampler_->start(); i < GradientSampler_->numMySamples(); i++ ) {
196 getGradient(*pointDual_,x,GradientSampler_->getMyPoint(i),tol);
197 sumDual_->axpy(GradientSampler_->getMyWeight(i),*
pointDual_);
198 ptgs.push_back(pointDual_->clone());
199 (ptgs.back())->
set(*pointDual_);
201 error = GradientSampler_->computeError(ptgs,x);
204 GradientSampler_->sumAll(*sumDual_,g);
205 gradient_->axpy(1.0,g);
207 GradientSampler_->setSamples();
212 hv.
zero(); pointDual_->zero(); sumDual_->zero();
213 for (
int i = 0; i < HessianSampler_->numMySamples(); i++ ) {
214 getHessVec(*pointDual_,v,x,HessianSampler_->getMyPoint(i),tol);
215 sumDual_->axpy(HessianSampler_->getMyWeight(i),*
pointDual_);
217 HessianSampler_->sumAll(*sumDual_,hv);
Provides the interface to evaluate objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
RiskNeutralObjective(const Teuchos::RCP< ParametrizedObjective< Real > > &pObj, const Teuchos::RCP< SampleGenerator< Real > > &vsampler, const Teuchos::RCP< SampleGenerator< Real > > &gsampler, const Teuchos::RCP< SampleGenerator< Real > > &hsampler, const bool storage=true)
Teuchos::RCP< SampleGenerator< Real > > GradientSampler_
Teuchos::RCP< Vector< Real > > sumDual_
Teuchos::RCP< SampleGenerator< Real > > ValueSampler_
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
void getHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
Teuchos::RCP< Vector< Real > > gradient_
virtual ~RiskNeutralObjective()
RiskNeutralObjective(const Teuchos::RCP< ParametrizedObjective< Real > > &pObj, const Teuchos::RCP< SampleGenerator< Real > > &vsampler, const Teuchos::RCP< SampleGenerator< Real > > &gsampler, const bool storage=true)
RiskNeutralObjective(const Teuchos::RCP< ParametrizedObjective< Real > > &pObj, const Teuchos::RCP< SampleGenerator< Real > > &sampler, const bool storage=true)
Teuchos::RCP< Vector< Real > > pointDual_
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
void getValue(Real &val, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply preconditioner to vector.
std::map< std::vector< Real >, Real > value_storage_
std::map< std::vector< Real >, Teuchos::RCP< Vector< Real > > > gradient_storage_
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
void getGradient(Vector< Real > &g, const Vector< Real > &x, const std::vector< Real > ¶m, Real &tol)
virtual Real value(const Vector< Real > &x, Real &tol)
Compute value.
virtual void set(const Vector &x)
Set where .
Teuchos::RCP< ParametrizedObjective< Real > > ParametrizedObjective_
Teuchos::RCP< SampleGenerator< Real > > HessianSampler_