PyArmadillo: linear algebra library for Python

[top] API Documentation for PyArmadillo 0.500

Preamble

This documentation assumes PyArmadillo was imported using
from pyarma import *

If PyArmadillo was imported using "import pyarma as pa"
then the pa. prefix must be added to each class and function;
for example: pa.mat

See the syntax conversion table for differences between Matlab and PyArmadillo

First time users: please see the short example program

If you discover any bugs or regressions, please report them

History of API additions

Please cite the following article if you use PyArmadillo in your research and/or software.
Citations are useful for the continued development and maintenance of the library.

Jason Rumengan, Terry Yue Zhuo, Conrad Sanderson.
PyArmadillo: a streamlined linear algebra library for Python.
Journal of Open Source Software, Vol. 6, No. 66, 2021.

Overview

matrix and cube classes
member functions & variables

generated vectors/matrices/cubes
functions of vectors/matrices/cubes

decompositions, factorisations, inverses and equation solvers

signal & image processing
statistics & clustering
miscellaneous

Matrix and Cube Classes

mat, fmat, cx_mat cx_fmat, umat, imat		dense matrix classes

cube, fcube, cx_cube cx_fcube, ucube, icube		dense cube class ("3D matrix")

operators		`+ − * @ / == != <= >= < > && \|\|`

Member Functions & Variables

attributes		.n_rows, .n_cols, .n_elem, .n_slices, ...
element access		element/object access via (), [] and .at()
element initialisation		set elements via initialiser lists

.zeros		set all elements to zero
.ones		set all elements to one
.eye		set elements along main diagonal to one and off-diagonal elements to zero
.randu / .randn		set all elements to random values

.fill		set all elements to specified value
.imbue		imbue (fill) with values provided by functor or lambda function
.clean		replace elements below a threshold with zeros

.replace		replace specific elements with a new value
.transform		transform each element via functor or lambda function
.for_each		apply a functor or lambda function to each element

.set_size		change size without keeping elements (fast)
.reshape		change size while keeping elements
.resize		change size while keeping elements and preserving layout
.copy_size		change size to be same as given object
.reset		change size to empty

submatrix views		read/write access to contiguous and non-contiguous submatrices
subcube views		read/write access to contiguous and non-contiguous subcubes
[diag]		read/write access to matrix diagonals

.set_imag / .set_real		set imaginary/real part
.insert_rows/cols/slices		insert vector/matrix at specified row/column
.shed_rows/cols/slices		remove specified rows/columns
.swap_rows/cols		swap specified rows or columns
.swap		swap contents with given object

iterators (matrices)		iterators and associated member functions for dense matrices and vectors
iterators (cubes)		iterators and associated member functions for cubes
iterators (submatrices)		iterators and associated member functions for submatrices

.as_col / .as_row		return flattened matrix as column or row vector
.t / .st		return matrix transpose
.i		return inverse of square matrix
.min / .max		return extremum value
.index_min / .index_max		return index of extremum value
.eval		force evaluation of delayed expression

.in_range		check whether given location or span is valid
.is_empty		check whether object is empty
.is_vec		check whether matrix is a vector
.is_sorted		check whether vector or matrix is sorted

.is_trimatu / .is_trimatl		check whether matrix is upper/lower triangular
.is_diagmat		check whether matrix is diagonal
.is_square		check whether matrix is square sized
.is_symmetric		check whether matrix is symmetric
.is_hermitian		check whether matrix is hermitian
.is_sympd		check whether matrix is symmetric/hermitian positive definite

.is_zero		check whether all elements are zero
.is_finite		check whether all elements are finite
.has_inf		check whether any element is +-Inf
.has_nan		check whether any element is NaN

.print		print object
.brief_print		print object in abridged form

.save/.load (matrices and cubes)		save/load matrices and cubes in files

Generated Vectors/Matrices/Cubes

linspace		generate vector with linearly spaced elements
logspace		generate vector with logarithmically spaced elements
regspace		generate vector with regularly spaced elements
randperm		generate vector with random permutation of a sequence of integers
eye		generate identity matrix
ones		generate object filled with ones
zeros		generate object filled with zeros
randu / randn		generate object with random values (uniform and normal distributions)
randg		generate object with random values (gamma distribution)
randi		generate object with random integer values in specified interval
toeplitz		generate Toeplitz matrix

Functions of Vectors/Matrices/Cubes

abs		obtain magnitude of each element
accu		accumulate (sum) all elements
affmul		affine matrix multiplication
all		check whether all elements are non-zero, or satisfy a relational condition
any		check whether any element is non-zero, or satisfies a relational condition
approx_equal		approximate equality
arg		phase angle of each element
as_scalar		convert 1x1 matrix to pure scalar
clamp		obtain clamped elements according to given limits
cond		condition number of matrix
conj		obtain complex conjugate of each element
cross		cross product
cumsum		cumulative sum
cumprod		cumulative product
det		determinant
diagmat		generate diagonal matrix from given matrix or vector
diagvec		extract specified diagonal
diff		differences between adjacent elements
dot / cdot / norm_dot		dot product
eps		obtain distance of each element to next largest floating point representation
expmat		matrix exponential
expmat_sym		matrix exponential of symmetric matrix
find		find indices of non-zero elements, or elements satisfying a relational condition
find_finite		find indices of finite elements
find_nonfinite		find indices of non-finite elements
find_unique		find indices of unique elements
fliplr / flipud		flip matrix left to right or upside down
imag / real		extract imaginary/real part
ind2sub		convert linear index to subscripts
index_min / index_max		indices of extremum values
inplace_trans		in-place transpose
intersect		find common elements in two vectors/matrices
join_rows / join_cols		concatenation of matrices
join_slices		concatenation of cubes
kron		Kronecker tensor product
log_det		log determinant
logmat		matrix logarithm
logmat_sympd		matrix logarithm of symmetric matrix
min / max		return extremum values
nonzeros		return non-zero values
norm		various norms of vectors and matrices
normalise		normalise vectors to unit p-norm
prod		product of elements
powmat		matrix power
rank		rank of matrix
rcond		reciprocal of condition number
repelem		replicate elements
repmat		replicate matrix in block-like fashion
reshape		change size while keeping elements
resize		change size while keeping elements and preserving layout
reverse		reverse order of elements
roots		roots of polynomial
shift		shift elements
shuffle		randomly shuffle elements
size		obtain dimensions of given object
sort		sort elements
sort_index		vector describing sorted order of elements
sqrtmat		square root of matrix
sqrtmat_sympd		square root of symmetric matrix
sum		sum of elements
sub2ind		convert subscripts to linear index
symmatu / symmatl		generate symmetric matrix from given matrix
trace		sum of diagonal elements
trans		transpose of matrix
trapz		trapezoidal numerical integration
trimatu / trimatl		copy upper/lower triangular part
trimatu_ind / trimatl_ind		obtain indices of upper/lower triangular part
unique		return unique elements
vectorise		flatten matrix into vector
misc functions		miscellaneous element-wise functions: exp, log, pow, sqrt, round, sign, ...
trig functions		trigonometric element-wise functions: cos, sin, ...

Decompositions, Factorisations, Inverses and Equation Solvers (Dense Matrices)

chol		Cholesky decomposition
eig_sym		eigen decomposition of dense symmetric/hermitian matrix
eig_gen		eigen decomposition of dense general square matrix
eig_pair		eigen decomposition for pair of general dense square matrices
hess		upper Hessenberg decomposition
inv		inverse of general square matrix
inv_sympd		inverse of symmetric positive definite matrix
lu		lower-upper decomposition
null		orthonormal basis of null space
orth		orthonormal basis of range space
pinv		pseudo-inverse
qr		QR decomposition
qr_econ		economical QR decomposition
qz		generalised Schur decomposition
schur		Schur decomposition
solve		solve systems of linear equations
svd		singular value decomposition
svd_econ		economical singular value decomposition
syl		Sylvester equation solver

Signal & Image Processing

conv		1D convolution
conv2		2D convolution
fft / ifft		1D fast Fourier transform and its inverse
fft2 / ifft2		2D fast Fourier transform and its inverse
interp1		1D interpolation
interp2		2D interpolation
polyfit		find polynomial coefficients for data fitting
polyval		evaluate polynomial

Statistics & Clustering

stats functions		mean, median, standard deviation, variance
cov		covariance
cor		correlation
hist		histogram of counts
histc		histogram of counts with user specified edges
quantile		quantiles of a dataset
princomp		principal component analysis (PCA)
normpdf		probability density function of normal distribution
log_normpdf		logarithm version of probability density function of normal distribution
normcdf		cumulative distribution function of normal distribution
mvnrnd		random vectors from multivariate normal distribution
chi2rnd		random numbers from chi-squared distribution
wishrnd		random matrix from Wishart distribution
iwishrnd		random matrix from inverse Wishart distribution
running_stat		running statistics of scalars (one dimensional process/signal)
running_stat_vec		running statistics of vectors (multi-dimensional process/signal)
kmeans		cluster data into disjoint sets

Miscellaneous

constants		pi, inf, NaN, speed of light, ...
wall_clock		timer for measuring number of elapsed seconds
libraries		list of external libraries used by PyArmadillo
Matlab/PyArmadillo syntax differences		examples of Matlab syntax and conceptually corresponding PyArmadillo syntax
example program		short example program
API additions		API stability and list of API additions

Matrix and Cube Classes

mat, fmat, cx_mat, cx_fmat, umat, imat

Classes for dense matrices, with elements stored in column-major ordering (i.e. column by column)

The name of each matrix class is based on the element type contained:

`mat`	=	`matrix containing double-precision floating point elements`
`fmat`	=	`matrix containing single-precision floating point elements`
`cx_mat`	=	`matrix containing complex double-precision floating point elements`
`cx_fmat`	=	`matrix containing complex single-precision floating point elements`
`umat`	=	`matrix containing unsigned integers`
`imat`	=	`matrix containing signed integers`

In this documentation the mat class is used for convenience; it is possible to use other classes instead, e.g. cx_mat

Functions which use LAPACK (generally matrix decompositions) are only valid for the following classes: mat, fmat, cx_mat, cx_fmat

Constructors:

`mat()`
`mat(n_rows, n_cols)`		(elements are initialised to zero)
`mat(n_rows, n_cols, fill_type)`		(elements are initialised)
`mat(list)`		(treated as a column vector)
`mat(list of lists)`
`mat(str)`
`cx_mat(mat,mat)`		(for constructing a complex matrix out of two real matrices)
`mat(NumPy array)`		(converts a NumPy array of the same type into a PyArmadillo matrix)
`mat(submatrix/diagonal view)`		(explicitly copies a submatrix/diagonal view into a matrix)
`mat(matrix)`		(converts matrices of one type to another)

When using the mat(n_rows, n_cols) constructor, by default the elements are initialised to zero; the elements can be explicitly initialised by specifying the fill_type, which is one of: fill.zeros, fill.ones, fill.eye, fill.randu, fill.randn, fill.none, with the following meanings:

`fill.zeros`	=	set all elements to 0
`fill.ones`	=	set all elements to 1
`fill.eye`	=	set the elements along the main diagonal to 1 and off-diagonal elements to 0
`fill.randu`	=	set each element to a random value from a uniform distribution in the [0,1] interval
`fill.randn`	=	set each element to a random value from a normal/Gaussian distribution with zero mean and unit variance
`fill.none`	=	do not initialise the elements

For the mat(string) constructor, the format is elements separated by spaces, and rows denoted by semicolons; for example, the 2x2 identity matrix can be created using "1 0; 0 1".

Matrices with size Nx1 and 1xN are interpreted as column vectors and row vectors, respectively

PyArmadillo matrices can be converted into NumPy arrays using numpy.array(mat)

Plotting of PyArmadillo matrices is possible using matplotlib through conversion to a NumPy array

Examples:

A = mat(5, 5, fill.randu) # initialise 5x5 matrix with uniformly distributed random values
x = A[1,2]                # get element in row 1, column 2 (indices from 0)

M = umat(5, 5, fill.ones) # initialise unsigned integer matrix with ones
N = mat(M)                # convert to double-precision floating point

B = A + A                 # addition
C = A * B                 # matrix multiplication
D = A @ B                 # element-wise multiplication

X = cx_mat(A,B)           # initialise complex matrix with A as real, B as imaginary

B.zeros()                 # fill B with zeros
B.set_size(10,10)         # set size of B without preserving data
B.ones(5,6)               # set size of B to 5x5, fill with ones

B.print("B:")             # print B with "B:" as header

See also:

cube, fcube, cx_cube, cx_fcube, ucube, icube

Classes for dense cubes, also known as "3D matrices" or 3rd order tensors

The name of each class is named based on the element type contained:

`cube`	=	`cube containing double-precision floating point elements`
`fcube`	=	`cube containing single-precision floating point elements`
`cx_cube`	=	`cube containing complex double-precision floating point elements`
`cx_fcube`	=	`cube containing complex single-precision floating point elements`
`ucube`	=	`cube containing unsigned integers`
`icube`	=	`cube containing signed integers`

In this documentation the cube class is used for convenience; it is possible to use other classes instead, e.g. cx_cube

Constructors:

`cube()`
`cube(n_rows, n_cols, n_slices)`		(elements are initialised to zero)
`cube(n_rows, n_cols, n_slices, fill_type)`		(elements are initialised)
`cube(list of lists of lists)`
`cx_cube(cube,cube)`		(for constructing a complex cube out of two real cubes)
`cube(NumPy array)`		(converts a NumPy array of the same type into a PyArmadillo cube)
`cube(subcube view)`		(explicitly copies a subcube view into a cube)
`cube(cube)`		(converts cubes of one type to another)

When using the cube(n_rows, n_cols, n_slices) constructor, by default the elements are initialised to zero; the elements can be explicitly initialised by specifying the fill_type, which is one of: fill.zeros, fill.ones, fill.randu, fill.randn, fill.none, with the following meanings:

`fill.zeros`	=	set all elements to 0
`fill.ones`	=	set all elements to 1
`fill.randu`	=	set each element to a random value from a uniform distribution in the [0,1] interval
`fill.randn`	=	set each element to a random value from a normal/Gaussian distribution with zero mean and unit variance
`fill.none`	=	do not initialise the elements

PyArmadillo cubes can be converted into NumPy arrays using numpy.array(cube)

Examples:

A = cube(5, 5, 5, fill.randu) # initialise 5x5x5 cube with uniformly distributed random values
x = A[1,2,3]                  # get element in row 1, column 2, slice 3 (indices from 0)

Q = ucube(5, 5, 5, fill.ones) # initialise unsigned integer cube with ones
R = cube(Q)                   # convert to double-precision floating point

B = A + A                     # addition
C = 2 * B                     # scalar multiplication
D = A @ B                     # element-wise multiplication

X = cx_cube(A,B)              # initialise complex cube with A as real, B as imaginary

B.zeros()                     # fill B with zeros
B.set_size(10,10,10)          # set size of B without preserving data
B.ones(5,6,7)                 # set size of B to 5x6x7, fill with ones

B.print("B:")                 # print B with "B:" as header

See also:

operators: + − * @ / == != <= >= < > && ||

Overloaded operators for matrix and cube classes

Operations:

`+`		addition of two objects
`−`		subtraction of one object from another or negation of an object

`*`		matrix multiplication of two objects; not applicable to the cube class unless multiplying by a scalar

`@`		element-wise multiplication of two objects (Schur product)
`/`		element-wise division of an object by another object or a scalar

`==`		element-wise equality evaluation of two objects; generates a matrix of type umat
`!=`		element-wise non-equality evaluation of two objects; generates a matrix of type umat

`>=`		element-wise "greater than or equal to" evaluation of two objects; generates a matrix of type umat/ucube
`<=`		element-wise "less than or equal to" evaluation of two objects; generates a matrix of type umat/ucube

`>`		element-wise "greater than" evaluation of two objects; generates a matrix of type umat/ucube
`<`		element-wise "less than" evaluation of two objects; generates a matrix of type umat/ucube

`&&`		element-wise logical AND evaluation of two objects; generates a matrix of type umat/ucube
`\|\|`		element-wise logical OR evaluation of two objects; generates a matrix of type umat/ucube

For element-wise relational and logical operations (i.e. ==, !=, >=, <=, >, <, &&, ||) each element in the generated object is either 0 or 1, depending on the result of the operation

The generated object is of type umat for matrix inputs and ucube for cube inputs

Caveat: operators involving equality comparison (i.e. ==, !=, >=, <=) are not recommended for matrices of type mat, due to the necessarily limited precision of the underlying element types; consider using approx_equal() instead

Broadcasting operations (as in Matlab/Octave) occur automatically in PyArmadillo
- broadcasting in matrices: apply a vector operation to each row or column of a matrix
- broadcasting in cubes: apply a matrix operation to each slice of a cube

A RuntimeError is thrown if incompatible object sizes are used

Examples:

A = mat(5, 10, fill.randu)
B = mat(5, 10, fill.randu)
C = mat(10, 5, fill.randu)

P = A + B       # addition
Q = A - B       # subtraction
R = -B          # negation
S = A / 123.0   # element-wise division by a scalar
T = A @ B       # element-wise multiplication
U = A * C       # matrix multiplication

AA = umat("1 2 3; 4 5 6; 7 8 9;")
BB = umat("3 2 1; 6 5 4; 9 8 7;")

# compare elements
ZZ = (AA >= BB)

# broadcasting
X = mat(6, 5, fill.ones)
v = mat(6, 1, fill.randu)

# in-place addition of v to each column vector of X
X += v

# generate Y by adding v to each column vector of X
Y = X + v

See also:
- approx_equal()
- powmat()
- any()
- all()
- accu()
- as_scalar()
- find()
- .replace()
- miscellaneous element-wise functions (exp, log, pow, sqrt, square, round, ...)
- floating point representation in Wikipedia
- floating point representation in MathWorld

Member Functions & Variables

attributes

`.n_rows`		number of rows; present in mat and cube
`.n_cols`		number of columns; present in mat and cube
`.n_slices`		number of slices; present in cube
`.n_elem`		total number of elements; present in mat and cube

The variables are unsigned integers

The variables are read-only; to change the size, use .set_size(), .copy_size(), .zeros(), .ones(), or .reset()

Examples:

X = mat(4,5)
print("X has " + str(X.n_cols) + " columns")

Y = cube(4,5,6)
print("Y has " + str(Y.n_slices) + " slices")

See also:
- .set_size()
- .copy_size()
- .zeros()
- .ones()
- .reset()
- size()

element/object access using []

Provide access to individual elements or objects stored in matrices

`[n]`		For mat and cube classes, access the n-th element/object under the assumption of a flat layout, with column-major ordering of data (i.e. column by column). A RuntimeError is thrown if the requested element is out of bounds.

`[i,j]`		For mat class, access the element/object stored at the i-th row and j-th column. A RuntimeError is thrown if the requested element is out of bounds.

`[i,j,k]`		For cube class, access the element/object stored at the i-th row, j-th column and k-th slice. A RuntimeError is thrown if the requested element is out of bounds.

Examples:

A = mat(10, 10, fill.randu)
A[9,9] = 123.0  # set element in row 9, column 9 to 123.0
x = A[9,9]      # get element in row 9, column 9
y = A[99]       # get 99th element

