Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/SimpleLongLongTest/cxx_main.cpp
Go to the documentation of this file.
1//@HEADER
2// ************************************************************************
3//
4// Epetra: 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
43#include "Epetra_Map.h"
44#include "Epetra_Time.h"
45#include "Epetra_CrsMatrix.h"
46#include "Epetra_Vector.h"
47#include "Epetra_Flops.h"
48#ifdef EPETRA_MPI
49#include "Epetra_MpiComm.h"
50#include "mpi.h"
51#else
52#include "Epetra_SerialComm.h"
53#endif
54#include "../epetra_test_err.h"
55#include "Epetra_Version.h"
56
57// prototypes
58
59int main(int argc, char *argv[])
60{
61 int ierr = 0;
62
63#ifdef EPETRA_MPI
64
65 // Initialize MPI
66
67 MPI_Init(&argc,&argv);
68 int rank; // My process ID
69
70 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
71 Epetra_MpiComm Comm( MPI_COMM_WORLD );
72
73#else
74
75 int rank = 0;
77
78#endif
79
80 bool verbose = false;
81
82 // Check if we should print results to standard out
83 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
84
85 int verbose_int = verbose ? 1 : 0;
86 Comm.Broadcast(&verbose_int, 1, 0);
87 verbose = verbose_int==1 ? true : false;
88
89
90 // char tmp;
91 // if (rank==0) cout << "Press any key to continue..."<< std::endl;
92 // if (rank==0) cin >> tmp;
93 // Comm.Barrier();
94
95 Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
96 int MyPID = Comm.MyPID();
97 int NumProc = Comm.NumProc();
98
99 if(verbose && MyPID==0)
100 cout << Epetra_Version() << std::endl << std::endl;
101
102 if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
103 << " is alive."<<endl;
104
105 //bool verbose1 = verbose; // unused
106
107 // Redefine verbose to only print on PE 0
108 if(verbose && rank!=0) verbose = false;
109
110 int NumMyEquations = 1;
111 long long NumGlobalEquations = NumProc;
112
113 // Get update list and number of local equations from newly created Map
114 long long* MyGlobalElementsLL = new long long[NumMyEquations];
115
116 MyGlobalElementsLL[0] = 2000000000+MyPID;
117
118 // Construct a Map that puts approximately the same Number of equations on each processor
119
120 Epetra_Map MapLL(NumGlobalEquations, NumMyEquations, MyGlobalElementsLL, 0, Comm);
121
122 EPETRA_TEST_ERR(MapLL.GlobalIndicesInt(),ierr);
123 EPETRA_TEST_ERR(!(MapLL.GlobalIndicesLongLong()),ierr);
124
125 // Create an integer vector NumNz that is used to build the Petra Matrix.
126 // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
127
128 int* NumNzLL = new int[NumMyEquations];
129 NumNzLL[0] = 0;
130
131 // Create int types meant to add to long long matrix for test of failure
132 int* MyIntGlobalElementsLL = new int[NumMyEquations];
133 MyIntGlobalElementsLL[0] = 20000+MyPID;
134
135 // Create a long long Epetra_Matrix
136 Epetra_CrsMatrix A_LL(Copy, MapLL, NumNzLL);
138 EPETRA_TEST_ERR(A_LL.IndicesAreLocal(),ierr);
139
140 // Insert values
141 double one = 1.0;
142#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
143 // Try to add ints which should fail and be caught as an int
144 try {
145 A_LL.InsertGlobalValues(MyIntGlobalElementsLL[0], 1, &one, MyIntGlobalElementsLL+0);
146 } catch(int i) {
147 EPETRA_TEST_ERR(!(i==-1),ierr);
148 }
149#endif
150
151 // Add long longs which should succeed
152 EPETRA_TEST_ERR(!(A_LL.InsertGlobalValues(MyGlobalElementsLL[0], 1, &one, MyGlobalElementsLL+0)==0),ierr);
153 EPETRA_TEST_ERR(!(A_LL.IndicesAreGlobal()),ierr);
154 EPETRA_TEST_ERR(!(A_LL.FillComplete(false)==0),ierr);
155 EPETRA_TEST_ERR(!(A_LL.IndicesAreLocal()),ierr);
156
157
158 // Get update list and number of local equations from newly created Map
159 int* MyGlobalElementsInt = new int[NumMyEquations];
160
161 MyGlobalElementsInt[0] = 2000+MyPID;
162
163 // Create an integer vector NumNz that is used to build the Petra Matrix.
164 // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
165
166 int* NumNzInt = new int[NumMyEquations];
167 NumNzInt[0] = 0;
168
169 // Create int types meant to add to long long matrix for test of failure
170 long long* MyLLGlobalElementsInt = new long long[NumMyEquations];
171 MyLLGlobalElementsInt[0] = 2000000000+MyPID;
172
173#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
174 // Construct a Map that puts approximately the same Number of equations on each processor
175
176 Epetra_Map MapInt(NumGlobalEquations, NumMyEquations, MyGlobalElementsInt, 0LL, Comm);
177
178 EPETRA_TEST_ERR(!(MapInt.GlobalIndicesInt()),ierr);
180
181 // Create a int Epetra_Matrix
182 Epetra_CrsMatrix A_Int(Copy, MapInt, NumNzInt);
183 EPETRA_TEST_ERR(A_Int.IndicesAreGlobal(),ierr);
184 EPETRA_TEST_ERR(A_Int.IndicesAreLocal(),ierr);
185
186 // Insert values
187 try {
188 A_Int.InsertGlobalValues(MyLLGlobalElementsInt[0], 1, &one, MyLLGlobalElementsInt+0);
189 } catch(int i) {
190 EPETRA_TEST_ERR(!(i==-1),ierr);
191 }
192 // Add long longs which should succeed
193 EPETRA_TEST_ERR(!(A_Int.InsertGlobalValues(MyGlobalElementsInt[0], 1, &one, MyGlobalElementsInt+0)==0),ierr);
194 EPETRA_TEST_ERR(!(A_Int.IndicesAreGlobal()),ierr);
195 EPETRA_TEST_ERR(!(A_Int.FillComplete(false)==0),ierr);
196 EPETRA_TEST_ERR(!(A_Int.IndicesAreLocal()),ierr);
197#endif
198
199 delete [] MyGlobalElementsLL;
200 delete [] NumNzLL;
201 delete [] MyIntGlobalElementsLL;
202 delete [] MyGlobalElementsInt;
203 delete [] NumNzInt;
204 delete [] MyLLGlobalElementsInt;
205
206#ifdef EPETRA_MPI
207 MPI_Finalize() ;
208#endif
209
210/* end main
211*/
212return ierr ;
213}
std::string Epetra_Version()
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
bool GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
int FillComplete(bool OptimizeDataStorage=true)
Signal that data entry is complete. Perform transformations to local index space.
bool IndicesAreLocal() const
If matrix indices has been transformed to local, this query returns true, otherwise it returns false.
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Insert a list of elements in a given global row of the matrix.
bool IndicesAreGlobal() const
If matrix indices has not been transformed to local, this query returns true, otherwise it returns fa...
Epetra_Map: A class for partitioning vectors and matrices.
Definition Epetra_Map.h:119
Epetra_MpiComm: The Epetra MPI Communication Class.
int Broadcast(double *MyVals, int Count, int Root) const
Epetra_MpiComm Broadcast 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[])