SG++-Doxygen-Documentation
Loading...
Searching...
No Matches
sgpp::datadriven::DataMatrixDistributed Class Reference

Class to represent a DataMatrix which is distributed on a process grid. More...

#include <DataMatrixDistributed.hpp>

Public Types

enum class  DTYPE : int { DENSE = 1 , TRIDIAG_COEFFICIENT = 501 , TRIDIAG_RHS = 502 , OUT_OF_CORE = 602 }
 Enum for the datatypes of scalapack matrices, see http://netlib.org/scalapack/slug/node73.html. More...
 
enum class  TRIANGULAR { LOWER , UPPER }
 Enum that specifies the lower or upper triangular part of a matrix. More...
 

Public Member Functions

void add (const DataMatrixDistributed &a)
 Adds the values from another DataMatrix to the current values.
 
void copyFrom (const DataMatrixDistributed &other)
 Copies all values of another distributed data matrix to this object, resizes this object to the size of the other matrix.
 
 DataMatrixDistributed ()
 Creates an empty DataMatrixDistributed object.
 
 DataMatrixDistributed (double *input, std::shared_ptr< BlacsProcessGrid > grid, size_t globalRows, size_t globalColumns, size_t rowBlockSize, size_t columnBlockSize, DTYPE dtype=DTYPE::DENSE)
 Creates an two-dimensional DataMatrixDistributed with the specified input.
 
 DataMatrixDistributed (std::shared_ptr< BlacsProcessGrid > grid, size_t globalRows, size_t globalColumns, size_t rowBlockSize, size_t columnBlockSize, double value=0.0, DTYPE dtype=DTYPE::DENSE)
 Creates an two-dimensional DataMatrixDistributed filled with a value.
 
void distribute (const double *matrix, int masterRow=0, int masterCol=0)
 Distribute the matrix from the master on the process grid according to the 2d block cyclic scheme used in scalapack.
 
double get (size_t row, size_t col) const
 Returns the value of the element at position [row,col].
 
size_t getColumnBlockSize () const
 
int * getDescriptor ()
 
const int * getDescriptor () const
 
size_t getGlobalCols () const
 Returns the number of columns of the DataMatrix.
 
size_t getGlobalRows () const
 Returns the number of rows of the DataMatrix.
 
size_t getLocalColumns () const
 
double * getLocalPointer ()
 
const double * getLocalPointer () const
 
size_t getLocalRows () const
 
std::shared_ptr< BlacsProcessGridgetProcessGrid () const
 
size_t getRowBlockSize () const
 
size_t globalToColumnProcessIndex (size_t globalColumnIndex) const
 Calculates the column process index from the global column index.
 
size_t globalToLocalColumnIndex (size_t globalColumnIndex) const
 Calculates the local column index from the globalColumnIndex.
 
size_t globalToLocalRowIndex (size_t globalRowIndex) const
 Calculates the local row index from the globalRowIndex.
 
size_t globalToRowProcessIndex (size_t globalRowIndex) const
 Calculates the row process index from the global row index.
 
bool isProcessMapped () const
 
size_t localToGlobalColumnIndex (size_t localColumnIndex) const
 Calculates the global column index from the local column index.
 
size_t localToGlobalRowIndex (size_t localRowIndex) const
 Calculates the global row index from the local row index.
 
void mult (const DataMatrixDistributed &b, DataMatrixDistributed &c, bool transposeA=false, bool transposeB=false, double alpha=1.0, double beta=0.0) const
 Uses PBLAS to multiply this matrix with matrix B:
 
void mult (const DataVectorDistributed &x, DataVectorDistributed &y, bool transpose=false, double alpha=1.0, double beta=0.0) const
 Multiplies the matrix with a vector x and stores the result in another vector y.
 
double operator() (size_t row, size_t col) const
 Returns the value of the element at position [row,col].
 
void printMatrix () const
 Prints the matrix on stdout on process 0.
 
void resize (size_t rows, size_t cols)
 Resizes the matrix to rows and cols, data is discarded.
 
void set (size_t row, size_t col, double value)
 Sets the element at global position [row,col] to value.
 
void setAll (double value)
 Set all entries of the matrix to one value.
 