A = cube(10, 10, 10, fill.randu)
A[9,9,9] = 123.0  # set element in row 9, column 9, slice 9 to 123.0
x = A[9,9,9]      # get element in row 9, column 9, slice 9
y = A[99]         # get 99th element

See also:

element initialisation

Set elements in matrices and cubes using lists

Examples:

A = mat([ [1, 3, 5],
          [2, 4, 6] ])

Q = cube([ [ [ 1, 3 , 5  ],
             [ 2, 4 , 6  ] ],
           [ [ 7, 9 , 11 ],
             [ 8, 10, 12 ] ] ])

See also:

.zeros()
.zeros( n_rows, n_cols )
.zeros( size(X) )

Set the elements of an object to zero, optionally first changing the size to specified dimensions

Examples:

A = mat(5,10)       
A.zeros()           
# or:  A = mat(5,10,fill.zeros)

B = mat()           
B.zeros(10,20)      

C = mat()
C.zeros( size(B) )

See also:
- zeros() (standalone function)
- .ones()
- .clean()
- .is_zero()
- .randu()
- .fill()
- .reset()
- .set_size()
- size()

.ones()
.ones( n_rows, n_cols )
.ones( size(X) )

Set all the elements of an object to one, optionally first changing the size to specified dimensions

Examples:

A = mat(5,10)
A.ones()
# or:  A = mat(5,10,fill.ones)

B = mat()
B.ones(10,20)

C = mat()
C.ones( size(B) )

See also:
- ones() (standalone function)
- .eye()
- .zeros()
- .fill()
- .imbue()
- .randu()
- size()

.eye()
.eye( n_rows, n_cols )
.eye( size(X) )

Member functions of matrix

Set the elements along the main diagonal to one and off-diagonal elements to zero, optionally first changing the size to specified dimensions

An identity matrix is generated when n_rows = n_cols

Examples:

A = mat(5,5)
A.eye()
# or:  A = mat(5,5,fill.eye)

B = mat()
B.eye(5,5)

C = mat()
C.eye( size(B) )

See also:
- eye() (standalone function)
- .ones()
- .diag()
- diagmat()
- diagvec()
- size()

.randu()
.randu( n_rows, n_cols )			(member function of matrix)
.randu( n_rows, n_cols, n_slice )			(member function of cube)
.randu( size(X) )

.randn()
.randn( n_rows, n_cols )			(member function of matrix)
.randn( size(X) )			(member function of cube)

Set all the elements to random values, optionally first changing the size to specified dimensions

.randu() uses a uniform distribution in the [0,1] interval

.randn() uses a normal/Gaussian distribution with zero mean and unit variance

To change the RNG seed, use pyarma_rng.set_seed(value) or pyarma_rng.set_seed_random() functions

Examples:

A = mat(5,10)
A.randu()
# or:  A = mat(5,10,fill.randu)

B = mat()
B.randu(10,20)

C = mat()
C.randu( size(B) )

pyarma_rng.set_seed_random()  # set the seed to a random value

See also:
- randu() & randn() (standalone functions)
- .fill()
- .ones()
- .zeros()
- size()
- uniform distribution in Wikipedia
- normal distribution in Wikipedia

.fill( value )

Sets the elements to a specified value

The type of value must match the type of elements used by the container object (e.g. for mat the type is double)

Examples:

Note: to explicitly set all elements to zero during object construction, use the following more compact form:

See also:

.imbue( functor )
.imbue( lambda_function )

Imbue (fill) with values provided by a functor or lambda function

For matrices, filling is done column-by-column (i.e. column 0 is filled, then column 1, ...)

Examples:

from random import uniform
  
A = mat(4, 5)
  
A.imbue( lambda: uniform(0, 1) )

See also:
- .fill()
- .transform()
- element access
- function object at Wikipedia
- lambda expressions in Python

.clean( threshold )

For objects with non-complex elements: each element with an absolute value ≤ threshold is replaced by zero

For objects with complex elements: for each element, each component (real and imaginary) with an absolute value ≤ threshold is replaced by zero

Can be used to sparsify a matrix, in the sense of zeroing values with small magnitudes

Examples:

A = mat()

A.randu(1000, 1000)

A[12,34] =  datum.eps
A[56,78] = -datum.eps

A.clean(datum.eps)

See also:
- .replace()
- .is_zero()
- .zeros()
- clamp()
- datum.eps

.replace( old_value, new_value )

For all elements equal to old_value, set them to new_value

The type of old_value and new_value must match the type of elements used by the container object (e.g. for mat the type is double)

Caveat: floats are approximations due to their necessarily limited precision

Examples:

A = mat(5, 6, fill.randu)

A[diag].fill(datum.nan)

A.replace(datum.nan, 0) # replace each NaN with 0

See also:

.transform( functor )
.transform( lambda_function )

Transform each element using a functor or lambda function

For dense matrices, transformation is done column-by-column for all elements

Examples:

A = mat(4, 5, fill.ones)

# add 123 to every element
A.transform( lambda val: val + 123.0 )

See also:
- .for_each()
- .replace()
- .imbue()
- .clean()
- .fill()
- element access
- overloaded operators
- miscellaneous element-wise functions (exp, log, pow, sqrt, square, round, ...)
- function object at Wikipedia
- lambda expressions in Python

.for_each( functor )
.for_each( lambda_function )

For each element, pass its value to a functor or lambda function

For dense matrices and fields, the processing is done column-by-column for all elements

Caveat: element values are constant

Examples:

# print each element in a dense matrix

A = mat(4, 5, fill.ones)

A.for_each( lambda val: print(val) )

See also:
- .transform()
- .replace()
- element access
- miscellaneous element-wise functions (exp, log, pow, sqrt, square, round, ...)
- function object at Wikipedia

.set_size( n_rows, n_cols )			(member function of matrix)
.set_size( n_rows, n_cols, n_slices )			(member function of cube)
.set_size( size(X) )

Change the size of an object, without explicitly preserving data and without initialising the elements

If you need to initialise the elements to zero while changing the size, use .zeros() instead

If you need to explicitly preserve data while changing the size, use .reshape() or .resize() instead;
caveat: .reshape() and .resize() are considerably slower than .set_size()

Examples:

A = mat()
A.set_size(5,10)       # or:  A = mat(5,10)

B = mat()
B.set_size( size(A) )  # or:  B = mat( size(A) )

See also:
- .reset()
- .copy_size()
- .reshape()
- .resize()
- .zeros()
- size()

.reshape( n_rows, n_cols )			(member function of matrix)
.reshape( n_rows, n_cols, n_slices )			(member function of cube)
.reshape( size(X) )

Recreate the object according to given size specifications, with the elements taken from the previous version of the object in a column-wise manner; the elements in the generated object are placed column-wise (i.e. the first column is filled up before filling the second column)

The layout of the elements in the recreated object will be different to the layout in the previous version of the object

If the total number of elements in the previous version of the object is less than the specified size, the extra elements in the recreated object are set to zero

If the total number of elements in the previous version of the object is greater than the specified size, only a subset of the elements is taken

Caveats:
- do not use .reshape() if you simply want to change the size without preserving data; use .set_size() instead, which is much faster
- to grow/shrink the object while preserving the elements as well as the layout of the elements, use .resize() instead
- to flatten a matrix into a vector, use vectorise() or .as_col() / .as_row() instead

Examples:

A = mat(4, 5, fill.randu)

A.reshape(5,4)

See also:
- .resize()
- .set_size()
- .copy_size()
- .zeros()
- .reset()
- .as_col() / .as_row()
- reshape() (standalone function)
- vectorise()
- size()

.resize( n_rows, n_cols )			(member function of matrix)
.resize( size(X) )			(member function of cube)

Recreate the object according to given size specifications, while preserving the elements as well as the layout of the elements

Can be used for growing or shrinking an object (i.e. adding/removing rows, and/or columns)

Caveat: do not use .resize() if you simply want to change the size without preserving data; use .set_size() instead, which is much faster

Examples:

A = mat(4, 5, fill.randu)

A.resize(7,6)

See also:
- .reshape()
- .set_size()
- .copy_size()
- .zeros()
- .reset()
- .insert_rows/cols/slices
- .shed_rows/cols/slices
- resize() (standalone function)
- vectorise()
- size()

.copy_size( A )

Set the size to be the same as object A

Examples:

A = mat(5, 6, fill.randu)
B = mat()
B.copy_size(A)

print(str(B.n_rows))
print(str(B.n_cols))

See also:
- .reset()
- .set_size()
- .reshape()
- .resize()
- .zeros()
- size()

.reset()

Reset the size to zero (the object will have no elements)

Examples:

See also:

submatrix views

Read/write access to submatrix views

Contiguous views for matrix X:

to access a submatrix:

[

]

[

, size(

) ]

[

, size(

) ]

[ Y is a matrix ]

to access a specific column:

[

]

to access a specific row:

[

]

to access a contiguous span of columns:

[

]

to access a contiguous span of rows:

[

]

to access a subset of columns/rows from the first column/row:

[

head_cols

]

[

head_rows

]

to access a subset of columns/rows until the last column/row:

[

tail_cols

]

[

tail_rows

]

Non-contiguous views for matrix X:

to access elements:

[

]

to access rows/columns:

[

: ,

]

[

, :

]

to access submatrices:

[

]

Accessing diagonals:

[diag]

Instances of first : last can be replaced by : to indicate the entire range

For functions requiring one or more vector of indices, e.g. X[vector_of_row_indices, vector_of_column_indices], each vector of indices must be of type umat

In the function X[vector_of_indices], elements specified in vector_of_indices are accessed. X is interpreted as one long vector, with column-by-column ordering of the elements of X. The vector_of_indices must evaluate to a vector of type umat (e.g. generated by the find() function). The aggregate set of the specified elements is treated as a column vector (i.e. the output of X[vector_of_indices] is always a column vector).

Examples:

A = mat(5, 10, fill.zeros)

# set submatrix of rows 0-2 and columns 1-3 to uniformly distributed random elements
A[ 0:2, 1:3       ] = mat(3,3, fill.randu) 
A[ 0,1, size(3,3) ] = mat(3,3, fill.randu)

# copy submatrix of rows 0-2 and columns 1-3
B = A[ 0:2, 1:3        ].eval()                         
D = A[ 0, 1, size(3,3) ].eval()

# set column 1 to uniformly distributed random elements
A[:,1] = mat(5,1,fill.randu)            

X = mat(5, 5, fill.randu)

# get all elements of X that are greater than 0.5
q = X[ find(X > 0.5) ]

# add 123 to all elements of X greater than 0.5
X[ find(X > 0.5) ] += 123.0

# set four specific elements of X to 1
indices = umat([ 2, 3, 6, 8 ])

X[indices] = mat(4,1,fill.ones)

See also:

subcube views and slices

A collection of member functions of the cube class that provide subcube views

Contiguous views for cube Q:

to access a subcube:

[

]

[

, size(

) ]

[

, size(

) ]

[ R is a cube ]

to access a slice:

[

]

to access a single slice, with the slice provided as a matrix:

[

single_slice

]

to access a column:

[

]

to access a row:

[

]

to access columns:

[

]

to access rows:

[

]

to access slices:

[

]

to access a subset of slices from the first slice:

[

head_slices

]

to access a subset of slices until the last slice:

[

tail_slices

]

to access tubes:

[

]

[

]

[

]

[

, size(

) ]

Non-contiguous views for cube Q:

to access elements:

[

]

to access slices:

[

: , : ,

]

Instances of first : last can be replaced by:
- : to indicate the entire range
- a : a to indicate a particular row, column or slice

An individual slice, accessed via Q[single_slice, slice_number], is an instance of the mat class (a reference to a matrix is provided)

All tube forms are variants of subcube forms, using first_slice = 0 and last_slice = Q.n_slices-1

The tube [row,col] form uses row = first_row = last_row, and col = first_col = last_col

In the function Q[vector_of_indices], elements specified in vector_of_indices are accessed. Q is interpreted as one long vector, with column-by-column ordering of the elements of X. The vector_of_indices must evaluate to a vector of type umat (e.g. generated by the find() function). The aggregate set of the specified elements is treated as a column vector (i.e. the output of Q[vector_of_indices] is always a column vector).

Examples:

A = cube(2, 3, 4, fill.randu)

B = A[single_slice, 1]  # the slice is accessed as a matrix

A[:,:,0] = mat(2,3,fill.randu)
A[:,:,0][1,2] = 99.0

A[ 0:1, 0:1, 1:2        ] = cube(2,2,2,fill.randu)
A[ 0, 0, 1, size(2,2,2) ] = cube(2,2,2,fill.randu)

# add 123 to all elements of A greater than 0.5
A[ find(A > 0.5) ] += 123.0

C = A[head_slices,2]  # get first two slices

A[head_slices, 2] += 123.0

See also:

X[diag]
X[diag, k]

Read/write access to a diagonal in a matrix

The argument k is optional; by default the main diagonal is accessed (k = 0)

For k > 0, the k-th super-diagonal is accessed (top-right corner)

For k < 0, the k-th sub-diagonal is accessed (bottom-left corner)

The diagonal is interpreted as a column vector within expressions

Examples:

X = mat(5, 5, fill.randu)

a = X[diag].eval()            # copy main diagonal
b = X[diag, 1].eval()         # copy first super-diagonal
c = X[diag, -2].eval()        # copy second sub-diagonal

X[diag] = mat(5,1,fill.randu) # set main diagonal to uniformly distributed random elements
X[diag] += 6                  # add 6 to all main diagonal elements
X[diag].ones()                # set all main diagonal elements to one

See also:

.set_imag( X )
.set_real( X )

Set the imaginary/real part of an object

X must have the same size as the recipient object

Examples:

A = mat(4, 5, fill.randu)
B = mat(4, 5, fill.randu)

C = cx_mat(4, 5, fill.zeros)

C.set_real(A)
C.set_imag(B)

Caveat: to directly construct a complex matrix out of two real matrices, the following code is faster:

See also:
- matrix constructors
- imag() / real()

.insert_rows( row_number, X ) .insert_rows( row_number, number_of_rows ) .insert_rows( row_number, number_of_rows, set_to_zero )

.insert_cols( col_number, X ) .insert_cols( col_number, number_of_cols ) .insert_cols( col_number, number_of_cols, set_to_zero )

.insert_slices( slice_number, X ) .insert_slices( slice_number, number_of_slices ) .insert_slices( slice_number, number_of_slices, set_to_zero )		(member functions of cube)

Functions with the X argument: insert a copy of X at the specified row/column
- if inserting rows, X must have the same number of columns as the recipient object
- if inserting columns, X must have the same number of rows as the recipient object
- if inserting slices, X must have the same number of rows and columns as the recipient object (i.e. all slices must have the same size)

Functions with the number_of_... argument:
- expand the object by creating new rows/columns
- by default, the new rows/columns are set to zero
- if set_to_zero is False, the memory used by the new rows/columns will not be initialised

Examples:

A = mat(5, 10, fill.randu)
B = mat(5,  2, fill.ones )

# at column 2, insert a copy of B;
# A will now have 12 columns
A.insert_cols(2, B)

# at column 1, insert 5 zeroed columns;
# B will now have 7 columns
B.insert_cols(1, 5)

See also:

.shed_row( row_number ) .shed_rows( first_row, last_row ) .shed_rows( vector_of_indices )

.shed_col( column_number ) .shed_cols( first_column, last_column ) .shed_cols( vector_of_indices )

.shed_slices( slice_number ) .shed_slices( first_slice, last_slice ) .shed_slices( vector_of_indices )		(member functions of cube)

Functions with single scalar argument: remove the specified row/column

Functions with two scalar arguments: remove the specified range of rows/columns

The vector_of_indices must evaluate to a vector of type umat; it contains the indices of rows/columns/slices to remove

Examples:

A = mat(5, 10, fill.randu)
B = mat(5, 10, fill.randu)

A.shed_row(2)
A.shed_cols(2,4)

indices = umat([4, 6, 8])
B.shed_cols(indices)

See also:

.swap_rows( row1, row2 )
.swap_cols( col1, col2 )

Swap the contents of specified rows or columns

Examples:

X = mat(5, 5, fill.randu)
X.swap_rows(0,4)

See also:
- fliplr() & flipud()
- .swap()

.swap( X )

Swap contents with object X

Examples:

A = mat(4, 5, fill.zeros)
B = mat(6, 7, fill.ones )

A.swap(B)

See also:
- .swap_rows() & .swap_cols()

iterators (dense matrices)

Iterators traverse over all elements within the specified range

All iterators are read-only

All ranges are inclusive

Functions:

iter( X )		iterator from the first to last element
iterator( X, i, j )		iterator from the i-th element to the j-th element

col_iter( X, i, j )		iterator from the first element of the i-th column to the last element of the j-th column

row_iter( X, i, j )		iterator from the first element of the i-th row to the last element of the j-th row

i and j are optional; by default the first and last indices are used

Examples:

X = mat(5, 6, fill.randu)

for i in X:                # this prints all elements
  print(i)

it = iter(X)               # this also prints all elements
for i in it:
  print(i)

col_it = col_iter(X, 1, 3) # start of column 1 and end of column 3
for i in col_it:           # print all elements of columns 1-3
  print(i)

See also:

iterators (cubes)

Iterators traverse over all elements within the specified range

All iterators are read-only

All ranges are inclusive

Functions:

iter( Q )		iterator from the first to last element
iterator( Q, i, j )		iterator from the i-th element to the j-th element

slice_iter( Q, i, j )		iterator from the first element of the i-th slice to the last element of the j-th slice

i and j are optional; by default the first and last indices are used

Examples:

X = cube(5, 6, 7, fill.randu)

for i in X:                # this prints all elements
  print(i)

it = iter(X)               # this also prints all elements
for i in it:
  print(i)

slice_it = slice_iter(X, 1, 3) # start of slice 1 and end of cslice 3
for i in slice_it:             # print all elements of slices 1-3
  print(i)

See also:

iterators (dense submatrices & subcubes)

iterators for dense submatrix and subcube views

All iterators are read-only

Examples:

X = mat(100, 200, fill.randu)

for val in X[40:60, 50:100]:
  print(val)

See also:

.as_col()
.as_row()

.as_col(): return a flattened version of the matrix as a column vector; flattening is done by concatenating all columns

.as_row(): return a flattened version of the matrix as a row vector; flattening is done by concatenating all rows

Caveat: concatenating columns is faster than concatenating rows

Examples:

X = mat(4, 5, fill.randu)
v = X.as_col()

See also:

.t()
.st()

For real (non-complex) matrix:
- .t() provides a transposed copy of the matrix
- .st() not applicable

For complex matrix:
- .t() provides a Hermitian transpose (i.e. the conjugate of the elements is taken during the transpose)
- .st() provides a transposed copy without taking the conjugate of the elements

Examples:

See also:

.i()

Provides an inverse of the matrix expression

If the matrix expression is not square sized, a RuntimeError is thrown

If the matrix expression appears to be singular, the output matrix is reset and a RuntimeError is thrown

Caveats:
- if matrix A is know to be symmetric positive definite, it is faster to use inv_sympd() instead
- to solve a system of linear equations, such as Z = inv(X)*Y, using solve() can be faster and/or more accurate

Examples:

A = mat(4, 4, fill.randu)

X = A.i()

