[Pascal] Optymalizacja kodu

Witam:

Mam za zadanie znalezdz liczby doskonale z podanego zakresu.
Liczba doskonala to taka ktorej dzielniki po zsumowaniu rowne sa liczbie pierwotnej (np. 1 ma dzielnik tylko 1 , oraz 6 gdyz dzieli sie tylko przez 3 2 1 a 3+2+1=6)
Nie biezemy pod uwage liczby tej samej.
Program napisalem, lecz potrzebuje optymalizacji kodu by dzialal jak najszybciej.
(Moj PSOR powiedzial ze to bedzie gdzies tam gdzie zaznaczylem komentarzami {Start} i {STOP} ale nie sugerujcie sie tym az tak bardzo ;) ...)

wszystkie sugestie mile widziane
(Byc moze takie dzielniki liczb doskonalych sa liczbami pierwszymi ? - Jezeli tak to nie bedzie problemu z napisaniem kodu) plz help.

var i,j,zakres,suma:integer;

begin
   readln(zakres);
   for i:=1 to zakres do begin
     suma:=1;
     for j:=2 to i div 2 do
        {START}
       if i mod j = 0 then
           suma:=suma+j;
         {STOP}  
     if suma=i then
       write(suma,' ');
   end;


end.

Dla dowolnego k>0, jesli 2^{k-1 jest liczba pierwsza to (2}(k-1))*(2^k-1) jest liczba doskonala. To powinno ci troche pomoc :).

Oczywiście że mi pomogło :]
Stokrotne dzięki - tylko tyle, że teraz podając zakres np 500 szuka mi do wysokich liczb (wyzszych conajmniej 15 razy) nie widze niestety błędu z tym związanego :/
Jeżeli mógłby ktoś mi jeszcze raz pomóc z tym zakresem :/

var i,j,zakres,suma:longint; tablica:array[1..9999] of longint;


function potega(liczba,tyle:longint):longint;
var liczba2,ile:longint;
begin
liczba2:=liczba;
if tyle < 1 then write('') else
for ile:=1 to tyle-1 do
liczba:=liczba*liczba2;
potega:=liczba;
end;


begin
   readln(zakres);

       for i:=1 to zakres  do tablica[i]:=i;
          for i:=2 to zakres  do
            for j:=i+1 to zakres  do
              if (j mod i = 0) then tablica[j]:=0;



   for j:=2 to zakres do
    for i:=2 to zakres div 2 do
      if tablica[i]<>0 then
         if ((potega(2,j)-1) = tablica[i])  then
             write(' ',potega(2,(j-1))*(potega(2,j)-1));
 end.

W pierwszym rozwiązaniu pętla ?j? to szukanie podzielników, szukać można tylko do pierwiastka z ?i?, znajdujesz w ten sposób dwa dzielniki ?j? i ?i/j?

W pierwszym rozwiązaniu pętla ?j? to szukanie podzielników, szukać można tylko do pierwiastka z ?i?, znajdujesz w ten sposób dwa dzielniki ?j? i ?i/j?

Próbowałem przed chwilą ale coś robie nie tak bo znajduje mi tylko 1 ...
Jeżeli mógłbyś pokazać która to konkretnie linijka i co za nią wpisać byłbym wdzięczny.

Tak i co do rozwiązania drugiego nadal zostaje kwestia zakresu :/

Pobaw się z czymś takim:

function pier(x:longint):boolean;
var i:longint;
begin
if (not odd(x)) or (x<1) then pier:=false
else if x<9 then pier:=true
else begin
i:=round(sqrt(x)) or 1;
while x mod i<>0 do
dec(i,2);
pier:= i=1
end
end;

Po półtorej godziny poprzez 6, 28, 496, 8128 doszedłem do liczby 33 550 336.

Heh sposobem nr2: doszedłem do tej cyfry w 50 sekund :]
Tylko zeby nie ten zwalony zakres :/
Ten zakres to cięzki orzech do zgryzienia :
Wpisuje 9000 i mi szuka ponad 50 mln :/ - a nawet zresztą nie wiem bo mi się nie chce czekać :P

edit:

To jest nie możliwe po następnej kompilacji - na tym samym kompie nie znajduje mi juz liczby 33 miliony cośtam ...
w ogóle pętla nie staje - robi to jak na pentiumie "jedynce"

Mam wersję jeszcze szybszą! ;-)

if zakres<=6 then write(6);
if zakres<=28 then writeln(28);
..... 496,
8128, 33550336, 8589869056, 137438691328, 2305843008139952128
znanych jest jeszcze 31 doskonałych

readln(zakres);
j:=2;
while potega(2,j-1)*(potega(2,j)-1)<zakres do
inc(j);
zakres:=j;
.............

Sposób drugi to już nie optymalizacja kodu wynikająca z prostego spostrzeżenia co do pierwiastka, to już zmiana metody, miej obie.
for j:=2 to round(sqrt(j)) do
if i mod j =0 then
suma:=suma + j + i div 2;

potega(2,j) liczby longint;
włącz kontrolę zakresów, "j" czyli zakres może być maksymalnie równe 31.

Ehh teraz to już kompletnie namieszane :] - Nie dość że piszę nieczytelne kody ;) to teraz na dodatek swojego wlasnego programu nie rozumiem po tym co napisales :P
Gdybyś tak jeszcze podał gdzie co podstawić <> zamienić to już by był luxus ;P
Spróbuje do tego dojść sam, ale jak możesz to napisz na priv message chociazby zmieniony kod - jakis ostatnio czym wiecej sie ucze tym mniej wiem :]

