Hledání prvočísel, Eratosthenovo síto

Zadání je pro obě části stejné a může znít takto: Najděte nejefektivnější algoritmus pro hledání prvočísel v zadaném limitu. Pro práci s číselným rozsahem využijte prvků indexovaného pole, kde jednotlivé klíče pole budou reprezentovat celočíselnou testovanou posloupnost a hodnoty pole budou určovat prvočíselnost oněch klíčů.

• Vstup:
  100
• Výstup:
  2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97

Dedukce - hledání hrubou silou

Pomocí brainstormingu zde můžeme hned na začátku přijít na definici prvočísla. Položme tedy studentům otázku: Co je to prvočíslo?

Je pravděpodobné, že v tomto případě dojdeme rychle ke správné odpovědi, jelikož se jedná o známou definici. Pro kontrolu ještě můžeme přidat další otázku: Kdo si myslí, že číslo 80 je prvočíslo?

Kdo se přihlásí, zjevně nepochopil definici. Je proto nutné vysvětlit, proč tomu tak není. Pokračujme v otázkách i dále: Jak spolehlivě zjistíme, že nějaké číslo je prvočíslo?

Někdo se vytasí s poučkami o ciferných součtech, sudosti a dalších poučkách pro jednotlivé dělitele, ale ty nepokryjí všechny případy. Je třeba diskusi směřovat k závěru, že se mně číslo N nepodaří vydělit beze zbytku čísly od 2 do N-1. Možná někdo už přijde s odmocninou (viz níže). Abychom to ale nasměrovali žádoucím směrem, můžeme vyslovit požadavek, že chceme nejjednodušší algoritmus z hlediska zápisu, resp. náročnosti na pochopení, takže návrh na odmocninu a podobné zlepšováky oceníme, ale prozatím vynecháme.

Po brainstormingu definice prvočísla pracují studenti samostatně na svých programech. Pravděpodobným výstupem této samostatné práce jsou řešení hrubou silou (více, či méně efektivní).

Efektivita se odvíjí od způsobu, jakým je při zjišťování testována dělitelnost beze zbytku:

• Všechna menší čísla než je zjišťované číslo N – velká zátěž
• Všechna menší čísla než je N/2 – velká až střední zátěž
• Všechna menší čísla než je odmocnina z N – malá zátěž

Další důvody negativního ovlivnění rychlosti:

• Při nalezení dělitele není cyklus ukončen a nadále testuje zbytek možných dělitelů
• Možný dělitel je testováni od největšího čísla po nejmenší

V dalším kroku je vhodné promítnout přiložený program 1_brute-force v preferovaném jazyce a provést diskusi. Pro začátek mohou studenti hledat odlišnosti v jejich a promítaném programu (ve skupinách, či dohromady). Algoritmus je velice málo efektivní a zahrnuje všechny nedostatky shrnuté výše.

#include <iostream>
using namespace std;

int main() {

  // Deklarace a inicializace
  const int vel = 100; // Inicializace limitu - do kolika cisel hledame
  bool pole[vel + 1]; // Deklarujeme velikost pole 
  pole[0] = pole[1] = 0; // 0 a 1 nejsou prvocisla

  // Pripravime si pole - projdeme ho
  for(int i = 2; i <= vel; i++) {

    // Naplnime ho hodnotou true
    pole[i] = 1;

  }

  // Projdeme pole
  for(int i = 2; i <= vel; i++) {

    // Kazde cislo budeme zkouset delit cisly mensimi nez je ono samotne a vetsi nez nebo rovno 2
    for(int k = i - 1; k >= 2; k--) {

      // Pokud najdeme delitele jakehokoliv delitele
      if(i % k == 0) {

        // Cislo neni prvocislo
        pole[i] = 0;
        
      }

    }

  }
  
  // Projdeme pole
  for(int i = 0; i <= vel; i++) {

    // Pokud je cislo prvocislem
    if(pole[i] == 1) {

      // Vypiseme ho
      cout << i << endl;

    }

  }

}

