Programowanie w języku Delphi » Gotowce

Silnia na przykładzie szybkiego mnożenia

Obliczanie silnii dowolonej liczby
Prezentuję algorytm szybkiego mnożenia długich liczb, a jako przykład użycia - wykożystałem obliczenie silnii. Algorytm narodził się w czasach ZX Spectrum+ - wtedy pisało się programy tylko w assemblerze i to tak, aby kod był jaknajszybszy.

Składnia programu - język basic - ostatni krzyk mody - ibasic. Pomijam jego zalety bo to nie dotyczy algorytmu.

na początku trzeba założyć dwie tabele o pewnej ilości bajtów
char - określa zmienną jednobajtową, MaxDigits - liczba znaków w tabelach. Dla 65536 znaków można na spokojnie obliczać silnię 10000
tabela A to wynik, a tabela B to coś w stylu reszty. Jeden bajt w tabeli to jedna cyfra 0-9

const MaxDigits=65536
char a[MaxDigits],b[MaxDigits]


następnym krokiem jest podanie liczby n! - tu będzie nią N
trzeba sprawdzić czy n=0 | 1
jeśli n=0 albo 1 - wynik =n
ewentualnie czyścimy tabelę A:
RtlZeroMemory(&A,MaxDigits)

wpisujemy liczbę 2 na początek wyniku (zeby nie mnożyć przez 1 :) )

a[0]=2
LenA=1


LenA - zmienna przechowująca ilość cyfr wyniku. Jest potrzebna aby przyspieszyć obliczanie
czyścimy obie tabele

wchodzimy w pętlę mnożenia dla liczb 3->N
timestart=GetTickCount()
for BC=3 to N
  mul=BC /* mul będzie modyfikowane */
 
  RtlZeroMemory(&B,LenA) /* czyścimy 'resztę' */
 
  do                 /* mnożenie długiej liczby a[] przez mul */
    if (mul&1)    /* jeśli mul jest nieparzyste */
      AddAtoB() /* dodaj wynik do 'reszty' */
      mul--        /* zmniejsz mul o 1 */
    else
      RLCA()      /* pomnóż wynik przez 2 */
      mul/=2      /* podziel mul przez 2 */
    endif
  until mul<2    /* powtarzaj if N>1 */
 
  AddBtoA()    /* dodaj 'resztę' do wyniku */
next BC
 
time=GetTickCount()-timestart

zmienna time ma ilość ms wykonywania się obliczeń - podziej ją przez 1000, wyświetl wynik do dwóch miejsc po przecinku
print "Time:"+using("#.##",time/1000.0),"s"

no i teraz trzeba wyświetlić wynik lub wpisać do pliku. Zaczynamy od ostatniego bajtu (na pozycji LenA-1) kończąc na pozycji 0
for x=LenA - 1 to 0 step -1
  print (albo write do pliku) chr$(a[x]+48)
next x


a teraz funkcje:
Mnożenie wyniku przez 2

sub RLCA
        Carry  = 0
        for HL=0 to LenA-1
                a[HL] = a[HL]*2 + Carry
                Carry = (a[HL]>9)
                if Carry
                        a[HL] = a[HL]-10
                endif
        next HL        
        a[HL]+=Carry
        LenA+=Carry
 
        return
endsub


dodanie wyniku do reszty

sub AddAtoB
        Carry=0        
        for HL=0 to LenA-1
                b[HL]= b[HL] + a[HL] + Carry
                Carry = (b[HL]>9)
                if Carry
                        b[HL] = b[HL]-10
                endif
        next HL
        b[HL]+=Carry
        LenA+=Carry        
        return

endsub

dodanie reszty do wyniku

sub AddBtoA
        Carry  = 0
        for HL=0 to LenA-1
                a[HL]= a[HL] + b[HL] + Carry
                Carry = (a[HL]>9)
                if Carry
                        a[HL] = a[HL]-10
                endif
        next HL
        a[HL]+=Carry
        LenA+=Carry        
        return
endsub


Carry to zmienna przechowująca przeniesienie; Carry=True jeśli ostatnia operacja dodawania/mnożenia dała wynik>9
HL i BC - zwykłe zmienne, ich nazwy wywodzą się z Z80 :)

gotowca można ściągnąć


for HL=0 to LenA-1 działa tak: zmienna HL przyjmuje wartość 0. Wykonują się kolejne rozkazy aż program napotka rozkaz NEXT HL, wtedy HL przyjmuje kolejną wartość (tutaj o 1 większą) i dopuki HL<(LenA-1) - pętla powraca do następnego rozkazu po FOR


Carry = (a[HL]>9) działa tak: jeśli wyrażenie po znaku równości będzie TRUE to zmienna CARRY przyjmie wartość 1 (TRUE), w przeciwnym razie przyjmie wartość 0 (FALSE)


a[HL] - jeśli HL=0 to chodzi o pierwszy bajt z tabeli A - HL jest jakby wskaźnikiem/pointerem do kolejnych bajtów w tabeli. Pierwszy bajt w obu tabelach jest bajtem (cyfrą wyniku lub reszty) najmniej znaczącym (jedności)


