Elemental 0.77 documentation

The DistMatrix class

«  The Grid class   ::   Contents   ::   Partitioning  »

The DistMatrix class

The DistMatrix<T,U,V> class is meant to provide a distributed-memory analogue of the Matrix<T> class. Similarly to PLAPACK, roughly ten different matrix distributions are provided and it is trivial (in the programmability sense) to redistribute from one to another: in PLAPACK, one would simply call PLA_Copy, whereas, in Elemental, it is handled through overloading the \(=\) operator.

Since it is crucial to know not only how many processes to distribute the data over, but which processes, and in what manner they should be decomposed into a logical two-dimensional grid, an instance of the Grid class must be passed into the constructor of the DistMatrix<T,U,V> class.

Note

Since the DistMatrix<T,U,V> class makes use of MPI for message passing, custom interfaces must be written for nonstandard datatypes. As of now, the following datatypes are fully supported for DistMatrix<T,U,V>: int, float, double, Complex<float>, and Complex<double>.

AbstractDistMatrix

This abstract class defines the list of member functions that are guaranteed to be available for all matrix distributions.

type class AbstractDistMatrix<T>

The most general case, where the underlying datatype T is only assumed to be a ring; that is, it supports multiplication and addition and has the appropriate identities.

Basic information

int Height() const

Return the height of the matrix.

int Width() const

Return the width of the matrix.

int LocalHeight() const

Return the local height of the matrix.

int LocalWidth() const

Return the local width of the matrix.

int LocalLDim() const

Return the local leading dimension of the matrix.

size_t AllocatedMemory() const

Return the number of entries of type T that we have locally allocated space for.

const elem::Grid& Grid() const

Return the grid that this distributed matrix is distributed over.

T* LocalBuffer(int iLocal=0, int jLocal=0 )

Return a pointer to the portion of the local buffer that stores entry (iLocal,jLocal).

const T* LockedLocalBuffer(int iLocal=0, int jLocal=0 ) const

Return a pointer to the portion of the local buffer that stores entry (iLocal,jLocal), but do not allow for the data to be modified through the returned pointer.

Matrix<T>& LocalMatrix()

Return a reference to the local matrix.

const Matrix<T>& LockedLocalMatrix() const

Return an unmodifiable reference to the local matrix.

I/O

void Print(const std::string msg="") const

Print the distributed matrix to standard output (std::cout).

void Print(std::ostream& os, const std::string msg="") const

Print the distributed matrix to the output stream os.

void Write(const std::string filename, const std::string msg="") const

Print the distributed matrix to the file named filename.

Distribution details

void FreeAlignments()

Free all alignment constaints.

bool ConstrainedColAlignment() const

Return whether or not the column alignment is constrained.

bool ConstrainedRowAlignment() const

Return whether or not the row alignment is constrained.

int ColAlignment() const

Return the alignment of the columns of the matrix.

int RowAlignment() const

Return the alignment of the rows of the matrix.

int ColShift() const

Return the first global row that our process owns.

int RowShift() const

Return the first global column that our process owns.

int ColStride() const

Return the number of rows between locally owned entries.

int RowStride() const

Return the number of columns between locally owned entries.

Entry manipulation

T Get(int i, int j) const

Return the (i,j) entry of the global matrix. This operation is collective.

void Set(int i, int j, T alpha)

Set the (i,j) entry of the global matrix to \(\alpha\). This operation is collective.

void Update(int i, int j, T alpha)

Add \(\alpha\) to the (i,j) entry of the global matrix. This operation is collective.

T GetLocal(int iLocal, int jLocal) const

Return the (iLocal,jLocal) entry of our local matrix.

void SetLocal(int iLocal, int jLocal, T alpha)

Set the (iLocal,jLocal) entry of our local matrix to \(\alpha\).

void UpdateLocal(int iLoca, int jLocal, T alpha)

Add \(\alpha\) to the (iLocal,jLocal) entry of our local matrix.

Note

Many of the following routines are only valid for complex datatypes.

typename Base<T>::type GetRealPart(int i, int j) const

Return the real part of the (i,j) entry of the global matrix. This operation is collective.

typename Base<T>::type GetImagPart(int i, int j) const

Return the imaginary part of the (i,j) entry of the global matrix. This operation is collective.

void SetRealPart(int i, int j, typename Base<T>::type alpha)

