bancaldo™

The smallest number expressible as the sum of a prime square, prime cube, and prime fourth power is 28. In fact, there are exactly four numbers below fifty that can be expressed in such a way:

28 = 2^2 + 2^3 + 2^4
33 = 3^2 + 2^3 + 2^4
49 = 5^2 + 2^3 + 2^4
47 = 2^2 + 3^3 + 2^4

How many numbers below fifty million can be expressed as the sum of a prime square, prime cube, and prime fourth power?

import time


def is_prime(num):
    if num == 2:
        return True
    if num % 2 == 0 or num <= 1:
        return False
    sqr = int(num ** 0.5) + 1
    for divisor in xrange(3, sqr, 2):
        if num % divisor == 0:
            return False
    return True

def euler87(limit):
    """a**2 + b**3 + c**4"""
    results = set()
    primes = [x for x in xrange(2, int(limit ** 0.5))
              if is_prime(x)]
    for p4 in primes:
        for p3 in primes:
            partial = p4**4 + p3 ** 3
            if partial >= limit:
                break
            for p2 in primes:
                tot = p2 ** 2 + partial
                if tot >= limit:
                    break
                results.add(tot)
    return len(results)


if __name__ == "__main__":
    start_time = time.time()
    print "euler 87: {} \nelapsed time: {}s".format(
        euler87(50000000), time.time() - start_time)

Questi i 3 requisiti per giungere alla soluzione corretta:

1 – Per “100 decimal digits” si intendono tutte le cifre comprese quelle a sinistra della virgola”
Nel caso della radice di 2 (1, 414…) si deve contare anche “1”.414.
2 – Per evitare arrotondamenti sulla ultima cifra decimale (rounding), utilizzare una precisione maggiore di 100
(102 è sufficiente)
3 – Non vanno sommati i quadrati perfetti (es. quando calcoliamo la radice di 4, non
dobbiamo sommare 2 al totale). SOLO gli IRRAZIONALI

from decimal import Decimal, getcontext
import time

getcontext().prec = 102
start = time.time()

tot = 0
numbers = [n for n in range(2, 100) if (n**0.5) % 1 != 0]

for num in numbers:
    tot += sum(Decimal(num).sqrt().as_tuple()[1][:100])
print "euler 80: %s\nelapsed time: %ss." % (tot, time.time() - start)

Qui il testo del problema 89

Si deve calcolare il totale dei caratteri risparmiati durante la semplificazione dei numeri romani.
Ci ho messo un po’ a capirlo e continuavo a calcolare il totale dei caratteri della lista corretta….

import time


def filter_romans(path):
    with open(path) as inputfile:
        romans = [num.strip() for num in inputfile]
    startchr = sum([len(r) for r in romans])

    for bad, correct in [("VIIII", "IX"), ("IIII", "IV"),
                         ("LXXXX", "XC"), ("XXXX", "XL"),
                         ("DCCCC", "CM"), ("CCCC", "CD"),]:
        for roman in romans:
            index = romans.index(roman)
            if bad in roman:
                romans[index] = roman.replace(bad, correct)
    return startchr - sum([len(r) for r in romans])


if __name__ == "__main__":
    start = time.time()
    result = filter_romans(r"p089_roman.txt")
    print "euler 89: %s\nElapsed time %ss." % (result, time.time() - start)

It is possible to write ten as the sum of primes in exactly
five different ways:

7 + 3
5 + 5
5 + 3 + 2
3 + 3 + 2 + 2
2 + 2 + 2 + 2 + 2

What is the first value which can be written as the sum of
primes in over five thousand different ways?

Analisi:

stesso metodo risolutivo visto nei problemi 31 e 76.

python:

import time

START = time.time()

def is_prime(num):
    if num <= 1: return False
    elif num == 2: return True
    elif num % 2 == 0: return False
    else:
        d = 3
        r = int(num**0.5)
        while d <= r:
            if num % d == 0: return False
            d += 2
        return True

PRIMES = [p for p in xrange(2, 500) if is_prime(p)]

def eu_77(limit):
    target = 11
    while True:
        qty = [0] * (target + 1)
        qty[0] = 1
        for p in PRIMES:
            for j in range(p, target + 1):
                qty[j] += qty[j - p]
        if qty[target] > limit:
            return target
        target += 1
 
print "Euler 77: {}".format(eu_77(5000))
print time.time() - START

