Algorytmy obliczeniowe
Przybli anie ilorazu wymiernego jego sko czonym rozwini ciem
m m m
X X
D `" 0 Ò! "{Ri}: Dm = D
"R 1Ò! Q = D = Dm "R H" X"R
i i i
i=0 i=0 i=0
" dokładno ilorazu  okre lona precyzj wyznaczenia liczby R0R1...Rm
standaryzacja: (ujemny dzielnik zmieni znaki D oraz X)
D := Dsgn D, X := X sgn D
m-1 m -1 -m -m
normalizacja:² d" D < ² Ò! ² d" d = D² < 1 & x = X²
n n+1
0 < q <1Ò! (1- q)(1+ q)(1+ q2)(1+ q4)...(1+ q2 ) = (1- q2 ) E" 1
procedura:
di =1- z Ò! di+1 = di (1+ z) = (1- z2 ) = di (2 - di )
zbie no procedury  kwadratowa
-s -2s
1- di = z < ² Ò!1- di+1 = z2 < ²
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 1
Algorytmy obliczeniowe
Przybli anie ilorazu sko czonym rozwini ciem w systemie dwójkowym
{02,11,00,0-1,0-2,...}U2 +{12,11,10, x-1, x-2, x-3,...}U2 = {02,01,10, x-1, x-2, x-3,...}U2
D>0 Ò! 2 + {1, x-1, x-2, x-3,...}U2 = {1, x-1, x-2, x-3,...}NB
Procedura
0o R0 :1- 2-s d" D1 = R0D <1, Q1 = R0 X (np. R0=2 m: 1 2 sd"2 mDd"1)
1o obliczaj Ri = 2 - Di oraz Qi+1 = QiRi a 1- 2-n < Di+1 = DiRi d"1
" liczba wiod cych jedynek liczb Di zostaje podwojona w ka dej iteracji
1- 2- p d" Di <1Ò!1- 2-2 p d" Di+1 <1
" 1- 2-s d" D1 < 1 Ò! wzgl dn dokÅ‚adno 2 n ilorazu Q H" X R0R1...
zapewnia m = +1 iteracji
îÅ‚log n / sÅ‚Å‚
2
Ò! pierwszy mno nik R0  wyznaczony z dokÅ‚adno ci s+log2s bitów
R0 = f (D)  z matrycy ROM o rozmiarze 2s×( s+log2s) bitów
" przy pieszenie mno enia u ycie krótszych mno ników Ri
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 2
Algorytmy obliczeniowe
Obliczanie odwrotno ci dzielnika
metoda iteracyjna Newtona-Raphsona
y y=f(xi)
x
xi+2 xi+1 xi
y=f (xi+1)(x xi+1)+f(xi+1) y=f (xi)(x xi)+f(xi)
kolejne przybli enia miejsca zerowego f(x) okre la równanie rekurencyjne
f (xi )
xi+1 = xi - ,
2
f (xi )
-1
w odniesieniu do funkcji f (x) = x - D przybiera posta
xi+1 = xi (2 - Dxi )
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 3
Algorytmy obliczeniowe
Zbie no metody mno enia przez odwrotno dzielnika
Niech x0 = a oraz Da=1-q Ò! x1 = a(1+ q) = a(1- q)-1(1- q2)
i-1 i
xi = a(1+ q)...(1+ q2 ) = a(1- q)-1(1- q2 )
i i-1 i i
xi+1 = a(1- q)-1(1- q2 ){2 - Da[(1+ q)...(1+ q2 )]} = a(1- q)-1(1- q2 )(1+ q2 )
i
| q |<1Ò! lim xi = lima(1- q)-1(1- q2 ) = a(1- q)-1 = D-1
i" i"
1
" dzielnik znormalizowany d"| D |< 1 Ò! zbie no , je eli |a|<2 i aD>0.
2
-1
" zbie no kwadratowa  je li ´ = D - xi, to
i
-1 2 2
´ = D - xi+1 = D-1 - (D-1 -´ )[2 - D(D-1 -´ )] = D´ < ´
i+1 i i i i
" szybko zbie no ci zale y od dokładno ci pierwszego przybli enia
j+1 j -1
1
(k -1)2- j d" D < k2- j Ò! optymalne x0 (k) = 2 [2 + k - ]
2
wada  mniejsza dokładno ni uzyskiwana w dzieleniu sekwencyjnym
niezb dna korekcja wyniku dodatkowe działania arytmetyczne
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 4
Algorytmy obliczeniowe
Obliczanie pierwiastka i odwrotno ci pierwiastka kwadratowego
1
Liczba pierwiastkowana jest znormalizowana d" A < 1.
4
f (xi )
metoda iteracyjna Newtona  równanie rekurencyjne: xi+1 = xi -
2
f (xi )
Obliczanie pierwiastka kwadratowego:
Je li f(x)=x2 A, to x=sqrt(A) ("A) i wtedy f 2 (x)=2x, wi c
1
A
xi2 - A
2
1
xi+1 = xi - = xi +
2xi 2 xi
wada: konieczno dzielenia
Obliczanie odwrotno pierwiastka kwadratowego:
2 3
2
f x = x- - A i wtedy f x = - x- oraz
2
-
2
1
xi - A
1
3 2
xi+ = xi - = xi - xi A
-
- xi
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 5
Algorytmy obliczeniowe
Obliczanie odwrotno ci pierwiastków wy szych stopni
k
Je eli f (x) = x-k - A, to x = A-
2
Poniewa wtedy f x = -kx-k- , wi c otrzymujemy
-k
xi - A
xi+ = xi - = xi k + - xikA
-k- k
- kxi
Obliczenia wymagaj wielokrotnego mno enia / pot gowania je li k>2.
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 6
Algorytmy obliczeniowe
CORDIC  algorytm Voldera (1)
Obrót wektora zaczepionego w punkcie (0,0) przestrzeni kartezja skiej
y
(xb=Rcos(Ä…+´), yb=Rsin (Ä…+´)
yb
(xa=RcosÄ…, ya=Rsin Ä…)
ya
´
Ä…
x
xb xa
Z to samo ci trygonometrycznych
cos(Ä…+²)=cosÄ…cos² sinÄ…sin²
sin(Ä…+²)=sinÄ…cos²+cosÄ…sin²
wynika, e:
xb=Rcos(Ä…+´)=RcosÄ… cos´ -RsinÄ… sin´ = xacos´ -yasin´
yb=Rsin(Ä…+´)=RsinÄ… cos´ +R cosÄ… sin´ = yacos´ +xasin´
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 7
Algorytmy obliczeniowe
CORDIC  algorytm Voldera (2)
J.Volder (1956, sterowanie samolotu B-58)
Podstawiaj c t=tg´ (cos 2´ =1+t2), otrzymamy dla k tów w I wiartce:
(1+ t2)-1R cos(Ä… + arctgt) = RcosÄ… - tRsinÄ…
(1+ t2)-1Rsin(Ä… + arctgt) = RsinÄ… + tRcosÄ…
W krokach iteracji, to wynikiem obrotu wektora (xi,yi) o k t arctgti jest:
(1+ ti2)-1 xi+1 = xi - ti yi , i=0,1,&
(1+ ti2)-1 yi+1 = yi + tixi, i=0,1,&
Obie współrz dne s jednakowo skalowane, wi c w pojedynczym kroku mo na
zignorowa wydłu enie wektora, dokonuj c korekty w ostatnim kroku oblicze
* * * * * *
xi+1 = xi - ti yi oraz yi+1 = yi + tixi
* *
gdzie x0 = x0 , y0 = y0 oraz
n-1 n-1
* *
xn = xn 1+ ti2 , yn = yn 1+ ti2 .
" "
i =0 i =0
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 8
Algorytmy obliczeniowe
CORDIC  algorytm Voldera (3)
Podobne to samo ci dotycz funkcji hiperbolicznych sinh i cosh (wzór Eulera)
cosh(Ä…+´)=coshÄ… cosh´-sinhÄ… sinh´
sinh(Ä…+´)=sinhÄ… cosh´+coshÄ… sinh´
gdzie (wzór Eulera: expix=cosx+isinx)
2sinh x = -2isin x = exp x - exp(-x)
2cosh x = 2cos x = exp x + exp(-x)
exp x = cosh x + sinh x
Warto ci t=tg´=Ä…2 n, mo na Å‚atwo tablicowa i wtedy wszystkie obliczenia
mo na wykona za pomoc dodawania, odejmowania i przesuni cia.
Zale nie od znaku k ta wyró nia si obliczenia
" w trybie obrotu (rotation mode), gdy k t jest dodatni,
" w trybie normowania (vectoring mode), gdy k t jest ujemny
 jego wynikiem jest obliczenie długo ci wektora
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 9
Algorytmy obliczeniowe
CORDIC  algorytm Voldera (4)
Uogólnienie dla funkcji trygonometrycznych i hiperbolicznych.
trzecia zmienna  zi odległo k towa wektora od osi [0,x):
xi +1 = xi + mÃi ti yi
yi+1 = yi -Ãi ti xi
zi +1 = zi +Ãi (1/ m) arctg m ti
gdzie ti=2 S(m,i)  przyj ta sekwencja iteracji przyrostów
tryb obrotu (zi0) tryb normowania (yi0)
Ãi=-sign zi Ãi=sign (xiyi)
2 2
xnK x0 + y0
xnK(x0cosz0 y0sinz0)
m=1 arctg2 k
ynK(y0cosz0 x0sinz0) yn0
zn0 znz0 arctg2(y0/x0)
2 2
xnK x0 - y0
xnK(x0coshz0 y0sinhz0)
m= 1 tanh 12 k
ynK(y0coshz0 x0sinhz0) yn0
zn0 znz0 tanh 12(y0/x0)
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 10
hiperboliczny
trygonometr.
Algorytmy obliczeniowe
CORDIC  realizacja układowa
x
y
LUT
z
Zalety algorytmu CORDIC
obliczanie funkcji elementarnych za pomoc prostych działa arytmetycznych
prosta implementacja układowa algorytmu (Cyrix, procesory DSP)
Wada  wolna zbie no , konieczno wykonania du ej liczby oblicze ,
wersja ulepszona CORDIC 2.
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 11
Algorytmy obliczeniowe
Inne metody obliczania warto ci funkcji przest pnych
tablica odniesie (look-up table)
" zapami tanie warto ci funkcji jednej zmiennej z dokładno ci do n bitów
 matryca ROM o rozmiarze n× 2n bitów (dla n > 23 rozmiar > 128 Mb)
rozwini cie w szereg Taylora
" ró ne algorytmy dla poszczególnych funkcji
" wolna zbie no szeregu Taylora (zale y silnie od warto ci argumentu)
rozwini cia funkcji przest pnych w postaci ułamków wymiernych
" powszechnie stosowane w implementacjach programowych
" mog by bardzo skuteczne w realizacjach sprz towych, je li w dyspozycji
s szybkie zmiennoprzecinkowe sumatory i układy mno ce
algorytmy oparte na przybli eniach wielomianowych z u yciem tablic odniesie
" domena argumentu funkcji f(x) podzielona na przedziały równej długo ci
" warto ci graniczne f (xi ) w punktach podziału xi s w tablicy odniesie
" warto ci wewn trz przedziałów obliczane na podstawie aproksymacji
wielomianowej p(x - xi ) funkcji f (x - xi )
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 12
Algorytmy obliczeniowe*
Procedura normalizacji addytywnej
F  obliczana funkcja, xi+1 = xi - bi
i-1 i-1
{xi = x0 -
"b } 0 Ò! {yi = F(x0 - xi ) = F("b ) F(x0 )},
j j
j=0 j=0
i-1 i-1 i-1
lub {xi = x0 - g(bj )} 0 Ò! {yi = F( g(bj )) =
" " "b F(x0 )},
j
j=0 j=0 j=0
" rozwini cie w szereg Taylora funkcji g(x) wokół punktu x0
d F(x) d2 F(x)
2 2
F(x0 - µ ) = F(x0) - µ |x +µ |x +o(µ )
0 0
d x d x2
m-1 m-1
" xm = x0 -
"b = µ Ò! bÅ‚ d przybli enia F(x0) przez ym = F("b )
j j
j=0 j=0
d F(x)
´ = F(x0) - F(x0 - µ ) H" µ |x
0
d x
" µ 0 Ò! ´ 0 (´ zale y od pochodnej funkcji F w punkcie x0)
" {bi}  tak dobra aby w m krokach uzyska ´<2 n
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 13
Algorytmy obliczeniowe*
Procedura normalizacji multyplikatywnej
F  obliczana funkcja, xi+1 = xibi
i-1 i-1
-1
{xi = x0 j} 1Ò! {yi = F(x0xi -1) = F(
"b "b ) F(x0)}
j
j=0 j=0
m-1 m-1
-1
" 1- µ = x0 j Ò! F(
"b "b ) = F(x0(1- µ )-1)
j
j=0 j=0
" µ 0 Ò! F(x0(1- µ)-1) H" F(x0 + µ x0)
" bÅ‚ d przybli enia F(x0) rozwini cie w szereg Taylora dla x=x0(1 µ) 1
d F(x)
´ = F(x0(1- µ )-1) - F(x0 ) E" F(x0 + µ x0 ) - F(x0 ) = µ x0 |x +o(µ )
0
d x
" µ 0 Ò! ´ 0 (´ zale y od x0 i pochodnej funkcji F w punkcie x0)
" {bi}  tak dobra aby w m krokach uzyska ´<2 n
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 14
Algorytmy obliczeniowe*
Funkcja wykładnicza
" bi = (1+ si 2-i ), si "{1,0,1} (prostsze mno enie)
i-1 i-1 i-1
{xi = x0 -
"lnb } 0 Ò! {exp(x0 - xi ) = exp("lnb ) = "b exp(x0 )},
j j j
j=0 j=0 j=0
m-1 m-1
" ln(1- 2- j ) d" x0 d" ln(1+ 2- j ) Ò! istnieje ci g {si} taki, e xm 0
" "
j=1 j=1
Å„Å‚0, gdy
xi - ln(1+ 2-i ) < 0
0 < x0 <1Ò! si = Ò! {xi} 0
òÅ‚
ół1, gdy xi - ln(1+ 2-i ) e" 0
i 1+2 i ln(1+2 i) 1 2 i ln(1 2 i)
0 10,00000 00000 0,10110 00110 0 
1 1,10000 00000 0,01100 11111 0,10000 00000  0,10110 00110
2 1,01000 00000 0,00111 00100 0,11000 00000  0,01001 00111
3 1,00100 00000 0,00011 11001 0,11100 00000  0,00100 01001
4 1,00010 00000 0,00001 11110 0,11110 00000  0,00010 00010
5 1,00001 00000 0,00001 00000 0,11111 00000  0,00001 00001
6 1,00000 10000 0,00000 10000 0,11111 10000  0,00000 10000
7 1,00000 01000 0,00000 01000 0,11111 11000  0,00000 01000
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 15
Algorytmy obliczeniowe*
Funkcja wykładnicza (2)
" zapami tanie ln(1Ä… 2-i ) z dokÅ‚adno ci 2 n  matryca ROM 2n×n
" xn < 2-n exp(x0)  exp(x0 xn) <2-n eln 2 = 2-n+1
" dla argumentu spoza przedziału (0,1)
x = (xlog e)ln 2 x0 = (x loge- ðÅ‚x logeûÅ‚) ln 2 oraz
x x loge x0
ðÅ‚x logeûÅ‚ln 2+ x0
e = e = e = 2ðÅ‚x logeûÅ‚ e
wystarczaj cym zakresem wykładnika jest 0 d" x0 < ln 2
reguła wyboru, nie wymagaj ca cz stego odejmowania
" rozwini cie ln(1+z) w szereg Taylora wokół warto ci z = si 2-i
1 1 1
ln(1+ si 2-i ) = si 2-i - (si )2 2-2i + (si )3 2-3i - (si )4 2-4i + ...
2 3 4
w i-tym kroku procedury przy odejmowaniu ln(1+ si 2-i )
jest zerowany bit o wadze 2-i zbie no jest liniowa
odejmowanie xi - ln(1+ 2-i ) potrzebne gdy i-ty bit przybli enia xi jest 1
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 16
Algorytmy obliczeniowe*
Funkcja logarytmiczna
procedura normalizacji multyplikatywnej
i-1 i-1
-1
{xi = x0 j} 1Ò! {yi = ln(x0xi -1) = ln(
"b "b ) ln(x0 )}
j
j=0 j=0
" bi = 1+ si 2-i, gdzie si "{1,0,1} umo liwiaj uproszczenie mno enia
m-1 m-1 m-1
j -1
0,29 d" (1- s 2- j ) d" (1+ 2- j ) d" 4,77,
"(1- 2- ) d" x0 = " "
j
j=0 j=0 j=0
" znormalizowane ułamki dodatnie nale do dziedziny
x
" x0 nie spełnia ograniczenia przeskalowanie ln(x0 2E ) = ln x0 + Ex ln 2
-
" si "{0,1} Ò!1d" x01 d" 4,77 prosta reguÅ‚a wyboru si
1- 2- j d" xi d" 1- 2-( j+1) Ò! (s = 1) '" ("i d" l d" j -1:sl = 0)
j
je li xi ma je"i wiod cych jedynek ( xid" j e" 1- 2- j ), to si = 1, wi c
x = xi (1+ 2- j ) e" (1- 2- j )(1+ 2- j ) = 1- 2-2 j
j+1
d ln x
-1
" |x = x0 bÅ‚ d przybli enia ln x jest równy bÅ‚ dowi µ=1 xm
0
d x
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 17
Algorytmy obliczeniowe*
Metody przy pieszania oblicze
" rozwini cie funkcji ln(1+x) w szereg Taylora dla x = si2-i daje
ln(1+ si 2-i ) = si 2-i + o(2-2i-1)
dla i > n/2 ln(1+ si 2-i ) = si 2-i z dokładno ci nie gorsz ni 2 n
" przy obliczaniu warto ci funkcji wykładniczej xi = 0,0...0zi zi+1...,
z dokładno ci 2 n dla m>i>n/2 oraz si=zi mamy
m-1 m-1
j
ym = yi z 2- j + o(2-2 j-1)].
"(1 + z 2- ) = yi[1 + "
j j
j=i j=i
n/2 ostatnich kroków mo na zast pi mno eniem ym = yi (1+ xi ).
" przy obliczaniu funkcji logarytmicznej 1- xi = 0,0...0z z ....
je"i j+1
z dokładno ci 2 n dla m>i>n/2 oraz si=zi mamy (szereg Taylora)
m-1 m-1
j
ym = yi - z 2- j = yi - (1- xi )
"ln(1+ z 2- ) H" yi - "
j j
j=i j=i
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 18
Algorytmy obliczeniowe*
Funkcje hiperboliczne (1)
x
1 1
sinh x = (e - e-x ), cosh x = (ex + e-x )
2 2
procedura normalizacji addytywnej
xi+1 = xi - g(bi ), yi+1 = yi bi
" na podstawie to samo ci 1+ x = 1- x2 exp(tgh-1x) dla |x| `" 1 mamy
m-1 m-1 m-1
1- (si 2-i )2 exp( tgh-1(si 2-i )).
"(1+ si 2-i ) = " "
i=1 i=1 i=1
-1
" je li bi =1+ si 2-i oraz g(bi ) = tgh (si 2-i ), to (s0=0)
m-1 m-1
-1
xm = x0 -
"tgh (si 2-i ) ym = y0"(1+ si 2-i ) = y0Km exp(x0 - xm )
i=1 i=0
m-1
- -1
" si "{1,1} Ò! Km = 1- 2-2i oraz Km1 E" K =1,205136 Ä… 2 40
"
>20
i=1
m-1
ym = y0K
"(1+ si 2-i ) E"y0 K exp(x0 ) .
i=0
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 19
Algorytmy obliczeniowe*
Funkcje hiperboliczne (2)
Niech ti+1 = ti (1- si 2-i )
" yi i ti s zbie ne przy xi 0,
1 1
" zbie ne s wi c zi = ( yi + ti ) i wi = ( yi - ti )
2 2
zi+1 = zi + si 2-i wi, wi+1 = wi + si 2-i zi
" granicami zbie no ci s
-1
1 1 1 1
0 0
zm K y0 ex + t0 ex = ( y0 + t0 ) cosh x0 + (y0 - t0 ) sinh x0,
2 2 2 2
-1 x0 1 x0 1
1 1
wm K y0 e - t0 e = (y0 + t0 ) sinh x0 + (y0 - t0 ) cosh x0.
2 2 2 2
-1
" y0 = t0 = K Ò! zm cosh x0 oraz wm sinh x0
" zbie no bezwzgl dna ci gu {xi} do zera nie jest jednostajna
-1 -1
1
" tgh (2-(i+1) ) < tgh (2-i ) Ò! niektóre kroki musz by powtórzone
2
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 20
Algorytmy obliczeniowe*
Funkcje trygonometryczne (1)
algorytmy oparte na wzorze Eulera
eix = cos x + isin x, i = -1
procedura normalizacji addytywnej
xi+1 = xi - g(bi ), yi+1 = yi bi
eix przybiera warto ci zespolone b zespolone: bi = 1+ i si 2-i, wtedy
i
m-1 m-1 m-1
ym = y0 1+ (si 2-i )2 expëÅ‚i
ìÅ‚ ÷Å‚
"(1+ i si 2-i ) = y0" "arctg si 2-i öÅ‚
i=0 i=0 íÅ‚ i=0 Å‚Å‚
m-1
g(bi ) = arctg(si 2-i ) Ò!
"arctg(s 2-i ) = x0 - xm (arctg 2 i w matrycy ROM)
i
i=0
" arctg jest funkcj monotonicznie rosn c xi+1 = xi - arctg(si 2-i ) 0
m-1 m-1
ym = y0 1+ (si 2-i )2 expëÅ‚i
ìÅ‚ ÷Å‚
" "arctg si 2-i öÅ‚ y0 Kmexp(i x0)
i=0 íÅ‚ i=0 Å‚Å‚
m-1
-i
1
" obszar zbie no ci metody  | x0 |d" Ém =
"arctg(2 ) , Ém>20 E" 1,743 > Ä„
2
i=0
1
u yteczna domena 0 d" x0 d" Ä„ jest zawarta w przedziale zbie no ci
2
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 21
Algorytmy obliczeniowe*
Funkcje trygonometryczne (2)
yi = wi + i zi reguły iteracji dla cz ci rzeczywistej wi i urojonej zi
wi+1 + i zi+1 = (wi + i zi )(1+ i si 2-i ) = (wi - si 2-i zi ) + i(zi + si 2-i wi )
" szybka zbie no metody, ale zmienne Km
Å„Å‚
1, gdy xi > arctg(2-i ),
ôÅ‚0, gdy | xi |d" arctg(2-i ),
si = Ò! Km = Km (s1,s2,...,sm )
òÅ‚
ôÅ‚1, gdy xi < - arctg(2-i ).
ół
1
" zbie no wolniejsza, ale stałe Km ( Km>16 E" 1,6468 > Ą )
2
m-1
1, gdy xi e" 0,
Å„Å‚
si = 1+ 2-2i
òÅ‚1, gdy x < 0. Ò! K = "
m
i=0
ół i
" y0 = 1/ K Ò! wm = cos x0 oraz zm = sin x0
m
" liniowa zbie no metody
" wyrazy rozwini cia arctg w szereg Taylora dla i>n/3 s rz du < o(2 n)
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 22
Algorytmy obliczeniowe*
Odwrotne funkcje trygonometryczne (1)
procedura normalizacji multyplikatywnej
xi+1 = xi bi, yi+1 = yi + g(bi ) .
bi = 1+ i si 2-i g(bi ) = arctg(si 2-i ) (jak w obliczaniu funkcji sin i cos)
yi  rzeczywiste, xi = ui + i vi  warto ci zespolone, przy tym
ui+1 = ui - si 2-i vi, vi+1 = vi + si 2-i ui, yi+1 = yi + arctg(si 2-i ).
po m krokach algorytmu mamy
m-1 m-1
xm = x0 1+ (si 2-i )2
"(1+ i si 2-i ) = x0 Km (s)exp(i¸ ), K (s) = "
m m
i=0 i=0
si "{1,1} Ò! Km (s) = Km i wtedy w m krokach iteracji
xm x0 Km exp(i¸ ) = (u0 + i v0 )Km[cos¸ + i sin¸ ]
m-1
ym = y0 + arctg(si 2-i ) = y0 +¸ y0 +¸ .
"
m
j=0
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 23
Algorytmy obliczeniowe*
Odwrotne funkcje trygonometryczne (2)
Wynika st d, e
- -
um Km1 = u0 cos¸ - v0 sin¸ , vm Km1 = v0 cos¸ + u0 sin¸ .
m m m m
Aby wyznaczy arcus tangens nale y zapewni zbie no um do zera.
Zbie no procedur jest zapewniona przy v0 = 1. Wówczas przy u0 = c = tg¸
oraz y0 = 0 warto ym d y do ¸ = arctg c.
Zapewnienie dobrych warunków zbie no ci dla pozostałych odwrotnych
funkcji trygonometrycznych jest znacznie trudniejsze.
I tak, dla obliczenia warto ci funkcji arcsin c nale y tak dobra wektor s, aby
um d yÅ‚o do sin ¸, podobnie dla funkcji arccos c nale y znale wektor s
zapewniaj cy zbie no vm do cos ¸.
© Janusz Biernat, 07-06-Obliczenia numer.doc, 27 grudnia 2004 NUM 24