Set the real part of the (i,j) entry of the global matrix to \(\alpha\).

void SetImagPart(int i, int j, typename Base<T>::type alpha)

Set the imaginary part of the (i,j) entry of the global matrix to \(\alpha\).

void UpdateRealPart(int i, int j, typename Base<T>::type alpha)

Add \(\alpha\) to the real part of the (i,j) entry of the global matrix.

void UpdateImagPart(int i, int j, typename Base<T>::type alpha)

Add \(\alpha\) to the imaginary part of the (i,j) entry of the global matrix.

typename Base<T>::type GetRealPartLocal(int iLocal, int jLocal) const

Return the real part of the (iLocal,jLocal) entry of our local matrix.

typename Base<T>::type GetLocalImagPart(int iLocal, int jLocal) const

Return the imaginary part of the (iLocal,jLocal) entry of our local matrix.

void SetLocalRealPart(int iLocal, int jLocal, typename Base<T>::type alpha)

Set the real part of the (iLocal,jLocal) entry of our local matrix.

void SetLocalImagPart(int iLocal, int jLocal, typename Base<T>::type alpha)

Set the imaginary part of the (iLocal,jLocal) entry of our local matrix.

void UpdateRealPartLocal(int iLocal, int jLocal, typename Base<T>::type alpha)

Add \(\alpha\) to the real part of the (iLocal,jLocal) entry of our local matrix.

void UpdateLocalImagPart(int iLocal, int jLocal, typename Base<T>::type alpha)

Add \(\alpha\) to the imaginary part of the (iLocal,jLocal) entry of our local matrix.

Viewing

bool Viewing() const

Return whether or not this matrix is viewing another.

bool LockedView() const

Return whether or not this matrix is viewing another in a manner that does not allow for modifying the viewed data.

Utilities

void Empty()

Resize the distributed matrix so that it is \(0 \times 0\) and free all allocated storage.

void ResizeTo(int height, int width)

Reconfigure the matrix so that it is height \(\times\) width.

void SetGrid(const elem::Grid& grid)

Clear the distributed matrix’s contents and reconfigure for the new process grid.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of AbstractDistMatrix<T>.

type class AbstractDistMatrix<R>

Used to denote that the underlying datatype R is real.

type class AbstractDistMatrix<Complex<R>>

Used to denote that the underlying datatype Complex<R> is complex with base type R.

type class AbstractDistMatrix<F>

Used to denote that the underlying datatype F is a field.

DistMatrix

type class DistMatrix<T, U, V>

This templated class for manipulating distributed matrices is only defined for the following choices of the column and row Distribution‘s, U and V (T is a ring in this case).

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,U,V>.

type class DistMatrix<double, U, V>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, U, V>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, U, V>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, U, V>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, U, V>

The underlying datatype F is a field.

[MC,MR]

This is by far the most important matrix distribution in Elemental, as the vast majority of parallel routines expect the input to be in this form. For a \(7 \times 7\) matrix distributed over a \(2 \times 3\) process grid, individual entries would be owned by the following processes (assuming the column and row alignments are both 0):

\[\begin{split}\left(\begin{array}{ccccccc} 0 & 2 & 4 & 0 & 2 & 4 & 0 \\ 1 & 3 & 5 & 1 & 3 & 5 & 1 \\ 0 & 2 & 4 & 0 & 2 & 4 & 0 \\ 1 & 3 & 5 & 1 & 3 & 5 & 1 \\ 0 & 2 & 4 & 0 & 2 & 4 & 0 \\ 1 & 3 & 5 & 1 & 3 & 5 & 1 \\ 0 & 2 & 4 & 0 & 2 & 4 & 0 \end{array}\right)\end{split}\]

Similarly, if the column alignment is kept at 0 and the row alignment is changed to 2 (meaning that the third process column owns the first column of the matrix), the individual entries would be owned as follows:

\[\begin{split}\left(\begin{array}{ccccccc} 4 & 0 & 2 & 4 & 0 & 2 & 4 \\ 5 & 1 & 3 & 5 & 1 & 3 & 5 \\ 4 & 0 & 2 & 4 & 0 & 2 & 4 \\ 5 & 1 & 3 & 5 & 1 & 3 & 5 \\ 4 & 0 & 2 & 4 & 0 & 2 & 4 \\ 5 & 1 & 3 & 5 & 1 & 3 & 5 \\ 4 & 0 & 2 & 4 & 0 & 2 & 4 \end{array}\right)\end{split}\]