Y = (A+A).i()

See also:
- inv()
- rcond()
- pinv()
- solve()
- spsolve()

.min()
.max()

Return the extremum value of any matrix and cube expression

For objects with complex numbers, absolute values are used for comparison

Examples:

A = mat(5, 5, fill.randu)

max_val = A.max()

See also:
- .index_min() & .index_max()
- min() & max() (standalone functions with extended functionality)
- clamp()

.index_min()
.index_max()

Return the linear index of the extremum value of any matrix expression

For objects with complex numbers, absolute values are used for comparison

The returned index is an unsigned integer

Examples:

A = mat(5, 5, fill.randu)

i = A.index_max()

max_val = A[i]

See also:
- .min() & .max()
- index_min() & index_max() (standalone functions with extended functionality)
- ind2sub()
- sort_index()
- find()
- element access

.eval()

Explicitly copies a submatrix/diagonal view into a matrix

Examples:

A = umat(4,4,fill.ones)

B = A[0:1, 2:3].eval()

# any change to B does not affect A
B[0,0] = 5

# therefore, B and A[0:1, 2:3] are not equal
(B == A[0:1, 2:3]).print()

See also:

.in_range( i )
.in_range( start : end )

.in_range( row, col )			(member of mat)
.in_range( start_row : end_row, start_col : end_col )			(member of mat)

.in_range( row, col , slice)			(member of cube)
.in_range( start_row : end_row, start_col : end_col, start_slice: end_slice )			(member of cube)

.in_range( first_row, first_col, size(X) ) (X is a matrix)			(member of mat)
.in_range( first_row, first_col, size(n_rows, n_cols) )			(member of mat)

.in_range( first_row, first_col, first_slice, size(Q) ) (Q is a cube)			(member of cube)
.in_range( first_row, first_col, first_slice, size(n_rows, n_cols, n_slices) )			(member of cube)

Returns True if the given location is currently valid

Returns False if the object is empty, the location is out of bounds, or the span is out of bounds

Instances of start : end can be replaced by:
- : to indicate the entire range
- a : a to indicate a particular row, column or slice

Examples:

A = mat(4, 5, fill.randu)

print(A.in_range(0,0))  # True
print(A.in_range(3,4))  # True
print(A.in_range(4,5))  # False

See also:

.is_empty()

Returns True if the object has no elements

Returns False if the object has one or more elements

Examples:

A = mat(5, 5, fill.randu)
print(A.is_empty())

A.reset()
print(A.is_empty())

See also:

.is_vec()
.is_colvec()
.is_rowvec()

.is_vec():
- returns True if the matrix can be interpreted as a vector (either column or row vector)
- returns False if the matrix does not have exactly one column or one row

.is_colvec():
- returns True if the matrix can be interpreted as a column vector
- returns False if the matrix does not have exactly one column

.is_rowvec():
- returns True if the matrix can be interpreted as a row vector
- returns False if the matrix does not have exactly one row

Caveat: do not assume that the vector has elements if these functions return True; it is possible to have an empty vector (e.g. 0x1)

Examples:

A = mat(1, 5, fill.randu)
B = mat(5, 1, fill.randu)
C = mat(5, 5, fill.randu)

print(A.is_vec())
print(B.is_vec())
print(C.is_vec())

See also:

.is_sorted()		(form 1)
.is_sorted( sort_direction )		(form 2)
.is_sorted( sort_direction, dim )		(form 3)

For forms 1 and 2:
- if the matrix can be interpreted as a vector, return a bool indicating whether the elements are sorted
- otherwise, form 3 is used with dim = 0

For form 3, return a bool indicating whether the elements are sorted in each column (dim = 0), or each row (dim = 1)

The sort_direction argument is optional; sort_direction is one of:

`"ascend"`	↦	elements are ascending; consecutive elements can be equal; this is the default operation
`"descend"`	↦	elements are descending; consecutive elements can be equal
`"strictascend"`	↦	elements are strictly ascending; consecutive elements cannot be equal
`"strictdescend"`	↦	elements are strictly descending; consecutive elements cannot be equal

For matrices with complex numbers, order is checked via absolute values

Examples:

a = mat(10, 1, fill.randu)
b = sort(a)

check1 = a.is_sorted()
check2 = b.is_sorted()

A = mat(10, 10, fill.randu)

# check whether each column is sorted in descending manner
print(A.is_sorted("descend"))

# check whether each row is sorted in ascending manner
print(A.is_sorted("ascend", 1))

See also:
- sort()
- sort_index()

.is_trimatu()
.is_trimatl()

.is_trimatu():
- return True if the matrix is upper triangular, i.e. the matrix is square sized and all elements below the main diagonal are zero; return False otherwise
- caveat: if this function returns True, do not assume that the matrix contains non-zero elements on or above the main diagonal

.is_trimatl():
- return True if the matrix is lower triangular, i.e. the matrix is square sized and all elements above the main diagonal are zero; return False otherwise
- caveat: if this function returns True, do not assume that the matrix contains non-zero elements on or below the main diagonal

Examples:

A = mat(5, 5, fill.randu)
B = trimatl(A)

print(A.is_trimatu())
print(B.is_trimatl())

See also:

.is_diagmat()

Return True if the matrix is diagonal, i.e. all elements outside of the main diagonal are zero

Return False otherwise

Caveat: if this function returns True, do not assume that the matrix contains non-zero elements on the main diagonal

Examples:

A = mat(5, 5, fill.randu)
B = diagmat(A)

print(A.is_diagmat())
print(B.is_diagmat())

See also:

.is_square()

Returns True if the matrix is square, i.e. number of rows is equal to the number of columns

Returns False if the matrix is not square

Examples:

A = mat(5, 5, fill.randu)
B = mat(6, 7, fill.randu)
 
print(A.is_square())
print(B.is_square())

See also:

.is_symmetric()
.is_symmetric( tol )

Returns True if the matrix is symmetric

Returns False if the matrix is not symmetric

The tol argument is optional; if tol is specified, the given matrix X is considered symmetric if norm(X - X.st(), "inf") / norm (X, "inf") ≤ tol

Examples:

A = mat(5, 5, fill.randu)
B = A.t() * A

print(A.is_symmetric())
print(B.is_symmetric())

See also:

.is_hermitian()
.is_hermitian( tol )

Returns True if the matrix is hermitian (self-adjoint)

Returns False if the matrix is not hermitian

The tol argument is optional; if tol is specified, the given matrix X is considered hermitian if norm(X - X.t(), "inf") / norm (X, "inf") ≤ tol

Examples:

A = cx_mat(5, 5, fill.randu)
B = A.t() * A

print(A.is_hermitian())
print(B.is_hermitian())

See also:

.is_sympd()
.is_sympd( tol )

Returns True if the matrix is symmetric/hermitian positive definite within the tolerance given by tol

Returns False otherwise

The tol argument is optional; if tol is not specified, by default tol = 100 * datum.eps * norm(X, "fro")

Examples:

A = mat(5, 5, fill.randu)

B = A.t() * A

print(A.is_sympd())
print(B.is_sympd())

See also:

.is_zero()
.is_zero( tolerance )

For objects with non-complex elements: return True if each element has an absolute value ≤ tolerance; return False otherwise

For objects with complex elements: return True if for each element, each component (real and imaginary) has an absolute value ≤ tolerance; return False otherwise

The argument tolerance is optional; by default tolerance = 0

Examples:

A = mat(5, 5, fill.zeros)

A[0,0] = datum.eps

print(A.is_zero())
print(A.is_zero(datum.eps))

See also:

.is_finite()

Returns True if all elements of the object are finite

Returns False if at least one of the elements of the object is non-finite (±infinity or NaN)

Examples:

A = mat(5, 5, fill.randu)
B = mat(5, 5, fill.randu)

B[0,0] = datum.inf

print(A.is_finite())
print(B.is_finite())

See also:

.has_inf()

Returns True if at least one of the elements of the object is ±infinity

Returns False otherwise

Examples:

A = mat(5, 5, fill.randu)
B = mat(5, 5, fill.randu)

B[0,0] = datum.inf

print(A.has_inf())
print(B.has_inf())

See also:

.has_nan()

Returns True if at least one of the elements of the object is NaN (not-a-number)

Returns False otherwise

Caveat: NaN is not equal to anything, even itself

Examples:

A = mat(5, 5, fill.randu)
B = mat(5, 5, fill.randu)

B[0,0] = datum.nan

print(A.has_nan())
print(B.has_nan())

See also:

.print()
.print( header )

Print the contents of an object with an optional header string

Examples:

A = mat(5, 5, fill.randu)

A.print()

# "A:" is the optional header line
A.print("A:")

# print a transposed version of A
A.t().print()

See also:

.brief_print()
.brief_print( header )

Similar to the .print() member function, with the difference that an abridged version of the object is printed

Examples:

A = mat(123, 456, fill.randu)

A.brief_print()

# possible output:
# 
# [matrix size: 123x456]
#    0.8402   0.7605   0.6218   ...   0.9744
#    0.3944   0.9848   0.0409   ...   0.7799
#    0.7831   0.9350   0.4140   ...   0.8835
#         :        :        :     :        :        
#    0.4954   0.1826   0.9848   ...   0.1918

See also:
- .print()

saving/loading matrices and cubes

.save( filename )
.save( filename, file_type )

.load( filename )
.load( filename, file_type )

Store/retrieve data in a file

On success, .save() and .load() return True

On failure, .save() and .load() return False; additionally, .load() resets the object so that it has no elements

file_type can be one of the following:

auto_detect		Used only by .load() only: attempt to automatically detect the file type as one of the formats described below; [ default operation for .load() ]
arma_binary		Numerical data stored in machine dependent binary format, with a simple header to speed up loading. The header indicates the type and size of matrix. [ default operation for .save() ]
arma_ascii		Numerical data stored in human readable text format, with a simple header to speed up loading. The header indicates the type and size of matrix.
raw_binary		Numerical data stored in machine dependent raw binary format, without a header. Matrices are loaded to have one column. The .reshape() function can be used to alter the size of the loaded matrix without losing data.
raw_ascii		Numerical data stored in raw ASCII format, without a header. The numbers are separated by whitespace. The number of columns must be the same in each row. Data which was saved in Matlab/Octave using the -ascii option can be read in PyArmadillo, except for complex numbers. Complex numbers are stored in standard C++ notation, which is a tuple surrounded by brackets: e.g. (1.23,4.56) indicates 1.24 + 4.56i.
csv_ascii		Numerical data stored in comma separated value (CSV) text format, without a header. Handles complex numbers stored in the compound form of 1.24+4.56i
pgm_binary		Image data stored in Portable Gray Map (PGM) format. Saving int, float or double matrices is a lossy operation, as each element is copied and converted to an 8 bit representation. As such the matrix should have values in the [0,255] interval, otherwise the resulting image may not display correctly.
hdf5_binary		Numerical data stored in portable HDF5 binary format. for saving, the default dataset name within the HDF5 file is "dataset" for loading, the order of operations is: (1) try loading a dataset named "dataset", (2) try loading a dataset named "value", (3) try loading the first available dataset

Caveat: for saving/loading HDF5 files, the hdf5.h header file and a suitable HDF5 library must be available at the time of PyArmadillo installation

Examples:

A = mat(5, 6, fill.randu)

# default save format is arma_binary
A.save("A.bin")

# save in raw_ascii format
A.save("A.txt", raw_ascii)

# save in CSV format without a header
A.save("A.csv", csv_ascii)

# save in HDF5 format
A.save("A.h5", hdf5_binary)

# automatically detect format type while loading
B = mat()
B.load("A.bin")

# force loading in arma_ascii format
C = mat()
C.load("A.txt", arma_ascii)

# example of testing for success
D = mat()
ok = D.load("A.bin")

if not ok:
  print("problem with loading")

See also:

HDF in Wikipedia
CSV in Wikipedia

Generated Matrices

linspace( start, end )
linspace( start, end, N )

Generate an Nx1 matrix (of type mat) the values of the elements are linearly spaced from start to (and including) end

The argument N is optional; by default N = 100

Caveat: for N = 1, the generated matrix will have a single element equal to end

Example:

See also:

logspace( A, B )
logspace( A, B, N )

Generate an Nx1 matrix (of type mat) the values of the elements are logarithmically spaced from 10^A to (and including) 10^B

The argument N is optional; by default N = 50

Example:

See also:
- linspace()
- regspace()

regspace( start, end )
regspace( start, delta, end )

Generate an Nx1 matrix (of type mat) with regularly spaced elements:
[ (start + 0*delta), (start + 1*delta), (start + 2*delta), ⋯, (start + M*delta) ]
where M = floor((end-start)/delta), so that (start + M*delta) ≤ end

Similar in operation to the Matlab/Octave colon operator, i.e. start:end and start:delta:end

If delta is not specified:
- delta = +1, if start ≤ end
- delta = −1, if start > end (caveat: this is different to Matlab/Octave)

An empty matrix is generated when one of the following conditions is true:
- start < end, and delta < 0
- start > end, and delta > 0
- delta = 0

Examples:

a = regspace(0,  9)       # 0,  1, ...,   9

b = regspace(2,  2,  10)  # 2,  4, ...,  10

c = regspace(0, -1, -10)  # 0, -1, ..., -10

Caveat: do not use regspace() to specify ranges for contiguous submatrix views; use first : last instead

See also:
- linspace()
- logspace()

randperm( N )
randperm( N, M )

Generate an Nx1 matrix (of type umat) with a random permutation of integers from 0 to N-1

The optional argument M indicates the number of elements to return, sampled without replacement from 0 to N-1

Examples:

See also:

eye( n_rows, n_cols )
eye( size(X) )

Generate a matrix (of type mat) with the elements along the main diagonal set to one and off-diagonal elements set to zero

An identity matrix is generated when n_rows = n_cols

Examples:

A = eye(5,5)  // or:  A = mat(5,5,fill.eye)

B = 123.0 * eye(5,5)

C = 123.0 * fmat(eye(5,5))

See also:
- .eye() (member function of mat)
- .diag()
- ones()
- diagmat()
- diagvec()
- size()

ones( n_elem )
ones( n_rows, n_cols )
ones( n_rows, n_cols, n_slices )
ones( size(X) )

Generate a vector, matrix (of type mat) or cube (of type cube) with all elements set to one

Examples:

v = ones(10)
u = umat(ones(11))
A = ones(5,6)
Q = ones(5,6,7)

B = 123.0 * ones(5,6)

See also:
- .ones() (member function of mat and cube)
- eye()
- linspace()
- regspace()
- randperm()
- zeros()
- randu() & randn()
- size()

zeros( n_elem )
zeros( n_rows, n_cols )
zeros( n_rows, n_cols, n_slices )
zeros( size(X) )

Generate a vector, matrix (of type mat) or cube (of type cube) with the elements set to zero

Examples:

v = zeros(10)
u = umat(zeros(11))
A = zeros(5,6)
Q = zeros(5,6,7)

See also:
- .zeros() (member function of mat and cube)
- .ones() (member function of mat and cube)
- ones()
- randu() & randn()
- size()

randu( )
randu( n_elem )
randu( n_rows, n_cols )
randu( n_rows, n_cols, n_slices )
randu( size(X) )

randn( )
randn( n_elem )
randn( n_rows, n_cols )
randn( n_rows, n_cols, n_slices )
randn( size(X) )

Generate a scalar, vector, matrix or cube (of types double, mat, or cube) with the elements set to random floating point values

randu() uses a uniform distribution in the [0,1] interval

randn() uses a normal/Gaussian distribution with zero mean and unit variance

To change the RNG seed, use pyarma_rng.set_seed(value) or pyarma_rng.set_seed_random() functions

Caveat: to generate a matrix with random integer values instead of floating point values, use randi() instead

Examples:

v = randu(5)
A = randu(5,6)
Q = randu(5,6,7)

pyarma_rng.set_seed_random()  # set the seed to a random value

See also:
- .randu() & .randn() (member functions)
- randi()
- randg()
- mvnrnd()
- normpdf()
- .imbue()
- ones()
- zeros()
- shuffle()
- size()
- uniform distribution in Wikipedia
- normal distribution in Wikipedia

randg( )
randg( distr_param(a,b) )

randg( n_elem )
randg( n_elem, distr_param(a,b) )

randg( n_rows, n_cols )
randg( n_rows, n_cols, distr_param(a,b) )

randg( n_rows, n_cols, n_slices )
randg( n_rows, n_cols, n_slices, distr_param(a,b) )

randg( size(X) )
randg( size(X), distr_param(a,b) )

Generate a scalar, vector, matrix or cube (of types double, mat, or cube) with the elements set to random values from a gamma distribution:
where a is the shape parameter and b is the scale parameter, with constraints a > 0 and b > 0

The default distribution parameters are a=1 and b=1

Usage:
- s = randg( )
- s = randg( distr_param(a,b) )
- v = randg( n_elem )
- v = randg( n_elem, distr_param(a,b) )
- X = randg( n_rows, n_cols )
- X = randg( n_rows, n_cols, distr_param(a,b) )
- Y = randg( size(X) )
- Y = randg( size(X), distr_param(a,b) )
- Q = randg( n_rows, n_cols, n_slices )
- Q = randg( n_rows, n_cols, n_slices, distr_param(a,b) )
- R = randg( size(Q) )
- R = randg( size(Q), distr_param(a,b) )

To change the RNG seed, use pyarma_rng.set_seed(value) or pyarma_rng.set_seed_random() functions

Examples:

v = randg(100, distr_param(2,1))    # generates a vector

X = randg(10, 10, distr_param(2,1)) # generates a matrix

See also:

randi( )
randi( distr_param(a,b) )

randi( n_elem )
randi( n_elem, distr_param(a,b) )

randi( n_rows, n_cols )
randi( n_rows, n_cols, distr_param(a,b) )

randi( n_rows, n_cols, n_slices )
randi( n_rows, n_cols, n_slices, distr_param(a,b) )

randi( size(X) )
randi( size(X), distr_param(a,b) )

Generate a scalar, vector, matrix or cube (of types int, imat, or icube) with the elements set to random integer values in the [a,b] interval

The default distribution parameters are a=0 and b=maximum_int

Usage:
- v = randi( )
- v = randi( distr_param(a,b) )
- v = randi( n_elem )
- v = randi( n_elem, distr_param(a,b) )
- X = randi( n_rows, n_cols )
- X = randi( n_rows, n_cols, distr_param(a,b) )
- Y = randi( size(X) )
- Y = randi( size(X), distr_param(a,b) )
- Q = randi( n_rows, n_cols, n_slices )
- Q = randi( n_rows, n_cols, n_slices, distr_param(a,b) )
- R = randi( size(Q) )
- R = randi( size(Q), distr_param(a,b) )

To change the RNG seed, use pyarma_rng.set_seed(value) or pyarma_rng.set_seed_random() functions

Caveat: to generate a continuous distribution with floating point values (ie. float or double), use randu() instead

Examples:

A = randi(5, 6)                         # generates a matrix

A = randi(6, 7, distr_param(-10, +20))  # generates a matrix

pyarma_rng.set_seed_random()            # set the seed to a random value

See also:
- randu() & randn()
- randg()
- randperm()
- .imbue()
- ones()
- zeros()
- shuffle()
- size()

