search

JuliaAlgebra / MultivariatePolynomials.jl

Multivariate polynomials interface
https://juliaalgebra.github.io/MultivariatePolynomials.jl/stable/
Other
135 stars 27 forks source link

oneunit(T) fails for some types #245

Closed KlausC closed 1 year ago

KlausC commented 1 year ago 
julia> using DynamicPolynomials
julia> using LinearAlgebra
julia> @polyvar x y
(x, y)
julia> M = [x y; 1 x]
2×2 Matrix{Term{true, Int64}}:
 x  y
 1  x

julia> oneunit(typeof(1 / x))
ERROR: MethodError: no method matching RationalPoly{Term{true, Int64}, PolyVar{true}}(::RationalPoly{Term{true, Int64}, Monomial{true}})
Closest candidates are:
  RationalPoly{NT, DT}(::Any, ::Any) where {NT<:(AbstractPolynomialLike), DT<:(AbstractPolynomialLike)} at ~/.julia/packages/MultivariatePolynomials/sWAOE/src/rational.jl:5
  RationalPoly{NT, DT}(::Bool) where {NT, DT} at ~/.julia/packages/MultivariatePolynomials/sWAOE/src/rational.jl:10
Stacktrace:
 [1] oneunit(#unused#::Type{RationalPoly{Term{true, Int64}, PolyVar{true}}})
   @ Base ./number.jl:358
 [2] top-level scope
   @ REPL[35]:1

julia> inv(M)
ERROR: MethodError: no method matching RationalPoly{Polynomial{true, Int64}, Term{true, Int64}}(::RationalPoly{Polynomial{true, Int64}, Term{true, Int64}})
Closest candidates are:
  RationalPoly{NT, DT}(::Any, ::Any) where {NT<:(AbstractPolynomialLike), DT<:(AbstractPolynomialLike)} at ~/.julia/packages/MultivariatePolynomials/sWAOE/src/rational.jl:5
  RationalPoly{NT, DT}(::Bool) where {NT, DT} at ~/.julia/packages/MultivariatePolynomials/sWAOE/src/rational.jl:10
Stacktrace:
 [1] oneunit(#unused#::Type{RationalPoly{Polynomial{true, Int64}, Term{true, Int64}}})
   @ Base ./number.jl:358
 [2] lutype(T::Type)
   @ LinearAlgebra ~/julia-1.8.5/share/julia/stdlib/v1.8/LinearAlgebra/src/lu.jl:204
 [3] lu(A::Matrix{RationalPoly{Polynomial{true, Int64}, Term{true, Int64}}}, pivot::RowMaximum; check::Bool)
   @ LinearAlgebra ~/julia-1.8.5/share/julia/stdlib/v1.8/LinearAlgebra/src/lu.jl:279
 [4] lu(A::Matrix{RationalPoly{Polynomial{true, Int64}, Term{true, Int64}}}, pivot::RowMaximum) (repeats 2 times)
   @ LinearAlgebra ~/julia-1.8.5/share/julia/stdlib/v1.8/LinearAlgebra/src/lu.jl:278
 [5] inv(A::Matrix{Term{true, Int64}})
   @ LinearAlgebra ~/julia-1.8.5/share/julia/stdlib/v1.8/LinearAlgebra/src/dense.jl:893
 [6] top-level scope
   @ REPL[34]:1
blegat commented 1 year ago

Thanks for reporting this, I'm suprised the Base method does not call convert instead of the constructor. It should be fixed by https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/pull/246 Not that inv will still not work because abs is not defined. You should use the Cramer rule. We could add an implementation in this package

KlausC commented 1 year ago

IMO Cramer's rule is not worth while. In my experience it is sufficient to define abs(p::RationalPoly) = float(!iszero(p)) or similar to allow for "pivotized" LU decompostion.

blegat commented 1 year ago

Yes, that could also work, we can't define this here since it would lead to confusing behavior for other use cases though