Support a view of a matrix as a right-hand side b

docs/src/api.md

@@ -62,4 +62,5 @@ Krylov.ktypeof
 Krylov.kzeros
 Krylov.kones
 Krylov.vector_to_matrix
+Krylov.matrix_to_vector
 ```

src/krylov_utils.jl

@@ -219,7 +219,13 @@ function ktypeof(v::S) where S <: AbstractVector
 end
 
 function ktypeof(v::S) where S <: SubArray
-  return ktypeof(v.parent)
+  vp = v.parent
+  if isa(vp, DenseMatrix)
+    M = typeof(vp)
+    return matrix_to_vector(M)  # view of a row or a column of a matrix
+  else
+    return ktypeof(vp)  # view of a vector
+  end
 end
 
 """
@@ -228,18 +234,36 @@ end
 Return the dense matrix storage type `M` related to the dense vector storage type `S`.
 """
 function vector_to_matrix(::Type{S}) where S <: DenseVector
-  V = hasproperty(S, :body) ? S.body : S
-  par = V.parameters
+  T = hasproperty(S, :body) ? S.body : S
+  par = T.parameters
   npar = length(par)
   (2 ≤ npar ≤ 3) || error("Type $S is not supported.")
   if npar == 2
-    M = V.name.wrapper{par[1], 2}
+    M = T.name.wrapper{par[1], 2}
   else
-    M = V.name.wrapper{par[1], 2, par[3]}
+    M = T.name.wrapper{par[1], 2, par[3]}
   end
   return M
 end
 
+"""
+    S = matrix_to_vector(M)
+
+Return the dense vector storage type `S` related to the dense matrix storage type `M`.
+"""
+function matrix_to_vector(::Type{M}) where M <: DenseMatrix
+  T = hasproperty(M, :body) ? M.body : M
+  par = T.parameters
+  npar = length(par)
+  (2 ≤ npar ≤ 3) || error("Type $M is not supported.")
+  if npar == 2
+    S = T.name.wrapper{par[1], 1}
+  else
+    S = T.name.wrapper{par[1], 1, par[3]}
+  end
+  return S
+end
+
 """
     v = kzeros(S, n)
 

test/gpu/amd.jl

@@ -68,10 +68,8 @@ include("gpu.jl")
     #   Krylov.@kref!(n, x, y, c, s)
     # end
 
-    @testset "vector_to_matrix" begin
-      S = ROCVector{FC}
-      M2 = Krylov.vector_to_matrix(S)
-      @test M2 == M
+    @testset "conversion -- $FC" begin
+      test_conversion(S, M)
     end
 
     ε = eps(T)

test/gpu/gpu.jl

@@ -45,3 +45,8 @@ function test_solver(S, M)
   solver = GmresSolver(n, n, memory, S)
   solve!(solver, A, b)  # Test that we don't have errors
 end
+
+function test_conversion(S, M)
+  @test Krylov.vector_to_matrix(S) == M
+  @test Krylov.vector_to_matrix(M) == S
+end

test/gpu/intel.jl

@@ -76,10 +76,8 @@ end
     #   Krylov.@kref!(n, x, y, c, s)
     # end
 
-    @testset "vector_to_matrix" begin
-      S = oneVector{FC}
-      M2 = Krylov.vector_to_matrix(S)
-      @test M2 == M
+    @testset "conversion -- $FC" begin
+      test_conversion(S, M)
     end
 
     ε = eps(T)

test/gpu/metal.jl

@@ -80,10 +80,8 @@ end
     #   Krylov.@kref!(n, x, y, c, s)
     # end
 
-    @testset "vector_to_matrix" begin
-      S = MtlVector{FC}
-      M2 = Krylov.vector_to_matrix(S)
-      @test M2 == M
+    @testset "conversion -- $FC" begin
+      test_conversion(S, M)
     end
 
     ε = eps(T)

test/gpu/nvidia.jl

@@ -94,6 +94,7 @@ include("gpu.jl")
 
   for FC in (Float32, Float64, ComplexF32, ComplexF64)
     S = CuVector{FC}
+    V = CuSparseVector{FC}
     M = CuMatrix{FC}
     T = real(FC)
     n = 10
@@ -144,10 +145,8 @@ include("gpu.jl")
       Krylov.@kref!(n, x, y, c, s)
     end
 
-    @testset "vector_to_matrix" begin
-      S = CuVector{FC}
-      M2 = Krylov.vector_to_matrix(S)
-      @test M2 == M
+    @testset "conversion -- $FC" begin
+      test_conversion(S, M)
     end
 
     ε = eps(T)
@@ -177,5 +176,22 @@ include("gpu.jl")
     @testset "solver -- $FC" begin
       test_solver(S, M)
     end
+
+    @testset "ktypeof -- $FC" begin
+        dv = S(rand(FC, 10))
+        b = view(dv, 4:8)
+        @test Krylov.ktypeof(dv) <: S
+        @test Krylov.ktypeof(b)  <: S
+
+        dm = M(rand(FC, 10, 10))
+        b = view(dm, :, 4:8)
+        @test Krylov.ktypeof(b) <: S
+
+        sv = V(sprand(FC, 10, 0.5))
+        b = view(sv, 4:8)
+        @test Krylov.ktypeof(sv) <: S
+        @test Krylov.ktypeof(b)  <: S
+      end
+    end
   end
 end

test/test_aux.jl

@@ -129,25 +129,21 @@
 
   @testset "ktypeof" begin
     # test ktypeof
-    a = rand(Float32, 10)
-    b = view(a, 4:8)
-    @test Krylov.ktypeof(a) == Vector{Float32}
-    @test Krylov.ktypeof(b) == Vector{Float32}
-
-    a = rand(Float64, 10)
-    b = view(a, 4:8)
-    @test Krylov.ktypeof(a) == Vector{Float64}
-    @test Krylov.ktypeof(b) == Vector{Float64}
-
-    a = sprand(Float32, 10, 0.5)
-    b = view(a, 4:8)
-    @test Krylov.ktypeof(a) == Vector{Float32}
-    @test Krylov.ktypeof(b) == Vector{Float32}
-
-    a = sprand(Float64, 10, 0.5)
-    b = view(a, 4:8)
-    @test Krylov.ktypeof(a) == Vector{Float64}
-    @test Krylov.ktypeof(b) == Vector{Float64}
+    for FC in (Float32, Float64, ComplexF32, ComplexF64)
+      dv = rand(FC, 10)
+      b = view(dv, 4:8)
+      @test Krylov.ktypeof(dv) == Vector{FC}
+      @test Krylov.ktypeof(b)  == Vector{FC}
+
+      dm = rand(FC, 10, 10)
+      b = view(dm, :, 4:8)
+      @test Krylov.ktypeof(b) == Vector{FC}
+
+      sv = sprand(FC, 10, 0.5)
+      b = view(sv, 4:8)
+      @test Krylov.ktypeof(sv) == Vector{FC}
+      @test Krylov.ktypeof(b)  == Vector{FC}
+    end
   end
 
   @testset "vector_to_matrix" begin
@@ -159,6 +155,15 @@
     end
   end
 
+  @testset "matrix_to_vector" begin
+    # test matrix_to_vector
+    for FC in (Float32, Float64, ComplexF32, ComplexF64)
+      M = Matrix{FC}
+      S = Krylov.matrix_to_vector(M)
+      @test S == Vector{FC}
+    end
+  end
+
   @testset "macros" begin
     # test macros
     for FC ∈ (Float16, Float32, Float64, ComplexF16, ComplexF32, ComplexF64)