This challenge is to produce the shortest code for the constant \$\pi^{1/\pi}\$. Your code must output the first \$n\$ consecutive digits of \$\pi^{1/\pi}\$, where \$n\$ is given in the input. Alternatively, your program may accept no input, and output digits indefinitely

This is code golf, so the shortest submission (in bytes) wins except that it must output the 1000 digits for \$n = 1000\$ in less than 10 seconds on a reasonable PC.

You may not use a built-in for \$\pi\$, the gamma function or any trigonometic functions.

If the input is \$n = 1000\$ then the output should be:


Note the 1000 decimal places includes the first \$1\$.

Output format

Your code can output in any format you wish.

This related question tackles a simpler variant but with the same restrictions.


If it's helpful, you can assume \$\pi^{1/\pi}\$ is irrational and in fact you can even assume it is a Normal number.

  • 9
    \$\begingroup\$ Banning a built-in for pi or trigonometric functions is kind of useless. For example, we can use \$\Gamma(1/2)^2\$ \$\endgroup\$ – Luis Mendo May 20 at 9:05
  • 4
    \$\begingroup\$ I'm curious, is this even possible? Could it not be that digit number 1000000 of pi^(1/pi) relies on digit number 1000099? Or even 9999999, in theory at least..? At some point you'll have ...4599999999999999999999999999999999......, but that could be 4600000000000000000000000000000000..... if your calculations included more digits. I believe it's impossible to say how many more digits is needed..? Or? \$\endgroup\$ – Stewie Griffin May 20 at 10:01
  • 11
    \$\begingroup\$ Why the restriction on how the digits are output? It adds absolutely nothing to the challenge. \$\endgroup\$ – Shaggy May 20 at 11:12
  • 2
    \$\begingroup\$ For the record: math.stackexchange.com/q/3232906/92515 \$\endgroup\$ – Stewie Griffin May 20 at 12:46
  • 3
    \$\begingroup\$ Is ConvertDegreesToRadians(180) interesting enough as a way to produce pi? \$\endgroup\$ – someone May 20 at 15:26

Python, 149 bytes

Saved one byte due to @H.PWiz.

while i:p=2*z-i//2*p//~i;i-=2
while j:l=q//j+q*l//z;j-=1
while k:e=z+e*l//k//p;k-=1

Try it online!

Input n = 1000 finishes in less than 0.5s.

\$\sqrt[\pi]{\pi}\$ is calculated as the \$e^\frac{\ln(\pi)}{\pi}\$.\$\pi\$ is calculated with the usual Euler-Leibniz: \$\pi=\sum_{n=0}\limits^{\infty}{\frac{n!}{(2n+1)!!}}\$.

\$e^x\$ and \$\ln(x)\$ are both computed using a Taylor series: \$e^x=\sum\limits_{n=0}^{\infty}{\frac{x^n}{n!}}\$ and \$\ln(1-x)=-\sum\limits_{n=1}^{\infty}{\frac{x^n}{n}}\$. Because the iteration for \$\ln(x)\$ converges only when \$|x|<1\$, this is instead calculated as \$\ln(\pi)=-\ln(\frac{1}{\pi})\$.

Given that Python uses Karasuba Multiplication, the overall runtime complexity is \$\mathcal{O}(n^{1+\log_2(3)})\$ - in other words, twice as many digits will take approximately 6 times as long.

Subquadratic Complexity

import sys
from gmpy2 import isqrt, mpz

def piks(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) >> 1
  p1, q1, t1 = piks(a, m)
  p2, q2, t2 = piks(m+1, b)
  return (p1*p2, q1*q2, q2*t1 + p1*t2)

n = int(sys.argv[1])-1
m = n*20//3
z = mpz(10)**n

# n / log(777924, 10)
pi_terms = mpz(n*0.16975227728583067)

pp, pq, pt = piks(0, pi_terms)
pq *= 3528

pi2m = (pq << m) // pt

a, b = 2 << m, 8
while a != b:
  a, b = (a + b) >> 1, isqrt(a*b)

mlog2_pi = (z << m) // a

a, b = 2*pi2m, 8
while a != b:
  a, b = (a + b) >> 1, isqrt(a*b)
logpi_pi = z * pi2m // a - mlog2_pi

mlog2 = mlog2_pi * pq // pt

d = e = (15044673 << m) // 10450451
pt //= z
while d:
  a, b = 2*e, 8
  while a != b:
    a, b = (a + b) >> 1, isqrt(a*b)
  lnx = (pq * e) // (pt * a) - mlog2
  d = e * (lnx - logpi_pi) // z
  e -= d

print(e * z >> m)

Try it online!

Input n = 20000 finishes in less than one second.