It is possible to write five as a sum in exactly six different ways:
4 + 1
3 + 2
3 + 1 + 1
2 + 2 + 1
2 + 1 + 1 + 1
1 + 1 + 1 + 1 + 1
How many different ways can one hundred be written as a sum of at
least two positive integers?

Analisi:
Stesso metodo del problema 31

python:

import time

START = time.time()

def eu_76(limit):
    qty = [0] * (limit + 1)
    qty[0] = 1
    for number in range(1,limit + 1):
        for i in range(number, limit + 1):
            qty[i] += qty[i - number]
    return qty[limit]

print "euler 76: {}".format(eu_76(100))
print time.time() - START

variante:

“Euler invented a generating function which gives rise to a recurrence equation in P(n)..”
La formula interessante é:

ne consegue che:

partitions = {}
def partition(n):
    if n < 0: return 0
    if n == 0: return 1
    if n not in partitions:
        # http://mathworld.wolfram.com/PartitionFunctionP.html
        partitions[n] = sum([(-1)**(k+1) * (partition(n - (k * (3 * k - 1) / 2))
            + partition(n - (k * (3 * k + 1) / 2))) for k in range(1, n+1)])
    return partitions[n]


print partition(100) - 1

It turns out that 12 cm is the smallest length of wire that can be bent
to form an integer sided right angle triangle in exactly one way,
but there are many more examples.

12 cm: (3,4,5)
24 cm: (6,8,10)
30 cm: (5,12,13)
36 cm: (9,12,15)
40 cm: (8,15,17)
48 cm: (12,16,20)

In contrast, some lengths of wire, like 20 cm, cannot be bent to form
an integer sided right angle triangle, and other lengths allow more than
one solution to be found; for example, using 120 cm it is possible to form
exactly three different integer sided right angle triangles.

120 cm: (30,40,50), (20,48,52), (24,45,51)

Given that L is the length of the wire, for how many values of
L <= 1,500,000 can exactly one integer sided right angle triangle be formed?

Analisi:

Da wikipedia (Pythagorean triple), sappiamo che, dati due interi
‘n’ ed ‘m’, con n > m, possiamo calcolare la lunghezza dei tre lati del
triangolo. E’ importante che ‘n’ ed ‘m’ siano coprimi tra di loro, cioè che
abbiano massimo comune divisore (gcd) = 1 e che m – n sia dispari.
Questo è necessario per avere delle triplette pitagoriche primitive (3,4,5=12).
In caso contrario otterremmo delle triplette non primitive (multiple):
120 = (30,40,50) , (20,48,52), (24,45,51)
come si nota sono tutte triplette i cui ‘lati’ sono multipli dei lati
di triplette primitive;
30, 40, 50 (3, 4, 5)*10
20,48, 52 (5, 12, 13)*4
24, 45, 51 (8, 15, 17)*3
Queste sono le triplette che NON vogliamo contare durante il nostro ciclo.

a = m^2 - n^2
b = 2*m*n
c = m^2 + n^2

python:

import time
from fractions import gcd

st = time.time()


def eu_75(limit):
    lengths = [0] * limit
    for m in range(1, int(limit ** 0.5), 2):
      for n in range(2, int(limit ** 0.5) - m, 2):
        if gcd(m, n) == 1:
          per = abs(n * n - m * m) + (2 * m * n) + (m * m + n * n) # (a) + (b) + (c)
          for p in range(per, limit, per):
            lengths[p] += 1
    return lengths

print "Answer to PE75 = ", eu_75(1500000).count(1)
print time.time() - st

The number 145 is well known for the property that the sum of the factorial
 of its digits is equal to 145:

1! + 4! + 5! = 1 + 24 + 120 = 145

Perhaps less well known is 169, in that it produces the longest chain of
numbers that link back to 169; it turns out that there are only three such
loops that exist:

169 -> 363601 -> 1454 -> 169
871 -> 45361 -> 871
872 -> 45362 -> 872

It is not difficult to prove that EVERY starting number will eventually
get stuck in a loop. For example,

69 -> 363600 -> 1454 -> 169 -> 363601 (-> 1454)
78 -> 45360 -> 871 -> 45361 (-> 871)
540 -> 145 (-> 145)

Starting with 69 produces a chain of five non-repeating terms,
but the longest non-repeating chain with a starting number below one million
is sixty terms. How many chains, with a starting number below one million,
contain exactly sixty non-repeating terms?

Analisi:

Non molto da dire: ho utilizzato math.factorial per calcolare il fattoriale
di ogni numero. Per calcolare la somma dei fattoriali di ogni cifra di un
numero ho utilizzato un minimo di memoizzazione, per non sforare dal solito
minuto di tempo.
Poi per ogni numero, utilizzo una lista alla quale aggiungo la nuova somma
fino a chè questa somma non è già presente in lista. In caso contrario, con
l’apposita funzione, restituisco la lunghezza di tale lista (catena).
Per il resto è un Brute Force dove conteggio i numeri che hanno una catena di
sessanta elementi.

python:

import time
from math import factorial as fact


START = time.time()

_MEMO = {}
COUNT = 0


def sum_fact(num):
    '''Return the factorization sum'''
    if num not in _MEMO:
        res = 0
        for digit in str(num):
            res += fact(int(digit))
        _MEMO[num] = res
    return _MEMO[num]


def get_chains(num):
    '''Return chain length'''
    chains = []
    while num not in chains:
        chains.append(num)
        num = sum_fact(num)
    return len(chains)

for i in xrange(1, 10 ** 6):
    ch = get_chains(i)
    if ch == 60:
        COUNT += 1

print "Project Euler 74: {}".format(COUNT)
print "elapsed time: {}sec".format(time.time() - START)

Consider the fraction, n/d, where n and d are positive integers.
If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
If we list the set of reduced proper fractions for d <= 8 in
ascending order of size, we get:

1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8,
2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8

It can be seen that there are 3 fractions between 1/3 and 1/2.

How many fractions lie between 1/3 and 1/2 in the sorted set of reduced proper fractions for d <= 12,000?

Analisi:

Con un approccio Brute Force si sfora sicuramente dal minuto limite.
Utilizzando la sequenza di Farey, si ottimizza e non poco, il tempo
di esecuzione.
Date due frazioni adiacenti a/b e c/d, la formula che ci interessa per calcolare la successiva frazione e/f è:

e = int((n + b) / d) * c - a
f = int((n + b) / d) * d - b

Avendo come punto di partenza 1/3, considerando n = 12000, dobbiamo
trovare la frazione adiacente ad 1/3. Ci viene in soccorso dal link precedente, la seconda formula che dice che:
due frazioni a/b e c/d sono adiacenti se:

b * c – a * d = 1

in questo caso 1/3 e 4000/11999

Fatto questo, con un normale while-loop continueremo a calcolare le successive frazioni adiacenti fino a che il risultato non sarà 1/2.

python:

import time

START = time.time()

def farey(n):
    a, b, c, d = 1, 3, 4000, 11999
    res = 0
 
    while c != 1 and d != 2:
        res += 1
        k = int((n + b) / d)
        e = k * c - a;
        f = k * d - b;
        a, b = c, d
        c, d = e, f
    return res
            
print farey(12000)
print time.time() - START

Consider the fraction, n/d, where n and d are positive integers.
If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
If we list the set of reduced proper fractions for d <= 8 in
ascending order of size, we get:
1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5,
5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
It can be seen that there are 21 elements in this set.
How many elements would be contained in the set of reduced proper
fractions for d <= 1,000,000?

Analisi:

Come facciamo a determinare il numero di frazioni proprie ridotte di un numero?
Una frazione propria ridotta è quella frazione che ha il numeratore minore del denominatore
(propria) e con gcd (massimo comune denominatore) = 1.
Significa anche che numeratore e denominatore sono relativamente primi.
Quindi, dato un numero ‘i’, tutti i numeratori relativamente primi ad ‘i’ e minori di ‘i’, con
‘i’ al denominatore, sono frazioni proprie ridotte.
Per calcolare il numero di interi compresi tra 1 ed ‘i’ che sono coprimi con ‘i’, si usa
la funzione di Eulero.
La formula che ci interessa è:

phi(i) = i * [(1 - 1/p1) * (1 - 1/p2) ...]

Un approccio Brute Force purtroppo sfora dal minuto richiesto, quindi
utilizzeremo un algoritmo simile al crivello di eratostene.
Iteriamo cioè sun una lista contenente tutti i numeri da 2 al limite
richiesto; quando il valore contenuto in quell’indice di lista è uguale
al valore dell’indice stesso, allora iteriamo nuovamente su tutti i valori
compresi tra quel numero ed il limite (con passo pari al numero stesso, in modo
che phi venga calcolata solo sui divisori del numero stesso).

questo è il codice:

import time

st = time.time()

