2022-07-03 / Bartłomiej Kurek
Project Euler #5: Smallest multiple

Problem 5: Smallest multiple

2520 is the smallest number that can be divided by each of the numbers
from 1 to 10 without any remainder.

What is the smallest positive number that is evenly divisible by all of
the numbers from 1 to 20?

Kod: repozytorium git.


Strategia

Podejście naiwne

W naiwnym podejściu moglibyśmy iterować liczby aż znajdziemy taką, która po prostu dzieli się bez reszty przez wszystkie liczby z przedziału [1, 20]. Możemy nawet skorzystać z wiedzy, iż 2520 jest najmniejszą liczbą podzielną przez wszystkie liczby ze zbioru [1, 10]. Kod wyglądałby następująco:

    while(!found && i < ULONG_MAX) {
        if(i % 20 == 0 && /* 10, 5, 4, 2 */
           i % 19 == 0 &&
           i % 18 == 0 && /* 9, 6, 3 */
           i % 17 == 0 &&
           i % 16 == 0 && /* 8, 4 */
           i % 15 == 0 && /* 5, 3 */
           i % 14 == 0 && /* 7, 2 */
           i % 13 == 0 &&
           i % 12 == 0 && /* 6, 4, 3, 2 */
           i % 11 == 0)
        {
            found = i;
        }

        i++;
    }

Jak okaże się w dalszej części artykułu - nie jest to rozwiązanie optymalne.

Faktoryzacja

Podejściem optymalnym będzie rozkład na czynniki pierwsze liczb z zadanego przedziału i obliczenie przedmiotowej najmniejszej wielokrotności. Dla przypomnienia - przykładowy rozkład liczby 60:

  60
  |  \
  2  (30)
  |    \
  2    (15)
  |      \
  3      (5)
  |        \
  5        (1)

Czynniki to (2, 2, 3, 5). Jeśli przemnożymy je przez siebie, otrzymamy właśnie liczbę 60.

Wejściem naszego algorytmu jest jednak górny zakres przedziału, którego wszystkie elementy mają być dzielnikami poszukiwanej liczby. To właśnie dla liczb tego zakresu będziemy wykonywać faktoryzacji.

Algorytm zatem będzie działał następująco:

  • dla każdej liczby z przedziału [1, n] będziemy przeprowadzać rozkład na czynniki pierwsze (factorize())
  • otrzymane czynniki będziemy składować w liście pod warunkiem, że dla kolejnej faktoryzowanej liczby otrzymamy więcej czynników niż znamy do tej pory (merge_factors())
  • w ten sposób uzyskamy wszystkie czynniki dla liczb w zadanym przedziale
  • mnożąc te wszystkie czynniki przez siebie wyznaczymy najmniejszą wielokrotność

Przykład: rozpiszmy krok po kroku wyznaczenie najmniejszej liczby, którą dzielą wszystkie liczby z przedziału [1, 7].

# n = 7
i == 1      -> ... skip,
i == 2      -> Factors: [2]         merge: (new factor: 2)
i == 3      -> Factors: [3]         merge: (new factor: 3)
i == 4      -> Factors: [2, 2]      merge: (new factor: 2)
i == 5      -> Factors: [5]         merge: (new factor: 5)
i == 6      -> Factors: [2, 3]      no merge: 2, 3 already known
i == 7      -> Factors: [7]         merge: (new factor: 7)

Pomijamy oczywiście przypadek dla dzielnika wynoszącego 1.
Uzyskamy czynniki: [2, 3, 2, 5, 7], zatem najmniejszą liczbą podzielną przez wszystkie liczby z przedziału [1, 7] będzie:

[1]> (* 2 3 2 5 7)
420

Implementacje

Python

#!/usr/bin/env python

# Inspired by:
# http://code.jasonbhill.com/python/project-euler-problem-5/

import sys


def factorize(n):
    factors = []
    while n % 2 == 0:
        factors.append(2)
        n /= 2
    f = 3
    while n != 1:
        while n % f == 0:
            factors.append(f)
            n /= f
        f += 2
    return factors


