Amesos Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Private Attributes | List of all members
Amesos_Klu Class Reference

Amesos_Klu: A serial, unblocked code ideal for getting started and for very sparse matrices, such as circuit matrces. More...

#include <Amesos_Klu.h>

Inheritance diagram for Amesos_Klu:
Inheritance graph
[legend]

Private Attributes

int SerialXlda_
 
Teuchos::RCP< Amesos_Klu_PimplPrivateKluData_
 
Teuchos::RCP
< Amesos_StandardIndex
StdIndex_
 
Teuchos::RCP
< Amesos_StandardIndex
StdIndexRange_
 
Teuchos::RCP
< Amesos_StandardIndex
StdIndexDomain_
 
std::vector< int > Ap
 Ap, Ai, Aval form the compressed row storage used by Klu Ai and Aval can point directly into a matrix if it is StorageOptimized(), hence they may either be in vector form or may be a pointer into Epetra_CrsMatrix internals. More...
 
std::vector< int > VecAi
 
std::vector< double > VecAval
 
double * Aval
 
int * Ai
 
int UseDataInPlace_
 1 if Problem_->GetOperator() is stored entirely on process 0 More...
 
long long numentries_
 Number of non-zero entries in Problem_->GetOperator() More...
 
long long NumGlobalElements_
 Number of rows and columns in the Problem_->GetOperator() More...
 
Epetra_RowMatrixRowMatrixA_
 Operator converted to a RowMatrix. More...
 
Epetra_CrsMatrixCrsMatrixA_
 Operator converted to a CrsMatrix. More...
 
Teuchos::RCP< Epetra_MapSerialMap_
 Points to a Serial Map (unused if UseDataInPlace_ == 1 ) More...
 
Teuchos::RCP< Epetra_CrsMatrixSerialCrsMatrixA_
 Points to a Serial Copy of A (unused if UseDataInPlace_==1) More...
 
Epetra_RowMatrixStdIndexMatrix_
 Points to a Contiguous Copy of A. More...
 
Epetra_RowMatrixSerialMatrix_
 Points to a Serial Copy of A. More...
 
bool TrustMe_
 If true, no checks are made and the matrix is assume to be distributed. More...
 
int NumVectors_
 Number of vectors in RHS and LHS. More...
 
double * SerialXBvalues_
 Pointer to the actual values in the serial version of X and B. More...
 
double * SerialBvalues_
 
Teuchos::RCP< Epetra_MultiVectorSerialB_
 Serial versions of the LHS and RHS (may point to the original vector if serial) More...
 
Teuchos::RCP< Epetra_MultiVectorSerialX_
 
Teuchos::RCP< Epetra_MultiVectorSerialXextract_
 Serial versions of the LHS and RHS (if necessary) More...
 
Teuchos::RCP< Epetra_MultiVectorSerialBextract_
 
bool UseTranspose_
 If true, the transpose of A is used. More...
 
const Epetra_LinearProblemProblem_
 Pointer to the linear system problem. More...
 
std::vector< int > ColIndicesV_
 Only used for RowMatrices to extract copies. More...
 
std::vector< double > RowValuesV_
 Only used for RowMatrices to extract copies. More...
 
Teuchos::RCP< Epetra_ImportImportToSerial_
 Importer to process 0. More...
 
Teuchos::RCP< Epetra_ImportImportRangeToSerial_
 
Teuchos::RCP< Epetra_ImportImportDomainToSerial_
 
int MtxRedistTime_
 Quick access ids for the individual timings. More...
 
int MtxConvTime_
 
int VecRedistTime_
 
int SymFactTime_
 
int NumFactTime_
 
int SolveTime_
 
int OverheadTime_
 
- Private Attributes inherited from Amesos_Control
double AddToDiag_
 Add this value to the diagonal. More...
 
bool refactorize_
 
double rcond_threshold_
 If error is greater than this value, perform symbolic and numeric factorization with full partial pivoting. More...
 
int ScaleMethod_
 
bool AddZeroToDiag_
 Adds zero to diagonal of redistributed matrix (some solvers choke on a matrix with a partly empty diag) More...
 