It should also be noted that this is the default distribution format for the DistMatrix<T,U,V> class, as DistMatrix<T> defaults to DistMatrix<T,MC,MR>.

type class DistMatrix<T>
type class DistMatrix<T, MC, MR>

The most general case, where the underlying datatype T is only assumed to be a ring.

Constructors

DistMatrix( const elem::Grid& grid=DefaultGrid() )

Create a \(0 \times 0\) distributed matrix over the specified grid.

DistMatrix( int height, int width, const elem::Grid& grid=DefaultGrid() )

Create a height \(\times\) width distributed matrix over the specified grid.

DistMatrix(int height, int width, bool constrainedColAlignment, bool constrainedRowAlignment, int colAlignment, int rowAlignment, const elem::Grid& grid)

Create a height \(\times\) width distributed matrix distributed over the specified process grid, but with the top-left entry owned by the colAlignment process row and the rowAlignment process column. Each of these alignments may be constrained to remain constant when redistributing data into this DistMatrix<T>.

DistMatrix(int height, int width, bool constrainedColAlignment, bool constrainedRowAlignment, int colAlignment, int rowAlignment, int ldim, const elem::Grid& grid)

Same as above, but the local leading dimension is also specified.

DistMatrix(int height, int width, int colAlignment, int rowAlignment, const T* buffer, int ldim, const elem::Grid& grid)

View a constant distributed matrix’s buffer; the buffer must correspond to the local portion of an elemental distributed matrix with the specified row and column alignments and leading dimension, ldim.

DistMatrix(int height, int width, int colAlignment, int rowAlignment, T* buffer, int ldim, const elem::Grid& grid)

Same as above, but the contents of the matrix are modifiable.

DistMatrix(const DistMatrix<T, U, V>& A)

Build a copy of the distributed matrix A, but force it to be in the [MC,MR] distribution.

Redistribution

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, MC, MR>& A)

If this matrix can be properly aligned with A, then perform a local copy, otherwise perform an mpi::SendRecv() permutation first.

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, MC, STAR>& A)

Perform a local (filtered) copy to form an [MC,MR ] distribution and then, if necessary, fix the alignment of the MC distribution via an mpi::SendRecv() within process columns.

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, STAR, MR>& A)

Perform a local (filtered) copy to form an [MC,MR ] distribution and then, if necessary, fix the alignment of the MR distribution via an mpi::SendRecv() within process rows.

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, MD, STAR>& A)

Since the [MD,STAR] distribution is defined such that its columns are distributed like a diagonal of an [MC,MR] distributed matrix, this operation is not very common.

Note

This redistribution routine is not yet implemented.

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, STAR, MD>& A)

Note

This redistribution routine is not yet implemented.

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, MR, MC>& A)

This routine serves to transpose the distribution of A[MR,MC] into the standard matrix distribution, A[MC,MR]. This redistribution is implemented with four different approaches: one for matrices that are taller than they are wide, one for matrices that are wider than they are tall, one for column vectors, and one for row vectors.

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, MR, STAR>& A)

This is similar to the above routine, but with each row of A being undistributed, and only one approach is needed: \(A[M_C,M_R] \leftarrow A[V_C,\star] \leftarrow A[V_R,\star] \leftarrow A[M_R,\star]\).

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, STAR, MC>& A)

This routine is dual to the \(A[M_C,M_R] \leftarrow A[M_R,\star]\) redistribution and is accomplished through the sequence: \(A[M_C,M_R] \leftarrow A[\star,V_R] \leftarrow A[\star,V_C] \leftarrow A[\star,M_C]\).

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, VC, STAR>& A)

Perform an mpi::AllToAll() within process rows in order to redistribute to the [MC,MR] distribution (an mpi::SendRecv() within process columns may be required for alignment).

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, STAR, VC>& A)

Accomplished through the sequence \(A[M_C,M_R] \leftarrow A[\star,V_R] \leftarrow A[\star,V_C]\).

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, VR, STAR>& A)

Accomplished through the sequence \(A[M_C,M_R] \leftarrow A[V_C,\star] \leftarrow A[V_R,\star]\).

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, STAR, VR>& A)

