ROL
function/operator/test_01.cpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
28 #include "ROL_Stream.hpp"
29 #include "Teuchos_GlobalMPISession.hpp"
30 
31 typedef double RealT;
32 
33 int main(int argc, char *argv[]) {
34 
35 
36 
37  typedef std::vector<RealT> vector;
38 
39  typedef ROL::StdVector<RealT> SV;
40 
41  typedef ROL::StdLinearOperator<RealT> StdLinearOperator;
42 
43 
44  Teuchos::GlobalMPISession mpiSession(&argc, &argv);
45 
46  // This little trick lets us print to std::cout only if a (dummy) command-line argument is provided.
47  int iprint = argc - 1;
48  ROL::Ptr<std::ostream> outStream;
49  ROL::nullstream bhs; // outputs nothing
50  if (iprint > 0)
51  outStream = ROL::makePtrFromRef(std::cout);
52  else
53  outStream = ROL::makePtrFromRef(bhs);
54 
55  // Save the format state of the original std::cout.
56  ROL::nullstream oldFormatState;
57  oldFormatState.copyfmt(std::cout);
58 
59  int errorFlag = 0;
60 
61  // *** Test body.
62 
63  try {
64 
65  ROL::Ptr<vector> a_ptr = ROL::makePtr<vector>(
66  std::initializer_list<RealT>{4.0,2.0,1.0,3.0});
67  ROL::Ptr<vector> ai_ptr = ROL::makePtr<vector>(
68  std::initializer_list<RealT>{3.0/10.0, -2.0/10.0, -1.0/10.0, 4.0/10.0});
69 
70  ROL::Ptr<vector> x1_ptr = ROL::makePtr<vector>(
71  std::initializer_list<RealT>{1.0,-1.0});
72  ROL::Ptr<vector> b1_ptr = ROL::makePtr<vector>(2);
73 
74  ROL::Ptr<vector> x2_ptr = ROL::makePtr<vector>(2);
75  ROL::Ptr<vector> b2_ptr = ROL::makePtr<vector>(
76  std::initializer_list<RealT>{3.0,-1.0});
77 
78  ROL::Ptr<vector> y3_ptr = ROL::makePtr<vector>(
79  std::initializer_list<RealT>{-2.0,1.0});
80  ROL::Ptr<vector> c3_ptr = ROL::makePtr<vector>(2);
81 
82  ROL::Ptr<vector> y4_ptr = ROL::makePtr<vector>(2);
83  ROL::Ptr<vector> c4_ptr = ROL::makePtr<vector>(
84  std::initializer_list<RealT>{-6.0,1.0});
85 
86  StdLinearOperator A(a_ptr);
87  StdLinearOperator Ai(ai_ptr);
88 
89  SV x1(x1_ptr); SV x2(x2_ptr); SV y3(y3_ptr); SV y4(y4_ptr);
90  SV b1(b1_ptr); SV b2(b2_ptr); SV c3(c3_ptr); SV c4(c4_ptr);
91 
92  RealT tol = ROL::ROL_EPSILON<RealT>();
93 
94  // Test 1
95  *outStream << "\nTest 1: Matrix multiplication" << std::endl;
96  A.apply(b1,x1,tol);
97  *outStream << "x = [" << (*x1_ptr)[0] << "," << (*x1_ptr)[1] << "]" << std::endl;
98  *outStream << "b = [" << (*b1_ptr)[0] << "," << (*b1_ptr)[1] << "]" << std::endl;
99  b1.axpy(-1.0,b2);
100 
101  RealT error1 = b1.norm();
102  errorFlag += error1 > tol;
103  *outStream << "Error = " << error1 << std::endl;
104 
105  // Test 2
106  *outStream << "\nTest 2: Linear solve" << std::endl;
107  A.applyInverse(*x2_ptr,*b2_ptr,tol);
108  *outStream << "x = [" << (*x2_ptr)[0] << "," << (*x2_ptr)[1] << "]" << std::endl;
109  *outStream << "b = [" << (*b2_ptr)[0] << "," << (*b2_ptr)[1] << "]" << std::endl;
110  x2.axpy(-1.0,x1);
111 
112  RealT error2 = x2.norm();
113  errorFlag += error2 > tol;
114  *outStream << "Error = " << error2 << std::endl;
115 
116  // Test 3
117  *outStream << "\nTest 3: Transposed matrix multiplication" << std::endl;
118  A.applyAdjoint(*c3_ptr,*y3_ptr,tol);
119  *outStream << "y = [" << (*y3_ptr)[0] << "," << (*y3_ptr)[1] << "]" << std::endl;
120  *outStream << "c = [" << (*c3_ptr)[0] << "," << (*c3_ptr)[1] << "]" << std::endl;
121  c3.axpy(-1.0,c4);
122 
123  RealT error3 = c3.norm();
124  errorFlag += error3 > tol;
125  *outStream << "Error = " << error3 << std::endl;
126 
127  // Test 4
128  *outStream << "\nTest 4: Linear solve with transpose" << std::endl;
129  A.applyAdjointInverse(y4,c4,tol);
130  *outStream << "y = [" << (*y4_ptr)[0] << "," << (*y4_ptr)[1] << "]" << std::endl;
131  *outStream << "c = [" << (*c4_ptr)[0] << "," << (*c4_ptr)[1] << "]" << std::endl;
132  y4.axpy(-1.0,y3);
133 
134  RealT error4 = y4.norm();
135  errorFlag += error4 > tol;
136  *outStream << "Error = " << error4 << std::endl;
137 
138  *outStream << "x1 = "; x1.print(*outStream);
139  Ai.applyInverse(b1,x1,tol);
140  *outStream << "b1 = "; b1.print(*outStream);
141  A.apply(b1,x1,tol);
142  *outStream << "b1 = "; b1.print(*outStream);
143  A.applyInverse(x1,b1,tol);
144  *outStream << "x1 = "; x1.print(*outStream);
145  Ai.apply(x1,b1,tol);
146  *outStream << "x1 = "; x1.print(*outStream);
147 
148 
149  }
150  catch (std::logic_error& err) {
151  *outStream << err.what() << "\n";
152  errorFlag = -1000;
153  }; // end try
154 
155  if (errorFlag != 0)
156  std::cout << "End Result: TEST FAILED\n";
157  else
158  std::cout << "End Result: TEST PASSED\n";
159 
160  // reset format state of std::cout
161  std::cout.copyfmt(oldFormatState);
162 
163  return 0;
164 
165 }
166 
Provides the std::vector implementation to apply a linear operator, which is a std::vector representa...
Defines a no-output stream class ROL::NullStream and a function makeStreamPtr which either wraps a re...
Provides the ROL::Vector interface for scalar values, to be used, for example, with scalar constraint...
basic_nullstream< char, char_traits< char >> nullstream
Definition: ROL_Stream.hpp:38
int main(int argc, char *argv[])