1 #include "NonLinearEquation.h"
2 
3 using NEWMAT::ColumnVector;
4 using NEWMAT::Matrix;
5 using NEWMAT::SymmetricMatrix;
6 
7 //------------------------------------------------------------------------
8 // P.J. Williams
9 // Sandia National Laboratories
10 // pwillia@sandia.gov
11 // Last modified 03/20/2003
12 //------------------------------------------------------------------------
13 
14 namespace OPTPP {
15 
16 // Constructors
NonLinearEquation()17 NonLinearEquation::NonLinearEquation():
18    NonLinearConstraint(), b_(0), ctype_(0){}
19 
NonLinearEquation(NLP * nlprob,int numconstraints)20 NonLinearEquation::NonLinearEquation(NLP* nlprob, int numconstraints ):
21    NonLinearConstraint(nlprob,numconstraints), b_(numconstraints),
22    ctype_(numconstraints)
23    {b_ = 0.0; ctype_.ReSize(numOfCons_); ctype_ = NLeqn; }
24 
NonLinearEquation(NLP * nlprob,const ColumnVector & b,int numconstraints)25 NonLinearEquation::NonLinearEquation(NLP* nlprob, const ColumnVector& b,
26                                  int numconstraints ):
27    NonLinearConstraint(nlprob, b, numconstraints), b_(b),
28    ctype_(nlprob->getDim())
29    {ctype_.ReSize(numOfCons_); ctype_ = NLeqn;}
30 
31 // Functions for computing various quantities
evalCFGH(const ColumnVector & xc)32 void NonLinearEquation::evalCFGH(const ColumnVector & xc) const
33 {
34   nlp_->evalC(xc);
35 }
36 
evalResidual(const ColumnVector & xc)37 ColumnVector NonLinearEquation::evalResidual(const ColumnVector& xc) const
38 {
39       int i, index;
40       ColumnVector resid( numOfCons_);
41       cvalue_ = nlp_->evalCF(xc);
42 
43       // 08/06/2001 PJW - Okay to only loop over nnzl.
44       // We assume nnzl = total constraints and nnzu = 0
45       // for the equality constrained case.
46       for( i = 1; i <= nnzl_; i++){
47           index = constraintMappingIndices_[i-1];
48 	  resid(i) = cvalue_(index) - b_(index);
49       }
50       return resid;
51 }
52 
evalGradient(const ColumnVector & xc)53 Matrix NonLinearEquation::evalGradient(const ColumnVector& xc) const
54 {
55       int j, index;
56       Matrix grad(numOfVars_, numOfCons_);
57       Matrix constraint_grad = nlp_->evalCG(xc);
58 
59       for( j = 1; j <= nnzl_; j++){
60           index = constraintMappingIndices_[j-1];
61 	  grad.Column(j) = constraint_grad.Column(index);
62       }
63       return grad;
64 }
65 
evalHessian(ColumnVector & xc)66 SymmetricMatrix NonLinearEquation::evalHessian(ColumnVector& xc) const
67 {
68      SymmetricMatrix hess, constraint_hess;
69      constraint_hess = nlp_->evalCH(xc);
70 
71      hess = constraint_hess;
72      return hess;
73 }
74 
evalHessian(ColumnVector & xc,int darg)75 OptppArray<SymmetricMatrix> NonLinearEquation::evalHessian(ColumnVector& xc,
76 							  int darg) const
77 {
78      int i, index;
79      OptppArray<SymmetricMatrix> hess(numOfCons_);
80      OptppArray<SymmetricMatrix> constraint_hess = nlp_->evalCH(xc,darg);
81 
82      for( i = 1; i <= nnzl_; i++){
83          index = constraintMappingIndices_[i-1];
84 	 hess[i-1] = constraint_hess[index-1];
85      }
86      return hess;
87 }
88 
amIFeasible(const ColumnVector & xc,double epsilon)89 bool NonLinearEquation::amIFeasible(const ColumnVector& xc, double epsilon) const
90 {
91      int i;
92      bool feasible = true;
93      ColumnVector residual = evalResidual(xc);
94      for(i = 1; i <= numOfCons_; i++){
95         if( (residual(i) < -epsilon) || (residual(i) > epsilon) ){
96            feasible = false;
97            break;
98         }
99      }
100      return feasible;
101 }
102 
103 } // namespace OPTPP
104