Any sparse matrix: SparseMatrix, ListMatrix<SparseVector>, or an alias object referring to a mutable sparse matrix.
SparseMatrix deploys an adaptive reallocation strategy similar to std::vector, reserving additional stock memory by every reallocation. This way, the amortized reallocation costs are proportional to the logarithm of the number of appended rows/columns.
The real result type is a temporary object of type pm::RepeatedRow or pm::RepeatedCol. It contains only a single copy of a vector, or just a reference to it, but makes it to appear as a sequence of identical row (column) vectors.
Columns in the only_rows case, rows in the only_cols case.
The real result type is a temporary alias object of type pm::VectorSlice. It inherits the const-ness from the original matrix.
This assertion is checked only in the DIMENSION_CHECKS compilation mode.
Random-access matrix: Matrix, block matrix built of random-access matrices, or a minor of a random-access matrix with random-access row and column index sets (sequence or series).
The concrete matrix type hidden behind the GenericMatrix.
The real result type is a masquerade reference pm::Transposed&, referring to the originator matrix. It inherits the const-ness from the latter.
Either a functor class satisfying the extended requirements for binary operations, or an untyped template of such a class wrapped into BuildBinary or BuildBinaryIt (preferably expressed via the convenience typedef from namespace polymake::operations.
Rows in the only_rows case, columns in the only_cols case.
Any combination of GenericSet<...,int> , Complement<...,int> , and a special identifier All (choosing all rows/columns) is allowed here.
The real return type is a masquerade reference const pm::Transposed<SparseMatrix2x2>& .
The real result type is a temporary object of type pm::DiagMatrix.
Any persistent dense class: Matrix, or ListMatrix built with dense row vectors.
Random-column-access matrix: Matrix, SparseMatrix, block matrix built of random-column-access matrices, or a minor of a random-column-access matrix with a random-access column index set (sequence or series).
The real result type is a temporary alias object of type pm::MatrixMinor. It inherits the const-ness from the originator matrix.
Random-row-access matrix: Matrix, SparseMatrix, block matrix built of random-row-access matrices, or a minor of a random-row-access matrix with a random-access row index set (sequence or series).
The real return type is a temporary object (expression template) implementing the lazy evaluation of the operation result. You should have it in mind when passing the result of a matrix expression as a GenericMatrix argument to a template function where it has a chance to be evaluated multiple times. See also prevent_lazy class.
The vector and matrix operators defined in the Polymake Template Library do take care of this, so this remark applies more likely to your own functions.
Any persistent matrix: Matrix, SparseMatrix, ListMatrix, or an alias object referring to a mutable (non-const) martix.
The real result type is a masquerade reference pm::SingleRow& or SingleCol&, pointing to the source vector.
This assertion is checked only in the INDEX_CHECKS compilation mode.
Any persistent class: Matrix, SparseMatrix, or ListMatrix.
Either a functor class satisfying the extended requirements for unary operations, or an untyped template of such a class wrapped into BuildUnary or BuildUnaryIt (preferably expressed via the convenience typedef from namespace polymake::operations.
The real result type is a temporary object of type pm::SameElementSparseMatrix. It contains only a single copy of an element, or just a reference to an element, but makes it to appear as a matrix of identical elements.
The real result type is a temporary alias object of type pm::RowChain or pm::ColChain. It is mutable only if both operands are mutable too.
The real return type is a temporary proxy object, forwarding the read/write access operations to the sparse matrix. Prior to the destruction, which happens after the completion of the expression envolving the random access, it checks whether the element became zero; if so, it is removed from the matrix.

Prerequisits

#include <Matrix.h>
#include <SparseMatrix.h>
#include <ListMatrix.h>
using namespace polymake;

Introduction

The Matrix class family implement matrices in the usual algebraic notion. It contains three persistent matrix types with different data storage organization:

template <typename ElementType> class Matrix;

holds the elements in a contiguous array.

template <typename ElementType> class SparseMatrix;

is a two-dimensional associative array with row and column indices as keys; elements equal to the default value (ElementType(), which is 0 for most numerical types) are not stored, but implicitly encoded by the gaps in the key set. Each row and column is organized as a balanced binary search (AVL) tree.

template <typename VectorType> class ListMatrix;

is a list of row vectors. Based on std::list. Can be parameterized with any persistent vector type.

The following brief analysis might help you when choosing the most efficient representation for a concrete case. r and c are number of rows and columns respectively, and n is the average number of non-zero elements in a row/column.

Performance comparison
Operation Matrix SparseMatrix ListMatrix
rowwise iteration O(r) with approx. equal overhead
columnwise iteration O(c) O(c) O(c), but with rather high overhead
random row/column access O(1) O(1)
random element access O(1) O(log(n))
append a row O(r c) + reallocation costs O(r) + reallocation costs + n O(log(n)) O(1)
append a column O(r c) + reallocation costs
O(c) + reallocation costs + n O(log(n))
See also the special class RestrictedSparseMatrix.
c*(appending costs for a VectorType)

Matrix and SparseMatrix keep the data attached to a smart pointer with REFCOUNT.

Generic type

template <typename _Matrix, typename ElementType=typename _Matrix::element_type> class GenericMatrix;

The GenericClass for all matrix classes. Defines the following types and constants:

element_type ElementType, the persistent element type
top_type _Matrix, the concrete matrix object behind the generic reference
generic_type GenericMatrix<_Matrix> itself
persistent_type Matrix<ElementType> for dense matrices,
SparseMatrix<ElementType> for sparse matrices.
row_type
col_type
Alias classes representing single matrix rows and columns respectively. They always belong to the GenericVector family.
bool is_sparse true if the element sequence in an arbitrary row or column can contain gaps expressing implicit zero elements.

Methods and friend functions applicable to GenericMatrix directly are scattered over the entire page. They are discernible on the prefix GenericMatrix or parameters of type GenericMatrix&. To use other methods via a generic reference (or pointer), you must first move to the concrete object using the top() method.

Construction

Matrix();
Create a matrix with 0 rows and 0 columns.
Matrix (int m, int n);
Create a matrix with m rows and n columns, initialize all elements to 0 (in sparse cases implicitly.)
Matrix (int m, int n, const ElementType& value);
Create a matrix with m rows and n columns, initialize all elements with value.
template <typename Iterator> Matrix (int m, int n, Iterator src);
Create a matrix with m rows and n columns, initialize the elements from a data sequence. src should iterate either over m*n scalar values, corresponding to the elements in the row order (the column index changes first,) or over m vectors of dimension n, corresponding to the matrix rows. In the sparse case, zero input elements are filtered out.
template <typename OtherMatrix> Matrix (const GenericMatrix<OtherMatrix, ElementType>&);
Create a matrix as a copy of another matrix of the same element type.
template <typename OtherMatrix, typename OtherElementType> explicit Matrix (const GenericMatrix<OtherMatrix, OtherElementType>&);
Copying with element conversion is also possible, but requires an explicit constructor call.
GenericMatrix convert_to<OtherElementType> (const GenericMatrix&);
Create a matrix with another element type, using explicit type conversion.

Modification

template <typename OtherMatrix> GenericMatrix::operator= (const GenericMatrix<OtherMatrix, ElementType>&);
Persistent matrix objects have after the assignment the same dimensions as the right hand side operand. Alias objects, such as matrix minor or block matrix, cannot be resized, thus must have the same dimensions as on the right hand side.
void std::swap(GenericMatrix&, GenericMatrix&);
Exchange the contents of two matrices in a most efficient way. If at least one non-persistent object is involved, the operands must have equal dimensions.
void Matrix::resize(int m, int n);
Extend or truncate to new dimensions (m rows, n columns.) Surviving elements keep their values, new elements are initialized with 0 (in sparse cases implicitly.)
SparseMatrix deploys an adaptive reallocation strategy similar to std::vector, reserving additional stock memory by every reallocation. A special case, looking at the first glance like a "no operation": M.resize(M.rows(), M.cols()) , gets rid of this extra allocated storage.
void Matrix::clear();
Truncate to 0x0 matrix.
template <typename Scalar> void GenericMatrix::fill(const Scalar& value);
Fill with given value without changing the dimensions. value can be of arbitrary type assignable to ElementType.
void GenericMatrix::remove0s();
Physically remove all zero elements from a sparse matrix that might have creeped in by some previous operation.

The following methods are defined for SparseMatrix objects only.

void squeeze(); template <typename RowIndexConsumer> void squeeze(RowIndexConsumer row_index_consumer); template <typename RowIndexConsumer, typename ColIndexConsumer> void squeeze(RowIndexConsumer row_index_consumer, ColIndexConsumer col_index_consumer);
Remove all empty (i.e., consisting entirely of implicit zeroes,) rows and columns, renumber the rest, and reduce the dimensions. If you need to know the exact mapping between the old and new row (column) indices, you can supply output iterators (e.g., back_inserter of a std::list.) They will get the old indices of non-empty rows (columns) assigned in the ascending order.
template <typename Iterator> void permute_rows(Iterator perm); template <typename Iterator> void permute_cols(Iterator perm); template <typename Iterator> void permute_inv_rows(Iterator inv_perm); template <typename Iterator> void permute_inv_cols(Iterator inv_perm);
Permute the rows(columns) of the matrix without copying the elements. These operations are nevetherless expensive, as they need to visit each element and adjust its indices.
An iterator argument should run over a sequence of integer indices. In the straight variant (perm), it has the same meaning as in the select function. In the inverted variant (inv_perm), it is an inverse permutation: the first element specifies the new position of the 0-th row (column) etc.
If the index sequence does not build a proper permutation, the results are undefined. It is the responsibility of the application context to assure this condition.

The following methods are defined for ListMatrix objects only.

typename Rows<ListMatrix>::iterator insert_row (typename Rows<ListMatrix>::iterator where, const GenericVector&);
Insert a new row at the given position: before the row pointed to by where.
void ListMatrix::delete_row (typename Rows<ListMatrix>::iterator where);
Delete the row pointed to by where. Be aware that the iterator becomes void after the removal unless rescued in good time (with postincrement or postdecrement).

Element access

int GenericMatrix::rows() const; int GenericMatrix::cols() const;
Tell the current matrix dimensions.

Sequential access

The matrix classes as such do not implement any STL container interfaces: they are simply not defined for multi-dimensional containers. However, a matrix can be disguised in three different ways, yielding STL-conforming containers:

pseudo-container class producing function description
Rows<_Matrix>
const Rows<_Matrix>
rows(GenericMatrix&);
rows(const GenericMatrix&);
A sequence of row vectors; each element is a GenericVector.
Cols<_Matrix>
const Cols<_Matrix>
cols(GenericMatrix&);
cols(const GenericMatrix&);
A sequence of column vectors; each element is a GenericVector.
ConcatRows<_Matrix>
const ConcatRows<_Matrix>
concat_rows(GenericMatrix&)
concat_rows(const GenericMatrix&);
A sequence of elements in rowwise order (the column index changes first). The ConcatRows object itself is a GenericVector of dimension rows()*cols().

SparseMatrix specialities

int Rows<SparseMatrix>::iterator::index(); int Cols<SparseMatrix>::iterator::index();
Tell the current row (column) number. The numbers start, as usual, with 0.
int SparseMatrix::row_type::iterator::index(); int SparseMatrix::col_type::iterator::index();
Tell the current column (row) number of the element pointed to by the given iterator. Remember that rows and columns of SparseMatrix are sparse containers, so the iterators over them visit only the (presumably non-zero) elements physically stored in the matrix.
SparseMatrix::col_type::iterator cross_direction(const SparseMatrix::row_type::iterator&); SparseMatrix::row_type::iterator cross_direction(const SparseMatrix::col_type::iterator&);
For an element of a sparse matrix pointed to by the given iterator over the row, create an iterator over the column pointing to the same element (and vice versa.)

All these methods are, indeed, defined for the corresponding const_iterator classes too.

Random access

ElementType& Matrix::operator() (int i, int j); const ElementType& Matrix::operator() (int i, int j) const;
Get the element in i-th row, j-th column. The indices must lie in the valid ranges.
ElementType& SparseMatrix::operator() (int i, int j); const ElementType& SparseMatrix() (int i, int j) const;
The const operator returns a reference to a static zero object if the specified element is an implicit 0.
The indices must lie in the valid ranges.
GenericVector Matrix::row(int i); GenericVector Matrix::operator[] (int i); GenericVector Matrix::col(int i);
Get the i-th row or column respectively. The const variants of these methods create immutable slice objects.
Rows and columns of a SparseMatrix offer the binary tree data access methods additionally to the GenericVector interface.
The index must lie in the valid range.
operator[] is a mere synonym for row(int), and is defined with the sole purpose to allow the C-style expressions M[i][j]. Please be aware that M(i,j) is a much more efficient way to access a single matrix element, as it doesn't involve creation of a temporary vector row object.

Operations

const GenericMatrix operator- (const GenericMatrix&); const GenericMatrix operator+ (const GenericMatrix&, const GenericMatrix&); const GenericMatrix operator- (const GenericMatrix&, const GenericMatrix&); const GenericMatrix operator* (const Scalar&, const GenericMatrix&); const GenericMatrix operator* (const GenericMatrix&, const Scalar&); const GenericMatrix operator/ (const GenericMatrix&, const Scalar&); const GenericMatrix operator* (const GenericMatrix&, const GenericMatrix&); const GenericMatrix operator* (const GenericMatrix&, const GenericVector&); const GenericMatrix operator* (const GenericVector&, const GenericMatrix&);
Do exactly what the mnemonics suggest. Every type combination is allowed, including mixing of sparse and dense vectors and matrices. ElementType of matrix and vector operands, and Scalar type can be different, provided the arithmetic operator with corresponding arguments exists.
The dimensions of the operands must match.
_Matrix& GenericMatrix::negate(); _Matrix& GenericMatrix::operator+= (const GenericMatrix&); _Matrix& GenericMatrix::operator-= (const GenericMatrix&); _Matrix& GenericMatrix::operator*= (const Scalar&); _Matrix& GenericMatrix::operator/= (const Scalar&);
Corresponding assignment versions of the arithmetic operators.
bool operator== (const GenericMatrix&, const GenericMatrix&); bool operator!= (const GenericMatrix&, const GenericMatrix&); bool operator< (const GenericMatrix&, const GenericMatrix&); bool operator> (const GenericMatrix&, const GenericMatrix&); bool operator<= (const GenericMatrix&, const GenericMatrix&); bool operator>= (const GenericMatrix&, const GenericMatrix&);
Compare two matrices lexicographically (iterating over elements in row order). The dimensions of the operands must match.
GenericMatrix::operator bool() const; bool GenericMatrix::operator! () const;
Look for a non-zero element.

Other constructions

GenericMatrix GenericMatrix::minor(const RowIndexSet& row_indices, const ColIndexSet& column_indices);
Select a matrix minor lying on the intersection of the given row and column subsets. The const variant of this method creates an immutable minor object.
The indices must lie in the valid ranges.
GenericVector GenericMatrix::diagonal (int i=0); GenericVector GenericMatrix::anti_diagonal (int i=0);
Select the diagonal or anti-diagonal of a matrix:
i =0 the main diagonal (anti-diagonal)
i>0 the i-th diagonal below the main
i<0 the (-i)-th diagonal above the main
The const variants of these methods create immutable slice objects.
GenericMatrix T (GenericMatrix&);
Transpose the matrix. The roles of Rows and Cols get swapped.
GenericMatrix operator/ (const GenericMatrix&, const GenericMatrix&); GenericMatrix operator/ (const GenericVector&, const GenericMatrix&); GenericMatrix operator/ (const GenericMatrix&, const GenericVector&); GenericMatrix operator/ (const GenericVector&, const GenericVector&);
Create a block matrix, virtually appending the rows of the second matrix after the last row of the first matrix. A vector argument is treated as a matrix with one row. Column dimensions of the operands must be equal.
ElementType of the operands must be identical, otherwise you will get a rather cryptic message from the compiler.
Sorry for a rather weird mnemonics...
GenericMatrix operator| (const GenericMatrix&, const GenericMatrix&); GenericMatrix operator| (const GenericVector&, const GenericMatrix&); GenericMatrix operator| (const GenericMatrix&, const GenericVector&);
Create a block matrix, virtually appending the columns of the second matrix after the last column of the first matrix. A vector argument is treated as a matrix with one column. Row dimensions of the operands must be equal.
ElementType of the operands must be identical, otherwise you will get a rather cryptic message from the compiler.
Matrix& Matrix::operator/= (const GenericMatrix&); Matrix& Matrix::operator/= (const GenericVector&); Matrix& Matrix::operator|= (const GenericMatrix&); Matrix& Matrix::operator|= (const GenericVector&);
Append rows or columns to the matrix. All constraints for the non-assigning operators apply here too.
const GenericMatrix diag(const GenericMatrix&, const GenericMatrix&); const GenericMatrix diag(const GenericMatrix&, const GenericVector&); const GenericMatrix diag(const GenericVector&, const GenericMatrix&);
Create a block-diagonal matrix. Vector arguments are treated as square diagonal matrices.
GenericMatrix vector2row(const GenericVector& v); GenericMatrix vector2col(const GenericVector& v);
Disguise v as a matrix with 1 row (column).
const GenericMatrix diag(const GenericVector& v);
Create a square diagonal matrix. The elements of v appear on the main diagonal.
const GenericMatrix repeat_row(const GenericVector& v, int n); const GenericMatrix repeat_col(const GenericVector& v, int n);
Create a matrix with n rows (columns), each equal to v.
const GenericMatrix same_element_matrix(const ElementType& x, int m, int n); const GenericMatrix ones_matrix<ElementType>(int m, int n); const GenericMatrix zero_matrix<ElementType>(int m, int n);
Create a dense matrix with m rows and n columns, with all elements equal to x, 1, or 0, respectively. Note the obligatory explicit specification of the template parameter ElementType, where it cannot be deduced from the function arguments.
const GenericMatrix same_element_sparse_matrix<ElementType>(const GenericIncidenceMatrix& M); const GenericMatrix same_element_sparse_matrix(const GenericIncidenceMatrix& M, const ElementType& x);
Create a sparse matrix with the same structure as M, with elements equal to x in true cells and 0 in false cells.
If omitted, x is taken equal to 1; note the obligatory explicit template parameter specification in this case.
const GenericMatrix unit_matrix<ElementType>(int n);
Create a unit matrix n×n.

Restricted sparse matrix

enum sparse2d::restriction_kind { only_rows, only_cols }; template <typename ElementType, sparse2d::restriction_kind restriction=only_rows> class RestrictedSparseMatrix;

A special class allowing for efficient incremental building of sparse matrices.

If you have to fill a sparse matrix element for element, but do not know one of its dimensions in advance, you would have to resize it repeatedly, which would severely impact the performance. This special class maintains the internal data only for one direction, and thus can be filled cheaper than SparseMatrix.

On the other hand, almost all matrix operations described above need information about both directions (rows and columns.) Therefore, RestrictedSparseMatrix is deliberately not included in the GenericMatrix class family. Its only puspose is to gather elements and finally hand them over to a SparseMatrix object (without copying!)

RestrictedSparseMatrix();
Create a matrix with 0 rows and 0 columns.
explicit RestrictedSparseMatrix(int n);
Create a matrix with n rows and 0 columns (or vice versa in only_cols case), all elements are implicit zeroes.
template <typename Iterator> RestrictedSparseMatrix(int n, Iterator src);
Create a matrix with n rows and 0 columns (or vice versa in only_cols case), initialize the rows (columns) from the input sequence. Each item supplied by src should be a GenericVector.
void RestrictedSparseMatrix::clear();
Truncate to a 0x0 matrix.
void std::swap(RestrictedSparseMatrix&, RestrictedSparseMatrix&);
Swap the contents of two matrices in a most efficient way.
int RestrictedSparseMatrix::rows() const; int RestrictedSparseMatrix::cols() const;
Tell the current matrix dimensions. For the unsupported direction, the value returned is the highest column (row) index ever seen in an element access or modification operation.
ElementType& RestrictedSparseMatrix::operator() (int i, int j); const ElementType& RestrictedSparseMatrix() (int i, int j) const;
Random access to an element. The const operator returns a reference to a static zero object if the specified element is an implicit 0.
The index in the main direction must lie in the valid range, the opposite index can cause the matrix to grow dynamically.
Rows<RestrictedSparseMatrix>& rows(RestrictedSparseMatrix&); const Rows<RestrictedSparseMatrix>& rows(const RestrictedSparseMatrix&); Cols<RestrictedSparseMatrix>& cols(RestrictedSparseMatrix&); const Cols<RestrictedSparseMatrix>& cols(const RestrictedSparseMatrix&);
Sequential access to the rows and columns. Defined only for the main direction.
GenericVector RestrictedSparseMatrix::row (int i); GenericVector RestrictedSparseMatrix::col (int i);
Random access to a row (column). Defined only for the main direction.
RestrictedSparseMatrix& RestrictedSparseMatrix::operator/= (const GenericMatrix&); RestrictedSparseMatrix& RestrictedSparseMatrix::operator/= (const GenericVector&);
Append rows to the matrix. Allowed in both only_rows and only_cols modes.
RestrictedSparseMatrix& RestrictedSparseMatrix::operator|= (const GenericMatrix&); RestrictedSparseMatrix& RestrictedSparseMatrix::operator|= (const GenericVector&);
Append columns to the matrix. Allowed in both only_rows and only_cols modes.
explicit SparseMatrix(RestrictedSparseMatrix<ElementType>&); SparseMatrix& SparseMatrix::operator=(RestrictedSparseMatrix<ElementType>&);
The final and alone destination of a RestrictedSparseMatrix object is to hand over its internal data to a new created SparseMatrix. Any further action on the RestrictedSparseMatrix object (except destruction, obviously) is prohibited: it would lead to the immediate program crash.

SparseMatrix2x2

template <typename ElementType> struct SparseMatrix2x2 { int i,j; ElementType a_ii, a_ij, a_ji, a_jj; };

A compact encoding of matrices of special structure, which occur as companion matrices in elimination, triangularization, etc. algorithms. Note that this class is not derived from GenericMatrix.

For two given integer indices i and j, it encodes a sparse square matrix A of arbitrary dimensions with 4 designated elements Ai,i, Ai,j, Aj,i, and Aj,j.

Other elements are implicit zeroes except for the rest of the main diagonal, which is filled with ones.

SparseMatrix2x2(int _i, int _j, const ElementType& _a_ii, const ElementType& _a_ij, const ElementType& _a_ji, const ElementType& _a_jj);
The only constructor for this class.
SparseMatrix2x2 T (const SparseMatrix2x2&); ElementType det(const SparseMatrix2x2&);
Elementary operations: transposition, determinant.
void GenericMatrix::multiply_from_left(const SparseMatrix2x2& U);
Multiply (*this) with U from the left, that is, apply U to the rows.
void GenericMatrix::multiply_from_right(const SparseMatrix2x2& U);
Multiply (*this) with U from the right, that is, apply U to the columns.