W poście z 1:13 podałem zmienioną pętlę "j" dla pierwszego podanego przez ciebie programu. Pętla nie kręci się do "i div 2" tylko do "sqrt(i)".
Można by chyba pomyśleć jeszcze tak, jeżeli 2 nie jest podzielnikiem "i" to nie będzie już innych parzystych podzielników i "j" mogło by przyrastać o 2

W poście z 1:13 podałem zmienioną pętlę "j" dla pierwszego podanego przez ciebie programu. Pętla nie kręci się do "i div 2" tylko do "sqrt(i)".
Można by chyba pomyśleć jeszcze tak, jeżeli 2 nie jest podzielnikiem "i" to nie będzie już innych parzystych podzielników i "j" mogło by przyrastać o 2

Hmm co do pętli i div 2 to myślę że właśnie do tej pętli ma się jednak kręcić
Sqrt(i) to pierwiastek tak ? W każdym razie to do pierwszego przykładu.
A co do tego żeby przyrastało co dwa to jest dobry pomysł myślę że było by szybciej conajmniej o 1/3 :]

Dzięki za wszelką pomoc :)

Hmm co do pętli i div 2 to myślę że właśnie do tej pętli ma się jednak kręcić
Sqrt(i) to pierwiastek tak ? W każdym razie to do pierwszego przykładu.
A co do tego żeby przyrastało co dwa to jest dobry pomysł myślę że było by szybciej conajmniej o 1/3 :]
</quote>
Właśnie, że NIE!!!!
sqrt to pierwiastek, połowa dzielników jest <=sqrt, a druga połowa >=sqrt.

Do drugiego rozwiązania:
Mam pomysł jak pomniejszyć tablicę 32 razy! Do ?tab[x]? wpisujesz ?x?, później ewentualnie zerujesz, Czyli albo tab[x]==x, albo tab[x]==0, 4 razy pomniejszysz ją deklarując jako boolewską, a jeszcze 8 razy jeżeli będzie bajtowa i wykorzystasz bity jako flagi logiczne. Inicjujemy wpisując 255 ($FF) wszystkie prawdziwe, zerowanie może wyglądać:
Bajt:=x div 8;
Bit:=x mod 8;
Tab[bajt]:=tab[bajt] and (not (1 shl bit));

Testowanie
Bajt:=x div 8;
Bit:=x mod 8;
If tab[bajt] and (1 shl bit)) <>0 then ......

Co do sqrt, sprawdź jakie dzielniki znajduje się szukając do sqrt(n) (każdemu odpowiada drugi) i ile to trwa w porównaniu z szukaniem aż do połowy, wypróbuj np. na liczbie 2^{31-1=2147483647. Warto wiedzieć, że 2}m == 1 shl m (tu ?^? to potęgowanie, jedynka przesunięta m razy w lewo). Spora różnica kiedy pętla musi się kręcić miliard razy zamiast tylko 46 tysięcy!
Szukając do pierwiastka, powstaje pewien kłopot, każdy znaleziony podzielnik to ?suma:=suma+ (j) + (i div j)?, a co jeżeli sprawdzana liczba jest kwadratem? Wtedy j==i div j, czyli jeden podzielnik, a nie dwa, zrobiłem to tak:
?p:=trunc(sqrt(i)); if i=pp then suma:=p+1 else suma:=1?, a pętla szukająca podzielników ?while j<p do?, ważne ?<? nie ?<=?.
Osobne traktowanie liczb nieparzystych poprawia trochę szybkość.
Jeszcze jedna rzecz, w pętli można sprawdzać ?if suma>i then break?, to zrywa pętlę wcześniej bo już wiadomo że liczba doskonała nie jest.
Co do wykorzystania twierdzenia 2^(k-1)(2</sup>k-1), to pamiętaj że to już nie optymalizacja kodu tylko wykorzystanie teorii liczb, a co jeśli zostaniesz poproszony o dowód?
No i zaprowadziło cię ono do dziwnego programu. Zrobiłem to inaczej.
Zmienna ?k? to nasz wykładnik, czyli
Begin
K:=1;
Repeat
.....
inc(k)
until k=32 {zakres longint nas ogranicza}
end.
W tej pętli musimy sprawdzić czy 2^m-1 jest liczbą pierwszą, np.: ?L:=$7FFFFFFF shr (31-k)? (po znaku $ jest liczba szesnastkowa, ?7? i siedem razy ?F?), można niby prościej ?L:=1 shl k-1? ale dla k==31 powstaje błąd.
No to piszemy while j<....., po tej pętli wiemy czy L jest pierwsze, jeżeli tak to mamy liczbę doskonałą równą ?L*2^(k-1)?. 16 linijek kodu.

Znalazłem fajnego unita, liczy na liczbach całkowitych do około 2000 cyfr (liczba ma 6400 bitów).
"Drugim sposobem" znalazłem 10 liczb doskonałych w dwie minuty osiemnaście sekund. Dziesiąta ma 37 cyfr. Po nie całej minucie program zaczął sprawdzać kolejną (2^89-1) ale wynik będzie w przyszłym tygodniu ;-)

Poszperałem w "starociach" czyli klasykach, znalazłem algorytm sprawdzania "pierwszości" opaty tylko na mnożeniu i dodawaniu, ale kogo to dziś obchodzi/