Base.nullspace for general matrices

Let’s say I would like to use the nullspace method from julia Base module for a matrix of type AbstractMatrix{MyNumberType}. What interface should I implement in order to make it work? Does such an interface technically exist (since it is not documented)?

Looks like satisfying the asbtract array interface is sufficient. I’m not sure what kind of number you are dealing with though.

julia> struct Foo{T} <: AbstractMatrix{T} end

julia> Base.getindex(::Foo{T}, _...) where T = zero(T)

julia> Base.size(::Foo) = (3,3)

julia> using LinearAlgebra

julia> nullspace(Foo{Int}())
3×3 Matrix{Float64}:
 1.0  0.0  0.0
 0.0  1.0  0.0
 0.0  0.0  1.0

nullspace uses the SVD, so you need a package like GenericLinearAlgebra.jl that supports SVD of your number type.

1 Like

GenericLinearAlgebra despite its name is not really generic. It assumes your numbers are some kind of floating point objects. For example define a 51x51 Hilbert matrix

julia> m=[BigInt(1)//(n+m) for n in 1:51, m in 1:51];

julia> size(GenericLinearAlgebra.nullspace(m))
(51, 1)

while the matrix is invertible. I have written a few routines for really generic linear algebra (valid over any field, and some routines valid over any ring when applicable, similar to det_bareiss) which are in the package GenLinearAlgebra. These includes nullspace. I do not claim to have written the best possible algorithms, and I do not know if there has been other similar efforts. Note that the algorithms valid over any field are often not good ones for floating-point numbers.

1 Like

The algorithms for exact arithmetic are completely different from the algorithms for inexact arithmetic in many cases; GenericLinearAlgebra is geared towards the latter.

In inexact arithmetic, nullspace uses an SVD (which is inherently inexact due to the Abel–Ruffini theorem), with a cutoff for the singular values. In exact arithmetic, you’d compute a nullspace basis using elimination.

Yes, that was my point. I do not know what is the best vocabulary to talk about this dichotomy:
exact/inexact, exact/approximate, complete field/general field ? You could possibly also do approximate computation with p-adic numbers.