# Computing binomial(n,k) without overflow

**URL:** <https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349>\
**Category:** General Usage\
**Created:** [January 27, 2024, 10:16pm UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349 "2024-01-27T22:16:46Z")\
**Posts on this page:** 15\
**Page:** 1

<div class="post-metadata">

**Author:** ![mike.ingold](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mike.ingold/32/203749_2.png) [@mike.ingold](https://discourse.julialang.org/u/mike.ingold)\
**Post date:** [January 27, 2024, 10:16pm UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/1 "2024-01-27T22:16:46Z")

</div>

The function `binomial(n,k)` throws an `OverflowError` if the standard `::Int64` arguments are large enough that the result overflows an `::Int64`. This can be modified to `binomial(big(b),k)` to promote the output to `BigInt`, e.g.:

```julia
julia> binomial(73, 24)
ERROR: OverflowError: binomial(73, 24) overflows
Stacktrace:
 [1] binomial(n::Int64, k::Int64)
   @ Base .\intfuncs.jl:1114
 [2] top-level scope
   @ REPL[32]:1

julia> binomial(big(73), 24)
11844267374132633700

```

I’m trying to write a function as generically as possible that internally computes a binomial, and occasionally the inputs are large enough to require a `BigInt` type. In those cases the output of the overall function winds up turning into `BigFloat`.

Does anybody have an elegant generic way to compute `binomial` where the types are promoted to `BigInt` i.f.f. a `BigInt` output would be required? A `try/catch` could work but seems probably suboptimal.

---

<div class="post-metadata">

**Author:** ![mike.ingold](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mike.ingold/32/203749_2.png) [@mike.ingold](https://discourse.julialang.org/u/mike.ingold)\
**Post date:** [January 27, 2024, 10:39pm UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/2 "2024-01-27T22:39:41Z")

</div>

Less abstractly, this is where the issue crops up for me.

> <https://github.com/mikeingold/LineIntegrals.jl/blob/ebb30955a9bea68360a83edc47afaa3aa6da1ba4/src/utils.jl#L14>

In the short term I just forced everything up to `BigInt`, but the function outputs don’t require `BigFloat` level of precision even in those cases.

---

<div class="post-metadata">

**Author:** ![jar1](https://avatars.discourse-cdn.com/v4/letter/j/c0e974/32.png) [@jar1](https://discourse.julialang.org/u/jar1)\
**Post date:** [January 27, 2024, 10:42pm UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/3 "2024-01-27T22:42:56Z")

</div>

If your function’s return type depends on its argument values rather than just its inputs, then it isn’t “type stable”, which has some performance downsides and is more complex to use.

In the case where the result is small enough to fit in Int64, it’s faster to do try/catch than to always use `big`. In the overflow case, it’s faster to use `big` in the first place.

A third option is to convert to BigInt after the trycatch, so the function’s return value is statically known but it has internal dynamic types.

```julia
binomial_trycatch(n,k) =
    try
        binomial(n,k)
    catch e
        e isa OverflowError && binomial(big(n), k)
    end

binomial_trycatch_then_big(n,k)::BigInt =
    try
        binomial(n,k)
    catch e
        e isa OverflowError && binomial(big(n), k)
    end

binomial_big(n,k) = binomial(big(n),k)

using BenchmarkTools
# Fits in Int64.
@btime binomial_trycatch(73, 23);
@btime binomial_trycatch_then_big(73, 23);
@btime binomial_big(73, 23);
# Doesn't fit in Int64.
@btime binomial_trycatch(73, 24);
@btime binomial_trycatch_then_big(73, 24);
@btime binomial_big(73, 24);

```

```julia
julia> @btime binomial_trycatch(73, 23);
  177.726 ns (0 allocations: 0 bytes)

julia> @btime binomial_trycatch_then_big(73, 23);
  205.040 ns (2 allocations: 40 bytes)

julia> @btime binomial_big(73, 23);
  399.848 ns (12 allocations: 160 bytes)

julia> @btime binomial_trycatch(73, 24);
  37.892 μs (14 allocations: 248 bytes)

julia> @btime binomial_trycatch_then_big(73, 24);
  34.696 μs (14 allocations: 248 bytes)

julia> @btime binomial_big(73, 24);
  399.811 ns (11 allocations: 152 bytes)

```

---

<div class="post-metadata">

**Author:** ![Benny](https://avatars.discourse-cdn.com/v4/letter/b/49beb7/32.png) [@Benny](https://discourse.julialang.org/u/Benny)\
**Post date:** [January 27, 2024, 11:33pm UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/4 "2024-01-27T23:33:10Z")

</div>

Not sure from that exactly where the overhead of overflowing is, whether it’s recomputing `binomial` or catching an error, but I think the only way to cut down on that any further is to define a separate `binomial` function that `break`s the while loop where [`Base.binomial`](https://github.com/JuliaLang/julia/blob/3120989f39bb7ef7863c4aab8ab1227cf71eec66/base/intfuncs.jl#L1092) would throw the error. Then, it resumes the algorithm with `BigInt`, possibly with a rewrite that picks up where the loop left off instead of starting over.

---

<div class="post-metadata">

**Author:** ![mike.ingold](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mike.ingold/32/203749_2.png) [@mike.ingold](https://discourse.julialang.org/u/mike.ingold)\
**Post date:** [January 27, 2024, 11:59pm UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/5 "2024-01-27T23:59:13Z")

</div>

If I set some threshold on `n` where this is likely to occur, is there any way to dispatch to a helper function based on that value?

In this situation the ‘n’ is essentially inherent to the `BezierCurve` type, but I don’t think it’s accessible as a type parameter.

---

<div class="post-metadata">

**Author:** ![jar1](https://avatars.discourse-cdn.com/v4/letter/j/c0e974/32.png) [@jar1](https://discourse.julialang.org/u/jar1)\
**Post date:** [January 28, 2024, 12:44am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/6 "2024-01-28T00:44:13Z")

</div>

Another option is to use `Int128`

```julia
julia> binomial_i128(n,k) = binomial(Int128(n),k)
binomial_i128 (generic function with 1 method)

julia> @btime binomial_i128(73, 23);
  50.938 ns (0 allocations: 0 bytes)

julia> @btime binomial_i128(73, 24);
  49.638 ns (0 allocations: 0 bytes)

```

---

<div class="post-metadata">

**Author:** ![mike.ingold](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mike.ingold/32/203749_2.png) [@mike.ingold](https://discourse.julialang.org/u/mike.ingold)\
**Post date:** [January 28, 2024, 1:36am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/7 "2024-01-28T01:36:19Z")

</div>

Ooh, good point. Looks like converting `n` to `Int128` extends the non-erroring input range into sufficiently large numbers with good performance and prevents the need for the wrapping function to return a `BigFloat`.

@jar1 I marked your earlier post as the solution since it actually answered the original question, but this Int128 tip was the solution I actually needed. Thanks!

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [January 28, 2024, 1:51am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/8 "2024-01-28T01:51:59Z")

</div>

The accuracy of Float64 should be enough for most purposes. So perhaps the following will be enough:

```julia
using SpecialFunctions

logbinomial(n,k) = logfactorial(n)-logfactorial(k)-logfactorial(n-k)
B(i,n) = t -> exp(logbinomial(n,i) + i*log(t) + (n-i)*log1p(-t))

```

As an example:

```julia
julia> using UnicodePlots

julia> lineplot(0,1,B(3,6))
            ┌────────────────────────────────────────┐        
        0.4 │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│ #198(x)
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⣀⠔⠒⠒⠢⣄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⢀⠎⠁⠀⠀⠀⠀⠀⠑⡄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⡰⠃⠀⠀⠀⠀⠀⠀⠀⠀⠘⢆⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⡸⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠈⢣⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
   f(x) │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⡜⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠈⢣⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⡜⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⢣⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⡜⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⢣⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⠀⠀⡜⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⢣⠀⠀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⠀⢀⡜⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⢣⡀⠀⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⠀⢀⠎⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠱⡄⠀⠀⠀⠀⠀⠀│        
            │⠀⠀⠀⠀⠀⡠⠃⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠈⢆⠀⠀⠀⠀⠀│        
          0 │⣀⣀⣀⠤⠊⠁⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⠦⣄⣀⣀│        
            └────────────────────────────────────────┘        
            ⠀0⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀1⠀        

```

If performance becomes an issue, there might be faster ways to get at the results without the anonymous functions.

---

<div class="post-metadata">

**Author:** ![PeterSimon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petersimon/32/25193_2.png) [@PeterSimon](https://discourse.julialang.org/u/PeterSimon)\
**Post date:** [January 28, 2024, 3:07am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/9 "2024-01-28T03:07:29Z")

</div>

The `Int128` solution of @Jar1 always returns a `Float64` for integer inputs. As pointed out by @Dan this is perfectly fine for @mike.ingold 's specific application because he’s using these coefficients in a floating point calculation. So another, faster way to get a floating point approximation of the binomial is

```julia
binomial_float(n, k) = prod((n+1-i)/i for i in 1:min(k, n-k); init=1.0)

julia> @btime binomial_i128(73, 23)
  51.368 ns (0 allocations: 0 bytes)
5.685248339583665e18

julia> @btime binomial_float(73, 23)
  19.138 ns (0 allocations: 0 bytes)
5.685248339583665e18

```

Edit: After looking over the source of `binomial`, when called with a non-integer first argument, I see it’s almost the same as above. So a better definition for `binomial_float` would be

```julia
binomial_float(n, k) = binomial(float(n), k)

```

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [January 28, 2024, 3:39am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/10 "2024-01-28T03:39:10Z")

</div>

> [@PeterSimon](#):
>
> `19.138 ns (0 allocations: 0 bytes)`

I think the timings here are due to compile-time optimization.

---

<div class="post-metadata">

**Author:** ![PeterSimon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petersimon/32/25193_2.png) [@PeterSimon](https://discourse.julialang.org/u/PeterSimon)\
**Post date:** [January 28, 2024, 3:41am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/11 "2024-01-28T03:41:13Z")

</div>

```julia
julia> n = 73; k = 23; @btime binomial_float($n, $k)
  20.060 ns (0 allocations: 0 bytes)
5.685248339583665e18

```

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [January 28, 2024, 3:44am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/12 "2024-01-28T03:44:15Z")

</div>

Okay… still surprises me how fast computers have become over the years.

---

<div class="post-metadata">

**Author:** ![PeterSimon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petersimon/32/25193_2.png) [@PeterSimon](https://discourse.julialang.org/u/PeterSimon)\
**Post date:** [January 28, 2024, 3:44am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/13 "2024-01-28T03:44:38Z")

</div>

Maybe SIMD?

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [January 28, 2024, 3:45am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/14 "2024-01-28T03:45:21Z")

</div>

My version also in the same ballpark:

```julia
julia> n = 73; k = 23; @btime exp(logbinomial($n, $k))
  37.290 ns (0 allocations: 0 bytes)
5.685248339583862e18

```

---

<div class="post-metadata">

**Author:** ![PeterSimon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petersimon/32/25193_2.png) [@PeterSimon](https://discourse.julialang.org/u/PeterSimon)\
**Post date:** [January 28, 2024, 3:52am UTC](https://discourse.julialang.org/t/computing-binomial-n-k-without-overflow/109349/15 "2024-01-28T03:52:11Z")

</div>

Edited: After re-starting Julia I’m getting results consistent with you:

```julia
julia> n=73; k=23; @btime exp(logbinomial($n, $k))
  57.259 ns (0 allocations: 0 bytes)
5.685248339583862e18

```
