Algorytm CORDIC – Odcinek 2

Posted on by Mateusz Kowalski

Poprzedni odcinek znajduje się pod adresem https://www.kowalskimateusz.pl/algorytm-cordic-odcinek-1/

Polecam wersje pdf tego artykułu jest lepiej sformatowana.

Algorytm CORDIC odcinek 2Pobierz wersje pdf

Następny odcinek znajduje się pod adresem https://www.kowalskimateusz.pl/algorytm-cordic-odcinek-3/

Uogólniony Algorytm CORDIC

Algorytm CORDIC możemy uogólnić dokładając parametr \mu. Parametr \mu może przyjąć
wartość \mu \in \nawiasw{-1, 0, 1}. Oczywiście gdy \mu = 1 uzyskujemy tryb rotacyjny
i wektorowy wspólnie nazwany (Circular). Po polsku to będzie coś w stylu okrągły, kolisty,
okrężny, okręgowy. Dla \mu = 0 otrzymujemy tryb liniowy. Ten tryb pozwala na obliczanie
mnożenia i dzielenia bez wykonywania mnożenia i dzielenia wprost. Gdy \mu = -1 otrzymujemy
tryb hiperboliczny. Jednak w trybie hiperbolicznym jest kilka uwag co do tego wzoru.

      \[\begin{tabular}{|c|c|c|} \hline  Tryb                        & Rotacyjny          & Wektorowy          \\ \hline  Okrężny         ($\mu = 1$) & omówiony           & omówiony           \\ \hline  Liniowy         ($\mu = 0$) & teraz omówimy      & teraz omówimy      \\ \hline  Hiperboliczny   ($\mu = -1$)& w następnej części & w następnej części \\ \hline  \end{tabular}\]
      \[\begin{cases} x_{k + 1} = x_{k} - \mu d_k y_{k}2^{-k}\\ y_{k + 1} = y_{k} + d_k x_{k}2^{-k}\\ \beta_{k + 1} = \beta_{k} - d_k \gamma_k \\ \end{cases}\]

Wyróżnić możemy zatem 6 trybów. Aby określić, o którym trybie mówimy trzeba podać dwa określenia.
Pierwsze określenie aby wskazać jakie \mu mamy na myśli. Drugie czy chodzi o tryb rotacyjny
czy tryb wektorowy.

W literaturze zwykle zamiast \beta używane jest oznaczenie z, ale ja zostanę przy \beta,
bo z kojarzy się z osią z, tym bardziej, że zmienne x i y oznaczają właśnie
współrzędne w układzie współrzędnych.
Należy zdawać sobie sprawę, że \beta oznacza kąt tylko w trybie okrężnym.

Pewien przesuw pionowy

Rozważmy punkt A = \nawias{A_x, A_y} oraz wektor pionowy \varphi, tzn. równoległy do
osi y. Niech ten wektor będzie umieszczony na prostej o równania x = 1 i niech określa
przesunięcie punkt A w sposób przedstawiony na rysunku. W wyniku tego przesunięcia
punktu A otrzymujemy obraz, czyli punkt A' = \nawias{A'_x, A'_y}.
Współrzędna x jest taka sama A'_x = A_x, natomiast A'_y = A_y + A_x\varphi. Zależność ta wynika
z trójkątów podobnych. Oczywiście tak zdefiniowany przesuw jest addytywny, tzn.
skutek wykonania kolejnych dwóch jest taki sam jak jednego będącego ich sumą.

Rendered by QuickLaTeX.com

Rozkład liczby na sumę ułamków

Rozważmy ciąg geometryczny, w przypadku skończonym lub szereg geometryczny, w przypadku
nieskończonym o ilorazie q = \frac{1}{2}, czyli
1, \frac{1}{2}, \frac{1}{4}, \frac{1}{8}, \ldots.
Ten ciąg oznaczamy przez \nawias{\gamma_k} o wyrazie ogólnym \gamma_k = 2^{-k}.
W tym przypadku, \gamma_k nie oznacza kąta.

Następnie na bazie tego ciągu rozważmy nowy ciąg \nawias{\varphi_k}, którego wyrazy
\varphi_k są wyrazami tego ciągu z tym, że niektórym zmieniono znak na przeciwny.

      \[\varphi_k = \epsilon_k \gamma_k, \quad \mbox{ gdzie } \quad \epsilon_k\in\nawiasw{-1, 1}\]

Okazuje się wówczas, że każdą liczbę z przedziału od
-2 do 2 jesteśmy w stanie uzyskać jako sumę skończoną
b = \varphi_0 + \varphi_1 + \ldots + \varphi_{n-1} lub jako granicę takiego ciągu
b = \varphi_0 + \varphi_1 + \ldots. Pozostawiam to bez dowodu. Przykładowo

      \[b = \frac{7}{5} = 1 + \frac{1}{2} - \frac{1}{4} + \frac{1}{8} + \frac{1}{16} - \frac{1}{32}                      - \ldots\]

W zasadzie ten przedział jest domknięty czyli włącznie z liczbami 2 i -2, gdyż 2
otrzymamy w granicy gdybyśmy zawsze brali z plusem każdy wyraz, jest to wówczas suma
szeregu geometrycznego. Analogicznie -2 to także suma szeregu geometrycznego,
więc każdą liczbę a\in \nawiask{-2, 2} możemy tak rozpisać.
Jednak lepiej przyjąć ten przedział bez 2 i -2, tzn. a\in\nawias{-2,2}. Zauważmy, że 2
to inaczej 1\cdot 2, a 1 ma skończone rozwinięcie, w zasadzie to jest nim tylko 1 wyraz.
Mnożenie przez 2 to przesuniecie bitów w lewo.

Dobór znaków

Ponownie trzeba zapytać jak dobierać znaki? Analogicznie jak w poprzednich trybach.
Zamiast sumować kolejne \varphi_k do momentu, aż uzyskamy wartość b.
Możemy z kroku na krok modyfikować tę wartość aż osiągnie ona 0.
Niech \beta_k oznaczymy pozostałą cześć liczby b do rozłożenia,
po k-tej iteracji. Z kroku na krok modyfikujemy \beta_k, mianowicie
\beta_{k+1} = \beta_k - \varphi_k.
Jeśli bieżące \beta_k>0, to następny składnik \gamma_k należy odjąć.
Jak \beta_k<0, to następny składnik \gamma_k należy dodać.

      \[\beta_{k+1} = \beta_k - \varphi_k \quad \to \quad  \begin{cases} \beta_{k + 1} = \beta_k - \gamma_k &  \mbox{ dla } \beta_k\geq 0\\ \beta_{k + 1} = \beta_k + \gamma_k &  \mbox{ dla } \beta_k<0\\ \end{cases} \quad \to \quad \beta_{k + 1} = \beta_k - \sgn(\beta_k)\gamma_k\]

Rzecz jasna \beta_0 = b. Załóżmy, że chcemy przemnożyć daną liczbę a przez
daną liczbę b=\frac{7}{5}, wówczas

      \[ab = a\cdot \frac{7}{5} = a\cdot 1 + a\cdot\frac{1}{2} - a\cdot\frac{1}{4} + a\cdot\frac{1}{8}                          + a\cdot\frac{1}{16} - a\cdot\frac{1}{32} - \ldots\]
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com

