ML  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Classes | Typedefs | Functions
MLAPI Namespace Reference

MLAPI: Default namespace for all MLAPI objects and functions. More...

Classes

class  BaseLinearCombination
 
class  BaseObject
 Basic class for MLAPI objects. More...
 
class  BaseOperator
 Base class for all MLAPI objects. More...
 
class  CompObject
 Class to count flops. More...
 
class  DistributedMatrix
 
class  EpetraBaseOperator
 Basic class to wrap MLAPI::InverseOperator into Epetra_Operator. More...
 
struct  StackEntry
 
class  InverseOperator
 InverseOperator: basic class to define smoother and coarse solvers. More...
 
class  LinearCombinationAdd
 
class  LinearCombinationMixed
 
class  LinearCombinationScaled
 
class  MultiVectorScaled
 
class  MultiVectorCombination
 
class  BaseOperatorTimesMultiVector
 
class  Residual
 
class  LoadBalanceInverseOperator
 
class  LoadBalanceOperator
 
class  MATLABStream
 Basic stream to save in a MATLAB-compatible file MLAPI objects. More...
 
class  MultiLevelAdaptiveSA
 Black-box multilevel adaptive smoothed aggregation preconditioner. More...
 
class  MultiLevelSA
 Black-box multilevel smoothed aggregation preconditioner. More...
 
class  DoubleVector
 
class  MultiVector
 Basic class for distributed double-precision vectors. More...
 
class  Operator
 Operator: basic class to define operators within MLAPI. More...
 
class  ML_Operator_Box
 Simple wrapper for ML_Operator struct. More...
 
class  Epetra_SerialMatrix
 
class  SerialMatrix
 
class  Space
 Specifies the number and distribution among processes of elements. More...
 
class  TimeObject
 Class to track time spent in an object. More...
 

Typedefs

typedef struct MLAPI::StackEntry Entry
 

Functions

void GetPtent (const Operator &A, Teuchos::ParameterList &List, const MultiVector &ThisNS, Operator &Ptent, MultiVector &NextNS)
 Builds the tentative prolongator using aggregation.
 
void GetPtent (const Operator &A, Teuchos::ParameterList &List, Operator &Ptent)
 Builds the tentative prolongator with default null space.
 
void GetPtent (const Epetra_RowMatrix &A, Teuchos::ParameterList &List, double *thisns, Teuchos::RCP< Epetra_CrsMatrix > &Ptent, Teuchos::RCP< Epetra_MultiVector > &NextNS, const int domainoffset=0)
 Builds the tentative prolongator using aggregation. More...
 
void GetPtent (const Epetra_RowMatrix &A, Teuchos::ParameterList &List, double *thisns, Teuchos::RCP< Epetra_CrsMatrix > &Ptent, const int domainoffset=0)
 Builds the tentative prolongator using aggregation. More...
 
int GetAggregates (Epetra_RowMatrix &A, Teuchos::ParameterList &List, double *thisns, Epetra_IntVector &aggrinfo)
 Call ML aggregation on A according to parameters supplied in List. Return aggregates in aggrinfo. More...
 
int GetGlobalAggregates (Epetra_RowMatrix &A, Teuchos::ParameterList &List, double *thisns, Epetra_IntVector &aggrinfo)
 Call ML aggregation on A according to parameters supplied in List. Return aggregates in aggrinfo. More...
 
std::ostream & operator<< (std::ostream &os, const BaseObject &obj)
 
void SetDefaults (Teuchos::ParameterList &List)
 Sets default values in input List.
 
double MaxEigAnorm (const Operator &Op, const bool DiagonalScaling=false)
 Computes the maximum eigenvalue of Op using the A-norm of the operator.
 
double MaxEigCG (const Operator &Op, const bool DiagonalScaling=false)
 Computes the maximum eigenvalue of Op using the CG method.
 
double MaxEigPowerMethod (const Operator &Op, const bool DiagonalScaling=false)
 Computes the maximum eigenvalue of Op using the power method.
 