void sub (const DataMatrixDistributed &a)
 Subtracts the values from another DataMatrix of the current values.
 
DataMatrix toLocalDataMatrix () const
 
void toLocalDataMatrix (DataMatrix &localMatrix) const
 
DataMatrix toLocalDataMatrixBroadcast () const
 
void toLocalDataMatrixBroadcast (DataMatrix &localMatrix) const
 
DataVectorDistributed toVector () const
 Convert this matrix to a vector.
 
DataMatrixDistributed transpose ()
 Transposes this matrix.
 

Static Public Member Functions

static void add (DataMatrixDistributed &c, const DataMatrixDistributed &a, bool transposeA=false, double beta=1.0, double alpha=1.0)
 Calculates sub(C):=beta*sub(C) + alpha*op(sub(A))
 
static DataMatrixDistributed fromSharedData (const double *input, std::shared_ptr< BlacsProcessGrid > grid, size_t globalRows, size_t globalColumns, size_t rowBlockSize, size_t columnBlockSize, DTYPE dtype=DTYPE::DENSE)
 Creates a distributed data matrix from data which is already shared (mirrored) on each process.
 
static void mult (const DataMatrixDistributed &a, const DataMatrixDistributed &b, DataMatrixDistributed &c, bool transposeA=false, bool transposeB=false, double alpha=1.0, double beta=0.0)
 Uses BLAS to multiply matrix A with matrix B:
 
static void mult (const DataMatrixDistributed &a, const DataVectorDistributed &x, DataVectorDistributed &y, bool transpose=false, double alpha=1.0, double beta=0.0)
 Multiplies matrix A with a vector x and stores the result in another vector y.
 
static void solveCholesky (const DataMatrixDistributed &l, DataVectorDistributed &b, TRIANGULAR uplo=TRIANGULAR::LOWER)
 Solves a linear system of equations Ax=b using a previously computed Cholesky decomposition A=LL^T.
 
static void sub (DataMatrixDistributed &c, const DataMatrixDistributed &a, bool transposeA=false, double beta=1.0, double alpha=1.0)
 Calculates sub(C):=beta*sub(C) - alpha*op(sub(A))
 
static void transpose (const DataMatrixDistributed &a, DataMatrixDistributed &c, double alpha=1.0, double beta=0.0)
 Transposes matrix A and stores the result in C: sub(C):=beta*sub(C) + alpha*sub(A)'.
 

Detailed Description

Class to represent a DataMatrix which is distributed on a process grid.

The class provides a wrapper for ScaLAPACK methods on the matrix. See http://netlib.org/scalapack/slug/node76.html for information regarding the distribution scheme.

Member Enumeration Documentation

◆ DTYPE

Enum for the datatypes of scalapack matrices, see http://netlib.org/scalapack/slug/node73.html.

Enumerator
DENSE 
TRIDIAG_COEFFICIENT 
TRIDIAG_RHS 
OUT_OF_CORE 

◆ TRIANGULAR

Enum that specifies the lower or upper triangular part of a matrix.

Enumerator
LOWER 
UPPER 

Constructor & Destructor Documentation

◆ DataMatrixDistributed() [1/3]

sgpp::datadriven::DataMatrixDistributed::DataMatrixDistributed ( )

Creates an empty DataMatrixDistributed object.

Warning: This object cannot be used.

Referenced by transpose().

◆ DataMatrixDistributed() [2/3]

sgpp::datadriven::DataMatrixDistributed::DataMatrixDistributed ( std::shared_ptr< BlacsProcessGrid grid,
size_t  globalRows,
size_t  globalColumns,
size_t  rowBlockSize,
size_t  columnBlockSize,
double  value = 0.0,
DataMatrixDistributed::DTYPE  dtype = DTYPE::DENSE 
)

Creates an two-dimensional DataMatrixDistributed filled with a value.

Parameters
gridBLACS process grid this matrix will be distributed on
globalRowsnumber of global rows of the matrix
globalColumnsnumber of global columns of the matrix
rowBlockSizeblock size in the row dimension
columnBlockSizeblock size in the column dimension
valueinitial value for all elements of the matrix, default 0
dtypedatatype of the matrix, default DENSE

References sgpp::datadriven::csrc_, sgpp::datadriven::ctxt_, sgpp::datadriven::dtype_, isProcessMapped(), sgpp::datadriven::lld_, sgpp::datadriven::m_, sgpp::datadriven::mb_, sgpp::datadriven::n_, sgpp::datadriven::nb_, and sgpp::datadriven::rsrc_.

◆ DataMatrixDistributed() [3/3]

sgpp::datadriven::DataMatrixDistributed::DataMatrixDistributed ( double *  input,
std::shared_ptr< BlacsProcessGrid grid,
size_t  globalRows,
size_t  globalColumns,
size_t  rowBlockSize,
size_t  columnBlockSize,
DataMatrixDistributed::DTYPE  dtype = DTYPE::DENSE 
)

Creates an two-dimensional DataMatrixDistributed with the specified input.

Call this matrix on all processes in the grid to ensure proper initialization.

Parameters
inputpointer to data that will be distributed
gridBLACS process grid this matrix will be distributed on
globalRowsnumber of global rows of the matrix
globalColumnsnumber of global columns of the matrix
rowBlockSizeblock size in the row dimension
columnBlockSizeblock size in the column dimension
dtypedatatype of the matrix, default DENSE

References distribute(), and isProcessMapped().

Member Function Documentation

◆ add() [1/2]

void sgpp::datadriven::DataMatrixDistributed::add ( const DataMatrixDistributed a)

Adds the values from another DataMatrix to the current values.

Modifies the current values.

Parameters
aThe DataMatrix which is added to the current values

References add().

Referenced by add(), and sub().

◆ add() [2/2]

void sgpp::datadriven::DataMatrixDistributed::add ( DataMatrixDistributed c,
const DataMatrixDistributed a,
bool  transposeA = false,
double  beta = 1.0,
double  alpha = 1.0 
)
static

Calculates sub(C):=beta*sub(C) + alpha*op(sub(A))

Parameters
cmatrix C
amatrix A
transposeAtranspose matrix A if true, default false
betascalar factor for matrix C, default 1.0
alphascalar factor for matrix A, default 1.0

References alpha, getDescriptor(), getGlobalCols(), getGlobalRows(), getLocalPointer(), isProcessMapped(), sgpp::datadriven::pblasNoTranspose, and sgpp::datadriven::pblasTranspose.

◆ copyFrom()

void sgpp::datadriven::DataMatrixDistributed::copyFrom ( const DataMatrixDistributed other)

Copies all values of another distributed data matrix to this object, resizes this object to the size of the other matrix.

References sgpp::datadriven::lld_, sgpp::datadriven::m_, and sgpp::datadriven::n_.

Referenced by sgpp::datadriven::DataVectorDistributed::copyFrom().

◆ distribute()

void sgpp::datadriven::DataMatrixDistributed::distribute ( const double *  matrix,
int  masterRow = 0,
int  masterCol = 0 
)

Distribute the matrix from the master on the process grid according to the 2d block cyclic scheme used in scalapack.

More information: http://www.netlib.org/scalapack/slug/node76.html

Parameters
matrixPointer to the local matrix, only relevant for the master process
masterRowrow coordinate of the master process, default 0
masterColcol coordinate of the master process, default 0

References getLocalPointer().

Referenced by DataMatrixDistributed(), sgpp::datadriven::DataVectorDistributed::distribute(), and sgpp::datadriven::DBMatDMSChol::solveParallel().

◆ fromSharedData()

DataMatrixDistributed sgpp::datadriven::DataMatrixDistributed::fromSharedData ( const double *  input,
std::shared_ptr< BlacsProcessGrid grid,
size_t  globalRows,
size_t  globalColumns,
size_t  rowBlockSize,
size_t  columnBlockSize,
DTYPE  dtype = DTYPE::DENSE 
)
static

◆ get()

double sgpp::datadriven::DataMatrixDistributed::get ( size_t  row,
size_t  col 
) const

Returns the value of the element at position [row,col].

Parameters
rowRow
colColumn
Returns
Value of the element

References getLocalPointer(), and isProcessMapped().

Referenced by sgpp::datadriven::DBMatOfflineChol::compute_inverse_parallel(), sgpp::datadriven::DBMatOfflineOrthoAdapt::decomposeMatrixParallel(), sgpp::datadriven::DataVectorDistributed::get(), and operator()().

