Skip to content

Commit

Permalink
error with NaN; allow method to pass through
Browse files Browse the repository at this point in the history
  • Loading branch information
jverzani committed Feb 5, 2025
1 parent 9b4af66 commit db44627
Show file tree
Hide file tree
Showing 5 changed files with 27 additions and 9 deletions.
2 changes: 1 addition & 1 deletion Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ name = "Polynomials"
uuid = "f27b6e38-b328-58d1-80ce-0feddd5e7a45"
license = "MIT"
author = "JuliaMath"
version = "4.0.14"
version = "4.0.15"

[deps]
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Expand Down
6 changes: 2 additions & 4 deletions src/polynomials/multroot.jl
Original file line number Diff line number Diff line change
Expand Up @@ -145,7 +145,6 @@ function pejorative_manifold(
ϕ = 100, # residual growth factor
kwargs...
) where {T,X}

S = float(T)
u = convert(PnPolynomial{S,X}, coeffs0(p))
nu₂ = norm(u, 2)
Expand All @@ -156,9 +155,9 @@ function pejorative_manifold(
atol = ρ2, rtol = zero(ρ2))
ρⱼ /= nu₂

hasnan(v) && throw(ArgumentError("NaN in reduced polynomial"))
# root approximations
zs = roots(v)

# recover multiplicities
ls = pejorative_manifold_multiplicities(Val(method),
u, v, w,
Expand Down Expand Up @@ -206,7 +205,6 @@ function pejorative_manifold_multiplicities(
end

end

ls

end
Expand All @@ -225,6 +223,7 @@ function pejorative_manifold_multiplicities(

dv = derivative(v)
ls = w.(zs) ./ dv.(zs)

ls = round.(Int, real.(ls))

return ls
Expand Down Expand Up @@ -258,7 +257,6 @@ function pejorative_root(p, zs::Vector{S}, ls;

## Solve WJ Δz = W(Gl(z) - a)
## using weights min(1/|aᵢ|), i ≠ 1

m,n = sum(ls), length(zs)
# storage
a = p[2:end]./p[1] # a ~ (p[n-1], p[n-2], ..., p[0])/p[n]
Expand Down
1 change: 0 additions & 1 deletion src/polynomials/ngcd.jl
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,6 @@ function ngcd(p::P, q::Q,
p′ = P′{R,X}(ps[nz:end])
q′ = P′{R,X}(qs[nz:end])
out = NGCD.ngcd(p′, q′, args...; kwargs...)

## convert to original polynomial type
𝑷 = Polynomials.constructorof(P){R,X}
u,v,w = convert.(𝑷, (out.u,out.v,out.w))
Expand Down
13 changes: 10 additions & 3 deletions src/rational-functions/common.jl
Original file line number Diff line number Diff line change
Expand Up @@ -427,15 +427,22 @@ end

## ---- zeros, poles, ...
"""
poles(pq::AbstractRationalFunction; method=:numerical, kwargs...)
poles(pq::AbstractRationalFunction;
method=:numerical, multroot_method=:direct, kwargs...)
For a rational function `p/q`, first reduces to normal form, then finds the roots and multiplicities of the resulting denominator.
* `method` is used to pass to `lowest_terms`
* `multroot_method` is passed to the method argument of `multroot`, which can be `:direct` (the faster default) or `:iterative` (the slower, and possibly more robust alternate)
"""
function poles(pq::AbstractRationalFunction; method=:numerical, kwargs...)
function poles(pq::AbstractRationalFunction;
method=:numerical, # for lowest_terms
multroot_method=:direct, # or :iterative
kwargs...)
pq′ = lowest_terms(pq; method=method, kwargs...)
den = denominator(pq′)
mr = Multroot.multroot(den)
mr = Multroot.multroot(den; method=multroot_method)
(zs=mr.values, multiplicities = mr.multiplicities)
end

Expand Down
14 changes: 14 additions & 0 deletions test/StandardBasis.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1027,6 +1027,20 @@ end
out = Polynomials.Multroot.multroot(Polynomials.Polynomial(pb))
@test out.values [1.0] && out.multiplicities == [2]

## Issue #592
p = Polynomial([NaN,NaN,NaN,NaN,NaN,NaN,NaN])
@test_throws ArgumentError roots(p//p)

## issue #593
numcoeffs = Complex{BigFloat}[-0.0 + 0.0im, -0.0 + 0.0im, -0.0 + 0.0im, -0.0 + 0.0im, -0.0 + 0.0im, -0.0 + 0.0im, -6.4095e+257 - 3.314e+243im, -5.6118e+244 - 3.9514e+230im, -4.0102e+236 - 2.0735e+222im, -3.0725e+223 - 2.1631e+209im, -1.0975e+215 - 5.6751e+200im, -7.2072e+201 - 5.0751e+187im, -1.7167e+193 - 8.877e+178im, -9.3936e+179 - 6.6134e+165im, -1.678e+171 - 8.6782e+156im, -7.3448e+157 - 5.1693e+143im, -1.0499e+149 - 5.43e+134im, -3.4468e+135 - 2.4259e+121im, -4.1052e+126 - 2.1231e+112im, -8.9846e+112 - 6.3225e+98im, -9.1725e+103 - 4.742e+89im, -1.0036e+90 - 7.0656e+75im, -8.9646e+80 - 4.6354e+66im]
dencoeffs = Complex{BigFloat}[0.0 + 0.0im, 0.0 + 0.0im, -0.0 + 0.0im, -0.0 + 0.0im, 0.0 + 0.0im, 0.0 + 0.0im, -3.0837e+258 + 0.0im, -3.5997e+245 - 9.309e+230im, -1.9292e+237 - 7.0217e+217im, -1.9709e+224 - 5.0951e+209im, -5.2803e+215 - 3.2932e+196im, -4.6214e+202 - 1.1952e+188im, -8.2587e+193 - 6.4383e+174im, -6.0237e+180 - 1.5577e+166im, -8.0763e+171 - 6.7125e+152im, -4.711e+158 - 1.2183e+144im, -5.0507e+149 - 3.937e+130im, -2.2106e+136 - 5.7173e+121im, -1.9752e+127 - 1.2316e+108im, -5.7628e+113 - 1.4904e+99im, -4.4151e+104 - 1.6057e+85im, -6.4393e+90 - 1.6651e+76im, -4.3167e+81 + 0.0im]
p = Polynomial(numcoeffs)
q = Polynomial(dencoeffs)
r = p//q

@test_throws ArgumentError poles(r)
out = poles(r; multroot_method=:iterative)
@test out.multiplicities == [3]
end

@testset "critical points" begin
Expand Down

0 comments on commit db44627

Please sign in to comment.