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

Epetra_LinearProblemRedistor: A class for redistributing an Epetra_LinearProblem object. More...

#include <Epetra_LinearProblemRedistor.h>

Private Member Functions

int GenerateRedistMap ()
 

Private Attributes

Epetra_LinearProblemOrigProblem_
 
int NumProc_
 
Epetra_LinearProblemRedistProblem_
 
Epetra_MapRedistMap_
 
Epetra_RowMatrixTransposerTransposer_
 
Epetra_ExportRedistExporter_
 
bool Replicate_
 
bool ConstructTranspose_
 
bool MakeDataContiguous_
 
bool MapGenerated_
 
bool RedistProblemCreated_
 
int * ptr_
 

Constructors/destructors

 Epetra_LinearProblemRedistor (Epetra_LinearProblem *OrigProblem, const Epetra_Map &RedistMap)
 Epetra_LinearProblemRedistor constructor using pre-defined layout. More...
 
 Epetra_LinearProblemRedistor (Epetra_LinearProblem *OrigProblem, int NumProc, bool Replicate)
 Epetra_LinearProblemRedistor constructor specifying number of processor and replication bool. More...
 
 Epetra_LinearProblemRedistor (const Epetra_LinearProblemRedistor &Source)
 Epetra_LinearProblemRedistor copy constructor. More...
 
virtual ~Epetra_LinearProblemRedistor ()
 Epetra_LinearProblemRedistor destructor. More...
 

Forward transformation methods

int CreateRedistProblem (const bool ConstructTranspose, const bool MakeDataContiguous, Epetra_LinearProblem *&RedistProblem)
 Generate a new Epetra_LinearProblem as a redistribution of the one passed into the constructor. More...
 
int UpdateRedistProblemValues (Epetra_LinearProblem *ProblemWithNewValues)
 Update the values of an already-redistributed problem. More...
 
int UpdateRedistRHS (Epetra_MultiVector *RHSWithNewValues)
 Update the values of an already-redistributed RHS. More...
 

Reverse transformation methods

int UpdateOriginalLHS (Epetra_MultiVector *LHS)
 Update LHS of original Linear Problem object. More...
 

Attribute accessor methods

const Epetra_MapRedistMap () const
 Returns const reference to the Epetra_Map that describes the layout of the RedistLinearProblem. More...
 
const Epetra_ExportRedistExporter () const
 Returns const reference to the Epetra_Export object used to redistribute the original linear problem. More...
 

Utility methods

int ExtractHbData (int &M, int &N, int &nz, int *&ptr, int *&ind, double *&val, int &Nrhs, double *&rhs, int &ldrhs, double *&lhs, int &ldlhs) const
 Extract the redistributed problem data in a form usable for other codes that require Harwell-Boeing format. More...
 

I/O methods

virtual void Print (std::ostream &os) const
 Print method. More...
 

Detailed Description

Epetra_LinearProblemRedistor: A class for redistributing an Epetra_LinearProblem object.

This class provides capabilities to redistribute an existing Epetra_LinearProblem object across a parallel distributed memory machine. All or part of a linear problem object can be redistributed. Reverse distributions, value updates and matrix transposition are also supported. Specification of the redistribution can be done by

  1. providing a target Epetra_Map object describing the new distribution, or
  2. by specifying the number of processors to use and stating whether or not the problem should be completely replicated on all processors.

Definition at line 66 of file Epetra_LinearProblemRedistor.h.

Constructor & Destructor Documentation

Epetra_LinearProblemRedistor::Epetra_LinearProblemRedistor ( Epetra_LinearProblem OrigProblem,
const Epetra_Map RedistMap 
)

Epetra_LinearProblemRedistor constructor using pre-defined layout.

Parameters
Problem(In) An existing Epetra_LinearProblem object. The Epetra_RowMatrix, the LHS and RHS pointers do not need to be defined before this constructor is called.
RedistMap(In) An Epetra_Map describing the target layout of the redistribution.
Returns
Pointer to a Epetra_LinearProblemRedistor object.

Definition at line 55 of file Epetra_LinearProblemRedistor.cpp.

Epetra_LinearProblemRedistor::Epetra_LinearProblemRedistor ( Epetra_LinearProblem OrigProblem,
int  NumProc,
bool  Replicate 
)

Epetra_LinearProblemRedistor constructor specifying number of processor and replication bool.