Perform an mpi::AllToAll() within process columns in order to redistribute to the [MC,MR] distribution (an mpi::SendRecv() within process rows may be required for alignment).

const DistMatrix<T, MC, MR>& operator=(const DistMatrix<T, STAR, STAR>& A)

Perform an mpi::AllGather() over the entire grid in order to give every process a full copy of A.

Diagonal manipulation

void GetDiagonal(DistMatrix<T, MD, STAR>& d, int offset=0 ) const

The \([M_D,\star]\) distribution is defined such that its columns are distributed like diagonals of the standard matrix distribution, [MC,MR]. Thus, d can be formed locally if the distribution can be aligned with that of the offset diagonal of \(A[M_C,M_R]\).

void GetDiagonal(DistMatrix<T, STAR, MD>& d, int offset=0 ) const

This is the same as above, but d is a row-vector instead of a column-vector.

void SetDiagonal(const DistMatrix<T, MD, STAR>& d, int offset=0 )

Same as DistMatrix<T>::GetDiagonal(), but in reverse.

void SetDiagonal(const DistMatrix<T, STAR, MD>& d, int offset=0 )

Same as DistMatrix<T>::GetDiagonal(), but in reverse.

Note

Many of the following routines are only valid for complex datatypes and are analogous to their general counterparts from above in the obvious manner.

void GetRealPartOfDiagonal(DistMatrix<typename Base<T>::type, MD, STAR>& d, int offset=0 ) const
void GetImagPartOfDiagonal(DistMatrix<typename Base<T>::type, MD, STAR>& d, int offset=0 ) const
void GetRealPartOfDiagonal(DistMatrix<typename Base<T>::type, STAR, MD>& d, int offset=0 ) const
void GetImagPartOfDiagonal(DistMatrix<typename Base<T>::type, STAR, MD>& d, int offset=0 ) const
void SetRealPartOfDiagonal(const DistMatrix<typename Base<T>::type, MD, STAR>& d, int offset=0 )
void SetImagPartOfDiagonal(const DistMatrix<typename Base<T>::type, MD, STAR>& d, int offset=0 )
void SetRealPartOfDiagonal(const DistMatrix<typename Base<T>::type, STAR, MD>& d, int offset=0 )
void SetImagPartOfDiagonal(const DistMatrix<typename Base<T>::type, STAR, MD>& d, int offset=0 )

Alignment

All of the following clear the distributed matrix’s contents and then reconfigure the alignments as described.

void Align(int colAlignment, int rowAlignment)

Specify the process row, colAlignment, and process column, rowAlignment, which own the top-left entry.

void AlignCols(int colAlignment)

Specify the process row which owns the top-left entry.

void AlignRows(int rowAlignment)

Specify the process column which owns the top-left entry.

void AlignWith(const DistMatrix<S, MC, MR>& A)

Force the alignments to match those of A.

void AlignWith(const DistMatrix<S, MC, STAR>& A)

Force the column alignment to match that of A.

void AlignWith(const DistMatrix<S, STAR, MR>& A)

Force the row alignment to match that of A.

void AlignWith(const DistMatrix<S, MR, MC>& A)

Force the column alignment to match the row alignment of A (and vice-versa).

void AlignWith(const DistMatrix<S, MR, STAR>& A)

Force the row alignment to match the column alignment of A.

void AlignWith(const DistMatrix<S, STAR, MC>& A)

Force the column alignment to match the row alignment of A.

void AlignWith(const DistMatrix<S, VC, STAR>& A)

Force the column alignment to be equal to that of A (modulo the number of process rows).

void AlignWith(const DistMatrix<S, STAR, VC>& A)

Force the column alignment to equal the row alignment of A (modulo the number of process rows).

void AlignWith(const DistMatrix<S, VR, STAR>& A)

Force the row alignment to equal the column alignment of A (modulo the number of process columns).

void AlignWith(const DistMatrix<S, STAR, VR>& A)

Force the row alignment to equal the row alignment of A (modulo the number of process columns).

void AlignColsWith(const DistMatrix<S, MC, MR>& A)

Force the column alignment to match that of A.

void AlignColsWith(const DistMatrix<S, MC, STAR>& A)

Force the column alignment to match that of A.

void AlignColsWith(const DistMatrix<S, MR, MC>& A)

Force the column alignment to match the row alignment of A.

void AlignColsWith(const DistMatrix<S, STAR, MC>& A)

Force the column alignment to match the row alignment of A.

void AlignColsWith(const DistMatrix<S, VC, STAR>& A)

Force the column alignment to match the column alignment of A (modulo the number of process rows).

void AlignColsWith(const DistMatrix<S, STAR, VC>& A)

Force the column alignment to match the row alignment of A (modulo the number of process rows).

void AlignRowsWith(const DistMatrix<S, MC, MR>& A)

Force the row alignment to match that of A.

void AlignRowsWith(const DistMatrix<S, STAR, MR>& A)

Force the row alignment to match that of A.

void AlignRowsWith(const DistMatrix<S, MR, MC>& A)

Force the row alignment to match the column alignment of A.

void AlignRowsWith(const DistMatrix<S, MR, STAR>& A)

Force the row alignment to match the column alignment of A.

void AlignRowsWith(const DistMatrix<S, VR, STAR>& A)

Force the row alignment to match the column alignment of A (modulo the number of process columns).

void AlignRowsWith(const DistMatrix<S, STAR, VR>& A)

Force the row alignment to match the row alignment of A (modulo the number of process columns).

Views

void View(DistMatrix<T, MC, MR>& A)

Reconfigure this matrix such that it is essentially a copy of the distributed matrix A, but the local data buffer simply points to the one from A.

void LockedView(const DistMatrix<T, MC, MR>& A)

Same as above, but this matrix is “locked”, meaning that it cannot change the data from A that it points to.

void View(DistMatrix<T, MC, MR>& A, int i, int j, int height, int width)

View a subset of A rather than the entire matrix. In particular, reconfigure this matrix to behave like the submatrix defined from the [i,i+height) rows and [j,j+width) columns of A.

void LockedView(const DistMatrix<T, MC, MR>& A, int i, int j, int height, int width)

Same as above, but this matrix is “locked”, meaning that it cannot change the data from A that it points to.

void View(int height, int width, int colAlignment, int rowAlignment, T* buffer, int ldim, const elem::Grid& grid)

Reconfigure this distributed matrix around an implicit [MC,MR] distributed matrix of the specified dimensions, alignments, local buffer, local leading dimension, and process grid.

void LockedView(int height, int width, int colAlignment, int rowAlignment, const T* buffer, int ldim, const elem::Grid& grid)

Same as above, but the resulting matrix is “locked”, meaning that it cannot modify the underlying local data.

Note

The following functions have strict requirements on the input matrices and must be used with care in PureRelease and HybridRelease modes.

void View1x2(DistMatrix<T, MC, MR>& AL, DistMatrix<T, MC, MR>& AR)

Recombine two adjacent submatrices to form \([A_L A_R]\).

void LockedView1x2(const DistMatrix<T, MC, MR>& AL, const DistMatrix<T, MC, MR>& AR)

Same as above, but the result is “locked” (the data is not modifiable).

void View2x1(DistMatrix<T, MC, MR>& AT, DistMatrix<T, MC, MR>& AB)

Recombine two adjacent submatrices to form \([A_T; A_B]\).

void LockedView2x1(const DistMatrix<T, MC, MR>& AT, const DistMatrix<T, MC, MR>& AB)

Same as above, but the result is “locked” (the data is not modifiable).

void View2x2(DistMatrix<T, MC, MR>& ATL, DistMatrix<T, MC, MR>& ATR, DistMatrix<T, MC, MR>& ABL, DistMatrix<T, MC, MR>& ABR)

Recombine four adjacent submatrices to form \([A_{TL} A_{TR}; A_{BL} A_{BR}]\).

void LockedView2x2(const DistMatrix<T, MC, MR>& ATL, const DistMatrix<T, MC, MR>& ATR, const DistMatrix<T, MC, MR>& ABL, const DistMatrix<T, MC, MR>& ABR)

Same as above, but the result is “locked” (the data is not modifiable).

Custom communication routines

The following routines primarily exist as a means of avoiding the poor memory bandwidth which results from packing or unpacking large amounts of data without a unit stride. PLAPACK noticed this issue and avoided the problem by carefully (conjugate-)transposing matrices in strategic places, usually before a gather or after a scatter, and we follow suit.

void SumScatterFrom(const DistMatrix<T, MC, STAR>& A)