double MaxEigAnasazi (const Operator &Op, const bool DiagonalScaling=false)
 Computes the maximum eigenvalue of Op using Anasazi.
 
void Eig (const Operator &Op, MultiVector &ER, MultiVector &EI)
 Computes eigenvalues and eigenvectors using LAPACK (w/ one process only).
 
void Eigs (const Operator &A, int NumEigenvalues, MultiVector &ER, MultiVector &EI)
 
void StackPop ()
 
void StackPrint ()
 
MultiVectorCombination operator+ (const MultiVector &x, const MultiVector &y)
 Creates a new MultiVector, defined as x + y.
 
LinearCombinationAdd operator+ (const BaseLinearCombination &left, const BaseLinearCombination &right)
 
LinearCombinationMixed operator+ (const BaseLinearCombination &left, const MultiVector &right)
 
LinearCombinationMixed operator+ (const MultiVector &left, const BaseLinearCombination &right)
 
MultiVectorCombination operator+ (const MultiVectorScaled &left, const MultiVectorScaled &right)
 Creates a new MultiVector, defined as alpha * x + beta * y.
 
Residual operator+ (const MultiVectorScaled &left, const BaseOperatorTimesMultiVector &right)
 Creates a new MultiVector, defined as alpha * x + A * y.
 
Residual operator+ (const MultiVector &left, const BaseOperatorTimesMultiVector &right)
 Creates a new MultiVector, defined as x + A * y.
 
MultiVectorCombination operator- (const MultiVector &x, const MultiVector &y)
 Creates a new MultiVector, defined as x - y.
 
LinearCombinationAdd operator- (const BaseLinearCombination &left, const BaseLinearCombination &right)
 
LinearCombinationMixed operator- (const BaseLinearCombination &left, const MultiVector &right)
 
LinearCombinationMixed operator- (const MultiVector &left, const BaseLinearCombination &right)
 
Residual operator- (const MultiVector &left, const BaseOperatorTimesMultiVector &right)
 Creates a new MultiVector, defined as x - A * y.
 
MultiVector operator+ (const MultiVector &x, const double alpha)
 Creates a new MultiVector, defined as x + alpha.
 
MultiVector operator+ (const double alpha, const MultiVector &x)
 Creates a new MultiVector, defined as alpha + x.
 
MultiVector operator- (const MultiVector &x, const double alpha)
 Creates a new MultiVector, defined as x - alpha.
 
MultiVector operator- (const double alpha, const MultiVector &x)
 Creates a new MultiVector, defined as alpha - y.
 
Operator operator+ (const Operator &A, const Operator &B)
 Creates a new Operator, defined as A + B.
 
Operator operator- (const Operator &A, const Operator &B)
 Creates a new Operator, defined as A - B.
 
Operator operator* (const Operator &A, const Operator &B)
 Creates a new Operator, defined as A * B.
 
Operator operator* (const Operator &A, const double alpha)
 Creates a new Operator, defined as A * alpha.
 
Operator operator* (const double alpha, const Operator &A)
 Creates a new Operator, defined as alpha * A.
 
Operator operator/ (const Operator &A, const double alpha)
 Creates a new Operator, defined as A / alpha.
 
MultiVector operator* (const MultiVector &x, const double alpha)
 Creates a new MultiVector, defined as x * alpha.
 
MultiVector operator* (const double alpha, const MultiVector &x)
 
MultiVector operator/ (const MultiVector &x, const double alpha)
 Creates a new MultiVector y, such that y = x / alpha.
 
BaseOperatorTimesMultiVector operator* (const BaseOperator &A, const MultiVector &x)
 Creates a new MultiVector y, such that y = A * x.
 
BaseOperatorTimesMultiVector operator* (const BaseOperator &A, const BaseLinearCombination &x)
 Creates a new MultiVector y, such that y = A * x (x is a BaseLinearCombination)
 
double operator* (const MultiVector &x, const MultiVector &y)
 Computes the dot product between the first vector in x and y.
 
double operator* (const MultiVector &x, const BaseLinearCombination &y)
 
double operator* (const BaseLinearCombination &x, const MultiVector &y)
 
