Sparse linear algebra
Introduction
This chapter deals with sparse linear algebra over commutative rings and fields.
Sparse linear algebra, that is, linear algebra with sparse matrices, plays an important role in various algorithms in algebraic number theory. For example, it is one of the key ingredients in the computation of class groups and discrete logarithms using index calculus methods.
Sparse rows
Building blocks for sparse matrices are sparse rows, which are modelled by objects of type SRow
. More precisely, the type is of parametrized form SRow{T}
, where T
is the element type of the base ring . For example, SRow{ZZRingElem}
is the type for sparse rows over the integers.
It is important to note that sparse rows do not have a fixed number of columns, that is, they represent elements of . In particular any two sparse rows over the same base ring can be added.
Creation
sparse_row
— Methodsparse_row(R::Ring, J::Vector{Tuple{Int, T}}) -> SRow{T}
Constructs the sparse row with , where . The elements must belong to the ring .
sparse_row
— Methodsparse_row(R::Ring, J::Vector{Tuple{Int, T}}) -> SRow{T}
Constructs the sparse row with , where . The elements must belong to the ring .
sparse_row(R::Ring, J::Vector{Tuple{Int, Int}}) -> SRow
Constructs the sparse row over with , where .
sparse_row
— Methodsparse_row(R::Ring, J::Vector{Int}, V::Vector{T}) -> SRow{T}
Constructs the sparse row over with , where and .
Basic operations
Rows support the usual operations:
+
,-
,==
- multiplication by scalars
div
,divexact
getindex
— Methodgetindex(A::SRow, j::Int) -> RingElem
Given a sparse row and an index return .
add_scaled_row
— Methodadd_scaled_row(A::SRow{T}, B::SRow{T}, c::T) -> SRow{T}
Returns the row .
add_scaled_row
— Methodadd_scaled_row(A::SRow{T}, B::SRow{T}, c::T) -> SRow{T}
Returns the row .
transform_row
— Methodtransform_row(A::SRow{T}, B::SRow{T}, i::Int, j::Int, a::T, b::T, c::T, d::T)
Returns the tuple .
length
— Methodlength(A::SRow)
Returns the number of nonzero entries of .
Change of base ring
change_base_ring
— Methodchange_base_ring(R::Ring, A::SRow) -> SRow
Create a new sparse row by coercing all elements into the ring .
Maximum, minimum and 2-norm
maximum
— Methodmaximum(A::SRow{T}) -> T
Returns the largest entry of .
maximum
— Methodmaximum(A::SRow{T}) -> T
Returns the largest entry of .
minimum
— Methodminimum(A::SRow{T}) -> T
Returns the smallest entry of .
minimum(A::NfRelOrdIdl) -> NfOrdIdl
minimum(A::NfRelOrdIdl) -> NfRelOrdIdl
Returns the ideal where is the maximal order of the coefficient ideals of .
minimum
— Methodminimum(A::SRow{T}) -> T
Returns the smallest entry of .
norm2
— Methodnorm2(A::SRow{T} -> T
Returns .
Functionality for integral sparse rows
lift
— Methodlift(A::SRow{zzModRingElem}) -> SRow{ZZRingElem}
Return the sparse row obtained by lifting all entries in .
mod!
— Methodmod!(A::SRow{ZZRingElem}, n::ZZRingElem) -> SRow{ZZRingElem}
Inplace reduction of all entries of modulo to the positive residue system.
mod_sym!
— Methodmod_sym!(A::SRow{ZZRingElem}, n::ZZRingElem) -> SRow{ZZRingElem}
Inplace reduction of all entries of modulo to the symmetric residue system.
mod_sym!
— Methodmod_sym!(A::SRow{ZZRingElem}, n::Integer) -> SRow{ZZRingElem}
Inplace reduction of all entries of modulo to the symmetric residue system.
maximum
— Methodmaximum(abs, A::SRow{ZZRingElem}) -> ZZRingElem
Returns the largest, in absolute value, entry of .
Sparse matrices
Let be a commutative ring. Sparse matrices with base ring are modelled by objects of type SMat
. More precisely, the type is of parametrized form SRow{T}
, where T
is the element type of the base ring. For example, SMat{ZZRingElem}
is the type for sparse matrices over the integers.
In contrast to sparse rows, sparse matrices have a fixed number of rows and columns, that is, they represent elements of the matrices space . Internally, sparse matrices are implemented as an array of sparse rows. As a consequence, unlike their dense counterparts, sparse matrices have a mutable number of rows and it is very performant to add additional rows.
Construction
sparse_matrix
— Methodsparse_matrix(R::Ring) -> SMat
Return an empty sparse matrix with base ring .
sparse_matrix
— Methodsparse_matrix(R::Ring, n::Int, m::Int) -> SMat
Return a sparse times zero matrix over .
Sparse matrices can also be created from dense matrices as well as from julia arrays:
sparse_matrix
— Methodsparse_matrix(A::MatElem; keepzrows::Bool = true)
Constructs the sparse matrix corresponding to the dense matrix . If keepzrows
is false, then the constructor will drop any zero row of .
sparse_matrix
— Methodsparse_matrix(R::Ring, A::Matrix{T}) -> SMat
Constructs the sparse matrix over corresponding to .
sparse_matrix
— Methodsparse_matrix(R::Ring, A::Matrix{T}) -> SMat
Constructs the sparse matrix over corresponding to .
The normal way however, is to add rows:
push!
— Methodpush!(A::SMat{T}, B::SRow{T}) where T
Appends the sparse row B
to A
.
Sparse matrices can also be concatenated to form larger ones:
vcat!
— Methodvcat!(A::SMat, B::SMat) -> SMat
Vertically joins and inplace, that is, the rows of are appended to .
vcat
— Methodvcat(A::SMat, B::SMat) -> SMat
Vertically joins and .
hcat!
— Methodhcat!(A::SMat, B::SMat) -> SMat
Horizontally concatenates and , inplace, changing .
hcat
— Methodhcat(A::SMat, B::SMat) -> SMat
Horizontally concatenates and .
(Normal julia is also supported)
There are special constructors:
identity_matrix
— Methodidentity_matrix(::Type{SMat}, R::Ring, n::Int)
Return a sparse times identity matrix over .
zero_matrix
— Methodzero_matrix(::Type{SMat}, R::Ring, n::Int)
Return a sparse times zero matrix over .
zero_matrix
— Methodzero_matrix(::Type{SMat}, R::Ring, n::Int, m::Int)
Return a sparse times zero matrix over .
Slices:
sub
— Methodsub(A::SMat, r::UnitRange, c::UnitRange) -> SMat
Return the submatrix of , where the rows correspond to and the columns correspond to .
Transpose:
transpose
— Methodtranspose(A::SMat) -> SMat
Returns the transpose of .
Elementary Properties
sparsity
— Methodsparsity(A::SMat) -> Float64
Return the sparsity of A
, that is, the number of zero-valued elements divided by the number of all elements.
density
— Methoddensity(A::SMat) -> Float64
Return the density of A
, that is, the number of nonzero-valued elements divided by the number of all elements.
nnz
— Methodnnz(A::SMat) -> Int
Return the number of non-zero entries of .
nrows
— Methodnrows(A::SMat) -> Int
Return the number of rows of .
ncols
— Methodncols(A::SMat) -> Int
Return the number of columns of .
isone
— Methodisone(A::SMat)
Tests if is an identity matrix.
iszero
— Methodiszero(A::SMat)
Tests if is a zero matrix.
isupper_triangular
— Methodisupper_triangular(A::SMat)
Returns true if and only if is upper (right) triangular.
maximum
— Methodmaximum(A::SMat{T}) -> T
Finds the largest entry of .
minimum
— Methodminimum(A::SMat{T}) -> T
Finds the smallest entry of .
maximum
— Methodmaximum(abs, A::SMat{ZZRingElem}) -> ZZRingElem
Finds the largest, in absolute value, entry of .
elementary_divisors
— Methodelementary_divisors(A::SMat{ZZRingElem}) -> Vector{ZZRingElem}
The elementary divisors of , i.e. the diagonal elements of the Smith normal form of .
solve_dixon_sf
— Methodsolve_dixon_sf(A::SMat{ZZRingElem}, b::SRow{ZZRingElem}, is_int::Bool = false) -> SRow{ZZRingElem}, ZZRingElem
solve_dixon_sf(A::SMat{ZZRingElem}, B::SMat{ZZRingElem}, is_int::Bool = false) -> SMat{ZZRingElem}, ZZRingElem
For a sparse square matrix of full rank and a sparse matrix (row), find a sparse matrix (row) and an integer s.th. holds. The algorithm is a Dixon-based linear p-adic lifting method. If \code{is_int} is given, then is assumed to be . In this case rational reconstruction is avoided.
hadamard_bound2
— Methodhadamard_bound2(A::SMat{T}) -> T
The square of the product of the norms of the rows of .
echelon_with_transform
— Methodechelon_with_transform(A::SMat{zzModRingElem}) -> SMat, SMat
Find a unimodular matrix and an upper-triangular s.th. holds.
reduce_full
— Methodreduce_full(A::SMat{ZZRingElem}, g::SRow{ZZRingElem},
trafo = Val{false}) -> SRow{ZZRingElem}, Vector{Int}
Reduces modulo and assumes that is upper triangular.
The second return value is the array of pivot elements of that changed.
If trafo
is set to Val{true}
, then additionally an array of transformations is returned.
hnf!
— Methodhnf!(A::SMat{ZZRingElem})
Inplace transform of into upper right Hermite normal form.
hnf
— Methodhnf(A::SMat{ZZRingElem}) -> SMat{ZZRingElem}
Return the upper right Hermite normal form of .
snf
— Methodsnf(A::SMat{ZZRingElem})
The Smith normal form (snf) of .
hnf_extend!
— Methodhnf_extend!(A::SMat{ZZRingElem}, b::SMat{ZZRingElem}, offset::Int = 0) -> SMat{ZZRingElem}
Given a matrix in HNF, extend this to get the HNF of the concatenation with .
is_diagonal
— Methodis_diagonal(A::SMat) -> Bool
True iff only the i-th entry in the i-th row is non-zero.
det
— Methoddet(A::SMat{ZZRingElem})
The determinant of using a modular algorithm. Uses the dense (zzModMatrix) determinant on for various primes .
det_mc
— Methoddet_mc(A::SMat{ZZRingElem})
Computes the determinant of using a LasVegas style algorithm, i.e. the result is not proven to be correct. Uses the dense (zzModMatrix) determinant on for various primes .
valence_mc
— Methodvalence_mc{T}(A::SMat{T}; extra_prime = 2, trans = Vector{SMatSLP_add_row{T}}()) -> T
Uses a Monte-Carlo algorithm to compute the valence of . The valence is the valence of the minimal polynomial of , thus the last non-zero coefficient, typically .
The valence is computed modulo various primes until the computation stabilises for extra_prime
many.
trans
, if given, is a SLP (straight-line-program) in GL(n, Z). Then the valence of trans
* is computed instead.
saturate
— Methodsaturate(A::SMat{ZZRingElem}) -> SMat{ZZRingElem}
Computes the saturation of , that is, a basis for $\mathbf{Q}\otimes M \meet \mathbf{Z}^n$, where is the row span of and the number of rows of .
Equivalently, return for an invertible rational matrix , such that is integral and the elementary divisors of are all trivial.
hnf_kannan_bachem
— Methodhnf_kannan_bachem(A::SMat{ZZRingElem}) -> SMat{ZZRingElem}
Compute the Hermite normal form of using the Kannan-Bachem algorithm.
diagonal_form
— Methoddiagonal_form(A::SMat{ZZRingElem}) -> SMat{ZZRingElem}
A matrix that is diagonal and obtained via unimodular row and column operations. Like a snf without the divisibility condition.
Manipulation/ Access
getindex
— Methodgetindex(A::SMat, i::Int, j::Int)
Given a sparse matrix , return the entry .
getindex
— Methodgetindex(A::SMat, i::Int) -> SRow
Given a sparse matrix and an index , return the -th row of .
setindex!
— Methodsetindex!(A::SMat, b::SRow, i::Int)
Given a sparse matrix , a sparse row and an index , set the -th row of equal to .
swap_rows!
— Methodswap_rows!(A::SMat{T}, i::Int, j::Int)
Swap the -th and -th row of inplace.
swap_cols!
— Methodswap_cols!(A::SMat, i::Int, j::Int)
Swap the -th and -th column of inplace.
scale_row!
— Methodscale_row!(A::SMat{T}, i::Int, c::T)
Multiply the -th row of by inplace.
add_scaled_col!
— Methodadd_scaled_col!(A::SMat{T}, i::Int, j::Int, c::T)
Add times the -th column to the -th column of inplace, that is, , where denote the columns of .
add_scaled_row!
— Methodadd_scaled_row!(A::SMat{T}, i::Int, j::Int, c::T)
Add times the -th row to the -th row of inplace, that is, , where denote the rows of .
transform_row!
— Methodtransform_row!(A::SMat{T}, i::Int, j::Int, a::T, b::T, c::T, d::T)
Applies the transformation to , where are the rows of .
diagonal
— Methoddiagonal(A::SMat) -> ZZRingElem[]
The diagonal elements of in an array.
reverse_rows!
— Methodreverse_rows!(A::SMat)
Inplace inversion of the rows of .
mod_sym!
— Methodmod_sym!(A::SMat{ZZRingElem}, n::ZZRingElem)
Inplace reduction of all entries of modulo to the symmetric residue system.
find_row_starting_with
— Methodfind_row_starting_with(A::SMat, p::Int) -> Int
Tries to find the index such that and for all . It is assumed that is upper triangular. If such an index does not exist, find the smallest index larger.
reduce
— Methodreduce(A::SMat{ZZRingElem}, g::SRow{ZZRingElem}, m::ZZRingElem) -> SRow{ZZRingElem}
Given an upper triangular matrix over the integers, a sparse row and an integer , this function reduces modulo and returns modulo with respect to the symmetric residue system.
reduce
— Methodreduce(A::SMat{ZZRingElem}, g::SRow{ZZRingElem}) -> SRow{ZZRingElem}
Given an upper triangular matrix over a field and a sparse row , this function reduces modulo .
reduce
— Methodreduce(A::SMat{T}, g::SRow{T}) -> SRow{T}
Given an upper triangular matrix over a field and a sparse row , this function reduces modulo .
rand_row
— Methodrand_row(A::SMat) -> SRow
Return a random row of the sparse matrix .
Changing of the ring:
map_entries
— Methodmap_entries(f, A::SMat) -> SMat
Given a sparse matrix and a callable object , this function will construct a new sparse matrix by applying to all elements of .
change_base_ring
— Methodchange_base_ring(R::Ring, A::SMat)
Create a new sparse matrix by coercing all elements into the ring .
Arithmetic
Matrices support the usual operations as well
+
,-
,==
div
,divexact
by scalars- multiplication by scalars
Various products:
mul
— Methodmul(A::SMat{T}, b::AbstractVector{T}) -> Vector{T}
Return the product as a dense vector.
mul
— Methodmul(A::SMat{T}, b::AbstractMatrix{T}) -> Matrix{T}
Return the product as a dense array.
mul
— Methodmul(A::SMat{T}, b::MatElem{T}) -> MatElem
Return the product as a dense matrix.
mul
— Methodmul(A::SRow, B::SMat) -> SRow
Return the product as a sparse row.
Other:
sparse
— Methodsparse(A::SMat) -> SparseMatrixCSC
The same matrix, but as a sparse matrix of julia type SparseMatrixCSC
.
ZZMatrix
— MethodZZMatrix(A::SMat{ZZRingElem})
The same matrix , but as an ZZMatrix
.
ZZMatrix
— MethodZZMatrix(A::SMat{T}) where {T <: Integer}
The same matrix , but as an ZZMatrix
. Requires a conversion from the base ring of to .
Array
— MethodArray(A::SMat{T}) -> Matrix{T}
The same matrix, but as a two-dimensional julia array.