44#include "Teuchos_GlobalMPISession.hpp"
46#include "Teuchos_GlobalMPISession.hpp"
52#include "ROL_ParameterList.hpp"
66 Real
value(
const V &x, Real &tol ) {
72 void hessVec(
V &hv,
const V &v,
const V &x, Real &tol ) {
79 ROL::Ptr<const std::vector<Real> > xp =
82 for(
size_t i=0; i<xp->size(); ++i ) {
83 outStream << (*xp)[i] << std::endl;
92int main(
int argc,
char *argv[]) {
94 typedef std::vector<RealT> vector;
101 typedef ROL::ParameterList PL;
103 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
105 int iprint = argc - 1;
106 ROL::Ptr<std::ostream> outStream;
109 outStream = ROL::makePtrFromRef(std::cout);
111 outStream = ROL::makePtrFromRef(bhs);
118 PL &iplist = parlist.sublist(
"Step").sublist(
"Primal Dual Interior Point");
119 PL &lblist = iplist.sublist(
"Barrier Objective");
122 RealT kappaD = 1.e-4;
123 bool useLinearDamping =
true;
125 lblist.set(
"Use Linear Damping", useLinearDamping);
126 lblist.set(
"Linear Damping Coefficient",kappaD);
127 lblist.set(
"Initial Barrier Parameter",mu);
129 RealT ninf = ROL::ROL_NINF<RealT>();
130 RealT inf = ROL::ROL_INF<RealT>();
133 int numTestVectors = 19;
135 ROL::Ptr<vector> x_ptr = ROL::makePtr<vector>(
dim, 0.0);
136 ROL::Ptr<vector> d_ptr = ROL::makePtr<vector>(
dim, 0.0);
137 ROL::Ptr<vector> v_ptr = ROL::makePtr<vector>(
dim, 0.0);
138 ROL::Ptr<vector> l_ptr = ROL::makePtr<vector>(
dim, 0.0);
139 ROL::Ptr<vector> u_ptr = ROL::makePtr<vector>(
dim, 0.0);
140 ROL::Ptr<vector> e0_ptr = ROL::makePtr<vector>(
dim, 0.0);
147 (*l_ptr)[0] = ninf; (*u_ptr)[0] = 5.0;
148 (*l_ptr)[1] = ninf; (*u_ptr)[1] = inf;
149 (*l_ptr)[2] = -5.0; (*u_ptr)[2] = inf;
150 (*l_ptr)[3] = -5.0; (*u_ptr)[3] = 5.0;
156 ROL::Ptr<V> x = ROL::makePtr<SV>( x_ptr );
157 ROL::Ptr<V> d = ROL::makePtr<SV>( d_ptr );
158 ROL::Ptr<V> v = ROL::makePtr<SV>( v_ptr );
159 ROL::Ptr<V> l = ROL::makePtr<SV>( l_ptr );
160 ROL::Ptr<V> u = ROL::makePtr<SV>( u_ptr );
162 ROL::Ptr<const V> maskL, maskU;
167 std::vector<RealT> values(numTestVectors);
168 std::vector<RealT> exact_values(numTestVectors);
170 std::vector<ROL::Ptr<V> > x_test;
172 for(
int i=0; i<numTestVectors; ++i) {
173 x_test.push_back(x->clone());
175 RealT fillValue = xmax*(2.0*t-1.0);
176 x_test[i]->applyUnary(ROL::Elementwise::Fill<RealT>(fillValue));
179 ROL::Ptr<OBJ> obj = ROL::makePtr<NullObjective<RealT>>();
180 ROL::Ptr<BND> bnd = ROL::makePtr<ROL::Bounds<RealT>>(l,u);
184 maskL = ipobj.getLowerMask();
185 maskU = ipobj.getUpperMask();
187 ROL::Ptr<const std::vector<RealT> > maskL_ptr =
dynamic_cast<const SV&
>(*maskL).getVector();
188 ROL::Ptr<const std::vector<RealT> > maskU_ptr =
dynamic_cast<const SV&
>(*maskU).getVector();
190 *outStream <<
"\nLower bound vector" << std::endl;
193 *outStream <<
"\nUpper bound vector" << std::endl;
196 *outStream <<
"\nLower mask vector" << std::endl;
199 *outStream <<
"\nUpper mask vector" << std::endl;
202 *outStream <<
"\nChecking Objective value" << std::endl;
204 RealT tol = std::sqrt(ROL::ROL_EPSILON<RealT>());
205 *outStream << std::setw(16) <<
"x[i], i=0,1,2,3"
206 << std::setw(20) <<
"Computed Objective"
207 << std::setw(20) <<
"Exact Objective" << std::endl;
209 RealT valueError(0.0);
211 for(
int i=0; i<numTestVectors; ++i) {
212 values[i] = ipobj.value(*(x_test[i]),tol);
217 RealT xval = x_test[i]->dot(e0);
220 for(
int j=0; j<
dim; ++j) {
221 if( (*maskL_ptr)[j] ) {
222 RealT diff = xval-(*l_ptr)[j];
223 exact_values[i] -= mu*std::log(diff);
225 if( useLinearDamping && !(*maskU_ptr)[j] ) {
226 exact_values[i] += mu*kappaD*diff;
230 if( (*maskU_ptr)[j] ) {
231 RealT diff = (*u_ptr)[j]-xval;
232 exact_values[i] -= mu*std::log(diff);
234 if(useLinearDamping && !(*maskL_ptr)[j] ) {
235 exact_values[i] += mu*kappaD*diff;
241 *outStream << std::setw(16) << xval
242 << std::setw(20) << values[i]
243 << std::setw(20) << exact_values[i] << std::endl;
244 RealT valDiff = exact_values[i] - values[i];
245 valueError += valDiff*valDiff;
248 if(valueError>ROL::ROL_EPSILON<RealT>()) {
252 *outStream <<
"\nPerforming finite difference checks" << std::endl;
254 ipobj.checkGradient(*x,*v,
true,*outStream); *outStream << std::endl;
255 ipobj.checkHessVec(*x,*d,
true,*outStream); *outStream << std::endl;
256 ipobj.checkHessSym(*x,*d,*v,
true,*outStream); *outStream << std::endl;
259 catch (std::logic_error& err) {
260 *outStream << err.what() << std::endl;
265 std::cout <<
"End Result: TEST FAILED" << std::endl;
267 std::cout <<
"End Result: TEST PASSED" << std::endl;
void gradient(V &g, const V &x, Real &tol)
Compute gradient.
void hessVec(V &hv, const V &v, const V &x, Real &tol)
Apply Hessian approximation to vector.
Real value(const V &x, Real &tol)
Compute value.
Provides the interface to apply upper and lower bound constraints.
Provides the interface to evaluate the Interior Pointy log barrier penalty function with upper and lo...
Provides the interface to evaluate objective functions.
Provides the ROL::Vector interface for scalar values, to be used, for example, with scalar constraint...
Defines the linear algebra or vector space interface.
virtual void zero()
Set to zero vector.
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,...
int main(int argc, char *argv[])
void printVector(const ROL::Vector< Real > &x, std::ostream &outStream)