Teuchos Package Browser (Single Doxygen Collection)  Version of the Day
numerics/test/DenseMatrix/cxx_main.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Teuchos: Common Tools Package
5 // Copyright (2004) 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 Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
45 #include "Teuchos_Version.hpp"
46 
47 #define OTYPE int
48 #define STYPE std::complex<double>
49 
50 template<typename TYPE>
51 int PrintTestResults(std::string, TYPE, TYPE, bool);
52 
53 int ReturnCodeCheck(std::string, int, int, bool);
54 
55 typedef double Real;
58 //typedef Teuchos::SerialDenseVector<OTYPE, STYPE> DVector;
59 
60 int main(int argc, char* argv[])
61 {
62 
63  int i, j;
64  bool verbose = 0;
65  if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
66 
67  if (verbose)
68  std::cout << Teuchos::Teuchos_Version() << std::endl << std::endl;
69 
70  int numberFailedTests = 0;
71  int returnCode = 0;
72  std::string testName = "";
73 
74 
75 
76  if (verbose) std::cout<<std::endl<<"********** CHECKING TEUCHOS SERIAL DENSE MATRIX **********"<<std::endl<<std::endl;
77 
78  // default constructor test
79  DMatrix DefConTest;
80  if (verbose) std::cout <<"default constructor -- construct empty matrix ";
81  if ( DefConTest.values()!=NULL || DefConTest.numCols()!=0 || DefConTest.numRows()!=0 ||DefConTest.stride()!=0 ||DefConTest.empty()!=true ) {
82  if (verbose) std::cout << "unsuccessful."<<std::endl;
83  numberFailedTests++;
84  } else {
85  if (verbose) std::cout << "successful."<<std::endl;
86  }
87 
88  // constructor 1 (matrix w/ dimension but empty)
89 
90  DMatrix Con1Test( 3, 4 );
91  if (verbose) std::cout <<"constructor 1 -- empty matrix with given dimensions ";
92  if ( Con1Test.numRows()!=3 || Con1Test.numCols()!=4 || Con1Test( 1, 2 )!=0.0 ) {
93  if (verbose) std::cout << "unsuccessful."<<std::endl;
94  numberFailedTests++;
95  } else {
96  if (verbose) std::cout << "successful."<<std::endl;
97  }
98 
99  // constructor 2 (from array) tests
100 
101  STYPE a[9];
102  for(i = 0; i < 9; i++)
103  {
104  a[i] = i;
105  }
106  DMatrix Con2Test1ExpRes;
107  Con2Test1ExpRes.shape(2, 3);
108  Con2Test1ExpRes(0, 0) = 0; Con2Test1ExpRes(0, 1) = 2; Con2Test1ExpRes(0, 2) = 4;
109  Con2Test1ExpRes(1, 0) = 1; Con2Test1ExpRes(1, 1) = 3; Con2Test1ExpRes(1, 2) = 5;
110 
111  DMatrix Con2Test1(Teuchos::Copy, a, 2, 2, 3);
112  numberFailedTests += PrintTestResults("constructor 2 -- construct matrix from array subrange", Con2Test1, Con2Test1ExpRes, verbose);
113 
114 
115  // constructor 3 (copy constructor)
116 
117  DMatrix Con3TestCopy( Con2Test1ExpRes );
118  if(verbose) std::cout <<"constructor 3 -- copy constructor ";
119  if ( Con3TestCopy != Con2Test1ExpRes ) {
120  if (verbose) std::cout << "unsuccessful."<<std::endl;
121  numberFailedTests++;
122  } else {
123  if (verbose) std::cout << "successful."<<std::endl;
124  }
125 
126  DMatrix Con3TestCopyTrans( Con2Test1ExpRes, Teuchos::TRANS );
127  if(verbose) std::cout <<"constructor 3 -- copy constructor (transposed) ";
128  if ( Con3TestCopyTrans(2, 0) != Con2Test1ExpRes(0, 2) ) {
129  if (verbose) std::cout << "unsuccessful."<<std::endl;
130  numberFailedTests++;
131  } else {
132  if (verbose) std::cout << "successful."<<std::endl;
133  }
134 
135  // constructor 4 (submatrix)
136 
137  DMatrix Con4TestOrig(Teuchos::Copy, a, 3, 3, 3);
138  DMatrix Con4TestSubmatrix;
139  Con4TestSubmatrix.shape(2, 2);
140  Con4TestSubmatrix(0, 0) = 4; Con4TestSubmatrix(0, 1) = 7;
141  Con4TestSubmatrix(1, 0) = 5; Con4TestSubmatrix(1, 1) = 8;
142  DMatrix Con4TestCopy1(Teuchos::Copy, Con4TestOrig, 2, 2, 1, 1);
143  numberFailedTests += PrintTestResults("constructor 4 -- submatrix copy", Con4TestCopy1, Con4TestSubmatrix, verbose);
144  DMatrix Con4TestCopy2(Teuchos::Copy, Con4TestOrig, 3, 3, 0, 0);
145  numberFailedTests += PrintTestResults("constructor 4 -- full matrix copy", Con4TestCopy2, Con4TestOrig, verbose);
146  DMatrix Con4TestView1(Teuchos::View, Con4TestOrig, 2, 2, 1, 1);
147  numberFailedTests += PrintTestResults("constructor 4 -- full matrix view", Con4TestView1, Con4TestSubmatrix, verbose);
148  DMatrix Con4TestView2(Teuchos::View, Con4TestOrig, 3, 3, 0, 0);
149  numberFailedTests += PrintTestResults("constructor 4 -- submatrix view", Con4TestView2, Con4TestOrig, verbose);
150 
151  // Norm Tests
152 
153  DMatrix AAA;
154  AAA.shape(3, 3);
155  AAA(0, 0) = 1; AAA(0, 1) = 2; AAA(0, 2) = 3;
156  AAA(1, 0) = 4; AAA(1, 1) = 5; AAA(1, 2) = 6;
157  AAA(2, 0) = 7; AAA(2, 1) = 8; AAA(2, 2) = 9;
158  DMatrix BBB;
159  numberFailedTests += PrintTestResults("normOne of a 3x3", AAA.normOne(), 18.0, verbose);
160  numberFailedTests += PrintTestResults("normInf of a 3x3", AAA.normInf(), 24.0, verbose);
162  numberFailedTests += PrintTestResults("normFrobenius of a 3x3", AAA.normFrobenius(), 3.0, verbose);
163  numberFailedTests += PrintTestResults("normOne of a 0x0", BBB.normOne(), 0.0, verbose);
164  numberFailedTests += PrintTestResults("normInf of a 0x0", BBB.normInf(), 0.0, verbose);
165  numberFailedTests += PrintTestResults("normFrobenius of a 0x0", BBB.normFrobenius(), 0.0, verbose);
166 
167  // multiply() -- dimensions tests
168 
169  DMatrix DimTest0x0A, DimTest0x0B, DimTest2x0, DimTest1x2, DimTest2x1, DimTest2x2A, DimTest2x2B,
170  DimTest3x3, DimTest0x2, DimTest0x0Result, DimTest1x1Result, DimTest2x0Result, DimTest1x2Result, DimTest2x1Result, DimTest2x2Result,
171  DimTest2x3Result, DimTest0x2Result, DimTest3x3Result;
172 
173  DimTest0x2.shape(0, 2);
174  DimTest2x0.shape(2, 0);
175  DimTest1x2.shape(1, 2);
176  DimTest2x1.shape(2, 1);
177  DimTest2x2A.shape(2, 2);
178  DimTest2x2B.shape(2, 2);
179  DimTest3x3.shape(3, 3);
180  DimTest0x2Result.shape(0, 2);
181  DimTest1x1Result.shape(1, 1);
182  DimTest2x0Result.shape(2, 0);
183  DimTest1x2Result.shape(1, 2);
184  DimTest2x1Result.shape(2, 1);
185  DimTest2x2Result.shape(2, 2);
186  DimTest2x3Result.shape(2, 3);
187  DimTest3x3Result.shape(3, 3);
188 
189  returnCode = DimTest2x2Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x2A, DimTest2x2B, 1);
190  testName = "multiply() -- dimensions -- compatible square matrices";
191  numberFailedTests += ReturnCodeCheck(testName, returnCode, 0, verbose);
192  returnCode = DimTest2x3Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x2A, DimTest3x3, 1);
193  testName = "multiply() -- dimensions -- incompatible square matrices";
194  numberFailedTests += ReturnCodeCheck(testName, returnCode, 1, verbose);
195  returnCode = DimTest1x1Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest1x2, DimTest2x1, 1);
196  testName = "multiply() -- dimensions -- compatible nonsquare matrices";
197  numberFailedTests += ReturnCodeCheck(testName, returnCode, 0, verbose);
198  returnCode = DimTest2x2Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x1, DimTest1x2, 1);
199  testName = "multiply() -- dimensions -- compatible nonsquare matrices";
200  numberFailedTests += ReturnCodeCheck(testName, returnCode, 0, verbose);
201  returnCode = DimTest2x1Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x1, DimTest2x1, 1);
202  testName = "multiply() -- dimensions -- incompatible nonsquare matrices";
203  numberFailedTests += ReturnCodeCheck(testName, returnCode, 1, verbose);
204  returnCode = DimTest1x2Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest1x2, DimTest1x2, 1);
205  testName = "multiply() -- dimensions -- incompatible nonsquare matrices";
206  numberFailedTests += ReturnCodeCheck(testName, returnCode, 1, verbose);
207  returnCode = DimTest2x2Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x0, DimTest2x2A, 1);
208  testName = "multiply() -- dimensions -- first operand bad numCols";
209  numberFailedTests += ReturnCodeCheck(testName, returnCode, 1, verbose);
210  returnCode = DimTest0x2Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x2A, DimTest0x2, 1);
211  testName = "multiply() -- dimensions -- second operand bad numRows";
212  numberFailedTests += ReturnCodeCheck(testName, returnCode, 1, verbose);
213  returnCode = DimTest2x2Result.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, DimTest2x2A, DimTest2x0, 1);
214  testName = "multiply() -- dimensions -- second operand bad numCols";
215  numberFailedTests += ReturnCodeCheck(testName, returnCode, 1, verbose);
216 
217  // multiply() -- multiplication results tests
218 
219  DMatrix MultTest2x2A, MultTest2x2B, MultTest3x3A, MultTest3x3B, MultTest2x2ATimes2x2B,
220  MultTest3x3ATimes3x3B, MultTest2x2BTimes2x2A, MultTest3x3BTimes3x3A, MultTest2x2ATimes2x2BExpResult, MultTest2x2BTimes2x2AExpResult,
221  MultTest3x3ATimes3x3BExpResult, MultTest3x3BTimes3x3AExpResult, MultTest2x3A, MultTest2x3B, MultTest3x2A, MultTest3x2B,
222  MultTest2x3ATimes3x2B, MultTest3x2ATimes2x3B, MultTest2x3BTimes3x2A, MultTest3x2BTimes2x3A, MultTest2x3ATimes3x2BExpResult,
223  MultTest3x2ATimes2x3BExpResult, MultTest2x3BTimes3x2AExpResult, MultTest3x2BTimes2x3AExpResult;
224 
225  MultTest2x2A.shape(2, 2);
226  MultTest2x2B.shape(2, 2);
227  MultTest3x3A.shape(3, 3);
228  MultTest3x3B.shape(3, 3);
229  MultTest2x2ATimes2x2B.shape(2, 2);
230  MultTest2x2BTimes2x2A.shape(2, 2);
231  MultTest3x3ATimes3x3B.shape(3, 3);
232  MultTest3x3BTimes3x3A.shape(3, 3);
233  MultTest2x2ATimes2x2BExpResult.shape(2, 2);
234  MultTest2x2BTimes2x2AExpResult.shape(2, 2);
235  MultTest3x3ATimes3x3BExpResult.shape(3, 3);
236  MultTest3x3BTimes3x3AExpResult.shape(3, 3);
237  MultTest2x3A.shape(2, 3);
238  MultTest2x3B.shape(2, 3);
239  MultTest3x2A.shape(3, 2);
240  MultTest3x2B.shape(3, 2);
241  MultTest2x3ATimes3x2B.shape(2, 2);
242  MultTest3x2ATimes2x3B.shape(3, 3);
243  MultTest2x3BTimes3x2A.shape(2, 2);
244  MultTest3x2BTimes2x3A.shape(3, 3);
245  MultTest2x3ATimes3x2BExpResult.shape(2, 2);
246  MultTest3x2ATimes2x3BExpResult.shape(3, 3);
247  MultTest2x3BTimes3x2AExpResult.shape(2, 2);
248  MultTest3x2BTimes2x3AExpResult.shape(3, 3);
249 
250  for(i = 0; i < 2; i++)
251  {
252  for(j = 0; j < 2; j++)
253  {
254  MultTest2x2A(i, j) = i + j;
255  MultTest2x2B(i, j) = (i * j) + 1;
256  }
257  }
258  for(i = 0; i < 3; i++)
259  {
260  for(j = 0; j < 3; j++)
261  {
262  MultTest3x3A(i, j) = i + j;
263  MultTest3x3B(i, j) = (i * j) + 1;
264  }
265  }
266 
267  MultTest2x2ATimes2x2BExpResult(0, 0) = 1; MultTest2x2ATimes2x2BExpResult(0, 1) = 2;
268  MultTest2x2ATimes2x2BExpResult(1, 0) = 3; MultTest2x2ATimes2x2BExpResult(1, 1) = 5;
269  MultTest2x2BTimes2x2AExpResult(0, 0) = 1; MultTest2x2BTimes2x2AExpResult(0, 1) = 3;
270  MultTest2x2BTimes2x2AExpResult(1, 0) = 2; MultTest2x2BTimes2x2AExpResult(1, 1) = 5;
271  MultTest3x3ATimes3x3BExpResult(0, 0) = 3; MultTest3x3ATimes3x3BExpResult(0, 1) = 8; MultTest3x3ATimes3x3BExpResult(0, 2) = 13;
272  MultTest3x3ATimes3x3BExpResult(1, 0) = 6; MultTest3x3ATimes3x3BExpResult(1, 1) = 14; MultTest3x3ATimes3x3BExpResult(1, 2) = 22;
273  MultTest3x3ATimes3x3BExpResult(2, 0) = 9; MultTest3x3ATimes3x3BExpResult(2, 1) = 20; MultTest3x3ATimes3x3BExpResult(2, 2) = 31;
274  MultTest3x3BTimes3x3AExpResult(0, 0) = 3; MultTest3x3BTimes3x3AExpResult(0, 1) = 6; MultTest3x3BTimes3x3AExpResult(0, 2) = 9;
275  MultTest3x3BTimes3x3AExpResult(1, 0) = 8; MultTest3x3BTimes3x3AExpResult(1, 1) = 14; MultTest3x3BTimes3x3AExpResult(1, 2) = 20;
276  MultTest3x3BTimes3x3AExpResult(2, 0) = 13; MultTest3x3BTimes3x3AExpResult(2, 1) = 22; MultTest3x3BTimes3x3AExpResult(2, 2) = 31;
277  MultTest2x3A(0, 0) = 1; MultTest2x3A(0, 1) = 2; MultTest2x3A(0, 2) = 3;
278  MultTest2x3A(1, 0) = 4; MultTest2x3A(1, 1) = 5; MultTest2x3A(1, 2) = 6;
279  MultTest3x2A(0, 0) = 1; MultTest3x2A(0, 1) = 2;
280  MultTest3x2A(1, 0) = 3; MultTest3x2A(1, 1) = 4;
281  MultTest3x2A(2, 0) = 5; MultTest3x2A(2, 1) = 6;
282  MultTest2x3B(0, 0) = 0; MultTest2x3B(0, 1) = 2; MultTest2x3B(0, 2) = 4;
283  MultTest2x3B(1, 0) = 6; MultTest2x3B(1, 1) = 8; MultTest2x3B(1, 2) = 10;
284  MultTest3x2B(0, 0) = 0; MultTest3x2B(0, 1) = 2;
285  MultTest3x2B(1, 0) = 4; MultTest3x2B(1, 1) = 6;
286  MultTest3x2B(2, 0) = 8; MultTest3x2B(2, 1) = 10;
287  MultTest2x3ATimes3x2BExpResult(0, 0) = 32; MultTest2x3ATimes3x2BExpResult(0, 1) = 44;
288  MultTest2x3ATimes3x2BExpResult(1, 0) = 68; MultTest2x3ATimes3x2BExpResult(1, 1) = 98;
289  MultTest3x2ATimes2x3BExpResult(0, 0) = 12; MultTest3x2ATimes2x3BExpResult(0, 1) = 18; MultTest3x2ATimes2x3BExpResult(0, 2) = 24;
290  MultTest3x2ATimes2x3BExpResult(1, 0) = 24; MultTest3x2ATimes2x3BExpResult(1, 1) = 38; MultTest3x2ATimes2x3BExpResult(1, 2) = 52;
291  MultTest3x2ATimes2x3BExpResult(2, 0) = 36; MultTest3x2ATimes2x3BExpResult(2, 1) = 58; MultTest3x2ATimes2x3BExpResult(2, 2) = 80;
292  MultTest2x3BTimes3x2AExpResult(0, 0) = 26; MultTest2x3BTimes3x2AExpResult(0, 1) = 32;
293  MultTest2x3BTimes3x2AExpResult(1, 0) = 80; MultTest2x3BTimes3x2AExpResult(1, 1) = 104;
294  MultTest3x2BTimes2x3AExpResult(0, 0) = 8; MultTest3x2BTimes2x3AExpResult(0, 1) = 10; MultTest3x2BTimes2x3AExpResult(0, 2) = 12;
295  MultTest3x2BTimes2x3AExpResult(1, 0) = 28; MultTest3x2BTimes2x3AExpResult(1, 1) = 38; MultTest3x2BTimes2x3AExpResult(1, 2) = 48;
296  MultTest3x2BTimes2x3AExpResult(2, 0) = 48; MultTest3x2BTimes2x3AExpResult(2, 1) = 66; MultTest3x2BTimes2x3AExpResult(2, 2) = 84;
297 
298  MultTest2x2ATimes2x2B.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest2x2A, MultTest2x2B, 1);
299  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 2x2 * 2x2", MultTest2x2ATimes2x2B, MultTest2x2ATimes2x2BExpResult, verbose);
300  MultTest2x2BTimes2x2A.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest2x2B, MultTest2x2A, 1);
301  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 2x2 * 2x2", MultTest2x2BTimes2x2A, MultTest2x2BTimes2x2AExpResult, verbose);
302  MultTest3x3ATimes3x3B.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest3x3A, MultTest3x3B, 1);
303  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 3x3 * 3x3", MultTest3x3ATimes3x3B, MultTest3x3ATimes3x3BExpResult, verbose);
304  MultTest3x3BTimes3x3A.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest3x3B, MultTest3x3A, 1);
305  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 3x3 * 3x3", MultTest3x3BTimes3x3A, MultTest3x3BTimes3x3AExpResult, verbose);
306  MultTest2x3ATimes3x2B.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest2x3A, MultTest3x2B, 1);
307  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 2x3 * 3x2", MultTest2x3ATimes3x2B, MultTest2x3ATimes3x2BExpResult, verbose);
308  MultTest2x3BTimes3x2A.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest2x3B, MultTest3x2A, 1);
309  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 2x3 * 3x2", MultTest2x3BTimes3x2A, MultTest2x3BTimes3x2AExpResult, verbose);
310  MultTest3x2ATimes2x3B.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest3x2A, MultTest2x3B, 1);
311  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 3x2 * 2x3", MultTest3x2ATimes2x3B, MultTest3x2ATimes2x3BExpResult, verbose);
312  MultTest3x2BTimes2x3A.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1, MultTest3x2B, MultTest2x3A, 1);
313  numberFailedTests += PrintTestResults("multiply() -- mult. results -- 3x2 * 2x3", MultTest3x2BTimes2x3A, MultTest3x2BTimes2x3AExpResult, verbose);
314 
315  DMatrix MultTestHugeA, MultTestHugeB, MultTestHugeATimesHugeBExpResult,
316  MultTestHugeATimesHugeB;
317 
318  const int hugeSize = 100;
319  MultTestHugeA.shape(hugeSize, hugeSize);
320  MultTestHugeB.shape(hugeSize, hugeSize);
321  MultTestHugeATimesHugeBExpResult.shape(hugeSize, hugeSize);
322  MultTestHugeATimesHugeB.shape(hugeSize, hugeSize);
323 
324  for(i = 0; i < hugeSize; i++)
325  {
326  for(j = 0; j < hugeSize; j++)
327  {
328  MultTestHugeA(i, j) = j;
329  MultTestHugeB(i, j) = i;
330  MultTestHugeATimesHugeBExpResult(i, j) = 328350;
331  }
332  }
333 
334  MultTestHugeATimesHugeB.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1.0,
335  MultTestHugeA, MultTestHugeB, 1.0);
336  numberFailedTests += PrintTestResults(
337  "multiply() -- mult. results -- huge * huge",
338  MultTestHugeATimesHugeB, MultTestHugeATimesHugeBExpResult, verbose);
339 
340  //
341  // Check scale methods.
342  //
343  DMatrix ScalTest( 8, 8 );
345  // Scale the entries by 8, it should be 8.
346  if (verbose) std::cout << "scale() -- scale matrix by some number ";
347  returnCode = ScalTest.scale( 8.0 );
348  if (ScalTest(2, 3) == 8.0) {
349  if (verbose) std::cout<< "successful." <<std::endl;
350  } else {
351  if (verbose) std::cout<< "unsuccessful." <<std::endl;
352  numberFailedTests++;
353  }
354  // Pointwise scale the entries by zero, they all should be zero.
355  DMatrix ScalTest2( 8, 8 );
356  if (verbose) std::cout << "scale() -- point-wise scale matrix ";
357  ScalTest.scale( ScalTest2 );
358  if (ScalTest.normOne() == 0.0) {
359  if (verbose) std::cout<< "successful." <<std::endl;
360  } else {
361  if (verbose) std::cout<< "unsuccessful." <<std::endl;
362  numberFailedTests++;
363  }
364  //
365  // Check set methods.
366  //
367  DMatrix CCC( 5, 5 );
368  // Randomize the entries in CCC.
369  testName = "random() -- enter random entries into matrix";
370  returnCode = CCC.random();
371  numberFailedTests += ReturnCodeCheck(testName, returnCode, 0, verbose);
372  // Set the entries of CCC to 1.0.
373  testName = "putScalar() -- set every entry of this matrix to 1.0";
374  returnCode = CCC.putScalar(Teuchos::ScalarTraits<STYPE>::one());
375  numberFailedTests += ReturnCodeCheck(testName, returnCode, 0, verbose);
376  // Check assignment operator.
377  DMatrix CCC2( 5, 5 );
378  CCC2.assign( CCC );
379  if (verbose) std::cout << "assign() -- copy the values of an input matrix ";
380  if ( CCC( 3, 4 ) == Teuchos::ScalarTraits<STYPE>::one() ) {
381  if (verbose) std::cout<< "successful" <<std::endl;
382  } else {
383  if (verbose) std::cout<< "unsuccessful" <<std::endl;
384  numberFailedTests++;
385  }
386  // Create a view into a submatrix of CCC
387  DMatrix CCCview( Teuchos::View, CCC, 3, 3 );
388  DMatrix CCCtest1( 2, 3 );
389  CCCtest1 = CCCview;
390  if (verbose) std::cout << "operator= -- small(empty) = large(view) ";
391  if (CCCtest1.numRows()==3 && CCCtest1.values()==CCC.values()) {
392  if (verbose) std::cout<< "successful" <<std::endl;
393  } else {
394  if (verbose) std::cout<< "unsuccessful" <<std::endl;
395  numberFailedTests++;
396  }
397  CCCtest1 = CCC;
398  if (verbose) std::cout << "operator= -- small(view) = large(copy) ";
399  if (CCCtest1.numRows()==5 && CCCtest1.values()!=CCC.values()) {
400  if (verbose) std::cout<< "successful"<<std::endl;
401  } else {
402  if (verbose) std::cout<< "unsuccessful"<<std::endl;
403  numberFailedTests++;
404  }
405  DMatrix CCCtest2( 2, 2 );
406  CCCtest2 = 3.0;
407  CCCtest1 = CCCtest2;
408  if (verbose) std::cout << "operator= -- large(copy) = small(copy) ";
409  if (CCCtest1.numRows()==2 ) {
410  if (verbose) std::cout<< "successful"<<std::endl;
411  } else {
412  if (verbose) std::cout<< "unsuccessful"<<std::endl;
413  numberFailedTests++;
414  }
415  CCCtest1 = CCCview;
416  if (verbose) std::cout << "operator= -- large(copy) = small(view) ";
417  if (CCCtest1.numRows()==3 && CCCtest1.stride()==5) {
418  if(verbose) std::cout<<"successful" <<std::endl;
419  } else {
420  if (verbose) std::cout<<"unsuccessful"<<std::endl;
421  numberFailedTests++;
422  }
423 
424  DMatrix CCCtest3( CCCview );
425  CCCtest1 += CCCtest3;
426  if (verbose) std::cout << "operator+= -- add two matrices of the same size, but different leading dimension ";
427  if (CCCtest1(1,1)==2.0) {
428  if(verbose) std::cout<<"successful" <<std::endl;
429  } else {
430  if (verbose) std::cout<<"unsuccessful"<<std::endl;
431  numberFailedTests++;
432  }
433  if (verbose) std::cout << "operator+= -- add two matrices of different size (nothing should change) ";
434  CCCtest1 += CCC;
435  if (CCCtest1(1,1)==2.0) {
436  if(verbose) std::cout<<"successful" <<std::endl;
437  } else {
438  if (verbose) std::cout<<"unsuccessful"<<std::endl;
439  numberFailedTests++;
440  }
441  //
442  // Check overloaded operators.
443  //
444  bool op_result;
445  MultTestHugeATimesHugeB.reshape(10, 10);
446  op_result = (MultTestHugeATimesHugeB == MultTestHugeATimesHugeBExpResult);
447  if (verbose) {
448  std::cout << "operator== -- results -- small == huge "<< (op_result == false ? "successful" : "failed" )<<std::endl;
449  }
450  op_result = (MultTestHugeATimesHugeB != MultTestHugeATimesHugeBExpResult);
451  if (verbose) {
452  std::cout << "operator!= -- results -- small != huge "<< (op_result == true ? "successful" : "failed" )<<std::endl;
453  std::cout << std::endl<< MultTestHugeATimesHugeB << std::endl;
454  //These won't work unless boundschecking is enabled.
455  //std::cout << MultTestHugeATimesHugeB(100, 1) << std::endl;
456  //std::cout << MultTestHugeATimesHugeB(1, 100) << std::endl;
457  }
458 
459 
460  if (verbose) std::cout<<std::endl<<"********** CHECKING TEUCHOS SERIAL DENSE VECTOR **********"<<std::endl<<std::endl;
461 
462  DVector DefConTestV;
463  if (verbose) std::cout <<"default constructor -- construct empty std::vector ";
464  if ( DefConTestV.values()!=NULL || DefConTestV.length()!=0 || DefConTestV.numRows()!=0 ||DefConTestV.stride()!=0 ) {
465  if (verbose) std::cout << "unsuccessful."<<std::endl;
466  numberFailedTests++;
467  } else {
468  if (verbose) std::cout << "successful."<<std::endl;
469  }
470 
471  // constructor 1 (matrix w/ dimension but empty)
472 
473  DVector Con1TestV( 3 );
474  if (verbose) std::cout <<"constructor 1 -- empty std::vector with given dimensions ";
475  if ( Con1TestV.length()!=3 || Con1TestV.numCols()!=1 || Con1TestV( 1 )!=0.0 ) {
476  if (verbose) std::cout << "unsuccessful."<<std::endl;
477  numberFailedTests++;
478  } else {
479  if (verbose) std::cout << "successful."<<std::endl;
480  }
481 
482  // constructor 2 (from array) tests
483 
484  DVector Con2Test1V(Teuchos::Copy, a, 4);
485  if (verbose) std::cout <<"constructor 2 -- construct std::vector from array subrange ";
486  if ( Con2Test1V.numRows()!=4 || Con2Test1V.numCols()!=1 || Con2Test1V[ 2 ]!=2.0 ) {
487  if (verbose) std::cout << "unsuccessful."<<std::endl;
488  numberFailedTests++;
489  } else {
490  if (verbose) std::cout << "successful."<<std::endl;
491  }
492 
493  // constructor 3 (copy constructor)
494 
495  DVector Con3TestCopyV( Con2Test1V );
496  if(verbose) std::cout <<"constructor 3 -- copy constructor ";
497  if ( Con3TestCopyV != Con2Test1V ) {
498  if (verbose) std::cout << "unsuccessful."<<std::endl;
499  numberFailedTests++;
500  } else {
501  if (verbose) std::cout << "successful."<<std::endl;
502  }
503 
504  // non-member helper function (construct vector view of matrix column)
505 
507  DVector ColViewTestV = Teuchos::getCol<OTYPE,STYPE>( Teuchos::View, AAA, col );
508  if (verbose) std::cout <<"non-method helper function -- construct vector view of second column of matrix ";
509  if ( ColViewTestV.normInf() != 1.0 || ColViewTestV.normOne() != 3.0 ) {
510  if (verbose) std::cout << "unsuccessful."<<std::endl;
511  numberFailedTests++;
512  } else {
513  if (verbose) std::cout << "successful."<<std::endl;
514  }
515 
516  // checking norms
517 
518  numberFailedTests += PrintTestResults("normOne of a 3x1 std::vector", Con2Test1V.normOne(), 6.0, verbose);
519  numberFailedTests += PrintTestResults("normInf of a 3x1 std::vector", Con2Test1V.normInf(), 3.0, verbose);
520  Con2Test1V = Teuchos::ScalarTraits<STYPE>::one();
521  numberFailedTests += PrintTestResults("normFrobenius of a 3x1 std::vector", Con2Test1V.normFrobenius(), 2.0, verbose);
522 
523  // check size/resize
524 
525  DVector SizeTestV1;
526  SizeTestV1.size( 5 );
527  if(verbose) std::cout <<"size() -- test ";
528  if (SizeTestV1( 4 )!= 0.0) {
529  if (verbose) std::cout << "unsuccessful."<<std::endl;
530  numberFailedTests++;
531  } else {
532  if (verbose) std::cout << "successful."<<std::endl;
533  }
534  SizeTestV1 = 2.0*Teuchos::ScalarTraits<STYPE>::one();
535  SizeTestV1.resize( 10 );
536  if(verbose) std::cout <<"resize() -- test small --> large ";
537  if (SizeTestV1[ 4 ]!= 2.0 || SizeTestV1[ 8 ]!=0.0 ) {
538  if (verbose) std::cout << "unsuccessful."<<std::endl;
539  numberFailedTests++;
540  } else {
541  if (verbose) std::cout << "successful."<<std::endl;
542  }
543  SizeTestV1.resize( 3 );
544  if(verbose) std::cout <<"resize() -- test large --> small ";
545  if (SizeTestV1( 2 )!= 2.0) {
546  if (verbose) std::cout << "unsuccessful."<<std::endl;
547  numberFailedTests++;
548  } else {
549  if (verbose) std::cout << "successful."<<std::endl;
550  }
551 
552  DVector OpEqTestV1( 10 ); OpEqTestV1 = 3.0*Teuchos::ScalarTraits<STYPE>::one();
553  DVector OpEqTestV2( Teuchos::View, OpEqTestV1.values(), 3 );
554  DVector OpEqTestV3( 2 );
555  OpEqTestV3 = OpEqTestV2;
556  if (verbose) std::cout << "operator= -- small(empty) = large(view) ";
557  if (OpEqTestV3.length()==3 && OpEqTestV3.values()==OpEqTestV2.values()) {
558  if (verbose) std::cout<< "successful"<<std::endl;
559  } else {
560  if (verbose) std::cout<< "unsuccessful"<<std::endl;
561  numberFailedTests++;
562  }
563  OpEqTestV3 = OpEqTestV1;
564  if (verbose) std::cout << "operator= -- small(view) = large(copy) ";
565  if (OpEqTestV3.length()==10 && OpEqTestV3.values()!=OpEqTestV1.values()) {
566  if (verbose) std::cout<< "successful"<<std::endl;
567  } else {
568  if (verbose) std::cout<< "unsuccessful"<<std::endl;
569  numberFailedTests++;
570  }
571  OpEqTestV3.size(5);
572  OpEqTestV3 = OpEqTestV1;
573  if (verbose) std::cout << "operator= -- small(copy) = large(copy) ";
574  if (OpEqTestV3.length()==10 && OpEqTestV3.values()!=OpEqTestV1.values() && OpEqTestV3[ 9 ]==3.0) {
575  if (verbose) std::cout<< "successful"<<std::endl;
576  } else {
577  if (verbose) std::cout<< "unsuccessful"<<std::endl;
578  numberFailedTests++;
579  }
580 
581  DVector OpSumTestV1( OpEqTestV2 );
582  OpSumTestV1 += OpEqTestV2;
583  if (verbose) std::cout << "operator+= -- add two vectors of the same size, but different leading dimension ";
584  if (OpSumTestV1( 1 )==6.0) {
585  if (verbose) std::cout<<"successful" <<std::endl;
586  } else {
587  if (verbose) std::cout<<"unsuccessful"<<std::endl;
588  numberFailedTests++;
589  }
590  if (verbose) std::cout << "operator+= -- add two vectors of different size (nothing should change) ";
591  OpSumTestV1 += OpEqTestV1;
592  if (OpSumTestV1( 1 )==6.0) {
593  if (verbose) std::cout<<"successful" <<std::endl;
594  } else {
595  if (verbose) std::cout<<"unsuccessful"<<std::endl;
596  numberFailedTests++;
597  }
598 
599  DVector OpCompTestV1( 5 );
600  OpCompTestV1 = 2.0*Teuchos::ScalarTraits<STYPE>::one();
601  if(verbose) std::cout <<"operator== -- test large == small ";
602  if (OpCompTestV1 == SizeTestV1) {
603  if (verbose) std::cout << "unsuccessful."<<std::endl;
604  numberFailedTests++;
605  } else {
606  if (verbose) std::cout << "successful."<<std::endl;
607  }
608  if(verbose) std::cout <<"operator!= -- test large != small ";
609  if (OpCompTestV1 != SizeTestV1) {
610  if (verbose) std::cout << "successful."<<std::endl;
611  } else {
612  if (verbose) std::cout << "successful."<<std::endl;
613  numberFailedTests++;
614  }
615 
616  DVector ColSetTestV( AAA.numRows() );
617  ColSetTestV.putScalar( 2.0 );
618  bool ret = Teuchos::setCol<OTYPE,STYPE>( ColSetTestV, col, AAA );
619  if (verbose) std::cout <<"non-method helper function -- set second column of matrix with vector ";
620  if ( ColViewTestV.normInf() != 2.0 || ColViewTestV.normOne() != 6.0 || ret == false ) {
621  if (verbose) std::cout << "unsuccessful."<<std::endl;
622  numberFailedTests++;
623  } else {
624  if (verbose) std::cout << "successful."<<std::endl;
625  }
626  //
627  // If a test failed output the number of failed tests.
628  //
629  if(numberFailedTests > 0)
630  {
631  if (verbose) {
632  std::cout << "Number of failed tests: " << numberFailedTests << std::endl;
633  std::cout << "End Result: TEST FAILED" << std::endl;
634  return -1;
635  }
636  }
637  if(numberFailedTests == 0)
638  std::cout << "End Result: TEST PASSED" << std::endl;
639 
640  return 0;
641 }
642 
643 template<typename TYPE>
644 int PrintTestResults(std::string testName, TYPE calculatedResult, TYPE expectedResult, bool verbose)
645 {
646  int result;
647  if(calculatedResult == expectedResult)
648  {
649  if(verbose) std::cout << testName << " successful." << std::endl;
650  result = 0;
651  }
652  else
653  {
654  if(verbose) std::cout << testName << " unsuccessful." << std::endl;
655  result = 1;
656  }
657  return result;
658 }
659 
660 int ReturnCodeCheck(std::string testName, int returnCode, int expectedResult, bool verbose)
661 {
662  int result;
663  if(expectedResult == 0)
664  {
665  if(returnCode == 0)
666  {
667  if(verbose) std::cout << testName << " test successful." << std::endl;
668  result = 0;
669  }
670  else
671  {
672  if(verbose) std::cout << testName << " test unsuccessful. Return code was " << returnCode << "." << std::endl;
673  result = 1;
674  }
675  }
676  else
677  {
678  if(returnCode != 0)
679  {
680  if(verbose) std::cout << testName << " test successful -- failed as expected." << std::endl;
681  result = 0;
682  }
683  else
684  {
685  if(verbose) std::cout << testName << " test unsuccessful -- did not fail as expected. Return code was " << returnCode << "." << std::endl;
686  result = 1;
687  }
688  }
689  return result;
690 }
ScalarTraits< ScalarType >::magnitudeType normFrobenius() const
Returns the Frobenius-norm of the matrix.
static T one()
Returns representation of one for this ordinal type.
int PrintTestResults(std::string, TYPE, TYPE, bool)
Non-member helper functions on the templated serial, dense matrix/vector classes. ...
Templated serial dense matrix class.
ScalarType * values() const
Data array access method.
Teuchos::SerialDenseVector< int, std::complex< Real > > DVector
int multiply(ETransp transa, ETransp transb, ScalarType alpha, const SerialDenseMatrix< OrdinalType, ScalarType > &A, const SerialDenseMatrix< OrdinalType, ScalarType > &B, ScalarType beta)
Multiply A * B and add them to this; this = beta * this + alpha*A*B.
int scale(const ScalarType alpha)
Scale this matrix by alpha; *this = alpha**this.
This class creates and provides basic support for dense vectors of templated type as a specialization...
This structure defines some basic traits for a scalar field type.
int main(int argc, char *argv[])
ScalarTraits< ScalarType >::magnitudeType normInf() const
Returns the Infinity-norm of the matrix.
OrdinalType numRows() const
Returns the row dimension of this matrix.
int putScalar(const ScalarType value=Teuchos::ScalarTraits< ScalarType >::zero())
Set all values in the matrix to a constant value.
Teuchos::SerialDenseMatrix< int, std::complex< Real > > DMatrix
std::string Teuchos_Version()
int ReturnCodeCheck(std::string, int, int, bool)
bool empty() const
Returns whether this matrix is empty.
ScalarTraits< ScalarType >::magnitudeType normOne() const
Returns the 1-norm of the matrix.
int size(OrdinalType length_in)
Size method for changing the size of a SerialDenseVector, initializing entries to zero...
int random()
Set all values in the matrix to be random numbers.
int reshape(OrdinalType numRows, OrdinalType numCols)
Reshaping method for changing the size of a SerialDenseMatrix, keeping the entries.
Templated serial dense vector class.
OrdinalType stride() const
Returns the stride between the columns of this matrix in memory.
int shape(OrdinalType numRows, OrdinalType numCols)
Shape method for changing the size of a SerialDenseMatrix, initializing entries to zero...
OrdinalType numCols() const
Returns the column dimension of this matrix.
OrdinalType length() const
Returns the length of this vector.
SerialDenseMatrix< OrdinalType, ScalarType > & assign(const SerialDenseMatrix< OrdinalType, ScalarType > &Source)
Copies values from one matrix to another.
static T one()
Returns representation of one for this scalar type.
This class creates and provides basic support for dense rectangular matrix of templated type...