22 #include "Teuchos_GlobalMPISession.hpp"
29 int main(
int argc,
char *argv[]) {
31 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
34 int iprint = argc - 1;
35 ROL::Ptr<std::ostream> outStream;
38 outStream = ROL::makePtrFromRef(std::cout);
40 outStream = ROL::makePtrFromRef(bhs);
48 ROL::Ptr<ROL::Objective<RealT>> obj;
49 ROL::Ptr<ROL::Constraint<RealT>> constr;
50 ROL::Ptr<ROL::Vector<RealT>> x;
51 ROL::Ptr<ROL::Vector<RealT>> sol;
63 RealT left = -1e0, right = 1e0;
64 ROL::Ptr<std::vector<RealT>> xtest_ptr = ROL::makePtr<std::vector<RealT>>(
dim, 0.0);
65 ROL::Ptr<std::vector<RealT>> g_ptr = ROL::makePtr<std::vector<RealT>>(
dim, 0.0);
66 ROL::Ptr<std::vector<RealT>> d_ptr = ROL::makePtr<std::vector<RealT>>(
dim, 0.0);
67 ROL::Ptr<std::vector<RealT>> v_ptr = ROL::makePtr<std::vector<RealT>>(
dim, 0.0);
68 ROL::Ptr<std::vector<RealT>> vc_ptr = ROL::makePtr<std::vector<RealT>>(nc, 0.0);
69 ROL::Ptr<std::vector<RealT>> vl_ptr = ROL::makePtr<std::vector<RealT>>(nc, 0.0);
77 for (
int i=0; i<
dim; i++) {
78 (*xtest_ptr)[i] = ( (
RealT)rand() / (
RealT)RAND_MAX ) * (right - left) + left;
79 (*d_ptr)[i] = ( (
RealT)rand() / (
RealT)RAND_MAX ) * (right - left) + left;
80 (*v_ptr)[i] = ( (
RealT)rand() / (
RealT)RAND_MAX ) * (right - left) + left;
83 for (
int i=0; i<nc; i++) {
84 (*vc_ptr)[i] = ( (
RealT)rand() / (
RealT)RAND_MAX ) * (right - left) + left;
85 (*vl_ptr)[i] = ( (
RealT)rand() / (
RealT)RAND_MAX ) * (right - left) + left;
95 constr->checkApplyJacobian(xtest, v, vc,
true, *outStream); *outStream <<
"\n";
96 constr->checkApplyAdjointJacobian(xtest, vl, vc, xtest,
true, *outStream); *outStream <<
"\n";
97 constr->checkApplyAdjointHessian(xtest, vl, d, xtest,
true, *outStream); *outStream <<
"\n";
98 nlls.
checkGradient(xtest, d,
true, *outStream); *outStream <<
"\n";
99 nlls.
checkHessVec(xtest, v,
true, *outStream); *outStream <<
"\n";
100 nlls.
checkHessSym(xtest, d, v,
true, *outStream); *outStream <<
"\n";
103 ROL::ParameterList parlist;
104 parlist.sublist(
"Step").sublist(
"Trust Region").set(
"Subproblem Solver",
"Truncated CG");
105 parlist.sublist(
"Status Test").set(
"Gradient Tolerance",1.e-10);
106 parlist.sublist(
"Status Test").set(
"Constraint Tolerance",1.e-10);
107 parlist.sublist(
"Status Test").set(
"Step Tolerance",1.e-18);
108 parlist.sublist(
"Status Test").set(
"Iteration Limit",100);
109 ROL::Ptr<ROL::StatusTest<RealT>> status = ROL::makePtr<ROL::StatusTest<RealT>>(parlist);
110 ROL::Ptr<ROL::Step<RealT>> step = ROL::makePtr<ROL::TrustRegionStep<RealT>>(parlist);
114 *outStream <<
"\nSOLVE USING FULL HESSIAN\n";
116 algo.run(*x, nlls,
true, *outStream);
118 *outStream <<
"\nSOLVE USING GAUSS-NEWTON HESSIAN\n";
120 algo.run(*x, gnnlls,
true, *outStream);
122 catch (std::logic_error& err) {
123 *outStream << err.what() <<
"\n";
128 std::cout <<
"End Result: TEST FAILED\n";
130 std::cout <<
"End Result: TEST PASSED\n";
Ptr< Constraint< Real > > getEqualityConstraint(void) const
Defines a no-output stream class ROL::NullStream and a function makeStreamPtr which either wraps a re...
virtual std::vector< std::vector< Real > > checkGradient(const Vector< Real > &x, const Vector< Real > &d, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Finite-difference gradient check.
Provides the ROL::Vector interface for scalar values, to be used, for example, with scalar constraint...
Provides an interface to run optimization algorithms.
Ptr< Objective< Real > > getObjective(void) const
Contains definitions for the equality constrained NLP from Nocedal/Wright, 2nd edition, page 574, example 18.2; note the typo in reversing the initial guess and the solution.
void set(const Vector< Real > &x)
Set where .
Ptr< Vector< Real > > getSolution(const int i=0) const
basic_nullstream< char, char_traits< char >> nullstream
int main(int argc, char *argv[])
virtual std::vector< std::vector< Real > > checkHessVec(const Vector< Real > &x, const Vector< Real > &v, const bool printToStream=true, std::ostream &outStream=std::cout, const int numSteps=ROL_NUM_CHECKDERIV_STEPS, const int order=1)
Finite-difference Hessian-applied-to-vector check.
Provides the interface to evaluate nonlinear least squares objective functions.
virtual std::vector< Real > checkHessSym(const Vector< Real > &x, const Vector< Real > &v, const Vector< Real > &w, const bool printToStream=true, std::ostream &outStream=std::cout)
Hessian symmetry check.
Ptr< Vector< Real > > getInitialGuess(void) const