Amesos Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Amesos_Umfpack.cpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Amesos: Direct Sparse Solver Package
5// Copyright (2004) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// This library is free software; you can redistribute it and/or modify
11// it under the terms of the GNU Lesser General Public License as
12// published by the Free Software Foundation; either version 2.1 of the
13// License, or (at your option) any later version.
14//
15// This library is distributed in the hope that it will be useful, but
16// WITHOUT ANY WARRANTY; without even the implied warranty of
17// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18// Lesser General Public License for more details.
19//
20// You should have received a copy of the GNU Lesser General Public
21// License along with this library; if not, write to the Free Software
22// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23// USA
24// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
25//
26// ***********************************************************************
27// @HEADER
28
29#include "Amesos_Umfpack.h"
30extern "C" {
31#include "umfpack.h"
32}
33#include "Epetra_Map.h"
34#include "Epetra_Import.h"
35#include "Epetra_Export.h"
36#include "Epetra_RowMatrix.h"
37#include "Epetra_CrsMatrix.h"
38#include "Epetra_Vector.h"
39#include "Epetra_Util.h"
40
41//
42// Hack to deal with Bug #1418 - circular dependencies in amesos, umfpack and amd libraries
43//
44#include "Amesos_Klu.h"
45extern "C" {
46#include "amd.h"
47}
48
49//=============================================================================
50Amesos_Umfpack::Amesos_Umfpack(const Epetra_LinearProblem &prob ) :
51 Symbolic(0),
52 Numeric(0),
53 SerialMatrix_(0),
54 UseTranspose_(false),
55 Problem_(&prob),
56 Rcond_(0.0),
57 RcondValidOnAllProcs_(true),
58 MtxConvTime_(-1),
59 MtxRedistTime_(-1),
60 VecRedistTime_(-1),
61 SymFactTime_(-1),
62 NumFactTime_(-1),
63 SolveTime_(-1),
64 OverheadTime_(-1)
65{
66
67 // MS // move declaration of Problem_ above because I need it
68 // MS // set up before calling Comm()
69 Teuchos::ParameterList ParamList ;
70 SetParameters( ParamList ) ;
71
72 //
73 // Hack to deal with Bug #1418 - circular dependencies in amesos, umfpack and amd libraries
74 // This causes the amd files to be pulled in from libamesos.a
75 //
76 if ( UseTranspose_ ) {
77 double control[3];
78 // This should never be called
79 Amesos_Klu Nothing(*Problem_);
80 amd_defaults(control);
81 }
82
83}
84
85//=============================================================================
87{
88 if (Symbolic) umfpack_di_free_symbolic (&Symbolic);
89 if (Numeric) umfpack_di_free_numeric (&Numeric);
90
91 // print out some information if required by the user
92 if ((verbose_ && PrintTiming_) || verbose_ == 2) PrintTiming();
93 if ((verbose_ && PrintStatus_) || verbose_ == 2) PrintStatus();
94
95}
96
97//=============================================================================
98// If FirstTime is true, then build SerialMap and ImportToSerial,
99// otherwise simply re-ship the matrix, so that the numerical values
100// are updated.
101int Amesos_Umfpack::ConvertToSerial(const bool FirstTime)
102{
103 ResetTimer(0);
104 ResetTimer(1);
105
106 const Epetra_Map &OriginalMap = Matrix()->RowMatrixRowMap() ;
107
108 NumGlobalElements_ = Matrix()->NumGlobalRows();
109 numentries_ = Matrix()->NumGlobalNonzeros();
110 assert (NumGlobalElements_ == Matrix()->NumGlobalCols());
111
112 int NumMyElements_ = 0 ;
113 if (MyPID_ == 0) NumMyElements_ = NumGlobalElements_;
114
115 IsLocal_ = ( OriginalMap.NumMyElements() ==
116 OriginalMap.NumGlobalElements() )?1:0;
117
118 // if ( AddZeroToDiag_ ) IsLocal_ = 0 ; // bug # Umfpack does not support AddZeroToDiag_
119
120 Comm().Broadcast( &IsLocal_, 1, 0 ) ;
121
122 // Convert Original Matrix to Serial (if it is not already)
123 //
124 if (IsLocal_== 1) {
126 }
127 else
128 {
129 if (FirstTime)
130 {
131 SerialMap_ = rcp(new Epetra_Map(NumGlobalElements_,NumMyElements_,
132 0,Comm()));
133
134 if (SerialMap_.get() == 0)
135 AMESOS_CHK_ERR(-1);
136
137 ImportToSerial_ = rcp(new Epetra_Import (SerialMap(),OriginalMap));
138
139 if (ImportToSerial_.get() == 0)
140 AMESOS_CHK_ERR(-1);
141 }
142
143 SerialCrsMatrixA_ = rcp(new Epetra_CrsMatrix(Copy,SerialMap(),0));
144
145 if (SerialCrsMatrixA_.get() == 0)
146 AMESOS_CHK_ERR(-1);
147
148 SerialCrsMatrix().Import(*Matrix(), Importer(),Insert);
149
150#if 0
151
152 I was not able to make this work - 11 Feb 2006
153
154 if (AddZeroToDiag_ ) {
155 int OriginalTracebackMode = SerialCrsMatrix().GetTracebackMode() ;
156 SerialCrsMatrix().SetTracebackMode( EPETRA_MIN( OriginalTracebackMode, 0) ) ; // ExportToSerial is called both by PerformSymbolicFactorization() and PerformNumericFactorization(). When called by the latter, the call to insertglobalvalues is both unnecessary and illegal. Fortunately, Epetra allows us to ignore the error message by setting the traceback mode to 0.
157
158 //
159 // Add 0.0 to each diagonal entry to avoid empty diagonal entries;
160 //
161 double zero = 0.0;
162 for ( int i = 0 ; i < SerialMap_->NumGlobalElements(); i++ )
163 if ( SerialCrsMatrix().LRID(i) >= 0 )
164 SerialCrsMatrix().InsertGlobalValues( i, 1, &zero, &i ) ;
165 SerialCrsMatrix().SetTracebackMode( OriginalTracebackMode ) ;
166 }
167#endif
168 SerialCrsMatrix().FillComplete();
170 assert( numentries_ == SerialMatrix_->NumGlobalNonzeros()); // This should be set to an assignment if AddToZeroDiag is non -zero
171 }
172
173
174 MtxRedistTime_ = AddTime("Total matrix redistribution time", MtxRedistTime_, 0);
175 OverheadTime_ = AddTime("Total Amesos overhead time", OverheadTime_, 1);
176
177 return(0);
178}
179
180//=============================================================================
182{
183 ResetTimer(0);
184 ResetTimer(1);
185
186 // Convert matrix to the form that Umfpack expects (Ap, Ai, Aval),
187 // only on processor 0. The matrix has already been assembled in
188 // SerialMatrix_; if only one processor is used, then SerialMatrix_
189 // points to the problem's matrix.
190
191 if (MyPID_ == 0)
192 {
193 Ap.resize( NumGlobalElements_+1 );
194 Ai.resize( EPETRA_MAX( NumGlobalElements_, numentries_) ) ;
195 Aval.resize( EPETRA_MAX( NumGlobalElements_, numentries_) ) ;
196
197 int NumEntries = SerialMatrix_->MaxNumEntries();
198
199 int NumEntriesThisRow;
200 int Ai_index = 0 ;
201 int MyRow;
202 for (MyRow = 0 ; MyRow < NumGlobalElements_; MyRow++)
203 {
204 int ierr;
205 Ap[MyRow] = Ai_index ;
206 ierr = SerialMatrix_->ExtractMyRowCopy(MyRow, NumEntries,
207 NumEntriesThisRow,
208 &Aval[Ai_index], &Ai[Ai_index]);
209 if (ierr)
210 AMESOS_CHK_ERR(-1);
211
212#if 1
213 // MS // added on 15-Mar-05 and KSS restored 8-Feb-06
214 if (AddToDiag_ != 0.0) {
215 for (int i = 0 ; i < NumEntriesThisRow ; ++i) {
216 if (Ai[Ai_index+i] == MyRow) {
217 Aval[Ai_index+i] += AddToDiag_;
218 break;
219 }
220 }
221 }
222#endif
223 Ai_index += NumEntriesThisRow;
224 }
225
226 Ap[MyRow] = Ai_index ;
227 }
228
229 MtxConvTime_ = AddTime("Total matrix conversion time", MtxConvTime_, 0);
230 OverheadTime_ = AddTime("Total Amesos overhead time", OverheadTime_, 1);
231
232 return 0;
233}
234
235//=============================================================================
236int Amesos_Umfpack::SetParameters( Teuchos::ParameterList &ParameterList )
237{
238 // ========================================= //
239 // retrive UMFPACK's parameters from list. //
240 // default values defined in the constructor //
241 // ========================================= //
242
243 // retrive general parameters
244 SetStatusParameters( ParameterList ) ;
245 SetControlParameters( ParameterList ) ;
246
247 return 0;
248}
249
250//=============================================================================
252{
253 // MS // no overhead time in this method
254 ResetTimer(0);
255
256 int symbolic_ok = 0;
257
258 double *Control = (double *) NULL, *Info = (double *) NULL;
259
260 if (Symbolic)
261 umfpack_di_free_symbolic (&Symbolic) ;
262 if (MyPID_== 0) {
263 int status = umfpack_di_symbolic (NumGlobalElements_, NumGlobalElements_, &Ap[0],
264 &Ai[0], &Aval[0],
265 &Symbolic, Control, Info) ;
266 symbolic_ok = (status == UMFPACK_OK);
267 }
268
269 // Communicate the state of the numeric factorization with everyone.
270 Comm().Broadcast(&symbolic_ok, 1, 0);
271
272 if (!symbolic_ok) {
273 if (MyPID_ == 0) {
275 } else {
277 }
278 }
279
280 SymFactTime_ = AddTime("Total symbolic factorization time", SymFactTime_, 0);
281
282 return 0;
283}
284
285//=============================================================================
287{
288 int numeric_ok = 0;
289
290 // MS // no overhead time in this method
291 ResetTimer(0);
292
293 RcondValidOnAllProcs_ = false ;
294 if (MyPID_ == 0) {
295 std::vector<double> Control(UMFPACK_CONTROL);
296 std::vector<double> Info(UMFPACK_INFO);
297 umfpack_di_defaults( &Control[0] ) ;
298 if (Numeric) umfpack_di_free_numeric (&Numeric) ;
299 int status = umfpack_di_numeric (&Ap[0],
300 &Ai[0],
301 &Aval[0],
302 Symbolic,
303 &Numeric,
304 &Control[0],
305 &Info[0]) ;
306 Rcond_ = Info[UMFPACK_RCOND];
307
308#if NOT_DEF
309 std::cout << " Rcond_ = " << Rcond_ << std::endl ;
310
311 int lnz1 = 1000 ;
312 int unz1 = 1000 ;
313 int n = 4;
314 int * Lp = (int *) malloc ((n+1) * sizeof (int)) ;
315 int * Lj = (int *) malloc (lnz1 * sizeof (int)) ;
316 double * Lx = (double *) malloc (lnz1 * sizeof (double)) ;
317 int * Up = (int *) malloc ((n+1) * sizeof (int)) ;
318 int * Ui = (int *) malloc (unz1 * sizeof (int)) ;
319 double * Ux = (double *) malloc (unz1 * sizeof (double)) ;
320 int * P = (int *) malloc (n * sizeof (int)) ;
321 int * Q = (int *) malloc (n * sizeof (int)) ;
322 double * Dx = (double *) NULL ; /* D vector not requested */
323 double * Rs = (double *) malloc (n * sizeof (double)) ;
324 if (!Lp || !Lj || !Lx || !Up || !Ui || !Ux || !P || !Q || !Rs)
325 {
326 assert( false ) ;
327 }
328 int do_recip;
329 status = umfpack_di_get_numeric (Lp, Lj, Lx, Up, Ui, Ux,
330 P, Q, Dx, &do_recip, Rs, Numeric) ;
331 if (status < 0)
332 {
333 assert( false ) ;
334 }
335
336 printf ("\nL (lower triangular factor of C): ") ;
337 (void) umfpack_di_report_matrix (n, n, Lp, Lj, Lx, 0, &Control[0]) ;
338 printf ("\nU (upper triangular factor of C): ") ;
339 (void) umfpack_di_report_matrix (n, n, Up, Ui, Ux, 1, &Control[0]) ;
340 printf ("\nP: ") ;
341 (void) umfpack_di_report_perm (n, P, &Control[0]) ;
342 printf ("\nQ: ") ;
343 (void) umfpack_di_report_perm (n, Q, &Control[0]) ;
344 printf ("\nScale factors: row i of A is to be ") ;
345
346#endif
347
348 numeric_ok = (status == UMFPACK_OK);
349 }
350
351 // Communicate the state of the numeric factorization with everyone.
352 Comm().Broadcast(&numeric_ok, 1, 0);
353
354 if (!numeric_ok) {
355 if (MyPID_ == 0) {
357 } else {
359 }
360 }
361
362 NumFactTime_ = AddTime("Total numeric factorization time", NumFactTime_, 0);
363
364 return 0;
365}
366
367//=============================================================================
369{
370 if ( !RcondValidOnAllProcs_ ) {
371 Comm().Broadcast( &Rcond_, 1, 0 ) ;
372 RcondValidOnAllProcs_ = true;
373 }
374 return(Rcond_);
375}
376
377//=============================================================================
379{
380 bool OK = true;
381
382 if ( GetProblem()->GetOperator()->OperatorRangeMap().NumGlobalPoints() !=
383 GetProblem()->GetOperator()->OperatorDomainMap().NumGlobalPoints() ) OK = false;
384 return OK;
385}
386
387
388//=============================================================================
390{
391 // MS // NOTE: If you change this method, also change
392 // MS // NumericFactorization(), because it performs part of the actions
393 // MS // of this method. This is to avoid to ship the matrix twice
394 // MS // (once for the symbolic factorization, and once for the numeric
395 // MS // factorization) when it is not necessary.
396
399
400 CreateTimer(Comm(), 2);
401 // MS // Initialize two timers:
402 // MS // timer 1: this will track all time spent in Amesos most
403 // MS // important functions, *including* UMFPACK functions
404 // MS // timer 2: this will track all time spent in this function
405 // MS // that is not due to UMFPACK calls, and therefore it
406 // MS // tracks down how much Amesos costs. The timer starts
407 // MS // and ends in *each* method, unless the method does not
408 // MS // perform any real operation. If a method calls another
409 // MS // method, the timer will be stopped before the called
410 // MS // method, then restared.
411 // MS // All the time of this timer goes into "overhead"
412
413 MyPID_ = Comm().MyPID();
414 NumProcs_ = Comm().NumProc();
415
418
420
422
424
425 return 0;
426}
427
428//=============================================================================
430{
432
433 if (IsSymbolicFactorizationOK_ == false)
434 {
435 // Call here what is needed, to avoid double shipping of the matrix
436 CreateTimer(Comm(), 2);
437
438 MyPID_ = Comm().MyPID();
439 NumProcs_ = Comm().NumProc();
440
443
445
447
449
451 }
452 else
453 {
454 // need to reshuffle and reconvert because entry values may have changed
457
459 }
460
462
464
465 return 0;
466}
467
468//=============================================================================
470{
471 // if necessary, perform numeric factorization.
472 // This may call SymbolicFactorization() as well.
475
476 ResetTimer(1);
477
478 Epetra_MultiVector* vecX = Problem_->GetLHS();
479 Epetra_MultiVector* vecB = Problem_->GetRHS();
480
481 if ((vecX == 0) || (vecB == 0))
482 AMESOS_CHK_ERR(-1);
483
484 int NumVectors = vecX->NumVectors();
485 if (NumVectors != vecB->NumVectors())
486 AMESOS_CHK_ERR(-1);
487
488 Epetra_MultiVector *SerialB, *SerialX;
489
490 // Extract Serial versions of X and B
491 //
492 double *SerialXvalues ;
493 double *SerialBvalues ;
494
495 Epetra_MultiVector* SerialXextract = 0;
496 Epetra_MultiVector* SerialBextract = 0;
497
498 // Copy B to the serial version of B
499 //
500 ResetTimer(0);
501
502 if (IsLocal_ == 1) {
503 SerialB = vecB ;
504 SerialX = vecX ;
505 } else {
506 assert (IsLocal_ == 0);
507 SerialXextract = new Epetra_MultiVector(SerialMap(),NumVectors);
508 SerialBextract = new Epetra_MultiVector(SerialMap(),NumVectors);
509
510 SerialBextract->Import(*vecB,Importer(),Insert);
511 SerialB = SerialBextract;
512 SerialX = SerialXextract;
513 }
514
515 VecRedistTime_ = AddTime("Total vector redistribution time", VecRedistTime_, 0);
516
517 // Call UMFPACK to perform the solve
518 // Note: UMFPACK uses a Compressed Column Storage instead of compressed row storage,
519 // Hence to compute A X = B, we ask UMFPACK to perform A^T X = B and vice versa
520
521 OverheadTime_ = AddTime("Total Amesos overhead time", OverheadTime_, 1);
522
523 ResetTimer(0);
524
525 int SerialBlda, SerialXlda ;
526 int UmfpackRequest = UseTranspose()?UMFPACK_A:UMFPACK_At ;
527 int status = 0;
528
529 if ( MyPID_ == 0 ) {
530 int ierr;
531 ierr = SerialB->ExtractView(&SerialBvalues, &SerialBlda);
532 assert (ierr == 0);
533 ierr = SerialX->ExtractView(&SerialXvalues, &SerialXlda);
534 assert (ierr == 0);
535 assert( SerialBlda == NumGlobalElements_ ) ;
536 assert( SerialXlda == NumGlobalElements_ ) ;
537
538 for ( int j =0 ; j < NumVectors; j++ ) {
539 double *Control = (double *) NULL, *Info = (double *) NULL ;
540
541 status = umfpack_di_solve (UmfpackRequest, &Ap[0],
542 &Ai[0], &Aval[0],
543 &SerialXvalues[j*SerialXlda],
544 &SerialBvalues[j*SerialBlda],
545 Numeric, Control, Info) ;
546 }
547 }
548
549 if (status) AMESOS_CHK_ERR(status);
550
551 SolveTime_ = AddTime("Total solve time", SolveTime_, 0);
552
553 // Copy X back to the original vector
554
555 ResetTimer(0);
556 ResetTimer(1);
557
558 if ( IsLocal_ == 0 ) {
559 vecX->Export(*SerialX, Importer(), Insert ) ;
560 if (SerialBextract) delete SerialBextract ;
561 if (SerialXextract) delete SerialXextract ;
562 }
563
564 VecRedistTime_ = AddTime("Total vector redistribution time", VecRedistTime_, 0);
565
567 {
568 Epetra_RowMatrix* Matrix =
569 dynamic_cast<Epetra_RowMatrix*>(Problem_->GetOperator());
570 ComputeTrueResidual(*Matrix, *vecX, *vecB, UseTranspose(), "Amesos_Umfpack");
571 }
572
574 ComputeVectorNorms(*vecX, *vecB, "Amesos_Umfpack");
575 }
576
577 NumSolve_++;
578
579 OverheadTime_ = AddTime("Total Amesos overhead time", OverheadTime_, 1); // Amesos overhead
580
581 return(0);
582}
583
584// ======================================================================
586{
587 if (MyPID_ != 0) return;
588
589 PrintLine();
590
591 std::cout << "Amesos_Umfpack : Matrix has " << NumGlobalElements_ << " rows"
592 << " and " << numentries_ << " nonzeros" << std::endl;
593 std::cout << "Amesos_Umfpack : Nonzero elements per row = "
594 << 1.0*numentries_/NumGlobalElements_ << std::endl;
595 std::cout << "Amesos_Umfpack : Percentage of nonzero elements = "
596 << 100.0*numentries_/(pow(double(NumGlobalElements_),double(2.0))) << std::endl;
597 std::cout << "Amesos_Umfpack : Use transpose = " << UseTranspose_ << std::endl;
598
599 PrintLine();
600
601 return;
602}
603
604// ======================================================================
606{
607 if (Problem_->GetOperator() == 0 || MyPID_ != 0)
608 return;
609
610 double ConTime = GetTime(MtxConvTime_);
611 double MatTime = GetTime(MtxRedistTime_);
612 double VecTime = GetTime(VecRedistTime_);
613 double SymTime = GetTime(SymFactTime_);
614 double NumTime = GetTime(NumFactTime_);
615 double SolTime = GetTime(SolveTime_);
616 double OveTime = GetTime(OverheadTime_);
617
619 SymTime /= NumSymbolicFact_;
620
621 if (NumNumericFact_)
622 NumTime /= NumNumericFact_;
623
624 if (NumSolve_)
625 SolTime /= NumSolve_;
626
627 std::string p = "Amesos_Umfpack : ";
628 PrintLine();
629
630 std::cout << p << "Time to convert matrix to Umfpack format = "
631 << ConTime << " (s)" << std::endl;
632 std::cout << p << "Time to redistribute matrix = "
633 << MatTime << " (s)" << std::endl;
634 std::cout << p << "Time to redistribute vectors = "
635 << VecTime << " (s)" << std::endl;
636 std::cout << p << "Number of symbolic factorizations = "
637 << NumSymbolicFact_ << std::endl;
638 std::cout << p << "Time for sym fact = "
639 << SymTime * NumSymbolicFact_ << " (s), avg = "
640 << SymTime << " (s)" << std::endl;
641 std::cout << p << "Number of numeric factorizations = "
642 << NumNumericFact_ << std::endl;
643 std::cout << p << "Time for num fact = "
644 << NumTime * NumNumericFact_ << " (s), avg = "
645 << NumTime << " (s)" << std::endl;
646 std::cout << p << "Number of solve phases = "
647 << NumSolve_ << std::endl;
648 std::cout << p << "Time for solve = "
649 << SolTime * NumSolve_ << " (s), avg = " << SolTime << " (s)" << std::endl;
650 double tt = SymTime * NumSymbolicFact_ + NumTime * NumNumericFact_ + SolTime * NumSolve_;
651 if (tt != 0)
652 {
653 std::cout << p << "Total time spent in Amesos = " << tt << " (s) " << std::endl;
654 std::cout << p << "Total time spent in the Amesos interface = " << OveTime << " (s)" << std::endl;
655 std::cout << p << "(the above time does not include UMFPACK time)" << std::endl;
656 std::cout << p << "Amesos interface time / total time = " << OveTime / tt << std::endl;
657 }
658
659 PrintLine();
660
661 return;
662}
const int NumericallySingularMatrixError
const int StructurallySingularMatrixError
#define AMESOS_CHK_ERR(a)
void SetControlParameters(const Teuchos::ParameterList &ParameterList)
bool AddZeroToDiag_
Adds zero to diagonal of redistributed matrix (some solvers choke on a matrix with a partly empty dia...
double AddToDiag_
Add this value to the diagonal.
Amesos_Klu: A serial, unblocked code ideal for getting started and for very sparse matrices,...
Definition Amesos_Klu.h:116
bool PrintTiming_
If true, prints timing information in the destructor.
bool ComputeVectorNorms_
If true, prints the norms of X and B in Solve().
int verbose_
Toggles the output level.
int NumSymbolicFact_
Number of symbolic factorization phases.
bool IsNumericFactorizationOK_
If true, NumericFactorization() has been successfully called.
int NumSolve_
Number of solves.
bool ComputeTrueResidual_
If true, computes the true residual in Solve().
bool PrintStatus_
If true, print additional information in the destructor.
void SetStatusParameters(const Teuchos::ParameterList &ParameterList)
bool IsSymbolicFactorizationOK_
If true, SymbolicFactorization() has been successfully called.
int NumNumericFact_
Number of numeric factorization phases.
double GetTime(const std::string what) const
Gets the cumulative time using the string.
void ResetTimer(const int timerID=0)
Resets the internally stored time object.
Definition Amesos_Time.h:74
int AddTime(const std::string what, int dataID, const int timerID=0)
Adds to field what the time elapsed since last call to ResetTimer().
Definition Amesos_Time.h:80
void CreateTimer(const Epetra_Comm &Comm, int size=1)
Initializes the Time object.
Definition Amesos_Time.h:64
bool UseTranspose_
If true, solve the problem with the transpose.
int IsLocal_
1 if Problem_->GetOperator() is stored entirely on process 0
bool UseTranspose() const
Returns the current UseTranspose setting.
const Epetra_Comm & Comm() const
Returns a pointer to the Epetra_Comm communicator associated with this operator.
Teuchos::RCP< Epetra_Import > ImportToSerial_
Importer from distributed to serial (all rows on process 0).
int PerformNumericFactorization()
void PrintTiming() const
Prints timing information.
int NumericFactorization()
Performs NumericFactorization on the matrix A.
bool MatrixShapeOK() const
Returns true if UMFPACK can handle this matrix shape.
const Epetra_LinearProblem * Problem_
Pointer to the linear problem to solve.
void * Numeric
Umfpack internal opaque object.
Amesos_Umfpack(const Epetra_LinearProblem &LinearProblem)
Amesos_Umfpack Constructor.
int PerformSymbolicFactorization()
const Epetra_CrsMatrix & SerialCrsMatrix() const
Teuchos::RCP< Epetra_CrsMatrix > SerialCrsMatrixA_
std::vector< double > Aval
const Epetra_Map & SerialMap() const
void PrintStatus() const
Prints information about the factorization and solution phases.
const Epetra_Import & Importer() const
~Amesos_Umfpack(void)
Amesos_Umfpack Destructor.
int numentries_
Number of non-zero entries in Problem_->GetOperator()
int SymbolicFactorization()
Performs SymbolicFactorization on the matrix A.
void * Symbolic
Umfpack internal opaque object.
std::vector< int > Ai
std::vector< int > Ap
Ap, Ai, Aval form the compressed row storage used by Umfpack.
int Solve()
Solves A X = B (or AT x = B)
Epetra_RowMatrix * Matrix()
Returns a pointer to the linear system matrix.
int NumGlobalElements_
Number of rows and columns in the Problem_->GetOperator()
Teuchos::RCP< Epetra_Map > SerialMap_
Points to a Serial Map (unused if IsLocal == 1 )
const Epetra_LinearProblem * GetProblem() const
Returns the Epetra_LinearProblem.
int ConvertToSerial(const bool FirstTime)
Converts matrix to a serial Epetra_CrsMatrix.
double Rcond_
Reciprocal condition number estimate.
double GetRcond() const
Returns an estimate of the reciprocal of the condition number.
int SetParameters(Teuchos::ParameterList &ParameterList)
Updates internal variables.
int MtxConvTime_
Quick access pointers to internal timer data.
Epetra_RowMatrix * SerialMatrix_
Points to a Serial Copy of A.
void ComputeTrueResidual(const Epetra_RowMatrix &Matrix, const Epetra_MultiVector &X, const Epetra_MultiVector &B, const bool UseTranspose, const std::string prefix) const
Computes the true residual, B - Matrix * X, and prints the results.
void ComputeVectorNorms(const Epetra_MultiVector &X, const Epetra_MultiVector &B, const std::string prefix) const
Computes the norms of X and B and print the results.
void PrintLine() const
Prints line on std::cout.