Sparse with default value "Inf"

General Usage
,
23

Try with

struct Weight <: Number
           w::Float64
end

For some reason you can find the following line in stdlib/SparseArrays/src/higherorderfns.jl:

@inline _iszero(x) = x == 0
@inline _iszero(x::Number) = Base.iszero(x)
1 Like
24

This is not heavily tested, lacks convenience features, bounds checks, and various optimizations, but may get you started:

struct DefaultArray{T,N} <: AbstractArray{T,N}
    default::T
    size::NTuple{N,Int}
    elements::Dict{NTuple{N,Int},T}
end

function DefaultArray(default::T, size::NTuple{N,Int})  where {T,N}
    DefaultArray(default, size, Dict{NTuple{N,Int},T}())
end

Base.size(A::DefaultArray) = A.size

Base.getindex(A::DefaultArray, I::Int...) = get(A.elements, I, A.default)

Base.setindex!(A::DefaultArray, v, I::Int...) = A.elements[I] = v

julia> A = DefaultArray(Inf, (2,2))
2×2 DefaultArray{Float64,2}:
 Inf  Inf
 Inf  Inf

julia> A[1,2] = 3.0
3.0

julia> A
2×2 DefaultArray{Float64,2}:
 Inf   3.0
 Inf  Inf

(Even after years of using Julia, it never ceases to amaze me how easy it is to get a fairly fast plug-in prototype in <10 LOC.)

3 Likes
25

You could perhaps store the inverse of all entries and then when you access a zero it will be infinite :thinking:

2 Likes
26

I’m not 100% sure why mapped arrays doesn’t work here. I might be missing something

using MappedArrays
using SparseArrays

function f(x)
    if x == 0.0
        return Inf
    elseif x === NaN
        return 0.0
    else
        return x
    end
end 

function finv(y)
    if y == Inf
        return 0.0
    elseif y === 0.0
        return NaN
    else
        return y
    end
end 
A = spzeros(3, 3)
mA = mappedarray(f, finv, A)

julia> mA
3×3 mappedarray(f, finv, ::SparseMatrixCSC{Float64,Int64}) with eltype Float64:
 Inf  Inf  Inf
 Inf  Inf  Inf
 Inf  Inf  Inf

julia> mA .= Inf
3×3 mappedarray(f, finv, ::SparseMatrixCSC{Float64,Int64}) with eltype Float64:
 Inf  Inf  Inf
 Inf  Inf  Inf
 Inf  Inf  Inf

julia> nnz(mA.data)
0

julia> mA .-= 0
3×3 mappedarray(f, finv, ::SparseMatrixCSC{Float64,Int64}) with eltype Float64:
 Inf  Inf  Inf
 Inf  Inf  Inf
 Inf  Inf  Inf

julia> nnz(mA.data)
0

julia> mA[1, 2] = 3; mA[1,1] = 0
julia> mA
3×3 mappedarray(f, finv, ::SparseMatrixCSC{Float64,Int64}) with eltype Float64:
   0.0    3.0  Inf
 Inf    Inf    Inf
 Inf    Inf    Inf

julia> nnz(mA.data)
2

Unless you are looking for something different.

1 Like
27

It doesn’t surprise me that Mathematica allows this, but I’m not aware of any numerical sparse matrix libraries (UMFPACK, SuiteSparse, ARPACK, etc.) that support it. It would be a large amount of work to support the full range of efficient operations on a sparse matrix type that deviates from the standard CSC representation. Work that Mathematica has done, or course, but it’s not like we can use that.

2 Likes
28

I perfectly understand @StefanKarpinski, via this thread i discovered that just a tiny bit of work on my part can turn into a working solution that does not imply augmenting julia itself. Altough via the 0 \to \infty monstrosity,i could see that a compressed array like structure similar to SparseArrays, without the linear algebra could be truly useful, it would also benefit Graph Algorithms.

at the moment i’m using @Tamas_Papp solution to move forward, if some more efficient sparse array like structure comes along i can also switch (since DefaultArray behaves like an Array :smiley: ), thank you man!

29

You are welcome.

If this is genuinely useful and otherwise not covered by any existing packages, please consider packaging it up and registering (all you need is bounds checks which I did not implement, and some tests).

1 Like
30

i already made a repo

obviously i will work on it alongside the real problem i’m solving

1 Like