ROL
ROL_SimConstraint_Def.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_CONSTRAINT_STATE_DEF_H
11 #define ROL_CONSTRAINT_STATE_DEF_H
12 
13 namespace ROL {
14 
15 template<typename Real>
16 SimConstraint<Real>::SimConstraint(const Ptr<Constraint_SimOpt<Real>> &con,
17  const Ptr<const Vector<Real>> &z,
18  bool inSolve) : con_(con), z_(z), inSolve_(inSolve), init_(false) {}
19 
20 template<typename Real>
21 void SimConstraint<Real>::update( const Vector<Real> &u, bool flag, int iter ) {
22  con_->update_1(u,flag,iter);
23  //con_->update(u,*z_,flag,iter);
24 }
25 
26 template<typename Real>
27 void SimConstraint<Real>::update( const Vector<Real> &u, UpdateType type, int iter ) {
28  if (inSolve_) con_->solve_update(u,*z_,type,iter);
29  else con_->update_1(u,type,iter);
30 }
31 
32 template<typename Real>
33 void SimConstraint<Real>::value(Vector<Real> &c,const Vector<Real> &u,Real &tol) {
34  con_->value(c,u,*z_,tol);
35 }
36 
37 template<typename Real>
38 void SimConstraint<Real>::applyJacobian(Vector<Real> &jv,const Vector<Real> &v,const Vector<Real> &u,Real &tol) {
39  con_->applyJacobian_1(jv,v,u,*z_,tol);
40 }
41 
42 template<typename Real>
43 void SimConstraint<Real>::applyAdjointJacobian(Vector<Real> &ajv,const Vector<Real> &v,const Vector<Real> &u,Real &tol) {
44  con_->applyAdjointJacobian_1(ajv,v,u,*z_,tol);
45 }
46 
47 template<typename Real>
48 void SimConstraint<Real>::applyAdjointHessian(Vector<Real> &ahwv,const Vector<Real> &w,const Vector<Real> &v,const Vector<Real> &u,Real &tol) {
49  con_->applyAdjointHessian_11(ahwv,w,v,u,*z_,tol);
50 }
51 
52 template<typename Real>
53 void SimConstraint<Real>::applyPreconditioner(Vector<Real> &pv,const Vector<Real> &v,const Vector<Real> &u,const Vector<Real> &g,Real &tol) {
54  if (!init_) {
55  ijv_ = u.clone();
56  init_ = true;
57  }
58  con_->applyInverseJacobian_1(*ijv_,v,u,*z_,tol);
59  con_->applyInverseAdjointJacobian_1(pv,ijv_->dual(),u,*z_,tol);
60 }
61 
62 template<typename Real>
63 void SimConstraint<Real>::setParameter(const std::vector<Real> &param) {
64  con_->setParameter(param);
66 }
67 
68 } // namespace ROL
69 
70 #endif
SimConstraint(const Ptr< Constraint_SimOpt< Real >> &con, const Ptr< const Vector< Real >> &z, bool inSolve=false)
void value(Vector< Real > &c, const Vector< Real > &u, Real &tol) override
Evaluate the constraint operator at .
void applyAdjointJacobian(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &u, Real &tol) override
Apply the adjoint of the the constraint Jacobian at , , to vector .
void update(const Vector< Real > &u, bool flag=true, int iter=-1) override
Update constraint functions. x is the optimization variable, flag = true if optimization variable is ...
virtual void setParameter(const std::vector< Real > &param)
void applyAdjointHessian(Vector< Real > &ahwv, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, Real &tol) override
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ...
void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &u, Real &tol) override
Apply the constraint Jacobian at , , to vector .
void setParameter(const std::vector< Real > &param) override
void applyPreconditioner(Vector< Real > &pv, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &g, Real &tol) override
Apply a constraint preconditioner at , , to vector . Ideally, this preconditioner satisfies the follo...