def euler_72(limit):
    phi = range(limit)
    for i in xrange(2, limit):
        if phi[i] == i: # sieve
            for j in xrange(i, limit, i):
                phi[j] *= 1 - 1.0/i
    return int(sum(phi) - 1)


print "Project Euler 72: {}".format(euler_72(10**6 + 1))
print "elapsed time: {}sec".format(time.time() - st)

Spiegazione:

partiamo con un esempio facile.
Cerchiamo il numero di frazioni proprie ridotte per n <= 4. Sappiamo che: phi(2) = 1 phi(3) = 2 (i primi hanno risultato p-1) phi(4) = 2 quindi la funzione deve dare come risultato 5 [sourcecode lang="text"] >>> euler_72(5) 5[/sourcecode] Come funziona? A. viene creata la lista phi [sourcecode lang="text"] >>> phi = range(5) >>> phi [0, 1, 2, 3, 4] >>> [/sourcecode] B. si itera sui numeri da 2 a 4 compreso quando phi[i] = i entra in azione il crivello. il primo numero è ovviamente 2, siccome phi[2] = 2, iniziamo con il secondo for-loop che itererà su numeri da 2 a 4 compreso, con step 2, quindi 2 e 4. al primo giro abbiamo che [sourcecode lang="text"]phi[j] = phi[j](1 - 1.0/i)[/sourcecode] siccome j=2, abbiamo che phi[2] = 2 quindi phi[2] = 2 * (1 - 1/2) = 1 in questo modo phi[2] = 1 (non più 2) ora al secondo giro avremo j=4 [sourcecode lang="text"]phi[j] = phi[j](1 - 1.0/i)[/sourcecode] phi[4] = 4 quindi phi[4] = 4 * (1 - 1/2) = 2 ora anche phi[4] ha cambiato il proprio valore. Finito il loop annidato, torniamo a quello esterno e da i = 2, passiamo ad i = 3 phi[3] non è stato ancora toccato infatti phi[3]=3 entriamo nel secondo loop dove itereremo da 3 a 5 con step = 3, quindi solo un "giro" (3). phi[3] = 3 quindi phi[3] = 3 * (1 - 1/3) = 2 usciamo dal loop annidato e torniamo in quello esterno: i = 4 phi[4] != 4 perchè modificato dal loop annidato durante il primo ciclo, quindi non succede nulla. La lista modificata sarà: [0, 1, 1.0, 2.0, 2.0] Ora basterà restituire la somma di tali valori, escluso il secondo (phi(1) = 1) che non ci interessa.

Consider the fraction, n/d, where n and d are positive integers.
If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
If we list the set of reduced proper fractions for d <= 8 in ascending order
of size, we get:
1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3,
 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
It can be seen that 2/5 is the fraction immediately to the left of 3/7.
By listing the set of reduced proper fractions for d <= 1,000,000
in ascending order of size, find the numerator of the fraction immediately
to the left of 3/7.

Analisi:
Sicuramente per avere un risultato della frazione tendente a 3/7, dovremo
partire comunque dal numero più grande al denominatore, cioè 10^6, dopodichè
trovare quel numeratore per il quale si abbia come risultato una frazione minore
di 3/7. Appena trovata la coppia di numeri con risultato minore e molto vicino a
3/7, controlleremo anche che soddisfi la condizione di avere gcd = 1 (maggiore
comune divisore).

python:

import time
from fractions import gcd
from fractions import Fraction as f

st = time.time()

def euler_71(limit):
    rate = round(float(3)/7, 6)
    for d in xrange(limit, 1, -1): # descending
        for n in xrange(limit/2, 1, -1): # descending
            if round(float(n)/d, 6) < rate:
                fract = f(n, d)
                if gcd(fract.numerator, fract.denominator) == 1:
                    return n

            
print euler_71(1000000)
print time.time() - st

In realtà esiste una soluzione molto più snella, terribilmente semplice:

3/7 = 0.428571

rappresentabile come:

428571/1000000

qual’è quella frazione, appena inferiore alla precedente e quindi prossima a 3/7?
Ecco appunto:

428570/1000000

vediamo se come gcd ci siamo?

python:

>>> from fractions import gcd
>>> from fractions import Fraction as f
>>> fr = f(428570, 1000000)
>>> gcd(fr.numerator, fr.denominator)
1

bingo!

Quindi un sciocca soluzione one-line, potrebbe essere:

>>> (round(float(3)/7, 6)*10**6)-1
428570.0