Mnożenie to suma liczb, które były odpowiednio poprzesuwane bitowo oraz możliwe, że
zmieniono im wartość na przeciwną. Zauważmy, że \beta_k oznacza także jaka cześć
liczby b pozostała jeszcze do przemnożenia z a.
Niech x_k oznacza wartość liczby a. W każdej iteracji x_k jest takie samo.
Ta zmienna w tym trybie jest w zasadzie zbędna, ale dla spójności z pozostałymi
trybami zostawmy ją. Niech zmienna y_k oznacza bieżąca wartość iloczynu po k-tej
iteracji.

    \[A'_y = A_y + A_x\varphi \quad \to \quad y_{k + 1} = y_{k} + a\varphi_k \quad \to \quad \]


    \[\quad \to \quad y_{k + 1} = y_{k} + \sgn(\beta_k)a\cdot 2^{-k}\]


Oczywiście y_0 = 0. Uzyskujemy zatem układ rekurencji

    \[\begin{cases}x_{k + 1} = x_{k} \\y_{k + 1} = y_{k} + \sgn(\beta_k)x_k\cdot 2^{-k} \\\beta_{k + 1} = \beta_k - \sgn(\beta_k)\gamma_k \\x_0 = a \\y_0 = 0 \\\beta_0 = b \\\end{cases}\]


Po n iteracjach wynik mnożenia jest w przybliżeniu równy y_n, a niekiedy nawet dokładnie równy

    \[a\cdot b \approx y_n\]


Dla ujemnego a wzór także da poprawne wyniki. W przypadku ujemnego b także będą poprawne.

Przykłady mnożenia

Na poniższych rysunkach zwizualizowano i rozpisano kolejne kroki pracy algorytmu CORDIC
podczas mnożenia, czyli w trybie liniowym, rotacyjnym.

  • Rysunek pierwszy to mnożenia a=\frac{3}{4} z b=\frac{5}{3}.
  • Rysunek drugi to mnożenie a=\frac{5}{4} z b=\frac{3}{2}.
  • Rysunek trzeci to mnożenie a=\frac{-3}{5} z b=\frac{-2}{3}.
  • Rysunek czwarty to mnożenie a=\frac{-5}{6} z b=\frac{9}{10}. 
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com

Zauważmy, że w drugim przypadku liczby wejściowe były tak dobrane, że już po drugiej iteracji
idealnie trafiliśmy w wynik, więc każda kolejna iteracja nic nie zmieniała.
Taka sytuacja zdarza się rzadko. W 4 przypadku dokładny wynik mnożenia to \frac{-3}{4},
ale mimo to, że liczby wejściowe nie miały w mianowniku potęgi dwójki, to algorytm w małej
liczbie kroków nie trafi dokładnie w wynik.

Algorytm dzielenia

Okazuje się, że możemy także dzielić, nie wykonując dzielenia wprost.
Wyjdźmy od otrzymanego układu rekurencji, który opisuje przesuwanie się po prostej x = a.
Zastanówmy się nad sposobem podzielenia liczby c przez a. W trybie mnożącym było tak, że liczbę a
mnożyliśmy przez b w wyniku otrzymaliśmy ab. Skutkowało to tym, że punkt A=\nawias{a, 0}
przesuwaliśmy pionowo, aż uzyskaliśmy punkt o współrzędnych A'=\nawias{a, a\cdot b}.
Należy zatem ten proces odwrócić. Podobnie jak w trybie okrężnym, w podtrybie
rotacyjnym obracaliśmy punkt A=\nawias{1, 0}, aż uzyskaliśmy obraz
A'=\nawias{\cos\varphi, \sin\varphi}, w podtrybie wektorowym odwrotnie punkt
A=\nawias{A_x, A_y} był obracany, aż uzyskaliśmy A'=\nawias{A'_x, 0}.

Mamy zatem punkt A=\nawias{a, c}, który trzeba przesunąć po linii pionowej, aż będzie
na osi x, czyli A'=\nawias{a, 0}. Każde kolejne przesunięcie jest o połowę mniejsze
od poprzedniego. Tak jak w trybie rotacyjnym każdy kolejny tangens kąta obrotu był
o połowę mniejszy od poprzedniego.

Można pomyśleć, że jeśli bieżący y_k jest większy od 0 to należy przesunąć w dół,
czyli \varphi_k<0. Jeśli bieżący y_k jest ujemny to należy przesunąć w górę \varphi_k>0.
Zauważmy jednak, iż x_k może być ujemne, wówczas \varphi_k było by źle określone.
Poprawiając będzie \varphi_k = \epsilon_k \gamma_k = -\sgn(x_k y_k)\gamma_k.

      \[\begin{cases} y_{k + 1} = y_k + \varphi_k x_k \\ y_k = c \end{cases}\]
      \[\begin{cases} x_{k + 1} = x_{k} \\ y_{k + 1} = y_{k} + x_0\varphi_k \\ \beta_{k + 1} = \beta_k -\varphi_k \\ x_0 = a \\ y_0 = 0 \\ \beta_0 = b \\ \end{cases}\]
Rendered by QuickLaTeX.com

Suma przesunięć \varphi_k to przesunięcie punktu A do punktu A'. Zatem szukane przesunięcie
A' do A jest z przeciwnym znakiem. Ten wynik to wynik dzielenia. Niech zatem \beta_k
oznacza aktualnie naliczony wynik po k iteracjach. Rzeczy jasna na początku wynosi 0, to
prowadzi nas do rekurencji.

    \[\begin{cases}\beta_{k + 1} = \beta_k - \varphi_k \\\beta_0 = 0\end{cases} \quad \to \quad \begin{cases}\beta_{k + 1} = \beta_k + \sgn(x_ky_k)\gamma_k \\\beta_0 = 0\end{cases}\]

Otrzymujemy układ rekurencji

    \[\begin{cases}x_{k + 1} = x_{k} \\y_{k + 1} = y_{k} - \sgn(x_k y_k)x_k\cdot 2^{-k} \\\beta_{k + 1} = \beta_k + \sgn(x_k y_k)\gamma_k \\x_0 = a \\y_0 = c \\\beta_0 = 0 \\\end{cases}\]

Przykłady dzielenia

Na poniższych rysunkach zwizualizowano i rozpisano kolejne kroki pracy algorytmu CORDIC
podczas dzielenia, czyli w trybie liniowym, wektorowym.

  • Rysunek pierwszy to dzielenie \frac{5}{4} przez \frac{3}{2}
  • Rysunek drugi to dzielenie \frac{3}{8} przez \frac{3}{4}.
  • Rysunek trzeci to dzielenie \frac{3}{5} przez \frac{3}{2}.
  • Rysunek czwarty to dzielenie \frac{5}{2} przez 4. 
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com

Zauważmy, iż w przypadku drugim liczby były tak dobrane, że już po drugim kroku pracy
algorytmu trafiliśmy idealnie w wynik. Kolejne iteracje nic nie zmieniały. Analogicznie
na rysunku czwartym, po czwartej nie ma zmian. Taka sytuacja zdarza się rzadko.

Co robić poza przedziałem zbieżności? – mnożenie

Algorytm CORDIC w trybie mnożenia działa poprawnie. Tylko jeśli b\in\nawias{-2,2}.
Ten algorytm można łatwo zmodyfikować aby działa dla dowolnego b.
W zasadzie to dla dowolnego, choć z góry ustalonego przedziału.
Należy zmodyfikować algorytm o dodatkowy parametr, nazwałem go r\in\ZZ.

      \[\begin{cases} x_{k + 1} = x_{k}  \\ y_{k + 1} = y_{k} + \sgn(\beta_k)x_k\cdot 2^{-k} \\ \beta_{k + 1} = \beta_k   - \sgn(\beta_k)\gamma_k \\ x_0 = a \\ y_0 = 0 \\ \beta_0 = b \\ \end{cases}\quad \to \quad \begin{cases} x_{k + 1} = x_{k}  \\ y_{k + 1} = y_{k} + \sgn(\beta_k)x_k\cdot 2^{-k+\R{r}} \\ \beta_{k + 1} = \beta_k   - \sgn(\beta_k)\gamma_k \\ x_0 = a \\ y_0 = 0 \\ \beta_0 = b \\ \end{cases}\]

Jeśli chcemy zwiększyć przedział dla b dwukrotnie, czyli do b\in\nawias{-4,4},
to należy przyjąć r = 1. Jeśli chcemy zwiększyć czterokrotnie to r = 2.
Jeśli chcemy zwiększyć ośmiokrotnie to r = 3, itd.
Jeśli nie chcemy modyfikować to r = 0. Aby zachować tę samą dokładność trzeba
będzie wykonać o r iteracji więcej. Czyli jeśli zwiększymy przedział 128 krotnie,
to trzeba wówczas wykonać o r=7 iteracji więcej, aby osiągnąć ten sam rząd
dokładności wyników co w przypadku bez rozszerzenia r=0. Można to wykorzystać np.
przy liczbach stałoprzecinkowych.

Natomiast jednak jeśli dysponujemy liczbami w postaci zmiennoprzecinkowej to nie
trzeba wprowadzać tego rozszerzenia. Wówczas przedział b\in\nawias{-2,2} nam wystarczy.
Liczba zmiennoprzecinkowa jest postaci a = Z\cdot M\cdot 2^W. Gdzie Z\in\nawiasw{-1,1}
jest znakiem liczby, M\in\left[1,2 \right) jest mantysą, natomiast W\in\ZZ
jest wykładnikiem. W praktyce jest jeszcze bias przy wykładniku, ale dla skupienia
uwagi nie komplikujmy, bo można go łatwo uwzględnić.

Mnożąc dwie liczby a_1 = Z_1\cdot M_1\cdot 2^{W_1}, a_2 = Z_2\cdot M_2\cdot 2^{W_2}
zmiennoprzecinkowe mamy a_1\cdot a_2 = Z_1\cdot Z_2\cdot M_1\cdot M_2\cdot 2^{W_1 + W_2}.
Widzimy zatem, że mnożenie liczb zmiennoprzecinkowych sprowadza się do mnożenia ich
mantys, gdzie każda z nich jest liczbą z przedziału \left[1, 2\right) oraz dodawania
wykładników. Wynik iloczynu znajduje się w przedziale \left[ 1, 4\right). Możliwe,
że wymagana będzie normalizacja. Przykładowo 3,9 = 1,8\cdot 2^1. Wówczas do sumy
wykładnika trzeba dodać jeszcze to 1. Do algorytmu nie musimy wprowadzać żadnych modyfikacji.
Nawet możemy do algorytmu podawać te mantysy ze znakiem, gdyż
Z\cdot M \in \left(-2, -1 \right] \cup \left[1, 2\right).

Pojawia się pytanie co w przypadku zera? Zero ma ściśle określona sekwencja
bitów. Ponadto wynik mnożenia przez 0 z góry jest znany. Można zatem wykryć czy któraś
z liczb nie jest 0 i wówczas dać od razu wynik bez używania algorytmu.

Co robić poza przedziałem zbieżności? – dzielenie

W przypadku dzielenia sytuacja jest podobna do mnożenia. Mamy analogiczne dwie opcje.
Jedna to rozszerzania zbioru do takiego jakich używamy liczb. Drugie to wykorzystanie
własności liczb zmiennoprzecinkowych. Są jednak pewne różnice.

Algorytm w trybie dzielenia działa poprawnie jeśli wynik jest w przedziale \nawias{-2, 2}.
Można rozszerzyć ten przedział stosując analogiczną modyfikację jak ta omówiona w trybie
mnożenia, mianowicie. Dokładając analogiczny parametr r. Podobnie jak przedział
chcemy przykładowo zwiększyć 32 razy, tzn. do \nawias{-64, 64}, to przyjmujemy r=5
i aby zachować ten sam rząd precyzji należy wykonać r dodatkowych iteracji.

      \[\begin{cases} x_{k + 1} = x_{k}  \\ y_{k + 1} = y_{k} - \sgn(x_k y_k)x_k\cdot 2^{-k} \\ \beta_{k + 1} = \beta_k + \sgn(x_k y_k)\gamma_k \\ x_0 = a \\ y_0 = c \\ \beta_0 = 0 \\ \end{cases}\quad \to \quad \begin{cases} x_{k + 1} = x_{k}  \\ y_{k + 1} = y_{k} - \sgn(x_k y_k)x_k\cdot 2^{-k+\R{r}} \\ \beta_{k + 1} = \beta_k + \sgn(x_k y_k)\gamma_k \\ x_0 = a \\ y_0 = c \\ \beta_0 = 0 \\ \end{cases}\]

W przypadku liczb zmiennoprzecinkowych nie trzeba modyfikować algorytmu. Wówczas warto
wykorzystać własności liczb zmiennoprzecinkowych. Należy rozpoznać czy
przypadkiem nie chcemy wykonać dzielenia przez zero. W tym celu sprawdzić należy czy
a\neq 0. Ponadto jeśli c=0, to bez wykonywania dzielenia od razu wiadomo jaki jest wynik.
W pozostałych przypadkach liczbę możemy przedstawić w postaci a = Z\cdot M\cdot 2^W.
Gdzie Z\in\nawiasw{-1, 1} jest znakiem liczby, M\in\left[1,2 \right) jest mantysą,
natomiast W\in\ZZ jest wykładnikiem.

Dzieląc dwie liczby a_1 = Z_1\cdot M_1\cdot 2^{W_1}, a_2 = Z_2\cdot M_2\cdot 2^{W_2}
zmiennoprzecinkowe mamy \frac{a_1}{a_2} = \frac{Z_1M_1}{Z_2M_2}\cdot 2^{W_1 - W_2}
Widzimy zatem, że dzielenie liczb zmiennoprzecinkowych sprowadza się do dzielenia ich
mantys, gdzie każda z nich jest liczbą z przedziału \left[1, 2\right) oraz odejmowania
wykładników. Wynik ilorazu znajduje się w przedziale \left[ \frac{1}{2}, 2\right).
Możliwe zatem, że wymagana będzie normalizacja. Przykładowo 0,7 = 1,4\cdot 2^{-1}.
Wówczas do różnicy wykładników trzeba odjąć jeszcze to 1. Do algorytmu nie musimy
wprowadzać żadnych modyfikacji. Nawet możemy do algorytmu podawać te mantysy z dołożonym
znakiem liczb a_1, a_2, gdyż Z\cdot M \in \left(-2, -1 \right] \cup \left[1, 2\right).
Gdy dzielna i dzielnik są w tej postaci to iloraz będzie liczbą z przedziału
\frac{a}{b}\in\left( -2, -\frac{1}{2} \right] \cup \left[ \frac{1}{2}, 2 \right),

Algorytm CORDIC – Odcinek 1

Posted on by Mateusz Kowalski

Polecam wersje pdf tego artykułu jest lepiej sformatowana.

Algorytm CORDIC Odcinek 1Pobierz wersje pdf

Postanowiłem napisać ten artykuł, gdyż materiały, z który to rozgryzałem nie były najlepsze. Często nie podawały wprost założeń jakie trzeba spełnić. Jak te ograniczenia obejść itp. W dodatku w języku polskim praktycznie nie ma materiałów na ten temat. Dlatego bardzo było by mi miło gdybyś napisał mi, czy i jak Ci ten materiał pomógł. Teraz do rzeczy.

