test/CrsRectMatrix/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 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 #include "Epetra_LocalMap.h"
00030 #include "Epetra_Map.h"
00031 #include "Epetra_Time.h"
00032 #include "Epetra_Vector.h"
00033 #include "Epetra_CrsMatrix.h"
00034 #include "Epetra_Flops.h"
00035 #include "Epetra_Export.h"
00036 #include "Epetra_Import.h"
00037 #include "Epetra_SerialDenseVector.h"
00038 #include "Epetra_IntSerialDenseVector.h"
00039 #ifdef EPETRA_MPI
00040 #include "Epetra_MpiComm.h"
00041 #include "mpi.h"
00042 #else
00043 #include "Epetra_SerialComm.h"
00044 #endif
00045 #include "../epetra_test_err.h"
00046 #include "Epetra_Version.h"
00047  
00048 int main(int argc, char *argv[])
00049 {
00050   int ierr = 0, i, j, forierr = 0;
00051   int ntrials = 1;
00052 
00053 #ifdef EPETRA_MPI
00054 
00055   // Initialize MPI
00056 
00057   MPI_Init(&argc,&argv);
00058   Epetra_MpiComm Comm( MPI_COMM_WORLD );
00059 
00060 #else
00061 
00062   Epetra_SerialComm Comm;
00063 
00064 #endif
00065 
00066   bool verbose = false;
00067 
00068   // Check if we should print results to standard out
00069   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00070 
00071 
00072 
00073   // char tmp; if (Comm.MyPID()==0) { cout << "Press any key to continue..."<< endl;  cin >> tmp;} Comm.Barrier();
00074 
00075   Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00076   int MyPID = Comm.MyPID();
00077   int NumProc = Comm.NumProc();
00078 
00079   if (verbose && Comm.MyPID()==0)
00080     cout << Epetra_Version() << endl << endl;
00081 
00082   if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00083               << " is alive."<<endl;
00084 
00085   bool verbose1 = verbose;
00086 
00087   // Redefine verbose to only print on PE 0
00088   if (verbose && Comm.MyPID()!=0) verbose = false;
00089 
00090   int NumMyEquations = 10000;
00091 
00092   int NumGlobalEquations = NumMyEquations*NumProc;
00093   int NumGlobalVariables = 2 * NumGlobalEquations+1;
00094 
00095   // Construct a Map that puts approximately the same Number of equations on each processor
00096 
00097   Epetra_Map RowMap(NumGlobalEquations, 0, Comm);
00098   Epetra_Map XMap(NumGlobalVariables, 0, Comm);
00099   Epetra_Map& YMap = RowMap;
00100   
00101   // Get update list and number of local equations from newly created Map
00102   int * MyGlobalElements = new int[RowMap.NumMyElements()];
00103   RowMap.MyGlobalElements(MyGlobalElements);
00104 
00105   // Get update list and number of local equations from newly created XMap
00106   int * XGlobalElements = new int[XMap.NumMyElements()];
00107   XMap.MyGlobalElements(XGlobalElements);
00108 
00109   // Get update list and number of local variables from newly created YMap
00110   int * YGlobalElements = new int[YMap.NumMyElements()];
00111   YMap.MyGlobalElements(YGlobalElements);
00112 
00113   // We need vectors to compute:
00114   // X = A^T*Y
00115   // AATY = A*A^T*Y = A*X
00116   //  and 
00117   // BY = B*Y
00118 
00119   Epetra_Vector Y(YMap);
00120   Epetra_Vector X(XMap);
00121   Epetra_Vector AATY(YMap);
00122   Epetra_Vector BY(YMap);
00123 
00124 
00125   // Fill Y Vector
00126   Y.Random();
00127   //Y.PutScalar(1.0);
00128 
00129   // To create A^T explicitly we need an assembly map that is two elements longer than
00130   // the XMap, because each processor will be making contributions to two rows beyond what
00131   // it will own.
00132   int ATAssemblyNumMyElements = 2*MyGlobalElements[NumMyEquations-1] + 2 - 2*MyGlobalElements[0] + 1;
00133   int * ATAssemblyGlobalElements = new int[ATAssemblyNumMyElements];
00134 
00135   for (i=0; i<ATAssemblyNumMyElements; i++) ATAssemblyGlobalElements[i] = 2*MyGlobalElements[0] + i;
00136   Epetra_Map ATAssemblyMap(-1, ATAssemblyNumMyElements, ATAssemblyGlobalElements, 0, Comm);
00137 
00138   // Create a Epetra_Matrix with the values of A
00139   // A is a simple 1D weighted average operator that mimics a restriction operator
00140   // that might be found in a multigrid code.
00141   // Also create A^T explicitly
00142 
00143   Epetra_CrsMatrix A(Copy, RowMap, 3);
00144   Epetra_CrsMatrix ATAssembly(Copy, ATAssemblyMap, 2);
00145   Epetra_CrsMatrix AT(Copy, XMap, 2);
00146   
00147   //cout << "ATAssemblyMap = "<< endl<< ATAssemblyMap << endl
00148   //     << "XMap = " << endl << XMap << endl
00149   //     << "RowMap = " << endl << RowMap << endl;
00150   // Add  rows one-at-a-time
00151   // Need some vectors to help
00152   // Off diagonal Values will always be -1
00153 
00154 
00155   double *Values = new double[3];
00156   int *Indices = new int[3];
00157   int NumEntries;
00158   /*
00159   Values[0] = 0.25;
00160   Values[1] = 0.5;
00161   Values[2] = 0.25;
00162   */
00163   Values[0] = 0.5;
00164   Values[1] = 0.25;
00165   Values[2] = 0.25;
00166   forierr = 0;
00167   for (i=0; i<NumMyEquations; i++)
00168     {
00169   /*
00170       Indices[0] = 2*MyGlobalElements[i];
00171       Indices[1] = 2*MyGlobalElements[i]+1;
00172       Indices[2] = 2*MyGlobalElements[i]+2;
00173    */
00174       Indices[0] = 2*MyGlobalElements[i]+1;
00175       Indices[1] = 2*MyGlobalElements[i]+2;
00176       Indices[2] = 2*MyGlobalElements[i];
00177       NumEntries = 3;
00178       forierr += !(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
00179       for (j=0; j<3; j++)
00180   forierr += !(ATAssembly.InsertGlobalValues(Indices[j],1, &(Values[j]), &(MyGlobalElements[i]))>=0);
00181     }
00182   EPETRA_TEST_ERR(forierr,ierr);
00183 
00184 
00185   EPETRA_TEST_ERR(!(ATAssembly.FillComplete()==0),ierr);
00186   // Gather AT values from ATAssembly matrix
00187   Epetra_Export Exporter(ATAssemblyMap, XMap);
00188   EPETRA_TEST_ERR(!(AT.Export(ATAssembly, Exporter, Add)==0),ierr);
00189 
00190   // Finish up
00191   EPETRA_TEST_ERR(!(A.FillComplete(XMap, YMap)==0),ierr);
00192   EPETRA_TEST_ERR(!(AT.FillComplete(YMap, XMap)==0),ierr);
00193 
00194 
00195   if (verbose1 && NumGlobalEquations<20) { 
00196     if (verbose) cout << "\n\n Matrix A\n" << endl;
00197     cout << A << endl;
00198     if (verbose) cout << " \n\n Matrix A Transpose\n" << endl;
00199     cout <<  AT << endl;
00200   }
00201 
00202 
00203   // Create a Epetra_Matrix containing B = A*A^T.
00204   // This matrix will be a square tridiagonal matrix.  We will use it to compare the results
00205   // of A*(A^T*X) using two methods: (1) with two calls to Multiply using A^T and then A and
00206   // (2) using B directly.
00207 
00208   Epetra_CrsMatrix B(Copy, RowMap, 3);
00209 
00210   Values[0] = 1.0/16.0;
00211   Values[1] = 3.0/8.0;
00212   Values[2] = 1.0/16.0;
00213   int Valstart;
00214   forierr = 0;
00215   for (i=0; i<NumMyEquations; i++)
00216     {
00217       if (MyGlobalElements[i] == 0) {
00218       Indices[0] = MyGlobalElements[i];
00219       Indices[1] = MyGlobalElements[i]+1;
00220       NumEntries = 2;
00221       Valstart = 1;
00222       }
00223       else {
00224   Indices[0] = MyGlobalElements[i]-1;
00225   Indices[1] = MyGlobalElements[i];
00226   Indices[2] = MyGlobalElements[i]+1;
00227   NumEntries = 3;
00228   Valstart = 0;
00229       }
00230       if (MyGlobalElements[i] == NumGlobalEquations-1) NumEntries--;
00231       forierr += !(B.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values+Valstart, Indices)==0);
00232     }
00233   EPETRA_TEST_ERR(forierr,ierr);
00234 
00235   // Finish up
00236   EPETRA_TEST_ERR(!(B.FillComplete()==0),ierr);
00237   if (verbose && NumGlobalEquations<20) cout << "\n\nMatrix B \n" << endl;
00238   if (verbose1 && NumGlobalEquations<20) cout << B << endl;
00239 
00240 
00241   Epetra_Flops counter;
00242   A.SetFlopCounter(counter);
00243   B.SetFlopCounter(A);
00244   Epetra_Time timer(Comm);
00245   for (i=0; i<ntrials; i++) {
00246     EPETRA_TEST_ERR(!(B.Multiply(false, Y, BY)==0),ierr); // Compute BY = B*Y
00247   }
00248   double elapsed_time = timer.ElapsedTime();
00249   double total_flops = B.Flops();
00250   counter.ResetFlops();
00251   
00252   double MFLOPs = total_flops/elapsed_time/1000000.0;
00253 
00254   if (verbose) cout << "\n\nTotal MFLOPs for B*Y = " << MFLOPs << endl<< endl;
00255   if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = B*Y \n";
00256   if (verbose1 && NumGlobalEquations<20) cout << BY << endl;
00257  
00258 
00260 
00261   timer.ResetStartTime();
00262   for (i=0; i<ntrials; i++) {
00263     EPETRA_TEST_ERR(!(A.Multiply(true, Y, X)==0),ierr); // Compute X = A^T*Y
00264   }
00265   elapsed_time = timer.ElapsedTime();
00266   total_flops = A.Flops();
00267   counter.ResetFlops();
00268   MFLOPs = total_flops/elapsed_time/1000000.0;
00269 
00270   if (verbose) cout << "\n\nTotal MFLOPs for A^T*Y using A and trans=true = " << MFLOPs << endl<< endl;
00271   if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = AT*Y \n";
00272   if (verbose1 && NumGlobalEquations<20) cout << X << endl;
00273 
00275 
00276   timer.ResetStartTime();
00277   EPETRA_TEST_ERR(!(A.Multiply(false, X, AATY)==0),ierr); // Compute AATY = A*X
00278   elapsed_time = timer.ElapsedTime();
00279   total_flops = A.Flops();
00280   MFLOPs = total_flops/elapsed_time/1000000.0;
00281   counter.ResetFlops();
00282   Epetra_Vector resid(YMap);
00283   resid.Update(1.0, BY, -1.0, AATY, 0.0);
00284   double residual;
00285   resid.Norm2(&residual);
00286 
00287   if (verbose) cout << "\n\nTotal MFLOPs for A*X using A and trans=false = " << MFLOPs << endl<< endl;
00288   if (verbose) cout << "Residual = " << residual << endl<< endl;
00289   if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = A*ATY \n";
00290   if (verbose1 && NumGlobalEquations<20) cout << AATY << endl;
00291 
00293 
00294   AT.SetFlopCounter(counter);
00295   timer.ResetStartTime();
00296   for (i=0; i<ntrials; i++) {
00297     EPETRA_TEST_ERR(!(AT.Multiply(false, Y, X)==0),ierr); // Compute X = A^T*Y
00298   }
00299   elapsed_time = timer.ElapsedTime();
00300   total_flops = AT.Flops();
00301   counter.ResetFlops();
00302   MFLOPs = total_flops/elapsed_time/1000000.0;
00303 
00304   if (verbose) cout << "\n\nTotal MFLOPs for A^T*Y using AT and trans=false = " << MFLOPs << endl<< endl;
00305   if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = AT*Y \n";
00306   if (verbose1 && NumGlobalEquations<20) cout << X << endl;
00307 
00309 
00310   timer.ResetStartTime();
00311   for (i=0; i<ntrials; i++) {
00312     EPETRA_TEST_ERR(!(AT.Multiply(true, X, AATY)==0),ierr); // Compute AATY = A*X
00313   }
00314   elapsed_time = timer.ElapsedTime();
00315   total_flops = AT.Flops();
00316   MFLOPs = total_flops/elapsed_time/1000000.0;
00317   counter.ResetFlops();
00318   resid.Update(1.0, BY, -1.0, AATY, 0.0);
00319   resid.Norm2(&residual);
00320 
00321   if (verbose) cout << "\n\nTotal MFLOPs for A*X using AT and trans=true = " << MFLOPs << endl<< endl;
00322   if (verbose) cout << "Residual = " << residual << endl<< endl;
00323   if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = A*ATY \n";
00324   if (verbose1 && NumGlobalEquations<20) cout <<AATY << endl;
00325 
00326 
00328   // Now test use of Epetra_LocalMap vectors: First test case of local replicated range vector
00329 
00330   {
00331     Epetra_CrsMatrix AL(Copy, RowMap, 3);
00332     for (i=0; i<NumMyEquations; i++)
00333       {
00334   forierr += !(A.ExtractGlobalRowCopy(MyGlobalElements[i], 3, NumEntries, Values, Indices)==0);
00335   forierr += !(AL.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
00336       }
00337     EPETRA_TEST_ERR(forierr,ierr);
00338 
00339     Epetra_LocalMap YLMap(NumGlobalEquations, 0, Comm);
00340     EPETRA_TEST_ERR(!(AL.FillComplete(XMap, YLMap)==0),ierr);
00341     AL.SetFlopCounter(A);
00342     Epetra_Vector YL(YLMap);
00343     Epetra_Vector ALX(YLMap);
00344   
00345     timer.ResetStartTime();
00346     for (i=0; i<ntrials; i++) {
00347       EPETRA_TEST_ERR(!(A.Multiply(false, X, Y)==0),ierr); // Compute Y= A*X
00348     }
00349     elapsed_time = timer.ElapsedTime();
00350     total_flops = A.Flops();
00351     counter.ResetFlops();
00352     MFLOPs = total_flops/elapsed_time/1000000.0;
00353 
00354 
00355 
00356     if (verbose) cout << "\n\nTotal MFLOPs for Y=A*X using global distributed Y = " << MFLOPs << endl<< endl;
00357     if (verbose && NumGlobalEquations<20) cout << "\n\nVector Y = A*X using distributed Y \n";
00358     if (verbose1 && NumGlobalEquations<20) cout << Y << endl;
00359     if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist Y range map\n";
00360     if (verbose1 && NumGlobalEquations<20) cout << A << endl;
00361 
00362     timer.ResetStartTime();
00363     for (i=0; i<ntrials; i++) {
00364       EPETRA_TEST_ERR(!(AL.Multiply(false, X, ALX)==0),ierr); // Compute YL= A*X
00365     }
00366     elapsed_time = timer.ElapsedTime();
00367     total_flops = AL.Flops();
00368     counter.ResetFlops();
00369     MFLOPs = total_flops/elapsed_time/1000000.0;
00370 
00371     if (verbose) cout << "\n\nTotal MFLOPs for Y=A*X using Local replicated Y = " << MFLOPs << endl<< endl;
00372     if (verbose && NumGlobalEquations<20) cout << "\n\nVector YL = AL*X using local replicated Y \n";
00373     if (verbose1 && NumGlobalEquations<20) cout << ALX << endl;
00374     if (verbose && NumGlobalEquations<20) cout << "\n\nA using local Y range map\n";
00375     if (verbose1 && NumGlobalEquations<20) cout << AL << endl;
00376 
00377     // Now gather Y values from the distributed Y and compare them to the local replicated Y values
00378     Epetra_Import g2limporter(YLMap, YMap); // Import from distributed Y map to local Y map
00379     EPETRA_TEST_ERR(!(YL.Import(Y, g2limporter, Insert)==0),ierr);
00380     if (verbose && NumGlobalEquations<20) cout << "\n\nVector YL = imported from distributed Y \n";
00381     if (verbose1 && NumGlobalEquations<20) cout << YL << endl;
00382     EPETRA_TEST_ERR(!(YL.Update(-1.0, ALX, 1.0)==0),ierr);
00383     EPETRA_TEST_ERR(!(YL.Norm2(&residual)==0),ierr);
00384     if (verbose) cout << "Residual = " << residual << endl<< endl;
00385       
00386 
00387     // 
00388     // Multiply by transpose
00389     //
00390 
00391     timer.ResetStartTime();
00392     for (i=0; i<ntrials; i++) {
00393       EPETRA_TEST_ERR(!(A.Multiply(true, Y, X)==0),ierr); // Compute X = A^TY
00394     }
00395     elapsed_time = timer.ElapsedTime();
00396     total_flops = A.Flops();
00397     counter.ResetFlops();
00398     MFLOPs = total_flops/elapsed_time/1000000.0;
00399 
00400 
00401 
00402     if (verbose) cout << "\n\nTotal MFLOPs for X=A^TY using global distributed Y = " << MFLOPs << endl<< endl;
00403     if (verbose && NumGlobalEquations<20) cout << "\n\nVector X using distributed Y \n";
00404     if (verbose1 && NumGlobalEquations<20) cout << Y << endl;
00405     if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist Y range map\n";
00406     if (verbose1 && NumGlobalEquations<20) cout << A << endl;
00407 
00408     Epetra_Vector X1(XMap);
00409 
00410     timer.ResetStartTime();
00411     for (i=0; i<ntrials; i++) {
00412       EPETRA_TEST_ERR(!(AL.Multiply(true, ALX, X1)==0),ierr); // Compute X1 = AL^T*Y
00413     }
00414     elapsed_time = timer.ElapsedTime();
00415     total_flops = AL.Flops();
00416     counter.ResetFlops();
00417     MFLOPs = total_flops/elapsed_time/1000000.0;
00418 
00419     if (verbose) cout << "\n\nTotal MFLOPs for X1=A^T*Y using Local replicated Y = " << MFLOPs << endl<< endl;
00420     if (verbose && NumGlobalEquations<20) cout << "\n\nVector X1 using local replicated Y \n";
00421     if (verbose1 && NumGlobalEquations<20) cout << X1 << endl;
00422 
00423     EPETRA_TEST_ERR(!(X1.Update(-1.0, X, 1.0)==0),ierr);
00424     EPETRA_TEST_ERR(!(X1.Norm2(&residual)==0),ierr);
00425     if (verbose) cout << "Residual = " << residual << endl<< endl;
00426   }
00427       
00429   // Finally test use of Epetra_LocalMap vectors using local replicated domain vector
00430 
00431   {
00432     Epetra_CrsMatrix AL(Copy, RowMap, 3);
00433     for (i=0; i<NumMyEquations; i++)
00434       {
00435   forierr += !(A.ExtractGlobalRowCopy(MyGlobalElements[i], 3, NumEntries, Values, Indices)==0);
00436   forierr += !(AL.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
00437       }
00438     EPETRA_TEST_ERR(forierr,ierr);
00439 
00440     Epetra_LocalMap XLMap(NumGlobalVariables, 0, Comm);
00441     EPETRA_TEST_ERR(!(AL.FillComplete(XLMap, YMap)==0),ierr);
00442     AL.SetFlopCounter(A);
00443     Epetra_Vector XL(XLMap);
00444     Epetra_Vector ALX(XLMap);
00445   
00446     timer.ResetStartTime();
00447     for (i=0; i<ntrials; i++) {
00448       EPETRA_TEST_ERR(!(A.Multiply(false, X, Y)==0),ierr); // Compute Y= A*X
00449     }
00450     elapsed_time = timer.ElapsedTime();
00451     total_flops = A.Flops();
00452     counter.ResetFlops();
00453     MFLOPs = total_flops/elapsed_time/1000000.0;
00454 
00455 
00456 
00457     if (verbose) cout << "\n\nTotal MFLOPs for Y=A*X using global distributed X = " << MFLOPs << endl<< endl;
00458     if (verbose && NumGlobalEquations<20) cout << "\n\nVector Y = A*X using distributed X \n";
00459     if (verbose1 && NumGlobalEquations<20) cout << Y << endl;
00460     //if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist X range map\n";
00461     //if (verbose1 && NumGlobalEquations<20) cout << A << endl;
00462 
00463     // Now gather X values from the distributed X 
00464     Epetra_Import g2limporter(XLMap, XMap); // Import from distributed X map to local X map
00465     EPETRA_TEST_ERR(!(XL.Import(X, g2limporter, Insert)==0),ierr);
00466     if (verbose && NumGlobalEquations<20) cout << "\n\nVector XL = imported from distributed X \n";
00467     if (verbose1 && NumGlobalEquations<20) cout << XL << endl;
00468     Epetra_Vector Y1(Y);
00469     timer.ResetStartTime();
00470     for (i=0; i<ntrials; i++) {
00471       EPETRA_TEST_ERR(!(AL.Multiply(false, XL, Y1)==0),ierr); // Compute Y1= AL*XL
00472     }
00473     elapsed_time = timer.ElapsedTime();
00474     total_flops = AL.Flops();
00475     counter.ResetFlops();
00476     MFLOPs = total_flops/elapsed_time/1000000.0;
00477 
00478     if (verbose) cout << "\n\nTotal MFLOPs for Y1=A*XL using Local replicated X = " << MFLOPs << endl<< endl;
00479     if (verbose && NumGlobalEquations<20) cout << "\n\nVector Y1 = AL*XL using local replicated X \n";
00480     if (verbose1 && NumGlobalEquations<20) cout << Y1 << endl;
00481     //if (verbose && NumGlobalEquations<20) cout << "\n\nA using local X domain map\n";
00482     //if (verbose1 && NumGlobalEquations<20) cout << AL << endl;
00483 
00484     EPETRA_TEST_ERR(!(Y1.Update(-1.0, Y, 1.0)==0),ierr);
00485     EPETRA_TEST_ERR(!(Y1.Norm2(&residual)==0),ierr);
00486     if (verbose) cout << "Residual = " << residual << endl<< endl;
00487       
00488 
00489     // 
00490     // Multiply by transpose
00491     //
00492 
00493     timer.ResetStartTime();
00494     for (i=0; i<ntrials; i++) {
00495       EPETRA_TEST_ERR(!(A.Multiply(true, Y, X)==0),ierr); // Compute X = A^TY
00496     }
00497     elapsed_time = timer.ElapsedTime();
00498     total_flops = A.Flops();
00499     counter.ResetFlops();
00500     MFLOPs = total_flops/elapsed_time/1000000.0;
00501 
00502 
00503 
00504     if (verbose) cout << "\n\nTotal MFLOPs for X=A^TY using global distributed X = " << MFLOPs << endl<< endl;
00505     if (verbose && NumGlobalEquations<20) cout << "\n\nVector X using distributed X \n";
00506     if (verbose1 && NumGlobalEquations<20) cout << X << endl;
00507     //if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist X domain map\n";
00508     //if (verbose1 && NumGlobalEquations<20) cout << A << endl;
00509 
00510     timer.ResetStartTime();
00511     for (i=0; i<ntrials; i++) {
00512       EPETRA_TEST_ERR(!(AL.Multiply(true, Y, XL)==0),ierr); // Compute XL = AL^T*Y
00513     }
00514     elapsed_time = timer.ElapsedTime();
00515     total_flops = AL.Flops();
00516     counter.ResetFlops();
00517     MFLOPs = total_flops/elapsed_time/1000000.0;
00518 
00519     if (verbose) cout << "\n\nTotal MFLOPs for XL=A^T*Y1 using Local replicated X = " << MFLOPs << endl<< endl;
00520     if (verbose && NumGlobalEquations<20) cout << "\n\nVector XL using local replicated X \n";
00521     if (verbose1 && NumGlobalEquations<20) cout << XL << endl;
00522 
00523     Epetra_Vector XL1(XLMap);
00524     EPETRA_TEST_ERR(!(XL1.Import(X, g2limporter, Insert)==0),ierr);
00525     EPETRA_TEST_ERR(!(XL1.Update(-1.0, XL, 1.0)==0),ierr);
00526     EPETRA_TEST_ERR(!(XL1.Norm2(&residual)==0),ierr);
00527     if (verbose) cout << "Residual = " << residual << endl<< endl;
00528   }     
00529   // Release all objects
00530   delete [] Values;
00531   delete [] Indices;
00532   delete [] MyGlobalElements;
00533   delete [] XGlobalElements;
00534   delete [] YGlobalElements;
00535   delete [] ATAssemblyGlobalElements;
00536 
00537 
00538 #ifdef EPETRA_MPI
00539   MPI_Finalize() ;
00540 #endif
00541 
00542 /* end main
00543 */
00544 return ierr ;
00545 }
00546 

Generated on Thu Sep 18 12:37:56 2008 for Epetra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1