6118. Funkcja σ [SIGMA_PL]

Zadanie:
https://pl.spoj.com/problems/SIGMA_PL

Skrócony opis problemu:

Dla każdego z $t$ testów składających się z liczby $n$ wypisać σ(n), czyli sumę dzielników $n$.

Rozwiązanie naiwne:

Sprawdzić wszystkie liczby z przedziału $\left<1;n\right>$ i zsumować te, które dzielą bez reszty $n$.

Ewentualnie można przyspieszyć lekko ten algorytm sprawdzając liczby z przedziału $\left<1; \frac{n}{2}\right>$, a na koniec do wyniku dodać $n$.

Złożoność: $O(t \cdot n)$, gdyż dla każdego z $t$ testów wykonujemy $\frac{n}{2}$ sprawdzeń.

Rozwiązanie sprytniejsze:

Możemy zauważyć, że jeśli $x|n$ to również $\frac{n}{x}|n$, jako że $n = x \cdot \frac{n}{x}$. Wystarczy więc sprawdzić dzielniki do $\sqrt n$ i dodawać do wyniku $x$ oraz $\frac{n}{x}$. Wyjątkiem jest $\sqrt n$, którego nie należy liczyć podwójnie.

Złożoność: $O\left(t \cdot \sqrt n\right)$, gdyż dla każdego z $t$ testów wykonujemy $\sqrt n$ sprawdzeń.

Rozwiązanie jeszcze sprytniejsze:

Możemy wykorzystać też fakt, że funkcja σ jest funkcją multiplikatywną (a więc $f(m \cdot n) = f(m) \cdot f(n)$ dla $m \perp n$, a więc jeśli $m$ i $n$ są względnie pierwsze). Drugi fakt, który wykorzystamy to wartość funkcji dla liczb pierwszych. Jeśli liczba jest pierwsza, to nie ma oczywiście żadnych dzielników oprócz 1 i siebie. Wartość funkcji dla liczby pierwszej $p$ to zatem $p+1$.

Jak wykorzystamy oba fakty? Otóż przedstawmy liczbę $x$ jako iloczyn liczb pierwszych (czyli rozbijmy ją na czynniki pierwsze). $x = p_1^{w_1} \cdot p_2^{w_2} \cdot \ldots \cdot p_n^{w_n}$. A zatem, że obliczymy wynik dla każdego $p_i^{w_i}$ i je wszystkie wymnożymy. Jak obliczyć $\sigma\left(p_i^{w_i}\right)$? Tutaj nie możemy skorzystać z multiplikatywności funkcji σ, gdyż $p$ i $p$ nie są względnie pierwsze. Możemy jednak ją bardzo prosto obliczyć jeśli ją sobie rozpiszemy i zauważymy jak wygląda. Otóż: $\sigma\left(p_i^{w_i}\right) = 1 + p_i + p_i^2 + p_i^3 + \ldots + p_i^{w_i} = \sum_{j=0}^{w_i} p_i^j$. A co nam to daje? Ano to, że przecież to jest ciąg geometryczny! Czyli możemy prosto obliczyć jego sumę. Otrzymujemy zatem: $\sigma\left(p_i^{w_i}\right) = \frac{p_i^{{w_i}+1}-1}{p_i-1}$.

Złożoność: $O\left(t \cdot \sqrt n\right)$, gdyż dla każdego testu sprawdzamy dzielniki do $\sqrt n$. Co prawda mamy pętlę while w tym forze, ale za to dzielimy $n$ przez kolejne dzielniki. Tak więc obie te rzeczy się uzupełniają - mamy sobie while'a w forze, ale za to sam while się wykonuje mniej razy. Okazuje się, że ten algorytm jest szybszy niż poprzedni. Jeśli ktoś ma dowód złożoności to zapraszam do podzielenia się w komentarzach.

Algorytm:

int main()
{
        int n, div;
        long long sigma, p;
        get(n);
        // obliczamy poniżej osobno wartość sigmy dla liczby pierwszej 2, żeby potem móc sprawdzać tylko nieparzyste liczby
        // dodatkowym przyspieszeniem tutaj są operacje na bitach (np. nie musimy dzielić przez 2, ale przesuwamy liczbę w prawo)
        // jeśli ktoś nie rozumie co dzieje się w 3 poniższych linijkach, to można je zignorować, albo wrócić do nich po przeczytaniu kolejnych
        for(p = 2; !(n&1); n >>= 1)
                p <<= 1;
        sigma = p-1;
        for(div = 3; div*div <= n; div += 2) // dla co drugiej liczby do sqrt(n) włącznie
                if(n%div == 0) // jeśli jest ona dzielnikiem
                {
                        // ogólnie tutaj i w pętli poniżej obliczamy p^(w+1), czyli licznik ze wzoru z wyjaśnienia powyżej; w to wykładnik, czyli ile razy n dzieli się przez p
                        p = div; // na początku p = div, a nie p = 1, bo mamy "w+1", a nie "w" w wykładniku
                        while(n%div == 0)
                                n /= div,
                                p *= div;
                        // tak jak pisałem - sigma(n) to iloczyn sigm z poszczególnych liczb pierwszych
                        // mnożymy więc wynik przez (p^(w+1)-1)/(div-1)
                        // nie chciało mi się korzystać z nawiasów poniżej ;-)
                        sigma *= --p/~-div;
                }
        if(n > 1) // jeśli coś nam zostało, to musi to być liczba pierwsza o wykładniku 1, bo jakby był on chociaż 2, to by nam nie wyszło z pętli
                sigma *= n+1; // sigma z liczby pierwszej p to p+1
        print(sigma);
}