Algorytm CORDIC pozwala na obliczenie wartości kilku funkcji  matematycznych bez wykonywania mnożenia

  • sinus f(x)=\sin x
  • kosinus f(x)=\cos x
  • arkus tangens f(x) = \arctg x
  • sinus hiberboliczny f(x) = \sinh x
  • kosinus hiberboliczny f(x) = \cosh x
  • area tangens hiberboliczny f(x) = \artgh x

Pozwala także na obliczenie mnożenia i dzielenia bez wykonywania go wprost. Istnieją też różne modyfikacje pozwalające ponoć liczyć f(x)= \arcsin x oraz f(x) = \arccos x, f(x)=\ln x i inne funkcje. Tego jednak nie sprawdzałem.

Algorytm CORDIC jest pewną alternatywą dla szeregów Taylora. Opiera się tylko o dodawanie, odejmowanie, przesuwanie bitowe i potrzebuje pamięci aby spamiętać kilkanaście liczb. W tym miejscu już widać, że pierwszym wymaganiem jest to aby liczby były w postaci bitowej. Co w praktyce w zasadzie zawsze jest spełnione. Nie ma znaczenia czy operujemy na liczbach stałoprzecinkowych czy zmiennoprzecinkowych. Są jednak pewne ograniczenia co do zbieżności.

      \[\begin{tabular}{|c|c|c|} \hline Nazwa funkcji       & Przedział teoretyczny zbieżności                       & Sugerowane użycie \\ \hline $\sin x, \cos x$    & $x\in\nawias{-99,7\st; 99,7\st}$                       & $x\in\nawias{-\frac{\pi}{2}, \frac{\pi}{2}}$ \\ \hline $\arctg x$          & $x\in \RR$                                             & $x\in \RR$ \\ \hline $a \cdot b$         & $a\in\RR, b\in\nawiask{-2,2}$                          & $a,b\in\nawias{-2, 2}$ \\ \hline $\frac{a}{b}$       & $a\in\RR, b\in \RR \wedge \frac{a}{b}\in\nawiask{-2,2}$& $a,b\in\left(-2,-1\right] \cup \left[1, 2\right)$ \\ \hline  $\sinh x, \cosh x$  & $x\in \nawias{-1,11; 1,11}$                            & $x\in\nawias{-1,11; 1,11}$ \\ \hline $\artgh x $         & $x\in \nawias{-1; 1}$                                  & $x\in\nawias{-1; 1}$ \\ \hline \end{tabular}\]

Możliwości obliczeniowe po połączeniu poszczególnych wariantów algorytmu, skorzystaniu z tożsamości matematycznych lub innych modyfikacjach.

      \[\begin{tabular}{|c|c|c|} \hline Nazwa funkcji       & Przedział rozszerzony         & W jaki sposób \\ \hline $\sin x, \cos x$    & $x\in \RR$                    & \begin{tabular}{@{}c@{}}trygonometryczne wzory \\                                                                                redukcyjne\end{tabular} \\ \hline $a \cdot b$         & $a,b\in\RR$                   & \begin{tabular}{@{}c@{}}wykorzystanie własności \\                                                                               liczby zmiennoprzecinkowej \end{tabular} \\ \hline $\frac{a}{b}$       & $a,b\in\RR \wedge b \neq 0$   & \begin{tabular}{@{}c@{}}wykorzystanie własności \\                                                                               liczby zmiennoprzecinkowej\end{tabular} \\ \hline  $\sinh x, \cosh x$  & $x\in \RR$                    & \begin{tabular}{@{}c@{}}rekurencyjne wykorzystanie \\                                                                               tożsamości i mnożenia  \\                                                                               $\sinh 2x = 2\sinh x \cosh x$ \\                                                                               $\cosh 2x = \cosh^2 x + \sinh^2 x$ \end{tabular}\\ \hline \end{tabular}\]

Dla funkcji sinus i kosinus szczególne kąty typu: 0\st, 90\st, 180\st, 270\st
wykryć i potraktować oddzielnie.

Algorytm oblicza funkcję sinus i kosinus jednocześnie. Jest zatem szczególnie przydatny podczas obliczań obrotów na płaszczyźnie lub w przestrzeni trójwymiarowej. Nie jestem pewny, ale prawdopodobnie w kartach graficznych może być stosowany.Ponadto może być bardzo użyteczny w sprzęcie, w którym nie mamy do dyspozycji układów
mnożących. Grupą szczególnie zainteresowanych będą programiści niskopoziomowi, elektronicy i mikroelektronicy.

Co ciekawe Algorytm CORDIC pozwala także na obliczanie mnożenia i dzielenia bez wykonywania mnożenia i dzielenia. Tu także są pewne ograniczenia, ale poprzez odpowiednie użycie można je obejść. Oczywiście dzięki temu możemy w sytuacjach koniecznych wykonać mnożenie nie wykonując mnożenia wprost.

W przypadku funkcji hiperbolicznych \sinh x oraz \cosh x przedział zbieżności jest ograniczony tylko do około \nawias{-1,11; 1,11}. Można podać skąd te wartości pochodzą, ale jest to skomplikowane i nieistotne, bo i tak będę proponował zmniejszenie tego przedziału do \nawiask{-1, 1}. Wykorzystując tożsamości matematyczne można to rozszerzyć na większy
przedział, potrzebne jednak będzie wykonywanie mnożenia, które może być realizowane inną wersją/instancją algorytmu, tą wykonującą mnożenie bez mnożenia.

Algorytm CORDIC jest algorytmem iteracyjnym. Cechuje się tym, że kolejna iteracja niekoniecznie musi dać dokładniejszy wyniki, niż poprzednia. Co więcej na poprawę “rekordu” możemy krótko lub długo czekać. Jednak w długiej perspektywie algorytm jest zbieżny. Problem jednak polega na tym, że aby uniknąć mnożenia trzeba z góry założyć ile tych iteracji wykonamy.

Istnieje możliwość obliczenia dodatkowych funkcji matematycznych pośrednio.

      \[\begin{tabular}{|c|c|c|} \hline Funkcja matematyczna                                             & Uwagi \\ \hline $\tg x = \frac{\sin x}{\cos x}$                                  & trzeba dzielić \\ \hline $\ctg x = \frac{\cos x}{\sin x}$                                 & trzeba dzielić \\ \hline $\arcctg x = \frac{\pi}{2} - \arctg x$                           & tylko odejmowanie  \\ \hline $\arctg_2(y,x) =                                                        \begin{cases}                                                     \arctg \frac{y}{x}         & \mbox{, } x > 0                \\    \arctg \frac{y}{x} + \pi   & \mbox{, } x < 0 \wedge y\geq 0 \\    \arctg \frac{y}{x} - \pi   & \mbox{, } x < 0 \wedge y < 0   \\    \frac{\pi}{2}              & \mbox{, } x = 0 \wedge y > 0   \\    -\frac{\pi}{2}             & \mbox{, } x = 0 \wedge y < 0   \\    \mbox{nieokreślone}        & \mbox{, } x = 0 \wedge y = 0   \\    \end{cases}$                                                     & dodawanie lub odejmowanie $\pi$ \\ \hline $\tgh x = \frac{\sinh x}{\cosh x} $                              & trzeba dzielić \\ \hline $\ctgh x = \frac{\cosh x}{\sinh x} $                             & trzeba dzielić \\ \hline $e^x = \sinh x + \cosh x $                                       & \begin{tabular}{@{}c@{}} tylko dodawanie, o ile mnożenia \\                                                                                             nie było w samym $\sinh x$ lub $\cosh x$                                                                     \end{tabular} \\ \hline $\sqrt{x^2 + y^2}$                                               & trzeba mnożyć  \\ \hline  $\sqrt{x^2 - y^2}$                                               & trzeba mnożyć  \\ \hline $\ln x = 2\artgh \frac{x-1}{x+1}$                                & \begin{tabular}{@{}c@{}} trzeba dzielić, zauważmy, że \\                                                                                             dla $x>0 \Leftrightarrow \frac{x-1}{x+1}\in\nawias{-1,1}$                                                                     \end{tabular} \\ \hline $\log_a x = \frac{\ln x}{\ln a}$                                 & trzeba dzielić \\ \hline $a^x = e^{x\ln a} $                                              & trzeba mnożyć  \\ \hline $\sqrt{x} = \sqrt{\nawias{x+\frac{1}{4}}^2 - \nawias{x-\frac{1}{4}}^2}$     & trzeba mnożyć  \\ \hline \end{tabular}\]