if Carry - działa taksamo jak if Carry=1


LenA+=Carry - jeśli wynik ma jedną cyfrę, załóżmy że chwilowy wynik=5 - pomnożenie cyfry 5 przez 2 da liczbę 10 (czyli 2 cyfry). i zmienna Carry zgodnie z założeniem będzie równa 1 [patrz linijka Carry = (a[HL]>9)].  Każdy rozkaz NEXT w basicu ZAWSZE najpierw zmienia zmienną pętlową a potem sprawdza warunek podany w FOR, więc po wyjściu z FOR-NEXT zmienna HL wskazuje na pierwszy wolny bajt w tabeli wyniku A (z lewej strony, tutaj na drógą cyfrę, której zresztą jeszcze nie ma). Pierwsza cyfra wyniku=0 ponieważ warunek (liczba>9) został spełniony - wpisano do tabeli 10 a potem odjęto 10 i ustawiono Carry dla następnego mnożenia przez 2 (popularne JEDEN DALEJ). Ale pętla przebiegała dla liczb 0-1 więc wykonała się tylko raz - teraz wystarczy dopisać jedynkę z lewej strony wyniku (b[HL]+=Carry) albo b[HL]=Carry i zwiększyć ważną liczbę cyfr wyniku o jeden: LenA+=Carry (jeśli Carry=0 to LenA nie zwiększy się)

czyli zmienna Carry pełni dwuznaczną rolę - ustawia się gdy dodawanie/mnożenie przekroczyło 9, a zarazem po wyjściu z pętli wiadomo czy trzeba dopisać cyfrę do wyniku


wyniki:
10000! - pentium4 2.4GHz - 91s, athlonXP 2.1GHz - 77s

9 komentarzy

adameczek4 2006-05-22 23:44

kurka , utrudniłem se życie :P <lol>

adameczek4 2006-05-22 23:35

Dla kompletnych alików, przedstawiam kod, w ktorym smai decydujecie do jakiego n silna ma byc obliczona, znaczy sie nadajecie taka możliwość programowi, poprzez dopisanie wiekszych wartość , liczby elementów naturalnych do tablicy.

Oto kod:

program Silnia_Rekurencyjna;

{$APPTYPE CONSOLE}

uses
  SysUtils;

   {tablica kolejnych liczb naturalnych od 1 do 15}
   
const A: Array[1..15] of Longint=(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15);

var
i ,n : Integer;                         //okreslenie zmiennych
silnia : Int64;                         //dla ogromnych wartości silni zmiencie typ na REAL

begin
silnia:=1;                             //nadanie wartosci poczatkowej
Writeln('Podaj liczbe naturalna: ');
Readln(n);                             //pobranie zmiennej n

for i:=1 to n do                       //wzrastanie potegowania [obliczania silni]
silnia:= silnia  * A[i];

Writeln('Silnia: ', silnia);           //wyswietlenie wyniku
Readln(n);

end.

Kodzik prosty, przyjemny i nie ma żadnych problemów z 0!  :P . Pozdro for all !

VsMaX 2006-03-02 21:57

Chyba dalo sie to prosciej zrobic :)
Ok zrobilem wlasnna silnie ponizej kod
--Cut Here--
function Silnia(liczba : integer; i : integer; n : integer):integer;
begin
  liczba := 1;
  for i := 1 to n do
  begin
  liczba := liczba * i;
  end;
  result := liczba;
end;

procedure TForm1.Button1Click(Sender: TObject);
var
d : integer;
begin
  d := Silnia(1,1,StrToInt(Edit1.Text));
  Memo1.Lines.Add(IntToStr(d))

end;
--Cut Here--
To caly kod :P oczywiscie zeby silnia dzialala trzeba wstawic na forme komponent memo, jeden edit i button :)

Tomuslaw 2004-11-13 21:13

" (..) jeśli n=0 albo 1 - wynik =n (...) "

Otóż nie, bo jeśli dla 0!=1 i dla 1!=1...

matys 2004-09-27 07:52

Wielkie dzieki za maila;]

sapero 2004-09-26 15:28

niestety ostatnio pisałem w C jakieś 6 lat temu na spectrumie, C tylko znam ze strony tłumaczenia z C na basica. Operacje na tablicach cyfr raczej są proste (tu jest tylko dodawannie, a mnożenie przez 2 to też dodawanie)

ale jeśli komuś ciężko się czyta taki kod, to mogę to przerobić na czysty algorytm.

A to że dział nie ten - też już wiem :) ale nie widzę możliwości przeniesienia tematu

matys 2004-09-25 17:59

Mam prosbe moglbys tez pokazac implementacje w C++?? Bo o delphi nie mam pojecia a chetnie kod bym przeczytal;]

ŁF 2004-09-20 20:21

to chyba nie ten dział; to pownno być w "z pogranicza"