ROL
ROL_SimulatedConstraint.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef ROL_SIMULATED_CONSTRAINT_H
11 #define ROL_SIMULATED_CONSTRAINT_H
12 
13 #include "ROL_SimulatedVector.hpp"
14 #include "ROL_RiskVector.hpp"
16 
17 namespace ROL {
18 
19 template <class Real>
20 class SimulatedConstraint : public Constraint<Real> {
21 private:
22  const ROL::Ptr<SampleGenerator<Real> > sampler_;
23  const ROL::Ptr<Constraint_SimOpt<Real> > pcon_;
24  const bool useWeights_;
25 
26 public:
27 
28  virtual ~SimulatedConstraint() {}
29 
30  SimulatedConstraint(const ROL::Ptr<SampleGenerator<Real> > & sampler,
31  const ROL::Ptr<Constraint_SimOpt<Real> > & pcon,
32  const bool useWeights = true)
33  : sampler_(sampler), pcon_(pcon), useWeights_(useWeights) {}
34 
35  void update( const Vector<Real> &x, bool flag = true, int iter = -1 ) {
36  pcon_->update(x,flag,iter);
37  }
38  void update( const Vector<Real> &x, UpdateType type, int iter = -1 ) {
39  pcon_->update(x,type,iter);
40  }
41 
43  const Vector<Real> &x,
44  Real &tol) {
45  c.zero();
46  SimulatedVector<Real> &pc = dynamic_cast<SimulatedVector<Real>&>(c);
47  const Vector_SimOpt<Real> &uz = dynamic_cast<const Vector_SimOpt<Real>&>(x);
48  ROL::Ptr<const Vector<Real> > uptr = uz.get_1();
49  ROL::Ptr<const Vector<Real> > zptr = uz.get_2();
50  try {
51  const RiskVector<Real> &rz = dynamic_cast<const RiskVector<Real>&>(*zptr);
52  zptr = rz.getVector();
53  }
54  catch (const std::bad_cast &e) {}
55  const SimulatedVector<Real> &pu = dynamic_cast<const SimulatedVector<Real>&>(*uptr);
56 
57  std::vector<Real> param;
58  Real weight(0), one(1);
59  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pu.numVectors(); ++i) {
60  param = sampler_->getMyPoint(static_cast<int>(i));
61  weight = sampler_->getMyWeight(static_cast<int>(i));
62  pcon_->setParameter(param);
63  pcon_->update(*(pu.get(i)), *zptr);
64  pcon_->value(*(pc.get(i)), *(pu.get(i)), *zptr, tol);
65  weight = (useWeights_) ? weight : one;
66  pc.get(i)->scale(weight);
67  }
68 
69  }
70 
71 
72  virtual void applyJacobian(Vector<Real> &jv,
73  const Vector<Real> &v,
74  const Vector<Real> &x,
75  Real &tol) {
76  jv.zero();
77  // cast jv
78  SimulatedVector<Real> &pjv = dynamic_cast<SimulatedVector<Real>&>(jv);
79  // split x
80  const Vector_SimOpt<Real> &xuz = dynamic_cast<const Vector_SimOpt<Real>&>(x);
81  ROL::Ptr<const Vector<Real> > xuptr = xuz.get_1();
82  ROL::Ptr<const Vector<Real> > xzptr = xuz.get_2();
83  try {
84  const RiskVector<Real> &rxz = dynamic_cast<const RiskVector<Real>&>(*xzptr);
85  xzptr = rxz.getVector();
86  }
87  catch (const std::bad_cast &e) {}
88  const SimulatedVector<Real> &pxu = dynamic_cast<const SimulatedVector<Real>&>(*xuptr);
89  // split v
90  const Vector_SimOpt<Real> &vuz = dynamic_cast<const Vector_SimOpt<Real>&>(v);
91  ROL::Ptr<const Vector<Real> > vuptr = vuz.get_1();
92  ROL::Ptr<const Vector<Real> > vzptr = vuz.get_2();
93  try {
94  const RiskVector<Real> &rvz = dynamic_cast<const RiskVector<Real>&>(*vzptr);
95  vzptr = rvz.getVector();
96  }
97  catch (const std::bad_cast &e) {}
98  const SimulatedVector<Real> &pvu = dynamic_cast<const SimulatedVector<Real>&>(*vuptr);
99 
100  std::vector<Real> param;
101  Real weight(0), one(1);
102  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pvu.numVectors(); ++i) {
103  param = sampler_->getMyPoint(static_cast<int>(i));
104  weight = sampler_->getMyWeight(static_cast<int>(i));
105  pcon_->setParameter(param);
106  Vector_SimOpt<Real> xi(ROL::constPtrCast<Vector<Real> >(pxu.get(i)), ROL::constPtrCast<Vector<Real> >(xzptr));
107  Vector_SimOpt<Real> vi(ROL::constPtrCast<Vector<Real> >(pvu.get(i)), ROL::constPtrCast<Vector<Real> >(vzptr));
108  pcon_->update(xi);
109  pcon_->applyJacobian(*(pjv.get(i)), vi, xi, tol);
110  weight = (useWeights_) ? weight : one;
111  pjv.get(i)->scale(weight);
112  }
113  }
114 
115 
117  const Vector<Real> &v,
118  const Vector<Real> &x,
119  Real &tol) {
120  ajv.zero();
121  // split ajv
122  Vector_SimOpt<Real> &ajvuz = dynamic_cast<Vector_SimOpt<Real>&>(ajv);
123  ROL::Ptr<Vector<Real> > ajvuptr = ajvuz.get_1();
124  ROL::Ptr<Vector<Real> > ajvzptr = ajvuz.get_2();
125  try {
126  RiskVector<Real> &rajvz = dynamic_cast<RiskVector<Real>&>(*ajvzptr);
127  ajvzptr = rajvz.getVector();
128  }
129  catch (const std::bad_cast &e) {}
130  SimulatedVector<Real> &pajvu = dynamic_cast<SimulatedVector<Real>&>(*ajvuptr);
131  // split x
132  const Vector_SimOpt<Real> &xuz = dynamic_cast<const Vector_SimOpt<Real>&>(x);
133  ROL::Ptr<const Vector<Real> > xuptr = xuz.get_1();
134  ROL::Ptr<const Vector<Real> > xzptr = xuz.get_2();
135  try {
136  const RiskVector<Real> &rxz = dynamic_cast<const RiskVector<Real>&>(*xzptr);
137  xzptr = rxz.getVector();
138  }
139  catch (const std::bad_cast &e) {}
140  const SimulatedVector<Real> &pxu = dynamic_cast<const SimulatedVector<Real>&>(*xuptr);
141  // cast v
142  const SimulatedVector<Real> &pv = dynamic_cast<const SimulatedVector<Real>&>(v);
143 
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();
148  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pv.numVectors(); ++i) {
149  param = sampler_->getMyPoint(static_cast<int>(i));
150  weight = sampler_->getMyWeight(static_cast<int>(i));
151  pcon_->setParameter(param);
152  Vector_SimOpt<Real> xi(ROL::constPtrCast<Vector<Real> >(pxu.get(i)), ROL::constPtrCast<Vector<Real> >(xzptr));
153  Vector_SimOpt<Real> ajvi(pajvu.get(i), tmp1);
154  pcon_->update(xi);
155  pcon_->applyAdjointJacobian(ajvi, *(pv.get(i)), xi, tol);
156  weight = (useWeights_) ? weight : one;
157  ajvi.scale(weight);
158  tmp2->plus(*tmp1);
159  }
160  sampler_->sumAll(*tmp2, *ajvzptr);
161 
162  }
163 
164 
165  virtual void applyAdjointHessian(Vector<Real> &ahuv,
166  const Vector<Real> &u,
167  const Vector<Real> &v,
168  const Vector<Real> &x,
169  Real &tol) {
170  ahuv.zero();
171  // split ahuv
172  Vector_SimOpt<Real> &ahuvuz = dynamic_cast<Vector_SimOpt<Real>&>(ahuv);
173  ROL::Ptr<Vector<Real> > ahuvuptr = ahuvuz.get_1();
174  ROL::Ptr<Vector<Real> > ahuvzptr = ahuvuz.get_2();
175  try {
176  RiskVector<Real> &rahuvz = dynamic_cast<RiskVector<Real>&>(*ahuvzptr);
177  ahuvzptr = rahuvz.getVector();
178  }
179  catch (const std::bad_cast &e) {}
180  SimulatedVector<Real> &pahuvu = dynamic_cast<SimulatedVector<Real>&>(*ahuvuptr);
181  // cast u
182  const SimulatedVector<Real> &pu = dynamic_cast<const SimulatedVector<Real>&>(u);
183  // split v
184  const Vector_SimOpt<Real> &vuz = dynamic_cast<const Vector_SimOpt<Real>&>(v);
185  ROL::Ptr<const Vector<Real> > vuptr = vuz.get_1();
186  ROL::Ptr<const Vector<Real> > vzptr = vuz.get_2();
187  try {
188  const RiskVector<Real> &rvz = dynamic_cast<const RiskVector<Real>&>(*vzptr);
189  vzptr = rvz.getVector();
190  }
191  catch (const std::bad_cast &e) {}
192  const SimulatedVector<Real> &pvu = dynamic_cast<const SimulatedVector<Real>&>(*vuptr);
193  // split x
194  const Vector_SimOpt<Real> &xuz = dynamic_cast<const Vector_SimOpt<Real>&>(x);
195  ROL::Ptr<const Vector<Real> > xuptr = xuz.get_1();
196  ROL::Ptr<const Vector<Real> > xzptr = xuz.get_2();
197  try {
198  const RiskVector<Real> &rxz = dynamic_cast<const RiskVector<Real>&>(*xzptr);
199  xzptr = rxz.getVector();
200  }
201  catch (const std::bad_cast &e) {}
202  const SimulatedVector<Real> &pxu = dynamic_cast<const SimulatedVector<Real>&>(*xuptr);
203 
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();
208  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pxu.numVectors(); ++i) {
209  param = sampler_->getMyPoint(static_cast<int>(i));
210  weight = sampler_->getMyWeight(static_cast<int>(i));
211  pcon_->setParameter(param);
212  Vector_SimOpt<Real> xi(ROL::constPtrCast<Vector<Real> >(pxu.get(i)), ROL::constPtrCast<Vector<Real> >(xzptr));
213  Vector_SimOpt<Real> vi(ROL::constPtrCast<Vector<Real> >(pvu.get(i)), ROL::constPtrCast<Vector<Real> >(vzptr));
214  Vector_SimOpt<Real> ahuvi(pahuvu.get(i), tmp1);
215  pcon_->update(xi);
216  pcon_->applyAdjointHessian(ahuvi, *(pu.get(i)), vi, xi, tol);
217  weight = (useWeights_) ? weight : one;
218  ahuvi.scale(weight);
219  tmp2->plus(*tmp1);
220  }
221  sampler_->sumAll(*tmp2, *ahuvzptr);
222 
223  }
224 
226  const Vector<Real> &v,
227  const Vector<Real> &x,
228  const Vector<Real> &g,
229  Real &tol) {
230  Pv.zero();
231  // cast Pv
232  SimulatedVector<Real> &ppv = dynamic_cast<SimulatedVector<Real>&>(Pv);
233  // split x
234  const Vector_SimOpt<Real> &xuz = dynamic_cast<const Vector_SimOpt<Real>&>(x);
235  ROL::Ptr<const Vector<Real> > xuptr = xuz.get_1();
236  ROL::Ptr<const Vector<Real> > xzptr = xuz.get_2();
237  try {
238  const RiskVector<Real> &rxz = dynamic_cast<const RiskVector<Real>&>(*xzptr);
239  xzptr = rxz.getVector();
240  }
241  catch (const std::bad_cast &e) {}
242  const SimulatedVector<Real> &pxu = dynamic_cast<const SimulatedVector<Real>&>(*xuptr);
243  // split g
244  const Vector_SimOpt<Real> &guz = dynamic_cast<const Vector_SimOpt<Real>&>(g);
245  ROL::Ptr<const Vector<Real> > guptr = guz.get_1();
246  ROL::Ptr<const Vector<Real> > gzptr = guz.get_2();
247  try {
248  const RiskVector<Real> &rgz = dynamic_cast<const RiskVector<Real>&>(*gzptr);
249  gzptr = rgz.getVector();
250  }
251  catch (const std::bad_cast &e) {}
252  const SimulatedVector<Real> &pgu = dynamic_cast<const SimulatedVector<Real>&>(*guptr);
253  // cast v
254  const SimulatedVector<Real> &pv = dynamic_cast<const SimulatedVector<Real>&>(v);
255 
256  std::vector<Real> param;
257  Real weight(0), one(1);
258  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pv.numVectors(); ++i) {
259  param = sampler_->getMyPoint(static_cast<int>(i));
260  weight = sampler_->getMyWeight(static_cast<int>(i));
261  pcon_->setParameter(param);
262  Vector_SimOpt<Real> xi(ROL::constPtrCast<Vector<Real> >(pxu.get(i)), ROL::constPtrCast<Vector<Real> >(xzptr));
263  Vector_SimOpt<Real> gi(ROL::constPtrCast<Vector<Real> >(pgu.get(i)), ROL::constPtrCast<Vector<Real> >(gzptr));
264  pcon_->update(xi);
265  pcon_->applyPreconditioner(*(ppv.get(i)), *(pv.get(i)), xi, gi, tol);
266  weight = (useWeights_) ? weight : one;
267  ppv.get(i)->scale(one/(weight*weight));
268  }
269 
270  }
271 
272 
273 }; // class SimulatedConstraint
274 
275 } // namespace ROL
276 
277 #endif
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.
Definition: ROL_Vector.hpp:133
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
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
Defines the constraint operator interface for simulation-based optimization.
Defines the general constraint operator interface.
ROL::Ptr< const Vector< Real > > get_1() const