Parameters
Problem(In) An existing Epetra_LinearProblem object. The Epetra_RowMatrix, the LHS and RHS pointers do not need to be defined before this constructor is called.
NumProc(In) Number of processors to use when redistributing the problem. Must be between 1 and the number of processor on the parallel machine.
Replicate(In) A bool that indicates if the linear problem should be fully replicated on all processors. If true, then a complete copy of the linear problem will be made on each processor. If false, then the problem will be roughly evenly spread across the total number of processors.
Returns
Pointer to a Epetra_LinearProblemRedistor object.

Definition at line 71 of file Epetra_LinearProblemRedistor.cpp.

Epetra_LinearProblemRedistor::Epetra_LinearProblemRedistor ( const Epetra_LinearProblemRedistor Source)

Epetra_LinearProblemRedistor copy constructor.

Definition at line 88 of file Epetra_LinearProblemRedistor.cpp.

Epetra_LinearProblemRedistor::~Epetra_LinearProblemRedistor ( )
virtual

Epetra_LinearProblemRedistor destructor.

Definition at line 105 of file Epetra_LinearProblemRedistor.cpp.

Member Function Documentation

int Epetra_LinearProblemRedistor::CreateRedistProblem ( const bool  ConstructTranspose,
const bool  MakeDataContiguous,
Epetra_LinearProblem *&  RedistProblem 
)

Generate a new Epetra_LinearProblem as a redistribution of the one passed into the constructor.

Constructs a new Epetra_LinearProblem that is a copy of the one passed in to the constructor. The new problem will have redistributed copies of the RowMatrix, LHS and RHS from the original problem. If any of these three objects are 0 pointers, then the corresponding pointer will be zero in the redistributed object.

The redistributed matrix will constructed as an Epetra_CrsMatrix. The LHS and RHS will be Epetra_MultiVector objects.

Two bools can be set when calling this method. The first, ConstructTranspose, will cause the Redistribute method to construct the transpose of the original row matrix. The second, MakeDataContiguous, forces the memory layout of the output matrix, RHS and LHS to be stored so that it is compatible with Fortran. In particular, the Epetra_CrsMatrix is stored so that value from row to row are contiguous, as are the indices. This is compatible with the Harwell-Boeing compressed row and compressed column format. The RHS and LHS are created so that there is a constant stride between the columns of the multivector. This is compatible with Fortran 2D array storage.

Parameters
ConstructTranspose(In) Causes the output matrix to be transposed. This feature can be used to support solvers that need the matrix to be stored in column format. This option has no impact on the LHS and RHS of the output problem.
MakeDataContiguous(In) Causes the output matrix, LHS and RHS to be stored in a form compatible with Fortran-style solvers. The output matrix will be compatible with the Harwell-Boeing compressed column format. The RHS and LHS will be stored such that the last value in column j of the multivector is stored next to the first value in column j+1.
RedistProblem(Out) The redistributed Epetra_LinearProblem. The RowMatrix, LHS and RHS that are generated as part of this problem will be destroyed when the Epetra_LinearProblemRedistor object is destroyed.
Returns
Integer error code, 0 if no errors, positive value if one or more of the LHS or RHS pointers were 0. Negative if some other fatal error occured.

Definition at line 178 of file Epetra_LinearProblemRedistor.cpp.

int Epetra_LinearProblemRedistor::UpdateRedistProblemValues ( Epetra_LinearProblem ProblemWithNewValues)

Update the values of an already-redistributed problem.

Updates the values of an already-redistributed problem. This method allows updating the redistributed problem without allocating new storage. All three objects in the RedistProblem will be updated, namely the Matrix, LHS and RHS. If the LHS or RHS are 0 pointers, they will be ignored.

Parameters
ProblemWithNewValues(In) The values from ProblemWithNewValues will be copied into the RedistProblem. The ProblemWithNewValues object must be identical in structure to the Epetra_LinearProblem object used to create this instance of Epetra_LinearProblemRedistor.
Returns
Integer error code, 0 if no errors, positive value if one or more of the input LHS or RHS pointers were 0. Negative if some other fatal error occured.

Definition at line 244 of file Epetra_LinearProblemRedistor.cpp.

int Epetra_LinearProblemRedistor::UpdateRedistRHS ( Epetra_MultiVector RHSWithNewValues)

Update the values of an already-redistributed RHS.

Updates the values of an already-redistributed RHS. This method allows updating the redistributed RHS without allocating new storage. This method updates only the RHS, and no other part of the RedistLinearProblem.

Parameters
RHSWithNewValues(In) The values from RHSWithNewValues will be copied into the RHS of the RedistProblem. The RHSWithNewValues object must be identical in structure to the Epetra_MultiVector object used to create this instance of Epetra_LinearProblemRedistor.
Returns
Integer error code, 0 if no errors.

Definition at line 344 of file Epetra_LinearProblemRedistor.cpp.