◆ getColumnBlockSize()

size_t sgpp::datadriven::DataMatrixDistributed::getColumnBlockSize ( ) const
Returns
the column block size.

Referenced by sgpp::datadriven::DBMatOfflineOrthoAdapt::decomposeMatrixParallel().

◆ getDescriptor() [1/2]

◆ getDescriptor() [2/2]

const int * sgpp::datadriven::DataMatrixDistributed::getDescriptor ( ) const
Returns
Const pointer to ScaLAPACK matrix descriptor

Referenced by python.uq.parameters.ParameterBuilder.GeneralParameterBuilder::new().

◆ getGlobalCols()

◆ getGlobalRows()

◆ getLocalColumns()

size_t sgpp::datadriven::DataMatrixDistributed::getLocalColumns ( ) const
Returns
number of columns assigned to the current process

Referenced by sgpp::datadriven::DataVectorDistributed::getLocalRows().

◆ getLocalPointer() [1/2]

◆ getLocalPointer() [2/2]

const double * sgpp::datadriven::DataMatrixDistributed::getLocalPointer ( ) const
Returns
const pointer to the local data of this process

◆ getLocalRows()

size_t sgpp::datadriven::DataMatrixDistributed::getLocalRows ( ) const
Returns
number of rows assigned to the current process

◆ getProcessGrid()

std::shared_ptr< BlacsProcessGrid > sgpp::datadriven::DataMatrixDistributed::getProcessGrid ( ) const

◆ getRowBlockSize()

size_t sgpp::datadriven::DataMatrixDistributed::getRowBlockSize ( ) const

◆ globalToColumnProcessIndex()

size_t sgpp::datadriven::DataMatrixDistributed::globalToColumnProcessIndex ( size_t  globalColumnIndex) const

Calculates the column process index from the global column index.

Parameters
globalColumnIndex

Referenced by sgpp::datadriven::DataVectorDistributed::globalToRowProcessIndex().

◆ globalToLocalColumnIndex()

size_t sgpp::datadriven::DataMatrixDistributed::globalToLocalColumnIndex ( size_t  globalColumnIndex) const

Calculates the local column index from the globalColumnIndex.

Parameters
globalColumnIndex

Referenced by sgpp::datadriven::DataVectorDistributed::globalToLocalRowIndex().

◆ globalToLocalRowIndex()

size_t sgpp::datadriven::DataMatrixDistributed::globalToLocalRowIndex ( size_t  globalRowIndex) const

Calculates the local row index from the globalRowIndex.

Parameters
globalRowIndex

◆ globalToRowProcessIndex()

size_t sgpp::datadriven::DataMatrixDistributed::globalToRowProcessIndex ( size_t  globalRowIndex) const

Calculates the row process index from the global row index.

Parameters
globalRowIndex

◆ isProcessMapped()

bool sgpp::datadriven::DataMatrixDistributed::isProcessMapped ( ) const

◆ localToGlobalColumnIndex()

size_t sgpp::datadriven::DataMatrixDistributed::localToGlobalColumnIndex ( size_t  localColumnIndex) const

Calculates the global column index from the local column index.

Parameters
localColumnIndex

Referenced by sgpp::datadriven::DataVectorDistributed::localToGlobalRowIndex().

◆ localToGlobalRowIndex()

size_t sgpp::datadriven::DataMatrixDistributed::localToGlobalRowIndex ( size_t  localRowIndex) const

Calculates the global row index from the local row index.

Parameters
localRowIndex

◆ mult() [1/4]

void sgpp::datadriven::DataMatrixDistributed::mult ( const DataMatrixDistributed a,
const DataMatrixDistributed b,
DataMatrixDistributed c,
bool  transposeA = false,
bool  transposeB = false,
double  alpha = 1.0,
double  beta = 0.0 
)
static

Uses BLAS to multiply matrix A with matrix B:

sub(C) := alpha*op(sub(A))*op(sub(B)) + beta*sub(C)

Parameters
[in]amatrix A
[in]bmatrix B
[in,out]cused as output matrix C, can also be added to the product using factor beta
[in]transposeAWhether or not to transpose this matrix
[in]transposeBWhether or not to transpose matrix B
[in]alphafactor alpha, default 1.0
[in]betafactor beta, default 0.0