Skupmy się teraz tylko nad zasadą działania w trybie rotacyjnym. Pozwalającym obliczać sinus i kosinus. Tę zasadę zmodyfikujemy nieco w trybie wektorowym, która pozwoli nam obliczyć \arctg x

Algorytm CORDIC w trybie rotacyjnym (rotation mode)

Obrót na płaszczyźnie o kąt \varphi przy nieruchomym układzie współrzędnych

Przypomnijmy, że obrót na płaszczyźnie można zrealizować przy pomocy macierzy obrotu.

      \[\begin{bmatrix}     \cos \varphi & -\sin \varphi \\     \sin \varphi &  \cos \varphi \\     \end{bmatrix}\]

Chcąc obrócić punkt A = \nawias{A_x, A_y} o kąt \varphi\in \RR w nieruchomym układzie współrzędnych, należy

      \[\begin{bmatrix}     A'_x \\     A'_y      \end{bmatrix} = \begin{bmatrix}     \cos \varphi & -\sin \varphi \\     \sin \varphi &  \cos \varphi \\     \end{bmatrix} \cdot \begin{bmatrix}     A_x \\     A_y      \end{bmatrix}\]
Rendered by QuickLaTeX.com

Zauważmy, że jeśli obracalibyśmy punkt o współrzędnych A = \nawias{1, 0}, to w wyniku obrotu otrzymamy punkt o współrzędnych A' = \nawias{\cos \varphi, \sin \varphi}.

    
      \[\begin{bmatrix}     \cos \varphi \\     \sin \varphi       \end{bmatrix} = \begin{bmatrix}     \cos \varphi & -\sin \varphi \\     \sin \varphi &  \cos \varphi \\     \end{bmatrix} \cdot \begin{bmatrix}     1 \\     0      \end{bmatrix}\]
Rendered by QuickLaTeX.com

Wykonując kolejno dwa obroty, najpierw o kąt \varphi, a potem o kąt \psi skutek będzie taki sam jak wykonanie obrotu od razu o kąt \varphi + \psi. Potwierdza to też rachunek macierzowy.

    
      \[\begin{bmatrix}     \cos \psi & -\sin \psi \\     \sin \psi &  \cos \psi \\     \end{bmatrix} \cdot \begin{bmatrix}     \cos \varphi & -\sin \varphi \\     \sin \varphi &  \cos \varphi \\     \end{bmatrix} = \begin{bmatrix}     \cos \nawias{\varphi + \psi} & -\sin \nawias{\varphi + \psi} \\     \sin \nawias{\varphi + \psi} &  \cos \nawias{\varphi + \psi} \\     \end{bmatrix}\]
Rendered by QuickLaTeX.com
Przekształcenie

Wróćmy do dowolnego punktu A = \nawias{A_x, A_y}. Efekt po obrocie zapiszmy w postaci układu równań.

      \[\begin{cases} A'_x  = A_x \cos \varphi - A_y \sin \varphi \\ A'_y  = A_x \sin \varphi + A_y \cos \varphi \\ \end{cases}\]

Jeżeli przyjmiemy, że wykonujemy obroty w ramach kątów
\varphi \in \nawias{-\frac{\pi}{2}, \frac{\pi}{2}}, to wówczas \cos \varphi \neq 0. W takim razie można w prawych stronach wyciągnąć przed nawias \cos \varphi, a właściwie to za nawias.

      \[\begin{cases} A'_x  = A_x \cos \varphi - A_y \sin \varphi \\ A'_y  = A_x \sin \varphi + A_y \cos \varphi \\ \end{cases} \dalej \begin{cases} A'_x  = \nawias{A_x - A_y\tg \varphi} \cos \varphi \\ A'_y  = \nawias{A_y + A_x\tg \varphi} \cos \varphi \\ \end{cases}\]

Zauważmy, że skoro \varphi \in \nawias{-\frac{\pi}{2}, \frac{\pi}{2}}, to zachodzi tożsamość

      \[\frac{1}{\sqrt{1 + \tg^2 \varphi}} = \frac{1}{\sqrt{ \frac{\cos^2 \varphi + \sin^2 \varphi}{\cos^2 \varphi}}} = \frac{1}{\frac{1}{|\cos \varphi|}} =  |\cos \varphi| = \cos \varphi\]


      \[\begin{cases} A'_x  = \nawias{A_x - A_y\tg \varphi} \cos \varphi \\ A'_y  = \nawias{A_y + A_x\tg \varphi} \cos \varphi \\ \end{cases} \dalej \begin{cases} A'_x  = \nawias{A_x - A_y\tg \varphi} \frac{1}{\sqrt{1 + \tg^2 \varphi}}  \\ A'_y  = \nawias{A_y + A_x\tg \varphi} \frac{1}{\sqrt{1 + \tg^2 \varphi}}  \\ \end{cases}\]
Obrót jako suma mniejszych obrotów

Jeżeli kąt \varphi zapiszemy jako sumę n innych kątów, tzn.

      \[\varphi = \varphi_0 + \varphi_1 + \varphi_2 + \ldots + \varphi_{n-1}\]

To wówczas można zapisać następujący układ rekurencji, którego każda iteracja odpowiada za obrót punktu (\tilde{x}_k, \tilde{y}_k) o kąt \varphi_k. %, czyli innymi słowy równania różnicowe

      \[\begin{cases} \tilde{x}_{k + 1} = \nawias{\tilde{x}_k - \tilde{y}_k\tg \varphi_k} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} \\  \tilde{y}_{k + 1} = \nawias{\tilde{y}_k + \tilde{x}_k\tg \varphi_k} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} \\  \tilde{x}_0 = A_x \\ \tilde{y}_0 = A_y \\ \end{cases}\]

Po n krokach otrzymujemy dobre przybliżenie obrazu punktu A czyli punkt A' \approx \nawias{\tilde{x}_n, \tilde{y}_n}. Zauważmy, że funkcja \arctg jest funkcją nieparzystą, tzn. \arctg \nawias{-x} = -\arctg x. Zatem jeśli znamy \arctg dla pewnego x to znamy także \arctg dla -x. Wystarczy zmienić znak \arctg na przeciwny. Zastanówmy się nad zbiorem takich kątów \gamma_k, które dodając lub odejmując mogłyby
zbliżyć się dowolnie blisko zadanemu kątowi \varphi przy odpowiednio dużej liczbie kroków.

