Fix indexing/slicing of LayoutQ (#42)

Project.toml

@@ -1,6 +1,6 @@
 name = "MatrixFactorizations"
 uuid = "a3b82374-2e81-5b9e-98ce-41277c0e4c87"
-version = "0.9.5"
+version = "0.9.6"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"

src/MatrixFactorizations.jl

@@ -109,6 +109,10 @@ if VERSION < v"1.10-"
         lmul!(Q, y)
     end
     Base.@propagate_inbounds layout_getindex(A::LayoutQ, I::CartesianIndex) = A[to_indices(A, (I,))...]
+    Base.@propagate_inbounds layout_getindex(A::LayoutQ, I::Int...) =
+        Base.invoke(Base.getindex, Tuple{AbstractQ, typeof.(I)...}, A, I...)
+    Base.@propagate_inbounds layout_getindex(A::LayoutQ, I::AbstractVector{Int}, J::AbstractVector{Int}) =
+        hcat((A[:, j][I] for j in J)...)
 
 
     (*)(Q::LayoutQ, P::LayoutQ) = mul(Q, P)
@@ -122,11 +126,8 @@ axes(Q::LayoutQ) = axes(Q, 1), axes(Q, 2)
 copy(Q::LayoutQ) = Q
 Base.@propagate_inbounds getindex(A::LayoutQ, I...) = layout_getindex(A, I...)
 # by default, fall back to AbstractQ  methods
-_layout_getindex(A::LayoutQ, I...) =
+layout_getindex(A::LayoutQ, I...) =
     Base.invoke(Base.getindex, Tuple{AbstractQ, typeof.(I)...}, A, I...)
-layout_getindex(A::LayoutQ, I...) = _layout_getindex(A, I...)
-# disambiguate with ArrayLayouts
-layout_getindex(A::LayoutQ, I::Int...) = _layout_getindex(A, I...)
 
 size(Q::LayoutQ, dim::Integer) = size(getfield(Q, :factors), dim == 2 ? 1 : dim)
 size(Q::LayoutQ) = size(Q, 1), size(Q, 2)

test/runtests.jl

@@ -74,6 +74,8 @@ rectangularQ(Q::LinearAlgebra.AbstractQ) = Matrix(Q) # convert(Array, Q)
                     @inferred qrunblocked(a)
                     q, r = qra.Q, qra.R
                     @test_throws ErrorException qra.Z
+                    @test q[1,1] ≈ Matrix(q)[1,1]
+                    @test q[1:2,1:2] ≈ Matrix(q)[1:2,1:2]
                     @test q'*squareQ(q) ≈ Matrix(I, a_1, a_1)
                     @test q*squareQ(q)' ≈ Matrix(I, a_1, a_1)
                     @test q'*Matrix(1.0I, a_1, a_1)' ≈ squareQ(q)'
@@ -103,6 +105,8 @@ rectangularQ(Q::LinearAlgebra.AbstractQ) = Matrix(Q) # convert(Array, Q)
                     @inferred qrunblocked(a[:, 1:n1], Val(false))
                     q,r   = qra.Q, qra.R
                     @test_throws ErrorException qra.Z
+                    @test q[1,1] ≈ Matrix(q)[1,1]
+                    @test q[1:2,1:2] ≈ Matrix(q)[1:2,1:2]
                     @test q'*squareQ(q) ≈ Matrix(I, a_1, a_1)
                     @test q'*rectangularQ(q) ≈ Matrix(I, a_1, n1)
                     @test q*r ≈ a[:, 1:n1]
@@ -127,6 +131,8 @@ rectangularQ(Q::LinearAlgebra.AbstractQ) = Matrix(Q) # convert(Array, Q)
 
                     qra = qrunblocked(a[:,1:n1], Val(false))
                     q, r = qra.Q, qra.R
+                    @test q[1,1] ≈ Matrix(q)[1,1]
+                    @test q[1:2,1:2] ≈ Matrix(q)[1:2,1:2]
                     @test rmul!(copy(squareQ(q)'), q) ≈ Matrix(I, n, n)
                     @test_throws DimensionMismatch rmul!(Matrix{eltya}(I, n+1, n+1),q)
                     @test rmul!(squareQ(q), adjoint(q)) ≈ Matrix(I, n, n)

test/test_rq.jl

@@ -92,6 +92,8 @@ const Our=MatrixFactorizations
                     @inferred rq(a)
                     q, r  = rqa.Q, rqa.R
                     @test_throws ErrorException rqa.Z
+                    @test q[1,1] ≈ Matrix(q)[1,1]
+                    @test q[1:2,1:2] ≈ Matrix(q)[1:2,1:2]
                     @test q'*q ≈ Matrix(I, a_1, a_1)
                     @test q*q' ≈ Matrix(I, a_1, a_1)
                     @test q'*Matrix(1.0I, a_1, a_1)' ≈ Matrix(q')