Zaawansowane algorytmy i struktury danych/Wykład 4

Abstrakt

Naturalna metoda dodawania dwóch wielomianów wymaga czasu $\mathrm{\Theta }(n)$, natomiast prosty algorytm mnożenia dwóch wielomianów stopnia $n$ wymaga czasu $\mathrm{\Theta }(n^2)$. W wykładzie tym pokażemy, jak z wykorzystaniem szybkiej transformaty Fouriera (STF), wykonać wszystkie podstawowe operacje na wielomianach w czasie większym niż $\mathrm{\Theta }(n)$ o czynnik polilogarytmiczny. Pokażemy jak dla wielomianów stopnia $n$:

• mnożyć je w czasie $O(n\log n)$,
• obliczać wielomian interpolacyjny w czasie $O(n{\log }^{2}n)$,
• obliczać wartość wielomianu w $n$ punktach w czasie $O(n{\log }^{2}n)$,
• dzielić wielomiany w czasie $O(n{\log }^{3}n)$.

Mnożenie wielomianów w punktach

Niech $A(x)=\sum _{i=0}^{n-1}{a}_i{x}^i$ i $B(x)=\sum _{i=0}^{n-1}{b}_i{x}^i$ będzą wielomianami stopnia $n$ nad ciałem $F$. Wielomiany te możemy jednoznaczne reprezentować poprzez ich wartości w $n$ punktach. Następujące twierdzenie zostało sformułowane w ramach wykładu z Metod Numerycznych.

Twierdzenie [Twierdzenie o interpolacji wielomianów]

Dla dowolnego zbioru $n$ par $X=\{(x_0,y_0),(x_1,y_1),\dots ,(x_{n-1},y_{n-1})\}$ takiego, że wszystkie wartości $x_i$ są parami różne, istnieje jedyny wielomian $C(x)$ stopnia $n$ taki, że $C(x_i)=y_i$ dla $i=0,1,\dots ,n-1.$

Niech $X$ będzie ustalonym zbiorem parami różnych punktów $x_0,\dots ,x_{2n-1}\in F$. Dla tego zbioru punktów możemy wyznaczyć zbiory wartości wielomianów:

