Trilinos_Util::CrsMatrixGallery Class Reference

#include <Trilinos_Util_CrsMatrixGallery.h>

Inheritance diagram for Trilinos_Util::CrsMatrixGallery:

Inheritance graph
[legend]
Collaboration diagram for Trilinos_Util::CrsMatrixGallery:

Collaboration graph
[legend]
List of all members.
Epetra_CrsMatrixGetMatrix ()
 Returns a pointer to the CrsMatrix.
Epetra_CrsMatrixGetMatrixRef ()
Epetra_MultiVectorGetExactSolution ()
 Returns a pointer to the exact solution.
Epetra_MultiVectorGetStartingSolution ()
 Returns a pointer to the starting solution (typically, for HB problems).
Epetra_MultiVectorGetRHS ()
 Returns a pointer to the rhs corresponding to the selected exact solution.
const Epetra_MapGetMap ()
 Returns a pointer the internally stored Map.
const Epetra_MapGetMapRef ()
Epetra_LinearProblemGetLinearProblem ()
 Returns a pointer to Epetra_LinearProblem.
void ComputeResidual (double *residual)
 Computes the 2-norm of the residual.
void ComputeDiffBetweenStartingAndExactSolutions (double *residual)
 Computes the 2-norm of the difference between the starting solution and the exact solution.
void PrintMatrixAndVectors (ostream &os)
 Print out matrix and vectors.
void PrintMatrixAndVectors ()
void GetCartesianCoordinates (double *&x, double *&y, double *&z)
 Get pointers to double vectors containing coordinates of points.
int WriteMatrix (const string &FileName, const bool UseSparse=true)
 Print matrix on file in MATLAB format.
ostream & operator<< (ostream &os, const Trilinos_Util::CrsMatrixGallery &G)
 Print out detailed information about the problem at hand.

Public Member Functions

 CrsMatrixGallery (const string name, const Epetra_Comm &comm)
 Triutils_Gallery Constructor.
 CrsMatrixGallery (const string name, const Epetra_Map &map)
 Creates an Triutils_Gallery object using a given map.
 ~CrsMatrixGallery ()
 Triutils_Gallery destructor.
int Set (const string parameter, const int value)
 Sets a gallery options using an interger value.
int Set (const string parameter, const string value)
 Sets a gallery options using a C++ string .
int Set (const string parameter, const double value)
 Sets a gallery options using an double value.
int Set (const string parameter, const Epetra_Vector &value)
 Sets a gallery options using an Epetra_Vector.
int Set (Trilinos_Util::CommandLineParser &CLP)
 Sets gallery options using values passed from the shell.

Protected Member Functions

void CreateMap ()
 Creates a map.
void CreateMatrix ()
 Creates the CrdMatrix.
void CreateExactSolution ()
 Creates the exact solution.
void CreateStartingSolution ()
 Creates the starting solution.
void CreateRHS ()
 Create the RHS corresponding to the desired exact solution.
void CreateEye ()
void CreateMatrixDiag ()
void CreateMatrixTriDiag ()
void CreateMatrixLaplace1d ()
void CreateMatrixLaplace1dNeumann ()
void CreateMatrixCrossStencil2d ()
void CreateMatrixCrossStencil2dVector ()
void CreateMatrixLaplace2d ()
void CreateMatrixLaplace2d_BC ()
void CreateMatrixLaplace2d_9pt ()
void CreateMatrixStretched2d ()
void CreateMatrixRecirc2d ()
void CreateMatrixRecirc2dDivFree ()
void CreateMatrixLaplace2dNeumann ()
void CreateMatrixUniFlow2d ()
void CreateMatrixLaplace3d ()
void CreateMatrixCrossStencil3d ()
void CreateMatrixCrossStencil3dVector ()
void CreateMatrixLehmer ()
void CreateMatrixMinij ()
void CreateMatrixRis ()
void CreateMatrixHilbert ()
void CreateMatrixJordblock ()
void CreateMatrixCauchy ()
void CreateMatrixFiedler ()
void CreateMatrixHanowa ()
void CreateMatrixKMS ()
void CreateMatrixParter ()
void CreateMatrixPei ()
void CreateMatrixOnes ()
void CreateMatrixVander ()
void ReadMatrix ()
void GetNeighboursCartesian2d (const int i, const int nx, const int ny, int &left, int &right, int &lower, int &upper)
void GetNeighboursCartesian3d (const int i, const int nx, const int ny, const int nz, int &left, int &right, int &lower, int &upper, int &below, int &above)
void ZeroOutData ()
void SetupCartesianGrid2D ()
void SetupCartesianGrid3D ()
void ExactSolQuadXY (double x, double y, double &u)
void ExactSolQuadXY (double x, double y, double &u, double &ux, double &uy, double &uxx, double &uyy)