Widzimy, że w każdej iteracji trzeba wykonać mnożenie \tilde{y}_k z \tg \varphi_k. Chcemy uniknąć mnożenia. Jedyne na co się godzimy to mnożenie przez potęgę 2 tzn. przez 2^m, gdzie m\in \ZZ. Mnożenie liczby a w systemie dwójkowym przez 2^m odpowiada przesuwaniu bitów liczby a w lewo lub w prawo. Jeśli m jest dodatnie to wówczas przesuwamy bity w lewo o m pozycji. Jeśli m jest ujemne wówczas przesuwamy bity w prawo o -m pozycji. Niech zatem to będzie taki zbiór kątów, dla których \tg \gamma_k będzie potęgą 2. Rozpatrzmy zbiór kątów dla, których \tg \gamma_k = 2^{-k}, gdzie k jest liczba naturalną. Matematycznie taki zbiór można zapisać tak:

      \[\nawiasw{\gamma_k \in \left( 0, \frac{\pi}{4} \right] :           \gamma_k = \arctg 2^{-k} \quad \wedge \quad k\in \NN }\]

Gdybyśmy z góry wcześniej wyznaczyli wartości tych kątów do pewnego niekoniecznie dużego n, np. n=40 i je zapamiętali, to wówczas zamiast wykonywać kosztowne mnożenie przesuwalibyśmy bity liczby.

Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com
      \[\begin{tabular}{|c|c|c|c|c|c|c|c|c|c|c|} \hline  $k$         & 0        & 1                  & 2                     & 3                     & 4                     & 5                     & 6                     & 7                     & 8                     & \ldots \\ \hline  $\gamma_k$ &$\arctg 1$ &$\arctg\frac{1}{2}$ &$\arctg\frac{1}{4}$    &$\arctg\frac{1}{8}$    &$\arctg\frac{1}{16}$   &$\arctg\frac{1}{32}$   &$\arctg\frac{1}{64}$   &$\arctg\frac{1}{128}$  &$\arctg\frac{1}{256}$  & \ldots \\ \hline  $\gamma_k$ &$45,0\st$  &$26,57\st$          &$14,04\st$             &$7,13\st$              &$3,58\st$              &$1,79\st$              &$0,9\st$               &$0,45\st$              &$0,22\st$              & \ldots \\ \hline  $\gamma_k$ & 0,785     & 0,464              & 0,245                 & 0,124                 & 0,062                 & 0,031                 & 0,016                 & 0,008                 & 0,004                 & \ldots \\ \hline  \end{tabular}\]
      \[\varphi_k = \epsilon_k\gamma_k, \quad \mbox{gdzie} \quad  \epsilon_k \in \nawiasw{-1, 1}\]

Należy się zastanowić czy każdy kąt z przedziału \varphi \in \nawias{-\frac{\pi}{2}, \frac{\pi}{2}} można przestawić jako sumę tych kątów, gdzie dany kąt może być wzięty z plusem lub z minusem. Przykładowo kąt

      \[74\st \approx \gamma_0 + \gamma_1 + \gamma_2 - \gamma_3 - \gamma_4 - \gamma_5 + \gamma_6\]


      \[74\st \approx 45,0\st + 26,57\st + 14,04\st - 7,13\st - 3,58\st - 1,79\st + 0,9\st\]

Dla każdego kąta można znaleźć taki ciąg plusów i minusów jaki należy dołożyć do kątów
\gamma_k aby ten ciąg dążył do zadanego kąta. Pozostawiam to bez dowodu.

Warto zaznaczyć, że w długiej perspektywie iteracji wyniki się poprawiają, lecz w krótkim zakresie parę kolejnych iteracji może dać gorszy wynik zanim trafimy na iterację, która wynik poprawi.

Rendered by QuickLaTeX.com

Pojawia się pytanie jak wybierać te znaki? Odpowiedź jest prostsza, niż się wydaje. Jak aktualnie naliczona suma kątów \varphi_0 + \varphi_1 + \ldots + \varphi_k jest mniejsza od szukanego kąta \varphi to następny kąt \gamma_{k + 1} trzeba dodać. Jak bieżąca suma jest większa od danego kąta \varphi to następny kąt \gamma_{k + 1} trzeba odjąć.

Tę zasadę można odwrócić i dany kąt \varphi z kroku na krok modyfikować, aż będzie dostatecznie blisko 0 lub gdy uznamy, że osiągnęliśmy maksymalną liczbę kroków n. Niech \beta_k oznacza kąt, po kroku k, jaki pozostał jeszcze do obrócenia, aby znaleźć się w punkcie A', będącym obrazem punkt A.

      \[\begin{cases} %\beta_{k + 1} = \beta_{k} - \underbrace{\sgn(\beta_k)}_{\epsilon_k} \gamma_k \\ \beta_{k + 1} = \beta_{k} - \varphi_k \\ \beta_0 = \varphi  \end{cases}\]

Gdzie po n krokach \beta_n \approx 0. Zauważmy także, że \varphi_k = \epsilon_k \gamma_k = \sgn(\beta_k) \gamma_k, a co za tym idzie

      \[\tg \varphi_k = \tg (\epsilon_k \gamma_k) = \epsilon_k\tg \gamma_k = \sgn(\beta_k)\tg \gamma_k = \sgn(\beta_k)2^{-k}\]
Rendered by QuickLaTeX.com

W ten sposób problem mnożenia \tilde{y}_k \cdot \tg \varphi_k oraz
\tilde{x}_k \cdot \tg \varphi_k został rozwiązany. Pozostała kwestia
mnożenia przez \frac{1}{\sqrt{1 + \tg^2 \varphi_k}}.

Mnożenie przez \frac{1}{\sqrt{1 + \tg^2 \varphi_k}}

Zaobserwujmy co się dzieje podczas dwóch kolejnych iteracji.

      \[\begin{cases} \tilde{x}_{1} = \nawias{\tilde{x}_0 - \tilde{y}_0\tg \varphi_0} \frac{1}{\sqrt{1 + \tg^2 \varphi_0}} \\  \tilde{y}_{1} = \nawias{\tilde{y}_0 + \tilde{x}_0\tg \varphi_0} \frac{1}{\sqrt{1 + \tg^2 \varphi_0}} \\  \end{cases} \dalej \begin{cases} \tilde{x}_{2} = \nawias{\tilde{x}_1 - \tilde{y}_1\tg \varphi_1} \frac{1}{\sqrt{1 + \tg^2 \varphi_1}} \\  \tilde{y}_{2} = \nawias{\tilde{y}_1 + \tilde{x}_1\tg \varphi_1} \frac{1}{\sqrt{1 + \tg^2 \varphi_1}} \\  \end{cases}\]


      \[\begin{cases} \tilde{x}_{2} = \nawiask{\nawias{\tilde{x}_0 - \tilde{y}_0\tg \varphi_0} - \nawias{\tilde{y}_0 + \tilde{x}_0\tg \varphi_0}\tg \varphi_1} \frac{1}{\sqrt{1 + \tg^2 \varphi_0}\sqrt{1 + \tg^2 \varphi_1}} \\  \tilde{y}_{2} = \nawiask{\nawias{\tilde{y}_0 + \tilde{x}_0\tg \varphi_0} + \nawias{\tilde{x}_0 - \tilde{y}_0\tg \varphi_0}\tg \varphi_1} \frac{1}{\sqrt{1 + \tg^2 \varphi_0}\sqrt{1 + \tg^2 \varphi_1}} \\  \end{cases}\]

Możemy zatem skupić się na rekurencji poniżej, w której po wykonaniu n iteracji otrzymamy przybliżenie obraz punkt A, czyli punkt A'.

      \[\begin{cases} x_{k + 1} = x_{k} - y_{k}\tg \varphi_{k} \\ y_{k + 1} = y_{k} + x_{k}\tg \varphi_{k} \\ x_{0} = A_x \\ y_{0} = A_y \\ \end{cases} \quad \to \quad  \begin{cases} %A'_x \approx  \tilde{x}_{n} = x_{n} \prod\limits_{k = 0}^{n - 1} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} \\  %A'_y \approx  \tilde{y}_{n} = y_{n} \prod\limits_{k = 0}^{n - 1} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} \\  \end{cases}\]

