5 Seconds to Find Pie
Python - 65535

http://ideone.com/knKRhn



from math import exp, log

def divnr(p, q):
  """
    Integer division p/q using Newton-Raphson Division.
    Assumes p &gt; q &gt; 0.
  """

  sp = p.bit_length()-1
  sq = q.bit_length()-1
  sr = sp - sq + 1

  s = []
  t = sr
  while t &gt; 15:
    s = [t] + s
    t = (t&gt;&gt;1) + 1
  # Base-case division
  r = (1 &lt;&lt; (t&lt;&lt;1)) / (q &gt;&gt; sq-t)

  for u in s:
    r = (r &lt;&lt; u-t+1) - (r*r * (q &gt;&gt; sq-u) &gt;&gt; (t&lt;&lt;1))
    t = u
  return (r * (p &gt;&gt; sq)) &gt;&gt; sr

def pibs(a, b):
  if a == b:
    if a == 0:
      return (1, 1, 1123)
    p = a*(a*(32*a-48)+22)-3
    q = a*a*a*24893568
    t = 21460*a+1123
    return (p, -q, p*t)
  m = (a+b) &gt;&gt; 1
  p1, q1, t1 = pibs(a, m)
  p2, q2, t2 = pibs(m+1, b)
  return (p1*p2, q1*q2, q2*t1 + p1*t2)

def ebs(a, b):
  if a == b:
    if a == 0:
      return (1, 1)
    return (1, a)
  m = (a+b) &gt;&gt; 1
  p1, q1 = ebs(a, m)
  p2, q2 = ebs(m+1, b)
  return (p1*q2+p2, q1*q2)

if __name__ == '__main__':
  n = input()

  pi_terms = int(n*0.16975227728583067)

  # 10^n == e^p
  p = n*2.3025850929940457

  # Lambert W_0(p/e) a la Newton
  k = log(p) - 1
  w = k - (k-1)/(k+1)
  while k &gt; w:
    k = w
    w -= (k - p*exp(-k-1))/(k+1)

  # InverseGamma(e^p) approximation
  e_terms = int(p / w)

  pp, pq, pt = pibs(0, pi_terms)
  ep, eq = ebs(0, e_terms)

  z = 10**n
  p = 3528*z*ep*abs(pq)
  q = eq*abs(pt)

  pie = divnr(p, q)
  print pie,




Ideone doesn't seem to have gmpy2 installed, which is unfortunate for at least two reasons. One, because it would make the calcuation a lot faster, and two, because it makes any formula requiring an arbitrary precision square root impractical.

The formula I use for π was listed by Ramanujan as Formula (39):



which converges at the rate of ~5.89 digits per term. To my knowledge, this is the fastest converging series of its kind that doesn't require the evaluation of an arbitrary precision square root. Formula (44) in the same paper (convergence rate ~7.98 digits per term) is most often referred to as the Ramanujan Formula.

The formula I use for e is the sum of inverse factorials. The number of terms required is calculated as Γ-1(10n), using an approximation I found on mathoverflow. The Lambert W0 component is found using Newton's Method.

The calculation of each of these summations is done via the Fast E-function Evalution (more generally known as binary-splitting), originally devised by Karatsuba. The method reduces a summation to n terms to a single rational value p/q. These two values are then multiplied to produce the final result.

Update:
Profiling revealed that over half of the time needed for the calculation was spent in the final division. Only the upper most log2(10n) bits of q are needed to obtain full precision, so I trim a few off beforehand. The code now fills the Ideone output buffer in 3.33s.

Update 2:
Since this is an optimization challenge, I decided to write my own division routine to combat CPython's slowness. The implementation of divnr above uses Newton-Raphson Division. The general idea is to calculate d = 1/q · 2n using Newton's Method, where n is the number of bits the result requires, and compute the result as p · d >> n. Runtime is now 2.87s - and this is without even chopping off bits before the calculation; it's unnecessary for this method.