int MatrixProperty_
 Set the matrix property. More...
 
int MaxProcesses_
 
bool Reindex_
 If true, the Amesos class should reindex the matrix to standard indexing (i.e. More...
 
- Private Attributes inherited from Amesos_Status
bool IsSymbolicFactorizationOK_
 If true, SymbolicFactorization() has been successfully called. More...
 
bool IsNumericFactorizationOK_
 If true, NumericFactorization() has been successfully called. More...
 
bool PrintTiming_
 If true, prints timing information in the destructor. More...
 
bool PrintStatus_
 If true, print additional information in the destructor. More...
 
bool ComputeVectorNorms_
 If true, prints the norms of X and B in Solve(). More...
 
bool ComputeTrueResidual_
 If true, computes the true residual in Solve(). More...
 
int verbose_
 Toggles the output level. More...
 
int debug_
 Sets the level of debug_ output. More...
 
int NumSymbolicFact_
 Number of symbolic factorization phases. More...
 
int NumNumericFact_
 Number of numeric factorization phases. More...
 
int NumSolve_
 Number of solves. More...
 
double Threshold_
 
int MyPID_
 
int NumProcs_
 
 Amesos_Klu (const Epetra_LinearProblem &LinearProblem)
 Amesos_Klu Constructor. More...
 
 ~Amesos_Klu (void)
 Amesos_Klu Destructor. More...
 
int SymbolicFactorization ()
 Performs SymbolicFactorization on the matrix A. More...
 
int NumericFactorization ()
 Performs NumericFactorization on the matrix A. More...
 
int Solve ()
 Solves A X = B (or AT x = B) More...
 
const Epetra_LinearProblemGetProblem () const
 Get a pointer to the Problem. More...
 
bool MatrixShapeOK () const
 Returns true if KLU can handle this matrix shape. More...
 
int SetUseTranspose (bool UseTranspose_in)
 SetUseTranpose(true) is more efficient in Amesos_Klu. More...
 
bool UseTranspose () const
 Returns the current UseTranspose setting. More...
 
const Epetra_CommComm () const
 Returns a pointer to the Epetra_Comm communicator associated with this operator. More...
 
int SetParameters (Teuchos::ParameterList &ParameterList)
 Updates internal variables. More...
 
int NumSymbolicFact () const
 Returns the number of symbolic factorizations performed by this object. More...
 
int NumNumericFact () const
 Returns the number of numeric factorizations performed by this object. More...
 
int NumSolve () const
 Returns the number of solves performed by this object. More...
 
void PrintTiming () const
 Prints timing information. More...
 
void PrintStatus () const
 Prints information about the factorization and solution phases. More...
 
void GetTiming (Teuchos::ParameterList &TimingParameterList) const
 Extracts timing information and places in parameter list. More...
 
int CreateLocalMatrixAndExporters ()
 
int ExportToSerial ()
 
int ConvertToKluCRS (bool firsttime)
 
int PerformSymbolicFactorization ()
 
int PerformNumericFactorization ()
 

Additional Inherited Members

- Public Member Functions inherited from Amesos_BaseSolver
virtual ~Amesos_BaseSolver ()
 Destructor. More...
 
virtual void setParameterList (Teuchos::RCP< Teuchos::ParameterList > const &paramList)
 Redefined from Teuchos::ParameterListAcceptor (Does Not Work) More...
 
virtual Teuchos::RCP
< Teuchos::ParameterList
getNonconstParameterList ()
 This is an empty stub. More...
 
virtual Teuchos::RCP
< Teuchos::ParameterList
unsetParameterList ()
 This is an empty stub. More...
 
- Public Member Functions inherited from Teuchos::ParameterListAcceptor
virtual RCP< const ParameterListgetParameterList () const
 
virtual RCP< const ParameterListgetValidParameters () const
 
- Private Member Functions inherited from Amesos_Time
 Amesos_Time ()
 Default constructor to create size timers. More...
 
virtual ~Amesos_Time ()
 Default destructor. More...
 
void CreateTimer (const Epetra_Comm &Comm, int size=1)
 Initializes the Time object. More...
 
void ResetTimer (const int timerID=0)
 Resets the internally stored time object. More...
 
int AddTime (const std::string what, int dataID, const int timerID=0)
 Adds to field what the time elapsed since last call to ResetTimer(). More...
 
double GetTime (const std::string what) const
 Gets the cumulative time using the string. More...
 
double GetTime (const int dataID) const
 Gets the cumulative time using the dataID. More...
 
void GetTiming (Teuchos::ParameterList &list) const
 Load up the current timing information into the parameter list. More...
 
- Private Member Functions inherited from Amesos_NoCopiable
 Amesos_NoCopiable ()
 Default constructor. More...
 
 ~Amesos_NoCopiable ()
 Default destructor. More...
 
- Private Member Functions inherited from Amesos_Utils
 Amesos_Utils ()
 Default constructor. More...
 
 ~Amesos_Utils ()
 Default destructor. More...
 
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. More...
 
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. More...
 
void PrintLine () const
 Prints line on std::cout. More...
 
void SetMaxProcesses (int &MaxProcesses, const Epetra_RowMatrix &A)
 
- Private Member Functions inherited from Amesos_Control
 Amesos_Control ()
 Default constructor. More...
 
 ~Amesos_Control ()
 Default destructor. More...
 
void SetControlParameters (const Teuchos::ParameterList &ParameterList)
 
- Private Member Functions inherited from Amesos_Status
 Amesos_Status ()
 Default constructor. More...
 
 ~Amesos_Status ()
 Default destructor. More...
 
void SetStatusParameters (const Teuchos::ParameterList &ParameterList)
 

Detailed Description

Amesos_Klu: A serial, unblocked code ideal for getting started and for very sparse matrices, such as circuit matrces.

Interface to UMFPACK.

Interface to KLU internal solver.

Class Amesos_Klu is an object-oriented wrapper for KLU. KLU, whose sources are distributed within Amesos, is a serial solver for sparse matrices. KLU will solve a linear system of equations: $ A X = B$, where A is an Epetra_RowMatrix and X and B are Epetra_MultiVector objects.

Amesos_Klu computes $ A^T X = B$ more efficiently than $ >A X = B$. The latter requires a matrix transpose – which costs both time and space.

KLU is Tim Davis' implementation of Gilbert-Peierl's left-looking sparse partial pivoting algorithm, with Eisenstat & Liu's symmetric pruning. Gilbert's version appears as [L,U,P]=lu(A) in MATLAB. It doesn't exploit dense matrix kernels, but it is the only sparse LU factorization algorithm known to be asymptotically optimal, in the sense that it takes time proportional to the number of floating-point operations. It is the precursor to SuperLU, thus the name ("clark Kent LU"). For very sparse matrices that do not suffer much fill-in (such as most circuit matrices when permuted properly) dense matrix kernels do not help, and the asymptotic run-time is of practical importance.

The klu_btf code first permutes the matrix to upper block triangular form (using two algorithms by Duff and Reid, MC13 and MC21, in the ACM Collected Algorithms). It then permutes each block via a symmetric minimum degree ordering (AMD, by Amestoy, Davis, and Duff). This ordering phase can be done just once for a sequence of matrices. Next, it factorizes each reordered block via the klu routine, which also attempts to preserve diagonal pivoting, but allows for partial pivoting if the diagonal is to small.

Date
Last updated on 24-May-05.

Definition at line 111 of file Amesos_Klu.h.

Constructor & Destructor Documentation

Amesos_Klu::Amesos_Klu ( const Epetra_LinearProblem LinearProblem)

Amesos_Klu Constructor.

Creates an Amesos_Klu instance, using an Epetra_LinearProblem, passing in an already-defined Epetra_LinearProblem object.

Note: The operator in LinearProblem must be an Epetra_RowMatrix.

Definition at line 89 of file Amesos_Klu.cpp.

Amesos_Klu::~Amesos_Klu ( void  )

Amesos_Klu Destructor.

Definition at line 111 of file Amesos_Klu.cpp.

Member Function Documentation

int Amesos_Klu::SymbolicFactorization ( )
virtual

Performs SymbolicFactorization on the matrix A.

In addition to performing symbolic factorization on the matrix A, the call to SymbolicFactorization() implies that no change will be made to the non-zero structure of the underlying matrix without a subsequent call to SymbolicFactorization().

<br >Preconditions:

<br >Postconditions:

Returns
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

Definition at line 617 of file Amesos_Klu.cpp.

int Amesos_Klu::NumericFactorization ( )
virtual

Performs NumericFactorization on the matrix A.

In addition to performing numeric factorization on the matrix A, the call to NumericFactorization() implies that no change will be made to the underlying matrix without a subsequent call to NumericFactorization().

<br >Preconditions:

  • GetProblem().GetOperator() != 0 (return -1)
  • MatrixShapeOk(GetProblem().GetOperator()) == true (return -6)
  • The non-zero structure of the matrix should not have changed since the last call to SymbolicFactorization(). (return -2 if the number of non-zeros changes) Other changes can have arbitrary consequences.
  • The distribution of the matrix should not have changed since the last call to SymbolicFactorization()
  • The matrix should be indexed from 0 to n-1, unless the parameter "Reindex" was set to "true" prior to the call to SymbolicFactorization(). (return -3 - if caught)
  • The paremeter "Reindex" should not be set to "true" except on CrsMatrices. (return -4)
  • The paremeter "Reindex" should not be set to "true" unless Amesos was built with EpetraExt, i.e. with –enable-epetraext on the configure line. (return -4)
  • Internal errors retur -5.

<br >Postconditions:

  • Numeric Factorization will be performed (or marked to be performed) allowing Solve() to be performed correctly despite a potential change in in the matrix values (though not in the non-zero structure).
Returns
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

Definition at line 681 of file Amesos_Klu.cpp.

int Amesos_Klu::Solve ( )
virtual

Solves A X = B (or AT x = B)

<br >Preconditions:

<br >Postconditions:

  • X will be set such that A X = B (or AT X = B), within the limits of the accuracy of the underlying solver.
Returns
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

Definition at line 729 of file Amesos_Klu.cpp.

const Epetra_LinearProblem* Amesos_Klu::GetProblem ( ) const
inlinevirtual

Get a pointer to the Problem.

Implements Amesos_BaseSolver.

Definition at line 147 of file Amesos_Klu.h.

bool Amesos_Klu::MatrixShapeOK ( ) const
virtual

Returns true if KLU can handle this matrix shape.

Returns true if the matrix shape is one that KLU can handle. KLU only works with square matrices.

Implements Amesos_BaseSolver.

Definition at line 603 of file Amesos_Klu.cpp.

int Amesos_Klu::SetUseTranspose ( bool  UseTranspose_in)
inlinevirtual

SetUseTranpose(true) is more efficient in Amesos_Klu.

If SetUseTranspose() is set to true, $ A^T X = B$ is computed.

Implements Amesos_BaseSolver.

Definition at line 160 of file Amesos_Klu.h.

bool Amesos_Klu::UseTranspose ( ) const
inlinevirtual

Returns the current UseTranspose setting.

Implements Amesos_BaseSolver.

Definition at line 162 of file Amesos_Klu.h.

const Epetra_Comm& Amesos_Klu::Comm ( ) const
inlinevirtual

Returns a pointer to the Epetra_Comm communicator associated with this operator.

Implements Amesos_BaseSolver.

Definition at line 164 of file Amesos_Klu.h.

int Amesos_Klu::SetParameters ( Teuchos::ParameterList ParameterList)
virtual

Updates internal variables.

  <br \>Preconditions:<ul>
  <li>None.</li>
  </ul>

  <br \>Postconditions:<ul> 
  <li>Internal variables controlling the factorization and solve will
  be updated and take effect on all subseuent calls to NumericFactorization() 
  and Solve().</li>
  <li>All parameters whose value are to differ from the default values must 

be included in ParameterList. Parameters not specified in ParameterList revert to their default values.

Returns
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

Definition at line 449 of file Amesos_Klu.cpp.

int Amesos_Klu::NumSymbolicFact ( ) const
inlinevirtual

Returns the number of symbolic factorizations performed by this object.

Implements Amesos_BaseSolver.

Definition at line 169 of file Amesos_Klu.h.

int Amesos_Klu::NumNumericFact ( ) const
inlinevirtual

Returns the number of numeric factorizations performed by this object.

Implements Amesos_BaseSolver.

Definition at line 172 of file Amesos_Klu.h.

int Amesos_Klu::NumSolve ( ) const
inlinevirtual

Returns the number of solves performed by this object.

Implements Amesos_BaseSolver.

Definition at line 175 of file Amesos_Klu.h.

void Amesos_Klu::PrintTiming ( ) const
virtual

Prints timing information.

Implements Amesos_BaseSolver.

Definition at line 949 of file Amesos_Klu.cpp.

void Amesos_Klu::PrintStatus ( ) const
virtual

Prints information about the factorization and solution phases.

Implements Amesos_BaseSolver.

Definition at line 926 of file Amesos_Klu.cpp.

void Amesos_Klu::GetTiming ( Teuchos::ParameterList TimingParameterList) const
inlinevirtual

Extracts timing information and places in parameter list.

Reimplemented from Amesos_BaseSolver.

Definition at line 184 of file Amesos_Klu.h.

int Amesos_Klu::CreateLocalMatrixAndExporters ( )
private

Definition at line 196 of file Amesos_Klu.cpp.

int Amesos_Klu::ExportToSerial ( )
private

Definition at line 119 of file Amesos_Klu.cpp.

int Amesos_Klu::ConvertToKluCRS ( bool  firsttime)
private

Definition at line 337 of file Amesos_Klu.cpp.

int Amesos_Klu::PerformSymbolicFactorization ( )
private

Definition at line 478 of file Amesos_Klu.cpp.

int Amesos_Klu::PerformNumericFactorization ( )
private

Definition at line 516 of file Amesos_Klu.cpp.

Member Data Documentation

int Amesos_Klu::SerialXlda_
private

Definition at line 258 of file Amesos_Klu.h.

Teuchos::RCP<Amesos_Klu_Pimpl> Amesos_Klu::PrivateKluData_
private

Definition at line 267 of file Amesos_Klu.h.

Teuchos::RCP<Amesos_StandardIndex> Amesos_Klu::StdIndex_
private

Definition at line 268 of file Amesos_Klu.h.

Teuchos::RCP<Amesos_StandardIndex> Amesos_Klu::StdIndexRange_
private

Definition at line 269 of file Amesos_Klu.h.

Teuchos::RCP<Amesos_StandardIndex> Amesos_Klu::StdIndexDomain_
private

Definition at line 270 of file Amesos_Klu.h.

std::vector<int> Amesos_Klu::Ap
private

Ap, Ai, Aval form the compressed row storage used by Klu Ai and Aval can point directly into a matrix if it is StorageOptimized(), hence they may either be in vector form or may be a pointer into Epetra_CrsMatrix internals.

Ap must always be constructed.

Definition at line 276 of file Amesos_Klu.h.

std::vector<int> Amesos_Klu::VecAi
private

Definition at line 277 of file Amesos_Klu.h.

std::vector<double> Amesos_Klu::VecAval
private

Definition at line 278 of file Amesos_Klu.h.

double* Amesos_Klu::Aval
private

Definition at line 279 of file Amesos_Klu.h.

int* Amesos_Klu::Ai
private

Definition at line 280 of file Amesos_Klu.h.

int Amesos_Klu::UseDataInPlace_
private

1 if Problem_->GetOperator() is stored entirely on process 0

Definition at line 283 of file Amesos_Klu.h.

long long Amesos_Klu::numentries_
private

Number of non-zero entries in Problem_->GetOperator()

Definition at line 285 of file Amesos_Klu.h.

long long Amesos_Klu::NumGlobalElements_
private

Number of rows and columns in the Problem_->GetOperator()

Definition at line 287 of file Amesos_Klu.h.

Epetra_RowMatrix* Amesos_Klu::RowMatrixA_
private

Operator converted to a RowMatrix.

Definition at line 290 of file Amesos_Klu.h.

Epetra_CrsMatrix* Amesos_Klu::CrsMatrixA_
private

Operator converted to a CrsMatrix.

Definition at line 292 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_Map> Amesos_Klu::SerialMap_
private

Points to a Serial Map (unused if UseDataInPlace_ == 1 )

Definition at line 302 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_CrsMatrix> Amesos_Klu::SerialCrsMatrixA_
private

Points to a Serial Copy of A (unused if UseDataInPlace_==1)

Definition at line 304 of file Amesos_Klu.h.

Epetra_RowMatrix* Amesos_Klu::StdIndexMatrix_
private

Points to a Contiguous Copy of A.

Definition at line 306 of file Amesos_Klu.h.

Epetra_RowMatrix* Amesos_Klu::SerialMatrix_
private

Points to a Serial Copy of A.

Definition at line 308 of file Amesos_Klu.h.

bool Amesos_Klu::TrustMe_
private

If true, no checks are made and the matrix is assume to be distributed.

Definition at line 314 of file Amesos_Klu.h.

int Amesos_Klu::NumVectors_
private

Number of vectors in RHS and LHS.

Definition at line 316 of file Amesos_Klu.h.

double* Amesos_Klu::SerialXBvalues_
private

Pointer to the actual values in the serial version of X and B.

Definition at line 318 of file Amesos_Klu.h.

double* Amesos_Klu::SerialBvalues_
private

Definition at line 319 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_MultiVector> Amesos_Klu::SerialB_
private

Serial versions of the LHS and RHS (may point to the original vector if serial)

Definition at line 321 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_MultiVector> Amesos_Klu::SerialX_
private

Definition at line 322 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_MultiVector> Amesos_Klu::SerialXextract_
private

Serial versions of the LHS and RHS (if necessary)

Definition at line 324 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_MultiVector> Amesos_Klu::SerialBextract_
private

Definition at line 325 of file Amesos_Klu.h.

bool Amesos_Klu::UseTranspose_
private

If true, the transpose of A is used.

Definition at line 328 of file Amesos_Klu.h.

const Epetra_LinearProblem* Amesos_Klu::Problem_
private

Pointer to the linear system problem.

Definition at line 330 of file Amesos_Klu.h.

std::vector<int> Amesos_Klu::ColIndicesV_
private

Only used for RowMatrices to extract copies.

Definition at line 333 of file Amesos_Klu.h.

std::vector<double> Amesos_Klu::RowValuesV_
private

Only used for RowMatrices to extract copies.

Definition at line 335 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_Import> Amesos_Klu::ImportToSerial_
private

Importer to process 0.

Definition at line 337 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_Import> Amesos_Klu::ImportRangeToSerial_
private

Definition at line 338 of file Amesos_Klu.h.

Teuchos::RCP<Epetra_Import> Amesos_Klu::ImportDomainToSerial_
private

Definition at line 339 of file Amesos_Klu.h.

int Amesos_Klu::MtxRedistTime_
private

Quick access ids for the individual timings.

Definition at line 341 of file Amesos_Klu.h.

int Amesos_Klu::MtxConvTime_
private

Definition at line 341 of file Amesos_Klu.h.

int Amesos_Klu::VecRedistTime_
private

Definition at line 341 of file Amesos_Klu.h.

int Amesos_Klu::SymFactTime_
private

Definition at line 342 of file Amesos_Klu.h.

int Amesos_Klu::NumFactTime_
private

Definition at line 342 of file Amesos_Klu.h.

int Amesos_Klu::SolveTime_
private

Definition at line 342 of file Amesos_Klu.h.

int Amesos_Klu::OverheadTime_
private

Definition at line 342 of file Amesos_Klu.h.


The documentation for this class was generated from the following files: