LU with empty matrices (#263)

* LU with empty matrices

* Add \ tests

Author: dlfivefifty

Project.toml

@@ -1,6 +1,6 @@
 name = "BandedMatrices"
 uuid = "aae01518-5342-5314-be14-df237901396f"
-version = "0.17.3"
+version = "0.17.4"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"

src/banded/BandedLU.jl

@@ -84,7 +84,11 @@ function lu!(A::BandedMatrix{T}, pivot::Union{Val{false}, Val{true}} = Val(true)
     end
     m= size(A,1)
     l,u = bandwidths(A) # l of the bands are ignored and overwritten
-    _, ipiv = LAPACK.gbtrf!(l, u-l, m, bandeddata(A))
+    if !iszero(m)
+        _, ipiv = LAPACK.gbtrf!(l, u-l, m, bandeddata(A))
+    else
+        ipiv = Int[]
+    end
     return BandedLU{T,typeof(A)}(A, ipiv, zero(BlasInt))
 end
 

src/banded/linalg.jl

@@ -25,7 +25,8 @@ function ldiv!(A::BandedLU{T,<:BandedMatrix}, B::StridedVecOrMat{T}) where {T<:B
     m = size(A.factors,1)
     l,u = bandwidths(A.factors)
     data = bandeddata(A.factors)
-    LAPACK.gbtrs!('N', l, u-l, m, data, A.ipiv, B)
+    iszero(m) || LAPACK.gbtrs!('N', l, u-l, m, data, A.ipiv, B)
+    B
 end
 
 function ldiv!(A::BandedLU{T}, B::AbstractVecOrMat{Complex{T}}) where T<:Real

test/test_bandedlu.jl

@@ -151,4 +151,13 @@ struct _foo <: Number end
         @test_throws BoundsError size(BLU, -1)
         @test_throws BoundsError size(BLU,  0)
     end
+
+    @testset "zero matrix" begin
+        for A in (BandedMatrix{Float64}(undef, 0,0, 1,1),
+                  BandedMatrix{Float64}(undef, 0,3, 1,1),
+                  BandedMatrix{Float64}(undef, 0,0, -1,-2))
+            @test lu(A).factors == zeros(size(A)...)
+            @test lu(A) \ zeros(0) == zeros(0)
+        end
+    end
 end