Simultaneously sum \(A[M_C,\star]\) within each process row and scatter the entries in each row to form the result in an \([M_C,M_R]\) distribution.

void SumScatterUpdate(T alpha, const DistMatrix<T, MC, STAR>& A)

Same as above, but add \(\alpha\) times the result onto the parent distributed matrix rather than simply assigning the result to it.

void SumScatterFrom(const DistMatrix<T, STAR, MR>& A)

Simultaenously sum \(A[\star,M_R]\) within each process column and scatter the entries in each column to form the result in an \([M_C,M_R]\) distribution.

void SumScatterUpdate(T alpha, const DistMatrix<T, STAR, MR>& A)

Same as above, but add \(\alpha\) times the result onto the parent distributed matrix rather than simply assigning the result to it.

void SumScatterFrom(const DistMatrix<T, STAR, STAR>& A)

Simultaneously sum \(A[\star,\star]\) over the entire process grid and scatter the entries in each row and column to form the result in an \([M_C,M_R]\) distribution.

void SumScatterUpdate(T alpha, const DistMatrix<T, STAR, STAR>& A)

Same as above, but add \(\alpha\) times the result onto the parent distributed matrix rather than simply assigning the result to it.

void AdjointFrom(const DistMatrix<T, STAR, MC>& A)

Set the parent matrix equal to the redistributed adjoint of \(A[\star,M_C]\); in particular, \((A[\star,M_C])^H = A^H[M_C,\star]\), so perform an \([M_C,M_R] \leftarrow [M_C,\star]\) redistribution on the adjoint of A, which typically just consists of locally copying (and conjugating) subsets of the data from \(A[\star,M_C]\).

void AdjointFrom(const DistMatrix<T, MR, STAR>& A)

This routine is the dual of the above routine, and performs an \([M_C,M_R] \leftarrow [\star,M_R]\) redistribution on the adjoint of A.

void TransposeFrom(const DistMatrix<T, STAR, MC>& A)

Same as the corresponding DistMatrix<T>::AdjointFrom(), but with no conjugation.

void TransposeFrom(const DistMatrix<T, MR, STAR>& A)

Same as the corresponding DistMatrix<T>::AdjointFrom(), but with no conjugation.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,MC,MR> = DistMatrix<T>.

type class DistMatrix<double>
type class DistMatrix<double, MC, MR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>>
type class DistMatrix<Complex<double>, MC, MR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R>
type class DistMatrix<R, MC, MR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>>
type class DistMatrix<Complex<R>, MC, MR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F>
type class DistMatrix<F, MC, MR>

The underlying datatype F is a field.

[MC,* ]

This distribution is often used as part of matrix-matrix multiplication. For a \(7 \times 7\) matrix distributed over a \(2 \times 3\) process grid, individual entries would be owned by the following processes (assuming the column alignment is 0):

\[\begin{split}\left(\begin{array}{ccccccc} \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} \\ \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} \\ \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} \\ \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} \\ \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} \\ \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} & \{1,3,5\} \\ \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} & \{0,2,4\} \end{array}\right)\end{split}\]
type class DistMatrix<T, MC, STAR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,MC,STAR>.

type class DistMatrix<double, MC, STAR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, MC, STAR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, MC, STAR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, MC, STAR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, MC, STAR>

The underlying datatype F is a field.

[* ,MR]

This distribution is also frequently used for matrix-matrix multiplication. For a \(7 \times 7\) matrix distributed over a \(2 \times 3\) process grid, individual entries would be owned by the following processes (assuming the row alignment is 0):

\[\begin{split}\left(\begin{array}{ccccccc} \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \\ \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \\ \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \\ \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \\ \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \\ \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \\ \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} & \{2,3\} & \{4,5\} & \{0,1\} \end{array}\right)\end{split}\]
type class DistMatrix<T, STAR, MR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,STAR,MR>.

type class DistMatrix<double, STAR, MR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, STAR, MR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, STAR, MR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, STAR, MR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, STAR, MR>

The underlying datatype F is a field.

[MR,MC]

This is essentially the transpose of the standard matrix distribution, [MC,MR]. For a \(7 \times 7\) matrix distributed over a \(2 \times 3\) process grid, individual entries would be owned by the following processes (assuming the column and row alignments are both 0):

