Fix ldiv for Julia v1.6 (#46)

* Fix ldiv for Julia v1.6

* Fix

* Fix

Author: blegat

src/chebyshev.jl

@@ -104,7 +104,12 @@ function SA.coeffs(
     ext = SA.coeffs(algebra_element(cfs, source), extended)
     return SA.SparseCoefficients(
         sub.monomials,
-        LinearAlgebra.UpperTriangular(A) \ ext,
+        #LinearAlgebra.UpperTriangular(A) \ ext, # Julia v1.6 converts `A` to the eltype of the `result` which is bad for JuMP
+        LinearAlgebra.ldiv!(
+            zeros(_promote_coef(eltype(ext), Chebyshev), size(A, 2)),
+            LinearAlgebra.UpperTriangular(A),
+            ext,
+        ),
     )
 end
 

src/monomial.jl

@@ -118,7 +118,9 @@ function Base.copy(basis::SubBasis)
     return typeof(basis)(copy(basis.monomials))
 end
 
-Base.:(==)(a::SubBasis, b::SubBasis) = a.monomials == b.monomials
+function Base.:(==)(a::SubBasis{B}, b::SubBasis{B}) where {B}
+    return a.monomials == b.monomials
+end
 
 function algebra_type(::Type{BT}) where {B,M,BT<:MonomialIndexedBasis{B,M}}
     return Algebra{BT,B,M}

src/orthogonal.jl

@@ -105,11 +105,8 @@ function univariate_orthogonal_basis(
     end
 end
 
-function explicit_basis_covering(
-    ::FullBasis{B,M},
-    monos::SubBasis{<:AbstractMonomial,M},
-) where {B<:AbstractMultipleOrthogonal,M}
-    to_add = collect(monos.monomials)
+function _covering(::FullBasis{B,M}, monos) where {B,M}
+    to_add = collect(monos)
     m = Set{M}(to_add)
     while !isempty(to_add)
         mono = pop!(to_add)
@@ -125,7 +122,22 @@ function explicit_basis_covering(
             end
         end
     end
-    return SubBasis{B}(MP.monomial_vector(collect(m)))
+    return collect(m)
+end
+
+function explicit_basis_covering(
+    full::FullBasis{BM,M},
+    monos::SubBasis{B,M},
+) where {BM<:AbstractMonomial,B<:AbstractMultipleOrthogonal,M}
+    full = FullBasis{B,M}()
+    return SubBasis{BM}(_covering(full, monos.monomials))
+end
+
+function explicit_basis_covering(
+    full::FullBasis{B,M},
+    monos::SubBasis{<:AbstractMonomial,M},
+) where {B<:AbstractMultipleOrthogonal,M}
+    return SubBasis{B}(_covering(full, monos.monomials))
 end
 
 function _scalar_product_function(::Type{<:AbstractMultipleOrthogonal}, i::Int) end

test/runtests.jl

@@ -21,9 +21,23 @@ function _test_basis(basis)
 end
 
 struct TypeA end
-struct TypeB end
 Base.zero(::Type{TypeA}) = TypeA()
+Base.iszero(::TypeA) = false
+LinearAlgebra.adjoint(::TypeA) = TypeA()
+struct TypeB end
+Base.zero(::Type{TypeB}) = TypeB()
+Base.iszero(::TypeB) = false
+LinearAlgebra.adjoint(::TypeB) = TypeB()
+
 Base.:*(::Float64, ::TypeA) = TypeB()
+Base.:*(::TypeA, ::Float64) = TypeB()
+Base.:*(::Float64, ::TypeB) = TypeB()
+Base.:*(::TypeB, ::Float64) = TypeB()
+Base.:/(::TypeA, ::Float64) = TypeB()
+Base.:/(::TypeB, ::Float64) = TypeB()
+Base.:+(::TypeB, ::TypeB) = TypeB()
+Base.:-(::TypeB, ::TypeB) = TypeB()
+Base.convert(::Type{TypeB}, ::TypeA) = TypeB()
 
 function api_test(B::Type{<:MB.AbstractMonomialIndexed}, degree)
     @polyvar x[1:2]
@@ -66,6 +80,15 @@ function api_test(B::Type{<:MB.AbstractMonomialIndexed}, degree)
         a = MB.algebra_element(ones(length(basis)), basis)
         _test_op(MB.implicit, a)
         @test SA.star(a) == a
+        if B == Chebyshev || B == ScaledMonomial
+            mono = explicit_basis_covering(
+                FullBasis{Monomial,eltype(basis.monomials)}(),
+                basis,
+            )
+            a = MB.algebra_element(fill(TypeA(), length(basis)), mono)
+            @test SA.coeffs(a, full_basis).values ==
+                  fill(TypeB(), length(basis))
+        end
     end
     mono = x[1]^2 * x[2]^3
     p = MB.Polynomial{B}(mono)