ROL
ROL_Constraint_SimOpt.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_CONSTRAINT_SIMOPT_H
45#define ROL_CONSTRAINT_SIMOPT_H
46
47#include "ROL_Constraint.hpp"
48#include "ROL_Vector_SimOpt.hpp"
49#include "ROL_Types.hpp"
50#include <iostream>
51
52namespace ROL {
53
54template <class Real>
55class Constraint_SimOpt;
56
57}
58
59#include "ROL_NonlinearLeastSquaresObjective.hpp"
60#include "ROL_SimConstraint.hpp"
61#include "ROL_Objective_FSsolver.hpp"
62#include "ROL_TypeU_TrustRegionAlgorithm.hpp"
63#include "ROL_TypeE_AugmentedLagrangianAlgorithm.hpp"
64
105namespace ROL {
106
107template <class Real>
108class Constraint_SimOpt : public Constraint<Real> {
109private:
110 // Additional vector storage for solve
111 Ptr<Vector<Real>> unew_;
112 Ptr<Vector<Real>> jv_;
113
114 // Default parameters for solve (backtracking Newton)
115 const Real DEFAULT_atol_;
116 const Real DEFAULT_rtol_;
117 const Real DEFAULT_stol_;
118 const Real DEFAULT_factor_;
119 const Real DEFAULT_decr_;
120 const int DEFAULT_maxit_;
121 const bool DEFAULT_print_;
122 const bool DEFAULT_zero_;
123 const int DEFAULT_solverType_;
124
125 // User-set parameters for solve (backtracking Newton)
126
127protected:
128 Real atol_;
129 Real rtol_;
130 Real stol_;
131 Real factor_;
132 Real decr_;
133 int maxit_;
134 bool print_;
135 bool zero_;
136 int solverType_;
137
138 // Flag to initialize vector storage in solve
139 bool firstSolve_;
140
141public:
142 Constraint_SimOpt()
143 : Constraint<Real>(),
144 unew_(nullPtr), jv_(nullPtr),
145 DEFAULT_atol_(1.e-4*std::sqrt(ROL_EPSILON<Real>())),
146 DEFAULT_rtol_(1.e0),
147 DEFAULT_stol_(std::sqrt(ROL_EPSILON<Real>())),
148 DEFAULT_factor_(0.5),
149 DEFAULT_decr_(1.e-4),
150 DEFAULT_maxit_(500),
151 DEFAULT_print_(false),
152 DEFAULT_zero_(false),
153 DEFAULT_solverType_(0),
154 atol_(DEFAULT_atol_), rtol_(DEFAULT_rtol_), stol_(DEFAULT_stol_), factor_(DEFAULT_factor_),
155 decr_(DEFAULT_decr_), maxit_(DEFAULT_maxit_), print_(DEFAULT_print_), zero_(DEFAULT_zero_),
156 solverType_(DEFAULT_solverType_), firstSolve_(true) {}
157
163 virtual void update( const Vector<Real> &u, const Vector<Real> &z, bool flag = true, int iter = -1 ) {
164 update_1(u,flag,iter);
165 update_2(z,flag,iter);
166 }
167 virtual void update( const Vector<Real> &u, const Vector<Real> &z, UpdateType type, int iter = -1 ) {
168 update_1(u,type,iter);
169 update_2(z,type,iter);
170 }
171
177 virtual void update_1( const Vector<Real> &u, bool flag = true, int iter = -1 ) {}
178 virtual void update_1( const Vector<Real> &u, UpdateType type, int iter = -1 ) {}
179
185 virtual void update_2( const Vector<Real> &z, bool flag = true, int iter = -1 ) {}
186 virtual void update_2( const Vector<Real> &z, UpdateType type, int iter = -1 ) {}
187
194 virtual void solve_update( const Vector<Real> &u, const Vector<Real> &z, UpdateType type, int iter = -1) {}
195
209 virtual void value(Vector<Real> &c,
210 const Vector<Real> &u,
211 const Vector<Real> &z,
212 Real &tol) = 0;
213
225 virtual void solve(Vector<Real> &c,
226 Vector<Real> &u,
227 const Vector<Real> &z,
228 Real &tol) {
229 if ( zero_ ) u.zero();
230 Ptr<std::ostream> stream = makeStreamPtr(std::cout, print_);
231 solve_update(u,z,UpdateType::Initial,0);
232 value(c,u,z,tol);
233 Real cnorm = c.norm();
234 const Real ctol = std::min(atol_, rtol_*cnorm);
235 if (solverType_==0 || solverType_==3 || solverType_==4) {
236 if ( firstSolve_ ) {
237 unew_ = u.clone();
238 jv_ = u.clone();
239 firstSolve_ = false;
240 }
241 const Real one(1);
242 Real alpha(1), tmp(0);
243 int cnt = 0;
244 *stream << std::endl;
245 *stream << " Default Constraint_SimOpt::solve" << std::endl;
246 *stream << " ";
247 *stream << std::setw(6) << std::left << "iter";
248 *stream << std::setw(15) << std::left << "rnorm";
249 *stream << std::setw(15) << std::left << "alpha";
250 *stream << std::endl;
251 for (cnt = 0; cnt < maxit_; ++cnt) {
252 // Compute Newton step
253 applyInverseJacobian_1(*jv_,c,u,z,tol);
254 unew_->set(u);
255 unew_->axpy(-alpha, *jv_);
256 solve_update(*unew_,z,UpdateType::Trial);
257 value(c,*unew_,z,tol);
258 tmp = c.norm();
259 // Perform backtracking line search
260 while ( tmp > (one-decr_*alpha)*cnorm &&
261 alpha > stol_ ) {
262 alpha *= factor_;
263 unew_->set(u);
264 unew_->axpy(-alpha,*jv_);
265 solve_update(*unew_,z,UpdateType::Trial);
266 value(c,*unew_,z,tol);
267 tmp = c.norm();
268 }
269 *stream << " ";
270 *stream << std::setw(6) << std::left << cnt;
271 *stream << std::scientific << std::setprecision(6);
272 *stream << std::setw(15) << std::left << tmp;
273 *stream << std::scientific << std::setprecision(6);
274 *stream << std::setw(15) << std::left << alpha;
275 *stream << std::endl;
276 // Update iterate
277 cnorm = tmp;
278 u.set(*unew_);
279 solve_update(u,z,UpdateType::Accept,cnt);
280 if (cnorm <= ctol) break; // = covers the case of identically zero residual
281 alpha = one;
282 }
283 }
284 if (solverType_==1 || (solverType_==3 && cnorm > ctol)) {
285 Ptr<Constraint<Real>> con = makePtr<SimConstraint<Real>>(makePtrFromRef(*this),makePtrFromRef(z),true);
286 Ptr<Objective<Real>> obj = makePtr<NonlinearLeastSquaresObjective<Real>>(con,u,c,true);
287 ParameterList parlist;
288 parlist.sublist("General").set("Output Level",1);
289 parlist.sublist("General").sublist("Krylov").set("Iteration Limit",100);
290 parlist.sublist("Step").sublist("Trust Region").set("Subproblem Solver","Truncated CG");
291 parlist.sublist("Status Test").set("Gradient Tolerance",ctol);
292 parlist.sublist("Status Test").set("Step Tolerance",stol_);
293 parlist.sublist("Status Test").set("Iteration Limit",maxit_);
294 Ptr<TypeU::Algorithm<Real>> algo = makePtr<TypeU::TrustRegionAlgorithm<Real>>(parlist);
295 algo->run(u,*obj,*stream);
296 value(c,u,z,tol);
297 }
298 if (solverType_==2 || (solverType_==4 && cnorm > ctol)) {
299 Ptr<Constraint<Real>> con = makePtr<SimConstraint<Real>>(makePtrFromRef(*this),makePtrFromRef(z),true);
300 Ptr<Objective<Real>> obj = makePtr<Objective_FSsolver<Real>>();
301 Ptr<Vector<Real>> l = c.dual().clone();
302 ParameterList parlist;
303 parlist.sublist("General").set("Output Level",1);
304 parlist.sublist("Status Test").set("Constraint Tolerance",ctol);
305 parlist.sublist("Status Test").set("Step Tolerance",stol_);
306 parlist.sublist("Status Test").set("Iteration Limit",maxit_);
307 Ptr<TypeE::Algorithm<Real>> algo = makePtr<TypeE::AugmentedLagrangianAlgorithm<Real>>(parlist);
308 algo->run(u,*obj,*con,*l,*stream);
309 value(c,u,z,tol);
310 }
311 if (solverType_ > 4 || solverType_ < 0) {
312 ROL_TEST_FOR_EXCEPTION(true, std::invalid_argument,
313 ">>> ERROR (ROL:Constraint_SimOpt:solve): Invalid solver type!");
314 }
315 }
316
337 virtual void setSolveParameters(ParameterList &parlist) {
338 ParameterList & list = parlist.sublist("SimOpt").sublist("Solve");
339 atol_ = list.get("Absolute Residual Tolerance", DEFAULT_atol_);
340 rtol_ = list.get("Relative Residual Tolerance", DEFAULT_rtol_);
341 maxit_ = list.get("Iteration Limit", DEFAULT_maxit_);
342 decr_ = list.get("Sufficient Decrease Tolerance", DEFAULT_decr_);
343 stol_ = list.get("Step Tolerance", DEFAULT_stol_);
344 factor_ = list.get("Backtracking Factor", DEFAULT_factor_);
345 print_ = list.get("Output Iteration History", DEFAULT_print_);
346 zero_ = list.get("Zero Initial Guess", DEFAULT_zero_);
347 solverType_ = list.get("Solver Type", DEFAULT_solverType_);
348 }
349
365 virtual void applyJacobian_1(Vector<Real> &jv,
366 const Vector<Real> &v,
367 const Vector<Real> &u,
368 const Vector<Real> &z,
369 Real &tol) {
370 Real ctol = std::sqrt(ROL_EPSILON<Real>());
371 // Compute step length
372 Real h = tol;
373 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
374 h = std::max(1.0,u.norm()/v.norm())*tol;
375 }
376 // Update state vector to u + hv
377 Ptr<Vector<Real>> unew = u.clone();
378 unew->set(u);
379 unew->axpy(h,v);
380 // Compute new constraint value
381 update(*unew,z,UpdateType::Temp);
382 value(jv,*unew,z,ctol);
383 // Compute current constraint value
384 Ptr<Vector<Real>> cold = jv.clone();
385 update(u,z,UpdateType::Temp);
386 value(*cold,u,z,ctol);
387 // Compute Newton quotient
388 jv.axpy(-1.0,*cold);
389 jv.scale(1.0/h);
390 }
391
392
408 virtual void applyJacobian_2(Vector<Real> &jv,
409 const Vector<Real> &v,
410 const Vector<Real> &u,
411 const Vector<Real> &z,
412 Real &tol) {
413 Real ctol = std::sqrt(ROL_EPSILON<Real>());
414 // Compute step length
415 Real h = tol;
416 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
417 h = std::max(1.0,u.norm()/v.norm())*tol;
418 }
419 // Update state vector to u + hv
420 Ptr<Vector<Real>> znew = z.clone();
421 znew->set(z);
422 znew->axpy(h,v);
423 // Compute new constraint value
424 update(u,*znew,UpdateType::Temp);
425 value(jv,u,*znew,ctol);
426 // Compute current constraint value
427 Ptr<Vector<Real>> cold = jv.clone();
428 update(u,z,UpdateType::Temp);
429 value(*cold,u,z,ctol);
430 // Compute Newton quotient
431 jv.axpy(-1.0,*cold);
432 jv.scale(1.0/h);
433 }
434
450 virtual void applyInverseJacobian_1(Vector<Real> &ijv,
451 const Vector<Real> &v,
452 const Vector<Real> &u,
453 const Vector<Real> &z,
454 Real &tol) {
455 ROL_TEST_FOR_EXCEPTION(true, std::logic_error,
456 "The method applyInverseJacobian_1 is used but not implemented!\n");
457 }
458
474 virtual void applyAdjointJacobian_1(Vector<Real> &ajv,
475 const Vector<Real> &v,
476 const Vector<Real> &u,
477 const Vector<Real> &z,
478 Real &tol) {
479 applyAdjointJacobian_1(ajv, v, u, z, v.dual(), tol);
480 }
481
482
500 virtual void applyAdjointJacobian_1(Vector<Real> &ajv,
501 const Vector<Real> &v,
502 const Vector<Real> &u,
503 const Vector<Real> &z,
504 const Vector<Real> &dualv,
505 Real &tol) {
506 Real ctol = std::sqrt(ROL_EPSILON<Real>());
507 Real h = tol;
508 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
509 h = std::max(1.0,u.norm()/v.norm())*tol;
510 }
511 Ptr<Vector<Real>> cold = dualv.clone();
512 Ptr<Vector<Real>> cnew = dualv.clone();
513 update(u,z,UpdateType::Temp);
514 value(*cold,u,z,ctol);
515 Ptr<Vector<Real>> unew = u.clone();
516 ajv.zero();
517 for (int i = 0; i < u.dimension(); i++) {
518 unew->set(u);
519 unew->axpy(h,*(u.basis(i)));
520 update(*unew,z,UpdateType::Temp);
521 value(*cnew,*unew,z,ctol);
522 cnew->axpy(-1.0,*cold);
523 cnew->scale(1.0/h);
524 ajv.axpy(cnew->dot(v),*((u.dual()).basis(i)));
525 }
526 update(u,z,UpdateType::Temp);
527 }
528
529
545 virtual void applyAdjointJacobian_2(Vector<Real> &ajv,
546 const Vector<Real> &v,
547 const Vector<Real> &u,
548 const Vector<Real> &z,
549 Real &tol) {
550 applyAdjointJacobian_2(ajv, v, u, z, v.dual(), tol);
551 }
552
553
571 virtual void applyAdjointJacobian_2(Vector<Real> &ajv,
572 const Vector<Real> &v,
573 const Vector<Real> &u,
574 const Vector<Real> &z,
575 const Vector<Real> &dualv,
576 Real &tol) {
577 Real ctol = std::sqrt(ROL_EPSILON<Real>());
578 Real h = tol;
579 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
580 h = std::max(1.0,u.norm()/v.norm())*tol;
581 }
582 Ptr<Vector<Real>> cold = dualv.clone();
583 Ptr<Vector<Real>> cnew = dualv.clone();
584 update(u,z,UpdateType::Temp);
585 value(*cold,u,z,ctol);
586 Ptr<Vector<Real>> znew = z.clone();
587 ajv.zero();
588 for (int i = 0; i < z.dimension(); i++) {
589 znew->set(z);
590 znew->axpy(h,*(z.basis(i)));
591 update(u,*znew,UpdateType::Temp);
592 value(*cnew,u,*znew,ctol);
593 cnew->axpy(-1.0,*cold);
594 cnew->scale(1.0/h);
595 ajv.axpy(cnew->dot(v),*((z.dual()).basis(i)));
596 }
597 update(u,z,UpdateType::Temp);
598 }
599
615 virtual void applyInverseAdjointJacobian_1(Vector<Real> &iajv,
616 const Vector<Real> &v,
617 const Vector<Real> &u,
618 const Vector<Real> &z,
619 Real &tol) {
620 ROL_TEST_FOR_EXCEPTION(true, std::logic_error,
621 "The method applyInverseAdjointJacobian_1 is used but not implemented!\n");
622 };
623
641 virtual void applyAdjointHessian_11(Vector<Real> &ahwv,
642 const Vector<Real> &w,
643 const Vector<Real> &v,
644 const Vector<Real> &u,
645 const Vector<Real> &z,
646 Real &tol) {
647 Real jtol = std::sqrt(ROL_EPSILON<Real>());
648 // Compute step size
649 Real h = tol;
650 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
651 h = std::max(1.0,u.norm()/v.norm())*tol;
652 }
653 // Evaluate Jacobian at new state
654 Ptr<Vector<Real>> unew = u.clone();
655 unew->set(u);
656 unew->axpy(h,v);
657 update(*unew,z,UpdateType::Temp);
658 applyAdjointJacobian_1(ahwv,w,*unew,z,jtol);
659 // Evaluate Jacobian at old state
660 Ptr<Vector<Real>> jv = ahwv.clone();
661 update(u,z,UpdateType::Temp);
662 applyAdjointJacobian_1(*jv,w,u,z,jtol);
663 // Compute Newton quotient
664 ahwv.axpy(-1.0,*jv);
665 ahwv.scale(1.0/h);
666 }
667
668
686 virtual void applyAdjointHessian_12(Vector<Real> &ahwv,
687 const Vector<Real> &w,
688 const Vector<Real> &v,
689 const Vector<Real> &u,
690 const Vector<Real> &z,
691 Real &tol) {
692 Real jtol = std::sqrt(ROL_EPSILON<Real>());
693 // Compute step size
694 Real h = tol;
695 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
696 h = std::max(1.0,u.norm()/v.norm())*tol;
697 }
698 // Evaluate Jacobian at new state
699 Ptr<Vector<Real>> unew = u.clone();
700 unew->set(u);
701 unew->axpy(h,v);
702 update(*unew,z,UpdateType::Temp);
703 applyAdjointJacobian_2(ahwv,w,*unew,z,jtol);
704 // Evaluate Jacobian at old state
705 Ptr<Vector<Real>> jv = ahwv.clone();
706 update(u,z,UpdateType::Temp);
707 applyAdjointJacobian_2(*jv,w,u,z,jtol);
708 // Compute Newton quotient
709 ahwv.axpy(-1.0,*jv);
710 ahwv.scale(1.0/h);
711 }
712
713
731 virtual void applyAdjointHessian_21(Vector<Real> &ahwv,
732 const Vector<Real> &w,
733 const Vector<Real> &v,
734 const Vector<Real> &u,
735 const Vector<Real> &z,
736 Real &tol) {
737 Real jtol = std::sqrt(ROL_EPSILON<Real>());
738 // Compute step size
739 Real h = tol;
740 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
741 h = std::max(1.0,u.norm()/v.norm())*tol;
742 }
743 // Evaluate Jacobian at new control
744 Ptr<Vector<Real>> znew = z.clone();
745 znew->set(z);
746 znew->axpy(h,v);
747 update(u,*znew,UpdateType::Temp);
748 applyAdjointJacobian_1(ahwv,w,u,*znew,jtol);
749 // Evaluate Jacobian at old control
750 Ptr<Vector<Real>> jv = ahwv.clone();
751 update(u,z,UpdateType::Temp);
752 applyAdjointJacobian_1(*jv,w,u,z,jtol);
753 // Compute Newton quotient
754 ahwv.axpy(-1.0,*jv);
755 ahwv.scale(1.0/h);
756 }
757
775 virtual void applyAdjointHessian_22(Vector<Real> &ahwv,
776 const Vector<Real> &w,
777 const Vector<Real> &v,
778 const Vector<Real> &u,
779 const Vector<Real> &z,
780 Real &tol) {
781 Real jtol = std::sqrt(ROL_EPSILON<Real>());
782 // Compute step size
783 Real h = tol;
784 if (v.norm() > std::sqrt(ROL_EPSILON<Real>())) {
785 h = std::max(1.0,u.norm()/v.norm())*tol;
786 }
787 // Evaluate Jacobian at new control
788 Ptr<Vector<Real>> znew = z.clone();
789 znew->set(z);
790 znew->axpy(h,v);
791 update(u,*znew,UpdateType::Temp);
792 applyAdjointJacobian_2(ahwv,w,u,*znew,jtol);
793 // Evaluate Jacobian at old control
794 Ptr<Vector<Real>> jv = ahwv.clone();
795 update(u,z,UpdateType::Temp);
796 applyAdjointJacobian_2(*jv,w,u,z,jtol);
797 // Compute Newton quotient
798 ahwv.axpy(-1.0,*jv);
799 ahwv.scale(1.0/h);
800}
801
840 virtual std::vector<Real> solveAugmentedSystem(Vector<Real> &v1,
841 Vector<Real> &v2,
842 const Vector<Real> &b1,
843 const Vector<Real> &b2,
844 const Vector<Real> &x,
845 Real &tol) {
846 return Constraint<Real>::solveAugmentedSystem(v1,v2,b1,b2,x,tol);
847 }
848
849
868 virtual void applyPreconditioner(Vector<Real> &pv,
869 const Vector<Real> &v,
870 const Vector<Real> &x,
871 const Vector<Real> &g,
872 Real &tol) {
873 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(x);
874 Ptr<Vector<Real>> ijv = (xs.get_1())->clone();
875
876 try {
877 applyInverseJacobian_1(*ijv, v, *(xs.get_1()), *(xs.get_2()), tol);
878 }
879 catch (const std::logic_error &e) {
880 Constraint<Real>::applyPreconditioner(pv, v, x, g, tol);
881 return;
882 }
883
884 const Vector_SimOpt<Real> &gs = dynamic_cast<const Vector_SimOpt<Real>&>(g);
885 Ptr<Vector<Real>> ijv_dual = (gs.get_1())->clone();
886 ijv_dual->set(ijv->dual());
887
888 try {
889 applyInverseAdjointJacobian_1(pv, *ijv_dual, *(xs.get_1()), *(xs.get_2()), tol);
890 }
891 catch (const std::logic_error &e) {
892 Constraint<Real>::applyPreconditioner(pv, v, x, g, tol);
893 return;
894 }
895
896 }
897
903 virtual void update( const Vector<Real> &x, bool flag = true, int iter = -1 ) {
904 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(
905 dynamic_cast<const Vector<Real>&>(x));
906 update(*(xs.get_1()),*(xs.get_2()),flag,iter);
907 }
908 virtual void update( const Vector<Real> &x, UpdateType type, int iter = -1 ) {
909 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(
910 dynamic_cast<const Vector<Real>&>(x));
911 update(*(xs.get_1()),*(xs.get_2()),type,iter);
912 }
913
914 virtual void value(Vector<Real> &c,
915 const Vector<Real> &x,
916 Real &tol) {
917 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(
918 dynamic_cast<const Vector<Real>&>(x));
919 value(c,*(xs.get_1()),*(xs.get_2()),tol);
920 }
921
922
923 virtual void applyJacobian(Vector<Real> &jv,
924 const Vector<Real> &v,
925 const Vector<Real> &x,
926 Real &tol) {
927 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(
928 dynamic_cast<const Vector<Real>&>(x));
929 const Vector_SimOpt<Real> &vs = dynamic_cast<const Vector_SimOpt<Real>&>(
930 dynamic_cast<const Vector<Real>&>(v));
931 applyJacobian_1(jv,*(vs.get_1()),*(xs.get_1()),*(xs.get_2()),tol);
932 Ptr<Vector<Real>> jv2 = jv.clone();
933 applyJacobian_2(*jv2,*(vs.get_2()),*(xs.get_1()),*(xs.get_2()),tol);
934 jv.plus(*jv2);
935 }
936
937
938 using Constraint<Real>::applyAdjointJacobian;
939 virtual void applyAdjointJacobian(Vector<Real> &ajv,
940 const Vector<Real> &v,
941 const Vector<Real> &x,
942 Real &tol) {
943 Vector_SimOpt<Real> &ajvs = dynamic_cast<Vector_SimOpt<Real>&>(
944 dynamic_cast<Vector<Real>&>(ajv));
945 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(
946 dynamic_cast<const Vector<Real>&>(x));
947 Ptr<Vector<Real>> ajv1 = (ajvs.get_1())->clone();
948 applyAdjointJacobian_1(*ajv1,v,*(xs.get_1()),*(xs.get_2()),tol);
949 ajvs.set_1(*ajv1);
950 Ptr<Vector<Real>> ajv2 = (ajvs.get_2())->clone();
951 applyAdjointJacobian_2(*ajv2,v,*(xs.get_1()),*(xs.get_2()),tol);
952 ajvs.set_2(*ajv2);
953 }
954
955
956 virtual void applyAdjointHessian(Vector<Real> &ahwv,
957 const Vector<Real> &w,
958 const Vector<Real> &v,
959 const Vector<Real> &x,
960 Real &tol) {
961 Vector_SimOpt<Real> &ahwvs = dynamic_cast<Vector_SimOpt<Real>&>(
962 dynamic_cast<Vector<Real>&>(ahwv));
963 const Vector_SimOpt<Real> &xs = dynamic_cast<const Vector_SimOpt<Real>&>(
964 dynamic_cast<const Vector<Real>&>(x));
965 const Vector_SimOpt<Real> &vs = dynamic_cast<const Vector_SimOpt<Real>&>(
966 dynamic_cast<const Vector<Real>&>(v));
967 // Block-row 1
968 Ptr<Vector<Real>> C11 = (ahwvs.get_1())->clone();
969 Ptr<Vector<Real>> C21 = (ahwvs.get_1())->clone();
970 applyAdjointHessian_11(*C11,w,*(vs.get_1()),*(xs.get_1()),*(xs.get_2()),tol);
971 applyAdjointHessian_21(*C21,w,*(vs.get_2()),*(xs.get_1()),*(xs.get_2()),tol);
972 C11->plus(*C21);
973 ahwvs.set_1(*C11);
974 // Block-row 2
975 Ptr<Vector<Real>> C12 = (ahwvs.get_2())->clone();
976 Ptr<Vector<Real>> C22 = (ahwvs.get_2())->clone();
977 applyAdjointHessian_12(*C12,w,*(vs.get_1()),*(xs.get_1()),*(xs.get_2()),tol);
978 applyAdjointHessian_22(*C22,w,*(vs.get_2()),*(xs.get_1()),*(xs.get_2()),tol);
979 C22->plus(*C12);
980 ahwvs.set_2(*C22);
981 }
982
983
984
985 virtual Real checkSolve(const Vector<Real> &u,
986 const Vector<Real> &z,
987 const Vector<Real> &c,
988 const bool printToStream = true,
989 std::ostream & outStream = std::cout) {
990 // Solve constraint for u.
991 Real tol = ROL_EPSILON<Real>();
992 Ptr<Vector<Real>> r = c.clone();
993 Ptr<Vector<Real>> s = u.clone();
994 solve(*r,*s,z,tol);
995 // Evaluate constraint residual at (u,z).
996 Ptr<Vector<Real>> cs = c.clone();
997 update(*s,z,UpdateType::Temp);
998 value(*cs,*s,z,tol);
999 // Output norm of residual.
1000 Real rnorm = r->norm();
1001 Real cnorm = cs->norm();
1002 if ( printToStream ) {
1003 std::stringstream hist;
1004 hist << std::scientific << std::setprecision(8);
1005 hist << "\nTest SimOpt solve at feasible (u,z):\n";
1006 hist << " Solver Residual = " << rnorm << "\n";
1007 hist << " ||c(u,z)|| = " << cnorm << "\n";
1008 outStream << hist.str();
1009 }
1010 return cnorm;
1011 }
1012
1013
1026 virtual Real checkAdjointConsistencyJacobian_1(const Vector<Real> &w,
1027 const Vector<Real> &v,
1028 const Vector<Real> &u,
1029 const Vector<Real> &z,
1030 const bool printToStream = true,
1031 std::ostream & outStream = std::cout) {
1032 return checkAdjointConsistencyJacobian_1(w, v, u, z, w.dual(), v.dual(), printToStream, outStream);
1033 }
1034
1035
1051 virtual Real checkAdjointConsistencyJacobian_1(const Vector<Real> &w,
1052 const Vector<Real> &v,
1053 const Vector<Real> &u,
1054 const Vector<Real> &z,
1055 const Vector<Real> &dualw,
1056 const Vector<Real> &dualv,
1057 const bool printToStream = true,
1058 std::ostream & outStream = std::cout) {
1059 Real tol = ROL_EPSILON<Real>();
1060 Ptr<Vector<Real>> Jv = dualw.clone();
1061 update(u,z,UpdateType::Temp);
1062 applyJacobian_1(*Jv,v,u,z,tol);
1063 //Real wJv = w.dot(Jv->dual());
1064 Real wJv = w.apply(*Jv);
1065 Ptr<Vector<Real>> Jw = dualv.clone();
1066 update(u,z,UpdateType::Temp);
1067 applyAdjointJacobian_1(*Jw,w,u,z,tol);
1068 //Real vJw = v.dot(Jw->dual());
1069 Real vJw = v.apply(*Jw);
1070 Real diff = std::abs(wJv-vJw);
1071 if ( printToStream ) {
1072 std::stringstream hist;
1073 hist << std::scientific << std::setprecision(8);
1074 hist << "\nTest SimOpt consistency of Jacobian_1 and its adjoint: \n |<w,Jv> - <adj(J)w,v>| = "
1075 << diff << "\n";
1076 hist << " |<w,Jv>| = " << std::abs(wJv) << "\n";
1077 hist << " Relative Error = " << diff / (std::abs(wJv)+ROL_UNDERFLOW<Real>()) << "\n";
1078 outStream << hist.str();
1079 }
1080 return diff;
1081 }
1082
1083
1096 virtual Real checkAdjointConsistencyJacobian_2(const Vector<Real> &w,
1097 const Vector<Real> &v,
1098 const Vector<Real> &u,
1099 const Vector<Real> &z,
1100 const bool printToStream = true,
1101 std::ostream & outStream = std::cout) {
1102 return checkAdjointConsistencyJacobian_2(w, v, u, z, w.dual(), v.dual(), printToStream, outStream);
1103 }
1104
1120 virtual Real checkAdjointConsistencyJacobian_2(const Vector<Real> &w,
1121 const Vector<Real> &v,
1122 const Vector<Real> &u,
1123 const Vector<Real> &z,
1124 const Vector<Real> &dualw,
1125 const Vector<Real> &dualv,
1126 const bool printToStream = true,
1127 std::ostream & outStream = std::cout) {
1128 Real tol = ROL_EPSILON<Real>();
1129 Ptr<Vector<Real>> Jv = dualw.clone();
1130 update(u,z,UpdateType::Temp);
1131 applyJacobian_2(*Jv,v,u,z,tol);
1132 //Real wJv = w.dot(Jv->dual());
1133 Real wJv = w.apply(*Jv);
1134 Ptr<Vector<Real>> Jw = dualv.clone();
1135 update(u,z,UpdateType::Temp);
1136 applyAdjointJacobian_2(*Jw,w,u,z,tol);
1137 //Real vJw = v.dot(Jw->dual());
1138 Real vJw = v.apply(*Jw);
1139 Real diff = std::abs(wJv-vJw);
1140 if ( printToStream ) {
1141 std::stringstream hist;
1142 hist << std::scientific << std::setprecision(8);
1143 hist << "\nTest SimOpt consistency of Jacobian_2 and its adjoint: \n |<w,Jv> - <adj(J)w,v>| = "
1144 << diff << "\n";
1145 hist << " |<w,Jv>| = " << std::abs(wJv) << "\n";
1146 hist << " Relative Error = " << diff / (std::abs(wJv)+ROL_UNDERFLOW<Real>()) << "\n";
1147 outStream << hist.str();
1148 }
1149 return diff;
1150 }
1151
1152 virtual Real checkInverseJacobian_1(const Vector<Real> &jv,
1153 const Vector<Real> &v,
1154 const Vector<Real> &u,
1155 const Vector<Real> &z,
1156 const bool printToStream = true,
1157 std::ostream & outStream = std::cout) {
1158 Real tol = ROL_EPSILON<Real>();
1159 Ptr<Vector<Real>> Jv = jv.clone();
1160 update(u,z,UpdateType::Temp);
1161 applyJacobian_1(*Jv,v,u,z,tol);
1162 Ptr<Vector<Real>> iJJv = u.clone();
1163 //update(u,z); // Does this update do anything?
1164 applyInverseJacobian_1(*iJJv,*Jv,u,z,tol);
1165 Ptr<Vector<Real>> diff = v.clone();
1166 diff->set(v);
1167 diff->axpy(-1.0,*iJJv);
1168 Real dnorm = diff->norm();
1169 Real vnorm = v.norm();
1170 if ( printToStream ) {
1171 std::stringstream hist;
1172 hist << std::scientific << std::setprecision(8);
1173 hist << "\nTest SimOpt consistency of inverse Jacobian_1: \n ||v-inv(J)Jv|| = "
1174 << dnorm << "\n";
1175 hist << " ||v|| = " << vnorm << "\n";
1176 hist << " Relative Error = " << dnorm / (vnorm+ROL_UNDERFLOW<Real>()) << "\n";
1177 outStream << hist.str();
1178 }
1179 return dnorm;
1180 }
1181
1182 virtual Real checkInverseAdjointJacobian_1(const Vector<Real> &jv,
1183 const Vector<Real> &v,
1184 const Vector<Real> &u,
1185 const Vector<Real> &z,
1186 const bool printToStream = true,
1187 std::ostream & outStream = std::cout) {
1188 Real tol = ROL_EPSILON<Real>();
1189 Ptr<Vector<Real>> Jv = jv.clone();
1190 update(u,z,UpdateType::Temp);
1191 applyAdjointJacobian_1(*Jv,v,u,z,tol);
1192 Ptr<Vector<Real>> iJJv = v.clone();
1193 //update(u,z);
1194 applyInverseAdjointJacobian_1(*iJJv,*Jv,u,z,tol);
1195 Ptr<Vector<Real>> diff = v.clone();
1196 diff->set(v);
1197 diff->axpy(-1.0,*iJJv);
1198 Real dnorm = diff->norm();
1199 Real vnorm = v.norm();
1200 if ( printToStream ) {
1201 std::stringstream hist;
1202 hist << std::scientific << std::setprecision(8);
1203 hist << "\nTest SimOpt consistency of inverse adjoint Jacobian_1: \n ||v-inv(adj(J))adj(J)v|| = "
1204 << dnorm << "\n";
1205 hist << " ||v|| = " << vnorm << "\n";
1206 hist << " Relative Error = " << dnorm / (vnorm+ROL_UNDERFLOW<Real>()) << "\n";
1207 outStream << hist.str();
1208 }
1209 return dnorm;
1210 }
1211
1212
1213
1214 std::vector<std::vector<Real>> checkApplyJacobian_1(const Vector<Real> &u,
1215 const Vector<Real> &z,
1216 const Vector<Real> &v,
1217 const Vector<Real> &jv,
1218 const bool printToStream = true,
1219 std::ostream & outStream = std::cout,
1220 const int numSteps = ROL_NUM_CHECKDERIV_STEPS,
1221 const int order = 1) {
1222 std::vector<Real> steps(numSteps);
1223 for(int i=0;i<numSteps;++i) {
1224 steps[i] = pow(10,-i);
1225 }
1226
1227 return checkApplyJacobian_1(u,z,v,jv,steps,printToStream,outStream,order);
1228 }
1229
1230
1231
1232
1233 std::vector<std::vector<Real>> checkApplyJacobian_1(const Vector<Real> &u,
1234 const Vector<Real> &z,
1235 const Vector<Real> &v,
1236 const Vector<Real> &jv,
1237 const std::vector<Real> &steps,
1238 const bool printToStream = true,
1239 std::ostream & outStream = std::cout,
1240 const int order = 1) {
1241
1242 ROL_TEST_FOR_EXCEPTION( order<1 || order>4, std::invalid_argument,
1243 "Error: finite difference order must be 1,2,3, or 4" );
1244
1245 Real one(1.0);
1246
1247 using Finite_Difference_Arrays::shifts;
1248 using Finite_Difference_Arrays::weights;
1249
1250 Real tol = std::sqrt(ROL_EPSILON<Real>());
1251
1252 int numSteps = steps.size();
1253 int numVals = 4;
1254 std::vector<Real> tmp(numVals);
1255 std::vector<std::vector<Real>> jvCheck(numSteps, tmp);
1256
1257 // Save the format state of the original outStream.
1258 nullstream oldFormatState;
1259 oldFormatState.copyfmt(outStream);
1260
1261 // Compute constraint value at x.
1262 Ptr<Vector<Real>> c = jv.clone();
1263 this->update(u,z,UpdateType::Temp);
1264 this->value(*c, u, z, tol);
1265
1266 // Compute (Jacobian at x) times (vector v).
1267 Ptr<Vector<Real>> Jv = jv.clone();
1268 this->applyJacobian_1(*Jv, v, u, z, tol);
1269 Real normJv = Jv->norm();
1270
1271 // Temporary vectors.
1272 Ptr<Vector<Real>> cdif = jv.clone();
1273 Ptr<Vector<Real>> cnew = jv.clone();
1274 Ptr<Vector<Real>> unew = u.clone();
1275
1276 for (int i=0; i<numSteps; i++) {
1277
1278 Real eta = steps[i];
1279
1280 unew->set(u);
1281
1282 cdif->set(*c);
1283 cdif->scale(weights[order-1][0]);
1284
1285 for(int j=0; j<order; ++j) {
1286
1287 unew->axpy(eta*shifts[order-1][j], v);
1288
1289 if( weights[order-1][j+1] != 0 ) {
1290 this->update(*unew,z,UpdateType::Temp);
1291 this->value(*cnew,*unew,z,tol);
1292 cdif->axpy(weights[order-1][j+1],*cnew);
1293 }
1294
1295 }
1296
1297 cdif->scale(one/eta);
1298
1299 // Compute norms of Jacobian-vector products, finite-difference approximations, and error.
1300 jvCheck[i][0] = eta;
1301 jvCheck[i][1] = normJv;
1302 jvCheck[i][2] = cdif->norm();
1303 cdif->axpy(-one, *Jv);
1304 jvCheck[i][3] = cdif->norm();
1305
1306 if (printToStream) {
1307 std::stringstream hist;
1308 if (i==0) {
1309 hist << std::right
1310 << std::setw(20) << "Step size"
1311 << std::setw(20) << "norm(Jac*vec)"
1312 << std::setw(20) << "norm(FD approx)"
1313 << std::setw(20) << "norm(abs error)"
1314 << "\n"
1315 << std::setw(20) << "---------"
1316 << std::setw(20) << "-------------"
1317 << std::setw(20) << "---------------"
1318 << std::setw(20) << "---------------"
1319 << "\n";
1320 }
1321 hist << std::scientific << std::setprecision(11) << std::right
1322 << std::setw(20) << jvCheck[i][0]
1323 << std::setw(20) << jvCheck[i][1]
1324 << std::setw(20) << jvCheck[i][2]
1325 << std::setw(20) << jvCheck[i][3]
1326 << "\n";
1327 outStream << hist.str();
1328 }
1329
1330 }
1331
1332 // Reset format state of outStream.
1333 outStream.copyfmt(oldFormatState);
1334
1335 return jvCheck;
1336 } // checkApplyJacobian
1337
1338
1339 std::vector<std::vector<Real>> checkApplyJacobian_2(const Vector<Real> &u,
1340 const Vector<Real> &z,
1341 const Vector<Real> &v,
1342 const Vector<Real> &jv,
1343 const bool printToStream = true,
1344 std::ostream & outStream = std::cout,
1345 const int numSteps = ROL_NUM_CHECKDERIV_STEPS,
1346 const int order = 1) {
1347 std::vector<Real> steps(numSteps);
1348 for(int i=0;i<numSteps;++i) {
1349 steps[i] = pow(10,-i);
1350 }
1351
1352 return checkApplyJacobian_2(u,z,v,jv,steps,printToStream,outStream,order);
1353 }
1354
1355
1356
1357
1358 std::vector<std::vector<Real>> checkApplyJacobian_2(const Vector<Real> &u,
1359 const Vector<Real> &z,
1360 const Vector<Real> &v,
1361 const Vector<Real> &jv,
1362 const std::vector<Real> &steps,
1363 const bool printToStream = true,
1364 std::ostream & outStream = std::cout,
1365 const int order = 1) {
1366
1367 ROL_TEST_FOR_EXCEPTION( order<1 || order>4, std::invalid_argument,
1368 "Error: finite difference order must be 1,2,3, or 4" );
1369
1370 Real one(1.0);
1371
1372 using Finite_Difference_Arrays::shifts;
1373 using Finite_Difference_Arrays::weights;
1374
1375 Real tol = std::sqrt(ROL_EPSILON<Real>());
1376
1377 int numSteps = steps.size();
1378 int numVals = 4;
1379 std::vector<Real> tmp(numVals);
1380 std::vector<std::vector<Real>> jvCheck(numSteps, tmp);
1381
1382 // Save the format state of the original outStream.
1383 nullstream oldFormatState;
1384 oldFormatState.copyfmt(outStream);
1385
1386 // Compute constraint value at x.
1387 Ptr<Vector<Real>> c = jv.clone();
1388 this->update(u,z,UpdateType::Temp);
1389 this->value(*c, u, z, tol);
1390
1391 // Compute (Jacobian at x) times (vector v).
1392 Ptr<Vector<Real>> Jv = jv.clone();
1393 this->applyJacobian_2(*Jv, v, u, z, tol);
1394 Real normJv = Jv->norm();
1395
1396 // Temporary vectors.
1397 Ptr<Vector<Real>> cdif = jv.clone();
1398 Ptr<Vector<Real>> cnew = jv.clone();
1399 Ptr<Vector<Real>> znew = z.clone();
1400
1401 for (int i=0; i<numSteps; i++) {
1402
1403 Real eta = steps[i];
1404
1405 znew->set(z);
1406
1407 cdif->set(*c);
1408 cdif->scale(weights[order-1][0]);
1409
1410 for(int j=0; j<order; ++j) {
1411
1412 znew->axpy(eta*shifts[order-1][j], v);
1413
1414 if( weights[order-1][j+1] != 0 ) {
1415 this->update(u,*znew,UpdateType::Temp);
1416 this->value(*cnew,u,*znew,tol);
1417 cdif->axpy(weights[order-1][j+1],*cnew);
1418 }
1419
1420 }
1421
1422 cdif->scale(one/eta);
1423
1424 // Compute norms of Jacobian-vector products, finite-difference approximations, and error.
1425 jvCheck[i][0] = eta;
1426 jvCheck[i][1] = normJv;
1427 jvCheck[i][2] = cdif->norm();
1428 cdif->axpy(-one, *Jv);
1429 jvCheck[i][3] = cdif->norm();
1430
1431 if (printToStream) {
1432 std::stringstream hist;
1433 if (i==0) {
1434 hist << std::right
1435 << std::setw(20) << "Step size"
1436 << std::setw(20) << "norm(Jac*vec)"
1437 << std::setw(20) << "norm(FD approx)"
1438 << std::setw(20) << "norm(abs error)"
1439 << "\n"
1440 << std::setw(20) << "---------"
1441 << std::setw(20) << "-------------"
1442 << std::setw(20) << "---------------"
1443 << std::setw(20) << "---------------"
1444 << "\n";
1445 }
1446 hist << std::scientific << std::setprecision(11) << std::right
1447 << std::setw(20) << jvCheck[i][0]
1448 << std::setw(20) << jvCheck[i][1]
1449 << std::setw(20) << jvCheck[i][2]
1450 << std::setw(20) << jvCheck[i][3]
1451 << "\n";
1452 outStream << hist.str();
1453 }
1454
1455 }
1456
1457 // Reset format state of outStream.
1458 outStream.copyfmt(oldFormatState);
1459
1460 return jvCheck;
1461 } // checkApplyJacobian
1462
1463
1464
1465 std::vector<std::vector<Real>> checkApplyAdjointHessian_11(const Vector<Real> &u,
1466 const Vector<Real> &z,
1467 const Vector<Real> &p,
1468 const Vector<Real> &v,
1469 const Vector<Real> &hv,
1470 const bool printToStream = true,
1471 std::ostream & outStream = std::cout,
1472 const int numSteps = ROL_NUM_CHECKDERIV_STEPS,
1473 const int order = 1 ) {
1474 std::vector<Real> steps(numSteps);
1475 for(int i=0;i<numSteps;++i) {
1476 steps[i] = pow(10,-i);
1477 }
1478
1479 return checkApplyAdjointHessian_11(u,z,p,v,hv,steps,printToStream,outStream,order);
1480
1481 }
1482
1483 std::vector<std::vector<Real>> checkApplyAdjointHessian_11(const Vector<Real> &u,
1484 const Vector<Real> &z,
1485 const Vector<Real> &p,
1486 const Vector<Real> &v,
1487 const Vector<Real> &hv,
1488 const std::vector<Real> &steps,
1489 const bool printToStream = true,
1490 std::ostream & outStream = std::cout,
1491 const int order = 1 ) {
1492 using Finite_Difference_Arrays::shifts;
1493 using Finite_Difference_Arrays::weights;
1494
1495 Real one(1.0);
1496
1497 Real tol = std::sqrt(ROL_EPSILON<Real>());
1498
1499 int numSteps = steps.size();
1500 int numVals = 4;
1501 std::vector<Real> tmp(numVals);
1502 std::vector<std::vector<Real>> ahpvCheck(numSteps, tmp);
1503
1504 // Temporary vectors.
1505 Ptr<Vector<Real>> AJdif = hv.clone();
1506 Ptr<Vector<Real>> AJp = hv.clone();
1507 Ptr<Vector<Real>> AHpv = hv.clone();
1508 Ptr<Vector<Real>> AJnew = hv.clone();
1509 Ptr<Vector<Real>> unew = u.clone();
1510
1511 // Save the format state of the original outStream.
1512 nullstream oldFormatState;
1513 oldFormatState.copyfmt(outStream);
1514
1515 // Apply adjoint Jacobian to p.
1516 this->update(u,z,UpdateType::Temp);
1517 this->applyAdjointJacobian_1(*AJp, p, u, z, tol);
1518
1519 // Apply adjoint Hessian at (u,z), in direction v, to p.
1520 this->applyAdjointHessian_11(*AHpv, p, v, u, z, tol);
1521 Real normAHpv = AHpv->norm();
1522
1523 for (int i=0; i<numSteps; i++) {
1524
1525 Real eta = steps[i];
1526
1527 // Apply adjoint Jacobian to p at (u+eta*v,z).
1528 unew->set(u);
1529
1530 AJdif->set(*AJp);
1531 AJdif->scale(weights[order-1][0]);
1532
1533 for(int j=0; j<order; ++j) {
1534
1535 unew->axpy(eta*shifts[order-1][j],v);
1536
1537 if( weights[order-1][j+1] != 0 ) {
1538 this->update(*unew,z,UpdateType::Temp);
1539 this->applyAdjointJacobian_1(*AJnew, p, *unew, z, tol);
1540 AJdif->axpy(weights[order-1][j+1],*AJnew);
1541 }
1542 }
1543
1544 AJdif->scale(one/eta);
1545
1546 // Compute norms of Jacobian-vector products, finite-difference approximations, and error.
1547 ahpvCheck[i][0] = eta;
1548 ahpvCheck[i][1] = normAHpv;
1549 ahpvCheck[i][2] = AJdif->norm();
1550 AJdif->axpy(-one, *AHpv);
1551 ahpvCheck[i][3] = AJdif->norm();
1552
1553 if (printToStream) {
1554 std::stringstream hist;
1555 if (i==0) {
1556 hist << std::right
1557 << std::setw(20) << "Step size"
1558 << std::setw(20) << "norm(adj(H)(u,v))"
1559 << std::setw(20) << "norm(FD approx)"
1560 << std::setw(20) << "norm(abs error)"
1561 << "\n"
1562 << std::setw(20) << "---------"
1563 << std::setw(20) << "-----------------"
1564 << std::setw(20) << "---------------"
1565 << std::setw(20) << "---------------"
1566 << "\n";
1567 }
1568 hist << std::scientific << std::setprecision(11) << std::right
1569 << std::setw(20) << ahpvCheck[i][0]
1570 << std::setw(20) << ahpvCheck[i][1]
1571 << std::setw(20) << ahpvCheck[i][2]
1572 << std::setw(20) << ahpvCheck[i][3]
1573 << "\n";
1574 outStream << hist.str();
1575 }
1576
1577 }
1578
1579 // Reset format state of outStream.
1580 outStream.copyfmt(oldFormatState);
1581
1582 return ahpvCheck;
1583 } // checkApplyAdjointHessian_11
1584
1588 std::vector<std::vector<Real>> checkApplyAdjointHessian_21(const Vector<Real> &u,
1589 const Vector<Real> &z,
1590 const Vector<Real> &p,
1591 const Vector<Real> &v,
1592 const Vector<Real> &hv,
1593 const bool printToStream = true,
1594 std::ostream & outStream = std::cout,
1595 const int numSteps = ROL_NUM_CHECKDERIV_STEPS,
1596 const int order = 1 ) {
1597 std::vector<Real> steps(numSteps);
1598 for(int i=0;i<numSteps;++i) {
1599 steps[i] = pow(10,-i);
1600 }
1601
1602 return checkApplyAdjointHessian_21(u,z,p,v,hv,steps,printToStream,outStream,order);
1603
1604 }
1605
1609 std::vector<std::vector<Real>> checkApplyAdjointHessian_21(const Vector<Real> &u,
1610 const Vector<Real> &z,
1611 const Vector<Real> &p,
1612 const Vector<Real> &v,
1613 const Vector<Real> &hv,
1614 const std::vector<Real> &steps,
1615 const bool printToStream = true,
1616 std::ostream & outStream = std::cout,
1617 const int order = 1 ) {
1618 using Finite_Difference_Arrays::shifts;
1619 using Finite_Difference_Arrays::weights;
1620
1621 Real one(1.0);
1622
1623 Real tol = std::sqrt(ROL_EPSILON<Real>());
1624
1625 int numSteps = steps.size();
1626 int numVals = 4;
1627 std::vector<Real> tmp(numVals);
1628 std::vector<std::vector<Real>> ahpvCheck(numSteps, tmp);
1629
1630 // Temporary vectors.
1631 Ptr<Vector<Real>> AJdif = hv.clone();
1632 Ptr<Vector<Real>> AJp = hv.clone();
1633 Ptr<Vector<Real>> AHpv = hv.clone();
1634 Ptr<Vector<Real>> AJnew = hv.clone();
1635 Ptr<Vector<Real>> znew = z.clone();
1636
1637 // Save the format state of the original outStream.
1638 nullstream oldFormatState;
1639 oldFormatState.copyfmt(outStream);
1640
1641 // Apply adjoint Jacobian to p.
1642 this->update(u,z,UpdateType::Temp);
1643 this->applyAdjointJacobian_1(*AJp, p, u, z, tol);
1644
1645 // Apply adjoint Hessian at (u,z), in direction v, to p.
1646 this->applyAdjointHessian_21(*AHpv, p, v, u, z, tol);
1647 Real normAHpv = AHpv->norm();
1648
1649 for (int i=0; i<numSteps; i++) {
1650
1651 Real eta = steps[i];
1652
1653 // Apply adjoint Jacobian to p at (u,z+eta*v).
1654 znew->set(z);
1655
1656 AJdif->set(*AJp);
1657 AJdif->scale(weights[order-1][0]);
1658
1659 for(int j=0; j<order; ++j) {
1660
1661 znew->axpy(eta*shifts[order-1][j],v);
1662
1663 if( weights[order-1][j+1] != 0 ) {
1664 this->update(u,*znew,UpdateType::Temp);
1665 this->applyAdjointJacobian_1(*AJnew, p, u, *znew, tol);
1666 AJdif->axpy(weights[order-1][j+1],*AJnew);
1667 }
1668 }
1669
1670 AJdif->scale(one/eta);
1671
1672 // Compute norms of Jacobian-vector products, finite-difference approximations, and error.
1673 ahpvCheck[i][0] = eta;
1674 ahpvCheck[i][1] = normAHpv;
1675 ahpvCheck[i][2] = AJdif->norm();
1676 AJdif->axpy(-one, *AHpv);
1677 ahpvCheck[i][3] = AJdif->norm();
1678
1679 if (printToStream) {
1680 std::stringstream hist;
1681 if (i==0) {
1682 hist << std::right
1683 << std::setw(20) << "Step size"
1684 << std::setw(20) << "norm(adj(H)(u,v))"
1685 << std::setw(20) << "norm(FD approx)"
1686 << std::setw(20) << "norm(abs error)"
1687 << "\n"
1688 << std::setw(20) << "---------"
1689 << std::setw(20) << "-----------------"
1690 << std::setw(20) << "---------------"
1691 << std::setw(20) << "---------------"
1692 << "\n";
1693 }
1694 hist << std::scientific << std::setprecision(11) << std::right
1695 << std::setw(20) << ahpvCheck[i][0]
1696 << std::setw(20) << ahpvCheck[i][1]
1697 << std::setw(20) << ahpvCheck[i][2]
1698 << std::setw(20) << ahpvCheck[i][3]
1699 << "\n";
1700 outStream << hist.str();
1701 }
1702
1703 }
1704
1705 // Reset format state of outStream.
1706 outStream.copyfmt(oldFormatState);
1707
1708 return ahpvCheck;
1709 } // checkApplyAdjointHessian_21
1710
1714 std::vector<std::vector<Real>> checkApplyAdjointHessian_12(const Vector<Real> &u,
1715 const Vector<Real> &z,
1716 const Vector<Real> &p,
1717 const Vector<Real> &v,
1718 const Vector<Real> &hv,
1719 const bool printToStream = true,
1720 std::ostream & outStream = std::cout,
1721 const int numSteps = ROL_NUM_CHECKDERIV_STEPS,
1722 const int order = 1 ) {
1723 std::vector<Real> steps(numSteps);
1724 for(int i=0;i<numSteps;++i) {
1725 steps[i] = pow(10,-i);
1726 }
1727
1728 return checkApplyAdjointHessian_12(u,z,p,v,hv,steps,printToStream,outStream,order);
1729
1730 }
1731
1732
1733 std::vector<std::vector<Real>> checkApplyAdjointHessian_12(const Vector<Real> &u,
1734 const Vector<Real> &z,
1735 const Vector<Real> &p,
1736 const Vector<Real> &v,
1737 const Vector<Real> &hv,
1738 const std::vector<Real> &steps,
1739 const bool printToStream = true,
1740 std::ostream & outStream = std::cout,
1741 const int order = 1 ) {
1742 using Finite_Difference_Arrays::shifts;
1743 using Finite_Difference_Arrays::weights;
1744
1745 Real one(1.0);
1746
1747 Real tol = std::sqrt(ROL_EPSILON<Real>());
1748
1749 int numSteps = steps.size();
1750 int numVals = 4;
1751 std::vector<Real> tmp(numVals);
1752 std::vector<std::vector<Real>> ahpvCheck(numSteps, tmp);
1753
1754 // Temporary vectors.
1755 Ptr<Vector<Real>> AJdif = hv.clone();
1756 Ptr<Vector<Real>> AJp = hv.clone();
1757 Ptr<Vector<Real>> AHpv = hv.clone();
1758 Ptr<Vector<Real>> AJnew = hv.clone();
1759 Ptr<Vector<Real>> unew = u.clone();
1760
1761 // Save the format state of the original outStream.
1762 nullstream oldFormatState;
1763 oldFormatState.copyfmt(outStream);
1764
1765 // Apply adjoint Jacobian to p.
1766 this->update(u,z,UpdateType::Temp);
1767 this->applyAdjointJacobian_2(*AJp, p, u, z, tol);
1768
1769 // Apply adjoint Hessian at (u,z), in direction v, to p.
1770 this->applyAdjointHessian_12(*AHpv, p, v, u, z, tol);
1771 Real normAHpv = AHpv->norm();
1772
1773 for (int i=0; i<numSteps; i++) {
1774
1775 Real eta = steps[i];
1776
1777 // Apply adjoint Jacobian to p at (u+eta*v,z).
1778 unew->set(u);
1779
1780 AJdif->set(*AJp);
1781 AJdif->scale(weights[order-1][0]);
1782
1783 for(int j=0; j<order; ++j) {
1784
1785 unew->axpy(eta*shifts[order-1][j],v);
1786
1787 if( weights[order-1][j+1] != 0 ) {
1788 this->update(*unew,z,UpdateType::Temp);
1789 this->applyAdjointJacobian_2(*AJnew, p, *unew, z, tol);
1790 AJdif->axpy(weights[order-1][j+1],*AJnew);
1791 }
1792 }
1793
1794 AJdif->scale(one/eta);
1795
1796 // Compute norms of Jacobian-vector products, finite-difference approximations, and error.
1797 ahpvCheck[i][0] = eta;
1798 ahpvCheck[i][1] = normAHpv;
1799 ahpvCheck[i][2] = AJdif->norm();
1800 AJdif->axpy(-one, *AHpv);
1801 ahpvCheck[i][3] = AJdif->norm();
1802
1803 if (printToStream) {
1804 std::stringstream hist;
1805 if (i==0) {
1806 hist << std::right
1807 << std::setw(20) << "Step size"
1808 << std::setw(20) << "norm(adj(H)(u,v))"
1809 << std::setw(20) << "norm(FD approx)"
1810 << std::setw(20) << "norm(abs error)"
1811 << "\n"
1812 << std::setw(20) << "---------"
1813 << std::setw(20) << "-----------------"
1814 << std::setw(20) << "---------------"
1815 << std::setw(20) << "---------------"
1816 << "\n";
1817 }
1818 hist << std::scientific << std::setprecision(11) << std::right
1819 << std::setw(20) << ahpvCheck[i][0]
1820 << std::setw(20) << ahpvCheck[i][1]
1821 << std::setw(20) << ahpvCheck[i][2]
1822 << std::setw(20) << ahpvCheck[i][3]
1823 << "\n";
1824 outStream << hist.str();
1825 }
1826
1827 }
1828
1829 // Reset format state of outStream.
1830 outStream.copyfmt(oldFormatState);
1831
1832 return ahpvCheck;
1833 } // checkApplyAdjointHessian_12
1834
1835 std::vector<std::vector<Real>> checkApplyAdjointHessian_22(const Vector<Real> &u,
1836 const Vector<Real> &z,
1837 const Vector<Real> &p,
1838 const Vector<Real> &v,
1839 const Vector<Real> &hv,
1840 const bool printToStream = true,
1841 std::ostream & outStream = std::cout,
1842 const int numSteps = ROL_NUM_CHECKDERIV_STEPS,
1843 const int order = 1 ) {
1844 std::vector<Real> steps(numSteps);
1845 for(int i=0;i<numSteps;++i) {
1846 steps[i] = pow(10,-i);
1847 }
1848
1849 return checkApplyAdjointHessian_22(u,z,p,v,hv,steps,printToStream,outStream,order);
1850
1851 }
1852
1853 std::vector<std::vector<Real>> checkApplyAdjointHessian_22(const Vector<Real> &u,
1854 const Vector<Real> &z,
1855 const Vector<Real> &p,
1856 const Vector<Real> &v,
1857 const Vector<Real> &hv,
1858 const std::vector<Real> &steps,
1859 const bool printToStream = true,
1860 std::ostream & outStream = std::cout,
1861 const int order = 1 ) {
1862 using Finite_Difference_Arrays::shifts;
1863 using Finite_Difference_Arrays::weights;
1864
1865 Real one(1.0);
1866
1867 Real tol = std::sqrt(ROL_EPSILON<Real>());
1868
1869 int numSteps = steps.size();
1870 int numVals = 4;
1871 std::vector<Real> tmp(numVals);
1872 std::vector<std::vector<Real>> ahpvCheck(numSteps, tmp);
1873
1874 // Temporary vectors.
1875 Ptr<Vector<Real>> AJdif = hv.clone();
1876 Ptr<Vector<Real>> AJp = hv.clone();
1877 Ptr<Vector<Real>> AHpv = hv.clone();
1878 Ptr<Vector<Real>> AJnew = hv.clone();
1879 Ptr<Vector<Real>> znew = z.clone();
1880
1881 // Save the format state of the original outStream.
1882 nullstream oldFormatState;
1883 oldFormatState.copyfmt(outStream);
1884
1885 // Apply adjoint Jacobian to p.
1886 this->update(u,z,UpdateType::Temp);
1887 this->applyAdjointJacobian_2(*AJp, p, u, z, tol);
1888
1889 // Apply adjoint Hessian at (u,z), in direction v, to p.
1890 this->applyAdjointHessian_22(*AHpv, p, v, u, z, tol);
1891 Real normAHpv = AHpv->norm();
1892
1893 for (int i=0; i<numSteps; i++) {
1894
1895 Real eta = steps[i];
1896
1897 // Apply adjoint Jacobian to p at (u,z+eta*v).
1898 znew->set(z);
1899
1900 AJdif->set(*AJp);
1901 AJdif->scale(weights[order-1][0]);
1902
1903 for(int j=0; j<order; ++j) {
1904
1905 znew->axpy(eta*shifts[order-1][j],v);
1906
1907 if( weights[order-1][j+1] != 0 ) {
1908 this->update(u,*znew,UpdateType::Temp);
1909 this->applyAdjointJacobian_2(*AJnew, p, u, *znew, tol);
1910 AJdif->axpy(weights[order-1][j+1],*AJnew);
1911 }
1912 }
1913
1914 AJdif->scale(one/eta);
1915
1916 // Compute norms of Jacobian-vector products, finite-difference approximations, and error.
1917 ahpvCheck[i][0] = eta;
1918 ahpvCheck[i][1] = normAHpv;
1919 ahpvCheck[i][2] = AJdif->norm();
1920 AJdif->axpy(-one, *AHpv);
1921 ahpvCheck[i][3] = AJdif->norm();
1922
1923 if (printToStream) {
1924 std::stringstream hist;
1925 if (i==0) {
1926 hist << std::right
1927 << std::setw(20) << "Step size"
1928 << std::setw(20) << "norm(adj(H)(u,v))"
1929 << std::setw(20) << "norm(FD approx)"
1930 << std::setw(20) << "norm(abs error)"
1931 << "\n"
1932 << std::setw(20) << "---------"
1933 << std::setw(20) << "-----------------"
1934 << std::setw(20) << "---------------"
1935 << std::setw(20) << "---------------"
1936 << "\n";
1937 }
1938 hist << std::scientific << std::setprecision(11) << std::right
1939 << std::setw(20) << ahpvCheck[i][0]
1940 << std::setw(20) << ahpvCheck[i][1]
1941 << std::setw(20) << ahpvCheck[i][2]
1942 << std::setw(20) << ahpvCheck[i][3]
1943 << "\n";
1944 outStream << hist.str();
1945 }
1946
1947 }
1948
1949 // Reset format state of outStream.
1950 outStream.copyfmt(oldFormatState);
1951
1952 return ahpvCheck;
1953 } // checkApplyAdjointHessian_22
1954
1955}; // class Constraint_SimOpt
1956
1957} // namespace ROL
1958
1959#endif
ROL_Constraint.hpp
ROL_NUM_CHECKDERIV_STEPS
#define ROL_NUM_CHECKDERIV_STEPS
ROL_NonlinearLeastSquaresObjective.hpp
ROL_Objective_FSsolver.hpp
type
Real type
Definition ROL_ScalarTraits.hpp:55
iter
int iter
Definition ROL_Secant.hpp:72
ROL_SimConstraint.hpp
ROL_TypeE_AugmentedLagrangianAlgorithm.hpp
ROL_TypeU_TrustRegionAlgorithm.hpp
ROL_Types.hpp
Contains definitions of custom data types in ROL.
cnorm
Real cnorm
Definition ROL_Types.hpp:152
value
Real value
Definition ROL_Types.hpp:149
flag
bool flag
Definition ROL_Types.hpp:156
ROL_Vector_SimOpt.hpp
Constraint_SimOpt
Defines the constraint operator interface for simulation-based optimization.
Constraint_SimOpt::checkInverseJacobian_1
virtual Real checkInverseJacobian_1(const Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const bool printToStream=true, std::ostream &outStream=std::cout)
Definition ROL_Constraint_SimOpt.hpp:1152
Constraint_SimOpt::checkApplyAdjointHessian_12
std::vector< std::vector< Real > > checkApplyAdjointHessian_12(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1733
Constraint_SimOpt::solveAugmentedSystem
virtual std::vector< Real > solveAugmentedSystem(Vector< Real > &v1, Vector< Real > &v2, const Vector< Real > &b1, const Vector< Real > &b2, const Vector< Real > &x, Real &tol)
Approximately solves the augmented system
Definition ROL_Constraint_SimOpt.hpp:840
Constraint_SimOpt::checkApplyAdjointHessian_21
std::vector< std::vector< Real > > checkApplyAdjointHessian_21(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
, , , ,
Definition ROL_Constraint_SimOpt.hpp:1588
Constraint_SimOpt::print_
bool print_
Definition ROL_Constraint_SimOpt.hpp:134
Constraint_SimOpt::applyJacobian
virtual void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Definition ROL_Constraint_SimOpt.hpp:923
Constraint_SimOpt::update
virtual 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 i...
Definition ROL_Constraint_SimOpt.hpp:903
Constraint_SimOpt::rtol_
Real rtol_
Definition ROL_Constraint_SimOpt.hpp:129
Constraint_SimOpt::value
virtual void value(Vector< Real > &c, const Vector< Real > &u, const Vector< Real > &z, Real &tol)=0
Evaluate the constraint operator at .
Constraint_SimOpt::applyAdjointHessian_22
virtual void applyAdjointHessian_22(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the optimization-space derivative of the adjoint of the constraint optimization-space Jacobian ...
Definition ROL_Constraint_SimOpt.hpp:775
Constraint_SimOpt::checkAdjointConsistencyJacobian_2
virtual Real checkAdjointConsistencyJacobian_2(const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &dualw, const Vector< Real > &dualv, const bool printToStream=true, std::ostream &outStream=std::cout)
Check the consistency of the Jacobian and its adjoint. This is the secondary interface,...
Definition ROL_Constraint_SimOpt.hpp:1120
Constraint_SimOpt::update
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Definition ROL_Constraint_SimOpt.hpp:908
Constraint_SimOpt::checkInverseAdjointJacobian_1
virtual Real checkInverseAdjointJacobian_1(const Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const bool printToStream=true, std::ostream &outStream=std::cout)
Definition ROL_Constraint_SimOpt.hpp:1182
Constraint_SimOpt::DEFAULT_print_
const bool DEFAULT_print_
Definition ROL_Constraint_SimOpt.hpp:121
Constraint_SimOpt::update
virtual void update(const Vector< Real > &u, const Vector< Real > &z, UpdateType type, int iter=-1)
Definition ROL_Constraint_SimOpt.hpp:167
Constraint_SimOpt::decr_
Real decr_
Definition ROL_Constraint_SimOpt.hpp:132
Constraint_SimOpt::applyJacobian_1
virtual void applyJacobian_1(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the partial constraint Jacobian at , , to the vector .
Definition ROL_Constraint_SimOpt.hpp:365
Constraint_SimOpt::DEFAULT_zero_
const bool DEFAULT_zero_
Definition ROL_Constraint_SimOpt.hpp:122
Constraint_SimOpt::checkApplyAdjointHessian_11
std::vector< std::vector< Real > > checkApplyAdjointHessian_11(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1483
Constraint_SimOpt::solve_update
virtual void solve_update(const Vector< Real > &u, const Vector< Real > &z, UpdateType type, int iter=-1)
Update SimOpt constraint during solve (disconnected from optimization updates).
Definition ROL_Constraint_SimOpt.hpp:194
Constraint_SimOpt::checkAdjointConsistencyJacobian_2
virtual Real checkAdjointConsistencyJacobian_2(const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const bool printToStream=true, std::ostream &outStream=std::cout)
Check the consistency of the Jacobian and its adjoint. This is the primary interface.
Definition ROL_Constraint_SimOpt.hpp:1096
Constraint_SimOpt::checkApplyJacobian_2
std::vector< std::vector< Real > > checkApplyJacobian_2(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &v, const Vector< Real > &jv, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1339
Constraint_SimOpt::checkApplyAdjointHessian_11
std::vector< std::vector< Real > > checkApplyAdjointHessian_11(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1465
Constraint_SimOpt::DEFAULT_atol_
const Real DEFAULT_atol_
Definition ROL_Constraint_SimOpt.hpp:115
Constraint_SimOpt::Constraint_SimOpt
Constraint_SimOpt()
Definition ROL_Constraint_SimOpt.hpp:142
Constraint_SimOpt::checkApplyJacobian_1
std::vector< std::vector< Real > > checkApplyJacobian_1(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &v, const Vector< Real > &jv, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1233
Constraint_SimOpt::applyAdjointJacobian_2
virtual void applyAdjointJacobian_2(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the adjoint of the partial constraint Jacobian at , , to vector . This is the primary interface...
Definition ROL_Constraint_SimOpt.hpp:545
Constraint_SimOpt::applyAdjointHessian_12
virtual void applyAdjointHessian_12(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the optimization-space derivative of the adjoint of the constraint simulation-space Jacobian at...
Definition ROL_Constraint_SimOpt.hpp:686
Constraint_SimOpt::applyAdjointHessian
virtual void applyAdjointHessian(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Definition ROL_Constraint_SimOpt.hpp:956
Constraint_SimOpt::update_1
virtual void update_1(const Vector< Real > &u, bool flag=true, int iter=-1)
Update constraint functions with respect to Sim variable. x is the optimization variable,...
Definition ROL_Constraint_SimOpt.hpp:177
Constraint_SimOpt::checkSolve
virtual Real checkSolve(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &c, const bool printToStream=true, std::ostream &outStream=std::cout)
Definition ROL_Constraint_SimOpt.hpp:985
Constraint_SimOpt::checkApplyAdjointHessian_22
std::vector< std::vector< Real > > checkApplyAdjointHessian_22(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1853
Constraint_SimOpt::applyInverseAdjointJacobian_1
virtual void applyInverseAdjointJacobian_1(Vector< Real > &iajv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the inverse of the adjoint of the partial constraint Jacobian at , , to the vector .
Definition ROL_Constraint_SimOpt.hpp:615
Constraint_SimOpt::checkApplyJacobian_2
std::vector< std::vector< Real > > checkApplyJacobian_2(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &v, const Vector< Real > &jv, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1358
Constraint_SimOpt::setSolveParameters
virtual void setSolveParameters(ParameterList &parlist)
Set solve parameters.
Definition ROL_Constraint_SimOpt.hpp:337
Constraint_SimOpt::update_2
virtual void update_2(const Vector< Real > &z, UpdateType type, int iter=-1)
Definition ROL_Constraint_SimOpt.hpp:186
Constraint_SimOpt::applyAdjointJacobian_1
virtual void applyAdjointJacobian_1(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &dualv, Real &tol)
Apply the adjoint of the partial constraint Jacobian at , , to the vector . This is the secondary int...
Definition ROL_Constraint_SimOpt.hpp:500
Constraint_SimOpt::zero_
bool zero_
Definition ROL_Constraint_SimOpt.hpp:135
Constraint_SimOpt::DEFAULT_decr_
const Real DEFAULT_decr_
Definition ROL_Constraint_SimOpt.hpp:119
Constraint_SimOpt::checkApplyAdjointHessian_12
std::vector< std::vector< Real > > checkApplyAdjointHessian_12(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
, , , ,
Definition ROL_Constraint_SimOpt.hpp:1714
Constraint_SimOpt::applyAdjointHessian_21
virtual void applyAdjointHessian_21(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the simulation-space derivative of the adjoint of the constraint optimization-space Jacobian at...
Definition ROL_Constraint_SimOpt.hpp:731
Constraint_SimOpt::unew_
Ptr< Vector< Real > > unew_
Definition ROL_Constraint_SimOpt.hpp:111
Constraint_SimOpt::DEFAULT_solverType_
const int DEFAULT_solverType_
Definition ROL_Constraint_SimOpt.hpp:123
Constraint_SimOpt::checkApplyAdjointHessian_22
std::vector< std::vector< Real > > checkApplyAdjointHessian_22(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1835
Constraint_SimOpt::applyPreconditioner
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 . In general, this preconditioner satisfies the fo...
Definition ROL_Constraint_SimOpt.hpp:868
Constraint_SimOpt::maxit_
int maxit_
Definition ROL_Constraint_SimOpt.hpp:133
Constraint_SimOpt::applyAdjointHessian_11
virtual void applyAdjointHessian_11(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the simulation-space derivative of the adjoint of the constraint simulation-space Jacobian at ...
Definition ROL_Constraint_SimOpt.hpp:641
Constraint_SimOpt::checkApplyAdjointHessian_21
std::vector< std::vector< Real > > checkApplyAdjointHessian_21(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &p, const Vector< Real > &v, const Vector< Real > &hv, const std::vector< Real > &steps, const bool printToStream=true, std::ostream &outStream=std::cout, const int order=1)
, , , ,
Definition ROL_Constraint_SimOpt.hpp:1609
Constraint_SimOpt::applyJacobian_2
virtual void applyJacobian_2(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the partial constraint Jacobian at , , to the vector .
Definition ROL_Constraint_SimOpt.hpp:408
Constraint_SimOpt::DEFAULT_factor_
const Real DEFAULT_factor_
Definition ROL_Constraint_SimOpt.hpp:118
Constraint_SimOpt::DEFAULT_maxit_
const int DEFAULT_maxit_
Definition ROL_Constraint_SimOpt.hpp:120
Constraint_SimOpt::checkAdjointConsistencyJacobian_1
virtual Real checkAdjointConsistencyJacobian_1(const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &dualw, const Vector< Real > &dualv, const bool printToStream=true, std::ostream &outStream=std::cout)
Check the consistency of the Jacobian and its adjoint. This is the secondary interface,...
Definition ROL_Constraint_SimOpt.hpp:1051
Constraint_SimOpt::DEFAULT_stol_
const Real DEFAULT_stol_
Definition ROL_Constraint_SimOpt.hpp:117
Constraint_SimOpt::update
virtual void update(const Vector< Real > &u, const Vector< Real > &z, bool flag=true, int iter=-1)
Update constraint functions. x is the optimization variable, flag = true if optimization variable i...
Definition ROL_Constraint_SimOpt.hpp:163
Constraint_SimOpt::DEFAULT_rtol_
const Real DEFAULT_rtol_
Definition ROL_Constraint_SimOpt.hpp:116
Constraint_SimOpt::update_2
virtual void update_2(const Vector< Real > &z, bool flag=true, int iter=-1)
Update constraint functions with respect to Opt variable. x is the optimization variable,...
Definition ROL_Constraint_SimOpt.hpp:185
Constraint_SimOpt::checkApplyJacobian_1
std::vector< std::vector< Real > > checkApplyJacobian_1(const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &v, const Vector< Real > &jv, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Definition ROL_Constraint_SimOpt.hpp:1214
Constraint_SimOpt::applyAdjointJacobian_2
virtual void applyAdjointJacobian_2(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const Vector< Real > &dualv, Real &tol)
Apply the adjoint of the partial constraint Jacobian at , , to vector . This is the secondary interfa...
Definition ROL_Constraint_SimOpt.hpp:571
Constraint_SimOpt::factor_
Real factor_
Definition ROL_Constraint_SimOpt.hpp:131
Constraint_SimOpt::jv_
Ptr< Vector< Real > > jv_
Definition ROL_Constraint_SimOpt.hpp:112
Constraint_SimOpt::applyInverseJacobian_1
virtual void applyInverseJacobian_1(Vector< Real > &ijv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the inverse partial constraint Jacobian at , , to the vector .
Definition ROL_Constraint_SimOpt.hpp:450
Constraint_SimOpt::atol_
Real atol_
Definition ROL_Constraint_SimOpt.hpp:128
Constraint_SimOpt::stol_
Real stol_
Definition ROL_Constraint_SimOpt.hpp:130
Constraint_SimOpt::update_1
virtual void update_1(const Vector< Real > &u, UpdateType type, int iter=-1)
Definition ROL_Constraint_SimOpt.hpp:178
Constraint_SimOpt::firstSolve_
bool firstSolve_
Definition ROL_Constraint_SimOpt.hpp:139
Constraint_SimOpt::checkAdjointConsistencyJacobian_1
virtual Real checkAdjointConsistencyJacobian_1(const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, const bool printToStream=true, std::ostream &outStream=std::cout)
Check the consistency of the Jacobian and its adjoint. This is the primary interface.
Definition ROL_Constraint_SimOpt.hpp:1026
Constraint_SimOpt::solve
virtual void solve(Vector< Real > &c, Vector< Real > &u, const Vector< Real > &z, Real &tol)
Given , solve for .
Definition ROL_Constraint_SimOpt.hpp:225
Constraint_SimOpt::value
virtual void value(Vector< Real > &c, const Vector< Real > &x, Real &tol)
Definition ROL_Constraint_SimOpt.hpp:914
Constraint_SimOpt::applyAdjointJacobian_1
virtual void applyAdjointJacobian_1(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z, Real &tol)
Apply the adjoint of the partial constraint Jacobian at , , to the vector . This is the primary inter...
Definition ROL_Constraint_SimOpt.hpp:474
Constraint_SimOpt::solverType_
int solverType_
Definition ROL_Constraint_SimOpt.hpp:136
Constraint_SimOpt::applyAdjointJacobian
virtual void applyAdjointJacobian(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Definition ROL_Constraint_SimOpt.hpp:939
ROL::ROL::Constraint_SimOpt
Definition ROL_Constraint_SimOpt.hpp:108
ROL::ROL::Constraint
Defines the general constraint operator interface.
Definition ROL_Constraint_SerialSimOpt.hpp:86
Generated by doxygen 1.12.0