00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029 #include "Epetra_Map.h"
00030 #include "Epetra_Time.h"
00031 #include "Epetra_CrsMatrix.h"
00032 #include "Epetra_Vector.h"
00033 #include "Epetra_SerialDenseVector.h"
00034 #include "Epetra_Flops.h"
00035 #ifdef EPETRA_MPI
00036 #include "Epetra_MpiComm.h"
00037 #include "mpi.h"
00038 #else
00039 #include "Epetra_SerialComm.h"
00040 #endif
00041 #include "Epetra_Version.h"
00042 #include "Trilinos_Util_CrsMatrixGallery.h"
00043 #include "EpetraExt_SubCopy_CrsMatrix.h"
00044 using namespace Trilinos_Util;
00045
00046
00047
00048
00049 int main(int argc, char *argv[])
00050 {
00051 int ierr = 0, i, forierr = 0;
00052 bool debug = false;
00053
00054 #ifdef EPETRA_MPI
00055
00056
00057
00058 MPI_Init(&argc,&argv);
00059 int size, rank;
00060
00061 MPI_Comm_size(MPI_COMM_WORLD, &size);
00062 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00063 Epetra_MpiComm Comm( MPI_COMM_WORLD );
00064
00065 #else
00066
00067 int size = 1;
00068 int rank = 0;
00069 Epetra_SerialComm Comm;
00070
00071 #endif
00072
00073 bool verbose = false;
00074
00075
00076 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00077
00078
00079
00080
00081
00082
00083 int MyPID = Comm.MyPID();
00084 int NumProc = Comm.NumProc();
00085
00086 if(verbose && MyPID==0)
00087 cout << Epetra_Version() << endl << endl;
00088
00089 if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
00090 << " is alive."<<endl;
00091
00092
00093 CrsMatrixGallery laplace_2d("laplace_2d", Comm);
00094 laplace_2d.Set("problem_size",Comm.NumProc()*Comm.NumProc()*4);
00095 Epetra_CrsMatrix * laplace_2d_matrix = laplace_2d.GetMatrix();
00096 if (verbose) cout << "Orig matrix = " << *laplace_2d_matrix << endl;
00097 const Epetra_Map * origRowMap = laplace_2d.GetMap();
00098 int * origGids = origRowMap->MyGlobalElements();
00099
00100
00101 Epetra_IntSerialDenseVector newRowMapGids(origRowMap->NumMyElements()/2 + 1);
00102 int numNewRowGids = 0;
00103 for (int i=0; i<origRowMap->NumMyElements(); i=i+2)
00104 newRowMapGids[numNewRowGids++] = origGids[i];
00105 Epetra_Map newRowMap(-1, numNewRowGids, newRowMapGids.Values(), 0, origRowMap->Comm());
00106
00107
00108
00109 EpetraExt::CrsMatrix_SubCopy subMatrixTransform(newRowMap);
00110
00111
00112 Epetra_CrsMatrix & subA = subMatrixTransform(*laplace_2d_matrix);
00113 if (verbose) cout << "Sub matrix (every other row/column) = " << subA << endl;
00114
00115
00116
00117 (*laplace_2d_matrix)[0][0] = 12.0;
00118 if (verbose) cout << "Orig matrix = " << *laplace_2d_matrix << endl;
00119 subMatrixTransform.fwd();
00120 assert(subA[0][0]==12.0);
00121 if (verbose) cout << "Sub matrix (every other row/column) = " << subA << endl;
00122
00123
00124
00125 subA[0][0] = 24.0;
00126 if (verbose) cout << "Sub matrix (every other row/column) = " << subA << endl;
00127 subMatrixTransform.rvs();
00128 assert((*laplace_2d_matrix)[0][0]==24.0);
00129 if (verbose) cout << "Orig matrix = " << *laplace_2d_matrix << endl;
00130
00131 if (Comm.MyPID()==0) cout << "EpetraExt::CrsMatrix_SubCopy tests passed." << endl;
00132 #ifdef EPETRA_MPI
00133 MPI_Finalize() ;
00134 #endif
00135
00136 return 0;
00137 }
00138
00139