21 #include "Teuchos_GlobalMPISession.hpp"
23 template<
typename Real>
26 void value(std::vector<Real> &c,
const std::vector<Real> &x, Real &tol) {
29 void applyJacobian(std::vector<Real> &jv,
const std::vector<Real> &v,
const std::vector<Real> &x, Real &tol) {
32 void applyAdjointJacobian(std::vector<Real> &ajv,
const std::vector<Real> &v,
const std::vector<Real> &x, Real &tol) {
40 int main(
int argc,
char *argv[]) {
42 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
46 int iprint = argc - 1;
47 ROL::Ptr<std::ostream> outStream;
50 outStream = ROL::makePtrFromRef(std::cout);
52 outStream = ROL::makePtrFromRef(bhs);
57 RealT tol = std::sqrt(ROL::ROL_EPSILON<RealT>());
59 ROL::Ptr<con2d<RealT>> con = ROL::makePtr<con2d<RealT>>();
61 ROL::Ptr<std::vector<RealT>> yptr = ROL::makePtr<std::vector<RealT>>(2);
62 (*yptr)[0] =
static_cast<RealT>(rand())/static_cast<RealT>(RAND_MAX);
63 (*yptr)[1] =
static_cast<RealT>(rand())/static_cast<RealT>(RAND_MAX);
68 ROL::Ptr<std::vector<RealT>> xptr = ROL::makePtr<std::vector<RealT>>(2);
69 (*xptr)[0] = (*yptr)[0];
70 (*xptr)[1] = (*yptr)[1];
73 ROL::Ptr<std::vector<RealT>> Pxptr = ROL::makePtr<std::vector<RealT>>(2,0.0);
79 ROL::Ptr<std::vector<RealT>> x0ptr = ROL::makePtr<std::vector<RealT>>(2);
80 (*x0ptr)[0] = ((*yptr)[0]-(*yptr)[1])/static_cast<RealT>(2);
81 (*x0ptr)[1] = -(*x0ptr)[0];
86 *outStream << std::setprecision(6) << std::scientific << std::endl;
87 *outStream <<
" x[0] = " << (*xptr)[0] <<
" x[1] = " << (*xptr)[1] << std::endl;
88 *outStream <<
" Px[0] = " << (*Pxptr)[0] <<
" Px[1] = " << (*Pxptr)[1] << std::endl;
89 *outStream <<
" x*[0] = " << (*x0ptr)[0] <<
" x*[1] = " << (*x0ptr)[1] << std::endl;
91 e0.
set(x0); e0.
axpy(static_cast<RealT>(-1),Px);
93 *outStream <<
" Error in Euclidean Projection: " << err << std::endl;
95 e0.
set(x); e0.
axpy(static_cast<RealT>(-1),x0);
96 *outStream <<
" ||x*-x||^2 = " << e0.
norm() << std::endl;
98 e0.
set(x); e0.
axpy(static_cast<RealT>(-1),Px);
99 *outStream <<
" ||Px-x||^2 = " << e0.
norm() << std::endl << std::endl;
101 errorFlag += (err > tol);
103 ROL::Ptr<std::vector<RealT>> dptr = ROL::makePtr<std::vector<RealT>>(2);
104 (*dptr)[0] =
static_cast<RealT>(1)+static_cast<RealT>(2)*
static_cast<RealT>(rand())/static_cast<RealT>(RAND_MAX);
105 (*dptr)[1] =
static_cast<RealT>(1)+static_cast<RealT>(5)*
static_cast<RealT>(rand())/static_cast<RealT>(RAND_MAX);
107 ROL::Ptr<std::vector<RealT>> x1ptr = ROL::makePtr<std::vector<RealT>>(2);
108 (*x1ptr)[0] = ((*dptr)[0]*(*yptr)[0]-(*dptr)[1]*(*yptr)[1])/((*dptr)[0]+(*dptr)[1]);
109 (*x1ptr)[1] = -(*x1ptr)[0];
112 ROL::Ptr<std::vector<RealT>> zptr = ROL::makePtr<std::vector<RealT>>(2);
113 (*zptr)[0] = (*yptr)[0];
114 (*zptr)[1] = (*yptr)[1];
117 ROL::Ptr<std::vector<RealT>> Pzptr = ROL::makePtr<std::vector<RealT>>(2,0.0);
123 ROL::Ptr<std::vector<RealT>> e1ptr = ROL::makePtr<std::vector<RealT>>(2);
126 *outStream << std::endl;
127 *outStream <<
" x[0] = " << (*zptr)[0] <<
" x[1] = " << (*zptr)[1] << std::endl;
128 *outStream <<
" Px[0] = " << (*Pzptr)[0] <<
" Px[1] = " << (*Pzptr)[1] << std::endl;
129 *outStream <<
" x*[0] = " << (*x1ptr)[0] <<
" x*[1] = " << (*x1ptr)[1] << std::endl;
131 e1.set(x1); e1.axpy(static_cast<RealT>(-1),Pz);
133 *outStream <<
" Error in Euclidean Projection: " << err << std::endl;
135 e1.set(z); e1.axpy(static_cast<RealT>(-1),x1);
136 *outStream <<
" ||x*-x||^2 = " << e1.norm() << std::endl;
138 e1.set(z); e1.axpy(static_cast<RealT>(-1),Pz);
139 *outStream <<
" ||Px-x||^2 = " << e1.norm() << std::endl << std::endl;
141 errorFlag += (err > tol);
144 catch (std::logic_error& err) {
145 *outStream << err.what() <<
"\n";
150 std::cout <<
"End Result: TEST FAILED\n";
152 std::cout <<
"End Result: TEST PASSED\n";
void axpy(const Real alpha, const Vector< Real > &x)
Compute where .
Defines the equality constraint operator interface for StdVectors.
virtual void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
Defines a no-output stream class ROL::NullStream and a function makeStreamPtr which either wraps a re...
Projects on to the null space of a linear constraint.
Real norm() const
Returns where .
Provides the ROL::Vector interface for scalar values, to be used, for example, with scalar constraint...
void applyJacobian(std::vector< Real > &jv, const std::vector< Real > &v, const std::vector< Real > &x, Real &tol)
Provides the std::vector implementation of the ROL::Vector interface that handles scalings in the inn...
void set(const Vector< Real > &x)
Set where .
basic_nullstream< char, char_traits< char >> nullstream
int main(int argc, char *argv[])
void applyAdjointJacobian(std::vector< Real > &ajv, const std::vector< Real > &v, const std::vector< Real > &x, Real &tol)
void value(std::vector< Real > &c, const std::vector< Real > &x, Real &tol)