Lo primero, no se entiende claramente qué es lo que quieres hacer. ¿Obtener la raiz de 2 con precisión arbitraria? ¿Implementar tu propia función para cálculo de raices? ¿Calcular raices de raices de raices.... hasta cualquier nivel de profundidad y siempre con precisión arbitraria?
Lo segundo, me temo que algunos de estos objetivos son imposibles. Las raices que no salen exactas, como bien sabes, son irracionales. Esto quiere decir que el número de decimales es infinito. Obviamente no vas a poder encontrar un resultado "exacto" y en algún momento habrá que truncar los decimales y perder el resto.
Dicho esto, vamos con algo de teoria
El tipo float
El tipo float de python equivale al tipo double de otros lenguajes. Es un formato estandarizado utilizado por todos los computadores y calculadoras. Utiliza 64 bits para guardar un número real. Ya puedes imaginarte que con 64 bits no hay para infinitos decimales. De hecho, de esos 64 bits sólo 52 se usan para guardar la parte decimal (en binario), siendo otros 11 para la parte exponente, y 1 para el signo.
En el caso particular de raiz de 2, que es 1.4142135... o en binario sería 1.011010100000100111100110011... se almacenarían solo los primeros 52 bits que van tras la coma. La precisión por tanto es de 2-52 que viene a ser del orden de 10-16 lo que significa que al mostrarlo en base 10 sólo tendrías correctos los primeros 16 decimales, aproximadamente.
Además de las limitaciones de espacio del formato float
, aparecen efectos secundarios extraños debido a que el número debe ser convertido a binario antes de ser almacenado. El problema es que números de aspecto "inocente" en base 10, como por ejemplo 0.2, al ser escritos en binario producen una cantidad infinita de decimales (en particular, 0.2 en binario sería 0.0011001100110011... que si lo observas bien es periódico puro).
Esto no es tan extraño como pueda parecer a primera vista. En base 10 un número racional como 1/3 también produce infinitos decimales (0.333333...) mientras que en otras bases podría no ser así (en base 3 sería simplemente 0.1).
El caso es que eso causa también que números "normales" como por ejemplo precios en euros tales como 15,20€ al pasarlos a binario ocupen infinitos decimales, y por tanto en un float
no quepan y acaben siendo truncados provocando errores de redondeo.
El tipo Decimal
El tipo decimal.Decimal
implementa los números no enteros de otra forma, que se basa en la base 10 en lugar de la base 2, con el objetivo de que casos como los 15,20€ sean representables de forma exacta y no den lugar a los errores de redondeo de los float
.
Esto los hace atractivos para operaciones bancarias en las que es frecuente operar con un número finito de decimales, pero se quiere precisión total y no errores de redondeo. Sin embargo no es muy útil para el caso de los números irracionales como raiz de 2, ya que estos tienen infinitos decimales en cualquier base que elijas, y por tanto van a tener siempre errores de redondeo.
La principal ventaja de los Decimal
es que al menos te permiten especificar cuántos decimales quieres mantener. El espacio que ocupa un Decimal
no está prefijado a 64 bits como en los float
, sino que se reserva más espacio según la precisión con que quieras trabajar. Pero esa precisión va a tener siempre un límite (aunque lo eliges tú).
Tu código
En tu código aparece la siguiente instrucción:
x = (Decimal(math.sqrt(Decimal(2+a))))
Esa instrucción no tiene mucho sentido, pues lo que hace es:
- Cuando
a
vale 0, Decimal(2+a)
será Decimal(2)
. Eso convierte el 2 (que sería un entero) al formato Decimal()
y por tanto lo deja como 2.00000..00
con un número finito de decimales, digamos 60
- Seguidamente le pasa ese dato a
math.sqrt()
. Ya que esta función sólo trabaja con float
, lo primero que hace es volver a convertir el 2.000....0 en un flotante, lo que ya podría introducir un error de redondeo si hubiera habido decimales. En este caso no los hay y el 2 queda correctamente representado. Seguidamente calcula su raiz cuadrada y ya aparecen aqui los problemas, ya que la raiz cuadrada de 2 en el formato float
está limitada a 52 bits de la parte fracción. O sea, se pierde precisión a partir del dígito 16 al pasarlo a base 10.
- El resultado de esa raiz se convierte de nuevo a
Decimal()
. Pero ya es demasiado tarde porque ese resultado ya estaba "mal" a partir del dígito 16. Los 60 digitos que estás almacenando en x
son sólo correctos hasta el 16.
La solución de @aeportugal
Como @aeportugal ha indicado, sería más correcto utilizar el método .sqrt()
de la clase Decimal
ya que esa trabajará directamente con el dato decimal sin las conversiones intermedias a float
que causaban la pérdida de precisión.
Y de hecho podemos verificar que la respuesta en este caso es correcta:
>>> import decimal
>>> decimal.getcontext().prec = 60
>>> decimal.Decimal(2).sqrt()
Decimal('1.41421356237309504880168872420969807856967187537694807317668')
El resultado mostrado coincide correctamente hasta el último dígito con la solución que produce Wolfram Alpha cuando le preguntas cuánto vale "sqrt(2) with 60 digits
", siendo su respuesta:
1.41421356237309504880168872420969807856967187537694807317668
El resto de tus cuestiones
Entonces. Tenemos una forma de calcular sqrt(2) exacto hasta el decimal 60. ¿Para qué querías esto exactamente?
Parece que quieres computar la relación de recurrencia:
x_nuevo = sqrt(2+x_previo)
lo cual se haría con un bucle como el siguiente:
>>> a = decimal.Decimal(0)
>>> for i in range(N):
... x = (a+2).sqrt()
... a = x
dices que para N hasta 28 funciona, pero luego no. Lo he probado con N=30 y me ha salido:
x: Decimal('1.99999999999999999785987069335328430453055355679597492519638')
¿Cómo sabes si ese es o no es el resultado correcto? ¿Contra qué lo has comparado? Yo no he sido capaz a encontrar el resultado con Wolfram Alpha. Dudo mucho que una calculadora común te de un resultado mejor...
i = 3
yvalorN = 0
elwhile
no se va a ejecutar.