\$\pi\$ can be computed in subquadratic time by use of Karatsuba splitting a.k.a. Fast E-function Evaluation, which reduces a summation to n terms to a single rational value p/q, splitting in binary descent. I've chosen to use Ramanujan #39, which is the fastest converging series of its kind that doesn't require an arbitrary precision square root, to my knowledge. Explicitly, this is computed as:

\$\large{\left.{3528}\middle/{\sum\limits_{n=0}^{\infty}\frac{(-1)^n (4n)! (1123+21460n)}{(n!)^4 14112^{2n}}}\right.}\$

For \$\ln(x)\$ and \$e^x\$, using the same technique wouldn't provide any benefit, because both are computed as a power series of an arbitrary precision variable. Fortunately, Gauss has gifted us with an elegant quadratically converging formula for \$\ln(x)\$, based on the Arithmetic-Geometric Mean:

\$\DeclareMathOperator{\AGM}{AGM}\large\ln(x)=\lim\limits_{n\rightarrow\infty}\frac{\pi x^n}{2n\AGM(x^n,4)}\$

or, in cases when large powers of x are inconvenient, this can also be computed as:

\$\large\ln(x)\approx\frac{\pi x 2^m}{2\AGM(x 2^m,4)}-m\log(2)\$

Conveniently, division by \$\pi\$ can be acheived simply by not multiplying through by \$\pi\$.

With the natural logarithm thusly defined, \$e^\frac{\ln(\pi)}{\pi}\$ can be computed via Newton's method on \$x_{n+1}=x_n-x_n(\ln(x_n)-\frac{\ln(\pi)}{\pi})\$. The overall complexity is then \$\mathcal{O}(n^*\log^3(n))\$, where \$n^*\$ will vary with the complexity of the multiplication algorithm GMP is using for any given bit length.

  • 1
    \$\begingroup\$ This is very nice! \$\endgroup\$ – Anush May 30 at 18:22
  • 3
    \$\begingroup\$ -1 \$\endgroup\$ – H.PWiz May 30 at 20:32
  • 1
    \$\begingroup\$ What a great update! \$\endgroup\$ – Anush Jun 12 at 13:49

Wolfram Language (Mathematica), 62 bytes

In my first try I used Zeta function but this one uses the imaginary part of ln(-1) for pi. (@someone)
Prints more than 4000 digits in the first 10 seconds,


Try it online!

9 bytes saved from @someone

  • \$\begingroup\$ I was hoping they might all be printed on the same line so the first 1000 digits look like the example I gave in the question. \$\endgroup\$ – Anush May 20 at 10:06
  • \$\begingroup\$ @Anush I can print digits for ever or give you any number of digits (but then it should stop). You must decide \$\endgroup\$ – J42161217 May 20 at 10:09
  • \$\begingroup\$ Your code should print digits forever. I just meant they should be printed horizontally and not vertically. \$\endgroup\$ – Anush May 20 at 10:11
  • 4
    \$\begingroup\$ This is the first time you are saying that. The question does not ask that. \$\endgroup\$ – J42161217 May 20 at 10:13
  • 13
    \$\begingroup\$ it is highly recommended that you don't change the rules of a challenge after even one answer is submitted. I would suggest you to use sandbox next time. \$\endgroup\$ – J42161217 May 20 at 10:19

Wolfram Language (Mathematica), 23 bytes


Try it online!

Taking a lot of inspiration from @J42161217 and displaying n digits instead of an infinite stream.

Note that I'm using ToString in Tio to suppress the trailing precision specifier in the output, which does not appear when this code is executed in Mathematica.

Maybe using 180° for pi would work; but maybe that's too close to being trigonometric and is disallowed.


AXIOM, 221 bytes

p(n:PI):String==(d:=digits(n+9);e:=10.^-digits();i:=s:=0;repeat(k:=m(i);k<e=>break;s:=s+k;i:=i+1);r:=concat split((s^(1/s))::String,char " ");digits(d);r.(1..(n+1)))

test and ungolf:

(3) -> p 1000
                                                             Type: String
           Time: 0.03 (IN) + 0.53 (EV) + 0.25 (OT) + 0.15 (GC) = 0.97 sec

-- Bailey-Borwein-Plouffe formula for pi
-- https://en.m.wikipedia.org/wiki/Bailey%E2%80%93Borwein%E2%80%93Plouffe_formula   
--       +oo
--      -----          --                                      --
--      \          1   |     4         2         1         1    |
--   pi= |      -------|  ------- - ------- - ------- - ------- |
--      /            k |   8*k+1     8*k+4     8*k+5     8*k+6  |
--      ----- k   16   --                                      --
--        0
      d:=digits(n+9); e:=10.^-digits(); i:=s:=0
      r:=concat split((s^(1/s))::String,char " ")

the function m calculate the term of the sum; the p() function loop, sum them in the variable s until the term is < than min float value (one can see as epsilon); p() function return one string of that number pi^(1/pi); it seems the required digits are returned in less than 1 second for input to p() function 1000.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.