Generic matrices

# Generic matrices

AbstractAlgebra.jl allows the creation of dense matrices over any computable commutative ring $R$. Generic matrices over a commutative ring are implemented in src/generic/Matrix.jl.

As well as implementing the entire Matrix interface, including the optional functionality, there are many additional generic algorithms implemented for matrix spaces. We describe this functionality below.

All of this generic functionality is part of the Generic submodule of AbstractAlgebra.jl. This is exported by default, so it is not necessary to qualify names of functions.

## Types and parent objects

Generic matrices in AbstractAlgebra.jl have type Generic.Mat{T} where T is the type of elements of the matrix. Internally, generic matrices are implemented using an object wrapping a Julia two dimensional array, though they are not themselves Julia arrays. See the file src/generic/GenericTypes.jl for details.

Parents of generic matrices (matrix spaces) have type Generic.MatSpace{T}.

The generic matrix types belong to the abstract type AbstractAlgebra.MatElem{T} and the matrix space parent types belong to AbstractAlgebra.MatSpace{T}. Note that both the concrete type of a matrix space parent object and the abstract class it belongs to have the name MatElem, therefore disambiguation is required to specify which is intended.

The dimensions and base ring $R$ of a generic matrix are stored in its parent object, however to allow creation of matrices without first creating the matrix space parent, generic matrices in Julia do not contain a reference to their parent. They contain the row and column numbers and the base ring on a per matrix basis. The parent object can then be reconstructed from this data on demand.

## Matrix space constructors

A matrix space in AbstractAlgebra.jl represents a collection of all matrices with given dimensions and base ring.

In order to construct matrices in AbstractAlgebra.jl, one can first constructs the matrix space itself. This is accomplished with the following constructor.

MatrixSpace(R::Ring, rows::Int, cols::Int; cache::Bool=true)

Construct the space of matrices with the given number of rows and columns over the given base ring. By default such matrix spaces are cached based on the base ring and numbers of rows and columns. If the optional named parameter cached is set to false, no caching occurs.

Here are some examples of creating matrix spaces and making use of the resulting parent objects to coerce various elements into the matrix space.

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)

A = S()
B = S(12)
C = S(R(11))

We also allow matrices over a given base ring to be constructed directly (see the Matrix interface).

## Matrix element constructors

In addition to coercing elements into a matrix space as above, we provide the following functions for constructing explicit matrices.

Also see the Matrix interface for a list of other ways to create matrices.

R[a b c...;...]

Create the matrix over the base ring $R$ consisting of the given rows (separated by semicolons). Each entry is coerced into $R$ automatically. Note that parentheses may be placed around individual entries if the lists would otherwise be ambiguous, e.g. R[1 2; 2 (-3)].

Beware that this syntax does not support the creation of column vectors. See the notation below for creating those.

R[a b c...]

Create the row vector with entries in $R$ consisting of the given entries (separated by spaces). Each entry is coerced into $R$ automatically. Note that parentheses may be placed around individual entries if the list would otherwise be ambiguous, e.g. R[1 2 (-3)].

R[a b c...]'

Create the column vector with entries in $R$ consisting of the given entries (separated by spaces). Each entry is coerced into $R$ automatically. Observe the dash that is used to transpose the row vector notation (for free) to turn it into a column vector. Note that parentheses may be placed around individual entries if the list would otherwise be ambiguous, e.g. R[1 2 (-3)]'.

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)

M = R[t + 1 1; t^2 0]
N = R[t + 1 2 t]
P = R[1 2 t]'

## Submatrices

In addition to the functionality described in the Matrix interface for taking submatrices of a matrix, the following function variant is also available.

sub(M::AbstractAlgebra.MatElem, r1::Int, c1::Int, r2::Int, c2::Int)

Return a copy of the submatrix of $M$ from $(r1, c1)$ to $(r2, c2)$ inclusive. Note that is the copy is modified, the original matrix is not.

source

Examples

M = ZZ[1 2 3; 2 3 4]

N = sub(M, 1, 1, 2, 2)

## Matrix functionality provided by AbstractAlgebra.jl

