Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
LL/memorytest_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 // This program tests the memory management system of the class CrsMatrix (memory leak, invalid free).
43 // It should be run using valgrind.
44 //
45 // Initially written to demonstrate bug #5499 and regressions introduced by the first patch.
46 
47 #include <vector>
48 
49 #include "Epetra_ConfigDefs.h"
50 #include "Epetra_Map.h"
51 #include "Epetra_CrsMatrix.h"
52 #ifdef EPETRA_MPI
53 #include "Epetra_MpiComm.h"
54 #include "mpi.h"
55 #else
56 #include "Epetra_SerialComm.h"
57 #endif
58 //#include "../epetra_test_err.h"
59 #include "Epetra_Version.h"
60 
61 int main(int argc, char *argv[])
62 {
63  int ierr = 0;
64 
65 #ifdef EPETRA_MPI
66 
67  // Initialize MPI
68 
69  MPI_Init(&argc, &argv);
70  int rank; // My process ID
71 
72  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
73  Epetra_MpiComm Comm( MPI_COMM_WORLD );
74 
75 #else
76 
77  int rank = 0;
78  Epetra_SerialComm Comm;
79 
80 #endif
81 
82  bool verbose = false;
83 
84  // Check if we should print results to standard out
85  if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
86 
87  int verbose_int = verbose ? 1 : 0;
88  Comm.Broadcast(&verbose_int, 1, 0);
89  verbose = verbose_int==1 ? true : false;
90 
91  Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
92  int MyPID = Comm.MyPID();
93  int NumProc = Comm.NumProc();
94 
95  if(verbose && MyPID==0)
96  std::cout << Epetra_Version() << std::endl << std::endl;
97 
98  if (verbose) std::cout << "Processor "<<MyPID<<" of "<< NumProc
99  << " is alive."<< std::endl;
100 
101  // unused: bool verbose1 = verbose;
102 
103  // Redefine verbose to only print on PE 0
104  if(verbose && rank!=0)
105  verbose = false;
106 
107  if (verbose) std::cout << "Test the memory management system of the class CrsMatrix (memory leak, invalid free)" << std::endl;
108 
109  //
110  // Test 1: code initially proposed to illustrate bug #5499
111  //
112 
113  if(Comm.NumProc() == 1) { // this is a sequential test
114 
115  if (verbose) std::cout << "* Using Copy, ColMap, Variable number of indices per row and Static profile (cf. bug #5499)." << std::endl;
116 
117  // Row Map
118  Epetra_Map RowMap(2LL, 0LL, Comm);
119 
120  // ColMap
121  std::vector<long long> colids(2);
122  colids[0]=0;
123  colids[1]=1;
124  Epetra_Map ColMap(-1LL, 2, &colids[0], 0LL, Comm);
125 
126  // NumEntriesPerRow
127  std::vector<int> NumEntriesPerRow(2);
128  NumEntriesPerRow[0]=2;
129  NumEntriesPerRow[1]=2;
130 
131  // Test
132  Epetra_CrsMatrix A(Copy, RowMap, ColMap, &NumEntriesPerRow[0], true);
133  // Bug #5499 shows up because InsertGlobalValues() is not called (CrsMatrix::Values_ not allocated but freed)
134  A.FillComplete();
135 
136  }
137 
138  //
139  // Test 1 Bis: same as Test1, but without ColMap and variable number of indices per row. Does not seems to matter
140  //
141 
142  if(Comm.NumProc() == 1) { // this is a sequential test
143 
144  if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Static profile" << std::endl;
145 
146  Epetra_Map RowMap(2LL, 0LL, Comm);
147 
148  // Test
149  Epetra_CrsMatrix A(Copy, RowMap, 1, true);
150  // Bug #5499 shows up because InsertGlobalValues() is not called (CrsMatrix::Values_ not allocated but freed)
151  A.FillComplete();
152 
153  }
154 
155  //
156  // Test 2: same as Test 1 Bis but with one call to InsertGlobalValues.
157  //
158 
159  if(Comm.NumProc() == 1) {
160 
161  if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Static profile + InsertGlobalValues()." << std::endl;
162 
163  Epetra_Map RowMap(2LL, 0LL, Comm);
164 
165  // Test
166  Epetra_CrsMatrix A(Copy, RowMap, 1, true);
167  std::vector<long long> Indices(1);
168  std::vector<double> Values(1);
169  Values[0] = 2;
170  Indices[0] = 0;
171 
172  A.InsertGlobalValues(0, 1, &Values[0], &Indices[0]); // Memory leak if CrsMatrix::Values not freed
173 
174  A.FillComplete();
175 
176  }
177 
178  //
179  // Test 3: check if the patch is not introducing some obvious regression
180  //
181 
182  if(Comm.NumProc() == 1) {
183 
184  if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Dynamic profile" << std::endl;
185 
186  Epetra_Map RowMap(2LL, 0LL, Comm);
187 
188  // Test
189  Epetra_CrsMatrix A(Copy, RowMap, 1, false);
190  A.FillComplete();
191 
192  }
193 
194  //
195  // Test 4: idem but with one call to InsertGlobalValues.
196  //
197 
198  if(Comm.NumProc() == 1) {
199 
200  if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Dynamic profile + InsertGlobalValues()." << std::endl;
201 
202  Epetra_Map RowMap(2LL, 0LL, Comm);
203 
204  // Test
205  Epetra_CrsMatrix A(Copy, RowMap, 1, false);
206  std::vector<long long> Indices(1);
207  std::vector<double> Values(1);
208  Values[0] = 2;
209  Indices[0] = 0;
210 
211  A.InsertGlobalValues(0, 1, &Values[0], &Indices[0]);
212  A.FillComplete();
213 
214  }
215 
216  if(Comm.NumProc() == 1) {
217 
218  if (verbose) std::cout << "* Using Copy, Static Graph()." << std::endl;
219 
220  Epetra_Map RowMap(1LL, 0LL, Comm);
221 
222  // Test
223  Epetra_CrsGraph G(Copy, RowMap, 1);
224  std::vector<long long> Indices(1);
225  Indices[0] = 0;
226  G.InsertGlobalIndices(0, 1, &Indices[0]);
227  G.FillComplete();
228 
229  Epetra_CrsMatrix A(Copy, G);
230  std::vector<double> Values(1);
231  Values[0] = 2;
232  A.ReplaceGlobalValues(0, 1, &Values[0], &Indices[0]);
233  A.FillComplete();
234  double norminf = A.NormInf();
235  if (verbose) std::cout << "** Inf Norm of Matrix = " << norminf << "." << std::endl;
236  std::cout << A << std::endl;
237 
238 
239  }
240 
241 
242  if(Comm.NumProc() == 1) {
243 
244  if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and static profile + InsertGlobalValues() for a single row." << std::endl;
245 
246  Epetra_Map RowMap(1LL, 0LL, Comm);
247 
248  // Test
249  Epetra_CrsMatrix A(Copy, RowMap, 1, true);
250  std::vector<long long> Indices(1);
251  std::vector<double> Values(1);
252  Values[0] = 2;
253  Indices[0] = 0;
254 
255  A.InsertGlobalValues(0, 1, &Values[0], &Indices[0]);
256  A.FillComplete();
257 
258  }
259 
260  /*
261  if (bool) {
262  if (verbose) std::cout << std::endl << "tests FAILED" << std::endl << std::endl;
263  }
264  else {*/
265  if (verbose) std::cout << std::endl << "tests PASSED" << std::endl << std::endl;
266  /* } */
267 
268 #ifdef EPETRA_MPI
269  MPI_Finalize();
270 #endif
271 
272  return ierr;
273 }
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
int FillComplete()
Tranform to local index space. Perform other operations to allow optimal matrix operations.
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.
double NormInf() const
Returns the infinity norm of the global matrix.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
int MyPID() const
Return my process ID.
Epetra_MpiComm: The Epetra MPI Communication Class.
std::string Epetra_Version()
int FillComplete(bool OptimizeDataStorage=true)
Signal that data entry is complete. Perform transformations to local index space. ...
int NumProc() const
Returns total number of processes (always returns 1 for SerialComm).
int InsertGlobalIndices(int GlobalRow, int NumIndices, int *Indices)
Enter a list of elements in a specified global row of the graph.
Epetra_SerialComm: The Epetra Serial Communication Class.
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
int main(int argc, char *argv[])
int Broadcast(double *MyVals, int Count, int Root) const
Epetra_SerialComm Broadcast function.
virtual int ReplaceGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Replace specified existing values with this list of entries for a given global row of the matrix...
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs. ...