def merge_factors(factors, new_factors, number):
    if not factors:
        factors.extend(new_factors)
        return

    uniq = set(new_factors)
    for f in uniq:
        if (new_factors.count(f) - factors.count(f)) > 0:
            factors.append(f)


def solution(n):
    factors = []
    for i in range(1, n + 1):
        merge_factors(factors, factorize(i), i)
    product = 1
    for f in factors:
        product *= f
    return product


if __name__ == "__main__":
    n = int(sys.argv[1]) if len(sys.argv) > 1 else 20

    import time

    start = time.perf_counter()
    product = solution(n)
    elapsed = (time.perf_counter() - start)
    print("{}, {:.12f}".format(product, elapsed))

Pod wskazanym w komentarzu w kodzie linkiem był kiedyś blog, który obecnie zdaje się nie działać.

Analiza kodu:

  • factorize(): dokonujemy rozkładu liczby na czynniki. Dopóki liczba jest podzielna przez 2, zbieramy dwójki jako czynnik i dzielimy liczbę przez 2. Następnie szukamy pozostałych czynników, aż do momentu, kiedy n będzie równe 1.
  • merge_factors(): uzupełniamy zebrane czynniki o nowe wystąpienia
  • solution(): funkcja, która dla każdej liczby z przedziału dzielników przeprowadza faktoryzację i zbieranie czynników, i która ostatecznie wyznacza oiszukiwaną minimalną liczbę

C (wersja "slow", nieoptymalna)

W wersji nieoptymalnej (naiwnej) moglibyśmy po prostu iterować kolejne liczby aż odnajdziemy pierwszą z nich, która jest podzielna przez wszystkie liczby z zadanego przedziału.

Dla przedmiotowego przedziału [1, 20] moglibyśmy zrealizować to w następujący sposób.

unsigned long long solution_slow(unsigned int unused)
{
    unsigned long long found = 0;
    unsigned int i = 2520;

    while(!found && i < UINT_MAX) {
        if(i % 20 == 0 && /* 10, 5, 4, 2 */
           i % 19 == 0 &&
           i % 18 == 0 && /* 9, 6, 3 */
           i % 17 == 0 &&
           i % 16 == 0 && /* 8, 4 */
           i % 15 == 0 && /* 5, 3 */
           i % 14 == 0 && /* 7, 2 */
           i % 13 == 0 &&
           i % 12 == 0 && /* 6, 4, 3, 2 */
           i % 11 == 0)
        {
            found = i;
            return found;
        }

        i++;
    }

    return found;
}

Po krótkich testach od razu zauważymy problem wydajności tego rozwiązania:

$ time ./slow.bin 10
2520

real    0m0.001s
user    0m0.001s
sys 0m0.000s

$ time ./slow.bin 20
232792560

real    0m0.971s
user    0m0.967s
sys 0m0.004s

Duży skok czasu. Dodatkowo, jeśli chcielibyśmy ten algorytm zaimplementować generycznie, zdalibyśmy się na dodatkową pętlę, która sprawdza kolejne dzielniki.
Byłoby to rozwiązanie jeszcze nieco wolniejsze:

unsigned long long solution_slow_generic(unsigned int n)
{
    unsigned long long found = 0;
    unsigned int i = 0;

    while(!found && i < UINT_MAX) {
        found = i;
        for(unsigned int f = 2; f <= n; f++) {
            if(i % f != 0) {
                found = 0;
                break;
            }
        }

        i++;
    }

    return found;
}

Testujemy:

$ time ./slow.bin 20
232792560

real    0m1.185s
user    0m1.185s
sys 0m0.000s

C (wersja "fast", faktoryzacja)

Zacznijmy od wylistowania potrzebnych funkcji:

  • factorize()
  • merge_factors()
  • solution_fast()

Faktoryzacja to prosty algorytm arytmetyczny, szczegóły widzieliśmy już wcześniej. Nieco większym wyzwaniem jest realizacja merge_factors(). Zamiast tablic wygodniej byłoby użyć jakiegoś kontenera typu vector. Taki właśnie zaimplementujemy we własnym zakresie.

Typ danych vector_t:

