Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
test/RowMatrix_LL/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_JadMatrix.h"
47 #include "Epetra_Vector.h"
48 #include "Epetra_Flops.h"
49 #ifdef EPETRA_MPI
50 #include "Epetra_MpiComm.h"
51 #include "mpi.h"
52 #else
53 #include "Epetra_SerialComm.h"
54 #endif
55 #include "../epetra_test_err.h"
56 #include "Epetra_Version.h"
57 #include <vector>
58 #include <algorithm>
59 #include <string>
60 
61 // prototypes
62 
63 int checkValues( double x, double y, string message = "", bool verbose = false) {
64  if (fabs((x-y)/x) > 0.01) {
65  return(1);
66  if (verbose) cout << "********** " << message << " check failed.********** " << endl;
67  }
68  else {
69  if (verbose) cout << message << " check OK." << endl;
70  return(0);
71  }
72 }
73 
74 int checkMultiVectors( Epetra_MultiVector & X, Epetra_MultiVector & Y, string message = "", bool verbose = false) {
75  int numVectors = X.NumVectors();
76  int length = Y.MyLength();
77  int badvalue = 0;
78  int globalbadvalue = 0;
79  for (int j=0; j<numVectors; j++)
80  for (int i=0; i< length; i++)
81  if (checkValues(X[j][i], Y[j][i])==1) badvalue = 1;
82  X.Map().Comm().MaxAll(&badvalue, &globalbadvalue, 1);
83 
84  if (verbose) {
85  if (globalbadvalue==0) cout << message << " check OK." << endl;
86  else cout << "********* " << message << " check failed.********** " << endl;
87  }
88  return(globalbadvalue);
89 }
90 
91 int check(Epetra_RowMatrix & A, Epetra_RowMatrix & B, bool verbose);
92 
94  Epetra_Vector & q, Epetra_Vector & z, Epetra_Vector & resid, bool verbose);
95 
96 int power_method(bool TransA, Epetra_RowMatrix& A,
97  Epetra_Vector& q,
98  Epetra_Vector& z0,
99  Epetra_Vector& resid,
100  double * lambda, int niters, double tolerance,
101  bool verbose);
102 
103 int main(int argc, char *argv[])
104 {
105  int ierr = 0, i, forierr = 0;
106 #ifdef EPETRA_MPI
107 
108  // Initialize MPI
109 
110  MPI_Init(&argc,&argv);
111  int rank; // My process ID
112 
113  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
114  Epetra_MpiComm Comm( MPI_COMM_WORLD );
115 
116 #else
117 
118  int rank = 0;
119  Epetra_SerialComm Comm;
120 
121 #endif
122 
123  bool verbose = false;
124 
125  // Check if we should print results to standard out
126  if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
127 
128  int verbose_int = verbose ? 1 : 0;
129  Comm.Broadcast(&verbose_int, 1, 0);
130  verbose = verbose_int==1 ? true : false;
131 
132 
133  // char tmp;
134  // if (rank==0) cout << "Press any key to continue..."<< endl;
135  // if (rank==0) cin >> tmp;
136  // Comm.Barrier();
137 
138  Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
139  int MyPID = Comm.MyPID();
140  int NumProc = Comm.NumProc();
141 
142  if(verbose && MyPID==0)
143  cout << Epetra_Version() << endl << endl;
144 
145  if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
146  << " is alive."<<endl;
147 
148  // Redefine verbose to only print on PE 0
149  if(verbose && rank!=0)
150  verbose = false;
151 
152  int NumMyEquations = 10000;
153  long long NumGlobalEquations = (NumMyEquations * NumProc) + EPETRA_MIN(NumProc,3);
154  if(MyPID < 3)
155  NumMyEquations++;
156 
157  // Construct a Map that puts approximately the same Number of equations on each processor
158 
159  Epetra_Map Map(NumGlobalEquations, NumMyEquations, 0LL, Comm);
160 
161  // Get update list and number of local equations from newly created Map
162  vector<long long> MyGlobalElements(Map.NumMyElements());
163  Map.MyGlobalElements(&MyGlobalElements[0]);
164 
165  // Create an integer vector NumNz that is used to build the Petra Matrix.
166  // NumNz[i] is the Number of OFF-DIAGONAL term for the ith global equation on this processor
167 
168  vector<int> NumNz(NumMyEquations);
169 
170  // We are building a tridiagonal matrix where each row has (-1 2 -1)
171  // So we need 2 off-diagonal terms (except for the first and last equation)
172 
173  for(i = 0; i < NumMyEquations; i++)
174  if((MyGlobalElements[i] == 0) || (MyGlobalElements[i] == NumGlobalEquations - 1))
175  NumNz[i] = 1;
176  else
177  NumNz[i] = 2;
178 
179  // Create a Epetra_Matrix
180 
181  Epetra_CrsMatrix A(Copy, Map, &NumNz[0]);
184 
185  // Add rows one-at-a-time
186  // Need some vectors to help
187  // Off diagonal Values will always be -1
188 
189 
190  vector<double> Values(2);
191  Values[0] = -1.0;
192  Values[1] = -1.0;
193  vector<long long> Indices(2);
194  double two = 2.0;
195  int NumEntries;
196 
197  forierr = 0;
198  for(i = 0; i < NumMyEquations; i++) {
199  if(MyGlobalElements[i] == 0) {
200  Indices[0] = 1;
201  NumEntries = 1;
202  }
203  else if (MyGlobalElements[i] == NumGlobalEquations-1) {
204  Indices[0] = NumGlobalEquations-2;
205  NumEntries = 1;
206  }
207  else {
208  Indices[0] = MyGlobalElements[i]-1;
209  Indices[1] = MyGlobalElements[i]+1;
210  NumEntries = 2;
211  }
212  forierr += !(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, &Values[0], &Indices[0])==0);
213  forierr += !(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])>0); // Put in the diagonal entry
214  }
215  EPETRA_TEST_ERR(forierr,ierr);
216 
217  // Finish up
218  A.FillComplete();
219  A.OptimizeStorage();
220 
221  Epetra_JadMatrix JadA(A);
222  Epetra_JadMatrix JadA1(A);
223  Epetra_JadMatrix JadA2(A);
224 
225  // Create vectors for Power method
226 
227  Epetra_Vector q(Map);
228  Epetra_Vector z(Map); z.Random();
229  Epetra_Vector resid(Map);
230 
231  Epetra_Flops flopcounter;
232  A.SetFlopCounter(flopcounter);
233  q.SetFlopCounter(A);
234  z.SetFlopCounter(A);
235  resid.SetFlopCounter(A);
236  JadA.SetFlopCounter(A);
237  JadA1.SetFlopCounter(A);
238  JadA2.SetFlopCounter(A);
239 
240 
241  if (verbose) cout << "=======================================" << endl
242  << "Testing Jad using CrsMatrix as input..." << endl
243  << "=======================================" << endl;
244 
245  A.ResetFlops();
246  powerMethodTests(A, JadA, Map, q, z, resid, verbose);
247 
248  // Increase diagonal dominance
249 
250  if (verbose) cout << "\n\nIncreasing the magnitude of first diagonal term and solving again\n\n"
251  << endl;
252 
253 
254  if (A.MyGlobalRow(0)) {
255  int numvals = A.NumGlobalEntries(0);
256  vector<double> Rowvals(numvals);
257  vector<long long> Rowinds(numvals);
258  A.ExtractGlobalRowCopy(0, numvals, numvals, &Rowvals[0], &Rowinds[0]); // Get A[0,0]
259 
260  for (i=0; i<numvals; i++) if (Rowinds[i] == 0) Rowvals[i] *= 10.0;
261 
262  A.ReplaceGlobalValues(0, numvals, &Rowvals[0], &Rowinds[0]);
263  }
264  JadA.UpdateValues(A);
265  A.ResetFlops();
266  powerMethodTests(A, JadA, Map, q, z, resid, verbose);
267 
268  if (verbose) cout << "================================================================" << endl
269  << "Testing Jad using Jad matrix as input matrix for construction..." << endl
270  << "================================================================" << endl;
271  JadA1.ResetFlops();
272  powerMethodTests(JadA1, JadA2, Map, q, z, resid, verbose);
273 
274 #ifdef EPETRA_MPI
275  MPI_Finalize() ;
276 #endif
277 
278 return ierr ;
279 }
280 
282  Epetra_Vector & q, Epetra_Vector & z, Epetra_Vector & resid, bool verbose) {
283 
284  // variable needed for iteration
285  double lambda = 0.0;
286  // int niters = 10000;
287  int niters = 300;
288  double tolerance = 1.0e-2;
289  int ierr = 0;
290 
292 
293  // Iterate
294 
295  Epetra_Time timer(Map.Comm());
296 
297  double startTime = timer.ElapsedTime();
298  EPETRA_TEST_ERR(power_method(false, A, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
299  double elapsed_time = timer.ElapsedTime() - startTime;
300  double total_flops = q.Flops();
301  double MFLOPs = total_flops/elapsed_time/1000000.0;
302  double lambdaref = lambda;
303  double flopsref = total_flops;
304 
305  if (verbose)
306  cout << "\n\nTotal MFLOPs for reference first solve = " << MFLOPs << endl
307  << "Total FLOPS = " <<total_flops <<endl<<endl;
308 
309  lambda = 0.0;
310  startTime = timer.ElapsedTime();
311  EPETRA_TEST_ERR(power_method(false, JadA, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
312  elapsed_time = timer.ElapsedTime() - startTime;
313  total_flops = q.Flops();
314  MFLOPs = total_flops/elapsed_time/1000000.0;
315 
316  if (verbose)
317  cout << "\n\nTotal MFLOPs for candidate first solve = " << MFLOPs << endl
318  << "Total FLOPS = " <<total_flops <<endl<<endl;
319 
320  EPETRA_TEST_ERR(checkValues(lambda,lambdaref," No-transpose Power Method result", verbose),ierr);
321  EPETRA_TEST_ERR(checkValues(total_flops,flopsref," No-transpose Power Method flop count", verbose),ierr);
322 
324 
325  // Solve transpose problem
326 
327  if (verbose) cout << "\n\nUsing transpose of matrix and solving again (should give same result).\n\n"
328  << endl;
329  // Iterate
330  lambda = 0.0;
331  startTime = timer.ElapsedTime();
332  EPETRA_TEST_ERR(power_method(true, A, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
333  elapsed_time = timer.ElapsedTime() - startTime;
334  total_flops = q.Flops();
335  MFLOPs = total_flops/elapsed_time/1000000.0;
336  lambdaref = lambda;
337  flopsref = total_flops;
338 
339  if (verbose)
340  cout << "\n\nTotal MFLOPs for reference transpose solve = " << MFLOPs << endl
341  << "Total FLOPS = " <<total_flops <<endl<<endl;
342 
343  lambda = 0.0;
344  startTime = timer.ElapsedTime();
345  EPETRA_TEST_ERR(power_method(true, JadA, q, z, resid, &lambda, niters, tolerance, verbose),ierr);
346  elapsed_time = timer.ElapsedTime() - startTime;
347  total_flops = q.Flops();
348  MFLOPs = total_flops/elapsed_time/1000000.0;
349 
350  if (verbose)
351  cout << "\n\nTotal MFLOPs for candidate transpose solve = " << MFLOPs << endl
352  << "Total FLOPS = " <<total_flops <<endl<<endl;
353 
354  EPETRA_TEST_ERR(checkValues(lambda,lambdaref,"Transpose Power Method result", verbose),ierr);
355  EPETRA_TEST_ERR(checkValues(total_flops,flopsref,"Transpose Power Method flop count", verbose),ierr);
356 
357  EPETRA_TEST_ERR(check(A, JadA, verbose),ierr);
358 
359  return(0);
360 }
361 int power_method(bool TransA, Epetra_RowMatrix& A, Epetra_Vector& q, Epetra_Vector& z0,
362  Epetra_Vector& resid, double* lambda, int niters, double tolerance, bool verbose)
363 {
364 
365  // Fill z with random Numbers
366  Epetra_Vector z(z0);
367 
368  // variable needed for iteration
369  double normz, residual;
370 
371  int ierr = 1;
372 
373  for(int iter = 0; iter < niters; iter++) {
374  z.Norm2(&normz); // Compute 2-norm of z
375  q.Scale(1.0/normz, z);
376  A.Multiply(TransA, q, z); // Compute z = A*q // SEGFAULT HAPPENS HERE
377  q.Dot(z, lambda); // Approximate maximum eigenvaluE
378  if(iter%100==0 || iter+1==niters) {
379  resid.Update(1.0, z, -(*lambda), q, 0.0); // Compute A*q - lambda*q
380  resid.Norm2(&residual);
381  if(verbose) cout << "Iter = " << iter << " Lambda = " << *lambda
382  << " Residual of A*q - lambda*q = " << residual << endl;
383  }
384  if(residual < tolerance) {
385  ierr = 0;
386  break;
387  }
388  }
389  return(ierr);
390 }
391 
392 int check(Epetra_RowMatrix& A, Epetra_RowMatrix & B, bool verbose) {
393 
394  int ierr = 0;
395  EPETRA_TEST_ERR(!A.Comm().NumProc()==B.Comm().NumProc(),ierr);
396  EPETRA_TEST_ERR(!A.Comm().MyPID()==B.Comm().MyPID(),ierr);
397  EPETRA_TEST_ERR(!A.Filled()==B.Filled(),ierr);
398  EPETRA_TEST_ERR(!A.HasNormInf()==B.HasNormInf(),ierr);
400  EPETRA_TEST_ERR(!A.Map().SameAs(B.Map()),ierr);
406  EPETRA_TEST_ERR(!A.NumMyCols()==B.NumMyCols(),ierr);
409  for (int i=0; i<A.NumMyRows(); i++) {
410  int nA, nB;
411  A.NumMyRowEntries(i,nA); B.NumMyRowEntries(i,nB);
412  EPETRA_TEST_ERR(!nA==nB,ierr);
413  }
414  EPETRA_TEST_ERR(!A.NumMyRows()==B.NumMyRows(),ierr);
420  EPETRA_TEST_ERR(!A.UseTranspose()==B.UseTranspose(),ierr);
421 
422  int NumVectors = 5;
423  { // No transpose case
424  Epetra_MultiVector X(A.OperatorDomainMap(), NumVectors);
425  Epetra_MultiVector YA1(A.OperatorRangeMap(), NumVectors);
426  Epetra_MultiVector YA2(YA1);
427  Epetra_MultiVector YB1(YA1);
428  Epetra_MultiVector YB2(YA1);
429  X.Random();
430 
431  bool transA = false;
432  A.SetUseTranspose(transA);
433  B.SetUseTranspose(transA);
434  A.Apply(X,YA1);
435  A.Multiply(transA, X, YA2);
436  EPETRA_TEST_ERR(checkMultiVectors(YA1,YA2,"A Multiply and A Apply", verbose),ierr);
437  B.Apply(X,YB1);
438  EPETRA_TEST_ERR(checkMultiVectors(YA1,YB1,"A Multiply and B Multiply", verbose),ierr);
439  B.Multiply(transA, X, YB2);
440  EPETRA_TEST_ERR(checkMultiVectors(YA1,YB2,"A Multiply and B Apply", verbose), ierr);
441 
442  }
443  {// transpose case
444  Epetra_MultiVector X(A.OperatorRangeMap(), NumVectors);
445  Epetra_MultiVector YA1(A.OperatorDomainMap(), NumVectors);
446  Epetra_MultiVector YA2(YA1);
447  Epetra_MultiVector YB1(YA1);
448  Epetra_MultiVector YB2(YA1);
449  X.Random();
450 
451  bool transA = true;
452  A.SetUseTranspose(transA);
453  B.SetUseTranspose(transA);
454  A.Apply(X,YA1);
455  A.Multiply(transA, X, YA2);
456  EPETRA_TEST_ERR(checkMultiVectors(YA1,YA2, "A Multiply and A Apply (transpose)", verbose),ierr);
457  B.Apply(X,YB1);
458  EPETRA_TEST_ERR(checkMultiVectors(YA1,YB1, "A Multiply and B Multiply (transpose)", verbose),ierr);
459  B.Multiply(transA, X,YB2);
460  EPETRA_TEST_ERR(checkMultiVectors(YA1,YB2, "A Multiply and B Apply (transpose)", verbose),ierr);
461 
462  }
463 
464  Epetra_Vector diagA(A.RowMatrixRowMap());
465  EPETRA_TEST_ERR(A.ExtractDiagonalCopy(diagA),ierr);
466  Epetra_Vector diagB(B.RowMatrixRowMap());
467  EPETRA_TEST_ERR(B.ExtractDiagonalCopy(diagB),ierr);
468  EPETRA_TEST_ERR(checkMultiVectors(diagA,diagB, "ExtractDiagonalCopy", verbose),ierr);
469 
470  Epetra_Vector rowA(A.RowMatrixRowMap());
471  EPETRA_TEST_ERR(A.InvRowSums(rowA),ierr);
472  Epetra_Vector rowB(B.RowMatrixRowMap());
473  EPETRA_TEST_ERR(B.InvRowSums(rowB),ierr)
474  EPETRA_TEST_ERR(checkMultiVectors(rowA,rowB, "InvRowSums", verbose),ierr);
475 
476  Epetra_Vector colA(A.RowMatrixColMap());
477  EPETRA_TEST_ERR(A.InvColSums(colA),ierr);
478  Epetra_Vector colB(B.RowMatrixColMap());
479  EPETRA_TEST_ERR(B.InvColSums(colB),ierr);
480  EPETRA_TEST_ERR(checkMultiVectors(colA,colB, "InvColSums", verbose),ierr);
481 
482  EPETRA_TEST_ERR(checkValues(A.NormInf(), B.NormInf(), "NormInf before scaling", verbose), ierr);
483  EPETRA_TEST_ERR(checkValues(A.NormOne(), B.NormOne(), "NormOne before scaling", verbose),ierr);
484 
485  EPETRA_TEST_ERR(A.RightScale(colA),ierr);
486  EPETRA_TEST_ERR(B.RightScale(colB),ierr);
487 
488 
489  EPETRA_TEST_ERR(A.LeftScale(rowA),ierr);
490  EPETRA_TEST_ERR(B.LeftScale(rowB),ierr);
491 
492 
493  EPETRA_TEST_ERR(checkValues(A.NormInf(), B.NormInf(), "NormInf after scaling", verbose), ierr);
494  EPETRA_TEST_ERR(checkValues(A.NormOne(), B.NormOne(), "NormOne after scaling", verbose),ierr);
495 
496  vector<double> valuesA(A.MaxNumEntries());
497  vector<int> indicesA(A.MaxNumEntries());
498  vector<double> valuesB(B.MaxNumEntries());
499  vector<int> indicesB(B.MaxNumEntries());
500  return(0);
501  for (int i=0; i<A.NumMyRows(); i++) {
502  int nA, nB;
503  EPETRA_TEST_ERR(A.ExtractMyRowCopy(i, A.MaxNumEntries(), nA, &valuesA[0], &indicesA[0]),ierr);
504  EPETRA_TEST_ERR(B.ExtractMyRowCopy(i, B.MaxNumEntries(), nB, &valuesB[0], &indicesB[0]),ierr);
505  EPETRA_TEST_ERR(!nA==nB,ierr);
506  for (int j=0; j<nA; j++) {
507  double curVal = valuesA[j];
508  int curIndex = indicesA[j];
509  bool notfound = true;
510  int jj = 0;
511  while (notfound && jj< nB) {
512  if (!checkValues(curVal, valuesB[jj])) notfound = false;
513  jj++;
514  }
515  EPETRA_TEST_ERR(notfound, ierr);
516  vector<int>::iterator p = find(indicesB.begin(),indicesB.end(),curIndex); // find curIndex in indicesB
517  EPETRA_TEST_ERR(p==indicesB.end(), ierr);
518  }
519 
520  }
521  if (verbose) cout << "RowMatrix Methods check OK" << endl;
522 
523  return (ierr);
524 }
Epetra_MultiVector: A class for constructing and using dense multi-vectors, vectors and matrices in p...
void SetFlopCounter(const Epetra_Flops &FlopCounter_in)
Set the internal Epetra_Flops() pointer.
bool MyGlobalRow(int GID) const
Returns true of GID is owned by the calling processor, otherwise it returns false.
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
int NumGlobalEntries(long long Row) const
Returns the current number of nonzero entries in specified global row on this processor.
virtual long long NumGlobalDiagonals64() const =0
int Random()
Set multi-vector values to random numbers.
virtual int RightScale(const Epetra_Vector &x)=0
Scales the Epetra_RowMatrix on the right with a Epetra_Vector x.
virtual const Epetra_Map & RowMatrixRowMap() const =0
Returns the Epetra_Map object associated with the rows of this matrix.
virtual int SetUseTranspose(bool UseTranspose)=0
If set true, transpose of this operator will be applied.
virtual double NormOne() const =0
Returns the one norm of the global matrix.
bool SameAs(const Epetra_BlockMap &Map) const
Returns true if this and Map are identical maps.
#define EPETRA_TEST_ERR(a, b)
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
double ElapsedTime(void) const
Epetra_Time elapsed time function.
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.
int MyLength() const
Returns the local vector length on the calling processor of vectors in the multi-vector.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
int power_method(Epetra_CrsMatrix &A, double &lambda, int niters, double tolerance, bool verbose)
virtual const Epetra_Map & OperatorDomainMap() const =0
Returns the Epetra_Map object associated with the domain of this operator.
Epetra_Vector: A class for constructing and using dense vectors on a parallel computer.
virtual int ExtractDiagonalCopy(Epetra_Vector &Diagonal) const =0
Returns a copy of the main diagonal in a user-provided vector.
#define EPETRA_MIN(x, y)
virtual int LeftScale(const Epetra_Vector &x)=0
Scales the Epetra_RowMatrix on the left with a Epetra_Vector x.
int MyPID() const
Return my process ID.
Epetra_MpiComm: The Epetra MPI Communication Class.
std::string Epetra_Version()
bool IndicesAreLocal() const
If matrix indices has been transformed to local, this query returns true, otherwise it returns false...
int check(Epetra_CrsGraph &A, int NumMyRows1, int NumGlobalRows1, int NumMyNonzeros1, int NumGlobalNonzeros1, int *MyGlobalElements, bool verbose)
virtual int InvRowSums(Epetra_Vector &x) const =0
Computes the sum of absolute values of the rows of the Epetra_RowMatrix, results returned in x...
int NumVectors() const
Returns the number of vectors in the multi-vector.
Epetra_JadMatrix: A class for constructing matrix objects optimized for common kernels.
virtual int MyPID() const =0
Return my process ID.
int FillComplete(bool OptimizeDataStorage=true)
Signal that data entry is complete. Perform transformations to local index space. ...
Epetra_Time: The Epetra Timing Class.
Definition: Epetra_Time.h:75
int Dot(const Epetra_MultiVector &A, double *Result) const
Computes dot product of each corresponding pair of vectors.
virtual bool LowerTriangular() const =0
If matrix is lower triangular in local index space, this query returns true, otherwise it returns fal...
virtual int NumMyCols() const =0
Returns the number of matrix columns owned by the calling processor.
virtual int MaxAll(double *PartialMaxs, double *GlobalMaxs, int Count) const =0
Epetra_Comm Global Max function.
int NumMyElements() const
Number of elements on the calling processor.
virtual int NumMyRowEntries(int MyRow, int &NumEntries) const =0
Returns the number of nonzero entries in MyRow.
virtual int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const =0
Returns the result of a Epetra_Operator applied to a Epetra_MultiVector X in Y.
virtual int MaxNumEntries() const =0
Returns the maximum of NumMyRowEntries() over all rows.
int powerMethodTests(Epetra_RowMatrix &A, Epetra_RowMatrix &JadA, Epetra_Map &Map, Epetra_Vector &q, Epetra_Vector &z, Epetra_Vector &resid, bool verbose)
virtual const Epetra_Map & OperatorRangeMap() const =0
Returns the Epetra_Map object associated with the range of this operator.
int Scale(double ScalarValue)
Scale the current values of a multi-vector, this = ScalarValue*this.
virtual const Epetra_Comm & Comm() const =0
Returns a pointer to the Epetra_Comm communicator associated with this operator.
virtual bool UseTranspose() const =0
Returns the current UseTranspose setting.
int ExtractGlobalRowCopy(int GlobalRow, int Length, int &NumEntries, double *Values, int *Indices) const
Returns a copy of the specified global row in user-provided arrays.
int OptimizeStorage()
Make consecutive row index sections contiguous, minimize internal storage used for constructing graph...
virtual const Epetra_BlockMap & Map() const =0
Returns a reference to the Epetra_BlockMap for this object.
virtual double NormInf() const =0
Returns the infinity norm of the global matrix.
virtual int NumMyRows() const =0
Returns the number of matrix rows owned by the calling processor.
virtual bool Filled() const =0
If FillComplete() has been called, this query returns true, otherwise it returns false.
int NumProc() const
Returns total number of processes (always returns 1 for SerialComm).
int UpdateValues(const Epetra_RowMatrix &Matrix, bool CheckStructure=false)
Update values using a matrix with identical structure.
virtual long long NumGlobalCols64() const =0
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
Epetra_SerialComm: The Epetra Serial Communication Class.
bool IndicesAreGlobal() const
If matrix indices has not been transformed to local, this query returns true, otherwise it returns fa...
Epetra_Flops: The Epetra Floating Point Operations Class.
Definition: Epetra_Flops.h:58
virtual long long NumGlobalNonzeros64() const =0
virtual int Multiply(bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const =0
Returns the result of a Epetra_RowMatrix multiplied by a Epetra_MultiVector X in Y.
void ResetFlops() const
Resets the number of floating point operations to zero for this multi-vector.
int Update(double ScalarA, const Epetra_MultiVector &A, double ScalarThis)
Update multi-vector values with scaled values of A, this = ScalarThis*this + ScalarA*A.
virtual int NumMyDiagonals() const =0
Returns the number of local nonzero diagonal entries, based on global row/column index comparisons...
virtual int NumProc() const =0
Returns total number of processes.
virtual bool HasNormInf() const =0
Returns true if the this object can provide an approximate Inf-norm, false otherwise.
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
int Norm2(double *Result) const
Compute 2-norm of each vector in multi-vector.
double Flops() const
Returns the number of floating point operations with this multi-vector.
virtual bool UpperTriangular() const =0
If matrix is upper triangular in local index space, this query returns true, otherwise it returns fal...
int main(int argc, char *argv[])
virtual int ExtractMyRowCopy(int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const =0
Returns a copy of the specified local row in user-provided arrays.
virtual const Epetra_Map & RowMatrixColMap() const =0
Returns the Epetra_Map object associated with the columns of this matrix.
virtual int InvColSums(Epetra_Vector &x) const =0
Computes the sum of absolute values of the columns of the Epetra_RowMatrix, results returned in x...
int Broadcast(double *MyVals, int Count, int Root) const
Epetra_SerialComm Broadcast function.
Epetra_RowMatrix: A pure virtual class for using real-valued double-precision row matrices...
int checkValues(double x, double y, string message="", bool verbose=false)
int checkMultiVectors(Epetra_MultiVector &X, Epetra_MultiVector &Y, string message="", bool verbose=false)
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 NumMyNonzeros() const =0
Returns the number of nonzero entries in the calling processor&#39;s portion of the matrix.
virtual long long NumGlobalRows64() const =0
const Epetra_BlockMap & Map() const
Returns the address of the Epetra_BlockMap for this multi-vector.