double operator* (const BaseLinearCombination &x, const BaseLinearCombination &y)
 
Operator Gallery (const std::string ProblemType, const Space &MySpace)
 Creates a matrix using the TRIUTILS gallery.
 
Operator GetShiftedLaplacian1D (const int NX, const double Factor=0.99)
 Creates a 1D shifted Laplacian.
 
Operator GetShiftedLaplacian2D (const int NX, const int NY, const double Factor=0.99, const bool RandomScale=false)
 Creates a 2D shifted Laplacian.
 
Operator ReadMatrix (const char *FileName)
 Reads a matrix in MATLAB format.
 
Operator GetRecirc2D (const int NX, const int NY, const double conv, const double diff)
 Creates a recirculation problem in 2D.
 
Teuchos::ParameterList ReadParameterList (const char *FileName)
 Populates a list from specified file.
 
void Krylov (const Operator &A, const MultiVector &LHS, const MultiVector &RHS, const BaseOperator &Prec, Teuchos::ParameterList &List)
 
MultiVector Duplicate (const MultiVector &y)
 Creates a new vector, x, such that x = y.
 
MultiVector Duplicate (const MultiVector &y, const int v)
 Creates a new vector, x, such that x = y(:,v)
 
MultiVector Extract (const MultiVector &y, const int v)
 Extracts a component from a vector.
 
MultiVector Redistribute (const MultiVector &y, const int NumEquations)
 Redistributes the entry of a vector as a multivector.
 
Operator GetRAP (const Operator &R, const Operator &A, const Operator &P)
 Performs a triple matrix-matrix product, res = R * A *P.
 
Operator GetTranspose (const Operator &A, const bool byrow=true)
 Returns a newly created transpose of A.
 
Operator GetIdentity (const Space &DomainSpace, const Space &RangeSpace)
 Returns the identity matrix.
 
MultiVector GetDiagonal (const Operator &A)
 Returns a vector containing the diagonal elements of A.
 
MultiVector GetDiagonal (const Operator &A, const int offset)
 Returns a vector containing the diagonal elements of A.
 
Operator GetDiagonal (const MultiVector &D)
 Returns a newly created operator, containing D on the diagonal.
 
Operator GetJacobiIterationOperator (const Operator &Amat, double Damping)
 Returns an operator defined as (I - Damping A).
 
Operator GetPtent1D (const MultiVector &D, const int offset=0)
 Returns a newly created operator, containing D on the diagonal.
 
int ML_Operator_Add2 (ML_Operator *A, ML_Operator *B, ML_Operator *C, int matrix_type, double scalarA, double scalarB)
 
void AnalyzeCheap (const Operator &A)
 Performs a cheap analysis of the properties of the input operator.
 
void PrintSparsity (const Operator &A, int NumPDEEquations=1)
 Prints on file the sparsity structure of input operator.
 
Operator GetScaledOperator (const Operator &A, const double alpha)
 Multiply A by a double value, alpha.
 
Operator Duplicate (const Operator &A)
 Duplicates a given operator.
 
void ReadSAMISMatrix (const char *filen, Operator &A, int &NumPDEEqns)
 Reads symmetric matrix from SAMIS binary format.
 
void ReadSAMISKernel (const char *myKerFileName, MultiVector &A, const int limKer=-1)
 Reads null space vectors from SAMIS binary format.
 
ML_Comm * GetML_Comm ()
 Returns a pointer to the ML_Comm object defined on MPI_COMM_WORLD.
 
Epetra_CommGetEpetra_Comm ()
 Returns a reference to the Epetra_Comm object defined on MPI_COMM_WORLD.
 
void Barrier ()
 Calls Mpi_Barrier() if MPI is enabled.
 
int GetMyPID ()
 Returns the ID of the calling process.
 
int GetNumProcs ()
 Returns the total number of processes in the computation.
 
int GetPrintLevel ()
 Retutns the level of output (always 0 if MyPID() != 0).
 
void SetPrintLevel (int Level)
 Sets the level of output (from 0 to 10, 0 being verbose).
 
