notebook.community

Edit and run

It is well known that if the square root of a natural number is not an integer, then it is irrational. The decimal expansion of such square roots is infinite without any repeating pattern at all.

The square root of two is 1.41421356237309504880..., and the digital sum of the first one hundred decimal digits is 475.

For the first one hundred natural numbers, find the total of the digital sums of the first one hundred decimal digits for all the irrational square roots.



In [1]:

    
from math import sqrt



In [2]:

    
N_DIGITS = 100  # This is the number of digits to sum.
N_DIGITS_SQUARED = 2 * (N_DIGITS + 5)



In [3]:

    
a = int(2 * 10 ** N_DIGITS_SQUARED)
a









    Out[3]:





2000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000L



In [4]:

    
sqrt(a)









    Out[4]:





1.414213562373095e+105



In [5]:

    
long(sqrt(a))









    Out[5]:





1414213562373095104787002847759798204858218373714391208861616932096043069664092695467347475431626131898368L



In [6]:

    
type(2 ** 1000)









    Out[6]:





long



In [7]:

    
MAX_ITERATIONS = 10
EPSILON = 1
TOLERANCE = 1

def newton(f, fprime, n):
    xold = 1
    xold = long(sqrt(n))
    for i in range(MAX_ITERATIONS):
        denominator = fprime(xold, n)
        assert abs(denominator) >= EPSILON, 'too small'
        xnew = xold - long(long(f(xold, n))/long(denominator))
        # print i, xnew
        if abs(xnew-xold) < TOLERANCE:
            return xnew
        xold = xnew
    assert False, 'failed to converge %d, %d' % (i, xnew)



In [8]:

    
def f(x, n):
    return x*x - n



In [9]:

    
def f_prime(x, n):
    return 2*x



In [10]:

    
print a
b = newton(f, f_prime, a)
len(str(b)), b, sum((int(c) for c in str(b)[:100])), sum((int(c) for c in str(b)[1:101]))









    



2000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000






    Out[10]:





(106,
 1414213562373095048801688724209698078569671875376948073176679737990732478462107038850387534327641572735014L,
 475,
 481)



In [11]:

    
def foo(n):
    total = 0
    for i in range(n + 1):
        j = int(sqrt(i))
        if j*j == i:
            continue
        # print i
        a = i * 10 ** N_DIGITS_SQUARED
        b = newton(f, f_prime, a)
        total += sum((int(c) for c in str(b)[:100]))
    return total



In [12]:

    
n = 2
%timeit foo(n)
foo(n)









    



10000 loops, best of 3: 134 us per loop






    Out[12]:





475



In [13]:

    
n = 100
%timeit foo(n)
foo(n)









    



100 loops, best of 3: 11.8 ms per loop






    Out[13]:





40886