That's very cool and close to what I need. Thanks. I'm using adjacency
lists (vectors of vectors) but I'm sure I can adapt this code.
On Sunday, November 15, 2015 at 6:12:19 PM UTC-8, Tim Holy wrote:
>
> Haven't tested timing carefully, but something like the following should
> be
> pretty decent:
>
> julia> immutable LowerTriangularIterator{A<:AbstractMatrix}
> data::A
> end
>
> julia> Base.start(iter::LowerTriangularIterator) = (1,1)
> start (generic function with 49 methods)
>
> julia> Base.done(iter::LowerTriangularIterator, state) = state[2] >
> size(iter.data, 2)
> done (generic function with 51 methods)
>
> julia> @inline function Base.next(iter::LowerTriangularIterator, state)
> item = iter.data[state[1], state[2]]
> newstate = ifelse(state[1] == size(iter.data, 1), (state[2]+1,
> state[2]+1), (state[1]+1, state[2]))
> item, newstate
> end
>
>
> To make it vectorize (i.e., SIMD), you might also need definitions like
> these:
>
> https://github.com/JuliaLang/julia/blob/aa24ee0d0f334e8d9ecbced64bd975a6cc5f240d/base/multidimensional.jl#L134-L150
>
>
> but perhaps vectorization isn't applicable to your problem anyway.
>
> --Tim
>
> On Sunday, November 15, 2015 05:06:23 PM Seth wrote:
> > Thank you, Tim - I'll definitely check it out. Right now I can't seem to
> > find an efficient way of doing it, but there are a few tricks in
> > Iterators.jl that might work.
> >
> > On Sunday, November 15, 2015 at 1:33:29 PM UTC-8, Tim Holy wrote:
> > > This may also be a case where writing an iterator that just visits the
> > > elements with i >= j might be what you want. You can study the
> > > implementation
> > > of CartesianIndex in Base for inspiration.
> > >
> > > --Tim
> > >
> > > On Saturday, November 14, 2015 02:18:45 PM hustf wrote:
> > > > I believe this work-in-progress type may be adapted? I wanted to
> make
> > >
> > > this
> > >
> > > > faster by using tuples or immutable fixed vectors to store the data
> or
> > > > whatnot, but I guess I'm currently stuck. Help wanted.
> > > >
> > > > Compared to looping unnecessarily over a full matrix, this speeds
> things
> > >
> > > up
> > >
> > > > by around 35%.
> > > >
> > > > immutable Maq7 <:AbstractMatrix{Int}
> > > > data::Array{Int64,2} # upper triangular except diagonal, stored as
> a
> > > > horizontal vector
> > > > width::Int64
> > > > function Maq7(mat::Matrix{Int64})
> > > > if issym(mat)
> > > > n=size(mat,1)
> > > > if n>1
> > > > v=mat[2:n,1]
> > > > for col = 2:(n-1)
> > > > append!(v,mat[(col+1):n,col])
> > > > end
> > > > new(v', size(mat,1))
> > > > else
> > > > error("Ouch, not > (1x1)")
> > > > end
> > > > else
> > > > error("Ouch, not square symmetric")
> > > > end
> > > > end
> > > > end
> > > > Base.size(A::Maq7)= A.width, A.width
> > > > Base.length(A::Maq7)= div(A.width*(A.width-1), 2)
> > > > Base.getindex(A::Maq7, i::Int) = Base.getindex(A.data, i)
> > > > function Base.getindex(A::Maq7, r::Int , c::Int )
> > > > if c > r
> > > > return getindex(A.data, div( (r-1)*(2*A.width -r) , 2 ) + c - r )
> > > > elseif r > c
> > > > return getindex(A.data, div( (c-1)*(2*A.width -c) , 2 ) + r - c )
> > > > else
> > > > return zero(eltype(A.data))
> > > > end
> > > > end
> > > > Base.linearindexing(::Type{Maq7}) = Base.LinearFast()
> > > >
> > > >
> > > >
> > > > Example:
> > > > julia> ms
> > > >
> > > > 4x4 Array{Int64,2}:
> > > > 7 1 2 3
> > > > 1 8 4 5
> > > > 2 4 9 6
> > > > 3 5 6 10
> > > >
> > > > julia> m=Maq7(ms)
> > > >
> > > > 4x4 Maq7:
> > > > 0 1 2 3
> > > > 1 0 4 5
> > > > 2 4 0 6
> > > > 3 5 6 0
> > > >
> > > > julia> collect(m)
> > > >
> > > > 6-element Array{Int64,1}:
> > > > 1
> > > > 2
> > > > 3
> > > > 4
> > > > 5
> > > > 6
> > > >
> > > > julia> length(m)
> > > > 6
>
>