Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
CrsMatrix/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_Map.h"
50#include "Epetra_CrsMatrix.h"
51#ifdef EPETRA_MPI
52#include "Epetra_MpiComm.h"
53#include "mpi.h"
54#else
55#include "Epetra_SerialComm.h"
56#endif
57//#include "../epetra_test_err.h"
58#include "Epetra_Version.h"
59
60int main(int argc, char *argv[])
61{
62 int ierr = 0;
63
64#ifdef EPETRA_MPI
65
66 // Initialize MPI
67
68 MPI_Init(&argc, &argv);
69 int rank; // My process ID
70
71 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
72 Epetra_MpiComm Comm( MPI_COMM_WORLD );
73
74#else
75
76 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 int verbose_int = verbose ? 1 : 0;
87 Comm.Broadcast(&verbose_int, 1, 0);
88 verbose = verbose_int==1 ? true : false;
89
90 Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
91 int MyPID = Comm.MyPID();
92 int NumProc = Comm.NumProc();
93
94 if(verbose && MyPID==0)
95 std::cout << Epetra_Version() << std::endl << std::endl;
96
97 if (verbose) std::cout << "Processor "<<MyPID<<" of "<< NumProc
98 << " is alive."<< std::endl;
99
100 // unused: bool verbose1 = verbose;
101
102 // Redefine verbose to only print on PE 0
103 if(verbose && rank!=0)
104 verbose = false;
105
106 if (verbose) std::cout << "Test the memory management system of the class CrsMatrix (memory leak, invalid free)" << std::endl;
107
108 //
109 // Test 1: code initially proposed to illustrate bug #5499
110 //
111
112 if(Comm.NumProc() == 1) { // this is a sequential test
113
114 if (verbose) std::cout << "* Using Copy, ColMap, Variable number of indices per row and Static profile (cf. bug #5499)." << std::endl;
115
116 // Row Map
117 Epetra_Map RowMap(2, 0, Comm);
118
119 // ColMap
120 std::vector<int> colids(2);
121 colids[0]=0;
122 colids[1]=1;
123 Epetra_Map ColMap(-1, 2, &colids[0], 0, Comm);
124
125 // NumEntriesPerRow
126 std::vector<int> NumEntriesPerRow(2);
127 NumEntriesPerRow[0]=2;
128 NumEntriesPerRow[1]=2;
129
130 // Test
131 Epetra_CrsMatrix A(Copy, RowMap, ColMap, &NumEntriesPerRow[0], true);
132 // Bug #5499 shows up because InsertGlobalValues() is not called (CrsMatrix::Values_ not allocated but freed)
133 A.FillComplete();
134
135 }
136
137 //
138 // Test 1 Bis: same as Test1, but without ColMap and variable number of indices per row. Does not seems to matter
139 //
140
141 if(Comm.NumProc() == 1) { // this is a sequential test
142
143 if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Static profile" << std::endl;
144
145 Epetra_Map RowMap(2, 0, Comm);
146
147 // Test
148 Epetra_CrsMatrix A(Copy, RowMap, 1, true);
149 // Bug #5499 shows up because InsertGlobalValues() is not called (CrsMatrix::Values_ not allocated but freed)
150 A.FillComplete();
151
152 }
153
154 //
155 // Test 2: same as Test 1 Bis but with one call to InsertGlobalValues.
156 //
157
158 if(Comm.NumProc() == 1) {
159
160 if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Static profile + InsertGlobalValues()." << std::endl;
161
162 Epetra_Map RowMap(2, 0, Comm);
163
164 // Test
165 Epetra_CrsMatrix A(Copy, RowMap, 1, true);
166 std::vector<int> Indices(1);
167 std::vector<double> Values(1);
168 Values[0] = 2;
169 Indices[0] = 0;
170
171 A.InsertGlobalValues(0, 1, &Values[0], &Indices[0]); // Memory leak if CrsMatrix::Values not freed
172
173 A.FillComplete();
174
175 }
176
177 //
178 // Test 3: check if the patch is not introducing some obvious regression
179 //
180
181 if(Comm.NumProc() == 1) {
182
183 if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Dynamic profile" << std::endl;
184
185 Epetra_Map RowMap(2, 0, Comm);
186
187 // Test
188 Epetra_CrsMatrix A(Copy, RowMap, 1, false);
189 A.FillComplete();
190
191 }
192
193 //
194 // Test 4: idem but with one call to InsertGlobalValues.
195 //
196
197 if(Comm.NumProc() == 1) {
198
199 if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and Dynamic profile + InsertGlobalValues()." << std::endl;
200
201 Epetra_Map RowMap(2, 0, Comm);
202
203 // Test
204 Epetra_CrsMatrix A(Copy, RowMap, 1, false);
205 std::vector<int> Indices(1);
206 std::vector<double> Values(1);
207 Values[0] = 2;
208 Indices[0] = 0;
209
210 A.InsertGlobalValues(0, 1, &Values[0], &Indices[0]);
211 A.FillComplete();
212
213 }
214
215 if(Comm.NumProc() == 1) {
216
217 if (verbose) std::cout << "* Using Copy, Static Graph()." << std::endl;
218
219 Epetra_Map RowMap(1, 0, Comm);
220
221 // Test
222 Epetra_CrsGraph G(Copy, RowMap, 1);
223 std::vector<int> Indices(1);
224 Indices[0] = 0;
225 G.InsertGlobalIndices(0, 1, &Indices[0]);
226 G.FillComplete();
227
229 std::vector<double> Values(1);
230 Values[0] = 2;
231 A.ReplaceGlobalValues(0, 1, &Values[0], &Indices[0]);
232 A.FillComplete();
233 double norminf = A.NormInf();
234 if (verbose) std::cout << "** Inf Norm of Matrix = " << norminf << "." << std::endl;
235 std::cout << A << std::endl;
236
237
238 }
239
240
241 if(Comm.NumProc() == 1) {
242
243 if (verbose) std::cout << "* Using Copy, Fixed number of indices per row and static profile + InsertGlobalValues() for a single row." << std::endl;
244
245 Epetra_Map RowMap(1, 0, Comm);
246
247 // Test
248 Epetra_CrsMatrix A(Copy, RowMap, 1, true);
249 std::vector<int> Indices(1);
250 std::vector<double> Values(1);
251 Values[0] = 2;
252 Indices[0] = 0;
253
254 A.InsertGlobalValues(0, 1, &Values[0], &Indices[0]);
255 A.FillComplete();
256
257 }
258
259 /*
260 if (bool) {
261 if (verbose) std::cout << std::endl << "tests FAILED" << std::endl << std::endl;
262 }
263 else {*/
264 if (verbose) std::cout << std::endl << "tests PASSED" << std::endl << std::endl;
265 /* } */
266
267#ifdef EPETRA_MPI
268 MPI_Finalize();
269#endif
270
271 return ierr;
272}
int main(int argc, char *argv[])
@ Copy
std::string Epetra_Version()
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.
int FillComplete()
Tranform to local index space. Perform other operations to allow optimal matrix operations.
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.
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.
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.
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.