// Deklarace a inicializace
const vel = 100; // Inicializace limitu - do kolika cisel hledame
var pole = []; // Deklarujeme pole
pole[0] = pole[1] = false; // 0 a 1 nejsou prvocisla

// Pripravime si pole - projdeme ho
for(var i = 2; i <= vel; i++) {

  // Naplnime ho hodnotou true
  pole[i] = true;

}

// Projdeme pole
for(var i = 2; i <= vel; i++) {

  // Kazde cislo budeme zkouset delit cisly mensimi nez je ono samotne a vetsi nez nebo rovno 2
  for(var k = i - 1; k >= 2; k--) {

    // Pokud najdeme delitele jakehokoliv delitele
    if(i % k == 0) {

      // Cislo neni prvocislo
      pole[i] = false;

    }      

  }

}

// Projdeme pole
for(var i = 0; i <= vel; i++) {

  // Pokud je cislo prvocislem
  if(pole[i] == true) {

    // Vypiseme ho
    console.log(i); 

  }

}

# Deklarace a inicializace
vel = 100 # Inicializace limitu - do kolika cisel hledame
pole = [] # Deklarujeme pole
pole.extend([False, False]) # Do pole pridame na indexy 0 a 1 false - nejsou prvocisla

# Pripravime si pole
for i in range(2, vel + 1):

  # Pridame pokazde novy prvek pole, ktery nastavime na True
  pole.append(True)

# Projdeme pole
for i in range(2, vel + 1):

  # Kazde cislo budeme zkouset delit cisly mensimi nez je ono samotne a vetsi nez nebo rovno 2
  for k in range(i - 1, 2, -1):

    # Pokud najdeme delitele jakehokoliv delitele
    if i % k == 0:

      # Cislo neni prvocislo
      pole[i] = False

# Projdeme pole
for i in range(0, vel + 1):

  # Pokud je cislo prvocislem
  if pole[i] == True:

    # Vypiseme ho
    print(i)

Je pravděpodobné, že studenti nedokáží přijít na minimální počet testovaných čísel, tedy odmocnina z testovaného čísla N.

Pokud přišel někdo z nich na limit N/2, položme otázku: Můžeme jít ještě níž?

V případě, že všichni testovali až do čísla N-1: Je nutné testovat všechna čísla?

Když žádná z těchto otázek nepodpoří navedení na správnou cestu, nic se neděje a prostě jim to vysvětlíme. Jestliže postupuji od 2 dále a poprvé se mi podaří N vydělit nějakým A, pak N / A má celočíselný výsledek, který označíme B, a toto B je logicky větší než A (kdyby nebylo, našli bychom ho dříve než A) a je to taky dělitel, jež je logicky menší než N. Odmocnina je pak jasná hranice, protože čím větší je A, tím menší je B, ale stále je A < B a tou odmocninou najdeme hranici, kdy A = B (i když jen teoretickou – často nevyjde celočíselně). K odbornému vysvětlení je možné provést důkaz:

Nechť N je součin libovolných přirozených čísel A a B, která jsou větší než 1 a zároveň menší než N:
N = A * B
1 < A ≤ B < N

Výraz A je větší nebo rovno B vynásobíme A a B, vzniknou dvě rovnice:
A2 ≤ A * B
A * B ≤ B2

Následně můžeme nerovnice sloučit a dosadit N:
A2 ≤ N ≤ B2

Jelikož počítáme v množině přirozených čísel, můžeme výsledný výraz odmocnit, z čehož vyplívá, že není možné, aby obě čísla byla zároveň větší, než je odmocnina z N:
A ≤ √N ≤ B

Nemělo by se ale zapomenout i na vysvětlení ukončení cyklu a proč je lepší testovat dělitelnost od dvojky výše. Takovou otázku můžeme položit opět studentům, protože odpověď je triviální. Jaká je pravděpodobnost, že číslo vydělím 2? 50 %. Jaká je pravděpodobnost, že číslo vydělím 3? 33 %.