toeplitz( A )
toeplitz( A, B )
circ_toeplitz( A )

toeplitz(): generate a Toeplitz matrix, with the first column specified by A, and (optionally) the first row specified by B

circ_toeplitz(): generate a circulant Toeplitz matrix

A and B must be vectors

Examples:

A = mat(5, 1, fill.randu)
X = toeplitz(A)
Y = circ_toeplitz(A)

See also:

Functions of Matrices and Cubes

abs( X )

Obtain the magnitude of each element

Examples:

A = mat(5, 5, fill.randu)
B = abs(A) 

X = cx_mat(5, 5, fill.randu)
Y = abs(X)

See also:

accu( X )

Accumulate (sum) all elements of a matrix or cube

Examples:

A = mat(5, 6, fill.randu)
B = mat(5, 6, fill.randu)

x = accu(A)

y = accu(A @ B)

# accu(A @ B) is a "multiply-and-accumulate" operation
# as operator @ performs element-wise multiplication

See also:
- sum()
- cumsum()
- trace()
- mean()
- dot()
- as_scalar()

affmul( A, B )

Multiply matrix A by an augmented form of B, where a row with ones is appended to B;
for example:

⎡ A₀₀ A₀₁ A₀₂ ⎤   ⎡ B₀ ⎤
⎢ A₁₀ A₁₁ A₁₂ ⎥ x ⎢ B₁ ⎥
⎣ A₂₀ A₂₁ A₂₂ ⎦   ⎣ 1 ⎦

A is typically an affine transformation matrix

The number of columns in A must be equal to number of rows in the augmented form of B (i.e. A.n_cols = B.n_rows+1)

B can be a vector or matrix

Examples:

A = mat(4, 4, fill.randu)
B = mat(3, 1, fill.randu)

C = affmul(A,B)

See also:

