Project Euler #3

project_euler_0003.md

from math import sqrt

def is_prime(N, primes):
    for p in primes:
        if p*p > N:
            break
        if N % p is 0:
            return False

    return True


def init_primes(max_n, primes=[2, 3, 5, 7]):
    for i in range(11, max_n+1, 2):
        if is_prime(i, primes):
            primes.append(i)

    return primes


primes = init_primes(int(sqrt(600851475143)))
for p in primes[::-1]:
    if 600851475143 % p is 0:
        print(p)
        break

output = 6857

另解：Fermat's factorization method

不適用於大量測資的情況。

原理

參 這裡 及 維基

另外，平方數的後兩位只有 22 種可能，證明請參 式子(9) ~ (13) 。

Code

from math import sqrt

def is_square_number(N):
    if (N % 100) not in [1, 4, 9, 16, 25, 36, 49, 64, 81, 0, \
                        21, 44, 69, 96, 56, 89, 24, 61, 41, 84, \
                        29, 76]:
        return False
    return (int(sqrt(N))**2 == N)

def factorize(N):
    if N < 0:
        return
    if N is 1:
        print(1)
    elif N % 2 is 0:
        print(2)
        factorize(N // 2)
    else:
        x = int(sqrt(N)) + 1
        end = N // 2
        while x <= end and is_square_number(x*x - N) is False:
            x = x + 1
        y = int(sqrt(x*x - N))

        a = x + y
        b = x - y

        if a is 1 or b is 1:
            print(a*b)
        else:
            factorize(a)
            factorize(b)

factorize(600851475143)