Conosciamo una proprietà fondamentale che lega seno e coseno, e cioè la somma tra il quadrato del seno di un angolo e il quadrato del coseno dello stesso angolo è uguale a 1. Se indichiamo con x il coseno e con y il seno, la proprietà si traduce in questa formula: x2+y2 = 1, che è l'equazione di una circonferenza con centro nell'origine e raggio uguale a 1, detta anche circonferenza goniometrica. Quella sulla quale si calcolano, appunto, seno e coseno.
Anche seno iperbolico e coseno iperbolico sono legati da una analoga proprietà, che, se indichiamo ancora con x il coseno e con y il seno, diventa questa: x2-y2 = 1, che è l'equazione di una iperbole.
Ecco una figurina:
In blu la circonferenza delle funzioni goniometriche normali (si dice circolari, in realtà), in rosso l'iperbole delle funzioni iperboliche.
L'algoritmo che abbiamo usato per calcolare le funzioni goniometriche si basava sulle formule di somma e sottrazione per seno e coseno. Ne esistono di analoghe anche per le funzioni iperboliche; eccole qua:
cosh(x±y) = cosh(x)cosh(y)±sinh(x)sinh(y),
sin(x±y) = sinh(x)cosh(y)±cosh(x)sinh(y).
Cambia un solo segno: se cambiamo quindi un solo segno nei nostri programmi, siamo in grado di calcolare in un attimo anche queste nuove funzioni iperboliche.
In realtà non è così, c'è un problema. Tutta la costruzione del ventaglio funzionava grazie a una proprietà: se immaginiamo un ventaglio infinito, ogni angolo è minore della somma di tutti quelli che lo seguono. Nel caso di un ventaglio finito, invece, ogni angolo è minore della somma di tutti quelli che lo seguono più l'ultimo angolo (insomma, l'angolo più piccolo viene contato due volte). Questa è la proprietà che viene usata per dimostrare che ogni angolo è approssimabile con un errore minore dell'ampiezza dell'angolo più piccolo. Per il ventaglio iperbolico questa proprietà non vale più.
Se proviamo ad avvicinarci a un angolo utilizzando il ventaglio iperbolico, l'approssimazione non funziona. Per esempio, ecco cosa succede durante l'avvicinamento al valore 0.549:
+ 0.5493061443 0.5493061443 - 0.2554128119 0.2938933325 + 0.1256572141 0.4195505466 + 0.0625815715 0.4821321181 + 0.0312601785 0.5133922966 + 0.0156262718 0.5290185683 + 0.0078126590 0.5368312273 + 0.0039062699 0.5407374971 + 0.0019531275 0.5426906246 + 0.0009765628 0.5436671874 + 0.0004882813 0.5441554687 + 0.0002441406 0.5443996093 + 0.0001220703 0.5445216797 + 0.0000610352 0.5445827148 + 0.0000305176 0.5446132324 + 0.0000152588 0.5446284912 + 0.0000076294 0.5446361206 + 0.0000038147 0.5446399353 + 0.0000019073 0.5446418426 + 0.0000009537 0.5446427963 + 0.0000004768 0.5446432731 + 0.0000002384 0.5446435116 + 0.0000001192 0.5446436308 + 0.0000000596 0.5446436904
In prima colonna i nuovi valori del ventaglio, ottenuti tramite l'arcotangente iperbolica, in seconda colonna il valore che vorremmo approssimare: si vede che il risultato è abbastanza distante da 0.549, l'errore è già sulla terza cifra decimale, ed è decisamente maggiore della più piccola ampiezza del ventaglio.
Come fare, quindi?
Dobbiamo truccare il ventaglio, inserendo qualche valore in più: duplichiamo cioè qualche settore. Per mantenere la stessa lunghezza dobbiamo quindi rinunciare un po' alla precisione. Ecco un ventaglio che funziona: l'asterisco sta a indicare le righe duplicate, in prima colonna i valori, in seconda la somma di tutti i valori successivi, in terza la differenza tra le due colonne più l'ultimo valore.
0.5493061443 0.5771220351 0.0278197054 0.2554128119 0.3217092232 0.0663002260 0.1256572141 0.1960520090 0.0703986096 0.0625815715 0.1334704376 0.0708926808 * 0.0625815715 0.0708888661 0.0083111093 0.0312601785 0.0396286876 0.0083723238 0.0156262718 0.0240024158 0.0083799588 0.0078126590 0.0161897569 0.0083809126 * 0.0078126590 0.0083770979 0.0005682537 0.0039062699 0.0044708281 0.0005683729 0.0019531275 0.0025177006 0.0005683878 0.0009765628 0.0015411378 0.0005683897 0.0004882813 0.0010528565 0.0005683899 * 0.0004882813 0.0005645752 0.0000801086 0.0002441406 0.0003204346 0.0000801086 0.0001220703 0.0001983643 0.0000801086 0.0000610352 0.0001373291 0.0000801086 * 0.0000610352 0.0000762939 0.0000190735 0.0000305176 0.0000457764 0.0000190735 0.0000152588 0.0000305176 0.0000190735 * 0.0000152588 0.0000152588 0.0000038147 0.0000076294 0.0000076294 0.0000038147 0.0000038147 0.0000038147 0.0000038147 * 0.0000038147
Ora l'avvicinamento a 0.549 funziona:
+ 0.5493061443 0.5493061443 - 0.2554128119 0.2938933325 + 0.1256572141 0.4195505466 + 0.0625815715 0.4821321181 + 0.0625815715 0.5447136895 + 0.0312601785 0.5759738680 - 0.0156262718 0.5603475963 - 0.0078126590 0.5525349373 - 0.0078126590 0.5447222784 + 0.0039062699 0.5486285482 + 0.0019531275 0.5505816757 - 0.0009765628 0.5496051129 - 0.0004882813 0.5491168316 - 0.0004882813 0.5486285503 + 0.0002441406 0.5488726910 + 0.0001220703 0.5489947613 + 0.0000610352 0.5490557964 - 0.0000610352 0.5489947613 + 0.0000305176 0.5490252789 - 0.0000152588 0.5490100201 - 0.0000152588 0.5489947613 + 0.0000076294 0.5490023907 - 0.0000038147 0.5489985760 + 0.0000038147 0.5490023907
Come si vede, l'errore con cui approssimiamo 0.549 è minore del valore dell'ultimo settore del ventaglio.
Ora, finalmente, possiamo aggiustare l'algoritmo per seno e coseno iperbolici. Dato che alcuni settori sono duplicati, non possiamo ad ogni passo fare una divisione per due come facevamo prima: dobbiamo fare attenzione agli indici duplicati.
Ecco l'algoritmo:
Z0 = T,
C0 = Q,
S0 = 0,
Zm+1 = Zm-DmBm,
Cm+1 = Cm+EmDmSm,
Sm+1 = Sm+EmDmCm.
Ora spiego i termini: D è uguale a +1 oppure -1 in modo da avere lo stesso segno di Z; i valori B sono le ampiezze del nuovo ventaglio, Q è la costante cosh(B0)…cosh(B23) e i valori E sono le tangenti iperboliche dei valori B, cioè sono potenze intere di 1/2, ma non sono in ordine progressivo perché alcune sono state duplicate, come abbiamo detto prima.
Ed ecco un programma (già che ci siamo, gli facciamo anche stampare il valore di Q, così vediamo quanto vale):
from math import atanh,cosh,sinh from math import atanh angoh=[] indici=[1,2,3,4,4,5,6,7,7,8,9,10,11,11,12,13,14,14,15,16,16,17,18,18]
# compila la tabella per il CORDIC
for i in indici:
angoh+=[atanh(1./2**i)]
#calcolo del numero magico Q
Q=1
for a in angoh:
Q*=cosh(a)
print "Q =",Q
def senoh(x):
z0=x
c0=Q
s0=0
for i,b in enumerate(angoh):
if z0>0:
d=1
else:
d=-1
z1 = z0-d*b
c1 = c0+d*s0/2**indici[i]
s1 = s0+d*c0/2**indici[i]
z0,c0,s0 = z1,c1,s1
return s0
v1=senoh(1)
v2=sinh(1)
print v1,v2,abs(v1-v2)
La lista indici contiene gli esponenti della frazione 1/2 usati per creare il ventaglio. Ed ecco l'output:
Q = 1.20753405668 1.17520248124 1.17520119364 1.28759919193e-06
Dato che sinh(x)+cosh(x) = ex, il calcolo dell'esponenziale è immediato (tanto più che il nostro algoritmo calcola seno e coseno iperbolico contemporaneamente). Siamo a un passo dal calcolo del logaritmo, finalmente.
6 commenti:
Meno male che avevi ancora poco da dire! (appena ho un po' di tempo libero torno a rileggermi tutta la serie)
Perché duplichiamo proprio quei valori e non altri?
@Knulp: la sto tirando in lungo :-)
@Anonimo: non lo so, ho visto che fanno così. Probabilmente la risposta è "perché così funziona" :-) L'importante è che sia valida la proprietà richiesta.
Con quelle duplicazioni nessun elemento della lista è superiore alla somma dei successivi (ho verificato) e mi pare che sia questo il 'guaio' che può impedire al ventaglio di convergere, perché non permette di tornare sempre 'indietro' a sufficienza, se ce ne fosse bisogno.
Almeno così mi viene da pensare.
Sì, è così, le duplicazioni fanno sì che quella proprietà sia valida. Non so se siano state trovate per tentativi o se ci sia qualche altro metodo.
Penso che ci siano più possibilità.
Forse quella scelta è in qualche modo ottimale oppure è stata effettuata 'mettendo le pezze' e verificando che tutto funzionasse :-)
Posta un commento