Old thread alert :-).
SubArrays have changed a lot in 0.4, and now you can create a SubArray with
irregular indexes. You just can't pass them directly to Lapack---only
SubArrays created with Range-type indexes have the uniform strides that Lapack
needs.
Your example exploits the fact that you can always copy the data to a new
array (which has uniform strides) and pass that to Lapack. You could imagine
doing that "automatically" behind the scenes, but then the scary situation
would be that users might inadvertently make such a copy 10^6 times when, had
the user realized it were necessary, it could be done just once. So SubArrays
with irregular indexing will simply throw an error if you try to pass them to
Lapack. To circumvent the error, just call copy on your SubArray and then call
matrix algebra routines.
Best,
--Tim
On Friday, March 27, 2015 08:00:45 AM René Hiemstra wrote:
> Hi Tim,
>
> I don't really understand the problem with passing the more flexible
> SubMatrix (or SubArray) to Lapack. Suppose you have something as follows,
>
> typealias Index Union(Int, Vector{Int}, RangeIndex, ...)
>
> type SubMatrix{T,I<:Index, J<:Index}
> data :: Matrix{T}
> rows :: I
> cols :: J
> end
>
> A = rand(10,10)
> B = SubMatrix(A,[1,3,5], [1,6,7,10])
> C = B.data[B.rows, B.cols]
>
> Then we can always pass C to Lapack, right?
>
> Rene
>
> On Friday, March 22, 2013 at 9:52:35 AM UTC+1, Tim Holy wrote:
> > I'm guessing it's because you've (properly) made it a type of
> > AbstractArray,
> > but you haven't provided a size(A, d) method. The generic show() method
> > for
> > AbstractArray relies on that being available.
> >
> > I'm guessing you should be able to at least see whether it's been properly
> > created using dump().
> >
> > --Tim
> >
> > On Thursday, March 21, 2013 05:51:24 PM dr wrote:
> > > Thanks for that pointer. I tried for a bit to get this to work but could
> > > not quite get there. My julia knowhow is a bit too basic to understand
> >
> > the
> >
> > > problem.
> > > I have attached my attempt at an "ArrayView".
> > >
> > > julia> require("ArrayView.jl")
> > >
> > > julia> x=reshape([1:100],10,10);
> > >
> > > julia> xx=view(x,[1,3,5],[2,5,7])
> > > ERROR: Error showing value of type
> >
> > AbstractArrayView{Int64,2,Array{Int64,2},(Array{Int64,1},Array{Int64,1})}:
> > > MethodError(size,(SYSTEM: show(lasterr) caused an error
> > > ERROR: no method
> >
> > size(AbstractArrayView{Int64,2,Array{Int64,2},(Array{Int64,1},Array{Int64,
> > 1}>
> > > )},) in length at abstractarray.jl:15
> > >
> > > in show_delim_array at show.jl:61
> > > in show at show.jl:85
> > > in show at show.jl:6
> > > in show at base.jl:83
> > > in error_show at repl.jl:22
> > > in error_show at repl.jl:67
> > > in anonymous at no file:57
> > > in with_output_color at error.jl:36
> > > in display_error at client.jl:55
> > >
> > > julia> xx[1,1]
> > > ERROR: no method
> >
> > ref(AbstractArrayView{Int64,2,Array{Int64,2},(Array{Int64,1},Array{Int64,1
> > })>
> > > },Int64,Int64)
> > >
> > > On Saturday, March 16, 2013 7:42:18 AM UTC-7, Tim Holy wrote:
> > > > A bit of archaeology here...if someone wants to add a "view" type that
> > > > supports irregular indexing, it seems that
> > > >
> > > > git show 49462addba9555810fce8938ccbedad8a558db05:j/tensor.j
> > > >
> > > > from the main julia directory will get you a version of SubArray that
> >
> > did
> >
> > > > this.
> > > >
> > > > Probably the main challenge will be coming up with a decent name for
> >
> > the
> >
> > > > type.
> > > >
> > > > --Tim
