Commit 406385c8 authored by Igor Konshin's avatar Igor Konshin
Browse files

Documentation: documented several header files

Added doxygen-style comments to common, partitione and solver headers
parent 831bfd7c
...@@ -108,12 +108,24 @@ ...@@ -108,12 +108,24 @@
#define INMOST_MPI_DATA_BIG_ENUM_TYPE INMOST_MPI_UNSIGNED #define INMOST_MPI_DATA_BIG_ENUM_TYPE INMOST_MPI_UNSIGNED
/// Cross-platform timer that return current time in seconds.
/// The timer is similar to MPI_Wtime() and omp_get_wtime() but is independent on both flags USE_MPI and USE_OMP.
///
/// double seconds = Timer(); @n
/// ...some code... @n
/// std::cout << "Time spent is " << seconds << " seconds." << std::endl;
long double Timer(); long double Timer();
namespace INMOST namespace INMOST
{ {
/// Types of errors may occur in INMOST.
/// All of these error are fatal ones.
/// If error is detected then "throw" exception is generated.
/// The names of the error type are very intuitive and self-explained ones.
/// Use "try{...} catch{...}" blocks to implement exception handling.
enum ErrorType enum ErrorType
{ {
/// The list of errors connected to mesh consistency.
Failure=100, Failure=100,
NoTagPosition, NoTagPosition,
WrongDataType, WrongDataType,
...@@ -136,7 +148,7 @@ namespace INMOST ...@@ -136,7 +148,7 @@ namespace INMOST
NoSpaceForMarker, NoSpaceForMarker,
ElementBelongsToNobody, ElementBelongsToNobody,
/// The list of general type errors.
BadFileName, BadFileName,
BadFile, BadFile,
CorruptedIerarchy, CorruptedIerarchy,
...@@ -148,7 +160,7 @@ namespace INMOST ...@@ -148,7 +160,7 @@ namespace INMOST
BadParameter, BadParameter,
TopologyCheckError, TopologyCheckError,
//solver /// The list of errors may occur in the Linear Solver.
ErrorInSolver = 400, ErrorInSolver = 400,
DataCorruptedInSolver, DataCorruptedInSolver,
DifferentCommunicatorInSolver, DifferentCommunicatorInSolver,
...@@ -158,11 +170,12 @@ namespace INMOST ...@@ -158,11 +170,12 @@ namespace INMOST
PrepareMatrixFirst, PrepareMatrixFirst,
CannotReusePreconditionerOfDifferentSize, CannotReusePreconditionerOfDifferentSize,
//partitioner /// The list of errors may occur in the Partitioner.
ErrorInPartitioner = 500, ErrorInPartitioner = 500,
UnknownWeightSize, UnknownWeightSize,
DistributionTagWasNotFilled, DistributionTagWasNotFilled,
/// The very tail of the errors list.
NotImplemented = 1000, NotImplemented = 1000,
Impossible Impossible
}; };
......
...@@ -4,7 +4,7 @@ ...@@ -4,7 +4,7 @@
#include "inmost.h" #include "inmost.h"
#if defined(USE_PARTITIONER) && !defined(USE_MESH) #if defined(USE_PARTITIONER) && !defined(USE_MESH)
#warning "USE_PARITIONER require USE_MESH" #warning "USE_PARTITIONER require USE_MESH"
#undef USE_PARTITIONER #undef USE_PARTITIONER
#endif #endif
...@@ -13,11 +13,30 @@ ...@@ -13,11 +13,30 @@
namespace INMOST namespace INMOST
{ {
/// Main class to modify or improve the mesh distribution for better load balancing.
class Partitioner class Partitioner
{ {
public: public:
enum Type{Zoltan_Parmetis, Zoltan_Scotch, Zoltan_PHG, Zoltan_RCB, Zoltan_RIB, Zoltan_HSFC, Parmetis, Inner_RCM}; /// Type of the Partitioner can be currently used in this version of INMOST.
enum Action{Partition, Repartition, Refine}; /// @see Zoltan: Parallel Partitioning, Load Balancing and Data-Management Services. http://www.cs.sandia.gov/Zoltan/
/// @see ParMETIS: Parallel Graph Partitioning and Fill-reducing Matrix Ordering. http://glaros.dtc.umn.edu/gkhome/metis/parmetis/overview
enum Type
{
Zoltan_Parmetis, ///< Parmetis partitioner with the Zoltan package interface.
Zoltan_Scotch, ///< Scotch partitioner with the Zoltan package interface.
Zoltan_PHG, ///< Zoltan topology-based method using Partitioning of HyperGraph.
Zoltan_RCB, ///< Zoltan geometry-based method using Recursive Coordinate Bisection.
Zoltan_RIB, ///< Zoltan geometry-based method using Recursive Inertial Bisection.
Zoltan_HSFC, ///< Zoltan geometry-based method using Hilbert Space-Filling Curve partitioning.
Parmetis, ///< Parmetis partitioner with the original interface.
Inner_RCM ///< Internal serial only partitioner based on the Reverse Cuthill–McKee algorithm ordering.
};
enum Action
{
Partition, ///< Partition "from scratch", not taking into account the current mesh distribution.
Repartition, ///< Repartition the existing partition but try to stay close to the current mesh distribution.
Refine ///< Refine the current partition assuming only small changes of mesh distribution.
};
private: private:
enum Type pt; enum Type pt;
enum Action pa; enum Action pa;
...@@ -25,18 +44,37 @@ namespace INMOST ...@@ -25,18 +44,37 @@ namespace INMOST
Tag weight_tag; Tag weight_tag;
Mesh * m; Mesh * m;
public: public:
/// Initialize the use of partitioner.
/// @param argc The number of arguments transmitted to the function main.
/// @param argv The pointer to arguments transmitted to the function main.
/// The shortest call to this function with the default solver parameters is the following: Initialize(NULL,NULL);
/// @see Partitioner::SetMethod
/// @see Partitioner::Finalize
static void Initialize(int * argc, char *** argv); static void Initialize(int * argc, char *** argv);
/// Finalize the use of partitioner.
/// @see Partitioner::Initialize
static void Finalize(); static void Finalize();
/// The default constructor of the partitioner for the specified mesh.
Partitioner(Mesh * m); Partitioner(Mesh * m);
Partitioner(const Partitioner & other); Partitioner(const Partitioner & other);
Partitioner & operator =(Partitioner const & other); Partitioner & operator =(Partitioner const & other);
~Partitioner(); ~Partitioner();
/// Evaluate the earlier specified partitioner.
/// @see Partitioner::SetMethod
void Evaluate(); void Evaluate();
/// Set the partitioner method to be used.
/// @param t The concrete Type of the partitioner from the selected package.
/// @param a The partitioner Action, the default is Repartition.
/// @see Partitioner::Evaluate
void SetMethod(enum Type t, enum Action a = Repartition); void SetMethod(enum Type t, enum Action a = Repartition);
/// Compute the specific weights for the selected partitioner.
void SetWeight(Tag weight); void SetWeight(Tag weight);
/// Reset the computed weights for the partitioner.
void ResetWeight(); void ResetWeight();
/// Get the Mesh pointer for the current partitioner.
Mesh * GetMesh(); Mesh * GetMesh();
/// Get the Tag of the computed weights for the current partitioner.
Tag GetWeight(); Tag GetWeight();
}; };
} }
......
...@@ -9,17 +9,26 @@ ...@@ -9,17 +9,26 @@
#if defined(USE_SOLVER) #if defined(USE_SOLVER)
namespace INMOST namespace INMOST
{ {
/// Main class to set and solve linear system.
/// Solver class is used to set the coefficient Matrix, the right-hand side Vector
/// and the initial guess Vector, construct the preconditioner and Solve
/// the linear system.
/// Formally, Solver class is independent of INMOST::Mesh class.
/// @see Solver::Matrix
/// @see Solver::Vector
/// @see Solver::Solve
class Solver class Solver
{ {
private: private:
static INMOST_MPI_Type RowEntryType; //prepared in Initialize static INMOST_MPI_Type RowEntryType; //prepared in Initialize
public: public:
/// Type of the Solver can be currently used in this version of INMOST.
enum Type enum Type
{ {
INNER_ILU2, INNER_ILU2, ///< inner Solver based on second order ILU factorization.
INNER_MLILUC, INNER_MLILUC, ///< inner Solver based on Saad multilevel ILU with pivoting.
PETSC, PETSC, ///< external Solver PETSc, @see http://www.mcs.anl.gov/petsc/.
ANI //sequential, fortran ANI ///< external Solver from ANI3D based on ILU2 (sequential Fortran version).
}; };
static INMOST_MPI_Type & GetRowEntryType() {return RowEntryType;} static INMOST_MPI_Type & GetRowEntryType() {return RowEntryType;}
...@@ -29,6 +38,7 @@ namespace INMOST ...@@ -29,6 +38,7 @@ namespace INMOST
class Matrix; class Matrix;
class Vector; class Vector;
/// Base class for low level operations with objects of Solver class.
class OrderInfo class OrderInfo
{ {
private: private:
...@@ -50,24 +60,42 @@ namespace INMOST ...@@ -50,24 +60,42 @@ namespace INMOST
int rank,size; int rank,size;
public: public:
void Clear(); void Clear();
/// Return true if Matrix data have already been specified.
bool & HaveMatrix() { return have_matrix; } bool & HaveMatrix() { return have_matrix; }
OrderInfo(); OrderInfo();
OrderInfo(const OrderInfo & other); OrderInfo(const OrderInfo & other);
OrderInfo & operator =(OrderInfo const & other); OrderInfo & operator =(OrderInfo const & other);
~OrderInfo(); ~OrderInfo();
/// Prepare parallel state of the Matrix with specified overlap size.
/// This state of the matrix can be used, for instance, to construct
/// the preconditioner for Additive Swartz method.
/// @param m Matrix to be expanded.
/// @param overlap Overlap size, viz. the number of overlap layers.
void PrepareMatrix(Matrix & m, INMOST_DATA_ENUM_TYPE overlap); void PrepareMatrix(Matrix & m, INMOST_DATA_ENUM_TYPE overlap);
/// Restore initial nonparallel state of the Matrix with no overlap.
void RestoreMatrix(Matrix & m); void RestoreMatrix(Matrix & m);
/// Prepare parallel state of the Vector.
void PrepareVector(Vector & v); void PrepareVector(Vector & v);
/// Restore initial nonparallel state of the Vector.
void RestoreVector(Vector & v); //Restore initial nonparallel state of the vector void RestoreVector(Vector & v); //Restore initial nonparallel state of the vector
/// Retrieve the processor number by binary search for the specified global index.
INMOST_DATA_ENUM_TYPE GetProcessor(INMOST_DATA_ENUM_TYPE gind) const; //retrive processor by binary search in global_to_proc INMOST_DATA_ENUM_TYPE GetProcessor(INMOST_DATA_ENUM_TYPE gind) const; //retrive processor by binary search in global_to_proc
void GetOverlapRegion(INMOST_DATA_ENUM_TYPE proc, INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const; void GetOverlapRegion(INMOST_DATA_ENUM_TYPE proc, INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const;
/// Get the local index region for the specified process.
void GetLocalRegion(INMOST_DATA_ENUM_TYPE proc, INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const; void GetLocalRegion(INMOST_DATA_ENUM_TYPE proc, INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const;
/// Get the local index region for the current process.
void GetVectorRegion(INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const {mbeg = local_vector_begin; mend = local_vector_end;} void GetVectorRegion(INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const {mbeg = local_vector_begin; mend = local_vector_end;}
/// Get the rank of the current communicator, i.e. the current process index.
INMOST_DATA_ENUM_TYPE GetRank() const {return rank;} INMOST_DATA_ENUM_TYPE GetRank() const {return rank;}
/// Get the size of the current communicator, i.e. the total number of processes used.
INMOST_DATA_ENUM_TYPE GetSize() const {return size;} INMOST_DATA_ENUM_TYPE GetSize() const {return size;}
/// Update the shared data in parallel vector.
void Update (Vector & x); // update parallel vector void Update (Vector & x); // update parallel vector
/// Sum shared values in parallel vector.
void Accumulate(Vector & x); // sum shared values in parallel vector void Accumulate(Vector & x); // sum shared values in parallel vector
/// Get the sum of num elements of real array on all processes.
void Integrate(INMOST_DATA_REAL_TYPE * inout, INMOST_DATA_ENUM_TYPE num); void Integrate(INMOST_DATA_REAL_TYPE * inout, INMOST_DATA_ENUM_TYPE num);
/// Get the communicator which the solver is associated with.
INMOST_MPI_Comm GetComm() {return comm;} INMOST_MPI_Comm GetComm() {return comm;}
// Access to arrays below allows to organize manual exchange // Access to arrays below allows to organize manual exchange
INMOST_MPI_Request * GetSendRequests() {assert(!send_requests.empty()); return &send_requests[0];} INMOST_MPI_Request * GetSendRequests() {assert(!send_requests.empty()); return &send_requests[0];}
...@@ -83,6 +111,10 @@ namespace INMOST ...@@ -83,6 +111,10 @@ namespace INMOST
//~ void EndSequentialCode() {for(int i = rank; i < size; i++) MPI_Barrier(comm);} //~ void EndSequentialCode() {for(int i = rank; i < size; i++) MPI_Barrier(comm);}
}; };
/// Distributed vector class.
/// This class can be used to store both local and distributed dense data of real type.
/// For example, to form the right-hand side or initial guess to the solution.
/// @see Solver::Solve
class Vector class Vector
{ {
public: public:
...@@ -95,47 +127,68 @@ namespace INMOST ...@@ -95,47 +127,68 @@ namespace INMOST
std::string name; std::string name;
bool is_parallel; bool is_parallel;
public: public:
/// Main constructor of the Vector class.
/// @param _name Name of the vector, empty string by default.
/// @param start Start of the local data interval.
/// @param end End of the local data interval.
/// @param _comm Communicator for parallel data exchanges, MPI_COMM_WORLD by default.
Vector(std::string _name = "", INMOST_DATA_ENUM_TYPE start = 0, INMOST_DATA_ENUM_TYPE end = 0, INMOST_MPI_Comm _comm = INMOST_MPI_COMM_WORLD); Vector(std::string _name = "", INMOST_DATA_ENUM_TYPE start = 0, INMOST_DATA_ENUM_TYPE end = 0, INMOST_MPI_Comm _comm = INMOST_MPI_COMM_WORLD);
Vector(const Vector & other); Vector(const Vector & other);
Vector & operator =(Vector const & other); Vector & operator =(Vector const & other);
~Vector(); ~Vector();
/// Return reference to i-th element of the vector.
INMOST_DATA_REAL_TYPE & operator [](INMOST_DATA_ENUM_TYPE i) {return data[i];} INMOST_DATA_REAL_TYPE & operator [](INMOST_DATA_ENUM_TYPE i) {return data[i];}
/// Return i-th element of the vector.
INMOST_DATA_REAL_TYPE operator [](INMOST_DATA_ENUM_TYPE i) const {return data[i];} INMOST_DATA_REAL_TYPE operator [](INMOST_DATA_ENUM_TYPE i) const {return data[i];}
/// Return the global size of the vector.
INMOST_DATA_ENUM_TYPE Size() const { return static_cast<INMOST_DATA_ENUM_TYPE>(data.size()); } INMOST_DATA_ENUM_TYPE Size() const { return static_cast<INMOST_DATA_ENUM_TYPE>(data.size()); }
iterator Begin() {return data.begin();} iterator Begin() {return data.begin();}
const_iterator Begin() const {return data.begin();} const_iterator Begin() const {return data.begin();}
iterator End() {return data.end();} iterator End() {return data.end();}
const_iterator End() const {return data.end();} const_iterator End() const {return data.end();}
bool Empty() const {return data.empty();} bool Empty() const {return data.empty();}
/// Set the start and the end of the distributed vector interval.
void SetInterval(INMOST_DATA_ENUM_TYPE start, INMOST_DATA_ENUM_TYPE end) {data.set_interval_beg(start); data.set_interval_end(end);} void SetInterval(INMOST_DATA_ENUM_TYPE start, INMOST_DATA_ENUM_TYPE end) {data.set_interval_beg(start); data.set_interval_end(end);}
/// Get the start and the end of the distributed vector interval.
void GetInterval(INMOST_DATA_ENUM_TYPE & start, INMOST_DATA_ENUM_TYPE & end) const {start = data.get_interval_beg(); end = data.get_interval_end();} void GetInterval(INMOST_DATA_ENUM_TYPE & start, INMOST_DATA_ENUM_TYPE & end) const {start = data.get_interval_beg(); end = data.get_interval_end();}
void ShiftInterval(INMOST_DATA_ENUM_TYPE shift) {data.shift_interval(shift);} void ShiftInterval(INMOST_DATA_ENUM_TYPE shift) {data.shift_interval(shift);}
/// Get the first index of the distributed vector interval.
INMOST_DATA_ENUM_TYPE GetFirstIndex() const {return data.get_interval_beg();} INMOST_DATA_ENUM_TYPE GetFirstIndex() const {return data.get_interval_beg();}
/// Get the communicator which the vector is associated with.
INMOST_MPI_Comm GetCommunicator() const {return comm;} INMOST_MPI_Comm GetCommunicator() const {return comm;}
/// Get the scalar product of the specified interval of the distributed vector.
INMOST_DATA_REAL_TYPE ScalarProd(Vector const & other, INMOST_DATA_ENUM_TYPE index_begin, INMOST_DATA_ENUM_TYPE index_end) const; INMOST_DATA_REAL_TYPE ScalarProd(Vector const & other, INMOST_DATA_ENUM_TYPE index_begin, INMOST_DATA_ENUM_TYPE index_end) const;
/// Save the distributed vector to a single data file using parallel MPI I/O.
void Save(std::string file); void Save(std::string file);
/// Load the vector from a single data file using the specified interval.
/// If interval is not specified, then it will be automatically constructed,
/// with the about equal block size (the last block may has larger dimension).
void Load(std::string file, INMOST_DATA_ENUM_TYPE mbeg = ENUMUNDEF, INMOST_DATA_ENUM_TYPE mend = ENUMUNDEF); void Load(std::string file, INMOST_DATA_ENUM_TYPE mbeg = ENUMUNDEF, INMOST_DATA_ENUM_TYPE mend = ENUMUNDEF);
bool & isParallel() {return is_parallel;} bool & isParallel() {return is_parallel;}
/// Get the vector name specified in the main constructor.
std::string GetName() {return name;} std::string GetName() {return name;}
/// Clear all data of the current vector.
void Clear() {data.clear();} void Clear() {data.clear();}
//~ friend class Solver; //~ friend class Solver;
}; };
/// Class to store the sparse matrix row.
class Row class Row
{ {
public: public:
/// Entry of the sparse matrix row.
typedef struct entry_s typedef struct entry_s
{ {
INMOST_DATA_ENUM_TYPE first; INMOST_DATA_ENUM_TYPE first; ///< the column number of the row element.
INMOST_DATA_REAL_TYPE second; INMOST_DATA_REAL_TYPE second; ///< the real value of the row element.
entry_s() :first(0), second(0.0) {} entry_s() :first(0), second(0.0) {}
entry_s(const entry_s & other) :first(other.first), second(other.second) {}//{std::cout << __FUNCTION__ << " " << first << " " << second << std::endl;} entry_s(const entry_s & other) :first(other.first), second(other.second) {}//{std::cout << __FUNCTION__ << " " << first << " " << second << std::endl;}
entry_s(INMOST_DATA_ENUM_TYPE first, INMOST_DATA_REAL_TYPE second):first(first),second(second){} entry_s(INMOST_DATA_ENUM_TYPE first, INMOST_DATA_REAL_TYPE second):first(first),second(second){}
...@@ -203,6 +256,7 @@ namespace INMOST ...@@ -203,6 +256,7 @@ namespace INMOST
} }
~Row() {} ~Row() {}
Row & operator = (Row const & other) { data = other.data; marker = other.marker; return *this; } Row & operator = (Row const & other) { data = other.data; marker = other.marker; return *this; }
/// The operator [] used to fill the sparse matrix row, but not to access individual elements of the row.
INMOST_DATA_REAL_TYPE & operator [](INMOST_DATA_ENUM_TYPE i) // use to fill matrix, not to access individual elements INMOST_DATA_REAL_TYPE & operator [](INMOST_DATA_ENUM_TYPE i) // use to fill matrix, not to access individual elements
{ {
//for sparse_data type //for sparse_data type
...@@ -219,6 +273,7 @@ namespace INMOST ...@@ -219,6 +273,7 @@ namespace INMOST
return data.back().second; return data.back().second;
} }
/// The operator [] used to access individual elements of the row.
INMOST_DATA_REAL_TYPE operator[](INMOST_DATA_ENUM_TYPE i) const INMOST_DATA_REAL_TYPE operator[](INMOST_DATA_ENUM_TYPE i) const
{ {
//for sparse data type //for sparse data type
...@@ -231,8 +286,10 @@ namespace INMOST ...@@ -231,8 +286,10 @@ namespace INMOST
return 1.0e20; return 1.0e20;
} }
//void Reserve(INMOST_DATA_ENUM_TYPE num) { data.reserve(num);} //void Reserve(INMOST_DATA_ENUM_TYPE num) { data.reserve(num);}
/// Clear all data of the current row.
void Clear() { data.clear(); } void Clear() { data.clear(); }
void Swap(Solver::Row & other) { data.swap(other.data); bool tmp = marker; marker = other.marker; other.marker = tmp; } void Swap(Solver::Row & other) { data.swap(other.data); bool tmp = marker; marker = other.marker; other.marker = tmp; }
/// The size of the sparse row, i.e. the total number of nonzero elements.
INMOST_DATA_ENUM_TYPE Size() const { return static_cast<INMOST_DATA_ENUM_TYPE>(data.size()); } INMOST_DATA_ENUM_TYPE Size() const { return static_cast<INMOST_DATA_ENUM_TYPE>(data.size()); }
INMOST_DATA_ENUM_TYPE & GetIndex(INMOST_DATA_ENUM_TYPE k) {assert(k < data.size()); return (data.begin()+k)->first;} INMOST_DATA_ENUM_TYPE & GetIndex(INMOST_DATA_ENUM_TYPE k) {assert(k < data.size()); return (data.begin()+k)->first;}
INMOST_DATA_REAL_TYPE & GetValue(INMOST_DATA_ENUM_TYPE k) {assert(k < data.size()); return (data.begin()+k)->second;} INMOST_DATA_REAL_TYPE & GetValue(INMOST_DATA_ENUM_TYPE k) {assert(k < data.size()); return (data.begin()+k)->second;}
...@@ -247,13 +304,20 @@ namespace INMOST ...@@ -247,13 +304,20 @@ namespace INMOST
reverse_iterator rEnd() { return data.rend(); } reverse_iterator rEnd() { return data.rend(); }
const_reverse_iterator rBegin() const { return data.rbegin(); } const_reverse_iterator rBegin() const { return data.rbegin(); }
const_reverse_iterator rEnd() const { return data.rend(); } const_reverse_iterator rEnd() const { return data.rend(); }
/// Return the scalar product of the current sparse row by a dense Vector.
INMOST_DATA_REAL_TYPE RowVec(Vector & x) const; // returns A(row) * x INMOST_DATA_REAL_TYPE RowVec(Vector & x) const; // returns A(row) * x
void MoveRow(Row & new_pos) {data = new_pos.data;} //here move constructor and std::move may be used in future void MoveRow(Row & new_pos) {data = new_pos.data;} //here move constructor and std::move may be used in future
/// Set the vector entries by zeroes.
void Zero() {for(iterator it = Begin(); it != End(); ++it) it->second = 0;} void Zero() {for(iterator it = Begin(); it != End(); ++it) it->second = 0;}
//! Use this function if you know that this index is not repeated in the row (and indexes are inserted in increasing order (for sparse_data)) /// Push specified element into sparse row.
/// This function should be used only if the index is not repeated in the row
/// and all indexes are inserted in increasing order.
void Push(INMOST_DATA_ENUM_TYPE ind, INMOST_DATA_REAL_TYPE val) {data.push_back(entry(ind,val));} void Push(INMOST_DATA_ENUM_TYPE ind, INMOST_DATA_REAL_TYPE val) {data.push_back(entry(ind,val));}
}; };
/// Class to store the distributed sparse matrix by compressed rows.
/// The format used to store sparse matrix is analogous to Compressed Row Storage format (CRS).
/// @see http://netlib.org/linalg/html_templates/node91.html
class Matrix class Matrix
{ {
public: public:
...@@ -266,23 +330,35 @@ namespace INMOST ...@@ -266,23 +330,35 @@ namespace INMOST
std::string name; std::string name;
bool is_parallel; bool is_parallel;
public: public:
/// Main constructor of the Matrix class.
/// @param _name Name of the matrix, empty string by default.
/// @param start Start of the local data interval.
/// @param end End of the local data interval.
/// @param _comm Communicator for parallel data exchanges, MPI_COMM_WORLD by default.
Matrix(std::string _name = "", INMOST_DATA_ENUM_TYPE start = 0, INMOST_DATA_ENUM_TYPE end = 0, INMOST_MPI_Comm _comm = INMOST_MPI_COMM_WORLD); Matrix(std::string _name = "", INMOST_DATA_ENUM_TYPE start = 0, INMOST_DATA_ENUM_TYPE end = 0, INMOST_MPI_Comm _comm = INMOST_MPI_COMM_WORLD);
Matrix(const Matrix & other); Matrix(const Matrix & other);
Matrix & operator =(Matrix const & other); Matrix & operator =(Matrix const & other);
~Matrix(); ~Matrix();
/// Return reference to i-th Row of the matrix.
Row & operator [](INMOST_DATA_ENUM_TYPE i) {return data[i];} Row & operator [](INMOST_DATA_ENUM_TYPE i) {return data[i];}
/// Return reference to i-th Row of the matrix.
const Row & operator [](INMOST_DATA_ENUM_TYPE i) const {return data[i];} const Row & operator [](INMOST_DATA_ENUM_TYPE i) const {return data[i];}
/// Return the total number of rows in the matrix.
INMOST_DATA_ENUM_TYPE Size() const { return static_cast<INMOST_DATA_ENUM_TYPE>(data.size()); } INMOST_DATA_ENUM_TYPE Size() const { return static_cast<INMOST_DATA_ENUM_TYPE>(data.size()); }
bool Empty() const {return data.empty();} bool Empty() const {return data.empty();}
iterator Begin() {return data.begin();} iterator Begin() {return data.begin();}
iterator End() {return data.end();} iterator End() {return data.end();}
const_iterator Begin() const {return data.begin();} const_iterator Begin() const {return data.begin();}
const_iterator End() const {return data.end();} const_iterator End() const {return data.end();}
/// Set the start and the end row numbers of the distributed matrix interval.
void SetInterval(INMOST_DATA_ENUM_TYPE start, INMOST_DATA_ENUM_TYPE end) {data.set_interval_beg(start); data.set_interval_end(end);} void SetInterval(INMOST_DATA_ENUM_TYPE start, INMOST_DATA_ENUM_TYPE end) {data.set_interval_beg(start); data.set_interval_end(end);}
/// Get the start and the end row numbers of the distributed matrix interval.
void GetInterval(INMOST_DATA_ENUM_TYPE & start, INMOST_DATA_ENUM_TYPE & end) const {start = data.get_interval_beg(); end = data.get_interval_end();} void GetInterval(INMOST_DATA_ENUM_TYPE & start, INMOST_DATA_ENUM_TYPE & end) const {start = data.get_interval_beg(); end = data.get_interval_end();}
void ShiftInterval(INMOST_DATA_ENUM_TYPE shift) {data.shift_interval(shift);} void ShiftInterval(INMOST_DATA_ENUM_TYPE shift) {data.shift_interval(shift);}
/// Get the first row index of the distributed matrix interval.
INMOST_DATA_ENUM_TYPE GetFirstIndex() const {return data.get_interval_beg();} INMOST_DATA_ENUM_TYPE GetFirstIndex() const {return data.get_interval_beg();}
/// Get the communicator which the matrix is associated with.
INMOST_MPI_Comm GetCommunicator() const {return comm;} INMOST_MPI_Comm GetCommunicator() const {return comm;}
void MoveRows(INMOST_DATA_ENUM_TYPE from, INMOST_DATA_ENUM_TYPE to, INMOST_DATA_ENUM_TYPE size); //for parallel void MoveRows(INMOST_DATA_ENUM_TYPE from, INMOST_DATA_ENUM_TYPE to, INMOST_DATA_ENUM_TYPE size); //for parallel
void Swap(Solver::Matrix & other) void Swap(Solver::Matrix & other)
...@@ -297,14 +373,24 @@ namespace INMOST ...@@ -297,14 +373,24 @@ namespace INMOST
other.is_parallel = ptmp; other.is_parallel = ptmp;
} }
void MatVec(INMOST_DATA_REAL_TYPE alpha, Solver::Vector & x, INMOST_DATA_REAL_TYPE beta, Solver::Vector & out) const; //y = alpha*A*x + beta * y /// Matrix-vector product of the form: y = alpha*A*x + beta * y.
/// @param y Input/output vector.
/// @see Solver::Vector::Zero
void MatVec(INMOST_DATA_REAL_TYPE alpha, Solver::Vector & x, INMOST_DATA_REAL_TYPE beta, Solver::Vector & y) const; //y = alpha*A*x + beta * y
/// Clear all data of the matrix.
void Clear() {for(Matrix::iterator it = Begin(); it != End(); ++it) it->Clear(); data.clear();} void Clear() {for(Matrix::iterator it = Begin(); it != End(); ++it) it->Clear(); data.clear();}
/// Load the matrix from a single data file in MTX format using the specified interval.
/// If interval is not specified, then it will be automatically constructed,
/// with the about equal block size (the last block may has larger dimension).
void Load(std::string file, INMOST_DATA_ENUM_TYPE beg = ENUMUNDEF, INMOST_DATA_ENUM_TYPE end = ENUMUNDEF); void Load(std::string file, INMOST_DATA_ENUM_TYPE beg = ENUMUNDEF, INMOST_DATA_ENUM_TYPE end = ENUMUNDEF);
/// Save the distributed matrix to a single data file in MTX format using parallel MPI I/O.
/// @see http://math.nist.gov/MatrixMarket/formats.html
void Save(std::string file); void Save(std::string file);
bool & isParallel() { return is_parallel; } bool & isParallel() { return is_parallel; }
/// Get the matrix name specified in the main constructor.
std::string GetName() {return name;} std::string GetName() {return name;}
//~ friend class Solver; //~ friend class Solver;
}; };
...@@ -330,25 +416,67 @@ namespace INMOST ...@@ -330,25 +416,67 @@ namespace INMOST
Solver(const Solver & other);// prohibit copy Solver(const Solver & other);// prohibit copy
Solver & operator =(Solver const & other); //prohibit assignment Solver & operator =(Solver const & other); //prohibit assignment
public: public:
/// Set the solver parameter of the integer type.
void SetParameterEnum(std::string name, INMOST_DATA_ENUM_TYPE value); void SetParameterEnum(std::string name, INMOST_DATA_ENUM_TYPE value);
/// Set the solver parameter of the real type.
void SetParameterReal(std::string name, INMOST_DATA_REAL_TYPE value); void SetParameterReal(std::string name, INMOST_DATA_REAL_TYPE value);
/// Get the used defined name of the Solver.
std::string GetName() {return name;} std::string GetName() {return name;}
/// Get the package Type.
Type GetPackage() const {return _pack;} Type GetPackage() const {return _pack;}
/// Return the number of iterations performed by the last solution.
/// @see Solver::Solve
INMOST_DATA_ENUM_TYPE Iterations(); INMOST_DATA_ENUM_TYPE Iterations();
/// Return the final residual achieved by the last solution.
/// @see Solver::Solve
INMOST_DATA_REAL_TYPE Residual(); INMOST_DATA_REAL_TYPE Residual();
/// Set the matrix and construct the preconditioner.
/// @param OldPreconditioner If this parameter is set to true,
/// then the previous preconditioner will be used,