10 #ifndef ROL_SIMULATED_CONSTRAINT_H
11 #define ROL_SIMULATED_CONSTRAINT_H
22 const ROL::Ptr<SampleGenerator<Real> >
sampler_;
23 const ROL::Ptr<Constraint_SimOpt<Real> >
pcon_;
32 const bool useWeights =
true)
36 pcon_->update(x,flag,iter);
39 pcon_->update(x,type,iter);
48 ROL::Ptr<const Vector<Real> > uptr = uz.get_1();
49 ROL::Ptr<const Vector<Real> > zptr = uz.get_2();
54 catch (
const std::bad_cast &e) {}
57 std::vector<Real> param;
58 Real weight(0), one(1);
60 param =
sampler_->getMyPoint(static_cast<int>(i));
61 weight =
sampler_->getMyWeight(static_cast<int>(i));
62 pcon_->setParameter(param);
64 pcon_->value(*(pc.
get(i)), *(pu.
get(i)), *zptr, tol);
66 pc.
get(i)->scale(weight);
81 ROL::Ptr<const Vector<Real> > xuptr = xuz.get_1();
82 ROL::Ptr<const Vector<Real> > xzptr = xuz.get_2();
87 catch (
const std::bad_cast &e) {}
91 ROL::Ptr<const Vector<Real> > vuptr = vuz.
get_1();
92 ROL::Ptr<const Vector<Real> > vzptr = vuz.
get_2();
97 catch (
const std::bad_cast &e) {}
100 std::vector<Real> param;
101 Real weight(0), one(1);
103 param =
sampler_->getMyPoint(static_cast<int>(i));
104 weight =
sampler_->getMyWeight(static_cast<int>(i));
105 pcon_->setParameter(param);
109 pcon_->applyJacobian(*(pjv.
get(i)), vi, xi, tol);
111 pjv.
get(i)->scale(weight);
123 ROL::Ptr<Vector<Real> > ajvuptr = ajvuz.
get_1();
124 ROL::Ptr<Vector<Real> > ajvzptr = ajvuz.
get_2();
129 catch (
const std::bad_cast &e) {}
133 ROL::Ptr<const Vector<Real> > xuptr = xuz.
get_1();
134 ROL::Ptr<const Vector<Real> > xzptr = xuz.
get_2();
139 catch (
const std::bad_cast &e) {}
144 std::vector<Real> param;
145 Real weight(0), one(1);
146 ROL::Ptr<Vector<Real> > tmp1 = ajvzptr->clone();
147 ROL::Ptr<Vector<Real> > tmp2 = ajvzptr->clone();
149 param =
sampler_->getMyPoint(static_cast<int>(i));
150 weight =
sampler_->getMyWeight(static_cast<int>(i));
151 pcon_->setParameter(param);
155 pcon_->applyAdjointJacobian(ajvi, *(pv.
get(i)), xi, tol);
173 ROL::Ptr<Vector<Real> > ahuvuptr = ahuvuz.
get_1();
174 ROL::Ptr<Vector<Real> > ahuvzptr = ahuvuz.
get_2();
179 catch (
const std::bad_cast &e) {}
185 ROL::Ptr<const Vector<Real> > vuptr = vuz.get_1();
186 ROL::Ptr<const Vector<Real> > vzptr = vuz.get_2();
191 catch (
const std::bad_cast &e) {}
195 ROL::Ptr<const Vector<Real> > xuptr = xuz.
get_1();
196 ROL::Ptr<const Vector<Real> > xzptr = xuz.
get_2();
201 catch (
const std::bad_cast &e) {}
204 std::vector<Real> param;
205 Real weight(0), one(1);
206 ROL::Ptr<Vector<Real> > tmp1 = ahuvzptr->clone();
207 ROL::Ptr<Vector<Real> > tmp2 = ahuvzptr->clone();
209 param =
sampler_->getMyPoint(static_cast<int>(i));
210 weight =
sampler_->getMyWeight(static_cast<int>(i));
211 pcon_->setParameter(param);
216 pcon_->applyAdjointHessian(ahuvi, *(pu.
get(i)), vi, xi, tol);
235 ROL::Ptr<const Vector<Real> > xuptr = xuz.get_1();
236 ROL::Ptr<const Vector<Real> > xzptr = xuz.get_2();
241 catch (
const std::bad_cast &e) {}
245 ROL::Ptr<const Vector<Real> > guptr = guz.
get_1();
246 ROL::Ptr<const Vector<Real> > gzptr = guz.
get_2();
251 catch (
const std::bad_cast &e) {}
256 std::vector<Real> param;
257 Real weight(0), one(1);
259 param =
sampler_->getMyPoint(static_cast<int>(i));
260 weight =
sampler_->getMyWeight(static_cast<int>(i));
261 pcon_->setParameter(param);
265 pcon_->applyPreconditioner(*(ppv.
get(i)), *(pv.
get(i)), xi, gi, tol);
267 ppv.
get(i)->scale(one/(weight*weight));
typename PV< Real >::size_type size_type
virtual void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the constraint Jacobian at , , to vector .
ROL::Ptr< const Vector< Real > > get_2() const
Defines the linear algebra or vector space interface for simulation-based optimization.
virtual void applyAdjointHessian(Vector< Real > &ahuv, const Vector< Real > &u, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ...
const ROL::Ptr< Constraint_SimOpt< Real > > pcon_
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
const ROL::Ptr< SampleGenerator< Real > > sampler_
Defines the linear algebra of a vector space on a generic partitioned vector where the individual vec...
ROL::Ptr< const Vector< Real > > get(size_type i) const
virtual void applyPreconditioner(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, const Vector< Real > &g, Real &tol)
Apply a constraint preconditioner at , , to vector . Ideally, this preconditioner satisfies the follo...
size_type numVectors() const
virtual void applyAdjointJacobian(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the adjoint of the the constraint Jacobian at , , to vector .
SimulatedConstraint(const ROL::Ptr< SampleGenerator< Real > > &sampler, const ROL::Ptr< Constraint_SimOpt< Real > > &pcon, const bool useWeights=true)
void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update constraint function.
void value(Vector< Real > &c, const Vector< Real > &x, Real &tol)
Evaluate the constraint operator at .
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update constraint functions. x is the optimization variable, flag = true if optimization variable is ...
Ptr< const Vector< Real > > getVector(void) const
virtual ~SimulatedConstraint()
Defines the constraint operator interface for simulation-based optimization.
Defines the general constraint operator interface.
ROL::Ptr< const Vector< Real > > get_1() const