int Epetra_LinearProblemRedistor::UpdateOriginalLHS ( Epetra_MultiVector LHS)

Update LHS of original Linear Problem object.

Copies the values from the LHS of the RedistProblem Object into the LHS passed in to the method. If the RedistProblem is replicated, the LHS will be computed as an average from all processor.

Parameters
LHS(Out) On exit, the values in LHS will contain the values from the current RedistLinearProblem LHS. If the GIDs of the RedistMap are not one-to-one, e.g., if the map is replicated, the output values for each GID will be an average of all values at that GID. If the RedistProblem is being solved redundantly in any fashion, this approach to computing the values of LHS should produce a valid answer no matter how the RedistProblem LHS is distributed.
Returns
Error code, returns 0 if no error.

Definition at line 333 of file Epetra_LinearProblemRedistor.cpp.

const Epetra_Map& Epetra_LinearProblemRedistor::RedistMap ( ) const
inline

Returns const reference to the Epetra_Map that describes the layout of the RedistLinearProblem.

The RedistMap object can be used to construct other Epetra_DistObject objects whose maps are compatible with the redistributed linear problem map.

Warning
This method must not be called until after CreateRedistProblem() is called.

Definition at line 197 of file Epetra_LinearProblemRedistor.h.

const Epetra_Export& Epetra_LinearProblemRedistor::RedistExporter ( ) const
inline

Returns const reference to the Epetra_Export object used to redistribute the original linear problem.

The RedistExporter object can be used to redistribute other Epetra_DistObject objects whose maps are compatible with the original linear problem map, or a reverse distribution for objects compatible with the RedistMap().

Warning
This method must not be called until after CreateRedistProblem() is called.

Definition at line 204 of file Epetra_LinearProblemRedistor.h.

int Epetra_LinearProblemRedistor::ExtractHbData ( int &  M,
int &  N,
int &  nz,
int *&  ptr,
int *&  ind,
double *&  val,
int &  Nrhs,
double *&  rhs,
int &  ldrhs,
double *&  lhs,
int &  ldlhs 
) const

Extract the redistributed problem data in a form usable for other codes that require Harwell-Boeing format.

This method extract data from the linear problem for use with other packages, such as SuperLU, that require the matrix, rhs and lhs in Harwell-Boeing format. Note that the arrays returned by this method are owned by the Epetra_LinearProblemRedistor class, and they will be deleted when the owning class is destroyed.

Definition at line 288 of file Epetra_LinearProblemRedistor.cpp.

virtual void Epetra_LinearProblemRedistor::Print ( std::ostream &  os) const
virtual

Print method.

int Epetra_LinearProblemRedistor::GenerateRedistMap ( )
private

Definition at line 123 of file Epetra_LinearProblemRedistor.cpp.

Member Data Documentation

Epetra_LinearProblem* Epetra_LinearProblemRedistor::OrigProblem_
private

Definition at line 230 of file Epetra_LinearProblemRedistor.h.

int Epetra_LinearProblemRedistor::NumProc_
private

Definition at line 231 of file Epetra_LinearProblemRedistor.h.

Epetra_LinearProblem* Epetra_LinearProblemRedistor::RedistProblem_
private

Definition at line 232 of file Epetra_LinearProblemRedistor.h.

Epetra_Map* Epetra_LinearProblemRedistor::RedistMap_
private

Definition at line 233 of file Epetra_LinearProblemRedistor.h.

Epetra_RowMatrixTransposer* Epetra_LinearProblemRedistor::Transposer_
private

Definition at line 234 of file Epetra_LinearProblemRedistor.h.

Epetra_Export* Epetra_LinearProblemRedistor::RedistExporter_
private

Definition at line 235 of file Epetra_LinearProblemRedistor.h.

bool Epetra_LinearProblemRedistor::Replicate_
private

Definition at line 237 of file Epetra_LinearProblemRedistor.h.

bool Epetra_LinearProblemRedistor::ConstructTranspose_
private

Definition at line 238 of file Epetra_LinearProblemRedistor.h.

bool Epetra_LinearProblemRedistor::MakeDataContiguous_
private

Definition at line 239 of file Epetra_LinearProblemRedistor.h.

bool Epetra_LinearProblemRedistor::MapGenerated_
private

Definition at line 240 of file Epetra_LinearProblemRedistor.h.

bool Epetra_LinearProblemRedistor::RedistProblemCreated_
private

Definition at line 241 of file Epetra_LinearProblemRedistor.h.

int* Epetra_LinearProblemRedistor::ptr_
mutableprivate

Definition at line 243 of file Epetra_LinearProblemRedistor.h.


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