Následně by měli studenti upravit svůj program podle diskutovaných výkonových zlepšení. V případě nutnosti ukážeme způsob, jak se ve zvoleném programovacím jazyce implementuje matematická knihovna a práci s ní.

Otázky do diskuse

1. Najdete rozdíly mezi vaším a promítaným programem?
2. Mají tyto odchylky vliv na celkový výkon programu?

Možné problémy

• Student nedokáže implementovat algoritmus definice prvočísla. – Explicitní navedení na libovolném prvočísle pomocí důkazu. Proč je číslo 7 prvočíslo?
• Student neví, jak do programu zahrnout čísla 0 a 1. – Vysvětlit všem, že se jedná o ukázku toho, že nastanou situace, kdy není možné algoritmus použít na všechnu množinu čísel a je nutné některé prvky upravit (v našem případě vyjmout) „natvrdo“.

Eratosthenovo síto

Po předchozí diskusi představíme Eratosthenovo síto. Je dobré použít přiložený soubor eratosthenovo-sito.gif, který znovu poslouží k dedukci a následné debatě nad způsobem fungování síta.

Princip Eratosthenovo síta je následující: Pro čísla z rozsahu 2 až N vezmeme nejmenší číslo (v prvním kroku 2). Toto číslo je zároveň prvočíslem. Odstraníme všechny jeho násobky v rozsahu. Poté pokračujeme dalším číslem (v druhém kroku 3) a provedeme stejnou operaci. Takto pokračujeme až do doby, kdy je další číslo odmocninou z N (důvod byl popsán v první fázi), všechna zbylá čísla jsou prvočísla.

Toto je možné vysvětlit také pomocí přiloženého souboru: eratosthenovo-sito.png.

Studenti mají za úkol samostatně implementovat Eratosthenovo síto. Pro rychlejší práci mohou použít již napsaný program z první fáze a pouze v něm opravit cyklus, který určuje prvočíselnost. Výsledný program je přiložen jako 2_eratosthenovo-sito.

#include <iostream>
#include <cmath>
using namespace std;


int main() {

    // Deklarace a inicializace
    const int vel = 100; // Inicializace limitu - do kolika cisel hledame
    bool pole[vel + 1]; // Deklarujeme velikost pole 
    pole[0] = pole[1] = 0; // 0 a 1 nejsou prvocisla

    // Pripravime si pole - projdeme ho
    for(int i = 2; i <= vel; i++) {

        // Naplnime ho hodnotou true
        pole[i] = 1;

    }

    // Projdeme kazdy prvek pole az do odmocniny z velikosti zaokrouhlene dolu
    for(int i = 2; i <= sqrt(vel); i++) {

        // Pokud najdeme prvocislo
        if(pole[i] == 1) {

            // Projdeme si pole vzdy o nasobky prvocisla
            for(int j = i * i; j <= vel; j += i) {

                // Kazdy nasobek je delitelny prvocislem
                pole[j] = 0;

            }

        }

    }

    // Projdeme pole
    for(int i = 0; i <= vel; i++) {

        // Pokud je cislo prvocislem
        if(pole[i] == 1) {

            // Vypiseme ho
            cout << i << endl;

        }

    }

}

// Deklarace a inicializace
const vel = 100; // Inicializace limitu - do kolika cisel hledame
var pole = []; // Deklarujeme velikost pole 
pole[0] = pole[1] = false; // 0 a 1 nejsou prvocisla

// Pripravime si pole - projdeme ho
for(var i = 2; i <= vel; i++) {

    // Naplnime ho hodnotou true
    pole[i] = true;

}

// Projdeme kazdy prvek pole az do odmocniny z velikosti zaokrouhlene dolu
for(var i = 2; i <= Math.sqrt(vel + 1); i++) {

    // Pokud najdeme prvocislo
    if(pole[i] == true) {

        // Projdeme si pole vzdy o nasobky prvocisla
        for(var j = i * i; j <= vel; j += i) {

            // Kazdy nasobek je delitelny prvocislem
            pole[j] = false;

        }

    }

}    