### Basic matrix functionality

As well as the Ring and Matrix interfaces, the following functions are provided to manipulate matrices and to set and retrieve entries and other basic data associated with the matrices.

rows(a::AbstractAlgebra.MatElem)

Return the number of rows of the given matrix.

source
cols(a::AbstractAlgebra.MatElem)

Return the number of columns of the given matrix.

source

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)

A = S([t + 1 t R(1); t^2 t t; R(-2) t + 2 t^2 + t + 1])
B = S([R(2) R(3) R(1); t t + 1 t + 2; R(-1) t^2 t^3])

r = rows(B)
c = cols(B)
M = A + B
N = 2 + A
M1 = deepcopy(A)
A != B
isone(one(S)) == true
V = A[1:2, :]
W = A^3
Z = divexact(2*A, 2)

### Powering

powers{T <: RingElement}(a::AbstractAlgebra.MatElem{T}, d::Int)

Return an array of matrices $M$ wher $M[i + 1] = a^i$ for $i = 0..d$

source

Examples

M = ZZ[1 2 3; 2 3 4; 4 5 5]

A = powers(M, 4)

### Gram matrix

gram(x::AbstractAlgebra.MatElem)

Return the Gram matrix of $x$, i.e. if $x$ is an $r\times c$ matrix return the $r\times r$ matrix whose entries $i, j$ are the dot products of the $i$-th and $j$-th rows, respectively.

source

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)

A = S([t + 1 t R(1); t^2 t t; R(-2) t + 2 t^2 + t + 1])

B = gram(A)

### Trace

trace(x::AbstractAlgebra.MatElem)

Return the trace of the matrix $a$, i.e. the sum of the diagonal elements. We require the matrix to be square.

source

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)

A = S([t + 1 t R(1); t^2 t t; R(-2) t + 2 t^2 + t + 1])

b = trace(A)

### Content

content(x::AbstractAlgebra.MatElem)

Return the content of the matrix $a$, i.e. the greatest common divisor of all its entries, assuming it exists.

source

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)

A = S([t + 1 t R(1); t^2 t t; R(-2) t + 2 t^2 + t + 1])

b = content(A)

### Permutation

*(P::Generic.perm, x::AbstractAlgebra.MatElem)

Apply the pemutation $P$ to the rows of the matrix $x$ and return the result.

source

Examples

R, t = PolynomialRing(QQ, "t")
S = MatrixSpace(R, 3, 3)
G = PermGroup(3)

A = S([t + 1 t R(1); t^2 t t; R(-2) t + 2 t^2 + t + 1])
P = G([1, 3, 2])

B = P*A

### LU factorisation

lufact{T <: FieldElement}(A::AbstractAlgebra.MatElem{T}, P = PermGroup(rows(A)))

Return a tuple $r, p, L, U$ consisting of the rank of $A$, a permutation $p$ of $A$ belonging to $P$, a lower triangular matrix $L$ and an upper triangular matrix $U$ such that $p(A) = LU$, where $p(A)$ stands for the matrix whose rows are the given permutation $p$ of the rows of $A$.

source
fflu{T <: RingElement}(A::AbstractAlgebra.MatElem{T}, P = PermGroup(rows(A)))

