I logaritmi come li fanno i computer — parte 1

Dopo aver calcolato un logaritmo a mano, e in seguito a un commento in cui mi si chiedeva se quello che ho usato io fosse il metodo migliore, ho approfondito l'argomento e mi sono chiesto: ma come fanno le calcolatrici e i computer a calcolare i logaritmi?

E ho scoperto un mondo.

Un primo metodo che sfrutta il formato binario utilizzato dai computer è suggerito da Knuth. Ecco una piccola funzione in python che calcola il logaritmo in base 10 di un numero compreso tra 1 e 10 (incluso 1 e escluso 10).

def log10(x):
  
 a=x
 b=0
 d=1
 
 for n in xrange(1,Nbit): 
  a*=a
  d/=2.
  if a>=10:
   b+=d
   a/=10 
 return b

Dato che python non è un linguaggio a basso livello, questo programma simula soltanto quello che in effetti potrebbe fare una cpu (la divisione per 2 può essere molto rapida se il numero è memorizzato in forma binaria, per esempio). Il programmino funziona settando la variabile Nbit al numero di bit che vogliamo usare per memorizzare il logaritmo in forma binaria.

(Perché esistono anche i numeri binari con la virgola, eh. Per esempio, 1.01001 significa 1 + 2^-2 + 2^-5)

Il programma prende in input il numero x di cui calcola il logaritmo in base 10 nella forma b₁/2 + b₂/4 + b₃/8 + …

Utilizza una successione di valori a_k con la seguente proprietà:

a_k 10^{2k(b₁/2+…+b_k/2k)} = x^2k, che effettivamente è poco comprensibile. Ma se si applica il logaritmo in base 10 a destra e a sinistra, si porta l'esponente 2^k davanti al logaritmo e lo si semplifica, si vede che funziona, a parte un piccolo errore dovuto alla presenza di a_k. Giocando col valore di Nbit possiamo ottenere la precisione che vogliamo.

Questo, però, non è il metodo utilizzato dai computer o dalle calcolatrici: è un metodo didattico, ma niente di più. Ci fa capire quanto sia facile far fare il calcolo a un computer, rispetto a quello che si dovrebbe fare per eseguire a mano i vari passaggi seguendo il sistema di Briggs.

Knuth propone un altro algoritmo, che ho fatto molta fatica a capire. Eccolo qua, senza spiegazioni:

from math import log10

def log_10_Knuth(x):
 # N.B: 1<=y<2
 y=0
 z=x/2
 k=1

 while x<>1:
  if x-z < 1:
   z/=2
   k+=1
  else:
   x-=z
   z=x/2**k
   y+=log10(float(2**k)/(2**k-1))

 return y

Due note: x deve essere un numero appartenente all'intervallo [1,2), se non lo è, bisogna prima ridurlo. Come si vede, la funzione utilizza un'altra funzione di python che permette di calcolare il logaritmo in base 10 direttamente. E allora a cosa serve questa funzione, verrebbe da chiedersi?

In effetti questa funzione (che è sempre una simulazione di ciò che avviene in realtà) utilizza una tabella ausiliaria, memorizzata da qualche parte. Se abbiamo la possibilità di memorizzare log(2), log(4/3), log(8/7), e così via fino a ottenere la precisione desiderata, allora possiamo calcolare tutti i logaritmi che vogliamo utilizzando soltanto operazioni di somma, sottrazione e shift.

Si riesce a capire come funziona?