void Init (USR_COMM comm=USR_COMM_WORLD)
 Initialize the MLAPI workspace.
 
void Finalize ()
 Destroys the MLAPI workspace.
 
std::string GetString (const int &x)
 
std::string GetString (const double &x)
 
int GetMatrixType ()
 

Detailed Description

MLAPI: Default namespace for all MLAPI objects and functions.

Function Documentation

int MLAPI::GetAggregates ( Epetra_RowMatrix A,
Teuchos::ParameterList List,
double *  thisns,
Epetra_IntVector aggrinfo 
)

Call ML aggregation on A according to parameters supplied in List. Return aggregates in aggrinfo.

   On input, map of aggrinfo has to map row map of A. On output, aggrinfo[i]
   contains number of aggregate the row belongs to, where aggregates are
   numbered starting from 0.
   Return value is the processor-local number of aggregates build.
   If aggrinfo[i] >= return-value, then i is a processor local row
   of a row that ML has detected to be on a Dirichlet BC.
Parameters
A(in): Matrix to be aggregated on
List(in): ParameterList containing ML options
thisns(in): nullspace
aggrinfo(out),:vector containing aggregation information
Note
Map of aggrinfo has to match rowmap of A on input.
Returns
returns processor-local number of aggregates
Author
Michael Gee (gee@lnm.mw.tum.de)
int MLAPI::GetGlobalAggregates ( Epetra_RowMatrix A,
Teuchos::ParameterList List,
double *  thisns,
Epetra_IntVector aggrinfo 
)

Call ML aggregation on A according to parameters supplied in List. Return aggregates in aggrinfo.

   On input, map of aggrinfo has to map row map of A. On output, aggrinfo[i]
   contains number of global aggregate the row belongs to, where aggregates are
   numbered starting from 0 globally.
   Return value is the processor-local number of aggregates build.
   If aggrinfo[i] < 0, then i is a processor local row
   that ML has detected to be on a Dirichlet BC.
   if aggrinfo[i] >= 0, then i is a processor local row and aggrinfo[i] is
   a global aggregate id.
Parameters
A(in): Matrix to be aggregated on
List(in): ParameterList containing ML options
thisns(in): nullspace
aggrinfo(out),:vector containing aggregation information in global numbering
Note
Map of aggrinfo has to match rowmap of A on input.
Returns
returns processor-local number of aggregates
Author
Michael Gee (gee@lnm.mw.tum.de)
void MLAPI::GetPtent ( const Epetra_RowMatrix A,
Teuchos::ParameterList List,
double *  thisns,
Teuchos::RCP< Epetra_CrsMatrix > &  Ptent,
Teuchos::RCP< Epetra_MultiVector > &  NextNS,
const int  domainoffset = 0 
)

Builds the tentative prolongator using aggregation.

Build Ptent and NextNS as usual but as Epetra objects.

Parameters
A(in): Matrix to be aggregated on
List(in): ParameterList containing ML options
thisns(in): nullspace in format ML accepts
Ptent(out),:Matrix containing tentative prolongator
NextNS(out): MultiVector containing next level nullspace.
domainoffset(in,optional): give an offset such that the domainmap of Ptent starts global numbering from domainoffset instead from zero. This is useful to create block matrices.
Author
Michael Gee (gee@lnm.mw.tum.de)
void MLAPI::GetPtent ( const Epetra_RowMatrix A,
Teuchos::ParameterList List,
double *  thisns,
Teuchos::RCP< Epetra_CrsMatrix > &  Ptent,
const int  domainoffset = 0 
)

Builds the tentative prolongator using aggregation.

Build Ptent and NextNS as usual but as Epetra objects.

Parameters
A(in): Matrix to be aggregated on
List(in): ParameterList containing ML options
thisns(in): nullspace in format ML accepts
Ptent(out),:Matrix containing tentative prolongator
domainoffset(in,optional): give an offset such that the domainmap of Ptent starts global numbering from domainoffset instead from zero. This is useful to create block matrices.
Author
Michael Gee (gee@lnm.mw.tum.de)