all( X )		(form 1)
all( X, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, return a 1x1 matrix (of type umat) with its element (0 or 1) indicating whether all elements are non-zero or satisfy a relational condition
- otherwise, form 2 is used with dim = 0

For form 2,
- dim = 0 returns a row vector (of type umat), with each element (0 or 1) indicating whether the corresponding column of X has all non-zero elements
- dim = 1 returns a column vector (of type umat), with each element (0 or 1) indicating whether the corresponding row of X has all non-zero elements

Relational operators can be used instead of X, e.g. A > 0.5

Examples:

V = mat(10,1, fill.randu)
X = mat(5, 5, fill.randu)

# status1 will contain 1 if vector V has all non-zero elements
status1 = all(V)

# status2 will contain 1 if vector V has all elements greater than 0.5
status2 = all(V > 0.5)

# status3 will be contain 1 if matrix X has all elements greater than 0.6;
# note the use of vectorise()
status3 = all(vectorise(X) > 0.6)

# generate a row vector indicating which columns of X have all elements greater than 0.7
uA = all(X > 0.7)

See also:

any( X )		(form 1)
any( X, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, return a 1x1 matrix (of type umat) with its element (0 or 1) indicating whether any element is non-zero or satisfies a relational condition
- otherwise, form 2 is used with dim = 0
For form 2,
- dim = 0 returns a row vector (of type umat), with each element (0 or 1) indicating whether the corresponding column of X has any non-zero elements
- dim = 1 returns a column vector (of type umat), with each element (0 or 1) indicating whether the corresponding row of X has any non-zero elements

Relational operators can be used instead of X, e.g. A > 0.9

Examples:

V = mat(10,1, fill.randu)
X = mat(5, 5, fill.randu)

# status1 will contain 1 if matrix V has any non-zero elements
status1 = any(V)

# status2 will contain 1 if matrix V has any elements greater than 0.5
status2 = any(V > 0.5)

# status3 will contain 1 if matrix X has any elements greater than 0.6;
# note the use of vectorise()
status3 = any(vectorise(X) > 0.6)

# generate a row vector indicating which columns of X have elements greater than 0.7
uA = any(X > 0.7)

See also:

approx_equal( A, B, method, tol )
approx_equal( A, B, method, abs_tol, rel_tol )

Return True if all corresponding elements in A and B are approximately equal

Return False if any of the corresponding elements in A and B are not approximately equal, or if A and B have different dimensions

The argument method controls how the approximate equality is determined; it is one of:

`"absdiff"`	↦	scalars x and y are considered equal if \|x − y\| ≤ tol
`"reldiff"`	↦	scalars x and y are considered equal if \|x − y\| / max( \|x\|, \|y\| ) ≤ tol
`"both"`	↦	scalars x and y are considered equal if \|x − y\| ≤ abs_tol or \|x − y\| / max( \|x\|, \|y\| ) ≤ rel_tol

Examples:

A = mat(5, 5, fill.randu)
B = A + 0.001

same1 = approx_equal(A, B, "absdiff", 0.002)

C = 1000 * mat(5,5, fill.randu)
D = C + 1

same2 = approx_equal(C, D, "reldiff", 0.1)

same3 = approx_equal(C, D, "both", 2, 0.1)

See also:

arg( X )

Obtain the phase angle (in radians) of each element

Non-complex elements are treated as complex elements with zero imaginary component

Examples:

A = cx_mat(5, 5, fill.randu)
B = arg(A)

See also:

as_scalar( expression )

Evaluate an expression that results in a 1x1 matrix, followed by converting the 1x1 matrix to a pure scalar

Examples:

r = mat(1, 5, fill.randu)
q = mat(5, 1, fill.randu)

X = mat(5, 5, fill.randu)

a = as_scalar(r*q)
b = as_scalar(r*X*q)
c = as_scalar(r*diagmat(X)*q)
d = as_scalar(r*inv(diagmat(X))*q)

See also:
- .as_col() / .as_row()
- vectorise()
- accu()
- trace()
- dot()
- norm()

clamp( X, min_val, max_val )

Create a copy of X with each element clamped to the [min_val, max_val] interval;
any value lower than min_val will be set to min_val, and any value higher than max_val will be set to max_val

Examples:

A = mat(5, 5, fill.randu )

B = clamp(A, 0.2,     0.8) 

C = clamp(A, A.min(), 0.8) 

D = clamp(A, 0.2, A.max())

See also:

cond( A )

Return the condition number of matrix A (the ratio of the largest singular value to the smallest)

Large condition numbers suggest that matrix A is nearly singular

The computation is based on singular value decomposition; if the decomposition fails, a RuntimeError is thrown

Caveat: the rcond() function is faster for providing an estimate of the reciprocal of the condition number

Examples:

See also:

conj( X )

Obtain the complex conjugate of each element in a complex matrix or cube

Examples:

X = cx_mat(5, 5, fill.randu)
Y = conj(X)

See also:

cross( A, B )

Calculate the cross product between A and B, under the assumption that A and B are 3 dimensional vectors

Examples:

a = mat(3, 1, fill.randu)
b = mat(3, 1, fill.randu)

c = cross(a,b)

See also:

cumsum( X )		(form 1)
cumsum( X, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, return a vector of the same orientation, containing the cumulative sum of elements
- otherwise, form 2 is used with dim = 0

For form 2, return a matrix containing the cumulative sum of elements in each column (dim = 0), or each row (dim = 1)

Examples:

A = mat(5, 5, fill.randu)
B = cumsum(A)
C = cumsum(A, 1)

x = mat(10, 1, fill.randu)
y = cumsum(x)

See also:
- cumprod()
- accu()
- sum()
- diff()

cumprod( X )		(form 1)
cumprod( X, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, return a vector of the same orientation, containing the cumulative product of elements
- otherwise, form 2 is used with dim = 0

For form 2, return a matrix containing the cumulative product of elements in each column (dim = 0), or each row (dim = 1)

Examples:

A = mat(5, 5, fill.randu)
B = cumprod(A)
C = cumprod(A, 1)

x = mat(10, 1, fill.randu)
y = cumprod(x)

See also:
- cumsum()
- prod()

det( A )

Determinant of square matrix A

If A is not square sized, a RuntimeError is thrown

Caveat: for large matrices log_det() is more precise than det()

Examples:

See also:

diagmat( V )
diagmat( V, k )

diagmat( X )
diagmat( X, k )

Generate a diagonal matrix from vector V or matrix X

Given vector V, generate a square matrix with the k-th diagonal containing a copy of the vector; all other elements are set to zero

Given matrix X, generate a matrix with the k-th diagonal containing a copy of the k-th diagonal of X; all other elements are set to zero

The argument k is optional; by default the main diagonal is used (k = 0)

For k > 0, the k-th super-diagonal is used (above main diagonal, towards top-right corner)

For k < 0, the k-th sub-diagonal is used (below main diagonal, towards bottom-left corner)

Examples:

A = mat(5, 5, fill.randu)
B = diagmat(A)
C = diagmat(A,1)

v = mat(5, 1, fill.randu)
D = diagmat(v)
E = diagmat(v,1)

See also:

diagvec( A )
diagvec( A, k )

Extract the k-th diagonal from matrix A

The argument k is optional; by default the main diagonal is extracted (k = 0)

For k > 0, the k-th super-diagonal is extracted (top-right corner)

For k < 0, the k-th sub-diagonal is extracted (bottom-left corner)

The extracted diagonal is interpreted as a column vector

Examples:

A = mat(5, 5, fill.randu)

d = diagvec(A)

See also:
- .diag()
- diagmat()
- trace()
- vectorise()

diff( X )		(form 1)
diff( X, k )		(form 2)
diff( X, k, dim )		(form 3)

For forms 1 and 2:
- if matrix X can be interpreted return a vector of the same orientation, containing the differences between consecutive elements
- otherwise, form 3 is used with dim = 0

For form 3, return a matrix containing the differences between consecutive elements in each column (dim = 0), or each row (dim = 1)

The optional argument k indicates that the differences are calculated recursively k times; by default k=1 is used

The resulting number of differences is n - k, where n is the number of elements; if n ≤ k, the number of differences is zero (i.e. an empty vector/matrix is returned)

Examples:

See also:

dot( A, B )
cdot( A, B )
norm_dot( A, B )

dot(A,B): dot product of A and B, treating A and B as vectors

cdot(A,B): as per dot(A,B), but the complex conjugate of A is used

norm_dot(A,B): normalised dot product; equivalent to dot(A,B) / (∥A∥·∥B∥)

Caveat: norm() is more robust for calculating the norm, as it handles underflows and overflows

Examples:

a = mat(10, 1, fill.randu)
b = mat(10, 1, fill.randu)

x = dot(a,b)

See also:
- norm()
- as_scalar()
- cross()
- conj()

eps( X )

Obtain the positive distance of the absolute value of each element of X to the next largest representable floating point number

X can be a scalar (e.g. double), vector or matrix

Examples:

See also:

B = expmat( A )
expmat( B, A )

Matrix exponential of general square matrix A

If A is not square sized, a RuntimeError is thrown

If the matrix exponential cannot be found:
- B = expmat(A) resets B
- expmat(B,A) resets B and returns False

Caveat: if matrix A is symmetric, using expmat_sym() is faster

Caveat: the matrix exponential operation is generally not the same as applying the exp() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = expmat(A)

See also:

B = expmat_sym( A )
expmat_sym( B, A )

Matrix exponential of symmetric/hermitian matrix A

The computation is based on eigen decomposition

If A is not square sized, a RuntimeError is thrown

If the matrix exponential cannot be found:
- B = expmat_sym(A) resets B
- expmat_sym(B,A) resets B and returns False

Caveat: the matrix exponential operation is generally not the same as applying the exp() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = A*A.t()   # make symmetric matrix

C = expmat_sym(B)

See also:

find( X )
find( X, k )
find( X, k, s )

Return a column vector containing the indices of elements of X that are non-zero or satisfy a relational condition

X is interpreted as a vector, with column-by-column ordering of the elements of X

Relational operators can be used instead of X, e.g. A > 0.5

If k = 0 (default), return the indices of all non-zero elements, otherwise return at most k of their indices

If s="first" (default), return at most the first k indices of the non-zero elements

If s="last", return at most the last k indices of the non-zero elements

Caveats:
- to clamp values to an interval, clamp() is more efficient
- to replace a specific value, .replace() is more efficient

Examples:

A = mat(5, 5, fill.randu)
B = mat(5, 5, fill.randu)

q1 = find(A > B)
q2 = find(A > 0.5)
q3 = find(A > 0.5, 3, "last")

# change elements of A greater than 0.5 to 1
A[ find(A > 0.5) ].ones()

See also:

find_finite( X )

Return a column vector containing the indices of elements of X that are finite (i.e. not ±Inf and not NaN)

X is interpreted as a vector, with column-by-column ordering of the elements of X

Examples:

A = mat(5, 5, fill.randu)

A[1,1] = datum.inf

# accumulate only finite elements
val = accu( A[ find_finite(A) ] )

See also:

find_nonfinite( X )

Return a column vector containing the indices of elements of X that are non-finite (i.e. ±Inf or NaN)

X is interpreted as a vector, with column-by-column ordering of the elements of X

Examples:

A = mat(5, 5, fill.randu)

A[1,1] = datum.inf
A[2,2] = datum.nan

# change non-finite elements to zero
A[ find_nonfinite(A) ].zeros()

Caveat: to replace instances of a specific non-finite value (e.g. nan or inf), it is more efficient to use .replace()

See also:

find_unique( X )
find_unique( X, ascending_indices )

Return a column vector containing the indices of unique elements of X

X is interpreted as a vector, with column-by-column ordering of the elements of X

The ascending_indices argument is optional; it is one of:

Examples:

A = mat([ [ 2, 2, 4 ], 
          [ 4, 6, 6 ] ])

indices = find_unique(A)

See also:

fliplr( X )
flipud( X )

fliplr(): generate a copy of matrix X, with the order of the columns reversed

flipud(): generate a copy of matrix X, with the order of the rows reversed

Examples:

A = mat(5, 5, fill.randu)

B = fliplr(A)
C = flipud(A)

See also:

imag( X )
real( X )

Extract the imaginary/real part of a complex matrix or cube

Examples:

C = cx_mat(5, 5, fill.randu)

A = imag(C)
B = real(C)

See also:

ind2sub( size(X), index )		(form 1)
ind2sub( size(X), vector_of_indices )		(form 2)

Convert a linear index, or a vector of indices, to subscript notation

The argument size(X) can be replaced with size(n_rows, n_cols)

A RuntimeError is thrown if an index is out of range

When only one index is given (form 1), the subscripts are returned in a vector of type umat

When a vector of indices (of type umat) is given (form 2), the corresponding subscripts are returned in each column of an m x n matrix of type umat; m=2 for matrix subscripts, while m=3 for cube subscripts

Examples:

M = mat(4, 5, fill.randu)

s = ind2sub( size(M), 6 )

print("row: " + str(s[0]))
print("col: " + str(s[1]))


indices = find(M > 0.5)
t       = ind2sub( size(M), indices )


Q = cube(2, 3, 4)

u = ind2sub( size(Q), 8 )

print("row: " + str(s[0]))
print("col: " + str(s[1]))
print("slice: " + str(s[2]))

See also:

index_min( M )
index_min( M, dim )
index_min( Q )
index_min( Q, dim )

index_max( M )
index_max( M, dim )
index_max( Q )
index_max( Q, dim )

(form 1)
(form 2)

For form 1:
- if matrix M can be interpreted as a vector, return the linear index of the extremum value; the returned index is a 1x1 matrix (of type umat)
- otherwise, form 2 is used with dim = 0

For matrix M and:
- dim = 0, return a row vector (of type umat), with each column containing the index of the extremum value in the corresponding column of M
- dim = 1, return a column vector (of type umat), with each row containing the index of the extremum value in the corresponding row of M

For cube Q, return a cube (of type ucube) containing the indices of extremum values of elements along dimension dim, where dim ∈ { 0, 1, 2 }

For each column or row, the index starts at zero

For objects with complex numbers, absolute values are used for comparison

Examples:

v = mat(10, 1, fill.randu)

i = index_max(v)
max_val_in_v = v[i]


M = mat(5, 6, fill.randu)

ii = index_max(M)
jj = index_max(M,1)

max_val_in_col_2 = M[ ii[2], 2 ]

max_val_in_row_4 = M[ 4, jj[4] ]

See also:
- min() & max()
- .index_min() & .index_max() (member functions)
- sort_index()
- find()

inplace_trans( X )
inplace_strans( X )

In-place / in-situ transpose of matrix X

For real (non-complex) matrix:
- inplace_trans() performs a normal transpose
- inplace_strans() not applicable

For complex matrix:
- inplace_trans() performs a Hermitian transpose (i.e. the conjugate of the elements is taken during the transpose)
- inplace_strans() provides a transposed copy without taking the conjugate of the elements

Examples:

X = mat(4,     5,     fill.randu)
Y = mat(20000, 30000, fill.randu)

inplace_trans(X)

inplace_trans(Y)

See also:

C = intersect( A, B )		(form 1)
intersect( C, iA, iB, A, B )		(form 2)

For form 1:
- return the unique elements common to both A and B, sorted in ascending order

For form 2,
- store in C the unique elements common to both A and B, sorted in ascending order
- store in iA and iB the indices of the unique elements, such that C = A.elem(iA) and C = B.elem(iB)
- iA and iB must have the type umat (i.e. the indices are stored as unsigned integers)

C is a column vector if either A or B is a matrix or column vector; C is a row vector if both A and B are row vectors

For matrices and vectors with complex numbers, ordering is via absolute values

Examples:

A = regspace(4, 1)  # 4, 3, 2, 1
B = regspace(3, 6)  # 3, 4, 5, 6

C = intersect(A,B)  # 3, 4

CC = mat()
iA = umat()
iB = umat()

intersect(CC, iA, iB, A, B)

See also:

join_rows( A, B )
join_rows( A, B, C )
join_rows( A, B, C, D )

join_cols( A, B )
join_cols( A, B, C )
join_cols( A, B, C, D )

join_horiz( A, B )
join_horiz( A, B, C )
join_horiz( A, B, C, D )

join_vert( A, B )
join_vert( A, B, C )
join_vert( A, B, C, D )

join_rows() and join_horiz(): horizontal concatenation; join the corresponding rows of the given matrices; the given matrices must have the same number of rows

join_cols() and join_vert(): vertical concatenation; join the corresponding columns of the given matrices; the given matrices must have the same number of columns

Examples:

A = mat(4, 5, fill.randu)
B = mat(4, 6, fill.randu)
C = mat(6, 5, fill.randu)

AB = join_rows(A,B)
AC = join_cols(A,C)

See also:

join_slices( cube C, cube D )
join_slices( mat M, mat N )

join_slices( mat M, cube C )
join_slices( cube C, mat M )

for two cubes C and D: join the slices of C with the slices of D; cubes C and D must have the same number of rows and columns (ie. all slices must have the same size)

for two matrices M and N: treat M and N as cube slices and join them to form a cube with 2 slices; matrices M and N must have the same number of rows and columns

for matrix M and cube C: treat M as a cube slice and join it with the slices of C; matrix M and cube C must have the same number of rows and columns

Examples:

C = cube(5, 10, 3, fill.randu)
D = cube(5, 10, 4, fill.randu)

E = join_slices(C,D)

M = mat(10, 20, fill.randu)
N = mat(10, 20, fill.randu)

Q = join_slices(M,N)

R = join_slices(Q,M)

S = join_slices(M,Q)

See also:

kron( A, B )

Kronecker tensor product

Given matrix A (with n rows and p columns) and matrix B (with m rows and q columns), generate a matrix (with nm rows and pq columns) that denotes the tensor product of A and B

Examples:

A = mat(4, 5, fill.randu)
B = mat(5, 4, fill.randu)

K = kron(A,B)

See also:

result = log_det( A )

Complex log determinant of square matrix A

If A is not square sized, a RuntimeError is thrown

If matrix A is real and the determinant is positive:
- the real part of the result is the log determinant
- the imaginary part is zero

If matrix A is real and the determinant is negative:
- the real part of the result is the log of the absolute value of the determinant
- the imaginary part is equal to datum.pi

Examples:

A = mat(5, 5, fill.randu)
result = log_det(A)

See also:

B = logmat( A )
logmat( B, A )

Complex matrix logarithm of general square matrix A

If A is not square sized, a RuntimeError is thrown

If the matrix logarithm cannot be found:
- B = logmat(A) resets B
- logmat(B,A) resets B and returns False

Caveat: if matrix A is symmetric positive definite, using logmat_sympd() is faster

Caveat: the matrix logarithm operation is generally not the same as applying the log() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = logmat(A)

See also:

B = logmat_sympd( A )
logmat_sympd( B, A )

Matrix logarithm of symmetric/hermitian positive definite matrix A

If A is not square sized, a RuntimeError is thrown

If the matrix logarithm cannot be found:
- B = logmat_sympd(A) resets B
- logmat_sympd(B,A) resets B and returns False

Caveat: the matrix logarithm operation is generally not the same as applying the log() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = A*A.t()   # make symmetric matrix

C = logmat_sympd(B)

See also:

min( X )
min( X, dim )
min( A, B )

max( X )
max( X, dim )
max( A, B )

(form 1)
(form 2)
(form 3)

Form 1:
- if matrix M can be interpreted as a vector, return the extremum value as a 1x1 matrix
- otherwise, form 2 is used with dim = 0

Form 2:
- for matrix M, return the extremum value for each column (dim = 0), or each row (dim = 1)
- for cube Q, return the extremum values of elements along dimension dim, where dim ∈ { 0, 1, 2 }

For two matrices/cubes A and B, return a matrix/cube containing element-wise extremum values

For objects with complex numbers, absolute values are used for comparison

Examples:

v = mat(10, 1, fill.randu)
x = max(v)

M = mat(10, 10, fill.randu)

a = max(M)
b = max(M,0) 
c = max(M,1)

# element-wise maximum
X = mat(5, 6, fill.randu)
Y = mat(5, 6, fill.randu)
Z = max(X,Y)

See also:
- .min() & .max() (member functions)
- index_min() & index_max()
- clamp()
- statistics functions
- ensmallen - library for finding minimum of arbitrary function

nonzeros( X )

Return a column vector containing the non-zero values of X

Examples:

B = mat(100, 100, fill.eye)
b = nonzeros(B)

See also:

norm( X )
norm( X, p )

Compute the p-norm of X, where X can be a vector or matrix

For vectors, p is an integer ≥1, or one of: "-inf", "inf", "fro"

For matrices, p is one of: 1, 2, "inf", "fro"; the calculated norm is the induced norm (not entrywise norm)

"-inf" is the minimum norm, "inf" is the maximum norm, while "fro" is the Frobenius norm

The argument p is optional; by default p = 2 is used

For vector norm with p = 2 and matrix norm with p = "fro", a robust algorithm is used to reduce the likelihood of underflows and overflows

To obtain the zero norm or Hamming norm (i.e. the number of non-zero elements), use this expression: accu(X != 0)

Examples:

q = mat(5, 1, fill.randu)

x = norm(q, 2)
y = norm(q, "inf")

See also:

normalise( X )		(form 1)
normalise( X, p )		(form 2)
normalise( X, p, dim )		(form 3)

For forms 1 and 2:
- if matrix X can be interpreted as a vector, return its normalised version (i.e. having unit p-norm)
- otherwise, form 3 is used with dim = 0

For form 3, return its normalised version, where each column (dim = 0) or row (dim = 1) has been normalised to have unit p-norm

The p argument is optional; by default p = 2 is used

Examples:

A = mat(10, 1, fill.randu)
B = normalise(A)
C = normalise(A, 1)

X = mat(5, 6, fill.randu)
Y = normalise(X)
Z = normalise(X, 2, 1)

See also:

prod( X )		(form 1)
prod( X, p )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, return the product of all elements as a 1x1 matrix
- otherwise, form 2 is used with dim = 0

For form 2, return the product of elements in each column (dim = 0), or each row (dim = 1)

Examples:

v = mat(10, 1, fill.randu)
x = prod(v)

M = mat(10, 10, fill.randu)

a = prod(M)
b = prod(M,0)
c = prod(M,1)

See also:
- cumprod()
- Schur product

B = powmat( A, n )
powmat( B, A, n )

Matrix power operation: raise square matrix A to the power of n, where n is an integer or floating point number

If n is a floating point number, the resultant matrix B always has complex elements

For n = 0, an identity matrix is generated

If A is not square sized, a RuntimeError is thrown

If the matrix power cannot be found:
- B = powmat(A) resets B
- powmat(B,A) resets B and returns False

Caveats:
- to find the inverse of a matrix, use inv() instead
- to solve a system of linear equations, use solve() instead
- to find the matrix square root, use sqrtmat() instead
- the matrix power operation is generally not the same as applying the pow() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = powmat(A, 4)     #     integer exponent

C = powmat(A, 4.56)  # non-integer exponent

See also:

rank( X )
rank( X, tolerance )

Return the rank of matrix X

Any singular values less than tolerance are treated as zero

The tolerance argument is optional; by default tolerance is max(m,n)*max_sv*datum.eps, where:
- m = number of rows and n = number of columns in X
- max_sv = maximal singular value of X
- datum.eps = difference between 1 and the least value greater than 1 that is representable

The computation is based on singular value decomposition; if the decomposition fails, a RuntimeError is thrown

Caveat: to confirm whether a matrix is singular, use rcond() or cond()

Examples:

See also:
- cond()
- rcond()
- svd()
- orth()
- datum.eps
- Rank in MathWorld
- Rank in Wikipedia

rcond( A )

Return the 1-norm estimate of the reciprocal of the condition number of square matrix A

Values close to 1 suggest that A is well-conditioned

Values close to 0 suggest that A is badly conditioned

If A is not square sized, a RuntimeError is thrown

Examples:

A = mat(5, 5, fill.randu)

r = rcond(A)

See also:
- cond()
- rank()
- det()
- inv()
- solve()

repelem( A, num_copies_per_row, num_copies_per_col )

Generate a matrix by replicating each element of matrix A

The generated matrix has the following size:

Examples:

A = mat(2, 3, fill.randu)

B = repelem(A, 4, 5)

See also:
- repmat()
- kron()
- reshape()
- resize()

repmat( A, num_copies_per_row, num_copies_per_col )

Generate a matrix by replicating matrix A in a block-like fashion

The generated matrix has the following size:

Caveat: to apply a vector operation on each row or column of a matrix, it is generally more efficient to use broadcasting

Examples:

A = mat(2, 3, fill.randu)

B = repmat(A, 4, 5)

See also:
- broadcasting
- repelem()
- kron()
- reshape()
- resize()

reshape( X, n_rows, n_cols ) (X is a vector or matrix)
reshape( X, size(Y) )

reshape( Q, n_rows, n_cols, n_slices ) (Q is a cube)
reshape( Q, size(R) )

Generate a matrix/cube with given size specifications, whose elements are taken from the given object in a column-wise manner; the elements in the generated object are placed column-wise (i.e. the first column is filled up before filling the second column)

The layout of the elements in the generated object will be different to the layout in the given object

If the total number of elements in the given object is less than the specified size, the remaining elements in the generated object are set to zero

If the total number of elements in the given object is greater than the specified size, only a subset of elements is taken from the given object

Caveats:
- do not use reshape() if you simply want to change the size without preserving data; use .set_size() instead, which is much faster
- to grow/shrink a matrix while preserving the elements as well as the layout of the elements, use resize() instead
- to flatten a matrix into a vector, use vectorise() or .as_col() / .as_row() instead

Examples:

A = mat(10, 5, fill.randu)

B = reshape(A, 5, 10)

See also:
- .reshape() (member function)
- .set_size()
- resize()
- vectorise()
- as_scalar()
- diagmat()
- repmat()
- repelem()
- size()

resize( X, n_rows, n_cols ) (X is a vector or matrix)
resize( X, size(Y) )

resize( Q, n_rows, n_cols, n_slices ) (Q is a cube)
resize( Q, size(R) )

Generate a matrix/cube with given size specifications, whose elements as well as the layout of the elements are taken from the given object

Caveat: do not use resize() if you simply want to change the size without preserving data; use .set_size() instead, which is much faster

Examples:

A = mat(4, 5, fill.randu)

B = resize(A, 7, 6)

See also:
- .resize() (member function of mat and cube)
- .set_size() (member function of mat and cube)
- reshape()
- vectorise()
- as_scalar()
- repmat()
- repelem()
- size()

reverse( X )		(form 1)
reverse( X, p )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, generate a copy of the vector with the order of elements reversed
- otherwise, form 2 is used with dim = 0

For form 2, generate a copy of the matrix with the order of elements reversed in each column (dim = 0), or each row (dim = 1)

Examples:

v = mat(123, 1, fill.randu)
y = reverse(v)

A = mat(4, 5, fill.randu)
B = reverse(A)
C = reverse(A,1)

See also:

R = roots( P )
roots( R, P )

Find the complex roots of a polynomial function represented via vector P and store them in column vector R

The polynomial function is modelled as:
where p_{_i} is the i-th polynomial coefficient in vector P

The computation is based on eigen decomposition; if the decomposition fails:
- R = roots(P) resets R
- roots(R,P) resets R and returns False

Examples:

P = mat(5, 1, fill.randu)
  
R = roots(P)

See also:

shift( X, N )		(form 1)
shift( X, N, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, generate a copy of the vector with the elements shifted by N positions in a circular manner
- otherwise, form 2 is used with dim = 0

For form 2, generate a copy of the matrix with the elements shifted by N positions in each column (dim = 0), or each row (dim = 1)

N can be positive or negative

Examples:

A = mat(4, 5, fill.randu)
B = shift(A, -1)
C = shift(A, +1)

See also:

shuffle( X, )		(form 1)
shuffle( X, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, generate a copy of the vector with the elements shuffled
- otherwise, form 2 is used with dim = 0

For form 2, generate a copy of the matrix with the elements shuffled in each column (dim = 0), or each row (dim = 1)

Examples:

A = mat(4, 5, fill.randu)
B = shuffle(A)

See also:

size( X )
size( n_rows, n_cols )
size( n_rows, n_cols, n_slices )

Obtain the dimensions of object X, or explicitly specify the dimensions

The dimensions can be used in conjunction with:
- object constructors
- functions for changing size
- generator functions
- submatrix views
- subcube views
- ind2sub() and sub2ind()

The dimensions support simple arithmetic operations; they can also be printed and compared for equality/inequality

Examples:

A = mat(5,6)

B = mat(size(A), fill.zeros)

C = mat()

C.randu(size(A))

D = mat(10,20, fill.ones)
D[3,4,size(C)] = C    # access submatrix of E

E = mat( size(A) + size(E) )

F = mat( size(A) * 2 )

print("size of A: " + str(size(A)))

is_same_size = (size(A) == size(E))

See also:
- attributes

sort( X )		(form 1)
sort( X, sort_direction )		(form 2)
sort( X, sort_direction, dim )		(form 3)

For forms 1 and 2:
- if matrix X can be interpreted as a vector, return a vector which is a sorted version of the input vector
- otherwise, form 3 is used with dim = 0

For form 3, return a matrix with the elements of the input matrix sorted in each column (dim = 0), or each row (dim = 1)

The sort_direction argument is optional; sort_direction is either "ascend" or "descend"; by default "ascend" is used

For matrices and vectors with complex numbers, sorting is via absolute values

Examples:

A = mat(10, 10, fill.randu)
B = sort(A)

See also:

sort_index( X )
sort_index( X, sort_direction )

stable_sort_index( X )
stable_sort_index( X, sort_direction )

Return a vector which describes the sorted order of the elements of X (i.e. it contains the indices of the elements of X)

X is interpreted as a vector, with column-by-column ordering of the elements of X

The sort_direction argument is optional; sort_direction is either "ascend" or "descend"; by default "ascend" is used

The stable_sort_index() variant preserves the relative order of elements with equivalent values

For matrices and vectors with complex numbers, sorting is via absolute values

Examples:

q = mat(10, 1, fill.randu)

indices = sort_index(q)

See also:

B = sqrtmat( A )
sqrtmat( B, A )

Complex square root of general square matrix A

If A is not square sized, a RuntimeError is thrown

If matrix A appears to be singular, an approximate square root is attempted; additionally, sqrtmat(B,A) returns False

Caveat: if matrix A is symmetric positive definite, using sqrtmat_sympd() is faster

Caveat: the square root of a matrix is generally not the same as applying the sqrt() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = sqrtmat(A)

See also:

B = sqrtmat_sympd( A )
sqrtmat_sympd( B, A )

Square root of symmetric/hermitian positive definite matrix A

If A is not square sized, a RuntimeError is thrown

If the square root cannot be found:
- B = sqrtmat_sympd(A) resets B
- sqrtmat_sympd(B,A) resets B and returns False

Caveat: the matrix square root operation is generally not the same as applying the sqrt() function to each element

Examples:

A = mat(5, 5, fill.randu)

B = A*A.t()   # make symmetric matrix

C = sqrtmat_sympd(B)

See also:

sum( X )		(form 1)
sum( X, dim )		(form 2)

For form 1:
- if matrix M can be interpreted as a vector, return the sum of all elements as a 1x1 matrix
- otherwise, form 2 is used with dim = 0

For form 2:
- for matrix M, return the sum of elements in each column (dim = 0), or each row (dim = 1)
- for cube Q, return the sums of elements along dimension dim, where dim ∈ { 0, 1, 2 }; for example, dim = 0 indicates the sum of elements in each column within each slice

Caveat: to get a sum of all the elements regardless of the object type (i.e. vector, or matrix), use accu() instead

Examples:

v = mat(10, 1, fill.randu)
x = sum(v)

M = mat(10, 10, fill.randu)

a = sum(M)
b = sum(M,0)
c = sum(M,1)

y = accu(M)   # find the overall sum regardless of object type

See also:
- accu()
- cumsum()
- trace()
- trapz()
- mean()
- as_scalar()

index	= sub2ind( size(M), row, col )	(M is a matrix)
indices	= sub2ind( size(M), matrix_of_subscripts )

index	= sub2ind( size(Q), row, col, slice)	(Q is a cube)
indices	= sub2ind( size(M), matrix_of_subscripts )

Convert subscripts to a linear index

The argument size(X) can be replaced with size(n_rows, n_cols) or size(n_rows, n_cols, n_slices)

For the matrix_of_subscripts argument, the subscripts will be stored in each column of an m x n matrix of type umat; m=2 for matrix subscripts

A RuntimeError is thrown if a subscript is out of range

Examples:

M = mat(4,5)
Q = cube(4,5,6)

i = sub2ind( size(M), 2, 3 )
j = sub2ind( size(Q), 2, 3, 4 )

See also:

symmatu( A )
symmatu( A, do_conj )

symmatl( A )
symmatl( A, do_conj )

symmatu(A): generate symmetric matrix from square matrix A, by reflecting the upper triangle to the lower triangle

symmatl(A): generate symmetric matrix from square matrix A, by reflecting the lower triangle to the upper triangle

If A is a complex matrix, the reflection uses the complex conjugate of the elements; to disable the complex conjugate, set do_conj to False

If A is non-square, a RuntimeError is thrown

Examples:

A = mat(5, 5, fill.randu)

B = symmatu(A)
C = symmatl(A)

See also:

trace( X )

Sum of the elements on the main diagonal of matrix X

Examples:

A = mat(5, 5, fill.randu)

x = trace(A)

See also:
- accu()
- as_scalar()
- .diag()
- diagvec()
- sum()

trans( A )
strans( A )

For real (non-complex) matrix:
- trans() provides a transposed copy of the matrix
- strans() not applicable

For complex matrix:
- trans() provides a Hermitian transpose (i.e. the conjugate of the elements is taken during the transpose)
- strans() provides a transposed copy without taking the conjugate of the elements

Examples:

A = mat(5, 10, fill.randu)

B = trans(A)
C = A.t()    # equivalent to trans(A), but more compact

See also:
- .t()
- inplace_trans()

trapz( X, Y )
trapz( X, Y, dim )

trapz( Y )
trapz( Y, dim )

Compute the trapezoidal integral of Y with respect to spacing in X, in each column (dim = 0) or each row (dim = 1) of Y

X must be a vector; its length must equal either the number of rows in Y (when dim = 0), or the number of columns in Y (when dim = 1)

If X is not specified, unit spacing is used

If the dim argument is not specified and Y can be interpreted as a row vector, dim is set to 1

Examples:

X = linspace(0, datum.pi, 1000)
Y = sin(X)

Z = trapz(X,Y)

See also:

trimatu( A )
trimatu( A, k )

trimatl( A )
trimatl( A, k )

Create a new matrix by copying either the upper or lower triangular part from square matrix A, and setting the remaining elements to zero
- trimatu() copies the upper triangular part
- trimatl() copies the lower triangular part

The argument k specifies the diagonal which inclusively delineates the boundary of the triangular part
- for k > 0, the k-th super-diagonal is used (above main diagonal, towards top-right corner)
- for k < 0, the k-th sub-diagonal is used (below main diagonal, towards bottom-left corner)

The argument k is optional; by default the main diagonal is used (k = 0)

If A is non-square, a RuntimeError is thrown

Examples:

A = mat(5, 5, fill.randu)

U  = trimatu(A)
L  = trimatl(A)

UU = trimatu(A,  1)  # omit the main diagonal
LL = trimatl(A, -1)  # omit the main diagonal

See also:

trimatu_ind( size(A) )
trimatu_ind( size(A), k )

trimatl_ind( size(A) )
trimatl_ind( size(A), k )

Return a column vector containing the indices of elements that form the upper or lower triangle part of matrix A
- trimatu_ind() refers to the upper triangular part
- trimatl_ind() refers to the lower triangular part

The argument k specifies the diagonal which inclusively delineates the boundary of the triangular part
- for k > 0, the k-th super-diagonal is used (above main diagonal, towards top-right corner)
- for k < 0, the k-th sub-diagonal is used (below main diagonal, towards bottom-left corner)

The argument k is optional; by default the main diagonal is used (k = 0)

The argument size(A) can be replaced with size(n_rows, n_cols)

Examples:

A = mat(5, 5, fill.randu)

upper_indices = trimatu_ind( size(A) )
lower_indices = trimatl_ind( size(A) )

# extract upper/lower triangle into vector
upper_part = A[upper_indices]
lower_part = A[lower_indices]

# obtain indices without the main diagonal
alt_upper_indices = trimatu_ind( size(A),  1)
alt_lower_indices = trimatl_ind( size(A), -1)

See also:

unique( A )

Return the unique elements of A, sorted in ascending order

If A is a vector, the output is also a vector with the same orientation (row or column) as A; if A is a matrix, the output is always a column vector

Examples:

X = mat([ [ 1, 2 ],
          [ 2, 3 ] ])

Y = unique(X)

See also:
- find()
- find_unique()
- sort()
- shuffle()
- nonzeros()
- intersect()

vectorise( X )
vectorise( X, dim )
vectorise( Q )

Generate a flattened version of matrix X or cube Q

The argument dim is optional; by default dim = 0 is used

For dim = 0, the elements are copied from X column-wise, resulting in a column vector; equivalent to concatenating all the columns of X

For dim = 1, the elements are copied from X row-wise, resulting in a row vector; equivalent to concatenating all the rows of X

Caveat: column-wise vectorisation is faster than row-wise vectorisation

Examples:

X = mat(4, 5, fill.randu)

v = vectorise(X)

See also:

miscellaneous element-wise functions:

exp

log

pow

floor

erf

tgamma

sign

exp2

log2

square

ceil

erfc

lgamma

exp10

log10

sqrt

round

expm1

log1p

trunc

trunc_exp

trunc_log

Apply a function to each element

Usage:

B = fn(A), where fn(A) is one of the functions below
A and B must have the same matrix type or cube type, such as mat or cube

exp(A) base-e exponential: e^x

exp2(A) base-2 exponential: 2^x

exp10(A) base-10 exponential: 10^x

expm1(A) compute exp(A)-1 accurately for values of A close to zero (only for floating point elements)

trunc_exp(A) base-e exponential, truncated to avoid infinity (only for floating point elements)

log(A) natural log: log_e x

log2(A) base-2 log: log₂ x

log10(A) base-10 log: log₁₀ x

log1p(A) compute log(1+A) accurately for values of A close to zero (only for floating point elements)

trunc_log(A) natural log, truncated to avoid ±infinity (only for floating point elements)

pow(A, p) raise to the power of p: x^p

square(A) square: x²

sqrt(A) square root: √x

floor(A) largest integral value that is not greater than the input value

ceil(A) smallest integral value that is not less than the input value

round(A) round to nearest integer, with halfway cases rounded away from zero

trunc(A) round to nearest integer, towards zero

erf(A) error function (only for floating point elements)

erfc(A) complementary error function (only for floating point elements)

tgamma(A) gamma function (only for floating point elements)

lgamma(A) natural log of the absolute value of gamma function (only for floating point elements)

sign(A)

signum function; for each element a in A, the corresponding element b in B is:

	⎧	−1	if a < 0
b =	⎨	0	if a = 0
	⎩	+1	if a > 0

if a is complex and non-zero, then b = a / abs(a)

Caveats:
- all of the above functions are applied element-wise, where each element is treated independently
- the element-wise functions exp(), log(), sqrt() and pow() have the corresponding functions expmat(), logmat(), sqrtmat() and powmat() which take into account matrix structure

Examples:

See also:

trigonometric element-wise functions (cos, sin, tan, ...)

For single argument functions, B = trig_fn(A), where trig_fn is applied to each element in A, with trig_fn as one of:
- cos, acos, cosh, acosh
- sin, asin, sinh, asinh
- tan, atan, tanh, atanh
- sinc, defined as sinc(x) = sin(πx) / (πx) for x ≠ 0, and sinc(x) = 1 for x = 0

For dual argument functions, apply the function to each tuple of two corresponding elements in X and Y:
- Z = atan2(Y, X)
- Z = hypot(X, Y)

Examples:
See also:

Decompositions, Factorisations, Inverses and Equation Solvers

R = chol( X )
R = chol( X, layout )

chol( R, X )
chol( R, X, layout )

Cholesky decomposition of matrix X

Matrix X must be symmetric/hermitian and positive-definite

By default, R is upper triangular, such that R.t()*R = X

The argument layout is optional; layout is either "upper" or "lower", which specifies whether R is upper or lower triangular

If the decomposition fails:
- R = chol(X) resets R
- chol(R,X) resets R and returns False

Examples:

X = mat(5, 5, fill.randu)
Y = X.t()*X

R1 = chol(Y)
R2 = chol(Y, "lower")

See also:

eigval, eigvec = eig_sym( X )

eig_sym( eigval, X )

eig_sym( eigval, eigvec, X )

Eigen decomposition of dense symmetric/hermitian matrix X

The eigenvalues and corresponding eigenvectors are stored in eigval and eigvec, respectively

The eigenvalues are in ascending order

The eigenvectors are stored as column vectors

If X is not square sized, a RuntimeError is thrown

If the decomposition fails:
- eigval, eigvec = eig_sym(X) resets eigval & eigvec
- eig_sym(eigval,X) resets eigval and returns False
- eig_sym(eigval,eigvec,X) resets eigval & eigvec and returns False

Examples:

# for matrices with real elements

A = mat(50, 50, fill.randu)
B = A.t()*A  # generate a symmetric matrix

eigval = mat()
eigvec = mat()

eig_sym(eigval, eigvec, B)


# for matrices with complex elements

C = cx_mat(50, 50, fill.randu)
D = C.t()*C

eigval2 = mat()
eigvec2 = cx_mat()

eig_sym(eigval2, eigvec2, D)

See also:

eigval, leigvec, reigvec = eig_gen( X )
eigval, leigvec, reigvec = eig_gen( X, bal )

eig_gen( eigval, X )
eig_gen( eigval, X, bal )

eig_gen( eigval, eigvec, X )
eig_gen( eigval, eigvec, X, bal )

eig_gen( eigval, leigvec, reigvec, X )
eig_gen( eigval, leigvec, reigvec, X, bal )

Eigen decomposition of dense general (non-symmetric/non-hermitian) square matrix X

The eigenvalues and corresponding right eigenvectors are stored in eigval and eigvec, respectively

If both left and right eigenvectors are requested they are stored in leigvec and reigvec, respectively

The eigenvectors are stored as column vectors

The bal argument is optional; bal is one of:

If X is not square sized, a RuntimeError is thrown

If the decomposition fails:
- eigval, leigvec, reigvec = eig_gen(X) resets eigval, leigvec & reigvec
- eig_gen(eigval,X) resets eigval and returns False
- eig_gen(eigval,eigvec,X) resets eigval & eigvec and returns False
- eig_gen(eigval,leigvec,reigvec,X) resets eigval, leigvec & reigvec and returns False

Examples:

A = mat(10, 10, fill.randu)

eigval = cx_mat()
eigvec = cx_mat()

eig_gen(eigval, eigvec, A)
eig_gen(eigval, eigvec, A, "balance")

See also:

eigval, leigvec, reigvec = eig_pair( A, B )

eig_pair( eigval, A, B )

eig_pair( eigval, eigvec, A, B )

eig_pair( eigval, leigvec, reigvec, A, B )

Eigen decomposition for pair of general dense square matrices A and B of the same size, such that A*eigvec = B*eigvec*diagmat(eigval)

The eigenvalues and corresponding right eigenvectors are stored in eigval and eigvec, respectively

If both left and right eigenvectors are requested they are stored in leigvec and reigvec, respectively

The eigenvectors are stored as column vectors

If A or B is not square sized, a RuntimeError is thrown

If the decomposition fails:
- eigval, leigvec, reigvec = eig_pair(A,B) resets eigval, leigvec & reigvec
- eig_pair(eigval,A,B) resets eigval and returns False
- eig_pair(eigval,eigvec,A,B) resets eigval & eigvec and returns False
- eig_pair(eigval,leigvec,reigvec,A,B) resets eigval, leigvec & reigvec and returns False

Examples:

A = mat(10, 10, fill.randu)
B = mat(10, 10, fill.randu)

eigval = cx_mat()
eigvec = cx_mat()

eig_pair(eigval, eigvec, A, B)

See also:

U, H = hess( X )

hess( H, X )

hess( U, H, X )

Upper Hessenberg decomposition of square matrix X, such that X = U*H*U.t()

U is a unitary matrix containing the Hessenberg vectors

H is a square matrix known as the upper Hessenberg matrix, with elements below the first subdiagonal set to zero

If X is not square sized, a RuntimeError is thrown

If the decomposition fails:
- U, H = hess(X) resets U & H
- hess(H,X) resets H and returns False
- hess(U,H,X) resets U & H and returns False

Caveat: in general, upper Hessenberg decomposition is not unique

Examples:

X = mat(20,20, fill.randu)

U = mat()
H = mat()

hess(U, H, X)

See also:

B = inv( A )
inv( B, A )

Inverse of general square matrix A

If A is not square sized, a RuntimeError is thrown

If A appears to be singular:
- B = inv(A) resets B
- inv(B,A) resets B and returns False

Caveats:
- if matrix A is know to be symmetric positive definite, using inv_sympd() is faster
- if matrix A is know to be diagonal, use inv( diagmat(A) )
- if matrix A is know to be triangular, use inv( trimatu(A) ) or inv( trimatl(A) )
- to solve a system of linear equations, such as Z = inv(X)*Y, using solve() can be faster and/or more accurate

Examples:

See also:

B = inv_sympd( A )
inv_sympd( B, A )

Inverse of symmetric/hermitian positive definite matrix A

If A is not square sized, a RuntimeError is thrown

If A appears to be singular or not positive definite:
- B = inv_sympd(A) resets B
- inv_sympd(B,A) resets B and returns False

Caveat: to solve a system of linear equations, such as Z = inv(X)*Y, using solve() can be faster and/or more accurate

Examples:

A = mat(5, 5, fill.randu)
B = A.t() * A
C = inv_sympd(B)

See also:

L, U, P = lu( X )
lu( L, U, P, X )
lu( L, U, X )

Lower-upper decomposition (with partial pivoting) of matrix X

The first form provides a lower-triangular matrix L, an upper-triangular matrix U, and a permutation matrix P, such that P.t()*L*U = X

The second form provides permuted L and U, such that L*U = X; note that in this case L is generally not lower-triangular

If the decomposition fails:
- L, U, P = lu(X) resets L, U & P
- lu(L,U,P,X) resets L, U & P and returns False
- lu(L,U,X) resets L & U and returns False

Examples:

A = mat(5, 5, fill.randu)

L = mat()
U = mat()
P = mat()

lu(L, U, P, A)

B = P.t()*L*U

See also:

B = null( A )
B = null( A, tolerance )

null( B, A )
null( B, A, tolerance )

Find the orthonormal basis of the null space of matrix A

The dimension of the range space is the number of singular values of A not greater than tolerance

The tolerance argument is optional; by default tolerance is max(m,n)*max_sv*datum.eps, where:
- m = number of rows and n = number of columns in A
- max_sv = maximal singular value of A
- datum.eps = difference between 1 and the least value greater than 1 that is representable

The computation is based on singular value decomposition; if the decomposition fails:
- B = null(A) resets B
- null(B,A) resets B and returns False

Examples:

A = mat(5, 6, fill.randu)

A[0,:].zeros()
A[:,0].zeros()

B = null(A)

See also:

B = orth( A )
B = orth( A, tolerance )

orth( B, A )
orth( B, A, tolerance )

Find the orthonormal basis of the range space of matrix A, so that B.t()*B ≈ eye(r,r), where r = rank(A)

The dimension of the range space is the number of singular values of A greater than tolerance

The tolerance argument is optional; by default tolerance is max(m,n)*max_sv*datum.eps, where:
- m = number of rows and n = number of columns in A
- max_sv = maximal singular value of A
- datum.eps = difference between 1 and the least value greater than 1 that is representable

The computation is based on singular value decomposition; if the decomposition fails:
- B = orth(A) resets B
- orth(B,A) resets B and returns False

Examples:

See also:

B = pinv( A )
B = pinv( A, tolerance )

pinv( B, A )
pinv( B, A, tolerance )

Moore-Penrose pseudo-inverse of matrix A

The computation is based on singular value decomposition

The tolerance argument is optional

The default tolerance is max(m,n)*max_sv*datum.eps, where:
- m = number of rows and n = number of columns in A
- max_sv = maximal singular value of A
- datum.eps = difference between 1 and the least value greater than 1 that is representable

Any singular values less than tolerance are treated as zero

If the decomposition fails:
- B = pinv(A) resets B
- pinv(B,A) resets B and returns False

Examples:

A = mat(4, 5, fill.randu)

B = pinv(A)        # use default tolerance

C = pinv(A, 0.01)  # set tolerance to 0.01

See also:

Q, R = qr( X )		(form 1)
qr( Q, R, X )

Q, R, P = qr( X, "vector" )		(form 2)
qr( Q, R, P, X, "vector" )

Q, R, P = qr( X, "matrix" )		(form 3)
qr( Q, R, P, X, "matrix" )

Decomposition of X into an orthogonal matrix Q and a right triangular matrix R, with an optional permutation matrix/vector P
- form 1: decomposition has the form Q*R = X
- form 2: P is permutation vector with class umat; decomposition has the form Q*R = X[:,P]
- form 3: P is permutation matrix with class umat; decomposition has the form Q*R = X*P

If P is specified, a column pivoting decomposition is used; the diagonal entries of R are ordered from largest to smallest magnitude

If the decomposition fails:
- Q, R = qr(X) resets Q & R
- Q, R, P = qr(X, "vector") and Q, R, P = qr(X, "matrix") reset Q, R & P
- qr(Q, R, X) resets Q & R and returns False
- qr(Q, R, P, X, "vector") and qr(Q, R, P, X, "matrix") reset Q, R & P and return False

Examples:

X = mat(5, 5, fill.randu)

Q = mat()
R = mat()

qr(Q, R, X)

P_vec = umat()
P_mat = umat()

qr(Q, R, P_vec, X, "vector")
qr(Q, R, P_mat, X, "matrix")

See also:

Q, R = qr_econ( X )
qr_econ( Q, R, X )

Economical decomposition of X (with size m x n) into an orthogonal matrix Q and a right triangular matrix R, such that Q*R = X

If m > n, only the first n rows of R and the first n columns of Q are calculated (i.e. the zero rows of R and the corresponding columns of Q are omitted)

If the decomposition fails:
- Q, R = qr_econ(X) resets Q & R
- qr_econ(Q,R,X) resets Q & R and returns False

Examples:

X = mat(6, 5, fill.randu)

Q = mat()
R = mat()

qr_econ(Q, R, X)

See also:

AA, BB, Q, Z = qz( A, B )
AA, BB, Q, Z = qz( A, B, select )

qz( AA, BB, Q, Z, A, B )
qz( AA, BB, Q, Z, A, B, select )

Generalised Schur decomposition for pair of general square matrices A and B of the same size,
such that A = Q.t()*AA*Z.t() and B = Q.t()*BB*Z.t()

The select argument is optional and specifies the ordering of the top left of the Schur form; it is one of the following:

`"none"`		no ordering (default operation)
`"lhp"`		left-half-plane: eigenvalues with real part < 0
`"rhp"`		right-half-plane: eigenvalues with real part > 0
`"iuc"`		inside-unit-circle: eigenvalues with absolute value < 1
`"ouc"`		outside-unit-circle: eigenvalues with absolute value > 1

The left and right Schur vectors are stored in Q and Z, respectively

In the complex-valued problem, the generalised eigenvalues are found in diagvec(AA) / diagvec(BB)

If A or B is not square sized, a RuntimeError is thrown

If the decomposition fails:
- AA, BB, Q, Z = qz(A,B) resets AA, BB, Q & Z
- qz(AA, BB, Q, Z, A, B) resets AA, BB, Q & Z and returns False

Examples:

A = mat(10, 10, fill.randu)
B = mat(10, 10, fill.randu)

AA = mat()
BB = mat()
Q = mat()
Z = mat() 

qz(AA, BB, Q, Z, A, B)

See also:

U, S = schur( X )

schur( S, X )

schur( U, S, X )

Schur decomposition of square matrix X, such that X = U*S*U.t()

U is a unitary matrix containing the Schur vectors

S is an upper triangular matrix, called the Schur form of X

If X is not square sized, a RuntimeError is thrown

If the decomposition fails:
- U, S = schur(X) resets U & S
- schur(S,X) resets S and returns False
- schur(U,S,X) resets U & S and returns False

Caveat: in general, Schur decomposition is not unique

Examples:

X = mat(20,20, fill.randu)

U = mat()
S = mat()

schur(U, S, X)

See also:

X = solve( A, B )
X = solve( A, B, settings )

solve( X, A, B )
solve( X, A, B, settings )

Solve a dense system of linear equations, A*X = B, where X is unknown; similar functionality to the \ operator in Matlab/Octave, i.e. X = A \ B

A can be square sized (critically determined system), or non-square (under/over-determined system)

B can be a vector or matrix

The number of rows in A and B must be the same

By default, matrix A is analysed to automatically determine whether it is a general matrix, band matrix, diagonal matrix, or symmetric/hermitian positive definite (SPD) matrix; based on the detected matrix structure, a specialised solver is used for faster execution

The settings argument is optional; it is one of the following, or a combination thereof:

`solve_opts.fast`		fast mode: disable determining solution quality via rcond, disable iterative refinement, disable equilibration
`solve_opts.refine`		apply iterative refinement to improve solution quality (matrix A must be square)
`solve_opts.equilibrate`		equilibrate the system before solving (matrix A must be square)
`solve_opts.likely_sympd`		indicate that matrix A is likely symmetric/hermitian positive definite
`solve_opts.allow_ugly`		keep solutions of systems that are singular to working precision
`solve_opts.no_approx`		do not find approximate solutions for rank deficient systems
`solve_opts.no_band`		do not use specialised solver for band matrices or diagonal matrices
`solve_opts.no_trimat`		do not use specialised solver for triangular matrices
`solve_opts.no_sympd`		do not use specialised solver for symmetric/hermitian positive definite matrices

the above settings can be combined using the + operator; for example: solve_opts.fast + solve_opts.no_approx

Caveat: using solve_opts.fast will speed up finding the solution, but for poorly conditioned systems the solution may have lower quality

Caveat: not all SPD matrices are automatically detected; to skip the analysis step and directly indicate that matrix A is likely SPD, use solve_opts.likely_sympd

If no solution is found:
- X = solve(A,B) resets X
- solve(X,A,B) resets X and returns False

Internal implementation details are available in the following paper:

Examples:

A = mat(5, 5, fill.randu)
b = mat(5, 1, fill.randu)
B = mat(5, 5, fill.randu)

x1 = solve(A, b)

x2 = mat()
status = solve(x2, A, b)

X1 = solve(A, B)

X2 = solve(A, B, solve_opts.fast)  # enable fast mode

See also:
- inv()
- pinv()
- rcond()
- roots()
- syl()
- linear system of equations in MathWorld
- system of linear equations in Wikipedia
- band matrix in Wikipedia
- definiteness of a matrix in Wikipedia
- positive definite matrix in MathWorld
- iterative refinement
- An Adaptive Solver for Systems of Linear Equations - paper detailing the internal implementation of solve()

mat U, mat s, mat V = svd( mat X )

cx_mat U, mat s, cx_mat V = svd( cx_mat X )

svd( mat s, X )

svd( mat U, mat s, mat V, mat X )

svd( cx_mat U, mat s, cx_mat V, cx_mat X )

Singular value decomposition of dense matrix X

s is a column vector (i.e. Nx1 matrix)

If X is square, it can be reconstructed using X = U*diagmat(s)*V.t()

The singular values are in descending order

If the decomposition fails, the output objects are reset and:
- U, s, V = svd(X) resets U, s & V
- svd(s,X) resets s and returns False
- svd(U,s,V,X) resets U, s, V and returns False

Examples:

X = mat(5, 5, fill.randu)

U = mat()
s = mat()
V = mat()

svd(U,s,V,X)

See also:

mat U, mat s, mat V = svd_econ( mat X )
mat U, mat s, mat V = svd_econ( mat X, mode )

cx_mat U, mat s, cx_mat V = svd_econ( cx_mat X )
cx_mat U, mat s, cx_mat V = svd_econ( cx_mat X, mode )

svd_econ( mat U, mat s, mat V, mat X )
svd_econ( mat U, mat s, mat V, mat X, mode )

svd_econ( cx_mat U, mat s, cx_mat V, cx_mat X )
svd_econ( cx_mat U, mat s, cx_mat V, cx_mat X, mode )

Economical singular value decomposition of dense matrix X

s is a column vector (i.e. Nx1 matrix)

The singular values are in descending order

The mode argument is optional; mode is one of:

If the decomposition fails:
- U, s, V = svd_econ(X) resets U, s & V
- svd_econ(U,s,v,X) resets U, s & V and returns False

Examples:

X = mat(4, 5, fill.randu)

U = mat()
s = mat()
V = mat()

svd_econ(U, s, V, X)

See also:

X = syl( A, B, C )
syl( X, A, B, C )

Solve the Sylvester equation, i.e. AX + XB + C = 0, where X is unknown

Matrices A, B and C must be square sized

If no solution is found:
- syl(A,B,C) resets X
- syl(X,A,B,C) resets X and returns False

Examples:

A = mat(5, 5, fill.randu)
B = mat(5, 5, fill.randu)
C = mat(5, 5, fill.randu)

X1 = syl(A, B, C)

X2 = mat()
syl(X2, A, B, C)

See also:
- solve()
- Sylvester equation in Wikipedia

Signal & Image Processing

conv( A, B )
conv( A, B, shape )

1D convolution of vectors A and B

The orientation of the result vector is the same as the orientation of A (i.e. either column or row vector)

The shape argument is optional; it is one of:

The convolution operation is also equivalent to FIR filtering

Examples:

A = mat(256, 1, fill.randu)

B = mat(16, 1, fill.randu)

C = conv(A, B)

D = conv(A, B, "same")

See also:

conv2( A, B )
conv2( A, B, shape )

2D convolution of matrices A and B

The shape argument is optional; it is one of:

The implementation of 2D convolution in this version is preliminary; it is not yet fully optimised

Examples:

A = mat(256, 256, fill.randu)

B = mat(16, 16, fill.randu)

C = conv2(A, B)

D = conv2(A, B, "same")

See also:

cx_mat Y = fft( X )
cx_mat Y = fft( X, n )

cx_mat Z = ifft( cx_mat Y )
cx_mat Z = ifft( cx_mat Y, n )

fft(): fast Fourier transform of a vector or matrix (real or complex)

ifft(): inverse fast Fourier transform of a vector or matrix (complex only)

If given a matrix, the transform is done on each column vector of the matrix

The optional n argument specifies the transform length:
- if n is larger than the length of the input vector, a zero-padded version of the vector is used
- if n is smaller than the length of the input vector, only the first n elements of the vector are used

If n is not specified, the transform length is the same as the length of the input vector

Caveat: the transform is fastest when the transform length is a power of 2, e.g. 64, 128, 256, 512, 1024, ...

The implementation of the transform in this version is preliminary; it is not yet fully optimised

Examples:

X = mat(100, 1, fill.randu)
   
Y = fft(X, 128)

See also:

cx_mat Y = fft2( X )
cx_mat Y = fft2( X, n_rows, n_cols )

cx_mat Z = ifft2( cx_mat Y )
cx_mat Z = ifft2( cx_mat Y, n_rows, n_cols )

fft2(): 2D fast Fourier transform of a matrix (real or complex)

ifft2(): 2D inverse fast Fourier transform of a matrix (complex only)

The optional arguments n_rows and n_cols specify the size of the transform; a truncated and/or zero-padded version of the input matrix is used

Caveat: the transform is fastest when both n_rows and n_cols are a power of 2, e.g. 64, 128, 256, 512, 1024, ...

The implementation of the transform in this version is preliminary; it is not yet fully optimised

Examples:

A = mat(100, 100, fill.randu)
   
B = fft2(A)
C = fft2(A, 128, 128)

See also:

interp1( X, Y, XI, YI )
interp1( X, Y, XI, YI, method )
interp1( X, Y, XI, YI, method, extrapolation_value )

1D data interpolation

Given a 1D function specified in vectors X and Y (where X specifies locations and Y specifies the corresponding values),
generate vector YI which contains interpolated values at locations XI

The method argument is optional; it is one of:

`"nearest"`	=	interpolate using single nearest neighbour
`"linear"`	=	linear interpolation between two nearest neighbours (default setting)
`"*nearest"`	=	as per `"nearest"`, but faster by assuming that X and XI are monotonically increasing
`"*linear"`	=	as per `"linear"`, but faster by assuming that X and XI are monotonically increasing

If a location in XI is outside the domain of X, the corresponding value in YI is set to extrapolation_value

The extrapolation_value argument is optional; by default it is datum.nan (not-a-number)

Examples:

x = linspace(0, 3, 20)
y = square(x)

xx = linspace(0, 3, 100)

yy = mat()

interp1(x, y, xx, yy)  # use linear interpolation by default

interp1(x, y, xx, yy, "*linear")  # faster than "linear"

interp1(x, y, xx, yy, "nearest")

See also:

interp2( X, Y, Z, XI, YI, ZI )
interp2( X, Y, Z, XI, YI, ZI, method )
interp2( X, Y, Z, XI, YI, ZI, method, extrapolation_value )

2D data interpolation

Given a 2D function specified by matrix Z with coordinates given by vectors X and Y,
generate matrix ZI which contains interpolated values at the coordinates given by vectors XI and YI

The vector pairs (X, Y) and (XI, YI) define 2D coordinates in a grid;
for example, X defines the horizontal coordinates and Y defines the corresponding vertical coordinates,
so that ( X(m), Y(n) ) is the 2D coordinate of element Z(n,m)

The length of vector X must be equal to the number of columns in matrix Z

The length of vector Y must be equal to the number of rows in matrix Z

Vectors X, Y, XI, YI must contain monotonically increasing values (e.g. 0.1, 0.2, 0.3, ...)

The method argument is optional; it is one of:

If a coordinate in the 2D grid specified by (XI, YI) is outside the domain of the 2D grid specified by (X, Y),
the corresponding value in ZI is set to extrapolation_value

The extrapolation_value argument is optional; by default it is datum.nan (not-a-number)

Examples:

Z = mat()

Z.load("input_image.pgm", pgm_binary)  # load an image in pgm format

X = regspace(1, Z.n_cols)  # X = horizontal spacing
Y = regspace(1, Z.n_rows)  # Y = vertical spacing

XI = regspace(X.min(), 1.0/2.0, X.max()) # magnify by approx 2
YI = regspace(Y.min(), 1.0/3.0, Y.max()) # magnify by approx 3

ZI = mat()

interp2(X, Y, Z, XI, YI, ZI)  # use linear interpolation by default

ZI.save("output_image.pgm", pgm_binary)

See also:

P = polyfit( X, Y, N )
polyfit( P, X, Y, N )

Given a 1D function specified in vectors X and Y (where X holds independent values and Y holds the corresponding dependent values),
model the function as a polynomial of order N and store the polynomial coefficients in column vector P

The given function is modelled as:
where p_{_i} is the i-th polynomial coefficient; the coefficients are selected to minimise the overall error of the fit (least squares)

The column vector P has N+1 coefficients

N must be smaller than the number of elements in X

If the polynomial coefficients cannot be found:
- P = polyfit( X, Y, N ) resets P
- polyfit( P, X, Y, N ) resets P and returns False

Examples:

x = linspace(0,4*datum.pi,100)
y = cos(x)

p = polyfit(x,y,10)

See also:

Y = polyval( P, X )

Given vector P of polynomial coefficients and vector X containing the independent values of a 1D function,
generate vector Y which contains the corresponding dependent values

For each x value in vector X, the corresponding y value in vector Y is generated using:
where p_{_i} is the i-th polynomial coefficient in vector P

P must contain polynomial coefficients in descending powers (e.g. generated by the polyfit() function)

Examples:

x1 = linspace(0,4*datum.pi,100)
y1 = cos(x1)
p1 = polyfit(x1,y1,10)

y2 = polyval(p1,x1)

See also:

Statistics & Clustering

mean, median, stddev, var, spread

mean( M ) mean( M, dim ) mean( Q ) mean( Q, dim )	(form 1) (form 2) (form 1) (form 2)	⎫ ⎪ ⎬ mean (average value) ⎪ ⎭
median( M ) median( M, dim )	(form 1) (form 2)	⎫ ⎪ ⎬ median ⎪ ⎭
stddev( M ) stddev( M, norm_type ) stddev( M, norm_type, dim )	(form 1a) (form 1b) (form 2)	⎫ ⎪ ⎬ standard deviation ⎪ ⎭
var( M ) var( M, norm_type ) var( M, norm_type, dim )	(form 1a) (form 1b) (form 2)	⎫ ⎪ ⎬ variance ⎪ ⎭
spread( M ) spread( M, dim )	(form 1) (form 2)	⎫ ⎪ ⎬ spread (difference between max and min) ⎪ ⎭

Form 1:
- if matrix M can be interpreted as a vector, return the statistic calculated using all the elements of the vector as a 1x1 matrix
- otherwise, form 2 is used with dim = 0

Form 2:
- for matrix M, find the statistic for each column (dim = 0), or each row (dim = 1)
- for cube Q, find the statistics of elements along dimension dim, where dim ∈ { 0, 1, 2 }

The norm_type argument is optional; by default norm_type = 0 is used

For the var() and stddev() functions:
- the default norm_type = 0 performs normalisation using N-1 (where N is the number of samples), providing the best unbiased estimator
- using norm_type = 1 performs normalisation using N, which provides the second moment around the mean

Caveat: to obtain statistics for integer matrices, convert to a matrix/vector with floating point values using the mat() or fmat() constructors

Examples:

A = mat(5, 5, fill.randu)

B = mean(A)
C = var(A)
m = mean(mean(A))

v = mat(5, 1, fill.randu)
x = var(v)

See also:
- cov()
- cor()
- diff()
- min() & max()

cov( X, Y )
cov( X, Y, norm_type )

cov( X )
cov( X, norm_type )

For two matrix arguments X and Y, if each row of X and Y is an observation and each column is a variable, the (i,j)-th entry of cov(X,Y) is the covariance between the i-th variable in X and the j-th variable in Y

For vector arguments, the type of vector is ignored and each element in the vector is treated as an observation

For matrices, X and Y must have the same dimensions

For vectors, X and Y must have the same number of elements

cov(X) is equivalent to cov(X, X)

The default norm_type = 0 performs normalisation using N-1 (where N is the number of observations), providing the best unbiased estimation of the covariance matrix (if the observations are from a normal distribution). Using norm_type = 1 causes normalisation to be done using N, which provides the second moment matrix of the observations about their mean

Examples:

X = mat(4, 5, fill.randu)
Y = mat(4, 5, fill.randu)

C = cov(X,Y)
D = cov(X,Y, 1)

See also:

cor( X, Y )
cor( X, Y, norm_type )

cor( X )
cor( X, norm_type )

For two matrix arguments X and Y, if each row of X and Y is an observation and each column is a variable, the (i,j)-th entry of cor(X,Y) is the correlation coefficient between the i-th variable in X and the j-th variable in Y

For vector arguments, the type of vector is ignored and each element in the vector is treated as an observation

For matrices, X and Y must have the same dimensions

For vectors, X and Y must have the same number of elements

cor(X) is equivalent to cor(X, X)

The default norm_type = 0 performs normalisation of the correlation matrix using N-1 (where N is the number of observations). Using norm_type = 1 causes normalisation to be done using N

Examples:

X = mat(4, 5, fill.randu)
Y = mat(4, 5, fill.randu)

R = cor(X,Y)
S = cor(X,Y, 1)

See also:

hist( X )		(form 1a)
hist( X, n_bins )		(form 1b)
hist( X, centers )		(form 1c)
hist( X, centers, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, produce an unsigned vector of the same orientation as X (i.e. umat) that represents a histogram of counts
- otherwise, form 2 is used with dim = 0 and the centers vector must be specified

For form 2, produce a umat matrix containing either column histogram counts (for dim = 0), or row histogram counts (for dim = 1)

The bin centers can be automatically determined from the data, with the number of bins specified via n_bins (default is 10); the range of the bins is determined by the range of the data

The bin centers can also be explicitly specified via the centers vector; the vector must contain monotonically increasing values (e.g. 0.1, 0.2, 0.3, ...)

Examples:

v = mat(1000, 1, fill.randn) # Gaussian distribution

h1 = hist(v, 11)
h2 = hist(v, linspace(-2,2,11))

See also:

histc( X, edges )		(form 1)
histc( X, edges, dim )		(form 2)

For form 1:
- if matrix X can be interpreted as a vector, produce an unsigned vector of the same orientation as X (i.e. umat) that contains the counts of the number of values that fall between the elements in the edges vector
- otherwise, form 2 is used with dim = 0

For form 2, produce a umat matrix containing either column histogram counts (for dim = 0), or row histogram counts (for dim = 1)

The edges vector must contain monotonically increasing values (e.g. 0.1, 0.2, 0.3, ...)

Examples:

v = mat(1000, 1, fill.randn)  # Gaussian distribution

h = histc(v, linspace(-2,2,11))

See also:

quantile( X, P )		(form 1)
quantile( X, P, dim )		(form 2)

For a dataset stored in vector V or matrix X, calculate the quantiles corresponding to the cumulative probability values in the given P vector

For form 1:
- if matrix X can be interpreted as a vector, produce a vector with the same orientation as X and the same length as P
- otherwise, form 2 is used with dim = 0

For form 2, produce a matrix with the quantiles for each column vector (dim = 0) or each row vector (dim = 0)

The P vector must contain values in the [0,1] interval (e.g. 0.00, 0.25, 0.50, 0.75, 1.00)

The algorithm for calculating the quantiles is based on Definition 5 in:
Rob J. Hyndman and Yanan Fan. Sample Quantiles in Statistical Packages. The American Statistician, Vol. 50, No. 4, pp. 361-365, 1996. http://doi.org/10.2307/2684934

Examples:

V = mat(1000, 1, fill.randn)  # Gaussian distribution

P = mat([ 0.25, 0.50, 0.75 ])

Q = quantile(V, P)

See also:

mat coeff, mat score, mat latent, mat tsquared = princomp( mat X )
cx_mat coeff, cx_mat score, mat latent, cx_mat tsquared = princomp( cx_mat X )

princomp( mat coeff, mat X )
princomp( cx_mat coeff, cx_mat X )

princomp( mat coeff, mat score, mat X )
princomp( cx_mat coeff, cx_mat score, cx_mat X )

princomp( mat coeff, mat score, vec latent, mat X )
princomp( cx_mat coeff, cx_mat score, vec latent, cx_mat X )

princomp( mat coeff, mat score, vec latent, vec tsquared, mat X )
princomp( cx_mat coeff, cx_mat score, vec latent, cx_vec tsquared, cx_mat X )

Principal component analysis of matrix X

Each row of X is an observation and each column is a variable

output objects:
- coeff: principal component coefficients
- score: projected data
- latent: vector containing eigenvalues of the covariance matrix of X
- tsquared: vector containing Hotteling's statistic for each sample

The computation is based on singular value decomposition

If the decomposition fails:
- coeff, score, latent, tsquared = princomp(X) resets coeff, score, latent, and tsquared
- remaining forms of princomp() reset all output matrices and return False

Examples:

A = mat(5, 4, fill.randu)

coeff = mat()
score = mat()
latent = mat()
tsquared = mat()

princomp(coeff, score, latent, tsquared, A)

See also:

normpdf( X )
normpdf( X, M, S )

For each scalar x in X, compute its probability density function according to a Gaussian (normal) distribution using the corresponding mean value m in M and the corresponding standard deviation value s in S:

X can be a scalar, vector, or matrix

M and S can jointly be either scalars, vectors, or matrices

If M and S are omitted, their values are assumed to be 0 and 1, respectively

Caveat: to reduce the incidence of numerical underflows, consider using log_normpdf()

Examples:

X = mat(10, 1, fill.randu)
M = mat(10, 1, fill.randu)
S = mat(10, 1, fill.randu)

P1 = normpdf(X)
P2 = normpdf(X,    M,    S   )
P3 = normpdf(1.23, M,    S   )
P4 = normpdf(X,    4.56, 7.89)
P5 = normpdf(1.23, 4.56, 7.89)

See also:

log_normpdf( X )
log_normpdf( X, M, S )

For each scalar x in X, compute the logarithm version of probability density function according to a Gaussian (normal) distribution using the corresponding mean value m in M and the corresponding standard deviation value s in S:



  
    
          ⎧    1   ⎧   (x-m)² ⎫ ⎫
    
    
      y =  log ⎪ ‒‒‒‒‒‒‒  exp ⎪ −0.5  ‒‒‒‒‒‒  ⎪ ⎪
    
    
          ⎩ s √(2π)   ⎩     s² ⎭ ⎭
    
  




  
    
            ⎧   (x-m)² ⎫
    
    
        =  −log(s √(2π))  +  ⎪ −0.5  ‒‒‒‒‒‒  ⎪
    
    
            ⎩     s² ⎭

X can be a scalar, vector, or matrix

M and S can jointly be either scalars, vectors, or matrices

If M and S are omitted, their values are assumed to be 0 and 1, respectively

Examples:

X = mat(10, 1, fill.randu)
M = mat(10, 1, fill.randu)
S = mat(10, 1, fill.randu)

P1 = log_normpdf(X)
P2 = log_normpdf(X,    M,    S   )
P3 = log_normpdf(1.23, M,    S   )
P4 = log_normpdf(X,    4.56, 7.89)
P5 = log_normpdf(1.23, 4.56, 7.89)

See also:
- normpdf()
- normal distribution in Wikipedia

normcdf( X )
normcdf( X, M, S )

For each scalar x in X, compute its cumulative distribution function according to a Gaussian (normal) distribution using the corresponding mean value m in M and the corresponding standard deviation value s in S

X can be a scalar, vector, or matrix

M and S can jointly be either scalars, vectors, or matrices

If M and S are omitted, their values are assumed to be 0 and 1, respectively

Examples:

X = mat(10, 1, fill.randu)
M = mat(10, 1, fill.randu)
S = mat(10, 1, fill.randu)

P1 = normcdf(X)
P2 = normcdf(X,    M,    S   )
P3 = normcdf(1.23, M,    S   )
P4 = normcdf(X,    4.56, 7.89)
P5 = normcdf(1.23, 4.56, 7.89)

See also:

X = mvnrnd( M, C )
X = mvnrnd( M, C, N )

mvnrnd( X, M, C )
mvnrnd( X, M, C, N )

Generate a matrix with random column vectors from a multivariate Gaussian (normal) distribution with parameters M and C:
- M is the mean; must be a matrix that can be interpreted as a column vector
- C is the covariance matrix; must be symmetric positive semi-definite (preferably positive definite)

N is the number of column vectors to generate; if N is omitted, it is assumed to be 1

If generating the random vectors fails:
- X = mvnrnd(M, C) and X = mvnrnd(M, C, N) reset X
- mvnrnd(X, M, C) and mvnrnd(X, M, C, N) reset X and return False

Examples:

M = mat(5, 1, fill.randu)

B = mat(5, 5, fill.randu)
C = B.t() * B

X = mvnrnd(M, C, 100)

See also:

mat X = chi2rnd( DF )
double s = chi2rnd( DF_scalar )
mat v = chi2rnd( DF_scalar, n_elem )
mat X = chi2rnd( DF_scalar, n_rows, n_cols )
mat Y = chi2rnd( DF_scalar, size(X) )

Generate a random scalar, vector or matrix with elements sampled from a chi-squared distribution with the degrees of freedom specified by parameter DF or DF_scalar

DF is a real vector or matrix, while DF_scalar is a scalar of type double

Each value in DF and DF_scalar must be greater than zero

For the chi2rnd(DF) form, the output vector/matrix has the same size and type as DF; each element in DF specifies a separate degree of freedom

For forms using DF_scalar, the output vector/matrix is of type mat

Examples:

from random import randint
X = chi2rnd(2, 5, 6)

A = mat(5, 6)
A.imbue(lambda: randint(1,10)) # imbue with random integers
B = chi2rnd(A)

See also:

W = wishrnd( S, df )
W = wishrnd( S, df, D )

wishrnd( W, S, df )
wishrnd( W, S, df, D )

Generate a random matrix sampled from the Wishart distribution with parameters S and df:
- S is a symmetric positive definite matrix (e.g. a covariance matrix)
- df is a scalar specifying the degrees of freedom; it can be a non-integer value

D is an optional argument; it specifies the Cholesky decomposition of S; if D is provided, S is ignored;
using D is more efficient if you need to use wishrnd() many times for the same S matrix

If generating the random matrix fails:
- W = wishrnd(S, df) and W = wishrnd(S, df, D) reset W
- wishrnd(W, S, df) and wishrnd(W, S, df, D) reset W and return False

Examples:

X = mat(5, 5, fill.randu)

S = X.t() * X

W = wishrnd(S, 6.7)

See also:

W = iwishrnd( T, df )
W = iwishrnd( T, df, Dinv )

iwishrnd( W, T, df )
iwishrnd( W, T, df, Dinv )

Generate a random matrix sampled from the inverse Wishart distribution with parameters T and df:
- T is a symmetric positive definite matrix
- df is a scalar specifying the degrees of freedom; it can be a non-integer value

Dinv is an optional argument; it specifies the Cholesky decomposition of the inverse of T; if Dinv is provided, T is ignored
using Dinv is more efficient if you need to use iwishrnd() many times for the same T matrix

If generating the random matrix fails:
- W = iwishrnd(T, df) and W = iwishrnd(T, df, Dinv) reset W
- iwishrnd(W, T, df) and iwishrnd(W, T, df, Dinv) reset W and return False

Examples:

X = mat(5, 5, fill.randu)

T = X.t() * X

W = iwishrnd(T, 6.7)

See also:

running_stat		(double-precision floating point)
frunning_stat		(single-precision floating point)
cx_running_stat		(complex double-precision floating point)
cx_frunning_stat		(complex single-precision floating point)

Classes for running statistics (online statistics) of scalars (one dimensional process/signal)

Useful if the storage of all samples (scalars) is impractical, or if the number of samples is not known in advance

In this documentation the running_stat class is used for convenience; it is possible to use other classes instead, e.g. cx_running_stat

For an instance of running_stat named as X, the member functions are:

X(scalar)		update the statistics using the given scalar
X.min()		current minimum value
X.max()		current maximum value
X.spread()		current range
X.mean()		current mean or average value
X.var() and X.var(norm_type)		current variance
X.stddev() and X.stddev(norm_type)		current standard deviation
X.reset()		reset all statistics and set the number of samples to zero
X.count()		current number of samples

The norm_type argument is optional; by default norm_type = 0 is used

For the .var() and .stddev() functions, the default norm_type = 0 performs normalisation using N-1 (where N is the number of samples so far), providing the best unbiased estimator; using norm_type = 1 causes normalisation to be done using N, which provides the second moment around the mean

The return type of .count() depends on the underlying form of type: it is either float or double

Examples:

from random import normalvariate
stats = running_stat()

for i in range(10000):
  sample = normalvariate(0, 1) # normal distribution
  stats(sample)

print("mean = " + str(stats.mean()))
print("var  = " + str(stats.var()))
print("min  = " + str(stats.min()))
print("max  = " + str(stats.max()))

See also:
- running_stat_vec (running statistics of vectors)
- statistics functions

running_stat_vec		(double-precision floating point)
frunning_stat_vec		(single-precision floating point)
cx_running_stat_vec		(complex double-precision floating point)
cx_frunning_stat_vec		(complex single-precision floating point)

running_stat_vec(calc_cov)		(double-precision floating point)
frunning_stat_vec(calc_cov)		(single-precision floating point)
cx_running_stat_vec(calc_cov)		(complex double-precision floating point)
cx_frunning_stat_vec(calc_cov)		(complex single-precision floating point)

Class for running statistics (online statistics) of vectors (multi-dimensional process/signal)

Useful if the storage of all samples (vectors) is impractical, or if the number of samples is not known in advance

This class is similar to running_stat, with the difference that vectors are processed instead of scalars

In this documentation the running_stat_vec class is used for convenience; it is possible to use other classes instead, e.g. cx_running_stat_vec

For an instance of running_stat_vec named as X, the member functions are:

X(vector)		update the statistics using the given vector
X.min()		vector of current minimum values
X.max()		vector of current maximum values
X.spread()		vector of current ranges
X.mean()		vector of current means
X.var() and X.var(norm_type)		vector of current variances
X.stddev() and X.stddev(norm_type)		vector of current standard deviations
X.cov() and X.cov(norm_type)		matrix of current covariances; valid if calc_cov=true during construction of running_stat_vec
X.reset()		reset all statistics and set the number of samples to zero
X.count()		current number of samples

The calc_cov argument is optional; by default calc_cov=False, indicating that the covariance matrix will not be calculated; to enable the covariance matrix, use calc_cov=True during construction; for example: X = running_stat_vec(True)

The norm_type argument is optional; by default norm_type = 0 is used

For the .var() and .stddev() functions, the default norm_type = 0 performs normalisation using N-1 (where N is the number of samples so far), providing the best unbiased estimator; using norm_type = 1 causes normalisation to be done using N, which provides the second moment around the mean

The return type of .count() depends on the underlying form of vec_type: it is either float or double

Examples:

stats = running_stat_vec()

for i in range(10000):
  sample = mat(5, 1, fill.randu)
  stats(sample)

stats.mean().print("mean = ")
stats.var().print("var  = ")
stats.max().print("max  = ")

more_stats = running_stat_vec(True)

for i in range(20):
  sample = mat(1, 3, fill.randu)
  
  sample[1] -= sample[0]
  sample[2] += sample[1]
  
  more_stats(sample)

more_stats.cov().print("covariance matrix = ")

sd = more_stats.stddev()

(more_stats.cov() / (sd.t() * sd)).print("correlations = ")

See also:
- running_stat (running statistics of scalars)
- statistics functions
- cov()
- cor()

kmeans( means, data, k, seed_mode, n_iter, print_mode )

Cluster given data into k disjoint sets

The means parameter is the output matrix for storing the resulting centroids of the sets, with each centroid stored as a column vector

The data parameter is the input data matrix, with each sample stored as a column vector

The k parameter indicates the number of centroids; the number of samples in the data matrix should be much larger than k

The seed_mode parameter specifies how the initial centroids are seeded; it is one of:

`keep_existing`		use the centroids specified in the means matrix as the starting point
`static_subset`		use a subset of the data vectors (repeatable)
`random_subset`		use a subset of the data vectors (random)
`static_spread`		use a maximally spread subset of data vectors (repeatable)
`random_spread`		use a maximally spread subset of data vectors (random start)

caveat:

static_spread

random_spread

static_subset

random_subset

The n_iter parameter specifies the number of clustering iterations; this is data dependent, but about 10 is typically sufficient

The print_mode parameter is either True or False, indicating whether progress is printed during clustering

If the clustering fails, the means matrix is reset and False is returned

Examples:

d = 5       # dimensionality
N = 10000   # number of vectors

data = mat(d, N, fill.randu)

means = mat()

status = kmeans(means, data, 2, random_subset, 10, True)

if not status:
  print("clustering failed")

means.print("means:")

See also:

Miscellaneous

constants (pi, inf, speed of light, ...)

`datum.pi`		π, the ratio of any circle's circumference to its diameter
`datum.inf`		∞, infinity
`datum.nan`		“not a number” (NaN) caveat: NaN is not equal to anything, even itself

`datum.e`		base of the natural logarithm
`datum.sqrt2`		square root of 2
`datum.sqrt2pi`		square root of 2π
`datum.eps`		machine epsilon: the difference between 1 and the value least greater than 1 that is representable (type and machine dependent)

`datum.log_min`		log of minimum non-zero value (type and machine dependent)
`datum.log_max`		log of maximum value (type and machine dependent)
`datum.euler`		Euler's constant, aka Euler-Mascheroni constant

`datum.gratio`		golden ratio
`datum.m_u`		atomic mass constant (in kg)
`datum.N_A`		Avogadro constant

`datum.k`		Boltzmann constant (in joules per kelvin)
`datum.k_evk`		Boltzmann constant (in eV/K)
`datum.a_0`		Bohr radius (in meters)

`datum.mu_B`		Bohr magneton
`datum.Z_0`		characteristic impedance of vacuum (in ohms)
`datum.G_0`		conductance quantum (in siemens)

`datum.k_e`		Coulomb's constant (in meters per farad)
`datum.eps_0`		electric constant (in farads per meter)
`datum.m_e`		electron mass (in kg)

`datum.eV`		electron volt (in joules)
`datum.ec`		elementary charge (in coulombs)
`datum.F`		Faraday constant (in coulombs)

`datum.alpha`		fine-structure constant
`datum.alpha_inv`		inverse fine-structure constant
`datum.K_J`		Josephson constant

`datum.mu_0`		magnetic constant (in henries per meter)
`datum.phi_0`		magnetic flux quantum (in webers)
`datum.R`		molar gas constant (in joules per mole kelvin)

`datum.G`		Newtonian constant of gravitation (in newton square meters per kilogram squared)
`datum.h`		Planck constant (in joule seconds)
`datum.h_bar`		Planck constant over 2 pi, aka reduced Planck constant (in joule seconds)

`datum.m_p`		proton mass (in kg)
`datum.R_inf`		Rydberg constant (in reciprocal meters)
`datum.c_0`		speed of light in vacuum (in meters per second)

`datum.sigma`		Stefan-Boltzmann constant
`datum.R_k`		von Klitzing constant (in ohms)
`datum.b`		Wien wavelength displacement law constant

The constants are stored in the datum class

Caveat: datum.nan is not equal to anything, even itself; to check whether a scalar x is finite, use math.isfinite(x)

The physical constants were mainly taken from NIST 2018 CODATA values, and some from WolframAlpha (as of 2009-06-23)

Examples:

print("2.0 * pi = " + str(2.0 * datum.pi))

print("speed of light = " + str(datum.c_0))

print("log_max for doubles = " + str(datum.log_max))

See also:
- .is_finite()
- .fill()
- NaN in Wikipedia
- physical constant in Wikipedia

wall_clock

Simple wall clock timer class for measuring the number of elapsed seconds

Examples:

A = mat(1000, 1000, fill.randu)
B = mat(1000, 1000, fill.randu)

timer = wall_clock()

timer.tic()
C = A * B
n = timer.toc()

print("elapsed: " + str(n))

libraries

Prints the libraries used by PyArmadillo and their location, as well as the libraries not used

Examples:

libraries()

# possible output:
#
#
# PyArmadillo library status on install:
#
# MKL is not used.
# OpenBLAS is used and is located at /usr/lib/x86_64-linux-gnu/libopenblas.so
# ATLAS is not used.
# BLAS is not used.
# LAPACK is used and is located at /usr/lib/x86_64-linux-gnu/liblapack.so
# FlexiBLAS is not used.
# HDF5 is used and is located at /usr/lib/x86_64-linux-gnu/hdf5/serial/libhdf5.so
# ARPACK is used and is located at /usr/lib/x86_64-linux-gnu/libarpack.so
# SuperLU is used and is located at /usr/lib/x86_64-linux-gnu/libsuperlu.so
# Accelerate is not used.

Examples of Matlab/Octave syntax and conceptually corresponding PyArmadillo syntax

Matlab/Octave	PyArmadillo	Notes

`A(1, 1)`	`A[0, 0]`	indexing in PyArmadillo starts at 0
`A(k, k)`	`A[k-1, k-1]`

`size(A,1)`	`A.n_rows`	read only
`size(A,2)`	`A.n_cols`
`size(Q,3)`	`A.n_slices`	Q is a cube (3D array)
`numel(A)`	`A.n_elem`

`A(:, k)`	`A[:, k]`	this is a conceptual example only; exact conversion from Matlab/Octave to PyArmadillo syntax will require taking into account that indexing starts at 0
`A(k, :)`	`A[k, :]`
`A(:, p:q)`	`A[:, p:q]`
`A(p:q, :)`	`A[p:q, :]`
`A(p:q, r:s)`	`A[p:q, r:s]`

`Q(p:q, r:s, t:u)`	`Q[p:q, r:s, t:u]`	Q is a cube (3D array)
`Q(:, :, t:u)`	`Q[:, :, t:u]`
`Q(:, :, k)`	`Q[single_slice, k]`

`A'`	`A.t() or trans(A)`	matrix transpose / Hermitian transpose (for complex matrices, the conjugate of each element is taken)

`A = zeros(size(A))`	`A.zeros()`
`A = ones(size(A))`	`A.ones()`
`A = zeros(k)`	`A = zeros(k,k)`
`A = ones(k)`	`A = ones(k,k)`

`C = complex(A,B)`	`C = cx_mat(A,B)`

`A * B`	`A * B`	matrix multiplication
`A .* B`	`A @ B`	element-wise multiplication
`A ./ B`	`A / B`	element-wise division
`A \ B`	`solve(A,B)`	conceptually similar to inv(A)*B, but more efficient
`A = A + 1;`	`A += 1`
`A = A - 1;`	`A -= 1`

`A = [ 1 2; 3 4; ]`	`A = mat([ [ 1, 2 ], [ 3, 4 ] ])`	element initialisation

`X = A(:)`	`X = vectorise(A)`
`X = [ A B ]`	`X = join_horiz(A,B)`
`X = [ A; B ]`	`X = join_vert(A,B)`

`A`	`A.print("A:")`

`save ‑ascii 'A.txt' A`	`A.save("A.txt", raw_ascii)`	Matlab/Octave matrices saved as ascii are readable by PyArmadillo (and vice-versa)
`load ‑ascii 'A.txt'`	`A.load("A.txt", raw_ascii)`

example program

Form 1:

from pyarma import *

A = mat(4, 5, fill.randu)
B = mat(4, 5, fill.randu)
  
C = A*B.t()

C.print("C:")

Form 2:

import pyarma as pa

A = pa.mat(4, 5, pa.fill.randu)
B = pa.mat(4, 5, pa.fill.randu)
  
C = A*B.t()

C.print("C:")

Both forms are functionally identical

See also the example program that comes with the PyArmadillo archive

History of API Additions, Changes and Deprecations

API Stability and Versioning:
- Each release of PyArmadillo has its public API (functions, classes, constants) described in the accompanying API documentation specific to that release.
- Each release of PyArmadillo has its full version specified as A.B.C, where A is a major version number, B is a minor version number, and C is a patch level (indicating bug fixes).
- Within a major version (e.g. 1), each minor version has a public API that strongly strives to be backwards compatible (at the source level) with the public API of preceding minor versions. For example, user code written for version 1.100 should work with version 1.200, 1.300, 1.400, etc. However, as later minor versions may have more features (API extensions) than preceding minor versions, user code specifically written for version 1.400 may not work with 1.300.
- An increase in the patch level, while the major and minor versions are retained, indicates modifications to the code and/or documentation which aim to fix bugs without altering the public API.
- We don't like changes to existing public API and strongly prefer not to break any user software. However, to allow evolution, we reserve the right to alter the public API in future major versions of PyArmadillo while remaining backwards compatible in as many cases as possible (e.g. major version 2 may have slightly different public API than major version 1).
- Caveat: any function, class, constant or other code not explicitly described in the public API documentation is considered as part of the underlying internal implementation details, and may be removed or changed without notice. (In other words, don't use internal functionality).

List of additions and changes for each version:
- Version 0.500:
  - instances of mat and cube are initialised to contain zero-valued elements by default
  - added standalone zeros(), ones(), randu(), randn(), eye()
  - added pyarma_rng.set_seed(value) and pyarma_rng.set_seed_random()
  - added extra forms for lu(), qr(), qr_econ(), qz(), svd_econ()
  - added subscripting for size objects
  - range() renamed to spread() to prevent conflicts with built-in range() in Python
  - for solve(), solve_opts_* flags renamed to solve_opts.* (eg. solve_opts_fast is now solve_opts.fast)
  - for mat and cube constructors, fill_* flags renamed to fill.* (eg. fill_zeros is now fill.zeros)
- Version 0.400:
  - initial public release, adapting dense matrices and cubes from Armadillo
- Versions 0.100 - 0.300:
  - internal development releases

`True`	=	the returned indices are sorted to be ascending (default setting)
`False`	=	the returned indices are in arbitrary order (faster operation)

`"balance"`	↦	diagonally scale and permute X to improve conditioning of the eigenvalues
`"nobalance"`	↦	do not balance X; this is the default operation

n_rows	=	num_copies_per_row	*	A.n_rows
n_cols	=	num_copies_per_col	*	A.n_cols

`"both"`	=	compute both left and right singular vectors (default operation)
`"left"`	=	compute only left singular vectors
`"right"`	=	compute only right singular vectors

`"full"`	=	return the full convolution (default setting), with the size equal to A.n_elem + B.n_elem - 1
`"same"`	=	return the central part of the convolution, with the same size as vector A

`"nearest"`	=	interpolate using nearest neighbours
`"linear"`	=	linear interpolation between nearest neighbours (default setting)

	1		⎧		(x-m)²	⎫
y =	‒‒‒‒‒‒‒	exp	⎪	−0.5	‒‒‒‒‒‒	⎪
	s √(2π)		⎩		s²	⎭

			⎧		(x-m)²	⎫
=	−log(s √(2π))	+	⎪	−0.5	‒‒‒‒‒‒	⎪
			⎩		s²	⎭

		x^a-1 exp( -x / b )
p(x\|a,b)	=
		b^a Γ(a)