Zauważmy iż \tg^2 \gamma_k = \tg^2 \nawias{-\gamma_k}, oznacza to, że wszystkie \tg \varphi_k są znane z góry, bo znane są z góry wartości kątów \gamma_k, tzn. \tg \varphi_k nie zależy od wybrania \varphi.
Rozpatrzmy funkcję K:\NN \to \RR, która w zasadzie jest ciągiem, którego pierwszy wyraz jest o indeksie 1. Ciąg jest dany wzorem

      \[K_n = \prod\limits_{k = 0}^{n - 1} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} =        \frac{1}{\sqrt{\prod\limits_{k = 0}^{n-1} \nawias{1 + \tg^2 \varphi_k }}} =        \frac{1}{\sqrt{\prod\limits_{k = 0}^{n-1} \nawias{1 + 2^{-2k} }}}\]

Można policzyć jednokrotnie wyrazy tego ciągu i zapamiętać.

      \[\begin{tabular}{|c|c|c|c|c|c|c|c|c|c|} \hline  $n$        & 1      & 2      & 3      & 4      & 5      & 6      & 7      & 8       & \ldots \\ \hline  $K_n$      & 0,7071 & 0,6325 & 0,6135 & 0,6088 & 0,6076 & 0,6074 & 0,6073 & 0,6073  & \ldots \\ \hline  \end{tabular}\]

Zamiast obliczać rekurencję z \tilde{x}_k i \tilde{y}_k możemy obliczać rekurencję z x_k oraz y_k, a na koniec raz pomnożyć przez K_{n}.
Wygląda na to, że jednego mnożenia trzeba będzie zatem dokonać.

Okazuje się, że jednak tego mnożenia także można uniknąć.
Zauważmy, że zamiast mnożyć na końcu to można by mnożyć na początku.
Wywnioskować to można spoglądając na rozpisane wcześniej 2 iteracje.

      \[\begin{matrix} \begin{cases} x_{k + 1} = x_{k} - y_{k}\tg \varphi_{k} \\ y_{k + 1} = y_{k} + x_{k}\tg \varphi_{k} \\ x_{0} = \tilde{x}_0 = A_x \\ y_{0} = \tilde{y}_0 = A_y \\ \end{cases}\\ \begin{cases} %x_{n} = \tilde{x}_{n} \prod\limits_{k = 0}^{n - 1} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} \\  %y_{n} = \tilde{y}_{n} \prod\limits_{k = 0}^{n - 1} \frac{1}{\sqrt{1 + \tg^2 \varphi_k}} \\  \tilde{x}_{n} = x_{n} K_{n} \\  \tilde{y}_{n} = y_{n} K_{n} \\  \end{cases} \end{matrix} \quad \to \quad  \begin{matrix} \begin{cases} x_{k + 1} = x_{k} - y_{k}\tg \varphi_{k} \\ y_{k + 1} = y_{k} + x_{k}\tg \varphi_{k} \\ x_{0} = \tilde{x}_0 = A_x K_{n} \\ y_{0} = \tilde{y}_0 = A_y K_{n} \\ \end{cases} \\ \begin{cases} \tilde{x}_{n} = x_{n} \\  \tilde{y}_{n} = y_{n} \\  \end{cases} \end{matrix}\]

Możemy przestać skupiać się na zmiennych akcentowanych tyldą \tilde{x}_k oraz \tilde{y}_k, a skupić na zmiennych nieakcentowanych tyldą, czyli x_k oraz y_k. Jak już było wspomniane \tg \varphi_k = \sgn(\beta_k)2^{-k}

      \[\begin{cases} x_{k + 1} = x_{k} - y_{k}\tg \varphi_{k} \\ y_{k + 1} = y_{k} + x_{k}\tg \varphi_{k} \\ \beta_{k + 1}     = \beta_{k}     - \sgn(\beta_k)\gamma_k \\ x_{0}     = A_x K_{n} \\ y_{0}     = A_y K_{n} \\ \beta_0 = \varphi\\ \end{cases}  \quad \to \quad  \begin{cases} x_{k + 1} = x_{k} - \sgn(\beta_k)y_{k}2^{-k} \\ y_{k + 1} = y_{k} + \sgn(\beta_k)x_{k}2^{-k} \\ \beta_{k + 1}     = \beta_{k}     - \sgn(\beta_k)\gamma_k \\ x_{0}     = A_x K_{n} \\ y_{0}     = A_y K_{n} \\ \beta_0 = \varphi\\ \end{cases}\]
życie algorytmu w celu obliczenia \sin \varphi oraz \cos \varphi

Teraz tak jak powiedziane było na początku. Jeśli punkt A=\nawias{A_x, A_y} będzie punktem A=\nawias{1, 0} to w wyniku otrzymamy punkt o współrzędnych A'=\nawias{\cos\varphi, \sin\varphi}. Decydując się na konkretną liczbę iteracji np. n=40 z góry wiadomo, że

      \[\begin{cases} x_{0} = K_{n} \\ y_{0} = 0 \\ \beta_0 = \varphi\\ \end{cases}.\]

Ograniczeniem jest to, że z góry na starcie musimy zdecydować ile zrobimy iteracji tego algorytmu. W ten sposób utworzyliśmy równania rekurencyjne, które pozwalają jednocześnie obliczyć wartość funkcji sinus i kosinus dla dowolnego kąta. Wzorami redukcyjnymi sprowadzamy dowolny kąt do kąta z przedziału \nawias{-\frac{\pi}{2}, \frac{\pi}{2}}, a następnie wykonujemy n iteracji poniższego układu. %równań.

      \[\begin{cases} x_{k + 1} = x_{k} - \sgn(\beta_k)y_{k}2^{-k} \\ y_{k + 1} = y_{k} + \sgn(\beta_k)x_{k}2^{-k} \\ \beta_{k + 1}   = \beta_{k}   - \sgn(\beta_k)\gamma_k \\ x_{0} = K_{n} \\ y_{0} = 0 \\ \beta_0 = \varphi\\ \end{cases}\]

Przybliżone wartości funkcji sinus i kosinus to

      \[\begin{cases} \cos \varphi \approx x_k  \\ \sin \varphi \approx y_k  \\ \end{cases}\]

Przykłady trybu rotacyjnego

Poniżej zwizualizowane i rozpisane są dwa przykłady. Jeden dotyczy obliczania
\sin 74\st oraz \cos 74\st, a drugi obliczenia \sin 55\st oraz \cos 55\st

Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com

Algorytm CORDIC w trybie wektorowym (vectoring mode)

W trybie rotacyjnym zadany punkt A=(1, 0) był obracany kolejno przez kąty \varphi_k, dla k = 0,1,\ldots ,n-1, aż osiągnięty zostanie punkt A'=(\cos x, \sin x). W trybie wektorowym zaczynam od dowolnego punktu A=(A_x, A_y), dla którego nie znamy kąta \varphi. Obracamy o kolejne kąty \varphi_k aż współrzędna y obrazu będzie w przybliżeniu równa A'_y \approx 0. Współrzędna x obrazu będzie wówczas równa A'_x \approx \sqrt{A_x^2 + A_y^2}, wynika to z twierdzenia Pitagorasa. Tryb wektorowy pozwala na obliczenie wartości funkcji \arctg, a dokładniej to \arctg \frac{A_y}{A_x}. Wyjdźmy od postaci

      \[\begin{cases} x_{k + 1} = x_{k} - y_{k}\tg \varphi_{k} \\ y_{k + 1} = y_{k} + x_{k}\tg \varphi_{k} \\ x_{0} = A_x K_{n} \\ y_{0} = A_y K_{n} \\ \end{cases}\]