Protected Attributes

const Epetra_Commcomm_
Epetra_CrsMatrixmatrix_
Epetra_MultiVectorExactSolution_
Epetra_MultiVectorStartingSolution_
Epetra_MultiVectorrhs_
Epetra_Mapmap_
Epetra_LinearProblemLinearProblem_
string name_
int NumGlobalElements_
int NumMyElements_
int * MyGlobalElements_
string MapType_
bool ContiguousMap_
std::vector< int > MapMap_
string ExactSolutionType_
string StartingSolutionType_
string ExpandType_
string RhsType_
int nx_
int ny_
int nz_
int mx_
int my_
int mz_
double lx_
double ly_
double lz_
int NumPDEEqns_
int NumVectors_
Epetra_VectorVectorA_
Epetra_VectorVectorB_
Epetra_VectorVectorC_
Epetra_VectorVectorD_
Epetra_VectorVectorE_
Epetra_VectorVectorF_
Epetra_VectorVectorG_
double a_
double b_
double c_
double d_
double e_
double f_
double g_
double alpha_
double beta_
double gamma_
double delta_
double conv_
double diff_
double source_
double epsilon_
string FileName_
string ErrorMsg
string OutputMsg
bool verbose_

Constructor & Destructor Documentation

Trilinos_Util::CrsMatrixGallery::CrsMatrixGallery ( const string  name,
const Epetra_Comm comm 
)

Triutils_Gallery Constructor.

Creates a Triutils_Gallery instance.

The first parameter is the name of the matrix. We refer to the Trilinos Tutorial for a detailed description of available matrices.

Note:
The matrix name can be empty (""), and set later using, for example, Set("matrix_name","laplace_2d");
An example of program using this class is reported below.

int main(int argc, char *argv[])
{

#ifdef HAVE_MPI
  MPI_Init(&argc,&argv);
  Epetra_MpiComm Comm (MPI_COMM_WORLD);
#else
  Epetra_SerialComm Comm;
#endif

  // create an Epetra matrix reading an H/B matrix
  Trilinos_Util_CrsMatrixGallery Gallery("hb", Comm);

  // set the name of the matrix
  Gallery.Set("matrix name", "bcsstk14.rsa");
  
  Epetra_CrsMatrix * A;
  Epetra_Vector * ExactSolution;
  Epetra_Vector * RHS;
  Epetra_Vector * StartingSolution;

  // at this point the matrix is read from file
  A = Gallery.GetMatrix();
  ExactSolution = Gallery.GetExactSolution();

  // at this point the RHS is allocated and filled
  RHS = Gallery.GetRHS();
  StartingSolution = Gallery.GetStartingSolution();
  
  // create linear problem
  Epetra_LinearProblem Problem(A,StartingSolution,RHS);
  // create AztecOO instance
  AztecOO Solver(Problem);

  Solver.SetAztecOption( AZ_precond, AZ_dom_decomp );  
  Solver.Iterate(1000,1E-9);

  // compute residual
  double residual;
  
  Gallery.ComputeResidual(&residual);
  if( Comm.MyPID()==0 ) cout << "||b-Ax||_2 = " << residual << endl;
  
  Gallery.ComputeDiffBetweenStartingAndExactSolutions(&residual);
  if( Comm.MyPID()==0 ) cout << "||x_exact - x||_2 = " << residual << endl;

 #ifdef HAVE_MPI
  MPI_Finalize() ;
#endif

return 0 ;
  } 