$X_A=\{(x_0,A(x_0)),(x_1,A(x_1),\dots ,(x_{2n-1},A(x_{2n-1}))\}$

$X_B=\{(x_0,B(x_0)),(x_1,B(x_1),\dots ,(x_{2n-1},B(x_{2n-1}))\}$

Niech $C$ będzie wynikiem mnożenia wielomianów $A$ i $B$, mamy wtedy

$C(x_i)=A(x_i)\cdot B(x_i)$.

Ponieważ stopień wielomianu $C$ jest nie większy niż $2n$ to z Twierdzenia o interpolacji zbiór wartości

$X_{A\times B}=\{(x_0,A(x_0)B(x_0)),(x_1,A(x_1)B(x_1),\dots ,(2x_{n-1},A(x_{2n-1})B(x_{2n-1}))\}$,

jednoznacznie wyznacza wielomian $A\times B$. Mając zbiory $X_A$ i $X_B$ możemy wyznaczyć zbiór $X_C$ w czasie $O(n)$. Procedura ta jest przedstawiona na następującym rysunku:

<flash>file=Zasd_fft1.swf|width=460|height=350</flash>

Jednak aby ostatecznie otrzymać szybszy algorytm niż algorytm naiwny musimy pokazać jak rozwiązać problem obliczania wartości wielomianu w $n$ punktach w czasie szybszym niż $\mathrm{\Theta }(n^2)$. Podobnie musimy umieć obliczać wielomian interpolacyjny dla danego zbioru punktów.

Szybka transformata Fouriera (STF)

Problem obliczania wartości wielomianu w $n$ punktach i problem jego interpolacji rozwiążemy wykorzystując szybką transformatę Fouriera. W poprzednim rozdziale nie zakładaliśmy nic na temat zbioru punktów $X$. Głównym pomysłem w konstrukcji algorytmu STF będzie właśnie wybór odpowiedniego zbioru punktów X tak, aby jak największa ilość wykonywanych obliczeń powtarzała się.

Założymy chcemy obliczyć wartości wielomianu $A(x)=\sum _{i=0}^{n-1}{a}_i{x}^i$ oraz $n$ jest parzyste. Jeżeli $n$ jest nieparzyste to dodajemy na początek $A(x)$ jednomian $0x^{n+1}$ co nie zmienia nam wyniku działania algorytmu. Punkty $X_n=\{x_0,x_1,\dots ,x_{n-1}\}$ zdefiniujemy w następujący sposób:

$x_i=e^{\frac{2\pi i}{n}}$.

Dla wielomianu $A(x)$ definiujemy dwa nowe wielomiany $A^{[0]}(x)$ i $A^{[1]}(x)$ poprzez wybranie do nich współczynników $A(x)$ o numerach odpowiednio parzystych i nieparzystych:

$A^{[0]}(x)=a_0+a_{2}x+a_4{x}^2+\dots +a_{n-2}{x}^{\frac{n}{2}-1}$,

$A^{[1]}(x)=a_1+a_{3}x+a_5{x}^2+\dots +a_{n-1}{x}^{\frac{n}{2}-1}$.

Wielomiany $A^{[0]}(x)$ oraz $A^{[1]}(x)$ są stopnia co najwyżej $\frac{n}{2}$. Co więcej zachodzi:

$A(x)=A^{[0]}(x^2)+x{A}^{[1]}(x^2)$      (1)

Widzimy teraz, że problem ewaluacji wielomianu $A(x)$ w punktach ${\omega }_n^0,{\omega }_n^1,\dots ,{\omega }_n^{n-1}$ sprowadza się do:

• ewaluacji wielomianów $A^{[0]}(x)$ i $A^{[1]}(x)$

w punktach

$X^{\prime }=\{x_0^2,x_1^2,\dots ,x_{n-1}^2\}$.

• a następnie obliczenie wartości $A(x)$wyniku zgodnie ze wzorem (1).

Zauważmy, że z definicji punktów $x_i$ mamy:

$x_i^2={\left(e^{\frac{2\pi i}{n}}\right)}^2=e^{\frac{2\pi i}{n/2}}$.

Możemy teraz zauważyć, że zachodzi $x_i^2=x_{i+\frac{n}{2}}^2$, a więc $X^{\prime }=X_{\frac{n}{2}}$. Udało nam się więc zredukować problem rozmiaru $n$ - obliczenia wartości wielomianu $A(x)$ stopnia $n$ w $n$do punktach, do dwóch problemów rozmiaru $\frac{n}{2}$ - obliczenia wartości wielomianów $A^{[0]}(x)$ i $A^{[1]}(x)$ stopnia $\frac{n}{2}$ w $\frac{n}{2}$ punktach. Możemy teraz zastosować tą technikę rekurencyjne otrzymując następujący algorytm.

Algorytm Algorytm Szybkiej Transformaty Fouriera

---

 STF($a=(a_0,\dots ,a_{n-1})$)
 if $n$ nieparzyste then
   dodaj wyraz $a_n$ do $a$
   zwiększ $n$
 if $n=1$ then return a
 ${\omega }_n=e^{\frac{2\pi i}{n}}$
 $\omega =1$
 $a^{[0]}=(a_0,a_2,\dots ,a_{n-2})$
 $a^{[1]}=(a_1,a_3,\dots ,a_{n-1})$
 $y^{[0]}=SFT(a^{[0]})$
 $y^{[1]}=SFT(a^{[1]})$
 for k=0 to $\frac{n}{2}-1$ do
   $y_k=y_k^{[0]}+\omega y_k^{[1]}$
   $y_{k+\frac{n}{2}}=y_k^{[0]}-\omega y_k^{[1]}$
   \omega = \omega \omega_n
 return y

Algorytm ten najpierw oblicza SFT wielomianów $A^{[0]}(x)$ i $A^{[1]}(x)$ a następnie łączy te wyniki w celu wyliczenia SFT dla wielomianu $A(x)$. Przeanalizujmy teraz wykonanie pętli. Zauważmy najpierw, że w $k$'tym kroku pętli mamy $\omega ={\omega }_n^k=e^{\frac{2\pi ik}{n}}=x_k.$. Czyli:

$y_k=y_k^{[0]}+x_k{y}_k^{[1]}=A^{[0]}(e^{\frac{2\pi ik}{n/2}})+x_k{A}^{[1]}(e^{\frac{2\pi ik}{n/2}})=$ $=A^{[0]}\left({\left(e^{\frac{2\pi ik}{n}}\right)}^2\right)+x_k{A}^{[1]}\left({\left(e^{\frac{2\pi ik}{n}}\right)}^2\right)=A^{[0]}(x_k^2)+x_k{A}^{[1]}(x_k^2)=A(x_k),$

oraz

$y_{k+\frac{n}{2}}=y_k^{[0]}-x_k{y}_k^{[1]}=A^{[0]}(e^{\frac{2\pi ik}{n/2}})-e^{\frac{2\pi ik}{n}}{A}^{[1]}(e^{\frac{2\pi ik}{n/2}})=$ $=A^{[0]}(e^{\frac{2\pi i(k+n/2)}{n/2}})+e^{\frac{2\pi ik+n/2}{n}}{A}^{[1]}(e^{\frac{2\pi i(k+n/2)}{n/2}})=A^{[0]}(x_{k+n/2}^2)+x_{k+n/2}^2{A}^{[1]}(x_{k+n/2}^2)=A(x_{k+\frac{n}{2}}).$

Gdzie w ostatniej równości skorzystaliśmy ze wzoru (1). Widzimy zatem, że algorytm poprawnie oblicza wartość STF dla wielomianu $A(x)$. Równanie rekurencyjne na czas działania procedury STF wygląda następująco:

$T(n)=2T(\frac{n}{2})+\mathrm{\Theta }(n)=\mathrm{\Theta }(n\log n)$.

Odwrotna transformata Fouriera

Aby zakończyć konstrukcję algorytmu dla szybkiego mnożenia wielomianów pozostaje nam pokazanie jak wykonać obliczyć wielomian interpolujący dla zbioru punktów $X_n$. Obliczenie wykonane w czasie szybkiej transformaty Fouriera możemy przedstawić w postaci macierzowej jako mnożenie macierzy przez wektor $(A(x_0),A(x_1),\dots ,A(x_{n-1}){)}^T=V_n(a_0,a_1,\dots ,a_{n-1}{)}^T$, gdzie $V_n=V(x_0,\dots ,x_{n-1})$ jest macierzą Vandermonde'a zawierającą potęgi $x_j$

$\left(\begin{array}{c}A(x_0)\\ A(x_1)\\ A(x_2)\\ ⋮\\ A(x_{n-1})\end{array}\right)=\left[\begin{array}{cccccc}{x}_0^0& x_0^1& x_0^2& x_0^3& \cdots & x_0^{n-1}\\ x_1^0& x_1^1& x_1^2& x_1^3& \cdots & x_1^{n-1}\\ x_2^0& x_2^1& x_2^2& x_2^3& \cdots & x_2^{n-1}\\ ⋮& ⋮& ⋮& ⋮& \ddots & ⋮\\ x_{n-1}^0& x_{n-1}^1& x_{n-1}^2& x_{n-1}^3& \cdots & x_{n-1}^{n-1}\\ \end{array}\right]\left(\begin{array}{c}{a}_0\\ a_1\\ a_2\\ ⋮\\ a_{n-1}\end{array}\right)$

Element macierzy $V(x_0,\dots ,x_{n-1})$ dany jest jako

$(V_n{)}_{j,k}=V(x_0,\dots ,x_{n-1}{)}_{j,k}=x_j^k$.

Korzystając z definicji zbioru $X_n$ otrzymujemy

$V(x_0,\dots ,x_{n-1}{)}_{j,k}={\left(e^{\frac{2\pi ij}{n}}\right)}^k=e^{\frac{2\pi ijk}{n}}.$

W celu wykonania operacji odwrotnej do SFT, czyli obliczenia wielomianu interpolacyjnego, musimy wykonać mnożenie $V_n^{-1}(A(x_0),A(x_1),\dots ,A(x_{n-1}){)}^T$.

{{lemat||| Niech macierz $W_n$ będzie zdefiniowana jako

${\left(W_n\right)}_{j,k}=\frac{1}{n}{e}^{\frac{-2\pi ijk}{n}},$

jest macierzą odwrotną do macierzy $V_n$. } Dowód

Pokażemy, że $V_n{W}_n=I$. Rozważmy pozycję $(j,k)$ macierzy $V_n{W}_n$:

${\left(V_n{W}_n\right)}_{j,k}=\sum _{l=0}^{n-1}{\left(V_n\right)}_{j,l}{\left(W_n\right)}_{l,k}=\sum _{l=0}^{n-1}{e}^{\frac{2\pi ijl}{n}}\frac{1}{n}{e}^{\frac{-2\pi ilk}{n}}=\sum _{l=0}^{n-1}\frac{1}{n}{e}^{\frac{2\pi l(j-k)}{n}}=$

Jeżeli $j=k$ to $e^{\frac{2\pi k(j-k)}{n}}=1$ i suma ta jest równa $1$. W przeciwnym przypadku możemy skorzystać ze wzoru na sumę szeregu geometrycznego:

$=\frac{1}{n}\frac{1-e^{\frac{2\pi n(j-k)}{n}}}{1-e^{\frac{2\pi (j-k)}{n}}}=\frac{1}{n}\frac{1-1^{(j-k)}}{1-e^{\frac{2\pi (j-k)}{n}}}=0.$

Czyli rzeczywiście $V_n{W}_n=I$. \qed

Porównując postać macierzy $V_n$ oraz macierzy $W_n$ widzimy, że w celu obliczenia transformaty odwrotnej możemy użyć Algorytmu Szybkiej Transformaty Fouriera, musimy tylko zamienić linijkę ${\omega }_m=e^{\frac{2\pi i}{n}}$ na ${\omega }_m=e^{-\frac{2\pi i}{n}}$ i podzielić otrzymany wynik przez $n$.

Dzielenie wielomianów

W tej części wykładu skupimy się na problemie dzielenia dwóch wielomianów. Niech $A(x)$ będzie wielomianem stopnia $m$ a $B(x)$ wielomianem stopnie $n$. W problemie dzielenia wielomianów, chcemy obliczyć dwa wielomiany $D(x)$ i $R(x)$ takie, że

$A(x)=B(x)D(x)+R(x),$

oraz stopień wielomianu $R(x)$ jest ostro mniejszy niż $n$. Wielomian $D(x)$ nazywamy {{def|wynikiem dzielenia|}, a wielomian $R(x)$ to {{def|reszta z dzielenia|}. Pierwszy pomysł jaki się od razu nasuwa, to spróbować policzyć odwrotność wielomianu $B(x)$ i przemnożenie przez tą odwrotność stronami tego równania. Niestety wielomiany nie mają niestety odwrotności będących wielomianami. Jednak nie pozostajemy tutaj zupełnie bezradni. Możemy rozszerzyć naszą dziedzinę obliczeń tak aby zagwarantować istnienie pewnych odwrotności.

Obliczenia będziemy wykonywać nad zbiorem szeregów formalnych $F[[x]]$ nad ciałem $F$, patrz Wykład z matematyki dyskretnej. Dla części elementów $F[[x]]$ istnieją odwrotności. Elementy te są postaci $1+xA(x)$ gdzie $A(x)\in F[[x]]$.

Ćwiczenie odwrotność_formalna

{{{3}}}