1d set/getindex for sparse

[tanmaykm]

ref #7023

This implements setindex! and getindex with 1d and boolean indices for AbstractSparseMatrix. An optimized implementation for SparseMatrixCSC can be done separately.

[ViralBShah]

I suspect these are going to be very slow, especially the setindex cases.

[tanmaykm]

Yes. This PR is not marked to fix the issue.
I'll either update this or follow this up with another one with an implementation on SparseMatrixCSC.

[ViralBShah]

How about doing it as part of this PR?

[IainNZ]

I don't understand when this could would be used?

[ViralBShah]

Just like you index a dense matrix with a vector, so can you with a sparse matrix. See the linked issue. It is also needed for completeness. Currently the following does not work:

julia> s = sprand(5,5,0.5);
julia> x = [5,6,7];
julia> s[x]
ERROR: no method similar(SparseMatrixCSC{Float64,Int64}, Type{Float64}, (Int64,))
 in getindex at abstractarray.jl:364

[IainNZ]

I more meant the "abstractness" of it - surely this would be a lot slower than a specialized method? This would like it would for any AbstractMatrix.

[tanmaykm]

The abstract methods would make it quicker to add new concrete types and optimize them incrementally.

[ViralBShah]

I agree. The abstract implementations are not useful except as an implementation guide. We can merge when the implementation of the remaining methods is done.

[tanmaykm]

Added implementation for SparseMatrixCSC and rebased.

[ViralBShah]

Tests in the gcc build are failing. I can't quite tell if it is related to this PR or not.

[tanmaykm]

It passes now. I had to rename a local variable to something other than sum2, but I wonder why.

[tanmaykm]

@ViralBShah Rebased and passing the tests.

[ViralBShah]

Could you remove the AbstractSparseMatrix methods? These are only likely to be slow for any sparse matrix data structure. Otherwise, please feel free to merge whenever you are ready.

[ViralBShah]

Also Cc:@mauro3

[ViralBShah]

Is it possible to reuse the binarysearch method in sparsematrix.jl?

[tanmaykm]

I don't think it is possible to reuse the binarysearch method here. This expects slightly different result than what binarysearch would give. Here ridx points to the index uptill row, and has a valid value irrespective of whether row was found in rowval[r1:r2]. And r1 is advanced till ridx so that subsequent searches in the same column are narrower.

[mauro3]

This sounds like the binary search which is implemented in sort.jl under the name of searchsortedfirst (which was not what I needed, thus I wrote my own). Also, note that you should probably use >>>, see http://googleresearch.blogspot.co.uk/2006/06/extra-extra-read-all-about-it-nearly.html

[tanmaykm]

Thanks, searchsortedfirst indeed can be used. I shall update the PR.

[ViralBShah]

I am curious if this has any performance impact. Worth checking for sure.

[tanmaykm]

Will check that. There are a few other places in sparsematrix.jl where searchsorted can be used. Shall replace them as well if there's no performance impact.

[tanmaykm]

There is a difference in performance. Here's a code snippet and some results for comparison: https://gist.github.com/tanmaykm/6871d2d014df4a92f40f

But it does not seem significant.

[tanmaykm]

Removed the abstract methods.

[tanmaykm]

Added changes to use existing search methods. No noticeable performance impact.

Before:

setindex!: elapsed time: 0.013420411599999998 seconds (363840.0 bytes allocated)
diag: elapsed time: 0.0003199585 seconds (8088.0 bytes allocated)

With this change:

setindex!: elapsed time: 0.012916476900000001 seconds (363840.0 bytes allocated)
diag: elapsed time: 0.0003488753 seconds (8088.0 bytes allocated)

Test code: https://gist.github.com/tanmaykm/e878ff0fe14f6b01b6dd

[ViralBShah]

Thanks for working through the various things.

[jiahao]

Any idea what's going on here?

$ make
...
linalg.jl
Warning: New definition 
    Ac_ldiv_B! could not show value of type Tuple at linalg/umfpack.jl:356
is ambiguous with: 
    Ac_ldiv_B! could not show value of type Tuple at linalg/umfpack.jl:351.
To fix, define 
    Ac_ldiv_B! could not show value of type Tuple
before the new definition.

[ViralBShah]

I don't get this. I saw this yesterday though, but it seemed to go away after a make clean.

[jiahao]

Ok. The other odd linalg-related thing is this test failure; it's new to me.

[mbauman]

I think it's an inference problem. The test works just fine at top-level.

[ViralBShah]

@jiahao The test failure seems to be unrelated to this PR. Should we have a new issue? I don't get the failure on my mac.

[mbauman]

As far as I'm aware, it's only happening in my dict-print PR: #7147. I cannot reproduce on master on my local machine.

[mauro3]

I've seen the Ac_ldiv_B! warning popping up at random as well.