MOOCHO (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DenseLinAlgPack_TestVectorClass.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
5 // Copyright (2003) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #include <iomanip>
43 #include <ostream>
44 #include <vector>
45 #include <typeinfo>
46 
52 
53 namespace {
54 
55 // This template function checks that iterator and subscriping access all
56 // give the same results.
57 template<class V, class I>
58 void test_access( V* _v, I begin, I end, std::ostream*out, bool* success ) {
59  using std::setw;
61 
62  V &v = *_v;
63  bool result;
64 
65  if(out)
66  *out << "\nbegin + i == v[i] == v(i+1), for i = 0,1,...,v.dim()\n";
67 
68  I itr;
69  int i;
70  if(out)
71  *out << "\n i, *(begin + i) == v[i] == v(i+1) == ? "
72  << "\n--, ------------ ------ ------ ------\n";
73  for( itr = begin, i = 0; itr != end; ++itr, ++i ) {
74  result = update_success( *itr == v[i] && v[i] == v(i+1), success );
75  if(out)
76  *out << setw(2) << i << ','
77  << setw(14) << *itr
78  << setw(11) << v[i]
79  << setw(11) << v(i+1)
80  << setw(11) << std::right << result << std::endl << std::left;
81  }
82  if(out) *out << std::endl;
83 }
84 
85 // This template function checks that a subregion creates the expected view.
86 // Here rng must be rng.full_range() == true.
87 template<class V, class VS>
88 void test_subregion_access( V* _v, VS* _vs, const DenseLinAlgPack::Range1D& rng
89  , std::ostream* out, bool* success )
90 {
91  using std::setw;
93 
94  bool result;
95  V &v = *_v;
96  VS &vs = *_vs;
97 
98  if(out)
99  *out << "\nv.begin() + i1 == vs.begin() + i2 == vs[i2] == vs(i2+1)"
100  << ", for i1 = lb-1,..,ub-1, for i2 = 0,..,rng.size()-1\n";
101 
102  typename V::const_iterator itr1;
103  typename VS::const_iterator itr2;
104  int i1, i2;
105  if(out)
106  *out << "\ni1, i2, *(v.begin() + i1) == *(vs.begin() + i2) == vs[i2] == vs(i2+1) == ? "
107  << "\n--, --, ----------------- ------------------ ------ -------- ------\n";
108  for( i1 = rng.lbound()-1, itr1 = v.begin() + i1, i2 = 0, itr2 = vs.begin();
109  i2 < rng.size();
110  ++i1, ++itr1, ++i2, ++itr2 )
111  {
112  result = update_success( *itr1 == *itr2 && *itr2 == vs[i2] && *itr2 == vs(i2+1), success );
113  if(out)
114  *out << setw(2) << i1 << ','
115  << setw(3) << i2 << ','
116  << setw(18) << *itr1
117  << setw(24) << *itr2
118  << setw(12) << vs[i2]
119  << setw(14) << vs(i2+1)
120  << setw(12) << std::right << result << std::endl << std::left;
121  }
122  if(out) *out << std::endl;
123 }
124 
125 // Print out a string for overlap
126 const char* overlap_str( DenseLinAlgPack::EOverLap overlap ) {
127  switch(overlap) {
129  return "NO_OVERLAP";
131  return "SOME_OVERLAP";
133  return "SAME_MEM";
134  }
135  return "Invalid value for EOverLap";
136 }
137 
138 } // end namespace
139 
141 {
142  using DenseLinAlgPack::comp;
145 
146  bool success = true;
147  bool result, result1, result2;
148 
149  if(out)
150  *out << "\n**********************************************"
151  << "\n*** Testing DVector and DVectorSlice classes ***"
152  << "\n**********************************************\n"
153  << std::boolalpha;
154 
155  try {
156 
157  if(out)
158  *out << "\nLet vvz[i-1] = i + 0.1*i, for i = 1...,n\n";
159 
160  // std::vector<> which is starndard for comparisons
161  const DVector::size_type n = 6;
162  std::vector<DVector::value_type> vvz(6);
163  {for(int i = 1; i <= n; ++i)
164  vvz[i-1] = i + 0.1 * i;
165  }
166 
167  if(out) *out << "\nLet alpha1 = 22.5\n";
168  const value_type alpha1 = 22.5;
169 
170  // ///////////////////////
171  // Test Constructors
172 
173  if(out)
174  *out << "\n***\n*** Testing constructors\n***\n";
175 
176  // DVectorSlice Constructors
177 
178  if(out) *out << "\nVectorSlice vs1\n";
179  DVectorSlice vs1;
180  if(out) *out << "vs1 =\n" << vs1;
181 
182  if(out) *out << "\nVectorSlice vs2(vvz.begin(),n)\n";
183  DVectorSlice vs2(&vvz[0],n);
184  if(out) *out << "vs2 =\n" << vs2;
185 
186  if(out) *out << "\nVectorSlice vs3(vvz.begin(),n,Range1D())\n";
187  DVectorSlice vs3(&vvz[0],n,Range1D());
188  if(out) *out << "vs3 =\n" << vs3;
189 
190  if(out) *out << "\nVectorSlice vs4(vs3,Range1D())\n";
191  DVectorSlice vs4(vs3,Range1D());
192  if(out) *out << "vs4 =\n" << vs4;
193 
194  // DVector Constructors
195 
196  if(out) *out << "\nVector v1\n";
197  DVector v1;
198  if(out) *out << "v1 =\n" << v1;
199 
200  if(out) *out << "\nVector v2(alpha1,n)\n";
201  DVector v2(alpha1,n);
202  if(out) *out << "v2 =\n" << v2;
203 
204  if(out) *out << "\nVector v3(vvz.begin(),n)\n";
205  DVector v3(&vvz[0],n);
206  if(out) *out << "v3 =\n" << v3;
207 
208  if(out) *out << "\nVector v4(vs4)\n";
209  DVector v4(vs4);
210  if(out) *out << "v4 =\n" << v4;
211 
212  // //////////////////////////////////////
213  // Test Binding Views (DVectorSlice)
214 
215  if(out)
216  *out << "\n***\n*** Testing DVectorSlice binding and "
217  "conversion from DVector -> DVectorSlice\n***\n"
218  << "\nvs1.bind(v2());\n";
219  vs1.bind(v2());
220  if(out) *out << "vs1 =\n" << vs1;
221 
222  // ///////////////////////
223  // Test DVector Resizing
224 
225  if(out)
226  *out << "\n***\n*** Testing DVector resizing\n***\n"
227  << "\nv2.free();\n";
228  v2.free();
229  if(out) *out << "v2.dim() == 0 : " << update_success( v2.dim() == 0, &success ) << std::endl;
230 
231  if(out)
232  *out << "\nv2.resize(n,2*alpha1);\n";
233  v2.resize(n,2*alpha1);
234  result1 = update_success( v2.dim() == n, &success );
235  result2 = update_success( comp(v2,2*alpha1), &success );
236  if(out)
237  *out << "( (v2.dim() -> " << v2.dim() << ") == n && "
238  << "(comp(v2,2*alpha1) -> " << result2 << ") ) : " << (result1 && result2) << std::endl
239  << "v2 =\n" << v2;
240 
241  // //////////////////////////////////////////////////////////
242  // Test Iterator Access, Subscriping and Reverse VectorSlices
243 
244  if(out) *out << "\n***\n*** Testing Iterator Access, Subscriping and Reverse VectorSlices\n***\n";
245 
246  if(out)
247  *out << "\nLet v == v3, begin = v.begin()";
248  test_access(&v3,v3.begin(),v3.end(),out,&success);
249 
250  if(out)
251  *out << "\nLet v == const_cast<const DVector&>(v3), begin = v.begin()";
252  test_access(&const_cast<const DVector&>(v3),const_cast<const DVector&>(v3).begin()
253  ,const_cast<const DVector&>(v3).end(),out,&success);
254 
255  if(out)
256  *out << "\nLet v == vs3, begin = v.begin()";
257  test_access(&vs3,vs3.begin(),vs3.end(),out,&success);
258 
259  if(out)
260  *out << "\nLet v == const_cast<const DVectorSlice&>(vs3), begin = v.begin()";
261  test_access(&const_cast<const DVectorSlice&>(vs3),const_cast<const DVectorSlice&>(vs3).begin()
262  ,const_cast<const DVectorSlice&>(vs3).end(),out,&success);
263 
264  if(out)
265  *out << "\nLet v == v3.rev(), begin = v3.rbegin()";
266  test_access(&v3.rev(),v3.rbegin(),v3.rend(),out,&success);
267 
268  if(out)
269  *out << "\nLet v == const_cast<const DVector&>(v3).rev(), begin = const_cast<const DVector&>(v3).rbegin()";
270  test_access(&const_cast<const DVector&>(v3).rev(),const_cast<const DVector&>(v3).rbegin()
271  ,const_cast<const DVector&>(v3).rend(),out,&success);
272 
273  if(out)
274  *out << "\nLet v == vs3.rev(), begin = vs3.rbegin()";
275  test_access(&vs3.rev(),vs3.rbegin(),vs3.rend(),out,&success);
276 
277  if(out)
278  *out << "\nLet v == const_cast<const DVectorSlice&>(vs3).rev()"
279  ", begin = const_cast<const DVectorSlice&>(vs3).rbegin()";
280  test_access(&const_cast<const DVectorSlice&>(vs3).rev(),const_cast<const DVectorSlice&>(vs3).rbegin()
281  ,const_cast<const DVectorSlice&>(vs3).rend(),out,&success);
282 
283 #ifdef LINALGPACK_CHECK_RANGE
284 
285  if(out) *out << "\n*** Test subscriping out of bounds\n";
286 
287  if(out) *out << "\nv3(0); (should throw std::out_of_range)\n";
288  try{
289  v3(0);
290  result = false;
291  }
292  catch(const std::out_of_range&) {
293  result = true;
294  }
295  if(out) *out << "v3(0) threw std::out_of_range : " << result << std::endl;
296  update_success( result, &success );
297 
298  if(out) *out << "\nv3(n+1); (should throw std::out_of_range)\n";
299  try{
300  v3(n+1);
301  result = false;
302  }
303  catch(const std::out_of_range&) {
304  result = true;
305  }
306  if(out) *out << "v3(n+1) threw std::out_of_range : " << result << std::endl;
307  update_success( result, &success );
308 
309  if(out) *out << "\nvs3(0); (should throw std::out_of_range)\n";
310  try{
311  vs3(0);
312  result = false;
313  }
314  catch(const std::out_of_range &) {
315  result = true;
316  }
317  if(out) *out << "vs3(0) threw std::out_of_range : " << result << std::endl;
318  update_success( result, &success );
319 
320  if(out) *out << "\nvs3(n+1); (should throw std::out_of_range)\n";
321  try{
322  vs3(n+1);
323  result = false;
324  }
325  catch(const std::out_of_range&) {
326  result = true;
327  }
328  if(out) *out << "vs3(n+1) threw std::out_of_range : " << result << std::endl;
329  update_success( result, &success );
330 
331 #endif
332 
333  // ////////////////////////////////
334  // Test Subregion Access
335 
336  if(out) *out << "\n***\n*** Testing Subregion Access\n***\n";
337 
338  // DVector Subregions
339  Range1D rng;
340 
341  if(out) *out << "\nv = v3, rng = [1,n/2], vs = v(rng)";
342  rng = Range1D(1,n/2);
343  test_subregion_access( &v3, &v3(rng), rng, out, &success );
344 
345  if(out) *out << "\nv = v3, rng = [n/3,2*n/3], vs = v(rng)";
346  rng = Range1D(n/3,2*n/3);
347  test_subregion_access( &v3, &v3(rng), rng, out, &success );
348 
349  if(out) *out << "\nv = v3, rng = [n/2,n], vs = v(rng)";
350  rng = Range1D(n/2,n);
351  test_subregion_access( &v3, &v3(rng), rng, out, &success );
352 
353  if(out) *out << "\nv = const_cast<const DVector&>(v3), rng = [n/2,n], vs = v(rng)";
354  rng = Range1D(n/2,n);
355  test_subregion_access( &v3, &const_cast<const DVector&>(v3)(rng), rng, out, &success );
356 
357  if(out) *out << "\nv = v3, rng = [1,n/2], vs = v(1,n/2)";
358  rng = Range1D(1,n/2);
359  test_subregion_access( &v3, &v3(1,n/2), rng, out, &success );
360 
361  if(out) *out << "\nv = const_cast<const DVector&>(v3), rng = [n/2,n], vs = v(n/2,n)";
362  rng = Range1D(n/2,n);
363  test_subregion_access( &v3, &const_cast<const DVector&>(v3)(n/2,n), rng, out, &success );
364 
365  // DVectorSlice Subregions
366 
367  if(out) *out << "\nv = vs3, rng = [1,n/2], vs = v(rng)";
368  rng = Range1D(1,n/2);
369  test_subregion_access( &vs3, &vs3(rng), rng, out, &success );
370 
371  if(out) *out << "\nv = vs3, rng = [n/3,2*n/3], vs = v(rng)";
372  rng = Range1D(n/3,2*n/3);
373  test_subregion_access( &vs3, &vs3(rng), rng, out, &success );
374 
375  if(out) *out << "\nv = vs3, rng = [n/2,n], vs = v(rng)";
376  rng = Range1D(n/2,n);
377  test_subregion_access( &vs3, &vs3(rng), rng, out, &success );
378 
379  if(out) *out << "\nv = const_cast<const DVectorSlice&>(vs3), rng = [n/2,n], vs = v(rng)";
380  rng = Range1D(n/2,n);
381  test_subregion_access( &vs3, &const_cast<const DVectorSlice&>(vs3)(rng), rng, out, &success );
382 
383  if(out) *out << "\nv = vs3, rng = [1,n/2], vs = v(1,n/2)";
384  rng = Range1D(1,n/2);
385  test_subregion_access( &vs3, &vs3(1,n/2), rng, out, &success );
386 
387  if(out) *out << "\nv = const_cast<const DVectorSlice&>(vs3), rng = [n/2,n], vs = v(n/2,n)";
388  rng = Range1D(n/2,n);
389  test_subregion_access( &vs3, &const_cast<const DVectorSlice&>(vs3)(n/2,n), rng, out, &success );
390 
391  // ///////////////////////
392  // Test Assignment
393 
394  if(out) *out << "\n***\n*** Testing assignment operators\n***\n";
395 
396  // DVector Assignment
397 
398  if(out) *out << "\nv1.resize(n); v1 = 0.0;\n";
399  v1.resize(n);
400  v1 = 0.0;
401  result = update_success( comp( v1, 0.0 ), &success );
402  if(out)
403  *out << "v1 =\n" << v1
404  << "v1 == 0.0 : " << result << std::endl;
405 
406  if(out) *out << "\nv1 = 0.0; v1 = vs3;\n";
407  v1 = 0.0;
408  v1 = vs3;
409  result = update_success( comp( v1, vs3 ), &success );
410  if(out)
411  *out << "v1 =\n" << v1
412  << "v1 == vs3 : " << result << std::endl;
413 
414  if(out) *out << "\nv1 = 0.0; v1 = v3;\n";
415  v1 = 0.0;
416  v1 = v3;
417  result = update_success( comp( v1, v3 ), &success );
418  if(out)
419  *out << "v1 =\n" << v1
420  << "v1 == v3 : " << result << std::endl;
421 
422  // DVectorSlice Assignment
423 
424  if(out) *out << "\nv1.resize(n); v1 = 0.0; vs1.bind(v1());\n";
425  v1.resize(n);
426  v1 = 0.0;
427  vs1.bind(v1());
428 
429  if(out) *out << "\nvs1 = alpha1;\n";
430  vs1 = alpha1;
431  result = update_success( comp( vs1, alpha1 ), &success );
432  if(out)
433  *out << "vs1 =\n" << v1
434  << "vs1 == alpha1 : " << result << std::endl;
435 
436  if(out) *out << "\nvs1 = 0.0; vs1 = vs3;\n";
437  vs1 = 0.0;
438  vs1 = vs3;
439  result = update_success( comp( vs1, vs3 ), &success );
440  if(out)
441  *out << "vs1 =\n" << vs1
442  << "vs1 == vs3 : " << result << std::endl;
443 
444  // ////////////////////////
445  // Test overlap()
446 
447  if(out) *out << "\n***\n*** Testing overlap\n***\n";
448 
449  EOverLap ovlap;
450 
451  // DVector overlap
452 
453  if(out) *out << "\n*** DVector overlap\n";
454 
455  if(out) *out << "(v1.overlap(v3) -> ";
456  ovlap = v1.overlap(v3);
457  result = update_success( ovlap == NO_OVERLAP, &success );
458  if(out) *out << overlap_str(ovlap) << ") == NO_OVERLAP : " << result << std::endl;
459 
460  if(out) *out << "(v3.overlap(v3) -> ";
461  ovlap = v3.overlap(v3);
462  result = update_success( ovlap == SAME_MEM, &success );
463  if(out) *out << overlap_str(ovlap) << ") == SAME_MEM : " << result << std::endl;
464 
465  if(out) *out << "(v3.overlap(v3(1,n-1)) -> ";
466  ovlap = v3.overlap(v3(1,n-1));
467  result = update_success( ovlap == SOME_OVERLAP, &success );
468  if(out) *out << overlap_str(ovlap) << ") == SOME_OVERLAP : " << result << std::endl;
469 
470  if(out) *out << "(v3.overlap(v3(2,n-1)) -> ";
471  ovlap = v3.overlap(v3(2,n-1));
472  result = update_success( ovlap == SOME_OVERLAP, &success );
473  if(out) *out << overlap_str(ovlap) << ") == SOME_OVERLAP : " << result << std::endl;
474 
475  // DVectorSlice overlap
476 
477  if(out) *out << "\n*** DVectorSlice overlap\n"
478  << "vs1.bind(v3());\n";
479 
480  vs1.bind(v3());
481 
482  if(out) *out << "(vs1.overlap(v2) -> ";
483  ovlap = vs1.overlap(v2);
484  result = update_success( ovlap == NO_OVERLAP, &success );
485  if(out) *out << overlap_str(ovlap) << ") == NO_OVERLAP : " << result << std::endl;
486 
487  if(out) *out << "(vs1.overlap(vs1) -> ";
488  ovlap = vs1.overlap(vs1);
489  result = update_success( ovlap == SAME_MEM, &success );
490  if(out) *out << overlap_str(ovlap) << ") == SAME_MEM : " << result << std::endl;
491 
492  if(out) *out << "(vs1(1,n/2).overlap(vs1(n/2+1,n)) -> ";
493  ovlap = vs1(1,n/2).overlap(vs1(n/2+1,n));
494  result = update_success( ovlap == NO_OVERLAP, &success );
495  if(out) *out << overlap_str(ovlap) << ") == NO_OVERLAP : " << result << std::endl;
496 
497  if(out) *out << "(vs1(1,n/2).overlap(vs1(n/3,2*n/3)) -> ";
498  ovlap = vs1(1,n/2).overlap(vs1(n/3,2*n/3));
499  result = update_success( ovlap == SOME_OVERLAP, &success );
500  if(out) *out << overlap_str(ovlap) << ") == SOME_OVERLAP : " << result << std::endl;
501 
502  if(out) *out << "(vs1(n/3,2*n/3).overlap(vs1(n/2+1,n)) -> ";
503  ovlap = vs1(n/3,2*n/3).overlap(vs1(n/2+1,n));
504  result = update_success( ovlap == SOME_OVERLAP, &success );
505  if(out) *out << overlap_str(ovlap) << ") == SOME_OVERLAP : " << result << std::endl;
506 
507  if(out) *out << "(vs1(n/3,2*n/3).overlap(vs1(n/3,2*n/3)) -> ";
508  ovlap = vs1(n/3,2*n/3).overlap(vs1(n/3,2*n/3));
509  result = update_success( ovlap == SAME_MEM, &success );
510  if(out) *out << overlap_str(ovlap) << ") == SAME_MEM : " << result << std::endl;
511 
512  } // end try
513  catch( const std::exception& excpt ) {
514  success = false;
515  if(out)
516  (*out) << "\nError, a standard exception was thrown: "
517  << typeName(excpt) << ": "
518  << excpt.what() << std::endl;
519  }
520  catch(...) {
521  success = false;
522  if(out)
523  (*out) << "\nError, an unknown exception was thrown\n";
524  }
525 
526  if(out) {
527  if(success)
528  (*out)
529  << "\n*** Congradulations, DVector and DVectorSlice seem to check out. ***\n";
530  else
531  (*out)
532  << "\n*** Oops, all of the tests for DVector and DVectorSlice "
533  "where not successful. ***\n";
534  }
535 
536  return success;
537 }
538 
std::string typeName(const T &t)
Index size() const
Return the size of the range (ubound() - lbound() + 1)
. One-based subregion index range class.
std::ostream * out
EOverLap
Enumeration for returning the amount of overlap between two objects.
Index lbound() const
Return lower bound of the range.
FortranTypes::f_dbl_prec value_type
Typedef for the value type of elements that is used for the library.
bool update_success(bool result_check, bool *success)
Helper function for updating a flag for if an operation returned false.
bool comp(const DVectorSlice &vs1, const DVectorSlice &vs2)
int n