Amesos_Merikos Class Reference

Amesos_Merikos: A parallel divide and conquer solver. More...

#include <Amesos_Merikos.h>

Inheritance diagram for Amesos_Merikos:

Inheritance graph
[legend]
Collaboration diagram for Amesos_Merikos:

Collaboration graph
[legend]
List of all members.

Public Member Functions

 Amesos_Merikos (const Epetra_LinearProblem &LinearProblem)
 Amesos_Merikos Constructor.
 ~Amesos_Merikos (void)
 Amesos_Merikos Destructor.
int RedistributeA ()
 Performs SymbolicFactorization on the matrix A.
int ConvertToScalapack ()
int PerformNumericFactorization ()
int SymbolicFactorization ()
 Performs SymbolicFactorization on the matrix A.
int NumericFactorization ()
 Performs NumericFactorization on the matrix A.
int LSolve ()
 Solves L X = B.
int USolve ()
 Solves U X = B.
int Solve ()
 Solves A X = B.
const Epetra_LinearProblemGetProblem () const
 Get a pointer to the Problem.
bool MatrixShapeOK () const
 Returns true if MERIKOS can handle this matrix shape.
int SetUseTranspose (bool UseTranspose)
 SetUseTranpose() controls whether to compute AX=B or ATX = B.
bool UseTranspose () const
 Returns the current UseTranspose setting.
const Epetra_CommComm () const
 Returns a pointer to the Epetra_Comm communicator associated with this matrix.
int SetParameters (Teuchos::ParameterList &ParameterList)
 Updates internal variables.
int NumSymbolicFact () const
 Returns the number of symbolic factorizations performed by this object.
int NumNumericFact () const
 Returns the number of numeric factorizations performed by this object.
int NumSolve () const
 Returns the number of solves performed by this object.
void PrintTiming ()
 Print timing information.
void PrintStatus ()
 Print information about the factorization and solution phases.

Protected Attributes

bool UseTranspose_
const Epetra_LinearProblemProblem_
Epetra_CrsMatrixL
Epetra_CrsMatrixU
bool PrintTiming_
bool PrintStatus_
bool ComputeVectorNorms_
bool ComputeTrueResidual_
int verbose_
int debug_
double ConTime_
double SymTime_
double NumTime_
double SolTime_
double VecTime_
double MatTime_
int NumSymbolicFact_
int NumNumericFact_
int NumSolve_
Epetra_TimeTime_
Epetra_MapScaLAPACK1DMap_
Epetra_CrsMatrixScaLAPACK1DMatrix_
Epetra_MapVectorMap_
std::vector< double > DenseA_
std::vector< int > Ipiv_
int NumOurRows_
int NumOurColumns_
bool TwoD_distribution_
int grid_nb_
int mypcol_
int myprow_
Epetra_CrsMatrixFatOut_
int nb_
int lda_
int iam_
int nprow_
int npcol_
int NumGlobalElements_
int m_per_p_

Detailed Description

Amesos_Merikos: A parallel divide and conquer solver.



Merikos partitions the rows of a matrix into two or more disjoint submatrices. i.e. if rows i and j are in different submatrices, A[i,j] == 0 == A[j,i]. Rows/columns not in any of the submatrices, i.e. the rows/columsn of the separator, are permuted to the bottom right.



Merikos factors each of the disjoint submatrices in parallel, (potentially by calling Amesos_Merikos() recursively), updating the rows and columns of the separator which belong to it and forming the schur complement of those rows and columns of the separator.



Merikos updates the trailing block of the matrix and then factors it.



Merikos is a Greek word for partial, reflecting the fact that Amesos_Merikos uses a series of partial LU factorizations, performed in parallel, to piece together the full LU decomposition.


Constructor & Destructor Documentation

Amesos_Merikos::Amesos_Merikos ( const Epetra_LinearProblem LinearProblem  ) 

Amesos_Merikos Constructor.

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

Amesos_Merikos::~Amesos_Merikos ( void   ) 

Amesos_Merikos Destructor.

Completely deletes an Amesos_Merikos object.


Member Function Documentation

int Amesos_Merikos::LSolve (  ) 

Solves L X = B.

     | L11   0   0  |  X1      B1
     | L21 L22   0  |  X2   =  B2
     | L31 L32 L33  |  X3   =  B3

    Foreach subblock of the matrix do:
      Note:  this will happen in parallel
      Lsolve() 
        i.e. L11.Solve(X1, B1) and L22.Solve(X2, B2) 
      Update the elements of B corresponding to the seperator,
        i.e. B3 = B3 - L31 X1 - L32 X2 
    Endfor
    Perform a solve on the trailing matrix:
      i.e. L33.LSolve(X3,B3) 

Returns:
Integer error code, set to 0 if successful.

bool Amesos_Merikos::MatrixShapeOK (  )  const [virtual]

Returns true if MERIKOS can handle this matrix shape.

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

Implements Amesos_BaseSolver.

int Amesos_Merikos::NumericFactorization (  )  [virtual]

Performs NumericFactorization on the matrix A.

    Static pivoting (i.e. scale and permute the matrix to
      produce a zero-free diagonal and to minimize the need for
      pivoting later).
    Partition the matrix 
    Redistribute the matrix to match the partitioning
    Foreach subblock of the matrix do:
      Note:  this will happen in parallel
      Create an instance of an Amesos solver object (must  
        support the Amesos_Component interface)
      Call PartialFactorization 
      Add the Schur Complement into the trailing block of the matrix.
    Endfor
    Create an Amesos instance for the trailing block of the matrix.
    Call SymbolicFactorization on the trailing block 
    Call NumericFactorization on the trailing block 

Returns:
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

int Amesos_Merikos::RedistributeA (  ) 

Performs SymbolicFactorization on the matrix A.

SymbolicFactorization() takes no action in Amesos_Merikos().

Returns:
Integer error code, set to 0 if successful.

int Amesos_Merikos::SetParameters ( Teuchos::ParameterList ParameterList  )  [virtual]

Updates internal variables.

Returns:
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

int Amesos_Merikos::Solve (  )  [virtual]

Solves A X = B.

     | L11     U12     U13  |  X1      B1
     | L21     L22     U23  |  X2   =  B2
     | L31     L32         A33  |  X3   =  B3

    Foreach subblock of the matrix do:
      Note:  this will happen in parallel
      Lsolve() 
        i.e. L11.Solve(X1, B1) and L22.Solve(X2, B2) 
      Update the elements of B corresponding to the seperator,
        i.e. B3 = B3 - L31 X1 - L32 X2 
    Endfor
    Perform a solve on the trailing matrix:
      i.e. A33.Solve(X3,B3)

    B = X ;
    Foreach subblock of the matrix do:
      Note:  this will happen in parallel
      Update the elements of B corresponding to this block
        i.e. B2 = B2 - U23 X3 ; B1 = B1 - U13 X3 
      Usolve() 
        i.e. U11.Solve(X1, B1) and U22.Solve(X2, B2) 
    Endfor

Returns:
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

int Amesos_Merikos::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.

int Amesos_Merikos::USolve (  ) 

Solves U X = B.

     | U11 U12 U13  |  X1      B1
     |   0 U22 U23  |  X2   =  B2
     |   0   0 U33  |  X3   =  B3

    Perform a solve on the trailing matrix:
      i.e. U33.USolve(X3,B3) 
    Foreach subblock of the matrix do:
      Note:  this will happen in parallel
      Update the elements of B corresponding to this block
        i.e. B2 = B2 - U23 X3 ; B1 = B1 - U13 X3 
      Usolve() 
        i.e. U11.Solve(X1, B1) and U22.Solve(X2, B2) 
    Endfor

Returns:
Integer error code, set to 0 if successful.


The documentation for this class was generated from the following file:
Generated on Wed May 12 21:39:57 2010 for Amesos by  doxygen 1.4.7