Attempting to add reltol and abstol to rank function

[miguelraz]

Attempting to address JuliaLang/LinearAlgebra.jl#15

Any critiques/comments/help/future advice is more than appreciated, willing to start from scratch.

1. Added tests and documentation.
2. Shout out to @jw3126 for helping out.
   I am clearly doing many things wrong here - I bungled two different commits, but don't know how to rebase my commits from linalg to the most recent master.

...Any adults willing to help a drowning kid are more than welcome.

[StefanKarpinski]

Not sure what the status of these branches is, but this got the unzip stuff mixed in with it.

[fredrikekre]

Looks like its only the last commit that is relevant for this PR?

Lets try this (if you use command line git):

1. git checkout rankfix make sure we are on the correct branch
2. git branch rankfix-backup make a backup of this branch
3. git reset --hard master reset this branch to master branch
4. git cherry-pick a795d23acf9c368e4df3f75c661f9ca8a0b327e6 to apply only the last commit to the newly reset branch
5. git push -f origin rankfix force-push to the PR branch (if your own repo is called origin)

[jw3126]

In the code svdvals(A) is called twice. I would be very impressed if the compiler could optimize this out (I don't think so).

[miguelraz]

@fredrikekre, man, I owe you a beer and a half. Thanks for bailing me out - I thought of the solution as I was getting up from bed (back up, take a copy of master, cherry pick what you need onto a new branch) but the command line escaped me. Thanks a zillion!
I am still trying to get the push into the branch I want - I seem to have reset master way too far back.

@jw3126 just tried to run this. Does this look optimized to you? My LLVM is not up to snuff.

julia> function rank(A::AbstractMatrix; abstol=0.0, reltol=0.0)
           thresh = abstol + reltol * maximum(svdvals(A))
           sum(svdvals(A) .> thresh )
       end
rank (generic function with 1 method)

julia> A = rand(3,3)
3×3 Array{Float64,2}:
 0.78815   0.636962  0.222599
 0.165679  0.573335  0.591731
 0.492888  0.491938  0.207964

julia> @code_llvm rank(A)

; Function Attrs: uwtable
define i64 @julia_rank_67840(i8**) #0 !dbg !5 {
top:
  %1 = call i8**** @jl_get_ptls_states() JuliaLang/julia#4
  %2 = alloca [8 x i8**], align 8
  %.sub = getelementptr inbounds [8 x i8**], [8 x i8**]* %2, i64 0, i64 0
  %3 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 3
  %4 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 2
  %5 = bitcast i8*** %3 to i8*
  call void @llvm.memset.p0i8.i32(i8* %5, i8 0, i32 40, i32 8, i1 false)
  %6 = bitcast [8 x i8**]* %2 to i64*
  store i64 12, i64* %6, align 8
  %7 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 1
  %8 = bitcast i8**** %1 to i64*
  %9 = load i64, i64* %8, align 8
  %10 = bitcast i8*** %7 to i64*
  store i64 %9, i64* %10, align 8
  store i8*** %.sub, i8**** %1, align 8
  store i8** null, i8*** %4, align 8
  %11 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 7
  %12 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 6
  %13 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 5
  %14 = getelementptr [8 x i8**], [8 x i8**]* %2, i64 0, i64 4
  store i8** inttoptr (i64 2147434912 to i8**), i8*** %3, align 8
  store i8** inttoptr (i64 2151328992 to i8**), i8*** %14, align 8
  store i8** inttoptr (i64 2151328992 to i8**), i8*** %13, align 8
  store i8** inttoptr (i64 2147434904 to i8**), i8*** %12, align 8
  store i8** %0, i8*** %11, align 8
  %15 = call i8** @jl_invoke(i8** inttoptr (i64 2221565456 to i8**), i8*** %3, i32 5)
  store i8** %15, i8*** %4, align 8
  %16 = bitcast i8** %15 to i64*
  %17 = load i64, i64* %16, align 16
  %18 = load i64, i64* %10, align 8
  store i64 %18, i64* %8, align 8
  ret i64 %17
}

[jw3126]

I can't decrypt that piece of llvm either, I would just time both versions. Note also that A can be any user defined matrix type with any user definded svdvals method. So even if the compiler could optimize for Matrix, it does not guarantee that it will also optimize out for MyMatrix. So I would only call once anyway.

[miguelraz]

Ach so @jw3126 - you I take it you suggest to call svdvals(A) once, store to array and then use that value twice?
I can patch that up.

[jw3126]

Yes.

[oscardssmith]

We now have this.