e-mail    Debatní kniha    Mapa stránek    Hlavní  
 perličky 
 

Prvočísla

Pavel Pitner

Ve škole žáci často řeší problém zjištění a výpisu všech prvočísel od dvou do zadaného N. Obvykle používají metodu, při níž se postupně berou všechna čísla X počínaje dvojkou a konče N a tato se dělí všemi čísly od 2 do odmocniny z X a zjišťují se zbytky po dělení. Je-li některý zbytek roven nule, pak číslo není prvočíslo, v opačném případě ano. Při řešení se mohou použít různé optimalizace (např. nedělíme číslem 2, případně dělíme jen prvočísly do odmocniny z X, atd.).}

Kromě tohoto algoritmu však existuje ještě jeden, který je nesrovnatelně efektivnější. Nazývá se Erathostenovo síto a postupujeme při něm takto:

  1. vytvoříme bitové pole tak, že pro každé přirozené číslo od 1 do $N$ vyhradíme jeden bit ($N$ je horní hranice výpisu prvočísel)
  2. index pole bude uvádět číslo, hodnota bitu rovná 0 bude znamenat, že číslo je prvočíslo, 1, že není prvočíslo
  3. na začátku pole vynulujeme tak, jako by všechna čísla byla prvočísly, číslo jedna označíme jako neprvočíslo
  4. procházíme prvky pole a vyhledáme první nulový bit (v prvním průchodu to bude u indexu 2, dále 3, 5, 7 atd.)
  5. nalezené číslo (označme je P) je prvočíslo, proto na jeho místě necháme v poli nulu
  6. nyní bereme všechny prvočíselné násobky P počínaje P (P.P, P.(první prvočíslo za P), P.(druhé prvočíslo za P), atd.) a na místa těchto čísel ukládáme jedničky (zcela jistě to nejsou prvočísla). Postupujeme tak dlouho, dokud hodnota násobku nepřesáhne zadané N. (Při realizaci algoritmu vyjde lépe brát všechny celočíselné (případně je liché) násobky počínaje P.)
  7. vrátíme se k bodu 4) a cyklus opakujeme od současného P do té doby, dokud není P větší než odmocnina z N
  8. pole je správně vyplněno, jsme hotovi a chceme-li prvočísla zobrazit, vypíšeme indexy těch prvků pole, kde je nula

Obrázky názorně ukazují, jak postupujeme:

  1. takto vypadá pole po inicializaci 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
  2. první najdeme číslo 2, označíme násobky dvou (2.2, 3.2, 4.2, atd.) 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 1 0 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 * ^ ^ ^ ^ ^ ^ ^ ^
  3. další najdeme trojku, označíme násobky tří, (začneme od 3.3), atd. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 1 0 0 1 0 1 0 1 1 1 0 1 0 1 1 1 0 1 0 * ^ ^
  4. dále čísla 5, 7, ...

Uvedený program používá částečně upravené právě popsané metody Erathostenova síta (např. v poli vůbec nejsou násobky čísla 2). Umožňuje výpočet a výpis prvočísel až do jednoho miliónu. Jeho činnost je patrná z komentovaného výpisu. Po přeložení do .EXE se spustí s parametrem udávajícím horní hranici výpočtu.

Program je vytvořen a odladěn v Turbo Pascalu 6.0. Rutina výpočtu je z důvodu rychlosti celá napsána v assembleru, vstupy a výstupy hodnot jsou ponechány v Pascalu. Výstup je možné nechat provádět na obrazovku nebo přesměrovat do souboru či na tiskárnu.

Doba samotného výpočtu (bez výpisu) je na AT 286 16 MHz pro prvočísla od 1 do 100000 takřka neměřitelná a do maxima, tj. 1000000 je asi 4 sekundy.

Pozn.: Program slouží spíše jako ukázka, co a jak se také dá naprogramovat v assembleru, než pro praktické použití.

Jestli se někomu podaří výpočet dále zrychlit (alespoň o 5%) nebo dokonce vymyslet lepší algoritmus, informuje mě prosím na adrese: Pavel Pitner, Nová 6, 679 61 Letovice

{PRVOCIS.PAS}

program Prvocisla;

const BitMask : array [ 0 .. 7 ] of byte =
                  ( 1, 2, 4, 8, 16, 32, 64, 128 );
                                          { bitová maska }
      Sirka = 9;        { šířka jednoho čísla při výpisu }
      PocetNaRadek = 8;  { počet čísel vypis. na 1 řádek }

var BitPole : array [ 0 .. 62500 ] of byte;
                             { bitové pole pro prvočísla }
    N, a, c : longint;
    xN, xSqrtN, b, i : word;

procedure Vypocti; assembler;
asm
        mov     ax,seg BitPole
        mov     es,ax
        mov     di,offset BitPole
        mov     cx,xN
        shr     cx,1
        inc     cx
        xor     ax,ax                        { vynuluj AX }
        rep     stosw                 { vynuluj pole bitů }
        mov     si,ax                   { vynuluj BX a SI }
        mov     bx,ax
        inc     ah                         { do AH uloz 1 }
        add     xN,offset BitPole
@0:     shr     al,1            { obnov pův. hodn. AL, SI }
        shr     si,1           { po přičítání dvojnásobku }
@1:     inc     al         { posun na další číslo tabulky }
        rol     ah,1
        adc     si,0
        cmp     si,xSqrtN          { jen do odmocniny z N }
        ja      @E
        test    ah,[si+offset BitPole]    { je to prvoč.? }
        jnz     @1                      { ne - hledej dál }
        and     al,7       { platné jsou jen dolní 3 bity }
        mov     cl,al         { AL, SI zkopíruj do CL, DI }
        mov     di,si
        add     di,offset BitPole
        stc               { vynásob AL, SI dvěma, protože }
        rcl     al,1          { budeme uvažovat jen liche }
        shl     si,1       { násobky nalezeného prvočísla }
@2:     add     di,si           { sečti vyšší části čísla }
        jc      @0            { došlo k přetečení - konec }
        add     cl,al           { sečti nižší části čísla }
        mov     bl,cl      { uprav případný přenos přes 7 }
        and     cl,7
        shr     bl,1
        shr     bl,1
        shr     bl,1
        add     di,bx            { přičti případný přenos }
        cmp     di,xN                    { jen čísla do N }
        ja      @0
        mov     bl,cl                        { číslo bitu }
        mov     ch,[bx+offset BitMask] { vezmi masku bitu }
        or      [di],ch   { ulož příznak 'není prvočíslo' }
        jmp     @2       { opakuj pro další lichý násobek }
@E:     sub     xN,offset BitPole { obnov původní hodnotu }
end;

begin
  writeln ( 'PRVOČÍSLA V2.0  (C) 1993 Pavel Pitner' );
  writeln ( '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' );
  if ParamCount <> 1 then begin
    writeln ( 'Zadej: PRVOCIS rozsah (od 2 do 1000000)' );
    Halt ( 1 );
  end;
  Val ( ParamStr ( 1 ), N, i );
  if ( i <> 0 ) or ( N < 2 ) or ( N > 1000000 ) then begin
    writeln ( 'Špatně zadaný rozsah!' );
    Halt ( 2 );
  end;
  xN := N div 16;           { převeď na index do tabulky }
  xSqrtN := trunc ( sqrt ( N ) / 16 );   { odmocnina z N }
  Vypocti;
  writeln ('Všechna prvočísla od 2 do ', N, ':'); writeln;
  a := 1; c := 1;
  for i := 0 to xN do begin     { pro všechny prvky pole }
    b := 1;                                 { maska bitu }
    repeat
      if ( BitPole [ i ] and b ) = 0 then begin
                                 { našli jsme prvočíslo? }
        if a = 1 then
          write ( 2 : Sirka )    { dvojku místo jedničky }
        else
          write ( a : Sirka );
        if c = PocetNaRadek then begin
          writeln; c := 1;
        end
        else
          inc ( c );
      end;
      a := a + 2;                    { další liché číslo }
      b := b * 2;                        { odroluj masku }
    until ( b = 256 ) or ( a > N );
  end;
  writeln;
end.



Pascal - hlavní
Překladače
Vlastní články
Převzaté články
Věci na stáhnutí
Odkazy k tématu
BP7 buglist
Chyba Run-time 200

BASIC