Rozwiązanie najsprytniejsze:

Możemy skorzystać z tzw. preprocessingu. Na początku programu tracimy dodatkowy czas, żeby coś obliczyć, ale obliczone dane przyspieszają nam pozostałą część programu. Stosuje się go jeśli jest dużo testów. Wtedy preprocessing jest opłacalny, gdyż jest wykonywany raz, a z jego wyników korzysta się wiele razy i zmniejsza złożoność dla każdego testu. Sumaryczny zysk jest wtedy większy niż koszt poniesiony na początku.

Tak więc preprocessingiem tutaj będzie sito Eratostenesa. Brzmi może i strasznie, ale idea jest bardzo prosta. Prosiłbym o zapoznanie się z tym opisem (jest podany obszerny przykład). Po zapuszczeniu sita Eratostenesa będziemy mieli tablicę liczb pierwszych do 31607. 31607 to największa liczba pierwsza mniejsza lub równa od pierwiastka maksymalnego $n$. Mając już ową tablicę wykorzystamy poprzedni algorytm. Różnica będzie taka, że nie będziemy sprawdzać co 2 liczbę, ale będziemy poruszać się bezpośrednio po liczbach pierwszych.

Złożoność: $O\left(\left(n \log n\right)\left(\log \log n\right) + t \cdot \pi\left(\sqrt n\right)\right)$. Pierwszy składnik sumy to złożoność sita Eratostenesa, a drugi oznacza, że dla każdego testu sprawdzamy każdą liczbę pierwszą do pierwiastka z $n$. $\pi\left(x\right)$ oznacza bowiem ilość liczb pierwszych do $x$ i jest szacowana przez $\frac{x}{\ln x}$. Tak więc złożoność tego algorytmu to: $O\left(\left(n \log n\right)\left(\log \log n\right) + t \cdot \frac{\sqrt n}{\ln \sqrt n}\right)$. Jak widzimy, dzielimy ów pierwiastek przez małą liczbę, ale koszt sita jest jednorazowy i niewielki, więc ten algorytm wychodzi nieznacznie szybszy.

Algorytm:

// będzie maksymalnie 3400 liczb pierwszych w przedziale 0..31607
int primes[3409], sieve[31609], count = 1; // count=1 bo wstawiamy 2 osobno do tablicy
void Eratostenes()
{
        int i, j;
        primes[0] = 2; sieve[0] = sieve[1] = 1; // wartości początkowe; 0 i 1 nie są tutaj niezbędne, ale jakbyście kopiowali ten kod...
        for(i = 4; i < 31607; i += 2) // na początku liczby parzyste znów, żeby później iść do drugą liczbę
                sieve[i] = 1;
        for(i = 3; i*i < 31607; i += 2) // dla nieparzystych liczb
                if(!sieve[i]) // jeśli jest liczbą pierwszą
                {
                        primes[count++] = i; // dodaj ją do tablicy liczb pierwszych
                        for(j = i*i; j < 31607; j += i) // dla każdej jej wielokrotności
                                sieve[j] = 1; // oznacz ją jako liczbę złożoną; 1 to liczba złożona, bo na początku ammy tablicę 0, a szkoda czasu odwracać wszystkie liczby z 0 na 1
                }
        for(j = i+!(i%2); j <= 31607; j += 2) // zaczynając od nieparzystej liczby, na której skończyliśmy
                if(!sieve[j]) // jeśli jest ona liczbą pierwszą
                        primes[count++] = j; // dodaj ją do tablicy
}

int main()
{
        int n, i;
        long long sigma, p;

        Eratostenes();
        get(n);
        for(sigma = 1, i = 0; i < count && primes[i]*primes[i] <= n; ++i) // dopóki nie skończą nam się liczby pierwsze i do sqrt(n)
                if(n%primes[i] == 0) // analogicznie do poprzedniego algorytmu
                {
                        p = primes[i];
                        while(n%primes[i] == 0)
                                n /= primes[i],
                                p *= primes[i];
                        sigma *= --p/~-primes[i];
                }
        print(sigma *= (n>1)*n+1); // skróciłem tutaj poprzedni zapis; taki fetysz
}

Optymalizacja:
Można zastosować niewielki trick - wczytać najpierw wszystkie liczby $n$, a następnie oszacować dla największego $n$ maksymalną liczbę pierwszą, która będzie dzielić to $n$. Można np. przyjąć $\sqrt n$. Dzięki temu dla małych testów nie będziemy tworzyć maksymalnego sita.