JuliaIntervals / IntervalRootFinding.jl

Library for finding the roots of a function using interval arithmetic
https://juliaintervals.github.io/IntervalRootFinding.jl/
Other
125 stars 26 forks source link

Unexpected answer to trigonometric polynomial #195

Closed josephcbradley closed 2 months ago

josephcbradley commented 9 months ago

I was looking at the following case and got unexpected behaviour:

julia> using IntervalRootFinding

julia> f(θ) = 2*cos(θ)^2 - 5*cos(θ)^2 + 3
f (generic function with 1 method)

julia> domain = 0..2π
[0, 6.28319]

julia> roots(f, domain)
17-element Vector{Root{Interval{Float64}}}:
 Root([6.28318, 6.28319], :unknown)
 Root([6.28318, 6.28319], :unknown)
 Root([3.14159, 3.1416], :unknown)
 Root([0, 7.63551e-08], :unknown)
 Root([6.25383e-07, 6.78054e-07], :unknown)
 Root([3.14159, 3.1416], :unknown)
 Root([3.14159, 3.1416], :unknown)
 Root([6.28318, 6.28319], :unknown)
 Root([3.14153, 3.14154], :unknown)
 Root([6.28318, 6.28319], :unknown)
 Root([6.28318, 6.28319], :unknown)
 Root([6.28318, 6.28319], :unknown)
 Root([3.14159, 3.1416], :unknown)
 Root([1.53912e-07, 1.6754e-07], :unknown)
 Root([3.14159, 3.1416], :unknown)
 Root([7.6355e-08, 8.38739e-08], :unknown)
 Root([3.10248e-07, 3.36623e-07], :unknown)

Roots.jl handles things as expected:

julia> using Roots

julia> find_zeros(f, (0, 2π))
3-element Vector{Float64}:
 0.0
 3.1415926437586927
 6.283185307179586

IntervalRootFinding.jl has actually found all of the zeros at the multiples of pi, but all with :unknown flags and lots of duplicates. Is there any reason for this?

Session info:

julia> versioninfo()
Julia Version 1.9.3
Commit bed2cd540a1 (2023-08-24 14:43 UTC)
Build Info:
  Official https://julialang.org/ release
Platform Info:
  OS: macOS (x86_64-apple-darwin22.4.0)
  CPU: 16 × Intel(R) Core(TM) i9-9880H CPU @ 2.30GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-14.0.6 (ORCJIT, skylake)
  Threads: 8 on 16 virtual cores
Environment:
  JULIA_NUM_THREADS = 8
  JULIA_IMAGE_THREADS = 1

(@v1.9) pkg> status
Status `~/.julia/environments/v1.9/Project.toml`
  [cbdf2221] AlgebraOfGraphics v0.6.16
  [6e4b80f9] BenchmarkTools v1.3.2
  [336ed68f] CSV v0.10.11
⌅ [052768ef] CUDA v4.4.1
  [324d7699] CategoricalArrays v0.10.8
  [aaaa29a8] Clustering v0.15.4
  [5ae59095] Colors v0.12.10
  [a93c6f00] DataFrames v1.6.1
  [1313f7d8] DataFramesMeta v0.14.0
  [864edb3b] DataStructures v0.18.15
  [8bb1440f] DelimitedFiles v1.9.1
  [77a26b50] DiffEqNoiseProcess v5.19.0
  [0c46a032] DifferentialEquations v7.11.0
  [31c24e10] Distributions v0.25.102
  [634d3b9d] DrWatson v2.13.0
  [48062228] FilePathsBase v0.9.21
  [e9fadd82] FilteredGraphs v0.1.0 `~/.julia/dev/FilteredGraphs`
⌅ [587475ba] Flux v0.13.17
  [38e38edf] GLM v1.9.0
  [e9467ef8] GLMakie v0.8.10
  [bd48cda9] GraphRecipes v0.5.12
  [86223c79] Graphs v1.9.0
  [cd3eb016] HTTP v1.10.0
  [f5425a1d] Hurst v0.1.4
  [916415d5] Images v0.26.0
  [d2bf35a9] IntervalRootFinding v0.5.11
  [98e50ef6] JuliaFormatter v1.0.37
  [2b0e0bc5] LanguageServer v4.4.0
  [2fda8390] LsqFit v0.13.0
  [945b72a4] MarketData v0.13.12
  [626554b9] MetaGraphs v0.7.2
  [46757867] NetworkLayout v0.4.5
  [6fd5a793] Octavian v0.3.27
  [6fe1bfb0] OffsetArrays v1.12.10
  [a15396b6] OnlineStats v1.6.3
  [e9f21f70] OpenAI v0.8.5
  [90014a1f] PDMats v0.11.21
  [9b87118b] PackageCompiler v2.1.9
  [98572fba] Parquet2 v0.2.19
  [b98c9c47] Pipe v1.3.0
  [14b8a8f1] PkgTemplates v0.7.45
  [f0f68f2c] PlotlyJS v0.18.10
  [91a5bcdd] Plots v1.39.0
  [c3e4b0f8] Pluto v0.19.29
  [27ebfcd6] Primes v0.5.4
  [c46f51b8] ProfileView v1.7.2
  [92933f4c] ProgressMeter v1.9.0
  [ce6b1742] RDatasets v0.7.7
  [295af30f] Revise v3.5.6
  [f2b01f46] Roots v2.0.20
  [fc66bc1b] SNAPDatasets v0.2.0
  [47aef6b3] SimpleWeightedGraphs v1.4.0
⌅ [2913bbd2] StatsBase v0.33.21
  [f3b207a7] StatsPlots v0.15.6
  [ab02a1b2] TableOperations v1.2.0
  [bd369af6] Tables v1.11.0
  [ac1d9e8a] ThreadsX v0.1.11
  [9e3dc215] TimeSeries v0.23.2
  [21ca0261] Transformers v0.2.8
  [e88e6eb3] Zygote v0.6.65
dpsanders commented 9 months ago

This is because the roots of this function are multiple roots: the function touches the line $y = 0$, rather than crossing it. Interval methods are then unable to prove that there is a unique root there.