diff --git a/Project.toml b/Project.toml index 7388759d..5b1df0e1 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "FiniteDifferences" uuid = "26cc04aa-876d-5657-8c51-4c34ba976000" -version = "0.9.3" +version = "0.9.4" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" diff --git a/src/methods.jl b/src/methods.jl index 20b56ea9..41485561 100644 --- a/src/methods.jl +++ b/src/methods.jl @@ -157,10 +157,12 @@ end # Compute coefficients for the method function _coefs(grid::AbstractVector{<:Real}, p::Integer, q::Integer) - C = [g^i for i in 0:(p - 1), g in grid] - x = zeros(Int, p) + # For high precision on the \ we use Rational, and to prevent overflows we use Int128 + # At the end we go to Float64 for fast floating point math (rather than rational math) + C = [Rational{Int128}(g^i) for i in 0:(p - 1), g in grid] + x = zeros(Rational{Int128}, p) x[q + 1] = factorial(q) - return C \ x + return Float64.(C \ x) end # Estimate the bound on the function value and its derivatives at a point diff --git a/test/methods.jl b/test/methods.jl index 4b46d33e..fe886626 100644 --- a/test/methods.jl +++ b/test/methods.jl @@ -37,6 +37,17 @@ using FiniteDifferences: Forward, Backward, Central, Nonstandard @test central_fdm(5, 1)(abs, 0.001) ≈ 1.0 end + @testset "Accuracy at high orders, with high adapt" begin + # Regression test against issues with precision during computation of _coeffs + # see https://github.com/JuliaDiff/FiniteDifferences.jl/issues/64 + + @test fdm(central_fdm(9, 5), exp, 1.0, adapt=4) ≈ exp(1) atol=1e-7 + + poly(x) = 4x^3 + 3x^2 + 2x + 1 + @test fdm(central_fdm(9, 3), poly, 1.0, adapt=4) ≈ 24 atol=1e-11 + end + + @testset "Printing FiniteDifferenceMethods" begin @test sprint(show, central_fdm(2, 1)) == """ FiniteDifferenceMethod: