Skip to content

Instantly share code, notes, and snippets.

@tompng
Last active September 15, 2025 12:33
Show Gist options
  • Select an option

  • Save tompng/4dee06ddb572bc308e25bcb62fb5439e to your computer and use it in GitHub Desktop.

Select an option

Save tompng/4dee06ddb572bc308e25bcb62fb5439e to your computer and use it in GitHub Desktop.
def exp_bs(x, prec)
# Taylor series of exp can be converted to continued fraction as:
# exp(x)-1 = x*(1+x/2*(1+x/3*(1+x/4*(1+x/5*(1+...)))))
# Binary splitting can be applied to this continued fraction as:
# (1 + a0/a1 + x/a1 * (1 + b0/b1 + x/b1 * rest))
# (1 + ((a0+x)*b1+x*b0)/(a1*b1) + (x*x)/(a1*b1) * rest)
x = BigDecimal(x)
step = (1..).bsearch { |k| Math.lgamma(k)[0] / Math.log(10) - k * x.exponent > prec }
ms = (1..step).map { [0, BigDecimal(it)] }
while ms.size > 1
xn = xn ? xn.mult(xn, prec) : x
ms = ms.each_slice(2).map do |a, b|
b ||= [0, BigDecimal(step)]
[
(a[0] * b[1]).add(xn * (b[0] + b[1]), prec),
a[1].mult(b[1], prec)
]
end
end
numerator, denominator = ms.first
1 + numerator.div(denominator, prec)
end
def exp(x, prec)
prec2 = prec + BigDecimal.double_fig + [x.exponent, 0].max
x = x.mult(1, prec2)
n = 2
ans = 1
while x != 0 do
xx = x.truncate(n)
x -= xx
t=Time.now
ans = exp_bs(xx, prec2).mult(ans, prec2)
n *= 2
p n => Time.now-t
end
ans.mult(1, prec)
end
# Requires NTT multiplication
exp(BigDecimal('0.'+'7'*100000), 100000)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment