Epetra Package Browser (Single Doxygen Collection) Development
test/IntSerialDense/cxx_main.cpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //               Epetra: Linear Algebra Services Package 
00005 //                 Copyright 2001 Sandia Corporation
00006 // 
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00038 // 
00039 // ************************************************************************
00040 //@HEADER
00041 
00042 
00043 // Epetra_IntSerialDense Test routine
00044 
00045 #include "Epetra_IntSerialDenseMatrix.h"
00046 #include "Epetra_IntSerialDenseVector.h"
00047 #include "../epetra_test_err.h"
00048 #include "Epetra_ConfigDefs.h"
00049 #include "Epetra_DataAccess.h"
00050 #include "Epetra_Version.h"
00051 #ifdef EPETRA_MPI
00052 #include "Epetra_MpiComm.h"
00053 #include <mpi.h>
00054 #else
00055 #include "Epetra_SerialComm.h"
00056 #endif
00057 
00058 // matrix-testing prototypes
00059 int matrixCoverage(bool verbose, bool debug);
00060 int matrixCtr(bool verbose, bool debug);
00061 int matrixCpyCtr(bool verbose, bool debug);
00062 int matrixAssignment(bool verbose, bool debug);
00063 int matrixExceptions(bool verbose, bool debug);
00064 // vector-testing prototypes
00065 int vectorCoverage(bool verbose, bool debug);
00066 int vectorCtr(bool verbose, bool debug);
00067 int vectorCpyCtr(bool verbose, bool debug);
00068 int vectorAssignment(bool verbose, bool debug);
00069 int vectorExceptions(bool verbose, bool debug);
00070 // helper function prototypes
00071 bool identicalSignatures(Epetra_IntSerialDenseMatrix& a, Epetra_IntSerialDenseMatrix& b, bool testLDA = true);
00072 bool seperateData(Epetra_IntSerialDenseMatrix& a, Epetra_IntSerialDenseMatrix& b);
00073 int* getRandArray(int length);
00074 int randomInt();
00075 void printArray(int* array, int length);
00076 void printMat(const char* name, Epetra_IntSerialDenseMatrix& matrix);
00077 void printHeading(const char* heading);
00078 
00079 int main(int argc, char *argv[]) {
00080 //============================
00081 // Check for verbose or debug
00082 
00083 bool verbose = false;
00084 bool debug = false;
00085 if(argc > 1) {
00086     if((argv[1][0] == '-') && (argv[1][1] == 'v'))
00087         verbose = true;
00088     if((argv[1][0] == '-') && (argv[1][1] == 'd')) {
00089         debug = true;
00090         verbose = true;
00091     }
00092 }
00093 
00094 //============================
00095 // Initialize Comm
00096 
00097 #ifdef EPETRA_MPI
00098   // Initialize MPI
00099   MPI_Init(&argc,&argv);
00100   Epetra_MpiComm Comm(MPI_COMM_WORLD);
00101 #else
00102   Epetra_SerialComm Comm;
00103 #endif
00104 
00105   if (verbose && Comm.MyPID()==0)
00106     cout << Epetra_Version() << endl << endl;
00107   
00108 //============================
00109 // other initial setup
00110 
00111   if (!verbose) Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00112 
00113   int ierr = 0;
00114   int returnierr = 0;
00115 
00116 //============================
00117 // Test vector
00118 
00119   ierr = vectorExceptions(verbose, debug);
00120   EPETRA_TEST_ERR(ierr, returnierr);
00121 
00122   ierr = vectorCtr(verbose, debug);
00123   EPETRA_TEST_ERR(ierr, returnierr);
00124   
00125   ierr = vectorCpyCtr(verbose, debug);
00126   EPETRA_TEST_ERR(ierr, returnierr);
00127   
00128   ierr = vectorAssignment(verbose, debug);
00129   EPETRA_TEST_ERR(ierr, returnierr);
00130 
00131   ierr = vectorCoverage(verbose, debug);
00132   EPETRA_TEST_ERR(ierr, returnierr);
00133 
00134 //============================
00135 // Test matrix
00136 
00137   ierr = matrixExceptions(verbose, debug);
00138   EPETRA_TEST_ERR(ierr, returnierr);
00139 
00140   ierr = matrixCtr(verbose, debug);
00141   EPETRA_TEST_ERR(ierr, returnierr);
00142 
00143   ierr = matrixCpyCtr(verbose, debug);
00144   EPETRA_TEST_ERR(ierr, returnierr);
00145 
00146   ierr = matrixAssignment(verbose, debug);
00147   EPETRA_TEST_ERR(ierr, returnierr);
00148 
00149   ierr = matrixCoverage(verbose, debug);
00150   EPETRA_TEST_ERR(ierr, returnierr);
00151 
00152 //============================
00153 // end of program cleanup
00154 
00155 #ifdef EPETRA_MPI
00156   MPI_Finalize();
00157 #endif
00158   
00159   return(returnierr);
00160 }
00161 
00162 //=========================================================================
00163 //=========================================================================
00164 // Matrix-testing functions
00165 //=========================================================================
00166 //=========================================================================
00167 // call functions we don't have tests for, for sake of code coverage.
00168 // (the functions are called, but their output isn't checked for correctness).
00169 int matrixCoverage(bool verbose, bool debug) {
00170   if(verbose) printHeading("Testing other matrix functions");
00171 
00172   int* rand1 = getRandArray(9);
00173   if(debug) printArray(rand1, 9);
00174   Epetra_IntSerialDenseMatrix m1(Copy, rand1, 3, 3, 3);
00175   if(debug) printMat("m1",m1);
00176   delete[] rand1;
00177 
00178   if(verbose) cout << "calling one norm" << endl;
00179   int onenorm = m1.OneNorm();
00180   if(debug) cout << "m1.OneNorm() = " << onenorm << endl;
00181 
00182   if(verbose) cout << "calling infinity norm" << endl;
00183   int infnorm = m1.InfNorm();
00184   if(debug) cout << "m1.InfNorm() = " << infnorm << endl;
00185 
00186   if(verbose) cout << "calling random" << endl;
00187   Epetra_IntSerialDenseMatrix m2(6,6);
00188   if(debug) printMat("m2 (before)",m2);
00189   m2.Random();
00190   if(debug) printMat("m2 (after)",m2);
00191 
00192   if(verbose) cout << "Checked OK." << endl;
00193 
00194   return(0);
00195 }
00196 
00197 //=========================================================================
00198 // test matrix default constructor, user constructor (copy & view),
00199 // () [] operators (read & write), shape (larger), reshape (smaller)
00200 int matrixCtr(bool verbose, bool debug) {
00201   const int m1rows = 5;
00202   const int m1cols = 4;
00203   const int m1arows = 4;
00204   const int m1acols = 6;
00205   const int m2rows = 2;
00206   const int m2cols = 7;
00207   const int m3rows = 8;
00208   const int m3cols = 3;
00209   const int m3Rrows = 5; // should be smaller than m3rows
00210   const int m3Rcols = 2; // should be smaller than m3cols
00211 
00212   int ierr = 0;
00213   int returnierr = 0;
00214   if(verbose) printHeading("Testing matrix constructors");
00215 
00216   if(verbose) cout << "default constructor" << endl;
00217   Epetra_IntSerialDenseMatrix m1;
00218   EPETRA_TEST_ERR(!(m1.CV() == Copy), ierr);
00219   if(verbose) cout << "shaping" << endl;
00220   m1.Shape(m1rows, m1cols);
00221   EPETRA_TEST_ERR(!(m1.M() == m1rows), ierr);
00222   for(int i = 0; i < m1rows; i++)
00223     for(int j = 0; j < m1cols; j++)
00224       EPETRA_TEST_ERR(!(m1(i,j) == 0), ierr);
00225   if(debug) printMat("m1",m1);
00226   returnierr += ierr;
00227   if(ierr == 0) 
00228      if(verbose) cout << "Checked OK." << endl;
00229   ierr = 0;
00230   if(verbose) cout << "\nmanually setting values" << endl;
00231   int* m1rand = getRandArray(m1rows * m1cols);
00232   for(int i = 0; i < m1rows; i++)
00233     for(int j = 0; j < m1cols; j++)
00234       m1(i,j) = m1rand[i*m1cols + j];
00235   for(int i = 0; i < m1rows; i++)
00236     for(int j = 0; j < m1cols; j++)
00237     EPETRA_TEST_ERR(!(m1[j][i] == m1rand[i*m1cols + j]), ierr);
00238   if(debug) {
00239     printArray(m1rand, m1rows * m1cols);
00240     printMat("m1",m1);
00241   }
00242   delete[] m1rand;
00243   returnierr += ierr;
00244   if(ierr == 0) 
00245      if(verbose) cout << "Checked OK." << endl;
00246   ierr = 0;
00247   
00248   if(verbose) cout << "\nshaped constructor" << endl;
00249   Epetra_IntSerialDenseMatrix m1a(m1arows, m1acols);
00250   EPETRA_TEST_ERR(!(m1a.M() == m1arows), ierr);
00251   EPETRA_TEST_ERR(!(m1a.N() == m1acols), ierr);
00252   EPETRA_TEST_ERR(!(m1a.LDA() == m1arows), ierr);
00253   EPETRA_TEST_ERR(!(m1a.CV() == Copy),ierr);
00254   if(debug) printMat("m1a", m1a);
00255   returnierr += ierr;
00256   if(ierr == 0) 
00257      if(verbose) cout << "Checked OK." << endl;
00258   ierr = 0;
00259   
00260   if(verbose) cout << "\nuser data constructor (view)" << endl;
00261   int* m2rand = getRandArray(m2rows * m2cols);
00262   if(debug) printArray(m2rand, m2rows * m2cols);
00263   Epetra_IntSerialDenseMatrix m2(View, m2rand, m2rows, m2rows, m2cols);
00264   EPETRA_TEST_ERR(!(m2.CV() == View), ierr);
00265   EPETRA_TEST_ERR(!(m2.M() == m2rows), ierr);
00266   EPETRA_TEST_ERR(!(m2.N() == m2cols), ierr);
00267   EPETRA_TEST_ERR(!(m2.LDA() == m2rows), ierr);
00268   EPETRA_TEST_ERR(!(m2.A() == m2rand), ierr);
00269   if(debug) printMat("m2",m2);
00270   returnierr += ierr;
00271   if(ierr == 0) 
00272      if(verbose) cout << "Checked OK." << endl;
00273   ierr = 0;
00274 
00275   if(verbose) cout << "\nchanging value, checking for correct behavior" << endl;
00276   int* m2randcopy = new int[m2rows * m2cols]; // backup of original values
00277   for(int i = 0; i < m2rows * m2cols; i++)
00278     m2randcopy[i] = m2rand[i];
00279   m2(0,4) = m2(0,4) + 1; // magic numbers for which element to change
00280   m2randcopy[4 * m2rows] = m2randcopy[4 * m2rows] + 1;
00281   for(int i = 0; i < m2rows * m2cols; i++)
00282     EPETRA_TEST_ERR(!(m2rand[i] == m2randcopy[i]), ierr); // m2rand should have updated correctly
00283   if(debug) {
00284     printArray(m2rand, m2rows * m2cols);
00285     printMat("m2",m2);
00286   }
00287   delete[] m2rand;
00288   delete[] m2randcopy;
00289   returnierr += ierr;
00290   if(ierr == 0) 
00291      if(verbose) cout << "Checked OK." << endl;
00292   ierr = 0;
00293 
00294   if(verbose) cout << "\nuser data constructor (copy)" << endl;
00295   int* m3rand = getRandArray(m3rows * m3cols);
00296   if(debug) printArray(m3rand, m3rows * m3cols);
00297   int* m3randcopy = new int[m3rows * m3cols];
00298   for(int i = 0; i < m3rows * m3cols; i++)
00299     m3randcopy[i] = m3rand[i];
00300   Epetra_IntSerialDenseMatrix m3(Copy, m3rand, m3rows, m3rows, m3cols);
00301   if(debug) printMat("m3",m3);
00302   if(verbose) cout << "checking to see if data initialized correctly" << endl;
00303   EPETRA_TEST_ERR(!(m3.CV() == Copy), ierr);
00304   EPETRA_TEST_ERR(!(m3.M() == m3rows), ierr);
00305   EPETRA_TEST_ERR(!(m3.N() == m3cols), ierr);
00306   EPETRA_TEST_ERR(!(m3.LDA() == m3rows), ierr);
00307   EPETRA_TEST_ERR(!(m3.A() != m3rand), ierr); // should not be a view
00308   for(int i = 0; i < m3rows; i++)
00309     for(int j = 0; j < m3cols; j++)
00310       EPETRA_TEST_ERR(!(m3[j][i] == m3rand[j * m3rows + i]), ierr); // data should be identical to user array
00311   returnierr += ierr;
00312   if(ierr == 0) 
00313      if(verbose) cout << "Checked OK." << endl;
00314   ierr = 0;
00315 
00316   if(verbose) cout << "\nmodifying entry" << endl;
00317   m3[1][5] = m3[1][5] + 3; // magic numbers for which element to change
00318   for(int i = 0; i < m3rows * m3cols; i++)
00319     EPETRA_TEST_ERR(!(m3rand[i] == m3randcopy[i]), ierr); // user array should be unchanged
00320   m3rand[13] = m3rand[13] + 3; // same magic modification performed to user's array
00321   for(int i = 0; i < m3rows; i++)
00322     for(int j = 0; j < m3cols; j++)
00323       EPETRA_TEST_ERR(!(m3[j][i] == m3rand[j * m3rows + i]), ierr); // should equal user array with same modification performed
00324   if(debug) {
00325     printArray(m3rand, m3rows * m3cols);
00326     printMat("m3",m3);
00327   }
00328   returnierr += ierr;
00329   if(ierr == 0) 
00330      if(verbose) cout << "Checked OK." << endl;
00331   ierr = 0;
00332   
00333   if(verbose) cout << "\nreshaping" << endl;
00334   m3.Reshape(m3Rrows, m3Rcols);
00335   EPETRA_TEST_ERR(!(m3.M() == m3Rrows), ierr);
00336   EPETRA_TEST_ERR(!(m3.N() == m3Rcols), ierr);
00337   EPETRA_TEST_ERR(!(m3.LDA() == m3Rrows), ierr);
00338   for(int i = 0; i < m3Rrows; i++)
00339     for(int j = 0; j < m3Rcols; j++)
00340       EPETRA_TEST_ERR(!(m3[j][i] == m3rand[j * m3rows + i]),ierr);
00341   if(debug) printMat("m3",m3);
00342 
00343   delete[] m3rand;
00344   delete[] m3randcopy;
00345   returnierr += ierr;
00346   if(ierr == 0) 
00347      if(verbose) cout << "Checked OK." << endl;
00348   ierr = 0;
00349   
00350   if(verbose) cout << "\nChecking pointer on zero-sized matrix" << endl;
00351   int* before = m3.A();
00352   if(verbose) cout << "Reshaping to (4,0)" << endl;
00353   if(debug) cout << "Before = " << before << endl;
00354   EPETRA_TEST_ERR(!(before != 0), ierr);
00355   m3.Reshape(4,0);
00356   int* after = m3.A();
00357   EPETRA_TEST_ERR(!(after == 0), ierr);
00358   if(debug) cout << "After = " << after << endl;
00359   m3.Shape(3,3);
00360   before = m3.A();
00361   if(verbose) cout << "Shaping to (0,3)" << endl;
00362   if(debug) cout << "Before = " << before << endl;
00363   EPETRA_TEST_ERR(!(before != 0), ierr);
00364   m3.Shape(0,3);
00365   after = m3.A();
00366   EPETRA_TEST_ERR(!(after == 0), ierr);
00367   if(debug) cout << "After = " << after << endl;
00368   returnierr += ierr;
00369   if(ierr == 0) 
00370      if(verbose) cout << "Checked OK." << endl;
00371   ierr = 0;
00372 
00373   return(returnierr);
00374 }
00375 
00376 //=========================================================================
00377 // test matrix copy constructor (copy & view)
00378 int matrixCpyCtr(bool verbose, bool debug) {
00379   const int m1rows = 5;
00380   const int m1cols = 4;
00381   const int m2rows = 2;
00382   const int m2cols = 6;
00383 
00384   int ierr = 0;
00385   int returnierr = 0;
00386   if(verbose) printHeading("Testing matrix copy constructors");
00387 
00388   if(verbose) cout << "checking copy constructor (view)" << endl;
00389   int* m1rand = getRandArray(m1rows * m1cols);
00390   if(debug) printArray(m1rand, m1rows * m1cols);
00391   Epetra_IntSerialDenseMatrix m1(View, m1rand, m1rows, m1rows, m1cols);
00392   if(debug) {
00393     cout << "original matrix:" << endl;
00394     printMat("m1",m1);
00395   }
00396   Epetra_IntSerialDenseMatrix m1clone(m1);
00397   if(debug) {
00398     cout << "clone matrix:" << endl;
00399     printMat("m1clone",m1clone);
00400   }
00401   if(verbose) cout << "making sure signatures match" << endl;
00402   EPETRA_TEST_ERR(!identicalSignatures(m1, m1clone), ierr);
00403   delete[] m1rand;
00404   returnierr += ierr;
00405   if(ierr == 0)
00406     if(verbose) cout << "Checked OK." << endl;
00407   ierr = 0;
00408   
00409   if(verbose) cout << "\nchecking copy constructor (copy)" << endl;
00410   int* m2rand = getRandArray(m2rows * m2cols);
00411   if(debug) printArray(m2rand, m2rows * m2cols);
00412   Epetra_IntSerialDenseMatrix m2(Copy, m2rand, m2rows, m2rows, m2cols);
00413   if(debug) {
00414     cout << "original matrix:" << endl;
00415     printMat("m2",m2);
00416   }
00417   Epetra_IntSerialDenseMatrix m2clone(m2);
00418   if(debug) {
00419     cout << "clone matrix:" << endl;
00420     printMat("m2clone",m2clone);
00421   }
00422   if(verbose) cout << "checking that signatures match" << endl;
00423   EPETRA_TEST_ERR(!identicalSignatures(m2, m2clone), ierr);
00424   returnierr += ierr;
00425   if(ierr == 0)
00426     if(verbose) cout << "Checked OK." << endl;
00427   ierr = 0;
00428 
00429   if(verbose) cout << "\nmodifying entry in m2, m2clone should be unchanged" << endl;
00430   EPETRA_TEST_ERR(!seperateData(m2, m2clone), ierr);
00431   if(debug) {
00432     printArray(m2rand, m2rows * m2cols);
00433     cout << "orig:" << endl;
00434     printMat("m2",m2);
00435     cout << "clone:" << endl;
00436     printMat("m2clone",m2clone);
00437   }
00438   delete[] m2rand;
00439   returnierr += ierr;
00440   if(ierr == 0)
00441     if(verbose) cout << "Checked OK." << endl;
00442   ierr = 0;
00443   
00444   return(returnierr);
00445 }
00446 
00447 //=========================================================================
00448 // test matrix error-reporting
00449 int matrixExceptions(bool verbose, bool debug) {
00450   int returnierr = 0;
00451   int ierr = 0;
00452   bool caught = false;
00453   Epetra_IntSerialDenseMatrix* matrix;
00454 
00455   if(verbose) printHeading("Testing matrix error-reporting.\nExpect error messages if EPETRA_NO_ERROR_REPORTS is not defined.");
00456   
00457   // invalid dimension to sized ctr (2 cases)
00458   try {
00459     caught = false;
00460     if(verbose) cout << "Checking Epetra_IntSerialDenseMatrix(-1, 6) - invalid rows";
00461     matrix = new Epetra_IntSerialDenseMatrix(-1, 6);
00462   }
00463   catch(int error) {
00464     caught = true;
00465     EPETRA_TEST_ERR(error != -1, returnierr);
00466     if(error == -1)
00467       if(verbose) cout << "Checked OK." << endl;
00468   }
00469   EPETRA_TEST_ERR(!caught, returnierr);
00470   try {
00471     caught = false;
00472     if(verbose) cout << "\nChecking Epetra_IntSerialDenseMatrix(3, -5) - invalid cols";
00473     matrix = new Epetra_IntSerialDenseMatrix(3, -5);
00474   }
00475   catch(int error) {
00476     caught = true;
00477     EPETRA_TEST_ERR(error != -1, returnierr);
00478     if(error == -1)
00479       if(verbose) cout << "Checked OK." << endl;
00480   }
00481   EPETRA_TEST_ERR(!caught, returnierr);
00482 
00483   // invalid dimension to user-data ctr (3 cases)
00484      int* rand2 = getRandArray(2);
00485   try {
00486     caught = false;
00487     if(verbose) cout << "\nChecking Epetra_IntSerialDenseMatrix(Copy, int*, -1, 2, 2) - invalid lda";
00488     matrix = new Epetra_IntSerialDenseMatrix(Copy, rand2, -1, 2, 2);
00489   }
00490   catch(int error) {
00491     caught = true;
00492     EPETRA_TEST_ERR(error != -1, returnierr);
00493     if(error == -1)
00494       if(verbose) cout << "Checked OK." << endl;
00495   }
00496   EPETRA_TEST_ERR(!caught, returnierr);
00497   try {
00498     caught = false;
00499     if(verbose) cout << "\nChecking Epetra_IntSerialDenseMatrix(Copy, int*, 3, -2, 3) - invalid rows";
00500     matrix = new Epetra_IntSerialDenseMatrix(Copy, rand2, 3, -2, 3);
00501   }
00502   catch(int error) {
00503     caught = true;
00504     EPETRA_TEST_ERR(error != -1, returnierr);
00505     if(error == -1)
00506       if(verbose) cout << "Checked OK." << endl;
00507   }
00508   EPETRA_TEST_ERR(!caught, returnierr);
00509   try {
00510     caught = false;
00511     if(verbose) cout << "\nChecking Epetra_IntSerialDenseMatrix(Copy, int*, 4, 4, -4) - invalid cols";
00512     matrix = new Epetra_IntSerialDenseMatrix(Copy, rand2, -4, 4, -4);
00513   }
00514   catch(int error) {
00515     caught = true;
00516     EPETRA_TEST_ERR(error != -1, returnierr);
00517     if(error == -1)
00518       if(verbose) cout << "Checked OK." << endl;
00519   }
00520   EPETRA_TEST_ERR(!caught, returnierr);
00521      delete [] rand2;
00522   
00523   // null pointer to user-data ctr
00524   try {
00525     caught = false;
00526     if(verbose) cout << "\nChecking Epetra_IntSerialDenseMatrix(Copy, 0, 5, 5, 5) - null pointer";
00527     matrix = new Epetra_IntSerialDenseMatrix(Copy, 0, 5, 5, 5);
00528   }
00529   catch(int error) {
00530     caught = true;
00531     EPETRA_TEST_ERR(error != -3, returnierr);
00532     if(error == -3)
00533       if(verbose) cout << "Checked OK." << endl;
00534   }
00535   EPETRA_TEST_ERR(!caught, returnierr);
00536 
00537   // invalid parameter to shape (2 cases)
00538   Epetra_IntSerialDenseMatrix m1;
00539   if(verbose) cout << "\nChecking Shape(-2, 2) - invalid rows" << endl;
00540   ierr = m1.Shape(-2, 2);
00541   EPETRA_TEST_ERR(!(ierr == -1), returnierr);
00542   if(ierr == -1)
00543     if(verbose) cout << "Checked OK." << endl;
00544   if(verbose) cout << "\nChecking Shape(3, -2) - invalid cols" << endl;
00545   ierr = m1.Shape(3, -2);
00546   EPETRA_TEST_ERR(!(ierr == -1), returnierr);
00547   if(ierr == -1)
00548     if(verbose) cout << "Checked OK." << endl;
00549   
00550   // invalid parameter to reshape (2 cases)
00551   m1.Shape(5, 5);
00552   if(verbose) cout << "\nChecking Reshape(-4, 3) - invalid rows" << endl;
00553   ierr = m1.Reshape(-4, 3);
00554   EPETRA_TEST_ERR(!(ierr == -1), returnierr);
00555   if(ierr == -1)
00556     if(verbose) cout << "Checked OK." << endl;
00557   if(verbose) cout << "\nChecking Reshape(4, -3) - invalid cols" << endl;
00558   ierr = m1.Reshape(4, -3);
00559   EPETRA_TEST_ERR(!(ierr == -1), returnierr);
00560   if(ierr == -1)
00561     if(verbose) cout << "Checked OK." << endl;
00562 
00563 #ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK // only test op() and op[] exceptions if macro is defined.
00564   // out of range index to op() & op[] (6 cases)
00565   int* rand16 = getRandArray(16);
00566   Epetra_IntSerialDenseMatrix m2(View, rand16, 4, 4, 4);
00567 
00568   // op() too high
00569   try {
00570     caught = false;
00571     if(verbose) cout << "\nChecking operator () - row index too high";
00572     ierr = m2(5, 2);
00573   }
00574   catch(int error) {
00575     caught = true;
00576     EPETRA_TEST_ERR(error != -1, returnierr);
00577     if(error == -1)
00578       if(verbose) cout << "Checked OK." << endl;
00579   }
00580   EPETRA_TEST_ERR(!caught, returnierr);
00581   try {
00582     caught = false;
00583     if(verbose) cout << "\nChecking operator () - col index too high";
00584     ierr = m2(3, 4);
00585   }
00586   catch(int error) {
00587     caught = true;
00588     EPETRA_TEST_ERR(error != -2, returnierr);
00589     if(error == -2)
00590       if(verbose) cout << "Checked OK." << endl;
00591   }
00592   EPETRA_TEST_ERR(!caught, returnierr);
00593 
00594   // op() too low
00595   try {
00596     caught = false;
00597     if(verbose) cout << "\nChecking operator () - row index too low";
00598     ierr = m2(-1, 0);
00599   }
00600   catch(int error) {
00601     caught = true;
00602     EPETRA_TEST_ERR(error != -1, returnierr);
00603     if(error == -1)
00604       if(verbose) cout << "Checked OK." << endl;
00605   }
00606   EPETRA_TEST_ERR(!caught, returnierr);
00607   try {
00608     caught = false;
00609     if(verbose) cout << "\nChecking operator () - col index too low";
00610     ierr = m2(0, -1);
00611   }
00612   catch(int error) {
00613     caught = true;
00614     EPETRA_TEST_ERR(error != -2, returnierr);
00615     if(error == -2)
00616       if(verbose) cout << "Checked OK." << endl;
00617   }
00618   EPETRA_TEST_ERR(!caught, returnierr);
00619 
00620   // op[] too high
00621   try { 
00622     caught = false;
00623     if(verbose) cout << "\nChecking operator [] - col index too high";
00624     ierr = m2[4][2];
00625   }
00626   catch(int error) {
00627     caught = true;
00628     EPETRA_TEST_ERR(error != -2, returnierr);
00629     if(error == -2)
00630       if(verbose) cout << "Checked OK." << endl;
00631   }
00632 
00633   // op[] too low
00634   try {
00635     caught = false;
00636     if(verbose) cout << "\nChecking operator [] - col index too low";
00637     ierr = m2[-1][0];
00638   }
00639   catch(int error) {
00640     caught = true;
00641     EPETRA_TEST_ERR(error != -2, returnierr);
00642     if(error == -2)
00643       if(verbose) cout << "Checked OK." << endl;
00644   }
00645   EPETRA_TEST_ERR(!caught, returnierr);
00646      delete [] rand16;
00647 #endif // end of HAVE_EPETRA_ARRAY_BOUNDS_CHECK conditional
00648   
00649   // ISDM = ISDV
00650   Epetra_IntSerialDenseMatrix m3;
00651   Epetra_IntSerialDenseVector v1;
00652   try {
00653     caught = false;
00654     if(verbose) cout << "\nChecking op = - assigning ISDV to ISDM";
00655     m3 = v1;
00656   }
00657   catch(int error) {
00658     caught = true;
00659     EPETRA_TEST_ERR(error != -5, returnierr);
00660     if(error == -5)
00661       if(verbose) cout << "Checked OK." << endl;
00662   }
00663   EPETRA_TEST_ERR(!caught, returnierr);
00664   
00665   return(returnierr);
00666 }
00667 
00668 //=========================================================================
00669 // test matrix operator= (copy & view)
00670 int matrixAssignment(bool verbose, bool debug) {
00671   int ierr = 0;
00672   int returnierr = 0;
00673   if(verbose) printHeading("Testing matrix operator=");
00674 
00675   // each section is in its own block so we can reuse variable names
00676   // lhs = left hand side, rhs = right hand side
00677   
00678   {
00679     // copy->copy (more space needed)
00680     // orig and dup should have same signature
00681     // modifying orig or dup should have no effect on the other
00682     if(verbose) cout << "Checking copy->copy (new alloc)" << endl;
00683     Epetra_IntSerialDenseMatrix lhs(2,2);
00684     int* rand1 = getRandArray(25);
00685     Epetra_IntSerialDenseMatrix rhs(Copy, rand1, 5, 5, 5);
00686     if(debug) {
00687       cout << "before assignment:" << endl;
00688       printMat("rhs",rhs);
00689       printMat("lhs",lhs);
00690     }
00691     lhs = rhs;
00692     if(debug) {
00693       cout << "after assignment:" << endl;
00694       printMat("rhs",rhs);
00695       printMat("lhs",lhs);
00696     }
00697     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
00698     EPETRA_TEST_ERR(!seperateData(rhs,lhs), ierr);
00699     delete[] rand1;
00700   }
00701   returnierr += ierr;
00702   if(ierr == 0)
00703     if(verbose) cout << "Checked OK." << endl;
00704   ierr = 0;
00705   {
00706     // copy->copy (have enough space)
00707     // orig and dup should have same signature
00708     // modifying orig or dup should have no effect on the other
00709     if(verbose) cout << "\nChecking copy->copy (no alloc)" << endl;
00710     int* rand1 = getRandArray(25);
00711     int* rand2 = getRandArray(20);
00712     Epetra_IntSerialDenseMatrix lhs(Copy, rand1, 5, 5, 5);
00713     Epetra_IntSerialDenseMatrix rhs(Copy, rand2, 4, 4, 5);
00714     int* origA = lhs.A();
00715     int origLDA = lhs.LDA();
00716     if(debug) {
00717       cout << "before assignment:" << endl;
00718       printMat("rhs",rhs);
00719       printMat("lhs",lhs);
00720     }
00721     lhs = rhs;
00722     if(debug) {
00723       cout << "after assignment:" << endl;
00724       printMat("rhs",rhs);
00725       printMat("lhs",lhs);
00726     }
00727     // in this case, instead of doing a "normal" LDA test in identSig,
00728     // we do our own. Since we had enough space already, A and LDA should
00729     // not have been changed by the assignment. (The extra parameter to
00730     // identicalSignatures tells it not to test LDA).
00731     EPETRA_TEST_ERR((lhs.A() != origA) || (lhs.LDA() != origLDA), ierr);
00732     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs,false), ierr);
00733     EPETRA_TEST_ERR(!seperateData(rhs,lhs), ierr);
00734      delete [] rand1;
00735      delete [] rand2;
00736   }
00737   returnierr += ierr;
00738   if(ierr == 0)
00739     if(verbose) cout << "Checked OK." << endl;
00740   ierr = 0;
00741   {
00742     // view->copy
00743     // orig and dup should have same signature
00744     // modifying orig or dup should have no effect on the other
00745     if(verbose) cout << "\nChecking view->copy" << endl;
00746     int* rand1 = getRandArray(25);
00747     int* rand2 = getRandArray(64);
00748     Epetra_IntSerialDenseMatrix lhs(View, rand1, 5, 5, 5);
00749     Epetra_IntSerialDenseMatrix rhs(Copy, rand2, 8, 8, 8);
00750     if(debug) {
00751       cout << "before assignment:" << endl;
00752       printMat("rhs",rhs);
00753       printMat("lhs",lhs);
00754     }
00755     lhs = rhs;
00756     if(debug) {
00757       cout << "after assignment:" << endl;
00758       printMat("rhs",rhs);
00759       printMat("lhs",lhs);
00760     }
00761     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
00762     EPETRA_TEST_ERR(!seperateData(rhs,lhs), ierr);
00763     delete[] rand1;
00764     delete[] rand2;
00765   }
00766   returnierr += ierr;
00767   if(ierr == 0)
00768     if(verbose) cout << "Checked OK." << endl;
00769   ierr = 0;
00770   {
00771     // copy->view
00772     // orig and dup should have same signature
00773     // modifying orig or dup should change the other
00774     if(verbose) cout << "\nChecking copy->view" << endl;
00775     int* rand1 = getRandArray(10);
00776     Epetra_IntSerialDenseMatrix lhs(4,4);
00777     Epetra_IntSerialDenseMatrix rhs(View, rand1, 2, 2, 5);
00778     if(debug) {
00779       cout << "before assignment:" << endl;
00780       printMat("rhs",rhs);
00781       printMat("lhs",lhs);
00782     }
00783     lhs = rhs;
00784     if(debug) {
00785       cout << "after assignment:" << endl;
00786       printMat("rhs",rhs);
00787       printMat("lhs",lhs);
00788     }
00789     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
00790     EPETRA_TEST_ERR(seperateData(rhs,lhs), ierr);
00791     delete[] rand1;
00792   }
00793   returnierr += ierr;
00794   if(ierr == 0)
00795     if(verbose) cout << "Checked OK." << endl;
00796   ierr = 0; 
00797   {
00798     // view->view
00799     // orig and dup should have same signature
00800     // modifying orig or dup should change the other
00801     if(verbose) cout << "\nChecking view->view" << endl;
00802     int* rand1 = getRandArray(9);
00803     int* rand2 = getRandArray(18);
00804     Epetra_IntSerialDenseMatrix lhs(View, rand1, 3, 3, 3);
00805     Epetra_IntSerialDenseMatrix rhs(View, rand2, 3, 3, 6);
00806     if(debug) {
00807       cout << "before assignment:" << endl;
00808       printMat("rhs",rhs);
00809       printMat("lhs",lhs);
00810     }
00811     lhs = rhs;
00812     if(debug) {
00813       cout << "after assignment:" << endl;
00814       printMat("rhs",rhs);
00815       printMat("lhs",lhs);
00816     }
00817     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
00818     EPETRA_TEST_ERR(seperateData(rhs,lhs), ierr);
00819     delete[] rand1;
00820     delete[] rand2;
00821   }
00822   returnierr += ierr;
00823   if(ierr == 0)
00824     if(verbose) cout << "Checked OK." << endl;
00825   ierr = 0;
00826   {
00827     // test MakeViewOf
00828     // orig and dup should have same signature except for CV_
00829     // modifying orig or dup should change the other
00830     if(verbose) cout << "\nChecking CrsGraph's usage of MakeViewOf" << endl;
00831     int* rand1 = getRandArray(10);
00832     int* rand2 = getRandArray(10);
00833     Epetra_IntSerialDenseMatrix lhs(Copy, rand1, 5, 5, 2);
00834     Epetra_IntSerialDenseMatrix rhs(Copy, rand2, 5, 5, 2);
00835     if(debug) {
00836       cout << "before assignment:" << endl;
00837       printMat("rhs",rhs);
00838       printMat("lhs",lhs);
00839     }
00840     lhs.MakeViewOf(rhs);
00841     if(debug) {
00842       cout << "after assignment:" << endl;
00843       printMat("rhs",rhs);
00844       printMat("lhs",lhs);
00845     }
00846     EPETRA_TEST_ERR(!(lhs.CV() == View), ierr);
00847     EPETRA_TEST_ERR(!(lhs.M() == rhs.M()), ierr);
00848     EPETRA_TEST_ERR(!(lhs.N() == rhs.N()), ierr);
00849     EPETRA_TEST_ERR(!(lhs.LDA() == rhs.LDA()), ierr);
00850     EPETRA_TEST_ERR(!(lhs.A() == rhs.A()), ierr);
00851     EPETRA_TEST_ERR(seperateData(rhs,lhs), ierr);
00852     delete[] rand1;
00853     delete[] rand2;
00854   }
00855   returnierr += ierr;
00856   if(ierr == 0)
00857     if(verbose) cout << "Checked OK." << endl;
00858   ierr = 0;
00859 
00860   return(returnierr);
00861 }
00862 
00863 //=========================================================================
00864 //=========================================================================
00865 // Vector-testing functions
00866 //=========================================================================
00867 //=========================================================================
00868 // call functions we don't have tests for, for sake of code coverage.
00869 // (the functions are called, but their output isn't checked for correctness).
00870 int vectorCoverage(bool verbose, bool debug) {
00871   if(verbose) printHeading("Testing other vector functions");
00872 
00873   int* rand1 = getRandArray(9);
00874   if(debug) printArray(rand1, 9);
00875   Epetra_IntSerialDenseVector v1(Copy, rand1, 9);
00876   if(debug) printMat("v1",v1);
00877   delete[] rand1;
00878 
00879   if(verbose) cout << "calling one norm" << endl;
00880   int onenorm = v1.OneNorm();
00881   if(debug) cout << "v1.OneNorm() = " << onenorm << endl;
00882 
00883   if(verbose) cout << "calling infinity norm" << endl;
00884   int infnorm = v1.InfNorm();
00885   if(debug) cout << "v1.InfNorm() = " << infnorm << endl;
00886 
00887   if(verbose) cout << "calling random" << endl;
00888   Epetra_IntSerialDenseVector v2(12);
00889   if(debug) printMat("v2 (before)",v2);
00890   v2.Random();
00891   if(debug) printMat("v2 (after)",v2);
00892 
00893   if(verbose) cout << "Checked OK (I think)." << endl;
00894 
00895   return(0);
00896 }
00897 //=========================================================================
00898 // test vector default constructor, user constructor (copy & view),
00899 // () [] operators (read & write), size (larger), resize (smaller)
00900 int vectorCtr(bool verbose, bool debug) {
00901   const int v1size = 5;
00902   const int v1asize = 15;
00903   const int v2size = 10;
00904   const int v3size = 8;
00905   const int v3resize = 5; // should be smaller than v3size
00906 
00907   int ierr = 0;
00908   int returnierr = 0;
00909   if(verbose) printHeading("Testing vector constructors");
00910 
00911   if(verbose) cout << "default constructor" << endl;
00912   Epetra_IntSerialDenseVector v1;
00913   EPETRA_TEST_ERR(!(v1.CV() == Copy), ierr);
00914   if(verbose) cout << "sizing" << endl;
00915   v1.Size(v1size);
00916   EPETRA_TEST_ERR(!(v1.Length() == v1size), ierr);
00917   for(int i = 0; i < v1size; i++) {
00918     EPETRA_TEST_ERR(!(v1(i) == 0), ierr);
00919   }
00920   if(debug) printMat("v1",v1);
00921   returnierr += ierr;
00922   if(ierr == 0) 
00923      if(verbose) cout << "Checked OK." << endl;
00924   ierr = 0;
00925   if(verbose) cout << "\nmanually setting values" << endl;
00926   int* v1rand = getRandArray(v1size);
00927   for(int i = 0; i < v1size; i++)
00928     v1(i) = v1rand[i];
00929   for(int i = 0; i < v1size; i++)
00930     EPETRA_TEST_ERR(!(v1[i] == v1rand[i]), ierr);
00931   if(debug) {
00932     printArray(v1rand, v1size);
00933     printMat("v1",v1);
00934   }
00935   delete[] v1rand;
00936   returnierr += ierr;
00937   if(ierr == 0) 
00938      if(verbose) cout << "Checked OK." << endl;
00939   ierr = 0;
00940 
00941   if(verbose) cout << "\nsized constructor" << endl;
00942   Epetra_IntSerialDenseVector v1a(v1asize);
00943   EPETRA_TEST_ERR(!(v1a.Length() == v1asize), ierr);
00944   EPETRA_TEST_ERR(!(v1a.CV() == Copy),ierr);
00945   if(debug) printMat("v1a", v1a);
00946   returnierr += ierr;
00947   if(ierr == 0) 
00948      if(verbose) cout << "Checked OK." << endl;
00949   ierr = 0;
00950 
00951   if(verbose) cout << "\nuser data constructor (view)" << endl;
00952   int* v2rand = getRandArray(v2size);
00953   if(debug) printArray(v2rand, v2size);
00954   Epetra_IntSerialDenseVector v2(View, v2rand, v2size);
00955   EPETRA_TEST_ERR(!(v2.CV() == View), ierr);
00956   EPETRA_TEST_ERR(!(v2.Length() == v2size), ierr);
00957   EPETRA_TEST_ERR(!(v2.Values() == v2rand), ierr);
00958   if(debug) printMat("v2",v2);
00959   returnierr += ierr;
00960   if(ierr == 0) 
00961      if(verbose) cout << "Checked OK." << endl;
00962   ierr = 0;
00963 
00964   if(verbose) cout << "\nchanging value, checking for correct behavior" << endl;
00965   int* v2randcopy = new int[v2size]; // backup of original values
00966   for(int i = 0; i < v2size; i++)
00967     v2randcopy[i] = v2rand[i];
00968   v2(4) = v2(4) + 1; // magic number for which element to change
00969   v2randcopy[4] = v2randcopy[4] +1;
00970   for(int i = 0; i < v2size; i++)
00971     EPETRA_TEST_ERR(!(v2rand[i] == v2randcopy[i]), ierr); // v2rand should have updated correctly
00972   if(debug) {
00973     printArray(v2rand, v2size);
00974     printMat("v2",v2);
00975   }
00976   delete[] v2rand;
00977   delete[] v2randcopy;
00978   returnierr += ierr;
00979   if(ierr == 0) 
00980      if(verbose) cout << "Checked OK." << endl;
00981   ierr = 0;
00982 
00983   if(verbose) cout << "\nuser data constructor (copy)" << endl;
00984   int* v3rand = getRandArray(v3size);
00985   if(debug) printArray(v3rand, v3size);
00986   int* v3randcopy = new int[v3size];
00987   for(int i = 0; i < v3size; i++)
00988     v3randcopy[i] = v3rand[i];
00989   Epetra_IntSerialDenseVector v3(Copy, v3rand, v3size);
00990   if(debug) printMat("v3",v3);
00991   if(verbose) cout << "checking to see if data initialized correctly" << endl;
00992   EPETRA_TEST_ERR(!(v3.CV() == Copy), ierr);
00993   EPETRA_TEST_ERR(!(v3.Length() == v3size), ierr);
00994   EPETRA_TEST_ERR(!(v3.Values() != v3rand), ierr); // should not be a view
00995   for(int i = 0; i < v3size; i++)
00996     EPETRA_TEST_ERR(!(v3[i] == v3rand[i]), ierr); // data should be identical to user array
00997   returnierr += ierr;
00998   if(ierr == 0) 
00999      if(verbose) cout << "Checked OK." << endl;
01000   ierr = 0;
01001 
01002   if(verbose) cout << "\nmodifying entry" << endl;
01003   v3[5] = v3[5] + 3; // magic number for which element to change
01004   for(int i = 0; i < v3size; i++)
01005     EPETRA_TEST_ERR(!(v3rand[i] == v3randcopy[i]), ierr); // user array should be unchanged
01006   v3rand[5] = v3rand[5] + 3; // same magic modification performed to user's array
01007   for(int i = 0; i < v3size; i++)
01008     EPETRA_TEST_ERR(!(v3[i] == v3rand[i]), ierr); // should equal user array with same modification performed
01009   if(debug) {
01010     printArray(v3rand, v3size);
01011     printMat("v3",v3);
01012   }
01013   returnierr += ierr;
01014   if(ierr == 0) 
01015      if(verbose) cout << "Checked OK." << endl;
01016   ierr = 0;
01017 
01018   if(verbose) cout << "\nresizing" << endl;
01019   v3.Resize(v3resize);
01020   EPETRA_TEST_ERR(!(v3.Length() == v3resize), ierr);
01021   for(int i = 0; i < v3resize; i++)
01022     EPETRA_TEST_ERR(!(v3[i] == v3rand[i]),ierr);
01023   if(debug) printMat("v3",v3);
01024   delete[] v3rand;
01025   delete[] v3randcopy;
01026   returnierr += ierr;
01027   if(ierr == 0) 
01028      if(verbose) cout << "Checked OK." << endl;
01029   ierr = 0;
01030 
01031   if(verbose) cout << "\nChecking pointer on zero-sized vector" << endl;
01032   int* before = v3.Values();
01033   if(verbose) cout << "Resizing to 0" << endl;
01034   if(debug) cout << "Before = " << before << endl;
01035   EPETRA_TEST_ERR(!(before != 0), ierr);
01036   v3.Resize(0);
01037   int* after = v3.Values();
01038   EPETRA_TEST_ERR(!(after == 0), ierr);
01039   if(debug) cout << "After = " << after << endl;
01040   v3.Size(3);
01041   before = v3.Values();
01042   if(verbose) cout << "Sizing to 0" << endl;
01043   if(debug) cout << "Before = " << before << endl;
01044   EPETRA_TEST_ERR(!(before != 0), ierr);
01045   v3.Size(0);
01046   after = v3.Values();
01047   EPETRA_TEST_ERR(!(after == 0), ierr);
01048   if(debug) cout << "After = " << after << endl;
01049   returnierr += ierr;
01050   if(ierr == 0) 
01051      if(verbose) cout << "Checked OK." << endl;
01052   ierr = 0;
01053 
01054   return(returnierr);
01055 }
01056 //=========================================================================
01057 // test vector copy constructor (copy & view)
01058 int vectorCpyCtr(bool verbose, bool debug) {
01059   const int v1size = 15;
01060   const int v2size = 12;
01061 
01062   int ierr = 0;
01063   int returnierr = 0;
01064   if(verbose) printHeading("Testing vector copy constructors");
01065 
01066   if(verbose) cout << "checking copy constructor (view)" << endl;
01067   int* v1rand = getRandArray(v1size);
01068   if(debug) printArray(v1rand, v1size);
01069   Epetra_IntSerialDenseVector v1(View, v1rand, v1size);
01070   if(debug) {
01071     cout << "original vector:" << endl;
01072     printMat("v1",v1);
01073   }
01074   Epetra_IntSerialDenseVector v1clone(v1);
01075   if(debug) {
01076     cout << "clone vector:" << endl;
01077     printMat("v1clone",v1clone);
01078   }
01079   if(verbose) cout << "making sure signatures match" << endl;
01080   EPETRA_TEST_ERR(!identicalSignatures(v1, v1clone), ierr);
01081   delete[] v1rand;
01082   returnierr += ierr;
01083   if(ierr == 0)
01084     if(verbose) cout << "Checked OK." << endl;
01085   ierr = 0;
01086 
01087   if(verbose) cout << "\nchecking copy constructor (copy)" << endl;
01088   int* v2rand = getRandArray(v2size);
01089   if(debug) printArray(v2rand, v2size);
01090   Epetra_IntSerialDenseVector v2(Copy, v2rand, v2size);
01091   if(debug) {
01092     cout << "original vector:" << endl;
01093     printMat("v2",v2);
01094   }
01095   Epetra_IntSerialDenseVector v2clone(v2);
01096   if(debug) {
01097     cout << "clone vector:" << endl;
01098     printMat("v2clone",v2clone);
01099   }
01100   if(verbose) cout << "checking that signatures match" << endl;
01101   EPETRA_TEST_ERR(!identicalSignatures(v2, v2clone), ierr);
01102   returnierr += ierr;
01103   if(ierr == 0)
01104     if(verbose) cout << "Checked OK." << endl;
01105   ierr = 0;
01106 
01107   if(verbose) cout << "\nmodifying entry in v2, v2clone should be unchanged" << endl;
01108   EPETRA_TEST_ERR(!seperateData(v2, v2clone), ierr);
01109   if(debug) {
01110     printArray(v2rand, v2size);
01111     cout << "orig:" << endl;
01112     printMat("v2",v2);
01113     cout << "clone:" << endl;
01114     printMat("v2clone",v2clone);
01115   }
01116   delete[] v2rand;
01117   returnierr += ierr;
01118   if(ierr == 0)
01119     if(verbose) cout << "Checked OK." << endl;
01120   ierr = 0;
01121 
01122   return(returnierr);
01123 }
01124 //=========================================================================
01125 // test vector operator= (copy & view)
01126 int vectorAssignment(bool verbose, bool debug) {
01127   int ierr = 0;
01128   int returnierr = 0;
01129   if(verbose) printHeading("Testing vector operator=");
01130 
01131   // each section is in its own block so we can reuse variable names
01132   // lhs = left hand side, rhs = right hand side
01133   
01134   {
01135     // copy->copy (more space needed)
01136     // orig and dup should have same signature
01137     // modifying orig or dup should have no effect on the other
01138     if(verbose) cout << "Checking copy->copy (new alloc)" << endl;
01139     Epetra_IntSerialDenseVector lhs(2);
01140     int* rand1 = getRandArray(10);
01141     Epetra_IntSerialDenseVector rhs(Copy, rand1, 10);
01142     if(debug) {
01143       cout << "before assignment:" << endl;
01144       printMat("rhs",rhs);
01145       printMat("lhs",lhs);
01146     }
01147     lhs = rhs;
01148     if(debug) {
01149       cout << "after assignment:" << endl;
01150       printMat("rhs",rhs);
01151       printMat("lhs",lhs);
01152     }
01153     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
01154     EPETRA_TEST_ERR(!seperateData(rhs,lhs), ierr);
01155     delete[] rand1;
01156   }
01157   returnierr += ierr;
01158   if(ierr == 0)
01159     if(verbose) cout << "Checked OK." << endl;
01160   ierr = 0;
01161   {
01162     // copy->copy (have enough space)
01163     // orig and dup should have same signature
01164     // modifying orig or dup should have no effect on the other
01165     if(verbose) cout << "\nChecking copy->copy (no alloc)" << endl;
01166     int* rand1 = getRandArray(20);
01167     int* rand2 = getRandArray(15);
01168     Epetra_IntSerialDenseVector lhs(Copy, rand1, 20);
01169     Epetra_IntSerialDenseVector rhs(Copy, rand2, 15);
01170     int* origA = lhs.A();
01171     int origLDA = lhs.LDA();
01172     if(debug) {
01173       cout << "before assignment:" << endl;
01174       printMat("rhs",rhs);
01175       printMat("lhs",lhs);
01176     }
01177     lhs = rhs;
01178     if(debug) {
01179       cout << "after assignment:" << endl;
01180       printMat("rhs",rhs);
01181       printMat("lhs",lhs);
01182     }
01183     // in this case, instead of doing a "normal" LDA test in identSig,
01184     // we do our own. Since we had enough space already, A and LDA should
01185     // not have been changed by the assignment. (The extra parameter to
01186     // identicalSignatures tells it not to test LDA).
01187     EPETRA_TEST_ERR((lhs.A() != origA) || (lhs.LDA() != origLDA), ierr);
01188     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs,false), ierr);
01189     EPETRA_TEST_ERR(!seperateData(rhs,lhs), ierr);
01190     delete[] rand1;
01191     delete[] rand2;
01192   }
01193   returnierr += ierr;
01194   if(ierr == 0)
01195     if(verbose) cout << "Checked OK." << endl;
01196   ierr = 0;
01197   {
01198     // view->copy
01199     // orig and dup should have same signature
01200     // modifying orig or dup should have no effect on the other
01201     if(verbose) cout << "\nChecking view->copy" << endl;
01202     int* rand1 = getRandArray(5);
01203     int* rand2 = getRandArray(8);
01204     Epetra_IntSerialDenseVector lhs(View, rand1, 5);
01205     Epetra_IntSerialDenseVector rhs(Copy, rand2, 8);
01206     if(debug) {
01207       cout << "before assignment:" << endl;
01208       printMat("rhs",rhs);
01209       printMat("lhs",lhs);
01210     }
01211     lhs = rhs;
01212     if(debug) {
01213       cout << "after assignment:" << endl;
01214       printMat("rhs",rhs);
01215       printMat("lhs",lhs);
01216     }
01217     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
01218     EPETRA_TEST_ERR(!seperateData(rhs,lhs), ierr);
01219     delete[] rand1;
01220     delete[] rand2;
01221   }
01222   returnierr += ierr;
01223   if(ierr == 0)
01224     if(verbose) cout << "Checked OK." << endl;
01225   ierr = 0;
01226   {
01227     // copy->view
01228     // orig and dup should have same signature
01229     // modifying orig or dup should change the other
01230     if(verbose) cout << "\nChecking copy->view" << endl;
01231     int* rand1 = getRandArray(10);
01232     Epetra_IntSerialDenseVector lhs(4);
01233     Epetra_IntSerialDenseVector rhs(View, rand1, 10);
01234     if(debug) {
01235       cout << "before assignment:" << endl;
01236       printMat("rhs",rhs);
01237       printMat("lhs",lhs);
01238     }
01239     lhs = rhs;
01240     if(debug) {
01241       cout << "after assignment:" << endl;
01242       printMat("rhs",rhs);
01243       printMat("lhs",lhs);
01244     }
01245     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
01246     EPETRA_TEST_ERR(seperateData(rhs,lhs), ierr);
01247     delete[] rand1;
01248   }
01249   returnierr += ierr;
01250   if(ierr == 0)
01251     if(verbose) cout << "Checked OK." << endl;
01252   ierr = 0;
01253   {
01254     // view->view
01255     // orig and dup should have same signature
01256     // modifying orig or dup should change the other
01257     if(verbose) cout << "\nChecking view->view" << endl;
01258     int* rand1 = getRandArray(9);
01259     int* rand2 = getRandArray(11);
01260     Epetra_IntSerialDenseVector lhs(View, rand1, 9);
01261     Epetra_IntSerialDenseVector rhs(View, rand2, 11);
01262     if(debug) {
01263       cout << "before assignment:" << endl;
01264       printMat("rhs",rhs);
01265       printMat("lhs",lhs);
01266     }
01267     lhs = rhs;
01268     if(debug) {
01269       cout << "after assignment:" << endl;
01270       printMat("rhs",rhs);
01271       printMat("lhs",lhs);
01272     }
01273     EPETRA_TEST_ERR(!identicalSignatures(rhs,lhs), ierr);
01274     EPETRA_TEST_ERR(seperateData(rhs,lhs), ierr);
01275     delete[] rand1;
01276     delete[] rand2;
01277   }
01278   returnierr += ierr;
01279   if(ierr == 0)
01280     if(verbose) cout << "Checked OK." << endl;
01281   ierr = 0;
01282   {
01283     // test MakeViewOf
01284     // orig and dup should have same signature except for CV_
01285     // modifying orig or dup should change the other
01286     if(verbose) cout << "\nChecking CrsGraph's usage of MakeViewOf" << endl;
01287     int* rand1 = getRandArray(10);
01288     int* rand2 = getRandArray(10);
01289     Epetra_IntSerialDenseVector lhs(Copy, rand1, 10);
01290     Epetra_IntSerialDenseVector rhs(Copy, rand2, 10);
01291     if(debug) {
01292       cout << "before assignment:" << endl;
01293       printMat("rhs",rhs);
01294       printMat("lhs",lhs);
01295     }
01296     lhs.MakeViewOf(rhs);
01297     if(debug) {
01298       cout << "after assignment:" << endl;
01299       printMat("rhs",rhs);
01300       printMat("lhs",lhs);
01301     }
01302     EPETRA_TEST_ERR(!(lhs.CV() == View), ierr);
01303     EPETRA_TEST_ERR(!(lhs.M() == rhs.M()), ierr);
01304     EPETRA_TEST_ERR(!(lhs.N() == rhs.N()), ierr);
01305     EPETRA_TEST_ERR(!(lhs.LDA() == rhs.LDA()), ierr);
01306     EPETRA_TEST_ERR(!(lhs.A() == rhs.A()), ierr);
01307     EPETRA_TEST_ERR(seperateData(rhs,lhs), ierr);
01308     delete[] rand1;
01309     delete[] rand2;
01310   }
01311   returnierr += ierr;
01312   if(ierr == 0)
01313     if(verbose) cout << "Checked OK." << endl;
01314   ierr = 0;
01315 
01316   return(returnierr);
01317 }
01318 
01319 //=========================================================================
01320 // test vector error-reporting
01321 int vectorExceptions(bool verbose, bool debug) {
01322   int returnierr = 0;
01323   int ierr = 0;
01324   bool caught = false;
01325   Epetra_IntSerialDenseVector* vector;
01326 
01327   if(verbose) printHeading("Testing vector error-reporting.\nExpect error messages if EPETRA_NO_ERROR_REPORTS is not defined.");
01328   
01329   try { // invalid dimension to sized ctr
01330     caught = false;
01331     if(verbose) cout << "Checking Epetra_IntSerialDenseVector(-1)";
01332     vector = new Epetra_IntSerialDenseVector(-1);
01333   }
01334   catch(int error) {
01335     caught = true;
01336     EPETRA_TEST_ERR(error != -1, returnierr);
01337     if(error == -1)
01338       if(verbose) cout << "Checked OK." << endl;
01339   }
01340   EPETRA_TEST_ERR(!caught, returnierr);
01341 
01342      int* rand2 = getRandArray(2);
01343   try { // invalid dimension to user-data ctr
01344     caught = false;
01345     if(verbose) cout << "\nChecking Epetra_IntSerialDenseVector(Copy, int*, -3)";
01346     vector = new Epetra_IntSerialDenseVector(Copy, rand2, -3);
01347   }
01348   catch(int error) {
01349     caught = true;
01350     EPETRA_TEST_ERR(error != -1, returnierr);
01351     if(error == -1)
01352       if(verbose) cout << "Checked OK." << endl;
01353   }
01354   EPETRA_TEST_ERR(!caught, returnierr);
01355      delete[] rand2;
01356 
01357   try { // null pointer to user-data ctr
01358     caught = false;
01359     if(verbose) cout << "\nChecking Epetra_IntSerialDenseVector(Copy, 0, 5)";
01360     vector = new Epetra_IntSerialDenseVector(Copy, 0, 5);
01361   }
01362   catch(int error) {
01363     caught = true;
01364     EPETRA_TEST_ERR(error != -3, returnierr);
01365     if(error == -3)
01366       if(verbose) cout << "Checked OK." << endl;
01367   }
01368   EPETRA_TEST_ERR(!caught, returnierr);
01369 
01370   // invalid parameter to size
01371   if(verbose) cout << "\nChecking Size(-2)" << endl;
01372   Epetra_IntSerialDenseVector v1;
01373   ierr = v1.Size(-2);
01374   EPETRA_TEST_ERR(!(ierr == -1), returnierr);
01375   if(ierr == -1)
01376     if(verbose) cout << "Checked OK." << endl;
01377 
01378   // invalid parameter to resize
01379   if(verbose) cout << "\nChecking Resize(-4)" << endl;
01380   v1.Size(5);
01381   ierr = v1.Resize(-4);
01382   EPETRA_TEST_ERR(!(ierr == -1), returnierr);
01383   if(ierr == -1)
01384     if(verbose) cout << "Checked OK." << endl;
01385 
01386 #ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK // only test op() and op[] exceptions if macro is defined.
01387   // out of range index to op() & op[]
01388   int* rand17 = getRandArray(17);
01389   Epetra_IntSerialDenseVector v2(View, rand17, 17);
01390   try { // op() too high
01391     caught = false;
01392     if(verbose) cout << "\nChecking operator () - index too high";
01393     ierr = v2(17);
01394   }
01395   catch(int error) {
01396     caught = true;
01397     EPETRA_TEST_ERR(error != -1, returnierr);
01398     if(error == -1)
01399       if(verbose) cout << "Checked OK." << endl;
01400   }
01401   EPETRA_TEST_ERR(!caught, returnierr);
01402   
01403   try { // op() too low
01404     caught = false;
01405     if(verbose) cout << "\nChecking operator () - index too low";
01406     ierr = v2(-1);
01407   }
01408   catch(int error) {
01409     caught = true;
01410     EPETRA_TEST_ERR(error != -1, returnierr);
01411     if(error == -1)
01412       if(verbose) cout << "Checked OK." << endl;
01413   }
01414   EPETRA_TEST_ERR(!caught, returnierr);
01415 
01416   try { // op[] too high
01417     caught = false;
01418     if(verbose) cout << "\nChecking operator [] - index too high";
01419     ierr = v2[17];
01420   }
01421   catch(int error) {
01422     caught = true;
01423     EPETRA_TEST_ERR(error != -1, returnierr);
01424     if(error == -1)
01425       if(verbose) cout << "Checked OK." << endl;
01426   }
01427   EPETRA_TEST_ERR(!caught, returnierr);
01428   
01429   try { // op[] too low
01430     caught = false;
01431     if(verbose) cout << "\nChecking operator [] - index too low";
01432     ierr = v2[-1];
01433   }
01434   catch(int error) {
01435     caught = true;
01436     EPETRA_TEST_ERR(error != -1, returnierr);
01437     if(error == -1)
01438       if(verbose) cout << "Checked OK." << endl;
01439   }
01440   EPETRA_TEST_ERR(!caught, returnierr);
01441   delete[] rand17;
01442 #endif // end of HAVE_EPETRA_ARRAY_BOUNDS_CHECK conditional
01443 
01444   // we don't need to check for ISDV = ISDM, as that is a compile-time error
01445   
01446   return(returnierr);
01447 }
01448 
01449 //=========================================================================
01450 //=========================================================================
01451 // helper functions
01452 //=========================================================================
01453 //=========================================================================
01454 // checks the signatures of two matrices
01455 bool identicalSignatures(Epetra_IntSerialDenseMatrix& a, Epetra_IntSerialDenseMatrix& b, bool testLDA) {
01456   /*cout << "M: " << a.M() << "  " << b.M() << endl;
01457   cout << "N: " << a.N() << "  " << b.N() << endl;
01458   cout << "LDA: " << a.LDA() << "  " << b.LDA() << endl;
01459   cout << "CV: " << a.CV() << "  " << b.CV() << endl;
01460   cout << "A: " << a.A() << "  " << b.A() << endl;*/
01461 
01462   if((a.M()  != b.M()  )|| // check properties first
01463      (a.N()  != b.N()  )||
01464      (a.CV() != b.CV() ))
01465     return(false);
01466 
01467   if(testLDA == true)      // if we are coming from op= c->c #2 (have enough space)
01468     if(a.LDA() != b.LDA()) // then we don't check LDA (but we do check it in the test function)
01469       return(false);
01470 
01471   if(a.CV() == View) { // if we're still here, we need to check the data
01472     if(a.A() != b.A()) // for a view, this just means checking the pointers
01473       return(false);   // for a copy, this means checking each element
01474   }
01475   else { // CV == Copy
01476     const int m = a.M();
01477     const int n = a.N();
01478     for(int i = 0; i < m; i++)
01479       for(int j = 0; j < n; j++) {
01480         if(a(i,j) != b(i,j))
01481           return(false);
01482       }
01483   }
01484 
01485   return(true); // if we're still here, signatures are identical
01486 }
01487 //=========================================================================
01488 // checks if two matrices are independent or not
01489 bool seperateData(Epetra_IntSerialDenseMatrix& a, Epetra_IntSerialDenseMatrix& b) {
01490   bool seperate;
01491 
01492   int r = EPETRA_MIN(a.M(),b.M()) / 2; // ensures (r,c) is valid
01493   int c = EPETRA_MIN(a.N(),b.N()) / 2; // in both matrices
01494 
01495   int orig_a = a(r,c);
01496   int new_value = a(r,c) + 1;
01497   if(b(r,c) == new_value) // there's a chance b could be independent, but
01498     new_value++;          // already have new_value in (r,c).
01499   
01500   a(r,c) = new_value;
01501   if(b(r,c) == new_value)
01502     seperate = false;
01503   else
01504     seperate = true;
01505 
01506   a(r,c) = orig_a; // undo change we made to a
01507 
01508   return(seperate);
01509 }
01510 //=========================================================================
01511 // returns a int* array of a given length, with random values on interval [0,100]
01512 int* getRandArray(int length) {
01513   int* array = new int[length];
01514   for(int i = 0; i < length; i++)
01515     array[i] = randomInt();
01516 
01517   return(array);
01518 }
01519 //=========================================================================
01520 // returns a random integer on the interval [0-maxint).
01521 // this is the same generator used in IntSerialDenseMatrix
01522 int randomInt() {
01523   const int maxint = 100;
01524 
01525   const double a = 16807.0;
01526   const double BigInt = 2147483647.0;
01527   double seed = rand(); // Use POSIX standard random function;
01528 
01529   seed = fmod(a * seed, BigInt); // fmod returns remainder of floating point division
01530                                  // (a * seed) - (floor(a * seed / BigInt) * BigInt)
01531   double randdouble = (seed / BigInt); // should be [0,1)
01532   int randint = int(randdouble * maxint);
01533   
01534   return(randint);
01535 }
01536 //=========================================================================
01537 // prints int* array with formatting
01538 void printArray(int* array, int length) {
01539   cout << "user array (size " << length << "): ";
01540   for(int i = 0; i < length; i++)
01541     cout << array[i] << "  ";
01542   cout << endl;
01543 }
01544 //=========================================================================
01545 // prints IntSerialDenseMatrix/Vector with formatting
01546 void printMat(const char* name, Epetra_IntSerialDenseMatrix& matrix) {
01547   //cout << "--------------------" << endl;
01548   cout << "*** " << name << " ***" << endl;
01549   cout << matrix;
01550   //cout << "--------------------" << endl;
01551 }
01552 
01553 //=========================================================================
01554 // prints section heading with spacers/formatting
01555 void printHeading(const char* heading) {
01556   cout << "\n==================================================================\n";
01557   cout << heading << endl;
01558   cout << "==================================================================\n";
01559 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines