Amesos  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Protected Attributes | List of all members
Amesos_Scalapack Class Reference

Amesos_Scalapack: A serial and parallel dense solver. For now, we implement only the unsymmetric ScaLAPACK solver. More...

#include <Amesos_Scalapack.h>

Inheritance diagram for Amesos_Scalapack:
Inheritance graph
[legend]
Collaboration diagram for Amesos_Scalapack:
Collaboration graph
[legend]

Public Member Functions

 Amesos_Scalapack (const Epetra_LinearProblem &LinearProblem)
 Amesos_Scalapack Constructor. More...
 
 ~Amesos_Scalapack (void)
 Amesos_Scalapack 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.
 
bool MatrixShapeOK () const
 Returns true if SCALAPACK can handle this matrix shape. More...
 
int SetUseTranspose (bool UseTranspose)
 SetUseTranpose(true) is more efficient in Amesos_Scalapack. More...
 
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. More...
 
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 () const
 Print timing information.
 
void PrintStatus () const
 Print information about the factorization and solution phases.
 
void GetTiming (Teuchos::ParameterList &TimingParameterList) const
 Extracts timing information from the current solver and places it in the parameter list.
 
- Public Member Functions inherited from Amesos_BaseSolver
virtual ~Amesos_BaseSolver ()
 Destructor.
 
virtual void setParameterList (Teuchos::RCP< Teuchos::ParameterList > const &paramList)
 Redefined from Teuchos::ParameterListAcceptor (Does Not Work)
 
virtual Teuchos::RCP
< Teuchos::ParameterList
getNonconstParameterList ()
 This is an empty stub.
 
virtual Teuchos::RCP
< Teuchos::ParameterList
unsetParameterList ()
 This is an empty stub.
 
- Public Member Functions inherited from Teuchos::ParameterListAcceptor
virtual RCP< const ParameterListgetParameterList () const
 
virtual RCP< const ParameterListgetValidParameters () const
 

Protected Attributes

int iam_
 
int NumGlobalElements_
 
int NumGlobalNonzeros_
 
int nprow_
 
int npcol_
 
int ictxt_
 
int m_per_p_
 
int DescA_ [10]
 
Epetra_MapScaLAPACK1DMap_
 
Epetra_CrsMatrixScaLAPACK1DMatrix_
 
Epetra_MapVectorMap_
 
std::vector< double > DenseA_
 
std::vector< int > Ipiv_
 
int NumOurRows_
 
int NumOurColumns_
 
bool UseTranspose_
 
const Epetra_LinearProblemProblem_
 
double ConTime_
 
double SymTime_
 
double NumTime_
 
double SolTime_
 
double VecTime_
 
double MatTime_
 
bool TwoD_distribution_
 
int grid_nb_
 
int mypcol_
 
int myprow_
 
Epetra_CrsMatrixFatOut_
 
int nb_
 
int lda_
 
Epetra_TimeTime_
 

Detailed Description

Amesos_Scalapack: A serial and parallel dense solver. For now, we implement only the unsymmetric ScaLAPACK solver.

 Amesos_Scalapack, an object-oriented wrapper for LAPACK and ScaLAPACK, 

will solve a linear systems of equations: A X = B using Epetra objects and the ScaLAPACK library, where A is an Epetra_RowMatrix and X and B are Epetra_MultiVector objects.



Amesos_Scalapack can be competitive for matrices that are not particularly sparse. ScaLAPACK solves matrices for which the fill-in is roughly 10% to 20% of the matrix size in time comparable to that achieve by other Amesos classes. Amesos_Scalapack scales well and hence its performance advantage will be largest when large number of processes are involved.



Amesos_Scalapack uses the ScaLAPACK functions PDGETRF and PDGETRS if more than one process is used. If only one process is used, Amesos_ScaLAPACK uses the LAPACK function PDGETRF and PDGETRS.



AmesosScaLAPACK uses full partial pivoting and will therefore provide answers that are at least as accurate as any direct sparse solver.



AmesosScalapack makes sense under the following circumstances:

Common control parameters :

Amesos_Scalapack supports the following parameters which are common to across multiple Amesos solvers:

Amesos_Scalapack supports the following parameters specific to Amesos_Scalapack.
Teuchos::ParameterList ScalapackParams = ParameterList.sublist("Scalapack") ;

Limitations:

None of the following limitations would be particularly difficult to remove.



The present implementation limits the number of right hand sides to the number of rows assigned to each process. i.e. nrhs < n/p.



The present implementation does not take advantage of symmetric or symmetric positive definite matrices, although ScaLAPACK has separate routines to take advantages of such matrices.

Constructor & Destructor Documentation

Amesos_Scalapack::Amesos_Scalapack ( const Epetra_LinearProblem LinearProblem)

Amesos_Scalapack Constructor.

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

Note: The operator in LinearProblem must be an Epetra_RowMatrix.

References SetParameters().

Amesos_Scalapack::~Amesos_Scalapack ( void  )

Member Function Documentation

bool Amesos_Scalapack::MatrixShapeOK ( ) const
virtual

Returns true if SCALAPACK can handle this matrix shape.

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

Implements Amesos_BaseSolver.

References GetProblem().

int Amesos_Scalapack::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().

preconditions:

postconditions:

  • nprow_, npcol_, DescA_
  • DenseA will be factored
  • Ipiv_ contains the pivots
Returns
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

References Comm(), Amesos_Status::debug_, Epetra_LinearProblem::GetOperator(), Epetra_Comm::MyPID(), Epetra_BlockMap::NumGlobalElements(), Epetra_RowMatrix::NumGlobalNonzeros(), Amesos_Status::NumNumericFact_, and Epetra_RowMatrix::RowMatrixRowMap().

int Amesos_Scalapack::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 subsequent 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.

References Amesos_Status::debug_, Teuchos::ParameterList::get(), Teuchos::ParameterList::isParameter(), Teuchos::ParameterList::isSublist(), and Teuchos::ParameterList::sublist().

Referenced by Amesos_Scalapack().

int Amesos_Scalapack::SetUseTranspose ( bool  UseTranspose)
inlinevirtual

SetUseTranpose(true) is more efficient in Amesos_Scalapack.

Implements Amesos_BaseSolver.

References UseTranspose().

int Amesos_Scalapack::Solve ( void  )
virtual

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

preconditions:

postconditions:

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

Implements Amesos_BaseSolver.

References Epetra_Comm::Broadcast(), Comm(), Amesos_Status::ComputeTrueResidual_, Amesos_Status::ComputeVectorNorms_, Amesos_Status::debug_, Epetra_Time::ElapsedTime(), Epetra_LinearProblem::GetLHS(), Epetra_LinearProblem::GetMatrix(), Epetra_LinearProblem::GetOperator(), Epetra_LinearProblem::GetRHS(), Insert, Epetra_RowMatrix::Multiply(), Epetra_Comm::MyPID(), Amesos_Status::NumSolve_, Epetra_Time::ResetStartTime(), Epetra_RowMatrix::RowMatrixRowMap(), UseTranspose(), and Amesos_Status::verbose_.

int Amesos_Scalapack::SymbolicFactorization ( )
virtual

Performs SymbolicFactorization on the matrix A.

There is no symbolic factorization phase in ScaLAPACK, as it operates only on dense matrices. Hence, Amesos_Scalapack::SymbolicFactorization() takes no action.

Returns
Integer error code, set to 0 if successful.

Implements Amesos_BaseSolver.

References Comm(), Amesos_Status::debug_, and Amesos_Status::NumSymbolicFact_.


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