References alpha, getDescriptor(), getGlobalRows(), getLocalPointer(), isProcessMapped(), sgpp::datadriven::pblasNoTranspose, and sgpp::datadriven::pblasTranspose.

◆ mult() [2/4]

void sgpp::datadriven::DataMatrixDistributed::mult ( const DataMatrixDistributed a,
const DataVectorDistributed x,
DataVectorDistributed y,
bool  transpose = false,
double  alpha = 1.0,
double  beta = 0.0 
)
static

Multiplies matrix A with a vector x and stores the result in another vector y.

sub(y) := alpha*sub(A)'*sub(x) + beta*sub(y)

Parameters
[in]amatrix to be multiplied
[in]xvector to be multiplied
[in,out]yvector in which the result should be stored
[in]transposetranspose if true, default false
[in]alphafactor alpha, default 1.0
[in]betafactor beta, default 0.0

References getDescriptor(), sgpp::datadriven::DataVectorDistributed::getDescriptor(), getGlobalCols(), getGlobalRows(), getLocalPointer(), sgpp::datadriven::DataVectorDistributed::getLocalPointer(), isProcessMapped(), sgpp::datadriven::DataVectorDistributed::isProcessMapped(), sgpp::datadriven::pblasNoTranspose, sgpp::datadriven::pblasTranspose, and transpose().

◆ mult() [3/4]

void sgpp::datadriven::DataMatrixDistributed::mult ( const DataMatrixDistributed b,
DataMatrixDistributed c,
bool  transposeA = false,
bool  transposeB = false,
double  alpha = 1.0,
double  beta = 0.0 
) const

Uses PBLAS to multiply this matrix with matrix B:

sub(C) := alpha*op(sub(this))*op(sub(B)) + beta*sub(C)

Parameters
[in]bmatrix B
[in,out]cused as output matrix C, can also be added to the product using factor beta
[in]transposeAWhether or not to transpose this matrix
[in]transposeBWhether or not to transpose matrix B
[in]alphafactor alpha, default 1.0
[in]betafactor beta, default 0.0

References alpha, and mult().

◆ mult() [4/4]

void sgpp::datadriven::DataMatrixDistributed::mult ( const DataVectorDistributed x,
DataVectorDistributed y,
bool  transpose = false,
double  alpha = 1.0,
double  beta = 0.0 
) const

Multiplies the matrix with a vector x and stores the result in another vector y.

sub(y) := alpha*sub(this)'*sub(x) + beta*sub(y)

Parameters
[in]xvector to be multiplied
[in,out]yvector in which the result should be stored
[in]transposetranspose if true, default false
[in]alphafactor alpha, default 1.0
[in]betafactor beta, default 0.0

References alpha, mult(), and transpose().

Referenced by sgpp::datadriven::DBMatOfflineOrthoAdapt::compute_inverse_parallel(), sgpp::datadriven::ModelFittingDensityEstimationOnOffParallel::computeResidual(), mult(), mult(), sgpp::datadriven::DBMatOnlineDE_SMW::smw_adapt_parallel(), sgpp::datadriven::DBMatDMS_SMW::solveParallel(), and sgpp::datadriven::DBMatDMSOrthoAdapt::solveParallel().

◆ operator()()

double sgpp::datadriven::DataMatrixDistributed::operator() ( size_t  row,
size_t  col 
) const

Returns the value of the element at position [row,col].

Only a const version exists as setting an element is not necessarily done in local memory.

Parameters
rowRow
colColumn
Returns
constant reference to the element

References get().

◆ printMatrix()

void sgpp::datadriven::DataMatrixDistributed::printMatrix ( ) const

Prints the matrix on stdout on process 0.

References toLocalDataMatrix(), and sgpp::base::DataMatrix::toString().

◆ resize()

void sgpp::datadriven::DataMatrixDistributed::resize ( size_t  rows,
size_t  cols 
)

Resizes the matrix to rows and cols, data is discarded.

Parameters
rows
cols

References getGlobalCols(), getGlobalRows(), isProcessMapped(), sgpp::datadriven::lld_, sgpp::datadriven::m_, and sgpp::datadriven::n_.

