Root calculation failure

[Fabi293]

The following code fails:

using Polynomials

s = Polynomial([0,1],:s)

if true
	r = (15s^14 + 1e-16s^15)//(1s + 14s^14)
else
	r = (15s^14 + 1e-16s^15)//(s)
end

roots(r)

In the if case (if true), I get the error:

ERROR: DimensionMismatch: second dimension of A, 14, does not match first dimension of B, 13
Stacktrace:
[1] (::LinearAlgebra.var"#throw_dimerr#251")(::Matrix{Float64}, nA::Int64, mB::Int64)
@ LinearAlgebra $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\diagonal.jl:267
[2] _muldiag_size_check
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\diagonal.jl:269 [inlined]
[3] _muldiag_size_check
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\diagonal.jl:284 [inlined]
[4] _mul_diag!
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\diagonal.jl:423 [inlined]
[5] _mul!
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\diagonal.jl:430 [inlined]
[6] _mul!
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\bidiag.jl:434 [inlined]
[7] mul!
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\matmul.jl:285 [inlined]
[8] mul!
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\matmul.jl:253 [inlined]
[9] *
@ $USER.julia\juliaup\julia-1.11.3+0.x64.w64.mingw32\share\julia\stdlib\v1.11\LinearAlgebra\src\matmul.jl:114 [inlined]
[10] pejorative_root(p::Vector{Float64}, zs::Vector{Float64}, ls::Vector{Int64}; τ::Float64, maxsteps::Int64, verbose::Bool)
@ Polynomials.Multroot $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:275
[11] pejorative_root
@ $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:252 [inlined]
[12] #pejorative_root#5
@ $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:247 [inlined]
[13] pejorative_root
@ $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:244 [inlined]
[14] multroot(p::Polynomial{Float64, :s}; verbose::Bool, kwargs::@kwargs{method::Symbol})
@ Polynomials.Multroot $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:115
[15] multroot
@ $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:101 [inlined]
[16] #roots#187
@ $USER.julia\packages\Polynomials\893kl\src\rational-functions\common.jl:463 [inlined]
[17] roots(pq::RationalFunction{Float64, :s, Polynomial{Float64, :s}})
@ Polynomials $USER.julia\packages\Polynomials\893kl\src\rational-functions\common.jl:456
[18] top-level scope
@ $PATH:11

and in the else case (if false):

ERROR: ArgumentError: NaN in reduced polynomial
Stacktrace:
[1] pejorative_manifold(p::Polynomial{Float64, :s}; method::Symbol, θ::Float64, ρ::Float64, ϕ::Int64, kwargs::@kwargs{})
@ Polynomials.Multroot $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:158
[2] pejorative_manifold
@ $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:140 [inlined]
[3] multroot(p::Polynomial{Float64, :s}; verbose::Bool, kwargs::@kwargs{method::Symbol})
@ Polynomials.Multroot $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:114
[4] multroot
@ $USER.julia\packages\Polynomials\893kl\src\polynomials\multroot.jl:101 [inlined]
[5] #roots#187
@ $USER.julia\packages\Polynomials\893kl\src\rational-functions\common.jl:463 [inlined]
[6] roots(pq::RationalFunction{Float64, :s, Polynomial{Float64, :s}})
@ Polynomials $USER.julia\packages\Polynomials\893kl\src\rational-functions\common.jl:456
[7] top-level scope
@ $PATH:11

I can imagine that it fails due to numerical errors, but the error messages look kinda weird if this is the case.

[jverzani]

Thanks. I'll try and fix this in a few days.

[jverzani]

Yeah, your instinct is correct.

For each case, the lowest_terms method if first called. This gives

Second case:

(NaN + Inf*s) // (0.0)

So clearly that is problematic. I'll add a check for this case, it seems easy to do

The first case gives

(1.0714285714285714*s^13 + 7.142857142857142e-18*s^14) // (0.07142857142857142 + 1.0*s^13)

The roots function does work for the numerator, but the code path calls the multroot function (for aniticipated multiplicities) and this causes the failure. I have to see if I can add some check here.