LAPACK.trttf! and LAPACK.tfttr! do not mutate arguments

[dmbates]

AFAICS these two functions copy to freshly allocated storage so the names should not have a trailing !

[andreasnoack]

I agree. Are you using these functions? I'm considering to remove them from base and into a package. They are not exported, tested or up to date with the development in the rest of the linear algebra code so it would require some work to get this code in the right shape and I prefer to use that time on other parts of the linear algebra code.

[tkelman]

Are we using packed format anywhere in base?

[andreasnoack]

Here https://github.com/JuliaLang/julia/blob/b8aa30dc3464053ebeb8aaad63cf61c3539282d0/base/linalg/rectfullpacked.jl

I added them years ago, but right now they are not really integrated with the rest of the linear algebra code. In theory, we could use them for the usual Cholesky, but the packed version is slightly slower than the usual version and I guess most users can afford the extra memory consumption.

2015-05-01 6:07 GMT-04:00 Tony Kelman notifications@github.com:

Are we using packed format anywhere in base?

— Reply to this email directly or view it on GitHub https://github.com/JuliaLang/julia/issues/11072#issuecomment-98090863.

[dmbates]

@andreasnoack I was looking at those routines because I was going to test the speed of the rectangular full packed Cholesky. This is not urgent at all. I agree that it would probably be better to create a package for them as they are a kind of a niche item. It happens that I have some relatively large symmetric, positive-definite matrices sparse matrices to decompose and in most cases, even after a fill-reducing permutation, the number of nonzeros is about 1/3 of the total possible. I wanted to compare the speed and memory use of the dense and sparse representations.

[tkelman]

Thanks for the link Andreas, that's indeed pretty vestigial and minimal right now. IIRC packed storage isn't as friendly on the cache and hasn't been well-optimized in openblas/MKL/Accelerate etc. I agree, like the Woodbury types, if it's not being used by anyone it should probably move to a package eventually.

@dmbates I came across a new BSD-licensed library we may want to have a look at making a package for in JuliaSparse, http://portal.nersc.gov/project/sparse/strumpack/ - not sure whether it'll be suited for your problems at all but figured it was worth mentioning somewhere.

[andreasnoack]

@dmbates Ok. I'll probably remove them and hopefully also make them more usable during the process.

@tkelman This is a special kind of packed storage that is more BLAS-3 friendly than the classical BLAS/LAPACK packed storage scheme. The STRUMPACK package looks interesting. I'll have a look.

[dmbates]

@tkelman Thanks for the pointer to STRUMPACK. Right now I have a rather specialized problem in that the structure of the matrices I am working with is derived from a bipartite graph. In my wanderings around the sparse matrix literature I think I saw some methods for dealing with this structure but I can't recall exactly where and apparently the Google is not strong in me today and I am unable to find it again.

The structure of the matrix that I wish to decompose is has two diagonal matrices, sizes N and n with N > n, that are the blocks on the diagonal. The off-diagonal block corresponds to a bipartite graph. Performing a Cholesky decomposition by blocks preserves the structure of the N by N diagonal block and the n by N off-diagonal block. That is, fill-in only occurs in the n by n block and that fill-in will be determined by the bipartite graph.

[ViralBShah]

Move it to LinearAlgebra.jl or a different package?

[tkelman]

Ideally wherever it goes could be registered and reasonably well-tested before moving it out of base, if anyone is actually using it. Also no hurry, doesn't need to be done until after 0.4 IMO.

[ViralBShah]

I feel like we can remove this without too much pain for 0.4.