\[\begin{split}\left(\begin{array}{ccccccc} 0 & 1 & 0 & 1 & 0 & 1 & 0 \\ 2 & 3 & 2 & 3 & 2 & 3 & 2 \\ 4 & 5 & 4 & 5 & 4 & 5 & 4 \\ 0 & 1 & 0 & 1 & 0 & 1 & 0 \\ 2 & 3 & 2 & 3 & 2 & 3 & 2 \\ 4 & 5 & 4 & 5 & 4 & 5 & 4 \\ 0 & 1 & 0 & 1 & 0 & 1 & 0 \end{array}\right)\end{split}\]
type class DistMatrix<T, MR, MC>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,MR,MC>.

type class DistMatrix<double, MR, MC>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, MR, MC>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, MR, MC>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, MR, MC>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, MR, MC>

The underlying datatype F is a field.

[MR,* ]

This is the transpose of the [* ,MR] distribution and is, like many of the previous distributions, useful for matrix-matrix multiplication. For a \(7 \times 7\) matrix distributed over a \(2 \times 3\) process grid, individual entries would be owned by the following processes (assuming the column alignment is 0):

\[\begin{split}\left(\begin{array}{ccccccc} \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} \\ \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} \\ \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} \\ \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} \\ \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} & \{2,3\} \\ \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} & \{4,5\} \\ \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} & \{0,1\} \end{array}\right)\end{split}\]
type class DistMatrix<T, MR, STAR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,MR,STAR>.

type class DistMatrix<double, MR, STAR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, MR, STAR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, MR, STAR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, MR, STAR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, MR, STAR>

The underlying datatype F is a field.

[* ,MC]

This is the transpose of the [MC,*] distribution and is, like many of the previous distributions, useful for matrix-matrix multiplication. For a \(7 \times 7\) matrix distributed over a \(2 \times 3\) process grid, individual entries would be owned by the following processes (assuming the column alignment is 0):

\[\begin{split}\left(\begin{array}{ccccccc} \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \\ \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \\ \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \\ \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \\ \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \\ \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \\ \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} & \{1,3,5\} & \{0,2,4\} \end{array}\right)\end{split}\]
type class DistMatrix<T, STAR, MC>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,STAR,MC>.

type class DistMatrix<double, STAR, MC>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, STAR, MC>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, STAR, MC>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, STAR, MC>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, STAR, MC>

The underlying datatype F is a field.

[MD,* ]

TODO, but not as high of a priority since the \([M_D,\star]\) distribution is not as crucial for end users as many other details that have not yet been documented.

type class DistMatrix<T, MD, STAR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,MD,STAR>.

type class DistMatrix<double, MD, STAR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, MD, STAR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, MD, STAR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, MD, STAR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, MD, STAR>

The underlying datatype F is a field.

[* ,MD]

TODO, but not as high of a priority since the \([\star,M_D]\) distribution is not as crucial for end users as many other details that have not yet been documented.

type class DistMatrix<T, STAR, MD>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,STAR,MD>.

type class DistMatrix<double, STAR, MD>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, STAR, MD>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, STAR, MD>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, STAR, MD>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, STAR, MD>

The underlying datatype F is a field.

[VC,* ]

This distribution makes use of a 1d distribution which uses a column-major ordering of the entire process grid. Since 1d distributions are useful for distributing vectors, and a column-major ordering is used, the distribution symbol is VC. Again using the simple \(2 \times 3\) process grid, with a zero column alignment, each entry of a \(7 \times 7\) matrix would be owned by the following sets of processes:

