41 #ifndef ROL_CONSTRAINT_PARTITIONED_DEF_H
42 #define ROL_CONSTRAINT_PARTITIONED_DEF_H
46 template<
typename Real>
50 : cvec_(cvec), offset_(offset),
51 scratch_(nullPtr), ncval_(0), initialized_(false) {
55 template<
typename Real>
57 std::vector<bool> isInequality,
59 : cvec_(cvec), isInequality_(isInequality), offset_(offset),
60 scratch_(nullPtr), ncval_(0), initialized_(false) {}
62 template<
typename Real>
67 template<
typename Real>
69 if (ind < 0 || ind > static_cast<int>(cvec_.size())) {
75 template<
typename Real>
77 const int ncon =
static_cast<int>(cvec_.size());
78 for (
int i = 0; i < ncon; ++i) {
79 cvec_[i]->update(getOpt(x),type,iter);
83 template<
typename Real>
85 const int ncon =
static_cast<int>(cvec_.size());
86 for (
int i = 0; i < ncon; ++i) {
87 cvec_[i]->update(getOpt(x),flag,iter);
91 template<
typename Real>
96 const int ncon =
static_cast<int>(cvec_.size());
98 for (
int i = 0; i < ncon; ++i) {
99 cvec_[i]->value(*cpv.
get(i), getOpt(x), tol);
100 if (isInequality_[i]) {
101 cpv.
get(i)->axpy(static_cast<Real>(-1),getSlack(x,cnt));
108 template<
typename Real>
116 const int ncon =
static_cast<int>(cvec_.size());
118 for (
int i = 0; i < ncon; ++i) {
119 cvec_[i]->applyJacobian(*jvpv.
get(i), getOpt(v), getOpt(x), tol);
120 if (isInequality_[i]) {
121 jvpv.
get(i)->axpy(static_cast<Real>(-1),getSlack(v,cnt));
127 template<
typename Real>
133 scratch_ = getOpt(ajv).clone();
140 const int ncon =
static_cast<int>(cvec_.size());
143 for (
int i = 0; i < ncon; ++i) {
145 cvec_[i]->applyAdjointJacobian(*scratch_, *vpv.
get(i), getOpt(x), tol);
146 getOpt(ajv).plus(*scratch_);
147 if (isInequality_[i]) {
148 getSlack(ajv,cnt).set(*vpv.
get(i));
149 getSlack(ajv,cnt).scale(static_cast<Real>(-1));
155 template<
typename Real>
162 scratch_ = getOpt(ahuv).clone();
169 const int ncon =
static_cast<int>(cvec_.size());
172 for (
int i = 0; i < ncon; ++i) {
173 Ptr<const Vector<Real>> ui = upv.
get(i);
175 cvec_[i]->applyAdjointHessian(*scratch_, *ui, getOpt(v), getOpt(x), tol);
176 getOpt(ahuv).plus(*scratch_);
177 if (isInequality_[i]) {
178 getSlack(ahuv,cnt).zero();
184 template<
typename Real>
195 const int ncon =
static_cast<int>(cvec_.size());
196 for (
int i = 0; i < ncon; ++i) {
197 cvec_[i]->applyPreconditioner(*pvpv.
get(i), *vpv.
get(i), getOpt(x), getOpt(g), tol);
201 template<
typename Real>
204 const int ncon =
static_cast<int>(cvec_.size());
205 for (
int i = 0; i < ncon; ++i) {
206 cvec_[i]->setParameter(param);
210 template<
typename Real>
215 catch (std::exception &e) {
220 template<
typename Real>
225 catch (std::exception &e) {
230 template<
typename Real>
235 template<
typename Real>
void applyAdjointHessian(Vector< Real > &ahuv, const Vector< Real > &u, const Vector< Real > &v, const Vector< Real > &x, Real &tol) override
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ...
void setParameter(const std::vector< Real > ¶m) override
void value(Vector< Real > &c, const Vector< Real > &x, Real &tol) override
Evaluate the constraint operator at .
ROL::Ptr< const Vector< Real > > get(size_type i) const
Defines the linear algebra of vector space on a generic partitioned vector.
Defines the linear algebra or vector space interface.
void applyAdjointJacobian(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &x, Real &tol) override
Apply the adjoint of the the constraint Jacobian at , , to vector .
virtual void applyPreconditioner(Vector< Real > &pv, const Vector< Real > &v, const Vector< Real > &x, const Vector< Real > &g, Real &tol) override
Apply a constraint preconditioner at , , to vector . Ideally, this preconditioner satisfies the follo...
void zero()
Set to zero vector.
Vector< Real > & getOpt(Vector< Real > &xs) const
Constraint_Partitioned(const std::vector< Ptr< Constraint< Real >>> &cvec, bool isInequality=false, int offset=0)
void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &x, Real &tol) override
Apply the constraint Jacobian at , , to vector .
virtual void setParameter(const std::vector< Real > ¶m)
int getNumberConstraintEvaluations(void) const
Vector< Real > & getSlack(Vector< Real > &xs, int ind) const
std::vector< bool > isInequality_
void update(const Vector< Real > &x, UpdateType type, int iter=-1) override
Update constraint function.
Defines the general constraint operator interface.
Ptr< Constraint< Real > > get(int ind=0) const