Class CommandLineParser can be used as well. In this case, one may decide to use the following:

  Trilinos_Util::CommandLineParser CLP(argc,argv);
  // set a problem with no matrix name
  Trilinos_Util::CrsMatrixGallery Gallery("", Comm);
  // read parameters and settings from the shell line
  G.Set(CLP);
  // continue with your code...

Parameters:
In comm - Epetra communicator

Trilinos_Util::CrsMatrixGallery::CrsMatrixGallery ( const string  name,
const Epetra_Map map 
)

Creates an Triutils_Gallery object using a given map.

Create a Triutils_Gallery object using an Epetra_Map. Problem size must match the elements in map.

Parameters:
In name - definition of the problem to be created.
In map - Epetra_Map

Trilinos_Util::CrsMatrixGallery::~CrsMatrixGallery (  ) 

Triutils_Gallery destructor.


Member Function Documentation

void Trilinos_Util::CrsMatrixGallery::ComputeDiffBetweenStartingAndExactSolutions ( double *  residual  ) 

Computes the 2-norm of the difference between the starting solution and the exact solution.

void Trilinos_Util::CrsMatrixGallery::ComputeResidual ( double *  residual  ) 

Computes the 2-norm of the residual.

void Trilinos_Util::CrsMatrixGallery::CreateExactSolution (  )  [protected]

Creates the exact solution.