struct vector_t
{
    int *array;
    unsigned int max_size;
    unsigned int size;
};

Nasz wektor będzie w sobie przechowywał ograniczenia rozmiaru oraz automatycznie śledził będzie liczbę elementów, które się w danym obiekcie wektora znajdują.
Dla tego wektora możemy zaimplementować zestaw funkcji, których szczegółowego kodu już nie będziemy tutaj przedstawiać. Ich implementację można znaleźć w repozytorium.

struct vector_t *vector_create(int max_size);
void vector_destroy(struct vector_t **vec);
void vector_push(struct vector_t *vec, int n);
void vector_cat(struct vector_t *dest, const struct vector_t *src);
int  vector_count(const struct vector_t *vec, int n);
void vector_print(const struct vector_t *vec);

Implementacja

Funkcja factorize(): analogicznie jak w Python. Czynniki dodajemy do wektora.

void factorize(int n, struct vector_t* vec)
{
    int f = 2;
    while(n % 2 == 0) {
        vector_push(vec, 2);
        n /= 2;
    }

    f = 3;
    while(n != 1) {
        while(n % f == 0) {
            vector_push(vec, f);
            n /= f;
        }
        f += 2;
    }
}

Funkcja merge_factors(): otrzymuje dwa wektory i do pierwszego z nich dodaje nowe wystąpienia poszczególnych czynników.

void merge_factors(struct vector_t *all_factors, const struct vector_t *new_factors)
{
    if(!all_factors->size) {
        vector_cat(all_factors, new_factors);
        return;
    }

    unsigned int i;
    for(i = 0; i < new_factors->size; i++) {
        /* they are sorted */
        if(i && (new_factors->array[i] == new_factors->array[i - 1])) {
            continue;
        }
        int value = new_factors->array[i];
        int diff = (vector_count(new_factors, value) - vector_count(all_factors, value));

        if(diff > 0) {
            vector_push(all_factors, value);
        }
    }
}

Funkcja solution_fast(): tworzymy wektor na wynikowe czynniki, dla liczb w zadanym przedziale wykonujemy rozkład i zbieramy czynniki. Na końcu wyznaczamy poszukiwaną liczbę mnożąc wszystkie czynniki.

unsigned long long solution_fast(unsigned int n)
{
    struct vector_t *all_factors = vector_create(128);

    unsigned int i;
    for(i = 1; i <= n; i++) {
        struct vector_t *factors = vector_create(32);
        factorize(i, factors);
        merge_factors(all_factors, factors);
        vector_destroy(&factors);
    }

    // vector_print(all_factors);
    unsigned long long product = 1;

    for(i = 0; i < all_factors->size; i++) {
        product *= all_factors->array[i];
    }

    return product;
}

Porównanie prędkości

Czas podany w sekundach:

|   (N) |         Python |            C (slow) |           C (fast) |
|-------+----------------+---------------------+--------------------|
|    10 | 0.000033322984 |  0.0000142990000000 | 0.0000167350000000 |
|    20 | 0.000051133014 |  1.1986012340000001 | 0.0000181550000000 |
|    30 | 0.000071857998 | 22.1204933890000000 | 0.0000204660000000 |
|    40 | 0.000096797012 | -                   | 0.0000215690000000 |
|    50 | 0.000124501996 | -                   | -                  |
|    60 | 0.000152404013 | -                   | -                  |
|    70 | 0.000183138996 | -                   | -                  |
|    80 | 0.000218809990 | -                   | -                  |
|    90 | 0.000259660010 | -                   | -                  |
|   100 | 0.000287114002 | -                   | -                  |
|   200 | 0.000785820012 | -                   | -                  |
|   500 | 0.003415011975 | -                   | -                  |
|  1000 | 0.012154298980 | -                   | -                  |
| 10000 | 0.781079200999 | -                   | -                  |
|       |                |                     |                    |

Uwagi:

  • Przykładowy wynik dla Python (10000).
  • W C z typem unsigned long long dojdziemy maksymalnie do przedziału [1, 47].
  • Wariant naiwny w C ("slow") zajmuje już 22 sekundy przy przedziale [1, 30].