\[\begin{split}\left(\begin{array}{ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 1 & 1 & 1 & 1 & 1 & 1 \\ 2 & 2 & 2 & 2 & 2 & 2 & 2 \\ 3 & 3 & 3 & 3 & 3 & 3 & 3 \\ 4 & 4 & 4 & 4 & 4 & 4 & 4 \\ 5 & 5 & 5 & 5 & 5 & 5 & 5 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{array}\right)\end{split}\]
type class DistMatrix<T, VC, STAR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,VC,STAR>.

type class DistMatrix<double, VC, STAR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, VC, STAR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, VC, STAR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, VC, STAR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, VC, STAR>

The underlying datatype F is a field.

[* ,VC]

This is the transpose of the above [VC,* ] distribution. On the standard \(2 \times 3\) process grid with a row alignment of zero, a \(7 \times 7\) matrix would be distributed as:

\[\begin{split}\left(\begin{array}{ccccccc} 0 & 1 & 2 & 3 & 4 & 5 & 0 \\ 0 & 1 & 2 & 3 & 4 & 5 & 0 \\ 0 & 1 & 2 & 3 & 4 & 5 & 0 \\ 0 & 1 & 2 & 3 & 4 & 5 & 0 \\ 0 & 1 & 2 & 3 & 4 & 5 & 0 \\ 0 & 1 & 2 & 3 & 4 & 5 & 0 \\ 0 & 1 & 2 & 3 & 4 & 5 & 0 \end{array}\right)\end{split}\]
type class DistMatrix<T, STAR, VC>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,STAR,VC>.

type class DistMatrix<double, STAR, VC>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, STAR, VC>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, STAR, VC>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, STAR, VC>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, STAR, VC>

The underlying datatype F is a field.

[VR,* ]

This distribution makes use of a 1d distribution which uses a row-major ordering of the entire process grid. Since 1d distributions are useful for distributing vectors, and a row-major ordering is used, the distribution symbol is VR. Again using the simple \(2 \times 3\) process grid, with a zero column alignment, each entry of a \(7 \times 7\) matrix would be owned by the following sets of processes:

\[\begin{split}\left(\begin{array}{ccccccc} 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 2 & 2 & 2 & 2 & 2 & 2 & 2 \\ 4 & 4 & 4 & 4 & 4 & 4 & 4 \\ 1 & 1 & 1 & 1 & 1 & 1 & 1 \\ 3 & 3 & 3 & 3 & 3 & 3 & 3 \\ 5 & 5 & 5 & 5 & 5 & 5 & 5 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{array}\right)\end{split}\]
type class DistMatrix<T, VR, STAR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,VR,STAR>.

type class DistMatrix<double, VR, STAR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, VR, STAR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, VR, STAR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, VR, STAR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, VR, STAR>

The underlying datatype F is a field.

[* ,VR]

This is the transpose of the above [VR,* ] distribution. On the standard \(2 \times 3\) process grid with a row alignment of zero, a \(7 \times 7\) matrix would be distributed as:

\[\begin{split}\left(\begin{array}{ccccccc} 0 & 2 & 4 & 1 & 3 & 5 & 0 \\ 0 & 2 & 4 & 1 & 3 & 5 & 0 \\ 0 & 2 & 4 & 1 & 3 & 5 & 0 \\ 0 & 2 & 4 & 1 & 3 & 5 & 0 \\ 0 & 2 & 4 & 1 & 3 & 5 & 0 \\ 0 & 2 & 4 & 1 & 3 & 5 & 0 \\ 0 & 2 & 4 & 1 & 3 & 5 & 0 \end{array}\right)\end{split}\]
type class DistMatrix<T, STAR, VR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,STAR,VR>.

type class DistMatrix<double, STAR, VR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, STAR, VR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, STAR, VR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, STAR, VR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, STAR, VR>

The underlying datatype F is a field.

[* ,* ]

This “distribution” actually redundantly stores every entry of the associated matrix on every process. Again using a \(2 \times 3\) process grid, the entries of a \(7 \times 7\) matrix would be owned by the following sets of processes:

\[\begin{split}\left(\begin{array}{ccccccc} \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \\ \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \\ \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \\ \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \\ \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \\ \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \\ \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} & \{0,1,...,5\} \end{array}\right)\end{split}\]
type class DistMatrix<T, STAR, STAR>

TODO: Add the member functions.

Special cases used in Elemental

This list of special cases is here to help clarify the notation used throughout Elemental’s source (as well as this documentation). These are all special cases of DistMatrix<T,STAR,STAR>.

type class DistMatrix<double, STAR, STAR>

The underlying datatype is the set of double-precision real numbers.

type class DistMatrix<Complex<double>, STAR, STAR>

The underlying datatype is the set of double-precision complex numbers.

type class DistMatrix<R, STAR, STAR>

The underlying datatype R is real.

type class DistMatrix<Complex<R>, STAR, STAR>

The underlying datatype Complex<R> is complex with base type R.

type class DistMatrix<F, STAR, STAR>

The underlying datatype F is a field.

«  The Grid class   ::   Contents   ::   Partitioning  »