Allow more general eltypes in sparse array multiplication #33205

pablosanjose · 2019-09-09T21:53:35Z

Replaces one(T) by 1 in five-argument mul! calls with SparseArray types, where T was before the target eltype. This allows eltypes that don't define one by relying on the fact that 1 is always a multiplicative identity for AbtractArrays as far as I can tell ~~(even for Bool eltypes)~~ EDIT: changed 1, 0 to true, false.

With this we can in particular use StaticArrays as eltypes of SparseMatrixCSC. It does not allow Array eltypes, since they don't define zero, which is used in mul!. Solving the issue also for eltypes that don't define zero is of much broader scope than the present PR due to the pervasive assumption throughout SparseArrays that zero(eltype(A)) for A::AbstractSparseArray is defined.

As a consequence of allowing general (e.g. matrix-like) eltypes in all product combinations involving a sparse vector or matrices (see tests) a range of changes throughout SparseArrays were needed, including in copy(::AbstractSparseMatrix), halfperm!, spmatmul and _At_or_Ac_mul_B.

This fix does not address a possible version of #33169 without any sparse factors (i.e. all changes here are strictly confined to the SparseArrays lib)

(This is my first PR to Julia stdlib. It might be a good idea to show this to Nanosoldier, and to have some careful review if possible)

Rebase