Zauważmy, że dla kątów \varphi \in \nawias{-\frac{\pi}{2}, \frac{\pi}{2}},
współrzędna x_k > 0 to o jaki kąt należy obrócić w danej iteracji zależy
od znaku y_k. Kąt obrotu w danej iteracji to \varphi_k. Jeżeli y_k jest dodatnie to należy obrócić o kąt \varphi_k ujemny. Jak y_k jest ujemne to należy obrócić o kąt \varphi_k dodatni. Wynika nam z tego 

      \[\varphi_k = \epsilon_k \gamma_k = -\sgn(y_k) \gamma_k\]

Natomiast kąt \beta_k to będzie suma wykonanych obrotów \varphi_k
z przeciwnym znakiem. Jest tak dla tego, że suma katów obrotu sprowadzenia punkt A=(A_x, A_y) do punkt A'=(A'_x, 0) jest kątem o przeciwnym zwrocie do kąta od dodatniej półosi osi x do promienia wodzącego A. Otrzymujemy następująca rekurencję

      \[\begin{cases}  \beta_{k + 1} = \beta_k - \varphi_k\\ \beta_0 = 0\\ \end{cases}\]
Rendered by QuickLaTeX.com

Układ rekurencji dla tego algorytmu jest następujący

      \[\begin{cases} x_{k + 1} = x_{k} - y_{k}\tg \varphi_{k} \\ y_{k + 1} = y_{k} + x_{k}\tg \varphi_{k} \\ \beta_{k + 1} = \beta_k - \varphi_k \\ x_{0} = A_x K_{n} \\ y_{0} = A_y K_{n} \\ \beta_0 = 0 \\ \end{cases} \quad \to \quad  \begin{cases} x_{k + 1} = x_{k} + \sgn(y_k) y_{k}\tg \gamma_k \\ y_{k + 1} = y_{k} - \sgn(y_k) x_{k}\tg \gamma_k \\ \beta_{k + 1} = \beta_k + \sgn(y_k)\gamma_k \\ x_{0} = A_x K_{n} \\ y_{0} = A_y K_{n} \\ \beta_0 = 0 \\ \end{cases} \quad \to \quad  \begin{cases} x_{k + 1} = x_{k} + \sgn(y_k) y_{k}2^{-k}\\ y_{k + 1} = y_{k} - \sgn(y_k) x_{k}2^{-k}\\ \beta_{k + 1} = \beta_k + \sgn(y_k)\gamma_k \\ x_{0} = A_x K_{n} \\ y_{0} = A_y K_{n} \\ \beta_0 = 0 \\ \end{cases}\]

Pojawia się jednak problem bo skoro A_x i A_y jest dowolne to musimy wykonać mnożenie A_x K_{n} oraz A_y K_{n}. W tym wypadku też jest na to sposób, aby uniknąć tego mnożenia. Mianowicie po prostu o nim zapomnijmy :-). Skutkować to będzie tylko tym, że współrzędna x obrazu będzie błędna, a dokładniej to A'=\nawias{\frac{1}{K_{n}}\sqrt{A_x^2 + A_y^2}, 0}. Zwykle przyjmujemy, że A_x = 1, ale nie jest to konieczne. Jeśli przyjmiemy A_x = 1, to wtedy będziemy obliczać \arctg A_y.
Po n iteracjach odczytujemy znalezioną wartość funkcji \arctg

      \[\arctg \frac{A_y}{A_x} \approx x_n \\\]

Przykłady trybu wektorowego

Poniżej zwizualizowane i rozpisane są dwa przykłady. Pierwszy dotyczy obliczania
\arctg 2, a drugi obliczenia \arctg 3.

Rendered by QuickLaTeX.com
Rendered by QuickLaTeX.com

Algorytm obliczania funkcji f(x, y) = \sqrt{x^2 + y^2}

W trybie wektorowym istnieje także możliwość obliczenia pierwiastka kwadratowego. Jest jednak pewne ale, mianowicie nie można tego wykorzystać do obliczenia funkcji g(x) = \sqrt{x}. No chyba, że znamy rozkład x na a i b postaci x = a^2 + b^2. Zobaczmy, że nawet, gdyby przyjąć x = 0 w funkcji f(x, y) = \sqrt{x^2 + y^2}, to f(0, y) = \sqrt{0+ y^2} = |y|. Wychodzi na to, że znając y możesz policzyć jego wartość bezwzględna, czyli żaden z tego pożytek. Co więcej wymagało by to mnożenia na starcie, bo w tym przypadku nie możemy zignorować mnożenia K_{n}. Zignorowanie skutkowało by tym, że końcowy wynik byłby podzielony przez K_{n}.

W przypadku obliczania funkcji f(x, y) = \sqrt{x^2 + y^2}, bez mnożenia się nie obejdziemy. Jeśli mnożenie jest dostępne to wtedy faktycznie możemy obliczać wartości funkcji f(x,y) = \sqrt{x^2 + y^2}.

Podsumowanie – tryb rotacyjny i tryb wektorowy

Zestawiając ze sobą omówione dwa tryby, po chwili zauważamy, iż zarówno tryb rotacyjny jak i tryb wektorowy możemy zapisać wspólnie.

      \[\begin{tabular}{|c|c|} \hline  Tryb rotacyjny & Tryb wektorowy \\ \hline  $ \begin{cases} x_{k + 1} = x_{k} - \sgn(\beta_k)y_{k}2^{-k} \\ y_{k + 1} = y_{k} + \sgn(\beta_k)x_{k}2^{-k} \\ \beta_{k + 1}   = \beta_{k}   - \sgn(\beta_k)\gamma_k \\ x_{0} = A_x K_{n} \\ y_{0} = A_y K_{n} \\ \beta_0 = \varphi\\ \end{cases} $ & $ \begin{cases} x_{k + 1} = x_{k} + \sgn(y_k) y_{k}2^{-k}\\ y_{k + 1} = y_{k} - \sgn(y_k) x_{k}2^{-k}\\ \beta_{k + 1} = \beta_k + \sgn(y_k)\gamma_k \\ x_{0} = A_x K_{n} \\ y_{0} = A_y K_{n} \\ \beta_0 = 0 \\ \end{cases} $ \\ \hline  \end{tabular}\]
      \[\begin{cases} x_{k + 1} = x_{k} - d_k y_{k}2^{-k}\\ y_{k + 1} = y_{k} + d_k x_{k}2^{-k}\\ \beta_{k + 1} = \beta_k - d_k \gamma_k \\ \end{cases}\]

Tryb rotacyjny

      \begin{alignat*}{9} & d_k               & \    =    \ & \sgn(\beta_k) \\ & x_{0}             & \    =    \ & K_{n} \\ & y_{0}             & \    =    \ & 0 \\ & \beta_0           & \    =    \ & \varphi \\ & \cos \varphi      & \ \approx \ & x_n \\ & \sin \varphi      & \ \approx \ & y_n \\ \end{alignat*}

Tryb wektorowy

      \begin{alignat*}{9} & d_k               & \    =    \ & -\sgn(y_k) \\ & x_{0}             & \    =    \ & 1 \\ & y_{0}             & \    =    \ & A_y \\ & \beta_0           & \    =    \ & 0 \\ & \arctg A_y        & \ \approx \ & \beta_n \\ & \brak{\arctg A_y} &             &         \\ \end{alignat*}

W następnej części omówiony zostanie tryb liniowy. Następna cześć jest dostępna tutaj

https://www.kowalskimateusz.pl/algorytm-cordic-odcinek-2/