/demo/manydigits.py

  1"""
  2This script calculates solutions to some of the problems from the
  3"Many Digits" competition:
  4http://www.cs.ru.nl/~milad/manydigits/problems.php
  5
  6Run with:
  7
  8    python manydigits.py
  9
 10"""
 11
 12from mpmath import *
 13from mpmath.libmp import to_fixed, bin_to_radix
 14
 15dps = 100
 16mp.dps = dps + 10
 17
 18def pr(x):
 19    """Return the first dps digits after the decimal point"""
 20    x = x._mpf_
 21    p = int(dps*3.33 + 10)
 22    t = to_fixed(x, p)
 23    d = bin_to_radix(t, p, 10, dps)
 24    s = str(d).zfill(dps)[-dps:]
 25    return s[:dps//2] + "\n" + s[dps//2:]
 26
 27print """
 28This script prints answers to a selection of the "Many Digits"
 29competition problems: http://www.cs.ru.nl/~milad/manydigits/problems.php
 30
 31The output for each problem is the first 100 digits after the
 32decimal point in the result.
 33"""
 34
 35print "C01: sin(tan(cos(1)))"
 36print pr(sin(tan(cos(1))))
 37print
 38
 39print "C02: sqrt(e/pi)"
 40print pr(sqrt(e/pi))
 41print
 42
 43print "C03: sin((e+1)^3)"
 44print pr(sin((e+1)**3))
 45print
 46
 47print "C04: exp(pi*sqrt(2011))"
 48mp.dps += 65
 49print pr(exp(pi*sqrt(2011)))
 50mp.dps -= 65
 51print
 52
 53print "C05: exp(exp(exp(1/2)))"
 54print pr(exp(exp(exp(0.5))))
 55print
 56
 57print "C06: arctanh(1-arctanh(1-arctanh(1-arctanh(1/pi))))"
 58print pr(atanh(1-atanh(1-atanh(1-atanh(1/pi)))))
 59print
 60
 61print "C07: pi^1000"
 62mp.dps += 505
 63print pr(pi**1000)
 64mp.dps -= 505
 65print
 66
 67print "C08: sin(6^(6^6))"
 68print pr(sin(6**(6**6)))
 69print
 70
 71print "C09: sin(10*arctan(tanh(pi*(2011^(1/2))/3)))"
 72mp.dps += 150
 73print pr(sin(10*atan(tanh(pi*sqrt(2011)/3))))
 74mp.dps -= 150
 75print
 76
 77print "C10: (7+2^(1/5)-5*(8^(1/5)))^(1/3) + 4^(1/5)-2^(1/5)"
 78a = mpf(1)/5
 79print pr(((7 + 2**a - 5*(8**a))**(mpf(1)/3) + 4**a - 2**a))
 80print
 81
 82print "C11: tan(2^(1/2))+arctanh(sin(1))"
 83print pr((tan(sqrt(2)) + atanh(sin(1))))
 84print
 85
 86print "C12: arcsin(1/e^2) + arcsinh(e^2)"
 87print pr(asin(1/exp(2)) + asinh(exp(2)))
 88print
 89
 90print "C17: S= -4*Zeta(2) - 2*Zeta(3) + 4*Zeta(2)*Zeta(3) + 2*Zeta(5)"
 91print pr(-4*zeta(2) - 2*zeta(3) + 4*zeta(2)*zeta(3) + 2*zeta(5))
 92print
 93
 94print "C18: Catalan G = Sum{i=0}{\infty}(-1)^i/(2i+1)^2"
 95print pr(catalan)
 96print
 97
 98print "C21: Equation exp(cos(x)) = x"
 99print pr(findroot(lambda x: exp(cos(x))-x, 1))
100print
101
102print "C22: J = integral(sin(sin(sin(x)))), x=0..1"
103print pr(quadts(lambda x: sin(sin(sin(x))), [0, 1]))
104print