copy(::Adjoint/Transpose)+ simplesmatrix fixes added news fix promote_op for sparsevec mul! added + to SimpleSMatrix (required for matprod) trying to make mat * vec work :-( _At_or_Ac_mul_B! fix to accept target eltype more simplesmatrix fixes yet more fixes to simplesmatrix more stringent tests

pablosanjose · 2019-09-09T21:56:01Z

(I see some whitespace fixes are needed... Can somebody briefly point out how to address that?)

fredrikekre · 2019-09-09T22:20:35Z

Can somebody briefly point out how to address that?

Remove trailing whitespaces.

andreasnoack · 2019-09-11T07:26:15Z

Seems reasonable to me. @tkf any thoughts on this? Should we use Booleans here?

pablosanjose · 2019-09-11T07:42:09Z

I think you are right, booleans are in fact better. I somehow though I saw true * SArray{Bool}(...) being promoted to SArray{Int}, but it is not. The universal one and zero is indeed true and false. I'll push an update with booleans

tkf · 2019-09-11T08:02:00Z

I don't think it matters too much here as we are not using α=false case (which matters when the inputs produce NaNs). But maybe it'd play nice with some custom eltype implementing GF(2) or something? I think it's good to use Booleans anyway for consistency. (Thank you @pablosanjose for the change!)

pablosanjose · 2019-09-11T10:02:27Z

Funny, there is a strange win32 failure that didn't get triggered before the last commit. I assume this is unrelated to this PR?
EDIT: a restart fixed it

pablosanjose · 2019-09-16T14:22:08Z

A friendly bump looking for some volunteer to review if possible. Thanks!

timholy

Very nice. Completely optional to get rid of more type parameters, this could basically be merged as is.

stdlib/SparseArrays/src/linalg.jl

stdlib/SparseArrays/src/sparsevector.jl

timholy · 2019-09-18T14:41:33Z

LGTM. I should allow someone else take a look, but if you don't hear soon let's merge.

stdlib/SparseArrays/src/sparsevector.jl

stdlib/SparseArrays/test/simplesmatrix.jl

Co-Authored-By: Fredrik Ekre <[email protected]>

pablosanjose · 2019-09-18T14:55:53Z

Many thanks to both for the review and approval!

pablosanjose · 2019-09-20T10:15:32Z

Could this be merged, perhaps?

KristofferC · 2020-01-15T09:04:32Z

This seems like it broke the tests for JWAS.jl:

ERROR: LoadError: LoadError: MethodError: no method matching zero(::Type{Any})
Closest candidates are:
  zero(::Type{Union{Missing, T}}) where T at missing.jl:105
  zero(::Type{Missing}) at missing.jl:103
  zero(::Type{LibGit2.GitHash}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.5/LibGit2/src/oid.jl:220
  ...
Stacktrace:
 [1] zero(::Type{Any}) at ./missing.jl:105
 [2] _zeros_eltypes at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.5/SparseArrays/src/higherorderfns.jl:204 [inlined]
 [3] _noshapecheck_map(::typeof(+), ::SparseArrays.SparseMatrixCSC{Any,Int64}, ::SparseArrays.SparseMatrixCSC{Float32,Int64}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.5/SparseArrays/src/higherorderfns.jl:159
 [4] map(::typeof(+), ::SparseArrays.SparseMatrixCSC{Any,Int64}, ::SparseArrays.SparseMatrixCSC{Float32,Int64}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.5/SparseArrays/src/higherorderfns.jl:1153
 [5] +(::SparseArrays.SparseMatrixCSC{Any,Int64}, ::SparseArrays.SparseMatrixCSC{Float32,Int64}) at /Users/sabae/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.5/SparseArrays/src/sparsematrix.jl:1648
 [6] addVinv(::JWAS.MME) at /Users/kristoffer/.julia/packages/JWAS/4ltS3/src/1.JWAS/src/buildMME/random_effects.jl:205
 [7] getMME(::JWAS.MME, ::DataFrame) at /Users/kristoffer/.julia/packages/JWAS/4ltS3/src/1.JWAS/src/buildMME/build_MME.jl:264
 [8] init_mixed_model_equations(::JWAS.MME, ::DataFrame, ::Bool) at /Users/kristoffer/.julia/packages/JWAS/4ltS3/src/1.JWAS/src/JWAS.jl:565
 [9] runMCMC(::JWAS.MME, ::DataFrame; heterogeneous_residuals::Bool, chain_length::Int64, starting_value::Bool, burnin::Int64, output_samples_frequency::Int64, output_samples_file::String, update_priors_frequency::Int64, methods::String, estimate_variance::Bool, Pi::Float64, estimatePi::Bool, estimateScale::Bool, single_step_analysis::Bool, pedigree::Bool, categorical_trait::Bool, missing_phenotypes::Bool, constraint::Bool, causal_structure::Bool, outputEBV::Bool, output_heritability::Bool, seed::Int64, printout_model_info::Bool, printout_frequency::Int64, big_memory::Bool, double_precision::Bool) at /Users/kristoffer/.julia/packages/JWAS/4ltS3/src/1.JWAS/src/JWAS.jl:234
 [10] top-level scope at /Users/kristoffer/.julia/packages/JWAS/4ltS3/test/test_BayesianAlphabet.jl:100
 [11] include(::String) at ./client.jl:439
 [12] top-level scope at /Users/kristoffer/.julia/packages/JWAS/4ltS3/test/runtests.jl:3
 [13] include(::String) at ./client.jl:439
 [14] top-level scope at none:6

pablosanjose · 2020-01-15T09:20:20Z

First of all, my apologies, this was my first merged PR and I feel bad!
I'm trying to understand the reason for the failure. It seems the stacktrace involves codepaths in files not touched by this PR, except [5], which involves + in sparsematrix.jl, but this is a function not modified here. Any further clue why that failure is related to the PR?

KristofferC · 2020-01-15T09:34:50Z

First of all, my apologies, this was my first merged PR and I feel bad!

Please don't! This happens all the time, this just either means that we lack some tests in Base or that the package does something "weird". Consider it a rite of passage ;)
I'll try figure out a bit better what is happening, I just posted it directly here for reference and so someone else can also look into it if they feel like it.

The reason for suspecting this PR is that it was bisected here: #34238 (comment). As you can see in that thread, we started with over 200 new breakages for 1.4 so breaking packages is definitely not a big deal.

KristofferC · 2020-01-15T10:17:46Z

So the difference is in:

julia> a = convert(SparseMatrixCSC{AbstractFloat, Int}, sparse(rand(Float32, 2,2)));

julia> b = sparse(rand(Float32, 2,2));

julia> typeof(a*b)
SparseMatrixCSC{Any,Int64}

On 1.3, this returned a SparseMatrixCSC{AbstractFloat,Int64}. I'm not sure what is "correct" here...

KristofferC · 2020-01-15T10:28:57Z

I made a PR to the package to fix it: reworkhow/JWAS.jl#59.

Not sure if we need to change anything in Julia.

pablosanjose · 2020-01-15T10:35:13Z

I see! With this PR the target eltype of the product a*b is computed with

julia> Base.promote_op(LinearAlgebra.matprod, AbstractFloat, Float32)
Any

while in 1.3 we merely did

julia> promote_type(AbstractFloat, Float32)
AbstractFloat

Probably the correct solution within Julia is to make Base.promote_op smarter, returning AbstractFloat in this case. What do you think?

pablosanjose · 2020-01-15T10:39:09Z

Uhm, that promote_op method goes into the return_type abyss where there be dragons. Such a fix is far beyond the red line for me.

KristofferC · 2020-01-15T10:39:34Z

Yeah, promote_op is known to be a bit of a wart (for that reason).

pablosanjose added 3 commits September 3, 2018 22:10

Merge pull request #1 from JuliaLang/master

5912bb6

Rebase

Merge remote-tracking branch 'origin/master'

cc454a4

ViralBShah added the sparse Sparse arrays label Sep 9, 2019

pablosanjose added 2 commits September 10, 2019 08:51

whitespace fix

bca5043

News

0fafdac

pablosanjose changed the title ~~Allow general eltypes that implement zero in sparse mul!~~ Allow more general eltypes in sparse array multiplication Sep 10, 2019

use true, false instead of 1, 0

e3b9539

timholy approved these changes Sep 18, 2019

View reviewed changes

stdlib/SparseArrays/src/linalg.jl Show resolved Hide resolved

stdlib/SparseArrays/src/sparsevector.jl Outdated Show resolved Hide resolved

pablosanjose added 2 commits September 18, 2019 16:36

removed redundant type parameters

a5eeb86

fix a missing 1 -> true change

80e059a

fredrikekre approved these changes Sep 18, 2019

View reviewed changes

stdlib/SparseArrays/src/sparsevector.jl Outdated Show resolved Hide resolved

stdlib/SparseArrays/test/simplesmatrix.jl Show resolved Hide resolved

pablosanjose and others added 2 commits September 18, 2019 16:53

fix a missing zero(Ty) -> false

5403167

Co-Authored-By: Fredrik Ekre <[email protected]>

add a missing licence header

fb6e144

Co-Authored-By: Fredrik Ekre <[email protected]>

fredrikekre merged commit 34d2b87 into JuliaLang:master Sep 20, 2019

ericphanson mentioned this pull request Dec 30, 2019

edit and update NEWS for 1.4 #34219

Merged

maleadt mentioned this pull request Jan 15, 2020

WIP: Backports for 1.4-RC1 #34238

Merged

28 tasks

blegat mentioned this pull request Jan 27, 2020

Fix tests for Julia v1.4-rc1 jump-dev/MathOptInterface.jl#1010

Merged

Uh oh!

Allow more general eltypes in sparse array multiplication #33205

Allow more general eltypes in sparse array multiplication #33205

Uh oh!

Conversation

pablosanjose commented Sep 9, 2019 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

pablosanjose commented Sep 9, 2019

Uh oh!

fredrikekre commented Sep 9, 2019

Uh oh!

andreasnoack commented Sep 11, 2019

Uh oh!

pablosanjose commented Sep 11, 2019

Uh oh!

tkf commented Sep 11, 2019

Uh oh!

pablosanjose commented Sep 11, 2019 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

pablosanjose commented Sep 16, 2019

Uh oh!

timholy left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

timholy commented Sep 18, 2019

Uh oh!

Uh oh!

Uh oh!

pablosanjose commented Sep 18, 2019

Uh oh!

pablosanjose commented Sep 20, 2019

Uh oh!

KristofferC commented Jan 15, 2020

Uh oh!

pablosanjose commented Jan 15, 2020 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

KristofferC commented Jan 15, 2020 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

KristofferC commented Jan 15, 2020 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

KristofferC commented Jan 15, 2020

Uh oh!

pablosanjose commented Jan 15, 2020

Uh oh!

pablosanjose commented Jan 15, 2020

Uh oh!

KristofferC commented Jan 15, 2020

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

7 participants

pablosanjose commented Sep 9, 2019 •

edited

Loading

pablosanjose commented Sep 11, 2019 •

edited

Loading

pablosanjose commented Jan 15, 2020 •

edited

Loading

KristofferC commented Jan 15, 2020 •

edited

Loading

KristofferC commented Jan 15, 2020 •

edited

Loading