void Trilinos_Util::CrsMatrixGallery::CreateEye (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMap (  )  [protected]

Creates a map.

Creates an Epetra_Map. Before calling this function, the problem size must have been specified.

CreateMap() allows some different maps. The type of map is set using Set("map",value). Value is a string, defined as:

void Trilinos_Util::CrsMatrixGallery::CreateMatrix (  )  [protected]

Creates the CrdMatrix.

void Trilinos_Util::CrsMatrixGallery::CreateMatrixCauchy (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixCrossStencil2d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixCrossStencil2dVector (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixCrossStencil3d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixCrossStencil3dVector (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixDiag (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixFiedler (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixHanowa (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixHilbert (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixJordblock (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixKMS (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace1d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace1dNeumann (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace2d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace2d_9pt (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace2d_BC (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace2dNeumann (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLaplace3d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixLehmer (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixMinij (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixOnes (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixParter (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixPei (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixRecirc2d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixRecirc2dDivFree (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixRis (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixStretched2d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixTriDiag (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixUniFlow2d (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateMatrixVander (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::CreateRHS (  )  [protected]

Create the RHS corresponding to the desired exact solution.

void Trilinos_Util::CrsMatrixGallery::CreateStartingSolution (  )  [protected]

Creates the starting solution.

void Trilinos_Util::CrsMatrixGallery::ExactSolQuadXY ( double  x,
double  y,
double &  u,
double &  ux,
double &  uy,
double &  uxx,
double &  uyy 
) [protected]

void Trilinos_Util::CrsMatrixGallery::ExactSolQuadXY ( double  x,
double  y,
double &  u 
) [protected]

void Trilinos_Util::CrsMatrixGallery::GetCartesianCoordinates ( double *&  x,
double *&  y,
double *&  z 
)

Get pointers to double vectors containing coordinates of points.

Epetra_MultiVector * Trilinos_Util::CrsMatrixGallery::GetExactSolution (  ) 

Returns a pointer to the exact solution.

Returns a pointer to the exact solution.

Some choices are available to define the exact solution, using Set("exact solution", value). value can be:

Epetra_LinearProblem * Trilinos_Util::CrsMatrixGallery::GetLinearProblem (  ) 

Returns a pointer to Epetra_LinearProblem.

const Epetra_Map * Trilinos_Util::CrsMatrixGallery::GetMap (  ) 

Returns a pointer the internally stored Map.

const Epetra_Map & Trilinos_Util::CrsMatrixGallery::GetMapRef (  ) 

Epetra_CrsMatrix * Trilinos_Util::CrsMatrixGallery::GetMatrix (  ) 

Returns a pointer to the CrsMatrix.

Epetra_CrsMatrix & Trilinos_Util::CrsMatrixGallery::GetMatrixRef (  ) 

void Trilinos_Util::CrsMatrixGallery::GetNeighboursCartesian2d ( const int  i,
const int  nx,
const int  ny,
int &  left,
int &  right,
int &  lower,
int &  upper 
) [protected]

void Trilinos_Util::CrsMatrixGallery::GetNeighboursCartesian3d ( const int  i,
const int  nx,
const int  ny,
const int  nz,
int &  left,
int &  right,
int &  lower,
int &  upper,
int &  below,
int &  above 
) [protected]

Epetra_MultiVector * Trilinos_Util::CrsMatrixGallery::GetRHS (  ) 

Returns a pointer to the rhs corresponding to the selected exact solution.

Epetra_MultiVector * Trilinos_Util::CrsMatrixGallery::GetStartingSolution (  ) 

Returns a pointer to the starting solution (typically, for HB problems).

Returns a pointer to the starting solution. This is typically used while reading a HB problem. However, the user can set a starting solution using Set("starting solution", "value"). Value can be

void Trilinos_Util::CrsMatrixGallery::PrintMatrixAndVectors (  ) 

void Trilinos_Util::CrsMatrixGallery::PrintMatrixAndVectors ( ostream &  os  ) 

Print out matrix and vectors.

void Trilinos_Util::CrsMatrixGallery::ReadMatrix (  )  [protected]

int Trilinos_Util::CrsMatrixGallery::Set ( Trilinos_Util::CommandLineParser CLP  ) 

Sets gallery options using values passed from the shell.

int Trilinos_Util::CrsMatrixGallery::Set ( const string  parameter,
const Epetra_Vector value 
)

Sets a gallery options using an Epetra_Vector.

Sets a gallery options using an Epetra_Vector. The Epetra_Vector is copied into internal structures, and freed by the destructor.

int Trilinos_Util::CrsMatrixGallery::Set ( const string  parameter,
const double  value 
)

Sets a gallery options using an double value.

int Trilinos_Util::CrsMatrixGallery::Set ( const string  parameter,
const string  value 
)

Sets a gallery options using a C++ string .

int Trilinos_Util::CrsMatrixGallery::Set ( const string  parameter,
const int  value 
)

Sets a gallery options using an interger value.

void Trilinos_Util::CrsMatrixGallery::SetupCartesianGrid2D (  )  [protected]

void Trilinos_Util::CrsMatrixGallery::SetupCartesianGrid3D (  )  [protected]

int Trilinos_Util::CrsMatrixGallery::WriteMatrix ( const string &  FileName,
const bool  UseSparse = true 
)

Print matrix on file in MATLAB format.

void Trilinos_Util::CrsMatrixGallery::ZeroOutData (  )  [protected]


Friends And Related Function Documentation

ostream& operator<< ( ostream &  os,
const Trilinos_Util::CrsMatrixGallery G 
) [friend]

Print out detailed information about the problem at hand.


Member Data Documentation

double Trilinos_Util::CrsMatrixGallery::a_ [protected]

double Trilinos_Util::CrsMatrixGallery::alpha_ [protected]

double Trilinos_Util::CrsMatrixGallery::b_ [protected]

double Trilinos_Util::CrsMatrixGallery::beta_ [protected]

double Trilinos_Util::CrsMatrixGallery::c_ [protected]

const Epetra_Comm* Trilinos_Util::CrsMatrixGallery::comm_ [protected]

bool Trilinos_Util::CrsMatrixGallery::ContiguousMap_ [protected]

double Trilinos_Util::CrsMatrixGallery::conv_ [protected]

double Trilinos_Util::CrsMatrixGallery::d_ [protected]

double Trilinos_Util::CrsMatrixGallery::delta_ [protected]

double Trilinos_Util::CrsMatrixGallery::diff_ [protected]

double Trilinos_Util::CrsMatrixGallery::e_ [protected]

double Trilinos_Util::CrsMatrixGallery::epsilon_ [protected]

string Trilinos_Util::CrsMatrixGallery::ErrorMsg [protected]

Epetra_MultiVector* Trilinos_Util::CrsMatrixGallery::ExactSolution_ [protected]

string Trilinos_Util::CrsMatrixGallery::ExactSolutionType_ [protected]

string Trilinos_Util::CrsMatrixGallery::ExpandType_ [protected]

double Trilinos_Util::CrsMatrixGallery::f_ [protected]

string Trilinos_Util::CrsMatrixGallery::FileName_ [protected]

double Trilinos_Util::CrsMatrixGallery::g_ [protected]

double Trilinos_Util::CrsMatrixGallery::gamma_ [protected]

Epetra_LinearProblem* Trilinos_Util::CrsMatrixGallery::LinearProblem_ [protected]

double Trilinos_Util::CrsMatrixGallery::lx_ [protected]

double Trilinos_Util::CrsMatrixGallery::ly_ [protected]

double Trilinos_Util::CrsMatrixGallery::lz_ [protected]

Epetra_Map* Trilinos_Util::CrsMatrixGallery::map_ [protected]

std::vector<int> Trilinos_Util::CrsMatrixGallery::MapMap_ [protected]

string Trilinos_Util::CrsMatrixGallery::MapType_ [protected]

Epetra_CrsMatrix* Trilinos_Util::CrsMatrixGallery::matrix_ [protected]

int Trilinos_Util::CrsMatrixGallery::mx_ [protected]

int Trilinos_Util::CrsMatrixGallery::my_ [protected]

int* Trilinos_Util::CrsMatrixGallery::MyGlobalElements_ [protected]

int Trilinos_Util::CrsMatrixGallery::mz_ [protected]

string Trilinos_Util::CrsMatrixGallery::name_ [protected]

int Trilinos_Util::CrsMatrixGallery::NumGlobalElements_ [protected]

int Trilinos_Util::CrsMatrixGallery::NumMyElements_ [protected]

int Trilinos_Util::CrsMatrixGallery::NumPDEEqns_ [protected]

int Trilinos_Util::CrsMatrixGallery::NumVectors_ [protected]

int Trilinos_Util::CrsMatrixGallery::nx_ [protected]

int Trilinos_Util::CrsMatrixGallery::ny_ [protected]

int Trilinos_Util::CrsMatrixGallery::nz_ [protected]

string Trilinos_Util::CrsMatrixGallery::OutputMsg [protected]

Epetra_MultiVector* Trilinos_Util::CrsMatrixGallery::rhs_ [protected]

string Trilinos_Util::CrsMatrixGallery::RhsType_ [protected]

double Trilinos_Util::CrsMatrixGallery::source_ [protected]

Epetra_MultiVector* Trilinos_Util::CrsMatrixGallery::StartingSolution_ [protected]

string Trilinos_Util::CrsMatrixGallery::StartingSolutionType_ [protected]

Epetra_Vector* Trilinos_Util::CrsMatrixGallery::VectorA_ [protected]

Epetra_Vector * Trilinos_Util::CrsMatrixGallery::VectorB_ [protected]

Epetra_Vector * Trilinos_Util::CrsMatrixGallery::VectorC_ [protected]

Epetra_Vector * Trilinos_Util::CrsMatrixGallery::VectorD_ [protected]

Epetra_Vector * Trilinos_Util::CrsMatrixGallery::VectorE_ [protected]

Epetra_Vector * Trilinos_Util::CrsMatrixGallery::VectorF_ [protected]

Epetra_Vector * Trilinos_Util::CrsMatrixGallery::VectorG_ [protected]

bool Trilinos_Util::CrsMatrixGallery::verbose_ [protected]


The documentation for this class was generated from the following files:
Generated on Wed May 12 21:30:34 2010 for TriUtils by  doxygen 1.4.7