Integer overflow in innocent looking linear interpolation. Should some warning be issued?

[hurak]

Overflow is encountered during and innocent looking linear interpolation/extrapolation:

Interpolate linearly two points (1,1) and (10000,10000) and then find the value at 10001 by linear extrapolation. The following code does it

julia> itpi = LinearInterpolation([1,10000],[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
Ratios.SimpleRatio{Int64}(99980001, 99980001)
Ratios.SimpleRatio{Int64}(999800010000, 99980001)

julia> convert(Float64,itpi(10001)) 773.9376918087789

The result is clearly incorrect because the correct answer is 10001. Looking at the displayed ratios, I guess that an integer overflow is what happened here. There are no problems with smaller numbers:
```julia
julia> itp = LinearInterpolation([1,2],[1,2]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
 Ratios.SimpleRatio{Int64}(1, 1)
 Ratios.SimpleRatio{Int64}(2, 1)

julia> convert(Float64,itp(3))
3.0

or with the floating point data:

julia> itpf = LinearInterpolation([1.0,10000.0],[1.0,10000.0]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Float64},), ::Vector{Float64}, Gridded(Linear())), Line()) with element type Float64:
     1.0
 10000.0

julia> convert(Float64,itpf(10001.0))
10001.0

I am hesitating if this should be regarded (and filed) as an Issue. In general the responsibility for the overflow problems is on the application programmer, right? But shouldn’t the interpolation function give some warning here? The parameters of the problem look innocent, don’t they?

[mkitti]

If there is a bug here, it is that overflow occurs for such a simple case. If the user supplies an integer type, we stick with that integer type for the computation to take advantage of any performance gains along with any potential issues such as overflow.

Let's distill the bug to be a failure to simplify the SimpleRatio which makes us prematurely vulnerable to overflow in this case.

[timholy]

Can someone check whether it would be fixed with https://github.com/timholy/Ratios.jl/pull/19#issuecomment-894860010?

[mkitti]

Can someone check whether it would be fixed with timholy/Ratios.jl#19 (comment)?

julia> using Interpolations, Ratios

julia> itpi = LinearInterpolation([1,10000],[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
 SimpleRatio{Int64}(99980001, 99980001)
 SimpleRatio{Int64}(999800010000, 99980001)

julia> itpi(10001)
SimpleRatio{Int64}(7736281631652211921, 9996000599960001)

julia> convert(Float64, itpi(10001))
773.9376918087789

julia> Base.:+(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num + y.num, x.den) : SimpleRatio(x.num*y.den + x.den*y.num, x.den*y.den)

julia> Base.:-(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num - y.num, x.den) : SimpleRatio(x.num*y.den - x.den*y.num, x.den*y.den)

julia> itpi = LinearInterpolation([1,10000],[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
 SimpleRatio{Int64}(9999, 9999)
 SimpleRatio{Int64}(99990000, 9999)

julia> convert(Float64, itpi(10001))
10001.0

[mkitti]

Within Interpolations.jl, my thought would be to simplify the ratio before construction:

julia> using Interpolations, Ratios

julia> Interpolations.ratio(num::Integer, denom::Integer) = SimpleRatio(Base.divgcd(promote(num, denom)...)...)

julia> itpi = LinearInterpolation([1,10000],[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
 SimpleRatio{Int64}(1, 1)
 SimpleRatio{Int64}(10000, 1)

julia> itpi(10001)
SimpleRatio{Int64}(999899990001, 99980001)

julia> convert(Float64, itpi(10001))
10001.0

This also combines well with the above fix for Ratios.jl:

julia> Base.:+(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num + y.num, x.den) : SimpleRatio(x.num*y.den + x.den*y.num, x.den*y.den)

julia> Base.:-(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num - y.num, x.den) : SimpleRatio(x.num*y.den - x.den*y.num, x.den*y.den)

julia> itpi(10001)
SimpleRatio{Int64}(99999999, 9999)

Perhaps a simplify method would be generally useful?

julia> simplify(r::SimpleRatio) = (d = gcd(r.num, r.den); SimpleRatio(r.num ÷ d, r.den ÷ d))
simplify (generic function with 1 method)

julia> itpi(10001)
SimpleRatio{Int64}(99999999, 9999)

julia> simplify(itpi(10001))
SimpleRatio{Int64}(10001, 1)

At some point, it seems we are patching SimpleRatio until it becomes Rational.

[timholy]

At some point, it seems we are patching SimpleRatio until it becomes Rational.

Yes, and that would be a huge problem.

julia> using BenchmarkTools

julia> @btime Base.divgcd(7736281631652211921, 9996000599960001)
  127.862 ns (0 allocations: 0 bytes)
(7736281631652211921, 9996000599960001)

vs

julia> using Interpolations

julia> using Interpolations.Ratios

julia> Base.:+(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num + y.num, x.den) : SimpleRatio(x.num*y.den + x.den*y.num, x.den*y.den)

julia> Base.:-(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num - y.num, x.den) : SimpleRatio(x.num*y.den - x.den*y.num, x.den*y.den)

julia> itpi = LinearInterpolation([1,10000],[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
 SimpleRatio{Int64}(9999, 9999)
 SimpleRatio{Int64}(99990000, 9999)

julia> @btime $itpi(10001)
  7.331 ns (0 allocations: 0 bytes)
SimpleRatio{Int64}(99999999, 9999)

Do you see why I suggested that strategy first?

[mkitti]

Do you see why I suggested that strategy first?

It seems the denominator check seems to address this specific issue and is not too onerous. Do you want me to create an issue over at Ratios.jl? Would ifelse be better than a ternary operator to avoid branching logic?

[mkitti]

Regarding a warning or error that overflow occurred, perhaps working on the composition of this package, Ratios.jl, and SaferIntegers.jl would be useful.

Here's my attempt which requires some hot patching at the moment:

julia> using SaferIntegers, Interpolations, Ratios

julia> import SaferIntegers: ndigits0z, baseint, SaferIntegers

julia> for S in (:SafeSigned, :SafeUnsigned)
         @eval begin
               ndigits0z(x::T) where T<:$S = ndigits0z(baseint(x)) # do not reconvert
           ndigits0z(x::T1, b::T2) where T1<:$S where T2<:SafeInteger = ndigits0z(baseint(x), baseint(b)) # do not reconvert
           ndigits0z(x::T1, b::T2) where T1<:$S where T2<:Integer = ndigits0z(baseint(x), b) # do not reconvert
         end
       end

julia> Base.float(x::SafeInteger) = Base.float(SaferIntegers.baseint(x))

julia> Base.:-(x::SimpleRatio{T}) where {T<:SafeSigned} = SimpleRatio(-x.num, x.den)

julia> itpi = LinearInterpolation(SafeInt[1,10000],SafeInt[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{SafeInt64},), ::Vector{SafeInt64}, Gridded(Linear())), Line()) with element type Float64:
 SimpleRatio{SafeInt64}(99980001, 99980001)
 SimpleRatio{SafeInt64}(999800010000, 99980001)

julia> itpi(10001)
ERROR: OverflowError: 999800010000 * 99980001 overflowed for type Int64
Stacktrace:
 [1] throw_overflowerr_binaryop(op::Symbol, x::Int64, y::Int64)
   @ Base.Checked .\checked.jl:154
 [2] checked_mul
   @ .\checked.jl:288 [inlined]
 [3] *
   @ ~\.julia\packages\SaferIntegers\nOUgY\src\arith_ops.jl:21 [inlined]
 [4] +
   @ ~\.julia\packages\Ratios\uRs4y\src\Ratios.jl:31 [inlined]
 [5] extrapolate_axis
   @ ~\.julia\dev\Interpolations\src\extrapolation\extrapolation.jl:168 [inlined]
 [6] extrapolate_value
   @ ~\.julia\dev\Interpolations\src\extrapolation\extrapolation.jl:162 [inlined]
 [7] (::Interpolations.Extrapolation{Float64, 1, Interpolations.GriddedInterpolation{Float64, 1, SafeInt64, Gridded{Linear{Throw{OnGrid}}}, Tuple{Vector{SafeInt64}}}, Gridded{Linear{Throw{OnGrid}}}, Line{Nothing}})(x::Int64)
   @ Interpolations ~\.julia\dev\Interpolations\src\extrapolation\extrapolation.jl:52
 [8] top-level scope
   @ REPL[7]:1

[timholy]

I think before making any decisions, all the options should be benchmarked, including against the (buggy) implementation we have now. SaferIntegers is appealing in principle, but what are the performance implications? Interpolations is an extraordinarily performance-sensitive package.

[mkitti]

The SaferIntegers approach is orthogonal to the other approaches. Once we provide support for it, the user can make a choice whether to use a SafeInteger to detect an overflow or not.

The gcd approach is probably better handled by Rational. Perhaps we should have a Preferences.jl mechanism to use Rational instead of SimpleRatio when performance is not the primary need of the user?

To summarize, I think we have three orthogonal approaches:

1. Modify Ratios.jl to check if SimpleRatios has the same denominator during addition and subtraction. The added overhead of that check should be benchmarked carefully. Based on my informal experiments and the above times, the overhead looks reasonable to me.

2. Modify SaferIntegers.jl and Ratios.jl to work together with Interpolations.jl. While benchmarks would be informative, I do not think this influences whether we should do this or not since the choice to use SaferIntegers.jl lies with the user.

3. Modify Interpolations.jl to use Preferences.jl to modify the implementation of Interpolations.ratio between SimpleRatio or Rational.

[JeffreySarnoff]

How would SaferIntegers.jl need to be modified to work with Interpolations.jl?

[mkitti]

How would SaferIntegers.jl need to be modified to work with Interpolations.jl?

Hi @JeffreySarnoff , adding Base.float for SafeInteger would help.

Base.float(x::SafeInteger) = Base.float(SaferIntegers.baseint(x))

See https://github.com/JeffreySarnoff/SaferIntegers.jl/issues/19

[JeffreySarnoff]

I incorporated this

Base.float(x::S) where S<:SafeSigned = Base.float(baseint(x))
Base.float(x::U) where U<:SafeUnsigned = Base.float(baseint(x))

into SaferIntegers/src/promote.jl, bumped the version and reregistered it (I still need to add tests)

[mkitti]

On approach #2, we are then here:

julia> using Interpolations, SaferIntegers

julia> itpi = LinearInterpolation(SafeInt[1,10000],SafeInt[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{SafeInt64},), ::Vector{SafeInt64}, Gridded(Linear())), Line()) with element type Float64:
 Ratios.SimpleRatio{SafeInt64}(99980001, 99980001)
 Ratios.SimpleRatio{SafeInt64}(999800010000, 99980001)

julia> itpi(10001)
ERROR: MethodError: no method matching -(::Ratios.SimpleRatio{SafeInt64})
...

We just need to add the missing method to Ratios now since SafeInt64 is not Signed, but SafeSigned.

julia> using Interpolations.Ratios

julia> Base.:-(x::SimpleRatio{T}) where {T<:SafeSigned} = SimpleRatio(-x.num, x.den)

julia> itpi(10001)
ERROR: OverflowError: 999800010000 * 99980001 overflowed for type Int64
Stacktrace:
 [1] throw_overflowerr_binaryop(op::Symbol, x::Int64, y::Int64)
   @ Base.Checked .\checked.jl:154
 [2] checked_mul
   @ .\checked.jl:288 [inlined]
 [3] *
   @ ~\.julia\packages\SaferIntegers\zbJSp\src\arith_ops.jl:21 [inlined]
 [4] +
   @ ~\.julia\packages\Ratios\fgO34\src\Ratios.jl:35 [inlined]
 [5] extrapolate_axis
   @ ~\.julia\dev\Interpolations\src\extrapolation\extrapolation.jl:168 [inlined]
 [6] extrapolate_value
   @ ~\.julia\dev\Interpolations\src\extrapolation\extrapolation.jl:162 [inlined]
 [7] (::Interpolations.Extrapolation{Float64, 1, Interpolations.GriddedInterpolation{Float64, 1, SafeInt64, Gridded{Linear{Throw{OnGrid}}}, Tuple{Vector{SafeInt64}}}, Gridded{Linear{Throw{OnGrid}}}, Line{Nothing}})(x::Int64)
   @ Interpolations ~\.julia\dev\Interpolations\src\extrapolation\extrapolation.jl:52
 [8] top-level scope
   @ REPL[7]:1

[JeffreySarnoff]

We just need to add the missing method to Ratios since SafeInt64 is not Signed, but SafeSigned

This structuring

SafeInt64 <: Integer
SafeInt64 <: SafeSigned
!(SafeInt64 <: Signed)

Was borrowed from an early version of RoundingIntegers.jl by @timholy . Improvements + amendations inside Julia have allowed that to be reworked

export RSigned, RUnsigned, RInteger
abstract type RSigned   <: Signed end
abstract type RUnsigned <: Unsigned end

const RInteger = Union{RSigned,RUnsigned}

which recovers the desired, dispatch simplifying, subtyping relationships.

I should apply this approach to SaferIntegers -- unless there is something about the development of SaferIntegers that precluded doing so years ago. Honestly, I do recall considering it ... and do not recall precisely why it did not serve. There is a structural distinction in the implementation of RoundingIntegers.jl when compared to SaferIntegers.jl: less complexity enures to its overarching logical schema.

[mkitti]

I should apply this approach to SaferIntegers

I gave it a shot: https://github.com/JeffreySarnoff/SaferIntegers.jl/pull/20

[mkitti]

You should be able to detect overflow now without modification to any of the packages.

julia> using SaferIntegers, Interpolations

julia> itpi = LinearInterpolation(SafeInt[1,10000],SafeInt[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{SafeInt64},), ::Vector{SafeInt64}, Gridded(Linear())), Line()) with element type Float64:
 SimpleRatio{SafeInt64}(99980001, 99980001)
 SimpleRatio{SafeInt64}(999800010000, 99980001)

julia> itpi(10001)
ERROR: OverflowError: 999800010000 * 99980001 overflowed for type Int64
Stacktrace:

[JeffreySarnoff]

yes, that works for me using Julia 17ish

[stevengj]

If the user supplies an integer type, we stick with that integer type for the computation to take advantage of any performance gains along with any potential issues such as overflow.

This doesn't make sense to me.

It seems like the interpolation inputs should essentially be promoted to the type of the output, which is a floating-point type for integer inputs.

Significant performance gains from keeping in Int64 form seem unlikely; floating-point arithmetic is fast (we no longer live in the 1970s–80s). And the danger of catastrophic overflow is real — correctness beats small performance gains.

[mkitti]

At the moment the interpolations in this package return a SimpleRatio from Ratios.jl rather than a float when given an integer. That is a design decision which preceeds me. Tim would be better able to comment on this than I.

Integer based output from integer input seems reasonable to me, although I can see how some may not expect it or specifically SimpleRatio. The user could select float output by giving float input, just as they will get unitful output from unitful input.

I agree there may be distinct expectations of this package. For this reason, I am considering using Preferences.jl to make Interpolations.ratio configurable at compile time.

[timholy]

Not just integers, either. The idea of using a ratio comes from

julia> 3//4 * 1.0
0.75

julia> 3//4 * 1.0f0
0.75f0

This allows the correct type to be returned from itp(x, y). It's much harder to do that if you pick a particular floating-point type.

That said, we should still do the den trick if the performance implications are not too bad. I forget whether this has been carefully benchmarked (I have not done it seriously myself). Detecting overflow is all fine and good, but avoiding it is even better.

[mkitti]

I was thinking that perhaps we should create a new ratio type for the common denominator detection.

[timholy]

It looks like Interpolations is the only consumer:

julia> using RegistryCompatTools

help?> held_back_by
search: held_back_by held_back_packages print_held_back

  held_back_by(pkgname::AbstractString, version::VersionNumber)

  Return a list of packages would hold back the given not-yet-registered version of pkgname.

  ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────

  held_back_by(pkgname::AbstractString, d=held_back_packages())

  Return a list of packages that are holding back `pkgname`. `d` can be obtained from [`held_back_packages`](@ref).

julia> held_back_by("Ratios", v"1.0")
1-element Vector{String}:
 "Interpolations"

And it looks like it has no performance impact (the following must be run while using Revise):

julia> using Interpolations, BenchmarkTools
[ Info: Precompiling Interpolations [a98d9a8b-a2ab-59e6-89dd-64a1c18fca59]

julia> itp1 = interpolate(rand(15), BSpline(Quadratic(Flat(OnGrid()))));

julia> itp2 = interpolate(rand(15, 12), BSpline(Quadratic(Periodic(OnGrid()))));

julia> itp2l = interpolate(rand(15, 12), BSpline(Linear()));

julia> using Ratios

julia> SimpleRatio(2, 5) + SimpleRatio(1, 5)
SimpleRatio{Int64}(15, 25)

julia> @benchmark itp1(x) setup=(x=3+3*rand())
BenchmarkTools.Trial: 10000 samples with 995 evaluations.
 Range (min … max):  27.977 ns … 705.963 ns  ┊ GC (min … max): 0.00% … 95.72%
 Time  (median):     28.592 ns               ┊ GC (median):    0.00%
 Time  (mean ± σ):   30.399 ns ±  17.046 ns  ┊ GC (mean ± σ):  1.45% ±  2.52%

   ▄█▆▂                             ▂▂ ▁▃▁                     ▁
  █████▇▆▅▄▅▃▃▄▄▅▄▅▅▆▅▂▃▃▄▃▃▃▃▂▅██▆████████▆▅▆▆▆▄▆▆▆▅▅▅▅▅▄▄▅▄▄ █
  28 ns         Histogram: log(frequency) by time      41.7 ns <

 Memory estimate: 32 bytes, allocs estimate: 2.

julia> @benchmark itp2(x, y) setup=(x=3+3*rand(); y=4+2*rand())
BenchmarkTools.Trial: 10000 samples with 979 evaluations.
 Range (min … max):  64.200 ns … 932.400 ns  ┊ GC (min … max): 0.00% … 91.49%
 Time  (median):     65.199 ns               ┊ GC (median):    0.00%
 Time  (mean ± σ):   67.493 ns ±  26.095 ns  ┊ GC (mean ± σ):  1.20% ±  2.91%

  ▂▃▇█▆▃▁                                  ▁▂▂▁▁▂▁             ▂
  ███████▇▇▅▆███▆▇▇█▇▇▅▅▆▃▄▄▁▃▄▁▁▁▄▃▁▁▁▁▆▇█████████▆▇▇▇▆▇▆▅▇▇▇ █
  64.2 ns       Histogram: log(frequency) by time      81.9 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

julia> @benchmark itp2l(x, y) setup=(x=3+3*rand(); y=4+2*rand())
BenchmarkTools.Trial: 10000 samples with 994 evaluations.
 Range (min … max):  33.037 ns … 896.178 ns  ┊ GC (min … max): 0.00% … 94.73%
 Time  (median):     34.229 ns               ┊ GC (median):    0.00%
 Time  (mean ± σ):   36.485 ns ±  25.571 ns  ┊ GC (mean ± σ):  2.18% ±  3.02%

  ▂▃▂█▇▅▃▆▄                                ▁▁▁▂▂▁              ▂
  █████████▇▆▆▆▇▇▅▆▅▆▆▆▆▄▆▄▃▁▄▃▁▁▁▁▄▃▁▃▅▇▇█████████▇▇▇▆▇▇▇▇▆▇▆ █
  33 ns         Histogram: log(frequency) by time      50.4 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

shell> git stash pop
On branch master
Your branch is up to date with 'origin/master'.

Changes not staged for commit:
  (use "git add <file>..." to update what will be committed)
  (use "git restore <file>..." to discard changes in working directory)
    modified:   src/Ratios.jl

no changes added to commit (use "git add" and/or "git commit -a")
Dropped refs/stash@{0} (da13fd30491a55a8f70425dbfb94122e543f35e7)

julia> SimpleRatio(2, 5) + SimpleRatio(1, 5)
SimpleRatio{Int64}(3, 5)

julia> @benchmark itp1(x) setup=(x=3+3*rand())
BenchmarkTools.Trial: 10000 samples with 995 evaluations.
 Range (min … max):  27.946 ns …  1.055 μs  ┊ GC (min … max): 0.00% … 96.76%
 Time  (median):     28.552 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   30.614 ns ± 26.329 ns  ┊ GC (mean ± σ):  2.25% ±  2.56%

   ▅█▆▂                            ▁▂▂▁▂▂▁                    ▁
  █████▆▅▅▃▅▅▄▄▄▄▄▄▃▄▅▆▆▆▅▃▄▄▄▂▇██▆███████▇▆▅▅▄▄▄▄▅▅▅▅▆▆▅▆▆▅▅ █
  27.9 ns      Histogram: log(frequency) by time      41.6 ns <

 Memory estimate: 32 bytes, allocs estimate: 2.

julia> @benchmark itp2(x, y) setup=(x=3+3*rand(); y=4+2*rand())
BenchmarkTools.Trial: 10000 samples with 979 evaluations.
 Range (min … max):  63.468 ns …  1.133 μs  ┊ GC (min … max): 0.00% … 92.65%
 Time  (median):     64.884 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   67.562 ns ± 32.278 ns  ┊ GC (mean ± σ):  1.48% ±  2.96%

  ▁▆█▇▃                  ▁▂▃▂▂                                ▂
  ██████▇▇▇██▇▇▇▇▆▆▅▅▄▄▄▇██████▇▇▇▇██▇▅▅▄▁▃▅▄▄▁▄▁▃▁▁▄▃▁▁▁▄▄▁▃ █
  63.5 ns      Histogram: log(frequency) by time      94.2 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

julia> @benchmark itp2l(x, y) setup=(x=3+3*rand(); y=4+2*rand())
BenchmarkTools.Trial: 10000 samples with 994 evaluations.
 Range (min … max):  33.057 ns …  1.052 μs  ┊ GC (min … max): 0.00% … 96.43%
 Time  (median):     34.306 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   36.839 ns ± 30.796 ns  ┊ GC (mean ± σ):  2.61% ±  3.04%

  ▂▃▂█▆▃▄▅▁                             ▁ ▁▂▁                 ▁
  █████████▆▅▆▆▆▅▅▆▇▅▅▆▆▅▅▅▃▄▅▄▃▄▃▃▃▅▆▇████████▆▅▆▄▅▆▆▆▆▇▆▆▄▅ █
  33.1 ns      Histogram: log(frequency) by time      51.9 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

shell> git diff
diff --git a/src/Ratios.jl b/src/Ratios.jl
index c7a4a98..75ffb4f 100644
--- a/src/Ratios.jl
+++ b/src/Ratios.jl
@@ -32,8 +32,10 @@ Rational{T}(r::SimpleRatio{S}) where {T<:Integer, S<:Integer} = convert(T, r.num
 /(x::Integer, y::SimpleRatio) = SimpleRatio(x*y.den, y.num)
 +(x::Integer, y::SimpleRatio) = SimpleRatio(x*y.den + y.num, y.den)
 -(x::Integer, y::SimpleRatio) = SimpleRatio(x*y.den - y.num, y.den)
-+(x::SimpleRatio, y::SimpleRatio) = SimpleRatio(x.num*y.den + x.den*y.num, x.den*y.den)
--(x::SimpleRatio, y::SimpleRatio) = SimpleRatio(x.num*y.den - x.den*y.num, x.den*y.den)
++(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num + y.num, x.den) :
+                                                     SimpleRatio(x.num*y.den + x.den*y.num, x.den*y.den)
+-(x::SimpleRatio, y::SimpleRatio) = x.den == y.den ? SimpleRatio(x.num - y.num, x.den) :
+                                                     SimpleRatio(x.num*y.den - x.den*y.num, x.den*y.den)
 ^(x::SimpleRatio, y::Integer) = SimpleRatio(x.num^y, x.den^y)

 -(x::SimpleRatio{T}) where {T<:Signed} = SimpleRatio(-x.num, x.den)

[timholy]

https://github.com/timholy/Ratios.jl/pull/24

[mkitti]

Just to summarize, with Ratios.jl 0.4.2 we now get the correct answer without the need for overflow detection:

julia> using Interpolations

julia> itpi = LinearInterpolation(Int[1,10000],Int[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Int64},), ::Vector{Int64}, Gridded(Linear())), Line()) with element type Float64:
 Ratios.SimpleRatio{Int64}(9999, 9999)
 Ratios.SimpleRatio{Int64}(99990000, 9999)

julia> convert(Float64, itpi(10001))
10001.0

julia> itpi(10001)
Ratios.SimpleRatio{Int64}(99999999, 9999)

With SaferIntegers 3.x, you can also still watch for overflow:

julia> using SaferIntegers, Interpolations

julia> itpi = LinearInterpolation(SafeInt[1,10000],SafeInt[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{SafeInt64},), ::Vector{SafeInt64}, Gridded(Linear())), Line()) with element type Float64:
 Ratios.SimpleRatio{SafeInt64}(9999, 9999)
 Ratios.SimpleRatio{SafeInt64}(99990000, 9999)

julia> itpi(10001)
Ratios.SimpleRatio{SafeInt64}(99999999, 9999)

The benchmark indicates the performance impact is either negligible or tolerable.

Also, if you prefer floating point, you can always still use floating point input:

julia> itpi = LinearInterpolation(Float64[1,10000], Float64[1,10000]; extrapolation_bc=Line())
2-element extrapolate(interpolate((::Vector{Float64},), ::Vector{Float64}, Gridded(Linear())), Line()) with element type Float64:
     1.0
 10000.0

julia> itpi(10001)
10001.0

I'm still interested in providing a mechanism to configure Interpolations.ratio so that the end user has some choice of the backend implementation. Alternatives to SimpleRatio from Ratios.jl include Float64 from Base, Rationals from Base, and FastRationals.jl.