Matrix Base Class

class fastmat.Matrix

Bases: object

Matrix Base Class

Description: The baseclass of all matrix classes in fastmat. It also serves as wrapper around the standard Numpy Array [1].

H

Return the hermitian transpose

(read-only)

T

Return the transpose of the matrix as fastmat class

(read-only)

__init__()

Initialize an instance of a fastmat matrix.

This is the baseclass for all fastmat matrices and serves as a wrapper to define a matrix based on a two dimensional ndarray. Any specialized matrix type in fastmat is derived from this base class and defines its own __init__.

Every __init__ routine allows the specification of arbitrary keyworded arguments, which are passed in **options. Each specialized __init__ routine processes the options it accepts and passes the rest on to the initialization routines in the base class to define the basic behaviour of the class.

Parameters
arrMatrixnumpy.ndarray

A 2d array representing a dense matrix to be cast as a fastmat matrix.

forceContiguousInputbool, optional

If set, the input array is forced to be contiguous in the style as specified by fortranStyle. If the input array already fulfils the requirement nothing is done.

Defaults to False.

widenInputDatatypebool, optional

If set, the data type of the input array is promoted to at least match the output data type of the operation. Just like the minType option this parameter controls the accumulator width, however dynamically according to the output data type in this case.

Defaults to False.

fortranStylebool, optional

Control the style of contiguousity to be enforced by forceConfiguousInput. If this option is set to True, Fortran-style ordering (contiguous along columns) is enforced, if False C-Style (contiguous along rows).

Defaults to True.

minTypebool, optional

Specify a minimum data type for the input array to a transform. The input array data type will be promoted to at least the data type specified in this option before performing the actual transforms. Using this option is strongly advised for cases where small data types of both input array and matrix could cause range overflows otherwise, as the output data type promotion rules do not consider avoiding accumulator overflows due to performance reasons.

Defaults to numpy.int8.

bypassAllowbool, optional

Allow bypassing the implemented fastmat.Matrix.forward() and fastmat.Matrix.backward() transforms with dense matrix-vector products if runtime estimates suggest this is faster than using the implemented transforms. This requires valid calibration data to be available for the class of the to-be-created instance itself and the fastmat.Matrix base class at the time the new instance is created. If no valid performance calibration data exists this parameter is ignored and the implemented transforms will be used always.

Defaults to the value set in the package-wide fastmat.flags options.

bypassAutoArraybool, optional

Prevents the automatic generation of a dense matrix representation that would be used for bypassing the implemented transforms in case the performance profiles suggest this would be faster, if set to True. This is heavily advised if the matrix is unfeasibly large for a dense representation and does not feature fast transforms.

Defaults to the value as set in the package-wide :py:class`fastmat.flags` if no nested matrix of this instance has set this option to False. If just one has, this parameter defaults to False. If the matrix instance would disregard this, a nested instances’ AutoArray function would be called implicitly through this instances’ dense array constructur although this is disabled for the particular ndested matrix.

array
backward()

Backward Transform

Calculate the backward transform A^mathrm{H}*x where H is the hermitian transpose. Dimension-checking is performed to ensure valid fast transforms as these may succeed even when dimensions do not match. To support both single- and multidimensional input vectors x, single dimensional input will be reshaped to (n, 1) before processing and flattened to (n) after completion. This allows the use of both vectors and arrays. The actual transform code gets called by the callbacks specified in funcPython and funcCython, depending on the state of self._cythonCall.

Warning

Do not override this method

Note

The returned ndarray object may own its data, may be a view into another ndarray and may even be identical to the input array.

Parameters
arrXnumpy.ndarray

The input data array of either 1d or 2d. 1d arrays will be reshaped to 2d during internal processing.

Returns
The result of the operation as np.ndarray with the
same number of dimensions as arrX.
bypassAllow
bypassAutoArray
colNormalized

Return a column normalized matrix for this instance

(read-only)

colNorms

Return the column norms for this matrix instance

(read-only)

complexity

Complexity

(read-only)

Return the computational complexity of all functionality implemented in the class itself, not including calls of external code.

conj

Return the conjugate of the matrix as fastmat class

(read-only)

content
dtype
estimateRuntime()

Estimate the runtime of this matrix instances’ transforms.

Parameters
numVectorsint

Estimate the runtime for processing this number of vectors.

Returns
A tuple containing float estimates on the runtime of the
fastmat.Matrix.forward() and the
fastmat.Matrix.backward() transform if valid performance
profiles are available to this matrix instance. If not, return
(NaN, NaN)
forward()

Forward

Calculate the forward transform A * x. Dimension-checking is performed to ensure valid fast transforms as these may succeed even when dimensions do not match. To support both single- and multidimensional input vectors x, single dimensional input will be reshaped to (n, 1) before processing and flattened to (n) after completion. This allows the use of both vectors and arrays. The actual transform code gets called by the callbacks specified in funcPython and funcCython, depending on the state of self._cythonCall.

Warning

Do not override this method!

Note

The returned ndarray object may own its data, may be a view into another ndarray and may even be identical to the input array.

Parameters
arrXnumpy.ndarray

The input data array of either 1d or 2d. 1d arrays will be reshaped to 2d during internal processing.

Returns
The result of the operation as np.ndarray with the
same number of dimensions as arrX.
fusedType
getArray()

Return a dense array representation of this matrix.

getCol()

Return a column by index.

Parameters
idxint

Index of the column to return.

Returns
1d-numpy.ndarray holding the specified column.
getColNormalized()

Return a column normalized version of this matrix as fastmat matrix.

getColNorms()

Return a column normalized version of this matrix as fastmat matrix.

getCols()

Return a set of columns by index.

Parameters
indicesint, slice or numpy.ndarray

If an integer is given, this is equal to the output of getCol`(indices). If a 1D vector or slice is given, a 2D :py:class:`numpy.ndarray() containing the columns as requested by indices is returned.

Returns
1D or 2D (depending on type of indices) numpy.ndarray
holding the specified column(s).
getComplexity()

Return a transform complexity estimate for this matrix instance.

Returns a tuple containing the complexity estimates for the fastmat.Matrix.forward() and fastmat.Matrix.backward() transforms (in that order).

getConj()

Return the conjugate of this matrix as fastmat matrix.

getGram()

Return the gramian of this matrix as fastmat matrix.

getH()

Return the hermitian transpose of this matrix as fastmat matrix.

getInverse()

Return the hermitian transpose of this matrix as fastmat matrix.

getLargestEigenValue()

Largest Singular Value

For a given matrix \(A \in \mathbb{C}^{n \times n}\), so \(A\) is square, we calculate the absolute value of the largest eigenvalue \(\lambda \in \mathbb{C}\). The eigenvalues obey the equation

\[A \cdot v= \lambda \cdot v,\]

where \(v\) is a non-zero vector.

Input matrix \(A\), parameter \(0 < \varepsilon \ll 1\) as a stopping criterion Output largest eigenvalue \(\sigma_{\rm max}( A)\)

Note

This algorithm performs well if the two largest eigenvalues are not very close to each other on a relative scale with respect to their absolute value. Otherwise it might get trouble converging properly.

>>> # import the packages
>>> import numpy.linalg as npl
>>> import numpy as np
>>> import fastmat as fm
>>>
>>> # define the matrices
>>> n = 5
>>> H = fm.Hadamard(n)
>>> D = fm.Diag(np.linspace
>>>         1, 2 ** n, 2 ** n))
>>>
>>> K1 = fm.Product(H, D)
>>> K2 = K1.array
>>>
>>> # calculate the eigenvalue
>>> x1 = K1.largestEigenValue
>>> x2 = npl.eigvals(K2)
>>> x2 = np.sort(np.abs(x2))[-1]
>>>
>>> # check if the solutions match
>>> print(x1 - x2)

We define a matrix-matrix product of a Hadamard matrix and a diagonal matrix. Then we also cast it into a numpy-array and use the integrated EVD. For demonstration, try to increase \(n\)>10`and see what happens.

getLargestEigenVec()
getLargestSingularValue()

Largest Singular Value

For a given matrix \(A \in \mathbb{C}^{n \times m}\), we calculate the largest singular value \(\sigma_{\rm max}( A) > 0\), which is the largest entry of the diagonal matrix \(\Sigma \in \mathbb{C}^{n \times m}\) in the decomposition

\[A = U \Sigma V^{\rm H},\]

where \(U\) and \(V\) are matrices of the appropriate dimensions. This is done via the so called power iteration of \(A^{\rm H} \cdot A\).

  • Input matrix \(A\), parameter \(0 < \varepsilon \ll 1\) as a stopping criterion

  • Output largest singular value \(\sigma_{\rm max}( A)\)

Note

This algorithm performs well if the two largest singular values are not very close to each other on a relative scale. Otherwise it might get trouble converging properly.

>>> # import packages
>>> import numpy.linalg as npl
>>> import numpy as np
>>> import fastmat
>>>
>>> # define involved matrices
>>> n = 5
>>> H = fm.Hadamard(n)
>>> F = fm.Fourier(2**n)
>>> K1 = fm.Kron(H, F)
>>> K2 = K1
>>>
>>> # calculate the largest SV
>>> # and a reference solution
>>> x1 = largestSingularValue(K1.largestSingularValue
>>> x2 = npl.svd(K2,compute_uv
>>> # check if they match
>>> print(x1-x2)

We define a Kronecker product of a Hadamard matrix and a Fourier matrix. Then we also cast it into a numpy-array and use the integrated SVD. For demonstration, try to increase \(n\) to >10 and see what happens.

Returns
The largest singular value
getLargestSingularVectors()
getPseudoInverse()

Return the hermitian transpose of this matrix as fastmat matrix.

getRow()

Return a row by index.

Parameters
idxint

Index of the row to return.

Returns
1d-numpy.ndarray holding the specified row.
getRowNormalized()

Return a column normalized version of this matrix as fastmat matrix.

getRowNorms()

Return a row normalized version of this matrix as fastmat matrix.

getRows()

Return a set of rows by index.

Parameters
indicesint, slice or numpy.ndarray

If an integer is given, this is equal to the output of getRow`(indices). If a 1D vector or slice is given, a 2D :py:class:`numpy.ndarray() containing the rows as requested by indices is returned.

Returns
1D or 2D (depending on type of indices) numpy.ndarray
holding the specified row(s).
getScipyLinearOperator()
getT()

Return the transpose of this matrix as fastmat matrix.

gram

Return the gram matrix for this fastmat class

(read-only)

inverse

Return the inverse

(read-only)

largestEV
largestEigenValue

Return the largest eigenvalue for this matrix instance

(read-only)

largestEigenVec

Return the vector corresponding to the largest eigen value

(read-only)

largestSV
largestSingularValue

Return the largestSingularValue for this matrix instance

(read-only)

largestSingularVectors

Return the vectors corresponding to the largest singular value

This property returns a tuple (u, v) of the first columns of U and V in the singular value decomposition of

\[A = U \Sigma V^{\rm H},\]

which means that the tuple contains the leading left and right singular vectors of the matrix

(read-only)

nbytes
nbytesReference
next()

Stop iteration as __iter__ redirected here. Python2-Style.

normalized
numCols
numM
numN
numRows
numpyType
profileBackward
profileForward
pseudoInverse

Return the moore penrose inverse

(read-only)

reference()

Return explicit array reference of this matrix instance.

Return an explicit representation of the matrix without using any fastmat code. Provides type checks and raises errors if the matrix type (self.dtype) cannot hold the reference data. This implementation is meant to provide a reference version for testing and MUST not use any fastmat code for its implementation.

Returns
The array representation of this matrix instance as 2d
np.ndarray.
rowNormalized

Return a column normalized matrix for this instance

(read-only)

rowNorms

Return the row norms for this matrix instance

(read-only)

scipyLinearOperator

Return a Representation as scipy’s linear Operator

This property allows to make use of all the powerfull algorithms provided by scipy, that allow passing a linear operator to them, like optimization routines, system solvers or decomposition algorithms.

(read-only)

shape
tag