Return a tuple $r, d, p, L, U$ consisting of the rank of $A$, a denominator $d$, a permutation $p$ of $A$ belonging to $P$, a lower triangular matrix $L$ and an upper triangular matrix $U$ such that $p(A) = LD^1U$, where $p(A)$ stands for the matrix whose rows are the given permutation $p$ of the rows of $A$ and such that $D$ is the diagonal matrix diag$(p_1, p_1p_2, \ldots, p_{n-2}p_{n-1}, p_{n-1}$ where the $p_i$ are the inverses of the diagonal entries of $U$. The denominator $d$ is set to $\pm \mbox{det}(S)$ where $S$ is an appropriate submatrix of $A$ ($S = A$ if $A$ is square) and the sign is decided by the parity of the permutation.

source

Examples

R, x = PolynomialRing(QQ, "x")
K, a = NumberField(x^3 + 3x + 1, "a")
S = MatrixSpace(K, 3, 3)

A = S([K(0) 2a + 3 a^2 + 1; a^2 - 2 a - 1 2a; a^2 - 2 a - 1 2a])

r, P, L, U = lufact(A)
r, d, P, L, U = fflu(A)

### Reduced row-echelon form

rref{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Returns a tuple $(r, d, A)$ consisting of the rank $r$ of $M$ and a denominator $d$ in the base ring of $M$ and a matrix $A$ such that $A/d$ is the reduced row echelon form of $M$. Note that the denominator is not usually minimal.

source
rref{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Returns a tuple $(r, d, A)$ consisting of the rank $r$ of $M$ and a denominator $d$ in the base ring of $M$ and a matrix $A$ such that $A/d$ is the reduced row echelon form of $M$. Note that the denominator is not usually minimal.

source
rref{T <: FieldElement}(M::AbstractAlgebra.MatElem{T})

Returns a tuple $(r, A)$ consisting of the rank $r$ of $M$ and a reduced row echelon form $A$ of $M$.

source
isrref{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Return true if $M$ is in reduced row echelon form, otherwise return false.

source
isrref{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Return true if $M$ is in reduced row echelon form, otherwise return false.

source
isrref(M::AbstractAlgebra.MatElem{T}) where {T <: FieldElement}

Return true if $M$ is in reduced row echelon form, otherwise return false.

source

Examples

R, x = PolynomialRing(QQ, "x")
K, a = NumberField(x^3 + 3x + 1, "a")
S = MatrixSpace(K, 3, 3)

M = S([K(0) 2a + 3 a^2 + 1; a^2 - 2 a - 1 2a; a^2 + 3a + 1 2a K(1)])

r, A = rref(M)
isrref(A)

R, x = PolynomialRing(ZZ, "x")
S = MatrixSpace(R, 3, 3)

M = S([R(0) 2x + 3 x^2 + 1; x^2 - 2 x - 1 2x; x^2 + 3x + 1 2x R(1)])

r, d, A = rref(M)
isrref(A)

### Hermite normal form

hnf{T <: RingElement}(A::Mat{T})

Return the upper right row Hermite normal form of $A$.

source
hnf_with_trafo{T <: RingElement}(A::Mat{T}) -> Mat{T}, Mat{T}

Return the tuple $H, U$ consisting of the upper right row Hermite normal form $H$ of $A$ together with invertible matrix $U$ such that $UA = H$.

source

### Determinant

det(M)

Matrix determinant.

Example

julia> M = [1 0; 2 2]
2×2 Array{Int64,2}:
1  0
2  2

julia> det(M)
2.0
source
det{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Return the determinant of the matrix $M$. We assume $M$ is square.

source
det(M)

Matrix determinant.

Example

julia> M = [1 0; 2 2]
2×2 Array{Int64,2}:
1  0
2  2

julia> det(M)
2.0
source
det{T <: FieldElement}(M::AbstractAlgebra.MatElem{T})

Return the determinant of the matrix $M$. We assume $M$ is square.

source
det{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Return the determinant of the matrix $M$. We assume $M$ is square.

source

Examples

R, x = PolynomialRing(QQ, "x")
K, a = NumberField(x^3 + 3x + 1, "a")
S = MatrixSpace(K, 3, 3)

A = S([K(0) 2a + 3 a^2 + 1; a^2 - 2 a - 1 2a; a^2 + 3a + 1 2a K(1)])

d = det(A)

### Rank

rank{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Return the rank of the matrix $M$.

source
rank{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Return the rank of the matrix $M$.

source
rank{T <: FieldElement}(M::AbstractAlgebra.MatElem{T})

Return the rank of the matrix $M$.

source

Examples

R, x = PolynomialRing(QQ, "x")
K, a = NumberField(x^3 + 3x + 1, "a")
S = MatrixSpace(K, 3, 3)

A = S([K(0) 2a + 3 a^2 + 1; a^2 - 2 a - 1 2a; a^2 + 3a + 1 2a K(1)])

d = rank(A)

### Linear solving

solve{T <: FieldElement}(M::AbstractAlgebra.MatElem{T}, b::AbstractAlgebra.MatElem{T})

Given a non-singular $n\times n$ matrix over a field and an $n\times m$ matrix over the same field, return $x$ an $n\times m$ matrix $x$ such that $Ax = b$. If $A$ is singular an exception is raised.

source
solve_rational{T <: RingElement}(M::AbstractAlgebra.MatElem{T}, b::AbstractAlgebra.MatElem{T})

Given a non-singular $n\times n$ matrix over a ring and an $n\times m$ matrix over the same ring, return a tuple $x, d$ consisting of an $n\times m$ matrix $x$ and a denominator $d$ such that $Ax = db$. The denominator will be the determinant of $A$ up to sign. If $A$ is singular an exception is raised.

source
solve_triu{T <: FieldElement}(U::AbstractAlgebra.MatElem{T}, b::AbstractAlgebra.MatElem{T}, unit=false)

Given a non-singular $n\times n$ matrix over a field which is upper triangular, and an $n\times m$ matrix over the same field, return an $n\times m$ matrix $x$ such that $Ax = b$. If $A$ is singular an exception is raised. If unit is true then $U$ is assumed to have ones on its diagonal, and the diagonal will not be read.

source

Examples

R, x = PolynomialRing(QQ, "x")
K, a = NumberField(x^3 + 3x + 1, "a")
S = MatrixSpace(K, 3, 3)
U = MatrixSpace(K, 3, 1)

A = S([K(0) 2a + 3 a^2 + 1; a^2 - 2 a - 1 2a; a^2 + 3a + 1 2a K(1)])
b = U([2a a + 1 (-a - 1)]')

x = solve(A, b)

A = S([a + 1 2a + 3 a^2 + 1; K(0) a^2 - 1 2a; K(0) K(0) a])
b = U([2a a + 1 (-a - 1)]')

x = solve_triu(A, b, false)

R, x = PolynomialRing(ZZ, "x")
S = MatrixSpace(R, 3, 3)
U = MatrixSpace(R, 3, 2)

A = S([R(0) 2x + 3 x^2 + 1; x^2 - 2 x - 1 2x; x^2 + 3x + 1 2x R(1)])
b = U([2x x + 1 (-x - 1); x + 1 (-x) x^2]')

x, d = solve_rational(A, b)

S = MatrixSpace(ZZ, 3, 3)
T = MatrixSpace(ZZ, 3, 1)

A = S([BigInt(2) 3 5; 1 4 7; 9 2 2])
B = T([BigInt(4), 5, 7])

X, d = solve_rational(A, B)

### Inverse

inv(M)

Matrix inverse. Computes matrix N such that M * N = I, where I is the identity matrix. Computed by solving the left-division N = M \ I.

Example

julia> M = [2 5; 1 3]
2×2 Array{Int64,2}:
2  5
1  3

julia> N = inv(M)
2×2 Array{Float64,2}:
3.0  -5.0
-1.0   2.0

julia> M*N == N*M == eye(2)
true
source
inv{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Given a non-singular $n\times n$ matrix over a ring the tuple $X, d$ consisting of an $n\times n$ matrix $X$ and a denominator $d$ such that $AX = dI_n$, where $I_n$ is the $n\times n$ identity matrix. The denominator will be the determinant of $A$ up to sign. If $A$ is singular an exception is raised.

source
inv(M)

Matrix inverse. Computes matrix N such that M * N = I, where I is the identity matrix. Computed by solving the left-division N = M \ I.

Example

julia> M = [2 5; 1 3]
2×2 Array{Int64,2}:
2  5
1  3

julia> N = inv(M)
2×2 Array{Float64,2}:
3.0  -5.0
-1.0   2.0

julia> M*N == N*M == eye(2)
true
source
inv{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Given a non-singular $n\times n$ matrix over a ring the tuple $X, d$ consisting of an $n\times n$ matrix $X$ and a denominator $d$ such that $AX = dI_n$, where $I_n$ is the $n\times n$ identity matrix. The denominator will be the determinant of $A$ up to sign. If $A$ is singular an exception is raised.

source
inv{T <: FieldElement}(M::AbstractAlgebra.MatElem{T})

Given a non-singular $n\times n$ matrix over a field, return an $n\times n$ matrix $X$ such that $AX = I_n$ where $I_n$ is the $n\times n$ identity matrix. If $A$ is singular an exception is raised.

source

Examples

R, x = PolynomialRing(QQ, "x")
K, a = NumberField(x^3 + 3x + 1, "a")
S = MatrixSpace(K, 3, 3)

A = S([K(0) 2a + 3 a^2 + 1; a^2 - 2 a - 1 2a; a^2 + 3a + 1 2a K(1)])

X = inv(A)

R, x = PolynomialRing(ZZ, "x")
S = MatrixSpace(R, 3, 3)

A = S([R(0) 2x + 3 x^2 + 1; x^2 - 2 x - 1 2x; x^2 + 3x + 1 2x R(1)])

X, d = inv(A)

### Nullspace

nullspace(M)

Basis for nullspace of M.

Example

julia> M = [1 0 0; 0 1 0; 0 0 0]
3×3 Array{Int64,2}:
1  0  0
0  1  0
0  0  0

julia> nullspace(M)
3×1 Array{Float64,2}:
0.0
0.0
1.0
source
nullspace{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Returns a tuple $(\nu, N)$ consisting of the nullity $\nu$ of $M$ and a basis $N$ (consisting of column vectors) for the right nullspace of $M$, i.e. such that $MN$ is the zero matrix. If $M$ is an $m\times n$ matrix $N$ will be an $n\times \nu$ matrix. Note that the nullspace is taken to be the vector space kernel over the fraction field of the base ring if the latter is not a field. In AbstractAlgebra we use the name kernel'' for a function to compute an integral kernel.

source
nullspace(M)

Basis for nullspace of M.

Example

julia> M = [1 0 0; 0 1 0; 0 0 0]
3×3 Array{Int64,2}:
1  0  0
0  1  0
0  0  0

julia> nullspace(M)
3×1 Array{Float64,2}:
0.0
0.0
1.0
source
nullspace{T <: RingElement}(M::AbstractAlgebra.MatElem{T})

Returns a tuple $(\nu, N)$ consisting of the nullity $\nu$ of $M$ and a basis $N$ (consisting of column vectors) for the right nullspace of $M$, i.e. such that $MN$ is the zero matrix. If $M$ is an $m\times n$ matrix $N$ will be an $n\times \nu$ matrix. Note that the nullspace is taken to be the vector space kernel over the fraction field of the base ring if the latter is not a field. In AbstractAlgebra we use the name kernel'' for a function to compute an integral kernel.

source
nullspace{T <: FieldElement}(M::AbstractAlgebra.MatElem{T})

Returns a tuple $(\nu, N)$ consisting of the nullity $\nu$ of $M$ and a basis $N$ (consisting of column vectors) for the right nullspace of $M$, i.e. such that $MN$ is the zero matrix. If $M$ is an $m\times n$ matrix $N$ will be an $n\times \nu$ matrix. Note that the nullspace is taken to be the vector space kernel over the fraction field of the base ring if the latter is not a field. In Nemo we use the name kernel'' for a function to compute an integral kernel.

source

Examples

R, x = PolynomialRing(ZZ, "x")
S = MatrixSpace(R, 4, 4)

M = S([-6*x^2+6*x+12 -12*x^2-21*x-15 -15*x^2+21*x+33 -21*x^2-9*x-9;
-8*x^2+8*x+16 -16*x^2+38*x-20 90*x^2-82*x-44 60*x^2+54*x-34;
-4*x^2+4*x+8 -8*x^2+13*x-10 35*x^2-31*x-14 22*x^2+21*x-15;
-10*x^2+10*x+20 -20*x^2+70*x-25 150*x^2-140*x-85 105*x^2+90*x-50])

n, N = nullspace(M)

### Hessenberg form

hessenberg(A::AbstractAlgebra.MatElem{T}) where {T <: RingElement}

Returns the Hessenberg form of $M$, i.e. an upper Hessenberg matrix which is similar to $M$. The upper Hessenberg form has nonzero entries above and on the diagonal and in the diagonal line immediately below the diagonal.

source
ishessenberg{T <: RingElement}(A::AbstractAlgebra.MatElem{T})

Returns true if $M$ is in Hessenberg form, otherwise returns false.

source

Examples

R = ResidueRing(ZZ, 7)
S = MatrixSpace(R, 4, 4)

M = S([R(1) R(2) R(4) R(3); R(2) R(5) R(1) R(0);
R(6) R(1) R(3) R(2); R(1) R(1) R(3) R(5)])

A = hessenberg(M)
ishessenberg(A) == true

### Characteristic polynomial

charpoly{T <: RingElement}(V::Ring, Y::AbstractAlgebra.MatElem{T})

Returns the characteristic polynomial $p$ of the matrix $M$. The polynomial ring $R$ of the resulting polynomial must be supplied and the matrix is assumed to be square.

source

Examples

R = ResidueRing(ZZ, 7)
S = MatrixSpace(R, 4, 4)
T, x = PolynomialRing(R, "x")

M = S([R(1) R(2) R(4) R(3); R(2) R(5) R(1) R(0);
R(6) R(1) R(3) R(2); R(1) R(1) R(3) R(5)])

A = charpoly(T, M)

### Minimal polynomial

minpoly{T <: RingElement}(S::Ring, M::AbstractAlgebra.MatElem{T}, charpoly_only = false)

Returns the minimal polynomial $p$ of the matrix $M$. The polynomial ring $R$ of the resulting polynomial must be supplied and the matrix must be square.

source
minpoly{T <: FieldElement}(S::Ring, M::AbstractAlgebra.MatElem{T}, charpoly_only = false)

Returns the minimal polynomial $p$ of the matrix $M$. The polynomial ring $R$ of the resulting polynomial must be supplied and the matrix must be square.

source
minpoly{T <: RingElement}(S::Ring, M::AbstractAlgebra.MatElem{T}, charpoly_only = false)

Returns the minimal polynomial $p$ of the matrix $M$. The polynomial ring $R$ of the resulting polynomial must be supplied and the matrix must be square.

source

Examples

R = GF(13)
T, y = PolynomialRing(R, "y")

M = R[7 6 1;
7 7 5;
8 12 5]

A = minpoly(T, M)

### Transforms

similarity!{T <: RingElement}(A::AbstractAlgebra.MatElem{T}, r::Int, d::T)

Applies a similarity transform to the $n\times n$ matrix $M$ in-place. Let $P$ be the $n\times n$ identity matrix that has had all zero entries of row $r$ replaced with $d$, then the transform applied is equivalent to $M = P^{-1}MP$. We require $M$ to be a square matrix. A similarity transform preserves the minimal and characteristic polynomials of a matrix.

source

Examples

R = ResidueRing(ZZ, 7)
S = MatrixSpace(R, 4, 4)

M = S([R(1) R(2) R(4) R(3); R(2) R(5) R(1) R(0);
R(6) R(1) R(3) R(2); R(1) R(1) R(3) R(5)])

similarity!(M, 1, R(3))

### (Weak) Popov form

AbstractAlgebra.jl provides algorithms for computing the (weak) Popov of a matrix with entries in a univariate polynomial ring over a field.

weak_popov{T <: PolyElem}(A::Mat{T})

Return the weak Popov form of $A$.

source
weak_popov_with_trafo{T <: PolyElem}(A::Mat{T})

Compute a tuple $(P, U)$ where $P$ is the weak Popov form of $A$ and $U$ is a transformation matrix so that $P = UA$.

source
popov{T <: PolyElem}(A::Mat{T})

Return the Popov form of $A$.

source
popov_with_trafo{T <: PolyElem}(A::Mat{T})

Compute a tuple $(P, U)$ where $P$ is the Popov form of $A$ and $U$ is a transformation matrix so that $P = UA$.

source