cxx_tmpl_main.cpp

Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //                    Teuchos: Common Tools Package
00005 //                 Copyright (2004) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 // Kris
00030 // 07.24.03 -- Initial checkin
00031 // 08.08.03 -- All test suites except for TRSM are finished.
00032 // 08.14.03 -- The test suite for TRSM is finished (Heidi).
00033 /*
00034 
00035 This test program is intended to check an experimental default type (e.g. mp_real) against an "officialy supported" control type (e.g. double). For each test, the program will generate the appropriate scalars and randomly-sized vectors and matrices of random numbers for the routine being tested. All of the input data for the experimental type is casted into the control type, so in theory both BLAS routines should get the same input data. Upon return, the experimental output data is casted back into the control type, and the results are compared; if they are equal (within a user-definable tolerance) the test is considered successful.
00036 
00037 The test routine for TRSM is still being developed; all of the others are more or less finalized.
00038 
00039 */
00040 
00041 #include "Teuchos_BLAS.hpp"
00042 #include "Teuchos_Time.hpp"
00043 #include "Teuchos_Version.hpp"
00044 #include "Teuchos_GlobalMPISession.hpp"
00045 
00046 using Teuchos::BLAS;
00047 using Teuchos::ScalarTraits;
00048 
00049 // SType1 and SType2 define the datatypes for which BLAS output will be compared.
00050 // SType2 should generally be a control datatype "officially" supported by the BLAS; SType1 should be the experimental type being checked.
00051 
00052 // Define the first scalar type
00053 #ifdef HAVE_TEUCHOS_COMPLEX
00054 #define SType1     std::complex<float>
00055 #else
00056 #define SType1     float
00057 #endif
00058 
00059 // Define the second scalar type
00060 #ifdef HAVE_TEUCHOS_COMPLEX
00061 #define SType2     std::complex<double>
00062 #else
00063 #define SType2     double
00064 #endif
00065 
00066 // Define the ordinal type
00067 #define OType    long int 
00068 
00069 // MVMIN/MAX define the minimum and maximum dimensions of generated matrices and vectors, respectively.
00070 #define MVMIN      2
00071 #define MVMAX      20
00072 // SCALARMAX defines the maximum positive value (with a little leeway) generated for matrix and std::vector elements and scalars:
00073 // random numbers in [-SCALARMAX, SCALARMAX] will be generated.
00074 // Set SCALARMAX to a floating-point value (e.g. 10.0) to enable floating-point random number generation, such that
00075 // random numbers in (-SCALARMAX - 1, SCALARMAX + 1) will be generated.
00076 // Large values of SCALARMAX may cause problems with SType2 = int, as large integer values will overflow floating-point types.
00077 #define SCALARMAX  10
00078 // These define the number of tests to be run for each individual BLAS routine.
00079 #define ROTGTESTS  5
00080 #define ASUMTESTS  5
00081 #define AXPYTESTS  5
00082 #define COPYTESTS  5
00083 #define DOTTESTS   5
00084 #define IAMAXTESTS 5
00085 #define NRM2TESTS  5
00086 #define SCALTESTS  5
00087 #define GEMVTESTS  5
00088 #define GERTESTS   5
00089 #define TRMVTESTS  5
00090 #define GEMMTESTS  5
00091 #define SYMMTESTS  5
00092 #define TRMMTESTS  5
00093 #define TRSMTESTS  5
00094 
00095 // Returns ScalarTraits<TYPE>::random() (the input parameters are ignored)
00096 template<typename TYPE>
00097 TYPE GetRandom(TYPE, TYPE);
00098 
00099 // Returns a random integer between the two input parameters, inclusive
00100 template<>
00101 int GetRandom(int, int);
00102 
00103 // Returns a random double between the two input parameters, plus or minus a random number between 0 and 1
00104 template<>
00105 double GetRandom(double, double);
00106 
00107 template<typename TYPE>
00108 void PrintVector(TYPE* Vector, OType Size, std::string Name, bool Matlab = 0);
00109 
00110 template<typename TYPE>
00111 void PrintMatrix(TYPE* Matrix, OType Rows, OType Columns, OType LDM, std::string Name, bool Matlab = 0);
00112 
00113 template<typename TYPE1, typename TYPE2>
00114 bool CompareScalars(TYPE1 Scalar1, TYPE2 Scalar2, typename ScalarTraits<TYPE2>::magnitudeType Tolerance ); 
00115 
00116 template<typename TYPE1, typename TYPE2>
00117 bool CompareVectors(TYPE1* Vector1, TYPE2* Vector2, OType Size, typename ScalarTraits<TYPE2>::magnitudeType Tolerance ); 
00118 
00119 template<typename TYPE1, typename TYPE2>
00120 bool CompareMatrices(TYPE1* Matrix1, TYPE2* Matrix2, OType Rows, OType Columns, OType LDM, typename ScalarTraits<TYPE2>::magnitudeType Tolerance ); 
00121 
00122 // For most types, this function is just a wrapper for static_cast(), but for mp_real/double, it calls mp::dble()
00123 // The second input parameter is not used; it is only needed to determine what type to convert *to*
00124 template<typename TYPE1, typename TYPE2>
00125 TYPE2 ConvertType(TYPE1, TYPE2);
00126 
00127 // These functions return a random character appropriate for the BLAS arguments that share their names (uses GetRandom())
00128 Teuchos::ESide RandomSIDE();
00129 Teuchos::EUplo RandomUPLO();
00130 Teuchos::ETransp RandomTRANS();
00131 Teuchos::EDiag RandomDIAG();
00132 
00133 int main(int argc, char *argv[])
00134 {
00135   Teuchos::GlobalMPISession mpiSession(&argc, &argv);
00136   bool verbose = 0;
00137   bool debug = 0;
00138   bool matlab = 0;
00139   bool InvalidCmdLineArgs = 0;
00140   OType i, j, k;
00141   for(i = 1; i < argc; i++)
00142     {
00143       if(argv[i][0] == '-')
00144   {
00145     switch(argv[i][1])
00146       {
00147       case 'v':
00148         if(!verbose)
00149     {
00150       verbose = 1;
00151     }
00152         else
00153     {
00154       InvalidCmdLineArgs = 1;
00155     }
00156         break;
00157       case 'd':
00158         if(!debug)
00159     {
00160       debug = 1;
00161     }
00162         else
00163     {
00164       InvalidCmdLineArgs = 1;
00165     }
00166         break;
00167       case 'm':
00168         if(!matlab)
00169     {
00170       matlab = 1;
00171     }
00172         else
00173     {
00174       InvalidCmdLineArgs = 1;
00175     }
00176         break;
00177       default:
00178         InvalidCmdLineArgs = 1;
00179         break;
00180       }
00181   }
00182     }
00183 
00184   if (verbose)
00185     std::cout << Teuchos::Teuchos_Version() << std::endl << std::endl;
00186 
00187   if(InvalidCmdLineArgs || (argc > 4))
00188     {
00189       std::cout << "Invalid command line arguments detected. Use the following flags:" << std::endl
00190      << "\t -v enables verbose mode (reports number of failed/successful tests)" << std::endl
00191      << "\t -d enables debug mode (same as verbose with output of each test, not recommended for large numbers of tests)" << std::endl
00192      << "\t -m enables matlab-style output; only has an effect if debug mode is enabled" << std::endl;
00193       return 1;
00194     }
00195   typedef ScalarTraits<SType1>::magnitudeType MType1;
00196   typedef ScalarTraits<SType2>::magnitudeType MType2; 
00197   BLAS<OType, SType1> SType1BLAS;
00198   BLAS<OType, SType2> SType2BLAS;
00199   SType1 SType1zero = ScalarTraits<SType1>::zero();
00200   SType1 SType1one = ScalarTraits<SType1>::one();
00201   SType2 SType2one = ScalarTraits<SType2>::one();
00202   SType1* SType1A;
00203   SType1* SType1B;
00204   SType1* SType1C;
00205   SType1* SType1x;
00206   SType1* SType1y;
00207   SType1 SType1alpha, SType1beta;
00208   SType2* SType2A;
00209   SType2* SType2B;
00210   SType2* SType2C;
00211   SType2* SType2x;
00212   SType2* SType2y; 
00213   SType2 SType2alpha, SType2beta;
00214   SType1 SType1ASUMresult, SType1DOTresult, SType1NRM2result, SType1SINresult;
00215   SType2 SType2ASUMresult, SType2DOTresult, SType2NRM2result, SType2SINresult;
00216   MType1 SType1COSresult;
00217   MType2 SType2COSresult;
00218   OType incx, incy;
00219   OType SType1IAMAXresult;
00220   OType SType2IAMAXresult;
00221   OType TotalTestCount = 1, GoodTestSubcount, GoodTestCount = 0, M, M2, N, N2, P, LDA, LDB, LDC, Mx, My;
00222   Teuchos::EUplo UPLO;
00223   Teuchos::ESide SIDE;
00224   Teuchos::ETransp TRANS, TRANSA, TRANSB;
00225   Teuchos::EDiag DIAG;
00226   SType2 convertTo = ScalarTraits<SType2>::zero();
00227   MType2 TOL = 1e-5*ScalarTraits<MType2>::one();
00228   
00229   std::srand(time(NULL));
00230 
00231   //--------------------------------------------------------------------------------
00232   // BEGIN LEVEL 1 BLAS TESTS
00233   //--------------------------------------------------------------------------------
00234   // Begin ROTG Tests
00235   //--------------------------------------------------------------------------------
00236   GoodTestSubcount = 0;
00237   for(i = 0; i < ROTGTESTS; i++)
00238     {
00239       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
00240       SType2alpha = ConvertType(SType1alpha, convertTo);
00241       SType1beta = GetRandom(-SCALARMAX, SCALARMAX);
00242       SType2beta = ConvertType(SType1beta, convertTo);
00243       SType1COSresult = ScalarTraits<MType1>::zero();
00244       SType2COSresult = ConvertType(SType1COSresult, ScalarTraits<MType2>::zero());
00245       SType1SINresult = ScalarTraits<SType1>::zero();
00246       SType2SINresult = ConvertType(SType1SINresult, convertTo);
00247       
00248       if(debug)
00249   {
00250     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00251     std::cout << "SType1alpha = "  << SType1alpha << std::endl;
00252     std::cout << "SType2alpha = " << SType2alpha << std::endl;
00253     std::cout << "SType1beta = "  << SType1beta << std::endl;
00254     std::cout << "SType2beta = " << SType2beta << std::endl;
00255   }
00256       TotalTestCount++;
00257       SType1BLAS.ROTG(&SType1alpha, &SType1beta, &SType1COSresult, &SType1SINresult);
00258       SType2BLAS.ROTG(&SType2alpha, &SType2beta, &SType2COSresult, &SType2SINresult);
00259       if(debug)
00260   {
00261     std::cout << "SType1 ROTG COS result: " << SType1COSresult << std::endl;
00262     std::cout << "SType2 ROTG COS result: " << SType2COSresult << std::endl;
00263     std::cout << "SType1 ROTG SIN result: " << SType1SINresult << std::endl;
00264     std::cout << "SType2 ROTG SIN result: " << SType2SINresult << std::endl;
00265   }
00266       GoodTestSubcount += ( CompareScalars(SType1COSresult, SType2COSresult, TOL) && 
00267           CompareScalars(SType1SINresult, SType2SINresult, TOL) );
00268     }
00269   GoodTestCount += GoodTestSubcount;
00270   if(verbose || debug) std::cout << "ROTG: " << GoodTestSubcount << " of " << ROTGTESTS << " tests were successful." << std::endl;
00271   if(debug) std::cout << std::endl;
00272   //--------------------------------------------------------------------------------
00273   // End ROTG Tests
00274   //--------------------------------------------------------------------------------
00275 
00276   //--------------------------------------------------------------------------------
00277   // Begin ASUM Tests
00278   //--------------------------------------------------------------------------------
00279   GoodTestSubcount = 0;
00280   ScalarTraits<int>::seedrandom(0);
00281   for(i = 0; i < ASUMTESTS; i++)
00282     {
00283       incx = GetRandom(1, SCALARMAX);
00284       M = GetRandom(MVMIN, MVMAX);
00285       M2 = M*incx;
00286       SType1x = new SType1[M2];
00287       SType2x = new SType2[M2];
00288       for(j = 0; j < M2; j++)
00289   {
00290     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00291     SType2x[j] = ConvertType(SType1x[j], convertTo);
00292   }
00293       if(debug)
00294   {
00295     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00296     PrintVector(SType1x, M2, "SType1x", matlab);
00297     PrintVector(SType2x, M2, "SType2x", matlab);
00298   }
00299       TotalTestCount++;
00300       SType1ASUMresult = SType1BLAS.ASUM(M, SType1x, incx);
00301       SType2ASUMresult = SType2BLAS.ASUM(M, SType2x, incx);
00302       if(debug)
00303   {
00304     std::cout << "SType1 ASUM result: " << SType1ASUMresult << std::endl;
00305     std::cout << "SType2 ASUM result: " << SType2ASUMresult << std::endl;
00306   }
00307       GoodTestSubcount += CompareScalars(SType1ASUMresult, SType2ASUMresult, TOL);
00308       delete [] SType1x;
00309       delete [] SType2x;
00310     }
00311   GoodTestCount += GoodTestSubcount;
00312   if(verbose || debug) std::cout << "ASUM: " << GoodTestSubcount << " of " << ASUMTESTS << " tests were successful." << std::endl;
00313   if(debug) std::cout << std::endl;
00314 
00315   //--------------------------------------------------------------------------------
00316   // End ASUM Tests
00317   //--------------------------------------------------------------------------------
00318 
00319   //--------------------------------------------------------------------------------
00320   // Begin AXPY Tests
00321   //--------------------------------------------------------------------------------
00322   GoodTestSubcount = 0;
00323   for(i = 0; i < AXPYTESTS; i++)
00324     {
00325       incx = GetRandom(1, MVMAX);
00326       incy = GetRandom(1, MVMAX);
00327       M = GetRandom(MVMIN, MVMAX);
00328       Mx = M*std::abs(incx);
00329       My = M*std::abs(incy);
00330       if (Mx == 0) { Mx = 1; }
00331       if (My == 0) { My = 1; }
00332       SType1x = new SType1[Mx];
00333       SType1y = new SType1[My];
00334       SType2x = new SType2[Mx];
00335       SType2y = new SType2[My]; 
00336       for(j = 0; j < Mx; j++)
00337   {
00338     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00339     SType2x[j] = ConvertType(SType1x[j], convertTo);
00340   }
00341       for(j = 0; j < My; j++)
00342   {
00343     SType1y[j] = GetRandom(-SCALARMAX, SCALARMAX);
00344     SType2y[j] = ConvertType(SType1y[j], convertTo);
00345   }
00346       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
00347       SType2alpha = ConvertType(SType1alpha, convertTo);
00348       if(debug)
00349   {
00350     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00351     std::cout << "SType1alpha = "  << SType1alpha << std::endl;
00352     std::cout << "SType2alpha = " << SType2alpha << std::endl;
00353     PrintVector(SType1x, Mx, "SType1x", matlab);
00354     PrintVector(SType1y, My, "SType1y_before_operation", matlab);
00355     PrintVector(SType2x, Mx, "SType2x", matlab);
00356     PrintVector(SType2y, My, "SType2y_before_operation",  matlab);
00357   }
00358       TotalTestCount++;
00359       SType1BLAS.AXPY(M, SType1alpha, SType1x, incx, SType1y, incy);
00360       SType2BLAS.AXPY(M, SType2alpha, SType2x, incx, SType2y, incy);
00361       if(debug)
00362   {
00363     PrintVector(SType1y, My, "SType1y_after_operation", matlab);
00364     PrintVector(SType2y, My, "SType2y_after_operation", matlab);
00365   }
00366       GoodTestSubcount += CompareVectors(SType1y, SType2y, My, TOL);
00367       delete [] SType1x;
00368       delete [] SType1y;
00369       delete [] SType2x;
00370       delete [] SType2y;
00371     }
00372   GoodTestCount += GoodTestSubcount;
00373   if(verbose || debug) std::cout << "AXPY: " << GoodTestSubcount << " of " << AXPYTESTS << " tests were successful." << std::endl;
00374   if(debug) std::cout << std::endl;
00375   //--------------------------------------------------------------------------------
00376   // End AXPY Tests
00377   //--------------------------------------------------------------------------------
00378 
00379   //--------------------------------------------------------------------------------
00380   // Begin COPY Tests
00381   //--------------------------------------------------------------------------------
00382   GoodTestSubcount = 0;
00383   for(i = 0; i < COPYTESTS; i++)
00384     {
00385       incx = GetRandom(1, MVMAX);
00386       incy = GetRandom(1, MVMAX);
00387       M = GetRandom(MVMIN, MVMAX);
00388       Mx = M*std::abs(incx);
00389       My = M*std::abs(incy);
00390       if (Mx == 0) { Mx = 1; }
00391       if (My == 0) { My = 1; }
00392       SType1x = new SType1[Mx];
00393       SType1y = new SType1[My];
00394       SType2x = new SType2[Mx];
00395       SType2y = new SType2[My]; 
00396       for(j = 0; j < Mx; j++)
00397   {
00398     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00399     SType2x[j] = ConvertType(SType1x[j], convertTo);
00400   }
00401       for(j = 0; j < My; j++)
00402   {
00403     SType1y[j] = GetRandom(-SCALARMAX, SCALARMAX);
00404     SType2y[j] = ConvertType(SType1y[j], convertTo);
00405   }
00406       if(debug)
00407   {
00408     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00409     PrintVector(SType1x, Mx, "SType1x", matlab);
00410     PrintVector(SType1y, My, "SType1y_before_operation", matlab);
00411     PrintVector(SType2x, Mx, "SType2x", matlab);
00412     PrintVector(SType2y, My, "SType2y_before_operation", matlab);
00413   }
00414       TotalTestCount++;
00415       SType1BLAS.COPY(M, SType1x, incx, SType1y, incy);
00416       SType2BLAS.COPY(M, SType2x, incx, SType2y, incy);
00417       if(debug)
00418   {
00419     PrintVector(SType1y, My, "SType1y_after_operation", matlab);
00420     PrintVector(SType2y, My, "SType2y_after_operation", matlab);
00421   }
00422       GoodTestSubcount += CompareVectors(SType1y, SType2y, My, TOL);
00423       delete [] SType1x;
00424       delete [] SType1y;
00425       delete [] SType2x;
00426       delete [] SType2y;
00427     }
00428    GoodTestCount += GoodTestSubcount; if(verbose || debug) std::cout << "COPY: " << GoodTestSubcount << " of " << COPYTESTS << " tests were successful." << std::endl;
00429   if(debug) std::cout << std::endl;
00430   //--------------------------------------------------------------------------------
00431   // End COPY Tests
00432   //--------------------------------------------------------------------------------
00433 
00434   //--------------------------------------------------------------------------------
00435   // Begin DOT Tests
00436   //--------------------------------------------------------------------------------
00437   GoodTestSubcount = 0;
00438   for(i = 0; i < DOTTESTS; i++)
00439     {
00440       incx = GetRandom(1, MVMAX);
00441       incy = GetRandom(1, MVMAX);
00442       M = GetRandom(MVMIN, MVMAX);
00443       Mx = M*std::abs(incx);
00444       My = M*std::abs(incy);
00445       if (Mx == 0) { Mx = 1; }
00446       if (My == 0) { My = 1; }
00447       SType1x = new SType1[Mx];
00448       SType1y = new SType1[My];
00449       SType2x = new SType2[Mx];
00450       SType2y = new SType2[My]; 
00451       for(j = 0; j < Mx; j++)
00452   {
00453     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00454     SType2x[j] = ConvertType(SType1x[j], convertTo);
00455   }
00456       for(j = 0; j < My; j++)
00457   {
00458     SType1y[j] = GetRandom(-SCALARMAX, SCALARMAX);
00459     SType2y[j] = ConvertType(SType1y[j], convertTo);
00460   }
00461       if(debug)
00462   {
00463     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00464     PrintVector(SType1x, Mx, "SType1x", matlab);
00465     PrintVector(SType1y, My, "SType1y", matlab);
00466     PrintVector(SType2x, Mx, "SType2x", matlab);
00467     PrintVector(SType2y, My, "SType2y", matlab);
00468   }
00469       TotalTestCount++;
00470       SType1DOTresult = SType1BLAS.DOT(M, SType1x, incx, SType1y, incy);
00471       SType2DOTresult = SType2BLAS.DOT(M, SType2x, incx, SType2y, incy);
00472       if(debug)
00473   {
00474     std::cout << "SType1 DOT result: " << SType1DOTresult << std::endl;
00475     std::cout << "SType2 DOT result: " << SType2DOTresult << std::endl;
00476   }
00477       GoodTestSubcount += CompareScalars(SType1DOTresult, SType2DOTresult, TOL);
00478       delete [] SType1x;
00479       delete [] SType1y;
00480       delete [] SType2x;
00481       delete [] SType2y;
00482     }
00483   GoodTestCount += GoodTestSubcount;
00484   if(verbose || debug) std::cout << "DOT: " << GoodTestSubcount << " of " << DOTTESTS << " tests were successful." << std::endl;
00485   if(debug) std::cout << std::endl;
00486   //--------------------------------------------------------------------------------
00487   // End DOT Tests
00488   //--------------------------------------------------------------------------------
00489 
00490   //--------------------------------------------------------------------------------
00491   // Begin NRM2 Tests
00492   //--------------------------------------------------------------------------------
00493   GoodTestSubcount = 0;
00494   for(i = 0; i < NRM2TESTS; i++)
00495     {
00496       incx = GetRandom(1, SCALARMAX);
00497       M = GetRandom(MVMIN, MVMAX);
00498       M2 = M*incx; 
00499       SType1x = new SType1[M2];
00500       SType2x = new SType2[M2];
00501       for(j = 0; j < M2; j++)
00502   {
00503     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00504     SType2x[j] = ConvertType(SType1x[j], convertTo);
00505   }
00506       if(debug)
00507   {
00508     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00509     PrintVector(SType1x, M2, "SType1x", matlab);
00510     PrintVector(SType2x, M2, "SType2x", matlab);
00511   }
00512       TotalTestCount++;
00513       SType1NRM2result = SType1BLAS.NRM2(M, SType1x, incx);
00514       SType2NRM2result = SType2BLAS.NRM2(M, SType2x, incx);
00515       if(debug)
00516   {
00517     std::cout << "SType1 NRM2 result: " << SType1NRM2result << std::endl;
00518     std::cout << "SType2 NRM2 result: " << SType2NRM2result << std::endl;
00519   }
00520       GoodTestSubcount += CompareScalars(SType1NRM2result, SType2NRM2result, TOL);
00521       delete [] SType1x;
00522       delete [] SType2x;
00523     }
00524    GoodTestCount += GoodTestSubcount; if(verbose || debug) std::cout << "NRM2: " << GoodTestSubcount << " of " << NRM2TESTS << " tests were successful." << std::endl;
00525   if(debug) std::cout << std::endl;
00526   //--------------------------------------------------------------------------------
00527   // End NRM2 Tests
00528   //--------------------------------------------------------------------------------
00529 
00530   //--------------------------------------------------------------------------------
00531   // Begin SCAL Tests
00532   //--------------------------------------------------------------------------------
00533   GoodTestSubcount = 0;
00534   for(i = 0; i < SCALTESTS; i++)
00535     {
00536       // These will only test for the case that the increment is > 0, the
00537       // templated case can handle when incx < 0, but the blas library doesn't 
00538       // seem to be able to on some machines.
00539       incx = GetRandom(1, SCALARMAX);
00540       M = GetRandom(MVMIN, MVMAX);
00541       M2 = M*incx;
00542       SType1x = new SType1[M2];
00543       SType2x = new SType2[M2];
00544       for(j = 0; j < M2; j++)
00545   {
00546     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00547     SType2x[j] = ConvertType(SType1x[j], convertTo);
00548   }
00549       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
00550       SType2alpha = ConvertType(SType1alpha, convertTo);
00551       if(debug)
00552   {
00553     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00554     std::cout << "SType1alpha = " << SType1alpha << std::endl;
00555     std::cout << "SType2alpha = " << SType2alpha << std::endl;
00556     PrintVector(SType1x, M2, "SType1x_before_operation", matlab);
00557     PrintVector(SType2x, M2, "SType2x_before_operation", matlab);
00558   }
00559       TotalTestCount++;
00560       SType1BLAS.SCAL(M, SType1alpha, SType1x, incx);
00561       SType2BLAS.SCAL(M, SType2alpha, SType2x, incx);
00562       if(debug)
00563   {
00564     PrintVector(SType1x, M2, "SType1x_after_operation", matlab);
00565     PrintVector(SType2x, M2, "SType2x_after_operation", matlab);
00566   }
00567       GoodTestSubcount += CompareVectors(SType1x, SType2x, M2, TOL);
00568       delete [] SType1x;
00569       delete [] SType2x;
00570     }
00571   GoodTestCount += GoodTestSubcount;
00572   if(verbose || debug) std::cout << "SCAL: " << GoodTestSubcount << " of " << SCALTESTS << " tests were successful." << std::endl;
00573   if(debug) std::cout << std::endl;
00574   //--------------------------------------------------------------------------------
00575   // End SCAL Tests
00576   //--------------------------------------------------------------------------------
00577 
00578   //--------------------------------------------------------------------------------
00579   // Begin IAMAX Tests
00580   //--------------------------------------------------------------------------------
00581   GoodTestSubcount = 0;
00582   for(i = 0; i < IAMAXTESTS; i++)
00583     {
00584       incx = GetRandom(1, SCALARMAX);
00585       M = GetRandom(MVMIN, MVMAX);
00586       M2 = M*incx;
00587       SType1x = new SType1[M2];
00588       SType2x = new SType2[M2];
00589       for(j = 0; j < M2; j++)
00590   {
00591     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00592     SType2x[j] = ConvertType(SType1x[j], convertTo);
00593   }
00594       if(debug)
00595   {
00596     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00597     PrintVector(SType1x, M2, "SType1x", matlab);
00598     PrintVector(SType2x, M2, "SType2x", matlab);
00599   }
00600       TotalTestCount++;
00601       SType1IAMAXresult = SType1BLAS.IAMAX(M, SType1x, incx);
00602       SType2IAMAXresult = SType2BLAS.IAMAX(M, SType2x, incx);
00603       if(debug)
00604   {
00605     std::cout << "SType1 IAMAX result: " << SType1IAMAXresult << std::endl;
00606     std::cout << "SType2 IAMAX result: " << SType2IAMAXresult << std::endl;
00607   }
00608       GoodTestSubcount += (SType1IAMAXresult == SType2IAMAXresult);
00609       delete [] SType1x;
00610       delete [] SType2x;
00611     }
00612   GoodTestCount += GoodTestSubcount;
00613   if(verbose || debug) std::cout << "IAMAX: " << GoodTestSubcount << " of " << IAMAXTESTS << " tests were successful." << std::endl;
00614   if(debug) std::cout << std::endl;
00615   //--------------------------------------------------------------------------------
00616   // End IAMAX Tests
00617   //--------------------------------------------------------------------------------
00618 
00619   //--------------------------------------------------------------------------------
00620   // BEGIN LEVEL 2 BLAS TESTS
00621   //--------------------------------------------------------------------------------
00622   // Begin GEMV Tests
00623   //--------------------------------------------------------------------------------
00624   GoodTestSubcount = 0;
00625   for(i = 0; i < GEMVTESTS; i++)
00626     {
00627       // The parameters used to construct the test problem are chosen to be
00628       // valid parameters, so the GEMV routine won't bomb out.
00629       incx = GetRandom(1, MVMAX);
00630       while (incx == 0) {
00631           incx = GetRandom(1, MVMAX);
00632       }
00633       incy = GetRandom(1, MVMAX);
00634       while (incy == 0) {
00635           incy = GetRandom(1, MVMAX);
00636       }   
00637       M = GetRandom(MVMIN, MVMAX);
00638       N = GetRandom(MVMIN, MVMAX);
00639 
00640       TRANS = RandomTRANS();
00641       if (Teuchos::ETranspChar[TRANS] == 'N') { 
00642         M2 = M*std::abs(incy);
00643         N2 = N*std::abs(incx);   
00644       } else {
00645   M2 = N*std::abs(incy);
00646   N2 = M*std::abs(incx);
00647       }
00648 
00649       LDA = GetRandom(MVMIN, MVMAX);
00650       while (LDA < M) {
00651           LDA = GetRandom(MVMIN, MVMAX);
00652       }   
00653 
00654       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
00655       SType1beta = GetRandom(-SCALARMAX, SCALARMAX);
00656       SType2alpha = ConvertType(SType1alpha, convertTo);
00657       SType2beta = ConvertType(SType1beta, convertTo);
00658 
00659       SType1A = new SType1[LDA * N];
00660       SType1x = new SType1[N2];
00661       SType1y = new SType1[M2];
00662       SType2A = new SType2[LDA * N];
00663       SType2x = new SType2[N2];
00664       SType2y = new SType2[M2]; 
00665 
00666       for(j = 0; j < LDA * N; j++)
00667   {
00668     SType1A[j] = GetRandom(-SCALARMAX, SCALARMAX);
00669     SType2A[j] = ConvertType(SType1A[j], convertTo);
00670   }
00671       for(j = 0; j < N2; j++)
00672   {
00673     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00674     SType2x[j] = ConvertType(SType1x[j], convertTo);
00675   }
00676       for(j = 0; j < M2; j++)
00677   {
00678     SType1y[j] = GetRandom(-SCALARMAX, SCALARMAX);
00679     SType2y[j] = ConvertType(SType1y[j], convertTo);
00680   }
00681       if(debug)
00682   {
00683     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00684     std::cout << "SType1alpha = " << SType1alpha << std::endl;
00685     std::cout << "SType2alpha = " << SType2alpha << std::endl;
00686     std::cout << "SType1beta = " << SType1beta << std::endl;
00687     std::cout << "SType2beta = " << SType2beta << std::endl;
00688     PrintMatrix(SType1A, M, N, LDA, "SType1A", matlab);
00689     PrintVector(SType1x, N2, "SType1x", matlab);
00690     PrintVector(SType1y, M2, "SType1y_before_operation", matlab);
00691     PrintMatrix(SType2A, M, N, LDA, "SType2A", matlab);
00692     PrintVector(SType2x, N2, "SType2x", matlab);
00693     PrintVector(SType2y, M2, "SType2y_before_operation", matlab);
00694   }
00695       TotalTestCount++;
00696       SType1BLAS.GEMV(TRANS, M, N, SType1alpha, SType1A, LDA, SType1x, incx, SType1beta, SType1y, incy);
00697       SType2BLAS.GEMV(TRANS, M, N, SType2alpha, SType2A, LDA, SType2x, incx, SType2beta, SType2y, incy);
00698       if(debug)
00699   {
00700     PrintVector(SType1y, M2, "SType1y_after_operation", matlab);
00701     PrintVector(SType2y, M2, "SType2y_after_operation", matlab);
00702   }
00703       GoodTestSubcount += CompareVectors(SType1y, SType2y, M2, TOL);
00704       delete [] SType1A;
00705       delete [] SType1x;
00706       delete [] SType1y;
00707       delete [] SType2A;
00708       delete [] SType2x;
00709       delete [] SType2y;
00710     }
00711   GoodTestCount += GoodTestSubcount;
00712   if(verbose || debug) std::cout << "GEMV: " << GoodTestSubcount << " of " << GEMVTESTS << " tests were successful." << std::endl;
00713   if(debug) std::cout << std::endl;
00714   //--------------------------------------------------------------------------------
00715   // End GEMV Tests
00716   //--------------------------------------------------------------------------------
00717 
00718   //--------------------------------------------------------------------------------
00719   // Begin TRMV Tests
00720   //--------------------------------------------------------------------------------
00721   GoodTestSubcount = 0;
00722   for(i = 0; i < TRMVTESTS; i++)
00723     {
00724       UPLO = RandomUPLO();
00725       TRANSA = RandomTRANS();
00726       // Since the entries are integers, we don't want to use the unit diagonal feature,
00727       // this creates ill-conditioned, nearly-singular matrices.
00728       //DIAG = RandomDIAG();  
00729       DIAG = Teuchos::NON_UNIT_DIAG;
00730 
00731       N = GetRandom(MVMIN, MVMAX);
00732       incx = GetRandom(1, MVMAX);
00733       while (incx == 0) {
00734           incx = GetRandom(1, MVMAX);
00735       }
00736       N2 = N*std::abs(incx);
00737       SType1x = new SType1[N2];
00738       SType2x = new SType2[N2];
00739 
00740       for(j = 0; j < N2; j++)
00741   {
00742     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00743     SType2x[j] = ConvertType(SType1x[j], convertTo);
00744   }
00745 
00746       LDA = GetRandom(MVMIN, MVMAX);
00747       while (LDA < N) {
00748   LDA = GetRandom(MVMIN, MVMAX);
00749       }
00750       SType1A = new SType1[LDA * N];
00751       SType2A = new SType2[LDA * N];
00752 
00753       for(j = 0; j < N; j++)
00754   {      
00755     if(Teuchos::EUploChar[UPLO] == 'U') {
00756       // The operator is upper triangular, make sure that the entries are
00757       // only in the upper triangular part of A and the diagonal is non-zero.
00758       for(k = 0; k < N; k++) 
00759       {
00760         if(k < j) {
00761     SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
00762         } else {
00763     SType1A[j*LDA+k] = SType1zero;
00764         }
00765         SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
00766         if(k == j) {
00767     if (Teuchos::EDiagChar[DIAG] == 'N') {
00768       SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
00769       while (SType1A[j*LDA+k] == SType1zero) {
00770         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
00771       }
00772       SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
00773     } else {
00774       SType1A[j*LDA+k] = SType1one;
00775       SType2A[j*LDA+k] = SType2one;
00776     }
00777         }     
00778       }
00779     } else {
00780       // The operator is lower triangular, make sure that the entries are
00781       // only in the lower triangular part of A and the diagonal is non-zero.
00782       for(k = 0; k < N; k++) 
00783         {
00784     if(k > j) {
00785       SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
00786     } else {
00787       SType1A[j*LDA+k] = SType1zero;
00788     }
00789     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
00790     if(k == j) {
00791       if (Teuchos::EDiagChar[DIAG] == 'N') {
00792         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
00793         while (SType1A[j*LDA+k] == SType1zero) {
00794           SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
00795         }
00796         SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
00797       } else {
00798         SType1A[j*LDA+k] = SType1one;
00799         SType2A[j*LDA+k] = SType2one;
00800       }
00801     }     
00802         } // end for(k=0 ...    
00803     } // end if(UPLO == 'U') ...
00804   } // end for(j=0 ...      for(j = 0; j < N*N; j++)
00805       
00806       if(debug)
00807   {
00808     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00809     PrintMatrix(SType1A, N, N, LDA,"SType1A", matlab);
00810     PrintVector(SType1x, N2, "SType1x_before_operation", matlab);
00811     PrintMatrix(SType2A, N, N, LDA, "SType2A", matlab);
00812     PrintVector(SType2x, N2, "SType2x_before_operation", matlab);
00813   }
00814       TotalTestCount++;
00815       SType1BLAS.TRMV(UPLO, TRANSA, DIAG, N, SType1A, LDA, SType1x, incx);
00816       SType2BLAS.TRMV(UPLO, TRANSA, DIAG, N, SType2A, LDA, SType2x, incx);
00817       if(debug)
00818   {
00819     PrintVector(SType1x, N2, "SType1x_after_operation", matlab);
00820     PrintVector(SType2x, N2, "SType2x_after_operation", matlab);
00821   }
00822       GoodTestSubcount += CompareVectors(SType1x, SType2x, N2, TOL);
00823       delete [] SType1A;
00824       delete [] SType1x;
00825       delete [] SType2A;
00826       delete [] SType2x;
00827     }
00828   GoodTestCount += GoodTestSubcount;
00829   if(verbose || debug) std::cout << "TRMV: " << GoodTestSubcount << " of " << TRMVTESTS << " tests were successful." << std::endl;
00830   if(debug) std::cout << std::endl;
00831   //--------------------------------------------------------------------------------
00832   // End TRMV Tests
00833   //--------------------------------------------------------------------------------
00834 
00835   //--------------------------------------------------------------------------------
00836   // Begin GER Tests
00837   //--------------------------------------------------------------------------------
00838   GoodTestSubcount = 0;
00839   for(i = 0; i < GERTESTS; i++)
00840     {
00841       incx = GetRandom(1, MVMAX);
00842       while (incx == 0) {
00843           incx = GetRandom(1, MVMAX);
00844       }   
00845       incy = GetRandom(1, MVMAX);
00846       while (incy == 0) {
00847           incy = GetRandom(1, MVMAX);
00848       }   
00849       M = GetRandom(MVMIN, MVMAX);
00850       N = GetRandom(MVMIN, MVMAX);
00851 
00852       M2 = M*std::abs(incx);
00853       N2 = N*std::abs(incy);   
00854 
00855       LDA = GetRandom(MVMIN, MVMAX);
00856       while (LDA < M) {
00857           LDA = GetRandom(MVMIN, MVMAX);
00858       }   
00859 
00860       SType1A = new SType1[LDA * N];
00861       SType1x = new SType1[M2];
00862       SType1y = new SType1[N2];
00863       SType2A = new SType2[LDA * N];
00864       SType2x = new SType2[M2];
00865       SType2y = new SType2[N2];
00866       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
00867       SType2alpha = ConvertType(SType1alpha, convertTo);
00868       for(j = 0; j < LDA * N; j++)
00869   {
00870     SType1A[j] = GetRandom(-SCALARMAX, SCALARMAX);
00871     SType2A[j] = ConvertType(SType1A[j], convertTo);
00872   }
00873       for(j = 0; j < M2; j++)
00874   {
00875     SType1x[j] = GetRandom(-SCALARMAX, SCALARMAX);
00876     SType2x[j] = ConvertType(SType1x[j], convertTo);
00877   }
00878       for(j = 0; j < N2; j++)
00879   {
00880     SType1y[j] = GetRandom(-SCALARMAX, SCALARMAX);
00881     SType2y[j] = ConvertType(SType1y[j], convertTo);
00882   }
00883       if(debug)
00884   {
00885     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00886     std::cout << "SType1alpha = " << SType1alpha << std::endl;
00887     std::cout << "SType2alpha = " << SType2alpha << std::endl;
00888     PrintMatrix(SType1A, M, N, LDA,"SType1A_before_operation", matlab);
00889     PrintVector(SType1x, M2, "SType1x", matlab);
00890     PrintVector(SType1y, N2, "SType1y", matlab);
00891     PrintMatrix(SType2A, M, N, LDA,"SType2A_before_operation", matlab);
00892     PrintVector(SType2x, M2, "SType2x", matlab);
00893     PrintVector(SType2y, N2, "SType2y", matlab);
00894   }
00895       TotalTestCount++;
00896       SType1BLAS.GER(M, N, SType1alpha, SType1x, incx, SType1y, incy, SType1A, LDA);
00897       SType2BLAS.GER(M, N, SType2alpha, SType2x, incx, SType2y, incy, SType2A, LDA);
00898       if(debug)
00899   {
00900     PrintMatrix(SType1A, M, N, LDA, "SType1A_after_operation", matlab);
00901     PrintMatrix(SType2A, M, N, LDA, "SType2A_after_operation", matlab);
00902   }
00903       GoodTestSubcount += CompareMatrices(SType1A, SType2A, M, N, LDA, TOL);
00904       delete [] SType1A;
00905       delete [] SType1x;
00906       delete [] SType1y;
00907       delete [] SType2A;
00908       delete [] SType2x;
00909       delete [] SType2y;
00910     }
00911   GoodTestCount += GoodTestSubcount;
00912   if(verbose || debug) std::cout << "GER: " << GoodTestSubcount << " of " << GERTESTS << " tests were successful." << std::endl;
00913   if(debug) std::cout << std::endl;
00914   //--------------------------------------------------------------------------------
00915   // End GER Tests
00916   //--------------------------------------------------------------------------------
00917 
00918   //--------------------------------------------------------------------------------
00919   // BEGIN LEVEL 3 BLAS TESTS
00920   //--------------------------------------------------------------------------------
00921   // Begin GEMM Tests
00922   //--------------------------------------------------------------------------------
00923   GoodTestSubcount = 0;
00924   for(i = 0; i < GEMMTESTS; i++)
00925     { 
00926       TRANSA = RandomTRANS();      
00927       TRANSB = RandomTRANS();
00928       M = GetRandom(MVMIN, MVMAX);
00929       N = GetRandom(MVMIN, MVMAX);
00930       P = GetRandom(MVMIN, MVMAX);
00931 
00932       if(debug) {
00933     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
00934       }
00935       LDA = GetRandom(MVMIN, MVMAX);
00936       if (Teuchos::ETranspChar[TRANSA] == 'N') {
00937   while (LDA < M) {  LDA = GetRandom(MVMIN, MVMAX); }
00938   SType1A = new SType1[LDA * P];
00939   SType2A = new SType2[LDA * P];
00940   for(j = 0; j < LDA * P; j++)
00941   {
00942     SType1A[j] = GetRandom(-SCALARMAX, SCALARMAX);
00943     SType2A[j] = ConvertType(SType1A[j], convertTo);
00944   }
00945   if (debug) {
00946   PrintMatrix(SType1A, M, P, LDA, "SType1A", matlab);
00947   PrintMatrix(SType2A, M, P, LDA, "SType2A", matlab);
00948   }
00949       } else {
00950   while (LDA < P) {  LDA = GetRandom(MVMIN, MVMAX); }
00951   SType1A = new SType1[LDA * M];
00952   SType2A = new SType2[LDA * M];
00953   for(j = 0; j < LDA * M; j++)
00954   {
00955     SType1A[j] = GetRandom(-SCALARMAX, SCALARMAX);
00956     SType2A[j] = ConvertType(SType1A[j], convertTo);
00957   }
00958   if (debug) {
00959   PrintMatrix(SType1A, P, M, LDA, "SType1A", matlab);
00960   PrintMatrix(SType2A, P, M, LDA, "SType2A", matlab);
00961   }
00962       }
00963 
00964       LDB = GetRandom(MVMIN, MVMAX);
00965       if (Teuchos::ETranspChar[TRANSB] == 'N') {
00966   while (LDB < P) {  LDB = GetRandom(MVMIN, MVMAX); }
00967   SType1B = new SType1[LDB * N];
00968   SType2B = new SType2[LDB * N];
00969   for(j = 0; j < LDB * N; j++)
00970   {
00971     SType1B[j] = GetRandom(-SCALARMAX, SCALARMAX);
00972     SType2B[j] = ConvertType(SType1B[j], convertTo);
00973   }
00974   if (debug) {
00975     PrintMatrix(SType1B, P, N, LDB,"SType1B", matlab);
00976     PrintMatrix(SType2B, P, N, LDB,"SType2B", matlab);
00977   }
00978       } else { 
00979   while (LDB < N) {  LDB = GetRandom(MVMIN, MVMAX); }
00980   SType1B = new SType1[LDB * P];
00981   SType2B = new SType2[LDB * P];
00982   for(j = 0; j < LDB * P; j++)
00983   {
00984     SType1B[j] = GetRandom(-SCALARMAX, SCALARMAX);
00985     SType2B[j] = ConvertType(SType1B[j], convertTo);
00986   } 
00987   if (debug) {
00988     PrintMatrix(SType1B, N, P, LDB,"SType1B", matlab);
00989     PrintMatrix(SType2B, N, P, LDB,"SType2B", matlab);
00990   }
00991       }
00992 
00993       LDC = GetRandom(MVMIN, MVMAX);
00994       while (LDC < M) {  LDC = GetRandom(MVMIN, MVMAX); }
00995       SType1C = new SType1[LDC * N];
00996       SType2C = new SType2[LDC * N];
00997       for(j = 0; j < LDC * N; j++) {
00998     SType1C[j] = GetRandom(-SCALARMAX, SCALARMAX);
00999     SType2C[j] = ConvertType(SType1C[j], convertTo);
01000       }
01001       if(debug)
01002   {
01003     PrintMatrix(SType1C, M, N, LDC, "SType1C_before_operation", matlab);
01004     PrintMatrix(SType2C, M, N, LDC, "SType2C_before_operation", matlab);
01005   }
01006   
01007       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
01008       SType1beta = GetRandom(-SCALARMAX, SCALARMAX);
01009       SType2alpha = ConvertType(SType1alpha, convertTo);
01010       SType2beta = ConvertType(SType1beta, convertTo);
01011 
01012       TotalTestCount++;
01013       SType1BLAS.GEMM(TRANSA, TRANSB, M, N, P, SType1alpha, SType1A, LDA, SType1B, LDB, SType1beta, SType1C, LDC);
01014       SType2BLAS.GEMM(TRANSA, TRANSB, M, N, P, SType2alpha, SType2A, LDA, SType2B, LDB, SType2beta, SType2C, LDC);
01015       if(debug)
01016   {
01017     PrintMatrix(SType1C, M, N, LDC, "SType1C_after_operation", matlab);
01018     PrintMatrix(SType2C, M, N, LDC, "SType2C_after_operation", matlab);
01019   }
01020       GoodTestSubcount += CompareMatrices(SType1C, SType2C, M, N, LDC, TOL);
01021       delete [] SType1A;
01022       delete [] SType1B;
01023       delete [] SType1C;
01024       delete [] SType2A;
01025       delete [] SType2B;
01026       delete [] SType2C;
01027     }
01028   GoodTestCount += GoodTestSubcount;
01029   if(verbose || debug) std::cout << "GEMM: " << GoodTestSubcount << " of " << GEMMTESTS << " tests were successful." << std::endl;
01030   if(debug) std::cout << std::endl;
01031   //--------------------------------------------------------------------------------
01032   // End GEMM Tests
01033   //--------------------------------------------------------------------------------
01034 
01035   //--------------------------------------------------------------------------------
01036   // Begin SYMM Tests
01037   //--------------------------------------------------------------------------------
01038   GoodTestSubcount = 0;
01039   for(i = 0; i < SYMMTESTS; i++)
01040     { 
01041       M = GetRandom(MVMIN, MVMAX);
01042       N = GetRandom(MVMIN, MVMAX);
01043       SIDE = RandomSIDE();
01044       UPLO = RandomUPLO();
01045 
01046       LDA = GetRandom(MVMIN, MVMAX);
01047       if(Teuchos::ESideChar[SIDE] == 'L') {
01048   while (LDA < M) { LDA = GetRandom(MVMIN, MVMAX); }
01049   SType1A = new SType1[LDA * M];
01050   SType2A = new SType2[LDA * M];
01051   for(j = 0; j < LDA * M; j++) {
01052     SType1A[j] = GetRandom(-SCALARMAX, SCALARMAX);
01053     SType2A[j] = ConvertType(SType1A[j], convertTo);
01054   }
01055       } else {
01056   while (LDA < N) { LDA = GetRandom(MVMIN, MVMAX); }
01057   SType1A = new SType1[LDA * N];
01058   SType2A = new SType2[LDA * N];
01059   for(j = 0; j < LDA * N; j++) {
01060     SType1A[j] = GetRandom(-SCALARMAX, SCALARMAX);
01061     SType2A[j] = ConvertType(SType1A[j], convertTo);
01062   }
01063       }
01064 
01065       LDB = GetRandom(MVMIN, MVMAX);
01066       while (LDB < M) {  LDB = GetRandom(MVMIN, MVMAX); }
01067       SType1B = new SType1[LDB * N];
01068       SType2B = new SType2[LDB * N];
01069       for(j = 0; j < LDB * N; j++) {
01070     SType1B[j] = GetRandom(-SCALARMAX, SCALARMAX);
01071     SType2B[j] = ConvertType(SType1B[j], convertTo);
01072       }
01073     
01074       LDC = GetRandom(MVMIN, MVMAX);
01075       while (LDC < M) {  LDC = GetRandom(MVMIN, MVMAX); }
01076       SType1C = new SType1[LDC * N];
01077       SType2C = new SType2[LDC * N];
01078       for(j = 0; j < LDC * N; j++) {
01079     SType1C[j] = GetRandom(-SCALARMAX, SCALARMAX);
01080     SType2C[j] = ConvertType(SType1C[j], convertTo);
01081       }
01082       
01083       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
01084       SType1beta = GetRandom(-SCALARMAX, SCALARMAX);
01085       SType2alpha = ConvertType(SType1alpha, convertTo);
01086       SType2beta = ConvertType(SType1beta, convertTo);
01087       if(debug)
01088   {
01089     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
01090     std::cout << "SType1alpha = " << SType1alpha << std::endl;
01091     std::cout << "SType2alpha = " << SType2alpha << std::endl;
01092     std::cout << "SType1beta = " << SType1beta << std::endl;
01093     std::cout << "SType2beta = " << SType2beta << std::endl;
01094     if (Teuchos::ESideChar[SIDE] == 'L') {
01095         PrintMatrix(SType1A, M, M, LDA,"SType1A", matlab);
01096         PrintMatrix(SType2A, M, M, LDA,"SType2A", matlab);
01097     } else {
01098       PrintMatrix(SType1A, N, N, LDA, "SType1A", matlab);
01099       PrintMatrix(SType2A, N, N, LDA, "SType2A", matlab);
01100     }
01101     PrintMatrix(SType1B, M, N, LDB,"SType1B", matlab);
01102     PrintMatrix(SType1C, M, N, LDC,"SType1C_before_operation", matlab);
01103     PrintMatrix(SType2B, M, N, LDB,"SType2B", matlab);
01104     PrintMatrix(SType2C, M, N, LDC,"SType2C_before_operation", matlab);
01105   }
01106       TotalTestCount++;
01107 
01108       SType1BLAS.SYMM(SIDE, UPLO, M, N, SType1alpha, SType1A, LDA, SType1B, LDB, SType1beta, SType1C, LDC);
01109       SType2BLAS.SYMM(SIDE, UPLO, M, N, SType2alpha, SType2A, LDA, SType2B, LDB, SType2beta, SType2C, LDC);
01110       if(debug)
01111   {
01112     PrintMatrix(SType1C, M, N, LDC,"SType1C_after_operation", matlab);
01113     PrintMatrix(SType2C, M, N, LDC,"SType2C_after_operation", matlab);
01114   }
01115       GoodTestSubcount += CompareMatrices(SType1C, SType2C, M, N, LDC, TOL);
01116 
01117       delete [] SType1A;
01118       delete [] SType1B;
01119       delete [] SType1C;
01120       delete [] SType2A;
01121       delete [] SType2B;
01122       delete [] SType2C;
01123     }
01124   GoodTestCount += GoodTestSubcount;
01125   if(verbose || debug) std::cout << "SYMM: " << GoodTestSubcount << " of " << SYMMTESTS << " tests were successful." << std::endl;
01126   if(debug) std::cout << std::endl;
01127   //--------------------------------------------------------------------------------
01128   // End SYMM Tests
01129   //--------------------------------------------------------------------------------
01130 
01131   //--------------------------------------------------------------------------------
01132   // Begin TRMM Tests
01133   //--------------------------------------------------------------------------------
01134   GoodTestSubcount = 0;
01135   for(i = 0; i < TRMMTESTS; i++)
01136     { 
01137       M = GetRandom(MVMIN, MVMAX);
01138       N = GetRandom(MVMIN, MVMAX);
01139 
01140       LDB = GetRandom(MVMIN, MVMAX);
01141       while (LDB < M) {
01142     LDB = GetRandom(MVMIN, MVMAX);
01143       }
01144 
01145       SType1B = new SType1[LDB * N];
01146       SType2B = new SType2[LDB * N];
01147 
01148       SIDE = RandomSIDE();
01149       UPLO = RandomUPLO();
01150       TRANSA = RandomTRANS();
01151       DIAG = RandomDIAG();
01152 
01153       if(Teuchos::ESideChar[SIDE] == 'L')  // The operator is on the left side
01154   {
01155           LDA = GetRandom(MVMIN, MVMAX);
01156           while (LDA < M) {
01157         LDA = GetRandom(MVMIN, MVMAX);
01158           }
01159 
01160     SType1A = new SType1[LDA * M];
01161     SType2A = new SType2[LDA * M];
01162 
01163     for(j = 0; j < M; j++)
01164       {      
01165         if(Teuchos::EUploChar[UPLO] == 'U') {
01166     // The operator is upper triangular, make sure that the entries are
01167     // only in the upper triangular part of A and the diagonal is non-zero.
01168     for(k = 0; k < M; k++) 
01169     {
01170         if(k < j) {
01171             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01172         } else {
01173       SType1A[j*LDA+k] = SType1zero;
01174         }
01175               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01176         if(k == j) {
01177       if (Teuchos::EDiagChar[DIAG] == 'N') {
01178                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01179           while (SType1A[j*LDA+k] == SType1zero) {
01180         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01181           }
01182                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);  
01183           } else {
01184                 SType1A[j*LDA+k] = SType1one;
01185                     SType2A[j*LDA+k] = SType2one;
01186           }     
01187         }
01188     }
01189         } else {
01190     // The operator is lower triangular, make sure that the entries are
01191     // only in the lower triangular part of A and the diagonal is non-zero.
01192     for(k = 0; k < M; k++) 
01193     {
01194         if(k > j) {
01195             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01196         } else {
01197       SType1A[j*LDA+k] = SType1zero;
01198         }
01199               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01200         if(k == j) {
01201       if (Teuchos::EDiagChar[DIAG] == 'N') {
01202                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01203           while (SType1A[j*LDA+k] == SType1zero) {
01204         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01205           }
01206                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01207           } else {
01208                 SType1A[j*LDA+k] = SType1one;
01209                     SType2A[j*LDA+k] = SType2one;
01210       }
01211         }     
01212     } // end for(k=0 ...    
01213         } // end if(UPLO == 'U') ...
01214       } // end for(j=0 ...
01215   } // if(SIDE == 'L') ...
01216       else // The operator is on the right side
01217   {
01218           LDA = GetRandom(MVMIN, MVMAX);
01219           while (LDA < N) {
01220         LDA = GetRandom(MVMIN, MVMAX);
01221           }
01222 
01223     SType1A = new SType1[LDA * N];
01224     SType2A = new SType2[LDA * N];
01225 
01226     for(j = 0; j < N; j++)
01227       {      
01228         if(Teuchos::EUploChar[UPLO] == 'U') {
01229     // The operator is upper triangular, make sure that the entries are
01230     // only in the upper triangular part of A and the diagonal is non-zero.
01231     for(k = 0; k < N; k++) 
01232     {
01233         if(k < j) {
01234             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01235         } else {
01236       SType1A[j*LDA+k] = SType1zero;
01237         }
01238               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01239         if(k == j) {
01240       if (Teuchos::EDiagChar[DIAG] == 'N') {
01241                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01242           while (SType1A[j*LDA+k] == SType1zero) {
01243         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01244           }
01245                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01246           } else {
01247                 SType1A[j*LDA+k] = SType1one;
01248                     SType2A[j*LDA+k] = SType2one;
01249       }
01250         }     
01251     }
01252         } else {
01253     // The operator is lower triangular, make sure that the entries are
01254     // only in the lower triangular part of A and the diagonal is non-zero.
01255     for(k = 0; k < N; k++) 
01256     {
01257         if(k > j) {
01258             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01259         } else {
01260       SType1A[j*LDA+k] = SType1zero;
01261         }
01262               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01263         if(k == j) {
01264       if (Teuchos::EDiagChar[DIAG] == 'N') {
01265                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01266           while (SType1A[j*LDA+k] == SType1zero) {
01267         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01268           }
01269                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01270           } else {
01271                 SType1A[j*LDA+k] = SType1one;
01272                     SType2A[j*LDA+k] = SType2one;
01273       }
01274         }     
01275     } // end for(k=0 ...    
01276         } // end if(UPLO == 'U') ...
01277       } // end for(j=0 ...
01278   } // end if(SIDE == 'L') ...
01279 
01280       // Fill in the right hand side block B.
01281       for(j = 0; j < N; j++) {
01282     for(k = 0; k < M; k++) {
01283       SType1B[j*LDB+k] = GetRandom(-SCALARMAX, SCALARMAX);
01284       SType2B[j*LDB+k] = ConvertType(SType1B[j*LDB+k], convertTo);
01285     }
01286       }
01287       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
01288       SType2alpha = ConvertType(SType1alpha, convertTo);
01289       if(debug)
01290   {
01291     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
01292     std::cout << "SType1alpha = " << SType1alpha << std::endl;
01293     std::cout << "SType2alpha = " << SType2alpha << std::endl;
01294           if(Teuchos::ESideChar[SIDE] == 'L') { 
01295       PrintMatrix(SType1A, M, M, LDA, "SType1A", matlab);
01296       PrintMatrix(SType2A, M, M, LDA, "SType2A", matlab);
01297     } else {
01298       PrintMatrix(SType1A, N, N, LDA, "SType1A", matlab);
01299       PrintMatrix(SType2A, N, N, LDA, "SType2A", matlab);
01300     }
01301     PrintMatrix(SType1B, M, N, LDB,"SType1B_before_operation", matlab);
01302     PrintMatrix(SType2B, M, N, LDB,"SType2B_before_operation", matlab);
01303   }
01304       TotalTestCount++;
01305       SType1BLAS.TRMM(SIDE, UPLO, TRANSA, DIAG, M, N, SType1alpha, SType1A, LDA, SType1B, LDB);
01306       SType2BLAS.TRMM(SIDE, UPLO, TRANSA, DIAG, M, N, SType2alpha, SType2A, LDA, SType2B, LDB);
01307       if(debug)
01308   {
01309     PrintMatrix(SType1B, M, N, LDB, "SType1B_after_operation", matlab);
01310     PrintMatrix(SType2B, M, N, LDB, "SType2B_after_operation", matlab);
01311   }
01312       GoodTestSubcount += CompareMatrices(SType1B, SType2B, M, N, LDB, TOL);
01313       delete [] SType1A;
01314       delete [] SType1B;
01315       delete [] SType2A;
01316       delete [] SType2B;
01317     }
01318   GoodTestCount += GoodTestSubcount;
01319   if(verbose || debug) std::cout << "TRMM: " << GoodTestSubcount << " of " << TRMMTESTS << " tests were successful." << std::endl;
01320   if(debug) std::cout << std::endl;
01321   //--------------------------------------------------------------------------------
01322   // End TRMM Tests
01323   //--------------------------------------------------------------------------------
01324 
01325   //--------------------------------------------------------------------------------
01326   // Begin TRSM Tests
01327   //--------------------------------------------------------------------------------
01328   GoodTestSubcount = 0;
01329   for(i = 0; i < TRSMTESTS; i++)
01330     { 
01331       M = GetRandom(MVMIN, MVMAX);
01332       N = GetRandom(MVMIN, MVMAX);
01333 
01334       LDB = GetRandom(MVMIN, MVMAX);
01335       while (LDB < M) {
01336     LDB = GetRandom(MVMIN, MVMAX);
01337       }
01338 
01339       SType1B = new SType1[LDB * N];
01340       SType2B = new SType2[LDB * N];
01341 
01342       SIDE = RandomSIDE();
01343       UPLO = RandomUPLO();
01344       TRANSA = RandomTRANS();
01345       // Since the entries are integers, we don't want to use the unit diagonal feature,
01346       // this creates ill-conditioned, nearly-singular matrices.
01347       //DIAG = RandomDIAG();  
01348       DIAG = Teuchos::NON_UNIT_DIAG;
01349 
01350       if(Teuchos::ESideChar[SIDE] == 'L')  // The operator is on the left side
01351   {
01352           LDA = GetRandom(MVMIN, MVMAX);
01353           while (LDA < M) {
01354         LDA = GetRandom(MVMIN, MVMAX);
01355           }
01356 
01357     SType1A = new SType1[LDA * M];
01358     SType2A = new SType2[LDA * M];
01359 
01360     for(j = 0; j < M; j++)
01361       {      
01362         if(Teuchos::EUploChar[UPLO] == 'U') {
01363     // The operator is upper triangular, make sure that the entries are
01364     // only in the upper triangular part of A and the diagonal is non-zero.
01365     for(k = 0; k < M; k++) 
01366     {
01367         if(k < j) {
01368             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01369         } else {
01370       SType1A[j*LDA+k] = SType1zero;
01371         }
01372               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01373         if(k == j) {
01374       if (Teuchos::EDiagChar[DIAG] == 'N') {
01375                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01376           while (SType1A[j*LDA+k] == SType1zero) {
01377         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01378           }
01379                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);  
01380           } else {
01381                 SType1A[j*LDA+k] = SType1one;
01382                     SType2A[j*LDA+k] = SType2one;
01383           }     
01384         }
01385     }
01386         } else {
01387     // The operator is lower triangular, make sure that the entries are
01388     // only in the lower triangular part of A and the diagonal is non-zero.
01389     for(k = 0; k < M; k++) 
01390     {
01391         if(k > j) {
01392             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01393         } else {
01394       SType1A[j*LDA+k] = SType1zero;
01395         }
01396               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01397         if(k == j) {
01398       if (Teuchos::EDiagChar[DIAG] == 'N') {
01399                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01400           while (SType1A[j*LDA+k] == SType1zero) {
01401         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01402           }
01403                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01404           } else {
01405                 SType1A[j*LDA+k] = SType1one;
01406                     SType2A[j*LDA+k] = SType2one;
01407       }
01408         }     
01409     } // end for(k=0 ...    
01410         } // end if(UPLO == 'U') ...
01411       } // end for(j=0 ...
01412   } // if(SIDE == 'L') ...
01413       else // The operator is on the right side
01414   {
01415           LDA = GetRandom(MVMIN, MVMAX);
01416           while (LDA < N) {
01417         LDA = GetRandom(MVMIN, MVMAX);
01418           }
01419 
01420     SType1A = new SType1[LDA * N];
01421     SType2A = new SType2[LDA * N];
01422 
01423     for(j = 0; j < N; j++)
01424       {      
01425         if(Teuchos::EUploChar[UPLO] == 'U') {
01426     // The operator is upper triangular, make sure that the entries are
01427     // only in the upper triangular part of A and the diagonal is non-zero.
01428     for(k = 0; k < N; k++) 
01429     {
01430         if(k < j) {
01431             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01432         } else {
01433       SType1A[j*LDA+k] = SType1zero;
01434         }
01435               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01436         if(k == j) {
01437       if (Teuchos::EDiagChar[DIAG] == 'N') {
01438                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01439           while (SType1A[j*LDA+k] == SType1zero) {
01440         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01441           }
01442                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01443           } else {
01444                 SType1A[j*LDA+k] = SType1one;
01445                     SType2A[j*LDA+k] = SType2one;
01446       }
01447         }     
01448     }
01449         } else {
01450     // The operator is lower triangular, make sure that the entries are
01451     // only in the lower triangular part of A and the diagonal is non-zero.
01452     for(k = 0; k < N; k++) 
01453     {
01454         if(k > j) {
01455             SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01456         } else {
01457       SType1A[j*LDA+k] = SType1zero;
01458         }
01459               SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01460         if(k == j) {
01461       if (Teuchos::EDiagChar[DIAG] == 'N') {
01462                 SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01463           while (SType1A[j*LDA+k] == SType1zero) {
01464         SType1A[j*LDA+k] = GetRandom(-SCALARMAX, SCALARMAX);
01465           }
01466                     SType2A[j*LDA+k] = ConvertType(SType1A[j*LDA+k], convertTo);
01467           } else {
01468                 SType1A[j*LDA+k] = SType1one;
01469                     SType2A[j*LDA+k] = SType2one;
01470       }
01471         }     
01472     } // end for(k=0 ...    
01473         } // end if(UPLO == 'U') ...
01474       } // end for(j=0 ...
01475   } // end if(SIDE == 'L') ...
01476 
01477       // Fill in the right hand side block B.
01478       for(j = 0; j < N; j++)
01479   {
01480     for(k = 0; k < M; k++) 
01481       {
01482       SType1B[j*LDB+k] = GetRandom(-SCALARMAX, SCALARMAX);
01483       SType2B[j*LDB+k] = ConvertType(SType1B[j*LDB+k], convertTo);
01484       }
01485   }
01486       
01487       SType1alpha = GetRandom(-SCALARMAX, SCALARMAX);
01488       SType2alpha = ConvertType(SType1alpha, convertTo);
01489       
01490       if(debug)
01491   {
01492     std::cout << "Test #" << TotalTestCount << " --" << std::endl;
01493     std::cout << Teuchos::ESideChar[SIDE] << "\t" 
01494          << Teuchos::EUploChar[UPLO] << "\t" 
01495          << Teuchos::ETranspChar[TRANSA] << "\t" 
01496          << Teuchos::EDiagChar[DIAG] << std::endl;
01497     std::cout << "M="<<M << "\t" << "N="<<N << "\t" << "LDA="<<LDA << "\t" << "LDB="<<LDB << std::endl;
01498     std::cout << "SType1alpha = " << SType1alpha << std::endl;
01499     std::cout << "SType2alpha = " << SType2alpha << std::endl;
01500     if (Teuchos::ESideChar[SIDE] == 'L') {
01501         PrintMatrix(SType1A, M, M, LDA, "SType1A", matlab);
01502         PrintMatrix(SType2A, M, M, LDA, "SType2A", matlab);
01503     } else {
01504         PrintMatrix(SType1A, N, N, LDA, "SType1A", matlab);
01505         PrintMatrix(SType2A, N, N, LDA, "SType2A", matlab);
01506     }
01507     PrintMatrix(SType1B, M, N, LDB, "SType1B_before_operation", matlab);
01508     PrintMatrix(SType2B, M, N, LDB, "SType2B_before_operation", matlab);
01509   }
01510       TotalTestCount++;
01511 
01512       SType1BLAS.TRSM(SIDE, UPLO, TRANSA, DIAG, M, N, SType1alpha, SType1A, LDA, SType1B, LDB);
01513       SType2BLAS.TRSM(SIDE, UPLO, TRANSA, DIAG, M, N, SType2alpha, SType2A, LDA, SType2B, LDB);
01514  
01515       if(debug)
01516   {
01517     PrintMatrix(SType1B, M, N, LDB, "SType1B_after_operation", matlab);
01518     PrintMatrix(SType2B, M, N, LDB, "SType2B_after_operation", matlab);
01519   }
01520 
01521       if (CompareMatrices(SType1B, SType2B, M, N, LDB, TOL)==0)
01522   std::cout << "FAILED TEST!!!!!!" << std::endl;
01523       GoodTestSubcount += CompareMatrices(SType1B, SType2B, M, N, LDB, TOL);
01524 
01525       delete [] SType1A;
01526       delete [] SType1B;
01527       delete [] SType2A;
01528       delete [] SType2B;
01529     }
01530   GoodTestCount += GoodTestSubcount; 
01531   if(verbose || debug) std::cout << "TRSM: " << GoodTestSubcount << " of " << TRSMTESTS << " tests were successful." << std::endl;
01532   if(debug) std::cout << std::endl;
01533   //--------------------------------------------------------------------------------
01534   // End TRSM Tests
01535   //--------------------------------------------------------------------------------
01536 
01537   if((((TotalTestCount - 1) - GoodTestCount) != 0) || (verbose) || (debug))
01538     {
01539       std::cout << GoodTestCount << " of " << (TotalTestCount - 1) << " total tests were successful." << std::endl;
01540     }
01541 
01542   if ((TotalTestCount-1) == GoodTestCount) {
01543     std::cout << "End Result: TEST PASSED" << std::endl;
01544     return 0;
01545   }
01546 
01547   std::cout << "End Result: TEST FAILED" << std::endl;
01548   return (TotalTestCount-GoodTestCount-1);
01549 }
01550 
01551 template<typename TYPE>
01552 TYPE GetRandom(TYPE Low, TYPE High)
01553 {
01554   return ((TYPE)((double)((1.0 * ScalarTraits<int>::random()) / RAND_MAX) * (High - Low + 1)) + Low);
01555 }
01556 
01557 template<>
01558 int GetRandom(int Low, int High)
01559 {
01560   return ((int)((double)((1.0 * ScalarTraits<int>::random()) / RAND_MAX) * (High - Low + 1)) + Low);
01561 }
01562 
01563 template<>
01564 double GetRandom(double Low, double High)
01565 {
01566   return (((double)((1.0 * ScalarTraits<int>::random()) / RAND_MAX) * (High - Low + 1)) + Low + ScalarTraits<double>::random());
01567 }
01568 
01569 template<typename TYPE>
01570 void PrintVector(TYPE* Vector, OType Size, std::string Name, bool Matlab)
01571 {
01572   std::cout << Name << " =" << std::endl;
01573   OType i;
01574   if(Matlab) std::cout << "[";
01575   for(i = 0; i < Size; i++)
01576     {
01577       std::cout << Vector[i] << " ";
01578     }
01579   if(Matlab) std::cout << "]";
01580   if(!Matlab)
01581     {
01582       std::cout << std::endl << std::endl;
01583     }
01584   else
01585     {
01586       std::cout << ";" << std::endl;
01587     }
01588 }
01589 
01590 template<typename TYPE>
01591 void PrintMatrix(TYPE* Matrix, OType Rows, OType Columns, OType LDM, std::string Name, bool Matlab)
01592 {
01593   if(!Matlab)
01594     {
01595       std::cout << Name << " =" << std::endl;
01596       OType i, j;
01597       for(i = 0; i < Rows; i++)
01598   {
01599           for(j = 0; j < Columns; j++)
01600       {
01601         std::cout << Matrix[i + (j * LDM)] << " ";
01602       }
01603     std::cout << std::endl;
01604   }
01605       std::cout << std::endl;
01606     }
01607   else
01608     {
01609       std::cout << Name << " = ";
01610       OType i, j;
01611       std::cout << "[";
01612       for(i = 0; i < Rows; i++)
01613         {
01614     std::cout << "[";
01615           for(j = 0; j < Columns; j++)
01616       {
01617         std::cout << Matrix[i + (j * LDM)] << " ";
01618       }
01619     std::cout << "];";
01620   }
01621       std::cout << "];" << std::endl;
01622     }
01623 }
01624 
01625 template<typename TYPE1, typename TYPE2>
01626 bool CompareScalars(TYPE1 Scalar1, TYPE2 Scalar2, typename ScalarTraits<TYPE2>::magnitudeType Tolerance)
01627 {
01628   TYPE2 convertTo = ScalarTraits<TYPE2>::zero();
01629   typename ScalarTraits<TYPE2>::magnitudeType temp = ScalarTraits<TYPE2>::magnitude(Scalar2);
01630   typename ScalarTraits<TYPE2>::magnitudeType temp2 = ScalarTraits<TYPE2>::magnitude(ConvertType(Scalar1, convertTo) - Scalar2);
01631   if (temp != ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero()) {
01632     temp2 /= temp;
01633   }
01634   return( temp2 < Tolerance );
01635 }
01636 
01637 
01638 /*  Function:  CompareVectors
01639     Purpose:   Compares the difference between two vectors using relative euclidean-norms, i.e. ||v_1-v_2||_2/||v_2||_2
01640 */
01641 template<typename TYPE1, typename TYPE2>
01642 bool CompareVectors(TYPE1* Vector1, TYPE2* Vector2, OType Size, typename ScalarTraits<TYPE2>::magnitudeType Tolerance)
01643 {
01644   TYPE2 convertTo = ScalarTraits<TYPE2>::zero();
01645   TYPE2 temp = ScalarTraits<TYPE2>::zero();
01646   typename ScalarTraits<TYPE2>::magnitudeType temp2 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01647   typename ScalarTraits<TYPE2>::magnitudeType temp3 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01648   typename ScalarTraits<TYPE2>::magnitudeType sum = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01649   typename ScalarTraits<TYPE2>::magnitudeType sum2 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01650   OType i;
01651   for(i = 0; i < Size; i++)
01652     {
01653       sum2 += ScalarTraits<TYPE2>::magnitude(ScalarTraits<TYPE2>::conjugate(Vector2[i])*Vector2[i]);
01654       temp = ConvertType(Vector1[i], convertTo) - Vector2[i];
01655       sum += ScalarTraits<TYPE2>::magnitude(ScalarTraits<TYPE2>::conjugate(temp)*temp);
01656     }
01657   temp2 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::squareroot(sum2);
01658   if (temp2 != ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero())
01659     temp3 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::squareroot(sum)/temp2;
01660   else
01661     temp3 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::squareroot(sum);
01662   if (temp3 > Tolerance )
01663     return false;
01664   else
01665     return true;
01666 }
01667 
01668 /*  Function:  CompareMatrices
01669     Purpose:   Compares the difference between two matrices using relative frobenius-norms, i.e. ||M_1-M_2||_F/||M_2||_F
01670 */
01671 template<typename TYPE1, typename TYPE2>
01672 bool CompareMatrices(TYPE1* Matrix1, TYPE2* Matrix2, OType Rows, OType Columns, OType LDM, typename ScalarTraits<TYPE2>::magnitudeType Tolerance)
01673 {
01674   TYPE2 convertTo = ScalarTraits<TYPE2>::zero();
01675   TYPE2 temp = ScalarTraits<TYPE2>::zero();
01676   typename ScalarTraits<TYPE2>::magnitudeType temp2 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01677   typename ScalarTraits<TYPE2>::magnitudeType temp3 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01678   typename ScalarTraits<TYPE2>::magnitudeType sum = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01679   typename ScalarTraits<TYPE2>::magnitudeType sum2 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero();
01680   OType i,j;
01681   for(j = 0; j < Columns; j++)
01682     {
01683       for(i = 0; i < Rows; i++)
01684   {
01685     sum2 = ScalarTraits<TYPE2>::magnitude(ScalarTraits<TYPE2>::conjugate(Matrix2[j*LDM + i])*Matrix2[j*LDM + i]);
01686     temp = ConvertType(Matrix1[j*LDM + i],convertTo) - Matrix2[j*LDM + i]; 
01687     sum = ScalarTraits<TYPE2>::magnitude(ScalarTraits<TYPE2>::conjugate(temp)*temp);
01688   }
01689     }
01690   temp2 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::squareroot(sum2);
01691   if (temp2 != ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::zero())
01692     temp3 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::squareroot(sum)/temp2;
01693   else
01694     temp3 = ScalarTraits<typename ScalarTraits<TYPE2>::magnitudeType>::squareroot(sum);
01695   if (temp3 > Tolerance)
01696     return false;
01697   else
01698     return true;
01699 }
01700 
01701 
01702 template<typename TYPE1, typename TYPE2>
01703 TYPE2 ConvertType(TYPE1 T1, TYPE2 T2)
01704 {
01705   return static_cast<TYPE2>(T1);
01706 }
01707 
01708 Teuchos::ESide RandomSIDE()
01709 {
01710   Teuchos::ESide result;
01711   int r = GetRandom(1, 2);
01712   if(r == 1)
01713     {
01714       result = Teuchos::LEFT_SIDE;
01715     }
01716   else
01717     {
01718       result = Teuchos::RIGHT_SIDE;
01719     }
01720   return result;
01721 }
01722 
01723 Teuchos::EUplo RandomUPLO()
01724 {
01725   Teuchos::EUplo result;
01726   int r = GetRandom(1, 2);
01727   if(r == 1)
01728     {
01729       result = Teuchos::UPPER_TRI;
01730     }
01731   else
01732     {
01733       result = Teuchos::LOWER_TRI;
01734     }
01735   return result;
01736 }
01737 
01738 Teuchos::ETransp RandomTRANS()
01739 {
01740   Teuchos::ETransp result;
01741   int r = GetRandom(1, 4);
01742   if(r == 1 || r == 2)
01743     {
01744       result = Teuchos::NO_TRANS;
01745     }
01746   else if(r == 3)
01747     {
01748       result = Teuchos::TRANS;
01749     }
01750   else
01751     {
01752       result = Teuchos::CONJ_TRANS;
01753     }
01754   return result;
01755 }
01756 
01757 Teuchos::EDiag RandomDIAG()
01758 {
01759   Teuchos::EDiag result;
01760   int r = GetRandom(1, 2);
01761   if(r == 1)
01762     {
01763       result = Teuchos::NON_UNIT_DIAG;
01764     }
01765   else
01766     {
01767       result = Teuchos::UNIT_DIAG;
01768     }
01769   return result;
01770 }
01771 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

Generated on Tue Oct 20 10:13:59 2009 for Teuchos Package Browser (Single Doxygen Collection) by  doxygen 1.6.1