58 int NumGlobalNonzeros1,
int* MyGlobalElements,
bool verbose);
63int main(
int argc,
char *argv[]) {
74 MPI_Init(&argc,&argv);
77 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
91 if(argv[1][0]==
'-' && argv[1][1]==
'v') {
102 int MyPID = Comm.
MyPID();
105 if(verbose && MyPID==0)
108 if(verbose) cout <<
"Processor "<<MyPID<<
" of "<< NumProc <<
" is alive." << endl;
110 bool verbose1 = verbose;
113 if(verbose && rank != 0)
116 int NumMyEquations = 5;
117 int NumGlobalEquations = NumMyEquations*NumProc+
EPETRA_MIN(NumProc,3);
132 int* NumNz =
new int[NumMyEquations];
137 for(i=0; i<NumMyEquations; i++)
138 if(MyGlobalElements[i]==0 || MyGlobalElements[i] == NumGlobalEquations-1)
153 Indices =
new int[2];
157 for(i = 0; i < NumMyEquations; i++) {
158 if(MyGlobalElements[i] == 0) {
162 else if(MyGlobalElements[i] == NumGlobalEquations-1) {
163 Indices[0] = NumGlobalEquations-2;
167 Indices[0] = MyGlobalElements[i]-1;
168 Indices[1] = MyGlobalElements[i]+1;
178 int gRID = A.
GRID(0);
180 std::vector<int> indices_vec(numIndices);
184 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
188 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
202 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
205 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
207 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
210 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
212 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
214 if(verbose) cout <<
"\n*****Testing variable entry constructor\n" << endl;
216 int NumMyNonzeros = 3 * NumMyEquations;
219 if(A.
LRID(NumGlobalEquations-1) >= 0)
222 EPETRA_TEST_ERR(
check(A, NumMyEquations, NumGlobalEquations, NumMyNonzeros, 3*NumGlobalEquations-2,
223 MyGlobalElements, verbose),ierr);
225 for(i = 0; i < NumMyEquations; i++)
228 for(i = 0; i < NumMyEquations; i++)
232 if(verbose) cout <<
"NumIndices function check OK" << endl;
238 if(verbose) cout <<
"\n*****Testing constant entry constructor\n" << endl;
244 for(i = 0; i < NumMyEquations; i++)
262 EPETRA_TEST_ERR(
check(AA, NumMyEquations, NumGlobalEquations, NumMyEquations, NumGlobalEquations,
263 MyGlobalElements, verbose),ierr);
268 for(i = 0; i < NumMyEquations; i++)
272 if(verbose) cout <<
"NumIndices function check OK" << endl;
276 if(verbose) cout <<
"\n*****Testing copy constructor\n" << endl;
281 EPETRA_TEST_ERR(
check(B, NumMyEquations, NumGlobalEquations, NumMyEquations, NumGlobalEquations,
282 MyGlobalElements, verbose),ierr);
285 for(i = 0; i < NumMyEquations; i++)
289 if(verbose) cout <<
"NumIndices function check OK" << endl;
293 if(verbose) cout <<
"\n*****Testing post construction modifications\n" << endl;
299 delete[] MyGlobalElements;
308 int NumMyElements1 = 4;
309 int NumMyEquations1 = NumMyElements1;
310 int NumGlobalEquations1 = NumMyEquations1*NumProc;
321 int* NumNz1 =
new int[NumMyEquations1];
326 for(i = 0; i < NumMyEquations1; i++)
327 if(MyGlobalElements1[i]==1 || MyGlobalElements1[i] == NumGlobalEquations1)
341 int* Indices1 =
new int[2];
345 for(i = 0; i < NumMyEquations1; i++) {
346 if(MyGlobalElements1[i]==1) {
350 else if(MyGlobalElements1[i] == NumGlobalEquations1) {
351 Indices1[0] = NumGlobalEquations1-1;
355 Indices1[0] = MyGlobalElements1[i]-1;
356 Indices1[1] = MyGlobalElements1[i]+1;
367 if(verbose) cout <<
"Print out tridiagonal matrix, each part on each processor. Index base is one.\n" << endl;
373 delete[] MyGlobalElements1;
381 if(verbose) cout <<
"\n*****Checking cpy ctr, op=, and reference counting." << endl;
384 if(verbose && (tempierr == 0)) cout <<
"Checked OK." << endl;
388 if(verbose) cout <<
"\n*****Checking shared-ownership tests." << endl;
391 if(verbose && (tempierr == 0)) cout <<
"Checked OK." << endl;
408 const int NumMyElements = 10;
409 const int IndexBase = 0;
410 Epetra_Map Map1(-1, NumMyElements, IndexBase, Comm);
412 const int NumIndicesPerRow = 5;
418 array1[0] = NumIndicesPerRow / 2;
419 array1[1] = array1[0] + 1;
421 for(
int i = 0; i < NumIndicesPerRow; i++)
424 int soleOutput, sharedOutput;
427 if(verbose) cout <<
"InsertMyIndices..." << endl;
432 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
435 if(verbose) cout <<
"RemoveMyIndices(#0)..." << endl;
440 if (ierr != 0) cout <<
"tests FAILED" << std::endl;
464 if(verbose) cout <<
"FillComplete..." << endl;
469 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
472 if(verbose) cout <<
"OptimizeStorage..." << endl;
477 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
480 if(verbose) cout <<
"RemoveMyIndices..." << endl;
485 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
488 if(verbose) cout <<
"RemoveMyIndices(#2)..." << endl;
493 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
496 if(verbose) cout <<
"FillComplete(#2)..." << endl;
501 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
510 int GlobalRow = SoleOwnerG.
GRID(0);
513 if(verbose) cout <<
"InsertGlobalIndices..." << endl;
518 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
521 if(verbose) cout <<
"RemoveGlobalIndices..." << endl;
526 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
529 if(verbose) cout <<
"RemoveGlobalIndices(#2)..." << endl;
534 if(verbose && ierr > 0) cout <<
"soleOutput = " << soleOutput <<
" sharedOutput = " << sharedOutput << endl;
552 const int NumIndicesPerRow = 10;
553 const int NumGlobalElements = 50;
554 const int IndexBase = 0;
555 Epetra_Map Map1(NumGlobalElements, IndexBase, Comm);
560 if(verbose) cout <<
"Graph1 created (def ctr). data addr = " << g1addr <<
" ref. count = " << g1count << endl;
569 if(verbose) cout <<
"Graph2 created (cpy ctr). data addr = " << g2addr <<
" ref. count = " << g2count << endl;
575 if(verbose) cout <<
"Graph2 destroyed. Graph1 data addr = " << g1newaddr <<
" ref. count = " << g1newcount << endl;
582 if(verbose) cout <<
"Graph3 created (op= before). data addr = " << g3addr <<
" ref. count = " << g3count << endl;
588 if(verbose) cout <<
"Graph3 set equal to Graph1 (op= after). data addr = " << g3addr <<
" ref. count = " << g3count << endl;
595 int NumGlobalNonzeros1,
int* MyGlobalElements,
bool verbose)
602 int NumGlobalIndices;
606 int* MyCopyIndices =
new int[MaxNumIndices];
607 int* GlobalCopyIndices =
new int[MaxNumIndices];
612 if(verbose) cout <<
"Number of local Rows = " << NumMyRows << endl;
617 if(verbose) cout <<
"Number of local Nonzero entries = " << NumMyNonzeros << endl;
622 if(verbose) cout <<
"Number of global Rows = " << NumGlobalRows << endl;
627 if(verbose) cout <<
"Number of global Nonzero entries = " << NumGlobalNonzeros << endl;
649 for(i = 0; i < NumMyRows; i++) {
653 forierr += !(NumGlobalIndices==NumMyIndices);
654 for(j = 1; j < NumMyIndices; j++)
EPETRA_TEST_ERR(!(MyViewIndices[j-1]<MyViewIndices[j]),ierr);
655 for(j = 0; j < NumGlobalIndices; j++) {
656 forierr += !(GlobalCopyIndices[j]==A.
GCID(MyViewIndices[j]));
657 forierr += !(A.
LCID(GlobalCopyIndices[j])==MyViewIndices[j]);
662 for(i = 0; i < NumMyRows; i++) {
666 forierr += !(NumGlobalIndices==NumMyIndices);
667 for(j = 1; j < NumMyIndices; j++)
669 for(j = 0; j < NumGlobalIndices; j++) {
670 forierr += !(GlobalCopyIndices[j]==A.
GCID(MyCopyIndices[j]));
671 forierr += !(A.
LCID(GlobalCopyIndices[j])==MyCopyIndices[j]);
677 delete[] MyCopyIndices;
678 delete[] GlobalCopyIndices;
680 if(verbose) cout <<
"Rows sorted check OK" << endl;
std::string Epetra_Version()
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
int MinMyGID() const
Returns the minimum global ID owned by this processor.
int NumMyElements() const
Number of elements on the calling processor.
int MaxMyGID() const
Returns the maximum global ID owned by this processor.
Epetra_Comm: The Epetra Communication Abstract Base Class.
Epetra_CrsGraphData: The Epetra CrsGraph Data Class.
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs.
int InsertGlobalIndices(int GlobalRow, int NumIndices, int *Indices)
Enter a list of elements in a specified global row of the graph.
bool Filled() const
If FillComplete() has been called, this query returns true, otherwise it returns false.
int ExtractMyRowView(int LocalRow, int &NumIndices, int *&Indices) const
Get a view of the elements in a specified local row of the graph.
int MaxNumIndices() const
Returns the maximum number of nonzero entries across all rows on this processor.
const Epetra_BlockMap & DomainMap() const
Returns the DomainMap associated with this graph.
int RemoveMyIndices(int LocalRow, int NumIndices, int *Indices)
Remove a list of elements from a specified local row of the graph.
bool IndicesAreLocal() const
If column indices are in local range, this query returns true, otherwise it returns false.
bool StorageOptimized() const
If OptimizeStorage() has been called, this query returns true, otherwise it returns false.
const Epetra_BlockMap & RowMap() const
Returns the RowMap associated with this graph.
int NumMyDiagonals() const
Returns the number of diagonal entries in the local graph, based on global row/column index compariso...
const Epetra_CrsGraphData * DataPtr() const
Returns a pointer to the CrsGraphData instance this CrsGraph uses.
int ReferenceCount() const
Returns the reference count of CrsGraphData.
int GCID(int LCID_in) const
Returns the global column index for give local column index, returns IndexBase-1 if we don't have thi...
bool LowerTriangular() const
If graph is lower triangular in local index space, this query returns true, otherwise it returns fals...
int ExtractGlobalRowView(int GlobalRow, int &NumIndices, int *&Indices) const
Get a view of the elements in a specified global row of the graph.
bool IndicesAreGlobal() const
If column indices are in global range, this query returns true, otherwise it returns false.
int RemoveGlobalIndices(int GlobalRow, int NumIndices, int *Indices)
Remove a list of elements from a specified global row of the graph.
int ExtractMyRowCopy(int LocalRow, int LenOfIndices, int &NumIndices, int *Indices) const
Extract a list of elements in a specified local row of the graph. Put into storage allocated by calli...
bool MyLRID(int LRID_in) const
Returns true if the LRID passed in belongs to the calling processor in this map, otherwise returns fa...
bool UpperTriangular() const
If graph is upper triangular in local index space, this query returns true, otherwise it returns fals...
int FillComplete()
Tranform to local index space. Perform other operations to allow optimal matrix operations.
int NumMyIndices(int Row) const
Returns the current number of nonzero entries in specified local row on this processor.
bool MyGRID(int GRID_in) const
Returns true if the GRID passed in belongs to the calling processor in this map, otherwise returns fa...
int GRID(int LRID_in) const
Returns the global row index for give local row index, returns IndexBase-1 if we don't have this loca...
int NumGlobalDiagonals() const
Returns the number of diagonal entries in the global graph, based on global row/column index comparis...
int InsertMyIndices(int LocalRow, int NumIndices, int *Indices)
Enter a list of elements in a specified local row of the graph.
int NumMyRows() const
Returns the number of matrix rows on this processor.
int NumGlobalNonzeros() const
Returns the number of indices in the global graph.
int NumGlobalRows() const
Returns the number of matrix rows in global matrix.
int ExtractGlobalRowCopy(int GlobalRow, int LenOfIndices, int &NumIndices, int *Indices) const
Extract a list of elements in a specified global row of the graph. Put into storage allocated by call...
int NumMyBlockDiagonals() const
Returns the number of Block diagonal entries in the local graph, based on global row/column index com...
int NumGlobalBlockDiagonals() const
Returns the number of Block diagonal entries in the global graph, based on global row/column index co...
bool MyGlobalRow(int GID) const
Returns true of GID is owned by the calling processor, otherwise it returns false.
int NumGlobalIndices(long long Row) const
Returns the current number of nonzero entries in specified global row on this processor.
const Epetra_BlockMap & RangeMap() const
Returns the RangeMap associated with this graph.
int NumMyNonzeros() const
Returns the number of indices in the local graph.
int OptimizeStorage()
Make consecutive row index sections contiguous, minimize internal storage used for constructing graph...
int LRID(int GRID_in) const
Returns the local row index for given global row index, returns -1 if no local row for this global ro...
bool NoDiagonal() const
If graph has no diagonal entries in global index space, this query returns true, otherwise it returns...
int LCID(int GCID_in) const
Returns the local column index for given global column index, returns -1 if no local column for this ...
Epetra_IntSerialDenseVector: A class for constructing and using dense vectors.
int * Values()
Returns pointer to the values in vector.
Epetra_Map: A class for partitioning vectors and matrices.
Epetra_MpiComm: The Epetra MPI Communication Class.
void Barrier() const
Epetra_MpiComm Barrier function.
int NumProc() const
Returns total number of processes.
int MyPID() const
Return my process ID.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
Epetra_SerialComm: The Epetra Serial Communication Class.
#define EPETRA_TEST_ERR(a, b)
int main(int argc, char *argv[])
int checkSharedOwnership(Epetra_Comm &Comm, bool verbose)
int check(Epetra_CrsGraph &A, int NumMyRows1, int NumGlobalRows1, int NumMyNonzeros1, int NumGlobalNonzeros1, int *MyGlobalElements, bool verbose)
int checkCopyAndAssignment(Epetra_Comm &Comm, bool verbose)