10 #include "Teuchos_GlobalMPISession.hpp"
12 #include "Teuchos_GlobalMPISession.hpp"
18 #include "ROL_ParameterList.hpp"
32 Real
value(
const V &x, Real &tol ) {
38 void hessVec(
V &hv,
const V &v,
const V &x, Real &tol ) {
45 ROL::Ptr<const std::vector<Real> > xp =
48 for(
size_t i=0; i<xp->size(); ++i ) {
49 outStream << (*xp)[i] << std::endl;
58 int main(
int argc,
char *argv[]) {
60 typedef std::vector<RealT> vector;
67 typedef ROL::ParameterList PL;
69 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
71 int iprint = argc - 1;
72 ROL::Ptr<std::ostream> outStream;
75 outStream = ROL::makePtrFromRef(std::cout);
77 outStream = ROL::makePtrFromRef(bhs);
84 PL &iplist = parlist.sublist(
"Step").sublist(
"Primal Dual Interior Point");
85 PL &lblist = iplist.sublist(
"Barrier Objective");
89 bool useLinearDamping =
true;
91 lblist.set(
"Use Linear Damping", useLinearDamping);
92 lblist.set(
"Linear Damping Coefficient",kappaD);
93 lblist.set(
"Initial Barrier Parameter",mu);
95 RealT ninf = ROL::ROL_NINF<RealT>();
96 RealT inf = ROL::ROL_INF<RealT>();
99 int numTestVectors = 19;
101 ROL::Ptr<vector> x_ptr = ROL::makePtr<vector>(
dim, 0.0);
102 ROL::Ptr<vector> d_ptr = ROL::makePtr<vector>(
dim, 0.0);
103 ROL::Ptr<vector> v_ptr = ROL::makePtr<vector>(
dim, 0.0);
104 ROL::Ptr<vector> l_ptr = ROL::makePtr<vector>(
dim, 0.0);
105 ROL::Ptr<vector> u_ptr = ROL::makePtr<vector>(
dim, 0.0);
106 ROL::Ptr<vector> e0_ptr = ROL::makePtr<vector>(
dim, 0.0);
113 (*l_ptr)[0] = ninf; (*u_ptr)[0] = 5.0;
114 (*l_ptr)[1] = ninf; (*u_ptr)[1] = inf;
115 (*l_ptr)[2] = -5.0; (*u_ptr)[2] = inf;
116 (*l_ptr)[3] = -5.0; (*u_ptr)[3] = 5.0;
122 ROL::Ptr<V> x = ROL::makePtr<SV>( x_ptr );
123 ROL::Ptr<V> d = ROL::makePtr<SV>( d_ptr );
124 ROL::Ptr<V> v = ROL::makePtr<SV>( v_ptr );
125 ROL::Ptr<V> l = ROL::makePtr<SV>( l_ptr );
126 ROL::Ptr<V> u = ROL::makePtr<SV>( u_ptr );
128 ROL::Ptr<const V> maskL, maskU;
133 std::vector<RealT> values(numTestVectors);
134 std::vector<RealT> exact_values(numTestVectors);
136 std::vector<ROL::Ptr<V> > x_test;
138 for(
int i=0; i<numTestVectors; ++i) {
139 x_test.push_back(x->clone());
140 RealT t =
static_cast<RealT>(i)/static_cast<RealT>(numTestVectors-1);
141 RealT fillValue = xmax*(2.0*t-1.0);
142 x_test[i]->applyUnary(ROL::Elementwise::Fill<RealT>(fillValue));
145 ROL::Ptr<OBJ> obj = ROL::makePtr<NullObjective<RealT>>();
146 ROL::Ptr<BND> bnd = ROL::makePtr<ROL::Bounds<RealT>>(l,u);
150 maskL = ipobj.getLowerMask();
151 maskU = ipobj.getUpperMask();
153 ROL::Ptr<const std::vector<RealT> > maskL_ptr =
dynamic_cast<const SV&
>(*maskL).getVector();
154 ROL::Ptr<const std::vector<RealT> > maskU_ptr =
dynamic_cast<const SV&
>(*maskU).getVector();
156 *outStream <<
"\nLower bound vector" << std::endl;
159 *outStream <<
"\nUpper bound vector" << std::endl;
162 *outStream <<
"\nLower mask vector" << std::endl;
165 *outStream <<
"\nUpper mask vector" << std::endl;
168 *outStream <<
"\nChecking Objective value" << std::endl;
170 RealT tol = std::sqrt(ROL::ROL_EPSILON<RealT>());
171 *outStream << std::setw(16) <<
"x[i], i=0,1,2,3"
172 << std::setw(20) <<
"Computed Objective"
173 << std::setw(20) <<
"Exact Objective" << std::endl;
175 RealT valueError(0.0);
177 for(
int i=0; i<numTestVectors; ++i) {
178 values[i] = ipobj.value(*(x_test[i]),tol);
183 RealT xval = x_test[i]->dot(e0);
186 for(
int j=0; j<
dim; ++j) {
187 if( (*maskL_ptr)[j] ) {
188 RealT diff = xval-(*l_ptr)[j];
189 exact_values[i] -= mu*std::log(diff);
191 if( useLinearDamping && !(*maskU_ptr)[j] ) {
192 exact_values[i] += mu*kappaD*diff;
196 if( (*maskU_ptr)[j] ) {
197 RealT diff = (*u_ptr)[j]-xval;
198 exact_values[i] -= mu*std::log(diff);
200 if(useLinearDamping && !(*maskL_ptr)[j] ) {
201 exact_values[i] += mu*kappaD*diff;
207 *outStream << std::setw(16) << xval
208 << std::setw(20) << values[i]
209 << std::setw(20) << exact_values[i] << std::endl;
210 RealT valDiff = exact_values[i] - values[i];
211 valueError += valDiff*valDiff;
214 if(valueError>ROL::ROL_EPSILON<RealT>()) {
218 *outStream <<
"\nPerforming finite difference checks" << std::endl;
220 ipobj.checkGradient(*x,*v,
true,*outStream); *outStream << std::endl;
221 ipobj.checkHessVec(*x,*d,
true,*outStream); *outStream << std::endl;
222 ipobj.checkHessSym(*x,*d,*v,
true,*outStream); *outStream << std::endl;
225 catch (std::logic_error& err) {
226 *outStream << err.what() << std::endl;
231 std::cout <<
"End Result: TEST FAILED" << std::endl;
233 std::cout <<
"End Result: TEST PASSED" << std::endl;
Provides the interface to evaluate objective functions.
void gradient(V &g, const V &x, Real &tol)
Compute gradient.
void RandomizeVector(Vector< Real > &x, const Real &lower=0.0, const Real &upper=1.0)
Fill a ROL::Vector with uniformly-distributed random numbers in the interval [lower,upper].
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
Real value(const V &x, Real &tol)
Compute value.
void printVector(const ROL::Vector< Real > &x, std::ostream &outStream)
Provides the interface to evaluate the Interior Pointy log barrier penalty function with upper and lo...
Provides the interface to apply upper and lower bound constraints.
basic_nullstream< char, char_traits< char >> nullstream
void hessVec(V &hv, const V &v, const V &x, Real &tol)
Apply Hessian approximation to vector.
int main(int argc, char *argv[])