EpetraExt Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/BTF/cxx_main.cpp
Go to the documentation of this file.
1//@HEADER
2// ***********************************************************************
3//
4// EpetraExt: Epetra Extended - Linear Algebra Services Package
5// Copyright (2011) Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ***********************************************************************
40//@HEADER
41
42// CrsGraph_BTF Test routine
43
45#include "EpetraExt_Version.h"
46
47#include "Epetra_Time.h"
48#ifdef EPETRA_MPI
49#include "Epetra_MpiComm.h"
50#include <mpi.h>
51#endif
52
53#include "Epetra_SerialComm.h"
54#include "Epetra_Map.h"
55#include "Epetra_CrsGraph.h"
57#include "../epetra_test_err.h"
58
59
60int main(int argc, char *argv[]) {
61
62 int i, ierr=0, returnierr=0;
63
64#ifdef EPETRA_MPI
65
66 // Initialize MPI
67
68 MPI_Init(&argc,&argv);
69 int size, rank; // Number of MPI processes, My process ID
70
71 MPI_Comm_size(MPI_COMM_WORLD, &size);
72 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
73
74#else
75
76 int size = 1; // Serial case (not using MPI)
77 int rank = 0;
78
79#endif
80
81 bool verbose = false;
82
83 // Check if we should print results to standard out
84 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
85
86
87#ifdef EPETRA_MPI
88 Epetra_MpiComm Comm(MPI_COMM_WORLD);
89#else
91#endif
92 if (!verbose) Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
93
94 int MyPID = Comm.MyPID();
95 int NumProc = Comm.NumProc();
96
97 if (verbose) {
98 cout << EpetraExt::EpetraExt_Version() << endl << endl;
99 cout << Comm << endl << flush;
100 }
101
102 Comm.Barrier();
103 bool verbose1 = verbose;
104
105 if (verbose) verbose = (MyPID==0);
106
107 int NumMyElements = 3;
108 int NumGlobalElements = NumMyElements;
109 int IndexBase = 0;
110
111 cout << "MyPID: " << MyPID << ", NumMyElements: " << NumMyElements << endl;
112
113 Epetra_Map Map( NumMyElements, 0, Comm );
114 cout << Map << endl;
115
116 Epetra_CrsGraph Graph( Copy, Map, 1 );
117
118 int index = 2;
119 Graph.InsertGlobalIndices( 0, 1, &index );
120 index = 0;
121 Graph.InsertGlobalIndices( 1, 1, &index );
122 index = 1;
123 Graph.InsertGlobalIndices( 2, 1, &index );
124
125 Graph.FillComplete();
126 cout << Graph << endl;
127
128 EpetraExt::CrsGraph_BTF BTFTransform;
129 Epetra_CrsGraph & NewGraph = BTFTransform( Graph );
130
131 cout << NewGraph << endl;
132
133#ifdef EPETRA_MPI
134 MPI_Finalize();
135#endif
136
137 return returnierr;
138}
139
Copy
Block Triangular Factorization (Reordering) of Epetra_CrsGraph.
int InsertGlobalIndices(int_type GlobalRow, int NumIndices, int_type *Indices)
void Barrier() const
int NumProc() const
int MyPID() const
static void SetTracebackMode(int TracebackModeValue)
std::string EpetraExt_Version()
int main(int argc, char *argv[])