BLAS/LAPACK code assumes all StridedMatrices are column-major

[stevengj]

As I understand it, one of the purposes of the DenseMatrix type from JuliaLang/julia#987, as discussed in JuliaLang/julia#6212, is to include things like row-major NumPy arrays. However, it looks like most of the BLAS/LAPACK code assumes, without checking, that StridedMatrix arguments are column-major (or subsets of column-major matrices), i.e. they assume that stride(A,1) == 1.

Probably these routines should have stride(A,1) != 1 && ...generic fallback... for now. Ideally, they could do

if stride(A,1) != 1
   if stride(A,2) == 1
       ...handle row-major case by treating A as if it were transposed....
   else
       ...generic fallback
   end
end

cc: @andreasnoack, @lindahua

[lindahua]

I don't have problems with this. Looks like it is going to be a lot of work though.

[stevengj]

I agree that it is a lot of work. Until it is done, however, it will be dangerous to define any new DenseMatrix subtypes, as a consequence of linear algebra dispatching to LAPACK for DenseMatrix{Float64} etc.

For 0.3, I see several options:

• Leave things as-is, and hope that not too many people define such subtypes.
• Leave things as-is, but change the manual to recommend (for now) that DenseMatrix subtypes have unit stride along the first dimension.
• Leave things as-is, but don't document DenseMatrix at all; leave it as an undocumented internal type that is not ready for public consumption.
• Change StridedMatrix to be Union(Matrix,SubArray{Matrix}), as before.
• Try to fix this in time for 0.3, which may delay the release still further (unless someone is very quick with a patch).

@JeffBezanson, @jiahao?

[lindahua]

A stopgap: inserting this to each function that wraps BLAS:

stride(A, 1) == 1 || error( ... )

[lindahua]

The array view types in ArrayViews.jl, which are subtypes of DenseMatrix, now work well with BLAS & LAPACK. Changing StridedMatrix to Union(Matrix, SubArray{Matrix}) would cause problems.

[jiahao]

I'm in favour of undocumented DenseMatrix for 0.3.

[andreasnoack]

We should probably add chkstride1 to all BLAS and LAPACK methods.

cc: @jakebolewski

[kshyatt]

Since we recently added Row* etc I'm pinging @andyferris on this - can we close?

[tkelman]

that only handled vectors, not strided matrices or blas and lapack interfaces, so no

[tkelman]

though strided matrix is a union recently, and I think many of the calls do have chkstride1 now, so it may have been fixed just not by the rowvector work.

[ViralBShah]

Closing this since it is unlikely to do much about row-major arrays - and a package would be a better place for such things.

[stevengj]

Note that this issue is not about supporting row-major arrays, it is about making sure that we call chkstride1 or similar in any LAPACK call that accepts a StridedArray but requires column-major (i.e. all LAPACK calls).

[ViralBShah]

I believe we do chkstride1 widely throughout blas.jl and lapack.jl. Perhaps we need a specific issue for places we don't do it. Let's reopen if necessary.