Repository: invenia/NamedDims.jl (julia)

Issue:
Support eltype of array differing from eltype of factorization
# Explanation

An error occurs when left-dividing (`\`) certain `NamedDimsArray`s. Even if this error didn't happen there would still be #163, but they seem independent.

The error occurs when the `Factorization` is `QRPivoted`, which [has 2 parameters](https://github.com/JuliaLang/julia/blob/fbe28e4779f3d2b8683b5e721d7aa312ec656ac2/stdlib/LinearAlgebra/src/qr.jl#L163) rather than the 1 [assumed](https://github.com/invenia/NamedDims.jl/blob/236d9207d9ba9436063e11809047d4cf9315e7f3/src/functions_linearalgebra.jl#L12). The `Factorization` chosen for the operation depends on the arrays. Necessary conditions for the error seem to be (1) non-square matrix, (2) first array has eltype other than `Float64`, e.g. `Int64` or `AbstractFloat`.

# Examples

`Int64`

```julia
julia> using NamedDims

julia> X = NamedDimsArray{(:a, :b)}([1 2; 3 4]);

julia> y = NamedDimsArray{(:a,)}([1, 2]);

julia> X \ y
ERROR: TypeError: in NamedFactorization, in F, expected F<:LinearAlgebra.Factorization{T}, got Type{LinearAlgebra.LU{Float64, Matrix{Float64}}}
Stacktrace:
 [1] lu(::NamedDimsArray{(:a, :b), Int64, 2, Matrix{Int64}}; kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ NamedDims ~/.julia/packages/NamedDims/V4xXA/src/functions_linearalgebra.jl:43
 [2] lu
   @ ~/.julia/packages/NamedDims/V4xXA/src/functions_linearalgebra.jl:42 [inlined]
 [3] \(A::NamedDimsArray{(:a, :b), Int64, 2, Matrix{Int64}}, B::NamedDimsArray{(:a,), Int64, 1, Vector{Int64}})
   @ LinearAlgebra /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/generic.jl:1136
 [4] top-level scope
   @ REPL[90]:1
```

`AbstractFloat` (about the simplest I could get starting from a real use case, without understanding how this works)

```julia
julia> X = NamedDimsArray{(:a, :b)}(AbstractFloat[24.16 1.39044; 25.16 1.39946; 25.3 1.38443]);

julia> y = NamedDimsArray{(:a,)}([28.94, 28.75, 28.7]);  # doesn't need to be AbstractFloat

julia> X \ y
ERROR: TypeError: in NamedFactorization, in F, expected F<:LinearAlgebra.Factorization{T}, got Type{LinearAlgebra.QRPivoted{Float64, Matrix{Float64}}}
Stacktrace:
 [1] qr(nda::NamedDimsArray{(:a, :b), AbstractFloat, 2, Matrix{AbstractFloat}}, args::Val{true}; kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ NamedDims ~/.julia/packages/NamedDims/V4xXA/src/functions_linearalgebra.jl:43
 [2] qr
   @ ~/.julia/packages/NamedDims/V4xXA/src/functions_linearalgebra.jl:42 [inlined]
 [3] \(A::NamedDimsArray{(:a, :b), AbstractFloat, 2, Matrix{AbstractFloat}}, B::NamedDimsArray{(:a,), Float64, 1, Vector{Float64}})
   @ LinearAlgebra /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/generic.jl:1138
 [4] top-level scope
   @ REPL[66]:1
```

# What to do

This is easy enough to avoid by using `Float64` eltype.

I'm not sure if there's a way to support an arbitrary number of `Factorization` parameters with one constructor.

Pull request:
closes #164

Patch:
diff --git a/Project.toml b/Project.toml
index dfb1dab..db72b44 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "NamedDims"
 uuid = "356022a1-0364-5f58-8944-0da4b18d706f"
 authors = ["Invenia Technical Computing Corporation"]
-version = "0.2.32"
+version = "0.2.33"
 
 [deps]
 AbstractFFTs = "621f4979-c628-5d54-868e-fcf4e3e8185c"
diff --git a/src/functions_linearalgebra.jl b/src/functions_linearalgebra.jl
index 09707a3..caa1c93 100644
--- a/src/functions_linearalgebra.jl
+++ b/src/functions_linearalgebra.jl
@@ -13,6 +13,10 @@ struct NamedFactorization{L, T, F<:Factorization{T}} <: Factorization{T}
     fact::F
 end
 
+function NamedFactorization{L}(fact::F) where {L, T, F <: Factorization{T}}
+    return NamedFactorization{L,T,F}(fact)
+end
+
 # Need parent to explicitly call `getfield` to avoid infinitely calling `getproperty`
 Base.parent(named::NamedFactorization) = getfield(named, :fact)
 Base.size(named::NamedFactorization) = size(parent(named))
@@ -39,8 +43,7 @@ end
 for func in (:lu, :lu!, :lq, :lq!, :svd, :svd!, :qr, :qr!)
     @eval begin
         function LinearAlgebra.$func(nda::NamedDimsArray{L, T}, args...; kwargs...) where {L, T}
-            fact = $func(parent(nda), args...; kwargs...)
-            return NamedFactorization{L, T, typeof(fact)}(fact)
+            return NamedFactorization{L}($func(parent(nda), args...; kwargs...))
         end
     end
 end


Test patch:
diff --git a/test/functions_linearalgebra.jl b/test/functions_linearalgebra.jl
index f7bfe32..454fece 100644
--- a/test/functions_linearalgebra.jl
+++ b/test/functions_linearalgebra.jl
@@ -135,3 +135,9 @@ end
         end
     end
 end
+
+@testset "#164 factorization eltype not same as input eltype" begin
+    # https://github.com/invenia/NamedDims.jl/issues/164
+    nda = NamedDimsArray{(:foo, :bar)}([1 2 3; 4 5 6; 7 8 9])  # Int eltype
+    @test qr(nda) isa NamedDims.NamedFactorization{(:foo, :bar), Float64}
+end