// Projdeme pole
for(var i = 0; i <= vel; i++) {

    // Pokud je cislo prvocislem
    if(pole[i] == true) {

        // Vypiseme ho
        console.log(i);

    }

}

# Importujeme matematickou knihovnu
import math

# Deklarace a inicializace
vel = 100 # Inicializace limitu - do kolika cisel hledame
pole = [] # Deklarujeme pole
pole.extend([False, False]) # Do pole pridame na indexy 0 a 1 false - nejsou prvocisla

# Pripravime si pole
for i in range(2, vel + 1):

    # Pridame pokazde novy prvek pole, ktery nastavime na True
    pole.append(True)

# Projdeme kazdy prvek pole az do odmocniny z velikosti zaokrouhlene dolu
for i in range(2, int(math.sqrt(vel + 1)) + 1):

    # Pokud najdeme prvocislo
    if pole[i] == True:

        # Projdeme si pole vzdy o nasobky prvocisla
        for k in range(i * i, vel + 1, i):

            # Kazdy nasobek je delitelny prvocislem
            pole[k] = False

# Projdeme pole
for i in range(vel + 1):

    # Pokud je cislo prvocislem
    if pole[i] == True:

        # Vypiseme ho
        print(i)

Pro zpestření výuky je součástí také spustitelný soubor 3_porovnani-vykonu.exe (samotný kód programu se nachází v souboru 3_porovnani-vykonu.cpp). Program testuje rychlost použitých programů pomocí knihovny chrono (bez vypisování prvočísel).

#include <iostream>
#include <cmath>
using namespace std;


int main() {


    // Priprava programu
    const int vel = 10000;
    bool pole[vel + 1];
    pole[0] = pole[1] = 0;
    for(int i = 2; i <= vel; i++) {
        pole[i] = 1;
    }

    // Zaznamenani casu zacatku
    chrono::steady_clock::time_point begin = chrono::steady_clock::now();

    // Pocitani
    for(int i = 2; i <= vel; i++) {
        for(int k = i - 1; k >= 2; k--) {
            if(i % k == 0) {
                pole[i] = 0;
            }          
        }
    }

    // Zaznamenani casu konce
    chrono::steady_clock::time_point end = chrono::steady_clock::now();

    // Vypsani rozdilu casu
    cout << "Cas nutny k vypoctu prvocisel v rozsahu <2, " << vel << "> pomoci prvniho programu = " << chrono::duration_cast(end - begin).count() << " mikrosekund" << endl;


    // Pripraveni pole znovu
    for(int i = 2; i <= vel; i++) {
        pole[i] = 1;
    }

    // Zaznamenani casu zacatku
    chrono::steady_clock::time_point begin2 = chrono::steady_clock::now();

    // Pocitani
    for(int i = 2; i <= sqrt(vel); i++) {
        if(pole[i] == true) {
            for(int k = i * i; k <= vel; k += i) {
                pole[k] = false;
            }
        }
    }

    // Zaznamenani casu konce
    chrono::steady_clock::time_point end2 = chrono::steady_clock::now(); 

    // Vypsani rozdilu casu   
    cout << "Cas nutny k vypoctu prvocisel v rozsahu <2, " << vel << "> pomoci druheho programu (Eratosthenovo sito) = " << chrono::duration_cast(end2 - begin2).count() << " mikrosekund" << endl;
    
    // Podil dvou vyslednych casu
    cout << "Druhy program je priblizne " << (end - begin)/(end2 - begin2) << "x rychlejsi" << endl;

    // Prevence uzavreni exe po vykonani mereni
    system("pause");

}

Otázky do diskuse

1. Jakým způsobem funguje Eratosthenovo síto?
2. Jaká je horní mez zkoumaných čísel?
3. Kolikrát rychlejší bude Eratosthenovo síto, než program promítaný v první fázi?

Možné problémy

• Se znalostmi z první části se nepředpokládají dílčí problémy.