ROL
ROL_SimulatedEqualityConstraint.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Rapid Optimization Library (ROL) Package
5 // Copyright (2014) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact lead developers:
38 // Drew Kouri (dpkouri@sandia.gov) and
39 // Denis Ridzal (dridzal@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
44 #ifndef ROL_SIMULATED_EQUALITY_CONSTRAINT_H
45 #define ROL_SIMULATED_EQUALITY_CONSTRAINT_H
46 
47 #include "ROL_SimulatedVector.hpp"
48 #include "ROL_RiskVector.hpp"
50 
51 namespace ROL {
52 
53 template <class Real>
55 private:
56  const Teuchos::RCP<SampleGenerator<Real> > sampler_;
57  const Teuchos::RCP<EqualityConstraint_SimOpt<Real> > pcon_;
58  const bool useWeights_;
59 
60 public:
61 
63 
64  SimulatedEqualityConstraint(const Teuchos::RCP<SampleGenerator<Real> > & sampler,
65  const Teuchos::RCP<EqualityConstraint_SimOpt<Real> > & pcon,
66  const bool useWeights = true)
67  : sampler_(sampler), pcon_(pcon), useWeights_(useWeights) {}
68 
70  const Vector<Real> &x,
71  Real &tol) {
72  c.zero();
73  SimulatedVector<Real> &pc = Teuchos::dyn_cast<SimulatedVector<Real> >(c);
74  const Vector_SimOpt<Real> &uz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(x);
75  Teuchos::RCP<const Vector<Real> > uptr = uz.get_1();
76  Teuchos::RCP<const Vector<Real> > zptr = uz.get_2();
77  try {
78  const RiskVector<Real> &rz = Teuchos::dyn_cast<const RiskVector<Real> >(*zptr);
79  zptr = rz.getVector();
80  }
81  catch (const std::bad_cast &e) {}
82  const SimulatedVector<Real> &pu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*uptr);
83 
84  std::vector<Real> param;
85  Real weight(0), one(1);
86  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pu.numVectors(); ++i) {
87  param = sampler_->getMyPoint(static_cast<int>(i));
88  weight = sampler_->getMyWeight(static_cast<int>(i));
89  pcon_->setParameter(param);
90  pcon_->value(*(pc.get(i)), *(pu.get(i)), *zptr, tol);
91  weight = (useWeights_) ? weight : one;
92  pc.get(i)->scale(weight);
93  }
94 
95  }
96 
97 
98  virtual void applyJacobian(Vector<Real> &jv,
99  const Vector<Real> &v,
100  const Vector<Real> &x,
101  Real &tol) {
102  jv.zero();
103  // cast jv
104  SimulatedVector<Real> &pjv = Teuchos::dyn_cast<SimulatedVector<Real> >(jv);
105  // split x
106  const Vector_SimOpt<Real> &xuz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(x);
107  Teuchos::RCP<const Vector<Real> > xuptr = xuz.get_1();
108  Teuchos::RCP<const Vector<Real> > xzptr = xuz.get_2();
109  try {
110  const RiskVector<Real> &rxz = Teuchos::dyn_cast<const RiskVector<Real> >(*xzptr);
111  xzptr = rxz.getVector();
112  }
113  catch (const std::bad_cast &e) {}
114  const SimulatedVector<Real> &pxu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*xuptr);
115  // split v
116  const Vector_SimOpt<Real> &vuz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(v);
117  Teuchos::RCP<const Vector<Real> > vuptr = vuz.get_1();
118  Teuchos::RCP<const Vector<Real> > vzptr = vuz.get_2();
119  try {
120  const RiskVector<Real> &rvz = Teuchos::dyn_cast<const RiskVector<Real> >(*vzptr);
121  vzptr = rvz.getVector();
122  }
123  catch (const std::bad_cast &e) {}
124  const SimulatedVector<Real> &pvu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*vuptr);
125 
126  std::vector<Real> param;
127  Real weight(0), one(1);
128  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pvu.numVectors(); ++i) {
129  param = sampler_->getMyPoint(static_cast<int>(i));
130  weight = sampler_->getMyWeight(static_cast<int>(i));
131  pcon_->setParameter(param);
132  Vector_SimOpt<Real> xi(Teuchos::rcp_const_cast<Vector<Real> >(pxu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(xzptr));
133  Vector_SimOpt<Real> vi(Teuchos::rcp_const_cast<Vector<Real> >(pvu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(vzptr));
134  pcon_->applyJacobian(*(pjv.get(i)), vi, xi, tol);
135  weight = (useWeights_) ? weight : one;
136  pjv.get(i)->scale(weight);
137  }
138  }
139 
140 
142  const Vector<Real> &v,
143  const Vector<Real> &x,
144  Real &tol) {
145  ajv.zero();
146  // split ajv
147  Vector_SimOpt<Real> &ajvuz = Teuchos::dyn_cast<Vector_SimOpt<Real> >(ajv);
148  Teuchos::RCP<Vector<Real> > ajvuptr = ajvuz.get_1();
149  Teuchos::RCP<Vector<Real> > ajvzptr = ajvuz.get_2();
150  try {
151  RiskVector<Real> &rajvz = Teuchos::dyn_cast<RiskVector<Real> >(*ajvzptr);
152  ajvzptr = rajvz.getVector();
153  }
154  catch (const std::bad_cast &e) {}
155  SimulatedVector<Real> &pajvu = Teuchos::dyn_cast<SimulatedVector<Real> >(*ajvuptr);
156  // split x
157  const Vector_SimOpt<Real> &xuz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(x);
158  Teuchos::RCP<const Vector<Real> > xuptr = xuz.get_1();
159  Teuchos::RCP<const Vector<Real> > xzptr = xuz.get_2();
160  try {
161  const RiskVector<Real> &rxz = Teuchos::dyn_cast<const RiskVector<Real> >(*xzptr);
162  xzptr = rxz.getVector();
163  }
164  catch (const std::bad_cast &e) {}
165  const SimulatedVector<Real> &pxu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*xuptr);
166  // cast v
167  const SimulatedVector<Real> &pv = Teuchos::dyn_cast<const SimulatedVector<Real> >(v);
168 
169  std::vector<Real> param;
170  Real weight(0), one(1);
171  Teuchos::RCP<Vector<Real> > tmp1 = ajvzptr->clone();
172  Teuchos::RCP<Vector<Real> > tmp2 = ajvzptr->clone();
173  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pv.numVectors(); ++i) {
174  param = sampler_->getMyPoint(static_cast<int>(i));
175  weight = sampler_->getMyWeight(static_cast<int>(i));
176  pcon_->setParameter(param);
177  Vector_SimOpt<Real> xi(Teuchos::rcp_const_cast<Vector<Real> >(pxu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(xzptr));
178  Vector_SimOpt<Real> ajvi(pajvu.get(i), tmp1);
179  pcon_->applyAdjointJacobian(ajvi, *(pv.get(i)), xi, tol);
180  weight = (useWeights_) ? weight : one;
181  ajvi.scale(weight);
182  tmp2->plus(*tmp1);
183  }
184  sampler_->sumAll(*tmp2, *ajvzptr);
185 
186  }
187 
188 
189  virtual void applyAdjointHessian(Vector<Real> &ahuv,
190  const Vector<Real> &u,
191  const Vector<Real> &v,
192  const Vector<Real> &x,
193  Real &tol) {
194  ahuv.zero();
195  // split ahuv
196  Vector_SimOpt<Real> &ahuvuz = Teuchos::dyn_cast<Vector_SimOpt<Real> >(ahuv);
197  Teuchos::RCP<Vector<Real> > ahuvuptr = ahuvuz.get_1();
198  Teuchos::RCP<Vector<Real> > ahuvzptr = ahuvuz.get_2();
199  try {
200  RiskVector<Real> &rahuvz = Teuchos::dyn_cast<RiskVector<Real> >(*ahuvzptr);
201  ahuvzptr = rahuvz.getVector();
202  }
203  catch (const std::bad_cast &e) {}
204  SimulatedVector<Real> &pahuvu = Teuchos::dyn_cast<SimulatedVector<Real> >(*ahuvuptr);
205  // cast u
206  const SimulatedVector<Real> &pu = Teuchos::dyn_cast<const SimulatedVector<Real> >(u);
207  // split v
208  const Vector_SimOpt<Real> &vuz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(v);
209  Teuchos::RCP<const Vector<Real> > vuptr = vuz.get_1();
210  Teuchos::RCP<const Vector<Real> > vzptr = vuz.get_2();
211  try {
212  const RiskVector<Real> &rvz = Teuchos::dyn_cast<const RiskVector<Real> >(*vzptr);
213  vzptr = rvz.getVector();
214  }
215  catch (const std::bad_cast &e) {}
216  const SimulatedVector<Real> &pvu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*vuptr);
217  // split x
218  const Vector_SimOpt<Real> &xuz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(x);
219  Teuchos::RCP<const Vector<Real> > xuptr = xuz.get_1();
220  Teuchos::RCP<const Vector<Real> > xzptr = xuz.get_2();
221  try {
222  const RiskVector<Real> &rxz = Teuchos::dyn_cast<const RiskVector<Real> >(*xzptr);
223  xzptr = rxz.getVector();
224  }
225  catch (const std::bad_cast &e) {}
226  const SimulatedVector<Real> &pxu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*xuptr);
227 
228  std::vector<Real> param;
229  Real weight(0), one(1);
230  Teuchos::RCP<Vector<Real> > tmp1 = ahuvzptr->clone();
231  Teuchos::RCP<Vector<Real> > tmp2 = ahuvzptr->clone();
232  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pxu.numVectors(); ++i) {
233  param = sampler_->getMyPoint(static_cast<int>(i));
234  weight = sampler_->getMyWeight(static_cast<int>(i));
235  pcon_->setParameter(param);
236  Vector_SimOpt<Real> xi(Teuchos::rcp_const_cast<Vector<Real> >(pxu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(xzptr));
237  Vector_SimOpt<Real> vi(Teuchos::rcp_const_cast<Vector<Real> >(pvu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(vzptr));
238  Vector_SimOpt<Real> ahuvi(pahuvu.get(i), tmp1);
239  pcon_->applyAdjointHessian(ahuvi, *(pu.get(i)), vi, xi, tol);
240  weight = (useWeights_) ? weight : one;
241  ahuvi.scale(weight);
242  tmp2->plus(*tmp1);
243  }
244  sampler_->sumAll(*tmp2, *ahuvzptr);
245 
246  }
247 
249  const Vector<Real> &v,
250  const Vector<Real> &x,
251  const Vector<Real> &g,
252  Real &tol) {
253  Pv.zero();
254  // cast Pv
255  SimulatedVector<Real> &ppv = Teuchos::dyn_cast<SimulatedVector<Real> >(Pv);
256  // split x
257  const Vector_SimOpt<Real> &xuz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(x);
258  Teuchos::RCP<const Vector<Real> > xuptr = xuz.get_1();
259  Teuchos::RCP<const Vector<Real> > xzptr = xuz.get_2();
260  try {
261  const RiskVector<Real> &rxz = Teuchos::dyn_cast<const RiskVector<Real> >(*xzptr);
262  xzptr = rxz.getVector();
263  }
264  catch (const std::bad_cast &e) {}
265  const SimulatedVector<Real> &pxu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*xuptr);
266  // split g
267  const Vector_SimOpt<Real> &guz = Teuchos::dyn_cast<const Vector_SimOpt<Real> >(g);
268  Teuchos::RCP<const Vector<Real> > guptr = guz.get_1();
269  Teuchos::RCP<const Vector<Real> > gzptr = guz.get_2();
270  try {
271  const RiskVector<Real> &rgz = Teuchos::dyn_cast<const RiskVector<Real> >(*gzptr);
272  gzptr = rgz.getVector();
273  }
274  catch (const std::bad_cast &e) {}
275  const SimulatedVector<Real> &pgu = Teuchos::dyn_cast<const SimulatedVector<Real> >(*guptr);
276  // cast v
277  const SimulatedVector<Real> &pv = Teuchos::dyn_cast<const SimulatedVector<Real> >(v);
278 
279  std::vector<Real> param;
280  Real weight(0), one(1);
281  for (typename std::vector<SimulatedVector<Real> >::size_type i=0; i<pv.numVectors(); ++i) {
282  param = sampler_->getMyPoint(static_cast<int>(i));
283  weight = sampler_->getMyWeight(static_cast<int>(i));
284  pcon_->setParameter(param);
285  Vector_SimOpt<Real> xi(Teuchos::rcp_const_cast<Vector<Real> >(pxu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(xzptr));
286  Vector_SimOpt<Real> gi(Teuchos::rcp_const_cast<Vector<Real> >(pgu.get(i)), Teuchos::rcp_const_cast<Vector<Real> >(gzptr));
287  pcon_->applyPreconditioner(*(ppv.get(i)), *(pv.get(i)), xi, gi, tol);
288  weight = (useWeights_) ? weight : one;
289  ppv.get(i)->scale(one/(weight*weight));
290  }
291 
292  }
293 
294 
295 }; // class SimulatedEqualityConstraint
296 
297 } // namespace ROL
298 
299 #endif
void value(Vector< Real > &c, const Vector< Real > &x, Real &tol)
Evaluate the constraint operator at .
const Teuchos::RCP< EqualityConstraint_SimOpt< Real > > pcon_
const Teuchos::RCP< SampleGenerator< Real > > sampler_
Defines the linear algebra or vector space interface for simulation-based optimization.
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
Teuchos::RCP< const Vector< Real > > get_2() const
virtual void zero()
Set to zero vector.
Definition: ROL_Vector.hpp:157
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:74
Defines the equality constraint operator interface for simulation-based optimization.
Teuchos::RCP< const Vector< Real > > getVector(void) const
Defines the linear algebra of a vector space on a generic partitioned vector where the individual vec...
virtual void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the constraint Jacobian at , , to vector .
Defines the equality constraint operator interface.
Teuchos::RCP< const Vector< Real > > get_1() 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 .
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 ...
SimulatedEqualityConstraint(const Teuchos::RCP< SampleGenerator< Real > > &sampler, const Teuchos::RCP< EqualityConstraint_SimOpt< Real > > &pcon, const bool useWeights=true)
Teuchos::RCP< const Vector< Real > > get(size_type i) const