Referenced by sgpp::datadriven::DataVectorDistributed::resize().

◆ set()

void sgpp::datadriven::DataMatrixDistributed::set ( size_t  row,
size_t  col,
double  value 
)

◆ setAll()

void sgpp::datadriven::DataMatrixDistributed::setAll ( double  value)

Set all entries of the matrix to one value.

Parameters
value

References isProcessMapped().

Referenced by sgpp::datadriven::DataVectorDistributed::setAll().

◆ solveCholesky()

void sgpp::datadriven::DataMatrixDistributed::solveCholesky ( const DataMatrixDistributed l,
DataVectorDistributed b,
DataMatrixDistributed::TRIANGULAR  uplo = TRIANGULAR::LOWER 
)
static

Solves a linear system of equations Ax=b using a previously computed Cholesky decomposition A=LL^T.

Parameters
[in]llower triangular matrix L of the Cholesky decomposition
[in,out]bvector b of the linear system, is overwritten with solution x
[in]uplostorage format of the triangular matrix, either upper or lower triangular

References getDescriptor(), getGlobalRows(), getLocalPointer(), isProcessMapped(), LOWER, sgpp::datadriven::lowerTriangular, and sgpp::datadriven::upperTriangular.

Referenced by sgpp::datadriven::DBMatDMSChol::solveParallel().

◆ sub() [1/2]

void sgpp::datadriven::DataMatrixDistributed::sub ( const DataMatrixDistributed a)

Subtracts the values from another DataMatrix of the current values.

Modifies the current values.

Parameters
aThe DataMatrix which is subtracted from the current values

References sub().

Referenced by sgpp::datadriven::DBMatOnlineDE_SMW::smw_adapt_parallel(), and sub().

◆ sub() [2/2]

void sgpp::datadriven::DataMatrixDistributed::sub ( DataMatrixDistributed c,
const DataMatrixDistributed a,
bool  transposeA = false,
double  beta = 1.0,
double  alpha = 1.0 
)
static

Calculates sub(C):=beta*sub(C) - alpha*op(sub(A))

Parameters
cmatrix C
amatrix A
transposeAtranspose matrix A if true, default false
betascalar factor for matrix C, default 1.0
alphascalar factor for matrix A, default 1.0

References add(), and alpha.

◆ toLocalDataMatrix() [1/2]

DataMatrix sgpp::datadriven::DataMatrixDistributed::toLocalDataMatrix ( ) const

◆ toLocalDataMatrix() [2/2]

void sgpp::datadriven::DataMatrixDistributed::toLocalDataMatrix ( DataMatrix localMatrix) const
Parameters
[out]localMatrixthe gathered DataMatrix as a normal, not distributed, DataMatrix. Result can only be used on the master process.

References isProcessMapped().

◆ toLocalDataMatrixBroadcast() [1/2]

DataMatrix sgpp::datadriven::DataMatrixDistributed::toLocalDataMatrixBroadcast ( ) const

◆ toLocalDataMatrixBroadcast() [2/2]

void sgpp::datadriven::DataMatrixDistributed::toLocalDataMatrixBroadcast ( DataMatrix localMatrix) const
Parameters
[out]localMatrixthe whole DataMatrix is broadcasted to all processes in the grid

References isProcessMapped().

◆ toVector()

DataVectorDistributed sgpp::datadriven::DataMatrixDistributed::toVector ( ) const

Convert this matrix to a vector.

This operation is only possible if either globalRows or globalColumns equals one.

Returns
A DataVectorDistributed object with the same size and data

References toLocalDataMatrix().

◆ transpose() [1/2]

◆ transpose() [2/2]

void sgpp::datadriven::DataMatrixDistributed::transpose ( const DataMatrixDistributed a,
DataMatrixDistributed c,
double  alpha = 1.0,
double  beta = 0.0 
)
static

Transposes matrix A and stores the result in C: sub(C):=beta*sub(C) + alpha*sub(A)'.

Parameters
[in]amatrix A to transpose
[out]cresult matrix C
[in]alphafactor for matrix A, default 1.0
[in]betafactor for matrix C, default 0.0

References alpha, getDescriptor(), getGlobalCols(), getGlobalRows(), getLocalPointer(), and isProcessMapped().


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