Articulo de referencia

algoritmos de raíz cuadrada

Los algoritmos de raíz cuadrada calculan la raíz cuadrada no negativa. S {\displaystyle {\sqrt {S}}} de un número real positivo S {\displaystyle S} . Dado que todas las raíces c...

Los algoritmos de raíz cuadrada calculan la raíz cuadrada no negativa.S{\displaystyle {\sqrt {S}}}de un número real positivoS{\displaystyle S}. Dado que todas las raíces cuadradas de números naturales , excepto las de cuadrados perfectos , son irracionales , [ 1 ] las raíces cuadradas generalmente solo se pueden calcular con una precisión finita: estos algoritmos suelen construir una serie de aproximaciones cada vez más precisas .

La mayoría de los métodos de cálculo de la raíz cuadrada son iterativos: después de elegir una estimación inicial adecuada deS{\displaystyle {\sqrt {S}}}Se realiza un refinamiento iterativo hasta que se cumple algún criterio de terminación. Un método de refinamiento es el de Herón , un caso particular del método de Newton . Si la división es mucho más costosa que la multiplicación, puede ser preferible calcular la raíz cuadrada inversa .

Existen otros métodos para calcular la raíz cuadrada dígito a dígito o mediante series de Taylor . Las aproximaciones racionales de las raíces cuadradas pueden calcularse mediante expansiones en fracciones continuas .

El método empleado depende de la precisión requerida, así como de las herramientas y la capacidad de cálculo disponibles. Los métodos pueden clasificarse, a grandes rasgos, en aquellos adecuados para el cálculo mental, aquellos que generalmente requieren al menos papel y lápiz, y aquellos que se implementan como programas para ejecutarse en una computadora electrónica digital u otro dispositivo informático. Los algoritmos pueden tener en cuenta la convergencia (cuántas iteraciones se requieren para lograr una precisión específica), la complejidad computacional de las operaciones individuales (por ejemplo, la división) o de las iteraciones, y la propagación de errores (la exactitud del resultado final).

Algunos métodos, como la división sintética con lápiz y papel y el desarrollo en serie, no requieren un valor inicial. En ciertas aplicaciones, se requiere una raíz cuadrada entera , que es la raíz cuadrada redondeada o truncada al entero más cercano (en este caso, se puede emplear un procedimiento modificado).

Historia

Los procedimientos para hallar raíces cuadradas (en particular la raíz cuadrada de 2 ) se conocen al menos desde el período de la antigua Babilonia en el siglo XVII a. C. Los matemáticos babilonios calculaban la raíz cuadrada de 2 con tres dígitos sexagesimales después del 1, pero no se sabe con exactitud cómo. Sabían cómo aproximar una hipotenusa usando a2+b2a+b22a{\displaystyle {\sqrt {a^{2}+b^{2}}}\approx a+{\frac {b^{2}}{2a}}} (dando por ejemplo4160+153600{\displaystyle {\frac {41}{60}}+{\frac {15}{3600}}}para la diagonal de una puerta cuya altura es4060{\displaystyle {\frac {40}{60}}}varillas y cuyo ancho es1060{\displaystyle {\frac {10}{60}}}varillas) y es posible que hayan utilizado un enfoque similar para encontrar la aproximación de2.{\displaystyle {\sqrt {2}}.}[ 2 ]

El método de Herón del siglo I en Egipto fue el primer algoritmo verificable para calcular la raíz cuadrada. [ 3 ]

Los métodos analíticos modernos comenzaron a desarrollarse tras la introducción del sistema de numeración arábiga en Europa occidental a principios del Renacimiento. [ 4 ]

Hoy en día, casi todos los dispositivos informáticos disponen de una función de raíz cuadrada rápida y precisa, ya sea como una construcción del lenguaje de programación , una función intrínseca del compilador o de una biblioteca, o como un operador de hardware, basado en uno de los procedimientos descritos.

Estimación inicial

Muchos algoritmos iterativos de raíz cuadrada requieren un valor semilla inicial . La semilla debe ser un número positivo distinto de cero; debe estar entre 1 yS{\displaystyle S}, el número cuya raíz cuadrada se desea, porque la raíz cuadrada debe estar en ese rango. Si la semilla está lejos de la raíz, el algoritmo requerirá más iteraciones. Si se inicializa conincógnita0=1{\displaystyle x_{0}=1}(oS{\displaystyle S}), entonces aproximadamente12|registro2S|{\displaystyle {\tfrac {1}{2}}\vert \log _{2}S\vert }Se desperdiciarán iteraciones solo para obtener el orden de magnitud de la raíz. Por lo tanto, es útil tener una estimación aproximada, que puede tener una precisión limitada pero es fácil de calcular. En general, cuanto mejor sea la estimación inicial, más rápida será la convergencia. Para el método de Newton, una semilla ligeramente mayor que la raíz convergerá un poco más rápido que una semilla ligeramente menor que la raíz.

En general, una estimación se realiza de acuerdo con un intervalo arbitrario que se sabe que contiene la raíz (como por ejemplo[incógnita0,S/incógnita0]{\displaystyle [x_{0},S/x_{0}]}). La estimación es un valor específico de una aproximación funcional aF(incógnita)=incógnita{\displaystyle f(x)={\sqrt {x}}}en el intervalo. Obtener una mejor estimación implica obtener límites más ajustados en el intervalo o encontrar una mejor aproximación funcional paraF(incógnita){\displaystyle f(x)}Esto último generalmente implica el uso de un polinomio de orden superior en la aproximación, aunque no todas las aproximaciones son polinómicas. Los métodos comunes de estimación incluyen escalar, lineal, hiperbólico y logarítmico. Generalmente se utiliza una base decimal para la estimación mental o con lápiz y papel. Una base binaria es más adecuada para estimaciones por computadora. En la estimación, el exponente y la mantisa se tratan por separado, como se expresaría el número en notación científica.

Estimaciones decimales

Normalmente el númeroS{\displaystyle S}se expresa en notación científica comoa×102norte{\displaystyle a\times 10^{2n}}dónde1a<100{\displaystyle 1\leq a<100}y n es un número entero, y el rango de posibles raíces cuadradas esa×10norte{\displaystyle {\sqrt {a}}\times 10^{n}}dónde1a<10{\displaystyle 1\leq {\sqrt {a}}<10}.

Estimaciones escalares

Los métodos escalares dividen el rango en intervalos, y la estimación en cada intervalo está representada por un único número escalar. Si el rango se considera como un único intervalo, la media aritmética (5,5) o la media geométrica (103.16{\displaystyle {\sqrt {10}}\approx 3.16}) veces10norte{\displaystyle 10^{n}}Son estimaciones plausibles. El error absoluto y relativo de estas variará. En general, una sola magnitud escalar será muy imprecisa. Las mejores estimaciones dividen el rango en dos o más intervalos, pero las estimaciones escalares tienen una precisión inherentemente baja.

Para dos intervalos, divididos geométricamente, la raíz cuadradaS=a×10norte{\displaystyle {\sqrt {S}}={\sqrt {a}}\times 10^{n}}puede estimarse como [ Nota 1 ]S{2×10nortesi a<10,6×10nortesi a10.{\displaystyle {\sqrt {S}}\approx {\begin{cases}2\times 10^{n}&{\text{si }}a<10,\\6\times 10^{n}&{\text{si }}a\geq 10.\end{cases}}}

Esta estimación tiene un error absoluto máximo de4×10norte{\displaystyle 4\times 10^{n}}ena=100{\displaystyle a=100}y error relativo máximo del 100% ena=1{\displaystyle a=1}.

Ejemplo

ParaS=125348{\displaystyle S=125348}considerado como12.5348×104{\displaystyle 12.5348\times 10^{4}}, la estimación esS6102=600{\displaystyle {\sqrt {S}}\approx 6\cdot 10^{2}=600}.

125348=354.0{\displaystyle {\sqrt {125348}}=354.0}, un error absoluto de 246 y un error relativo de casi el 70%.

Estimaciones lineales

Una mejor estimación, y el método estándar utilizado, es una aproximación lineal a la función.y=incógnita2{\displaystyle y=x^{2}}sobre un pequeño arco. Si, como se indicó anteriormente, se factorizan las potencias de la base del número S y se reduce el intervalo a [ 1, 100 ] , se puede usar una línea secante que abarque el arco, o una línea tangente en algún punto a lo largo del arco como aproximación, pero una línea de regresión de mínimos cuadrados que interseque el arco será más precisa.

Una línea de regresión por mínimos cuadrados minimiza la diferencia promedio entre la estimación y el valor de la función. Su ecuación esy=8.7incógnita10{\displaystyle y=8.7x-10}. Reordenando,incógnita=0,115y+1.15{\displaystyle x=0.115y+1.15}. Redondeando los coeficientes para facilitar el cálculo, S(a/10+1.2)10norte{\displaystyle {\sqrt {S}}\approx (a/10+1.2)\cdot 10^{n}}

Esa es la mejor estimación promedio que se puede lograr con una aproximación lineal de una sola pieza de la función.y=incógnita2{\displaystyle y=x^{2}}en el intervalo [ 1, 100 ] . Tiene un error absoluto máximo de 1,2 en a = 100 y un error relativo máximo del 30 % en S = 1 y 10. [ Nota 2 ]

Para dividir por 10, se resta uno al exponente de a , o figurativamente se mueve la coma decimal un dígito a la izquierda. Para esta formulación, cualquier constante aditiva 1 más un pequeño incremento dará una estimación satisfactoria, por lo que recordar el número exacto no es una carga. La aproximación (redondeada o no) usando una sola línea que abarca el rango [ 1, 100 ] tiene menos de un dígito significativo de precisión; el error relativo es mayor que 1/2 2 , por lo que se proporcionan menos de 2 bits de información. La precisión está severamente limitada porque el rango es de dos órdenes de magnitud, bastante grande para este tipo de estimación.

Se puede obtener una estimación mucho mejor mediante una aproximación lineal por partes: múltiples segmentos de línea, cada uno aproximando algún subarco del original. Cuantos más segmentos de línea se utilicen, mejor será la aproximación. La forma más común es usar líneas tangentes; las decisiones críticas son cómo dividir el arco y dónde colocar los puntos de tangencia. Una forma eficaz de dividir el arco de y = 1 a y = 100 es geométricamente: para dos intervalos, los límites de los intervalos son la raíz cuadrada de los límites del intervalo original, 1 × 100, es decir [1, 2 100 ] y [ 2 100 ,100] . Para tres intervalos, los límites son las raíces cúbicas de 100: [1, 3 100 ], [ 3 100 ,( 3 100 ) 2 ] , y [( 3 100 ) 2 ,100] , etc. Para dos intervalos, 2 100 = 10 , un número muy conveniente. Las líneas tangentes son fáciles de derivar y se encuentran en incógnita=110{\displaystyle x={\sqrt {1{\sqrt {10}}}}}yincógnita=1010{\displaystyle x={\sqrt {10{\sqrt {10}}}}}Sus ecuaciones son:y=3.56incógnita3.16{\displaystyle y=3.56x-3.16} yy=11.2incógnita31.6{\displaystyle y=11.2x-31.6}Al invertir, las raíces cuadradas son:incógnita=0,28y+0,89{\displaystyle x=0.28y+0.89}yincógnita=0,089y+2.8{\displaystyle x=.089y+2.8}. Por lo tanto, paraS=a102norte{\displaystyle S=a\cdot 10^{2n}}: S{(0,28a+0,89)10nortesi a<10,(0,089a+2.8)10nortesi a10.{\displaystyle {\sqrt {S}}\approx {\begin{cases}(0.28a+0.89)\cdot 10^{n}&{\text{si }}a<10,\\(.089a+2.8)\cdot 10^{n}&{\text{si }}a\geq 10.\end{cases}}}

Los errores absolutos máximos se producen en los puntos más altos de los intervalos, en a = 10 y 100, y son 0,54 y 1,7 respectivamente. Los errores relativos máximos se encuentran en los extremos de los intervalos, en a = 1, 10 y 100, y son del 17 % en ambos casos. El 17 % o 0,17 es mayor que 1/10, por lo que el método arroja una precisión inferior a una cifra decimal.

Estimaciones hiperbólicas

En algunos casos, las estimaciones hiperbólicas pueden ser eficaces, porque una hipérbola también es una curva convexa y puede estar mejor situada a lo largo de un arco de y = que una línea. Las estimaciones hiperbólicas son computacionalmente más complejas, porque necesariamente requieren una división flotante. Una aproximación hiperbólica casi óptima a en el intervalo [ 1 , 100 ] esy=190/(10incógnita)20{\displaystyle y=190/(10-x)-20} . Transponiendo, la raíz cuadrada esincógnita=10190/(y+20){\displaystyle x=10-190/(y+20)} . Por lo tanto, paraS=a102norte{\displaystyle S=a\cdot 10^{2n}}: S(10190a+20)10norte{\displaystyle {\sqrt {S}}\approx \left(10-{\frac {190}{a+20}}\right)\cdot 10^{n}}

La división debe ser precisa solo hasta un dígito decimal, ya que la estimación general solo tiene esa precisión y puede hacerse mentalmente. Esta estimación hiperbólica es mejor en promedio que las estimaciones escalares o lineales. Tiene un error absoluto máximo de 1,58 en a = 100 y un error relativo máximo en a = 10 , donde la estimación de 3,67 es un 16,0 % mayor que la raíz de 3,16. Si en cambio se realizaran iteraciones de Newton-Raphson comenzando con una estimación de 10, se necesitarían dos iteraciones para llegar a 3,66, que coincide con la estimación hiperbólica. Para un caso más típico como 75, la estimación hiperbólica de 8,00 es solo un 7,6 % menor, y se requerirían 5 iteraciones de Newton-Raphson comenzando en 75 para obtener un resultado más preciso.

Estimaciones aritméticas

Un método análogo a la aproximación lineal por partes, pero que utiliza únicamente aritmética en lugar de ecuaciones algebraicas, emplea las tablas de multiplicar en orden inverso: la raíz cuadrada de un número entre 1 y 100 se encuentra entre 1 y 10. Por lo tanto, si sabemos que 25 es un cuadrado perfecto (5 × 5) y 36 es un cuadrado perfecto (6 × 6), entonces la raíz cuadrada de un número mayor o igual a 25 pero menor que 36 comienza con un 5. Lo mismo ocurre con los números entre otros cuadrados. Este método proporciona el primer dígito correcto, pero no es preciso al número entero más cercano: el primer dígito de la raíz cuadrada de 35, por ejemplo, es 5, pero la raíz cuadrada de 35 es casi 6.

Una mejor manera es dividir el rango en intervalos a la mitad entre los cuadrados. Así, cualquier número entre 25 y la mitad de 36, que es 30.5, se estima en 5; cualquier número mayor que 30.5 hasta 36, ​​se estima en 6. [ Nota 3 ] El procedimiento solo requiere un poco de aritmética para encontrar un número límite en el medio de dos productos de la tabla de multiplicar. Aquí hay una tabla de referencia de esos límites:

La operación final consiste en multiplicar la estimación k por la potencia de diez dividida por 2, de modo que paraS=a102norte{\displaystyle S=a\cdot 10^{2n}}, Sk10norte{\displaystyle {\sqrt {S}}\approx k\cdot 10^{n}}

El método proporciona implícitamente una cifra significativa de precisión, ya que redondea al mejor primer dígito.

El método se puede extender a 3 cifras significativas en la mayoría de los casos, interpolando entre los cuadrados más cercanos que delimitan el operando. Sik2a<(k+1)2{\displaystyle k^{2}\leq a<(k+1)^{2}}, entoncesa{\displaystyle {\sqrt {a}}}es aproximadamente k más una fracción, la diferencia entre a y k 2 dividida por la diferencia entre los dos cuadrados:

ak+R{\displaystyle {\sqrt {a}}\approx k+R}dóndeR=ak2(k+1)2k2=ak22k+1{\displaystyle R={\frac {ak^{2}}{(k+1)^{2}-k^{2}}}={\frac {ak^{2}}{2k+1}}}

La operación final, como se indicó anteriormente, consiste en multiplicar el resultado por la potencia de diez dividido entre 2; S=a10norte(k+R)10norte{\displaystyle {\sqrt {S}}={\sqrt {a}}\cdot 10^{n}\approx (k+R)\cdot 10^{n}}

k es un dígito decimal y R es una fracción que debe convertirse a decimal. Generalmente tiene un solo dígito en el numerador y uno o dos dígitos en el denominador, por lo que la conversión a decimal se puede realizar mentalmente.

Ejemplo

Encuentra la raíz cuadrada de 75.

75=751020{\displaystyle 75=75\cdot 10^{2\cdot 0}} , entonces a es 75 y n es 0. De las tablas de multiplicar, la raíz cuadrada de la mantisa debe ser 8 coma algo porque a está entre 8 × 8 = 64 y 9 × 9 = 81, entonces k es 8; algo es la representación decimal de R . En la fracción R , el numerador es 75 k 2 = 11 , y el denominador es 81 k 2 = 2 k + 1 = 17 . 11/17 es un poco menor que 12/18 = 2/3 = 0.67, así que adivinamos 0.66 (está bien adivinar aquí, el error es muy pequeño). La estimación final es 8 + 0.66 = 8.66 .

√75 redondeado a tres cifras significativas es 8,66, por lo que la estimación es precisa hasta tres cifras significativas. No todas las estimaciones realizadas con este método serán tan exactas, pero se aproximarán bastante.

Estimaciones binarias

Cuando se trabaja en el sistema numérico binario (como lo hacen internamente las computadoras), al expresar S comoa×22norte{\displaystyle a\times 2^{2n}}dónde0.12a<102{\displaystyle 0.1_{2}\leq a<10_{2}}, la raíz cuadradaS=a×2norte{\displaystyle {\sqrt {S}}={\sqrt {a}}\times 2^{n}}puede estimarse como S(0,485+0,485a)2norte{\displaystyle {\sqrt {S}}\approx (0.485+0.485a)\cdot 2^{n}}

que es la línea de regresión de mínimos cuadrados con coeficientes de 3 cifras significativas.a{\displaystyle {\sqrt {a}}}tiene un error absoluto máximo de 0,0408 ena=2{\displaystyle a=2}y un error relativo máximo del 3,0% ena=1{\displaystyle a=1}Una estimación redondeada computacionalmente conveniente (porque los coeficientes son potencias de 2) es:

S(0,5+0,5a)2norte{\displaystyle {\sqrt {S}}\approx (0.5+0.5a)\cdot 2^{n}}[ Nota 4 ]

que tiene un error absoluto máximo de 0,086 en 2 y un error relativo máximo del 6,1% en a = 0,5 y a = 2,0 .

ParaS=125348=111101001101001002=1.11101001101001002×216,{\displaystyle S=125348=1\;1110\;1001\;1010\;0100_{2}=1.1110\;1001\;1010\;0100_{2}\times 2^{16}\,,}la aproximación binaria daS(0,5+0,5a)28=1.01110100110100102×1000000002=1.456×256=372.8.{\displaystyle {\sqrt {S}}\approx (0.5+0.5a)\cdot 2^{8}=1.0111\;0100\;1101\;0010_{2}\times 1\;0000\;0000_{2}=1.456\times 256=372.8.}125348=354.0{\displaystyle {\sqrt {125348}}=354.0}Por lo tanto, la estimación tiene un error absoluto de 19 y un error relativo de 5,3%. El error relativo es un poco menor que 1/2 4 , por lo que la estimación es buena hasta 4+ bits.

Una estimación para un buen valor de 8 bits se puede obtener mediante una búsqueda en la tabla de los 8 bits superiores de un , recordando que el bit superior es implícito en la mayoría de las representaciones de punto flotante, y el bit inferior de los 8 debe redondearse. La tabla es de 256 bytes de valores de raíz cuadrada de 8 bits precalculados. Por ejemplo, para el índice 11101101 2 que representa 1.8515625 10 , la entrada es 10101110 2 que representa 1.359375 10 , la raíz cuadrada de 1.8515625 10 con una precisión de 8 bits (2+ dígitos decimales).

El método de Herón

El primer algoritmo explícito para aproximar S  {\displaystyle \ {\sqrt {S~}}\ }es conocido como el método de Herón , en honor al matemático griego del siglo I, Herón de Alejandría, quien describió el método en su obra Métrica del año 60 d . C. [ 3 ] Este método también se denomina método babilónico (que no debe confundirse con el método babilónico para aproximar hipotenusas ), aunque no hay evidencia de que el método fuera conocido por los babilonios .

Dado un número real positivoS{\displaystyle S}Sea x 0 > 0 cualquier estimación inicial positiva . El método de Herón consiste en calcular iterativamente incógnitanorte+1=12(incógnitanorte+Sincógnitanorte),{\displaystyle x_{n+1}={\frac {1}{2}}\left(x_{n}+{\frac {S}{x_{n}}}\right),} hasta que se alcance la precisión deseada. La secuencia ( incógnita0, incógnita1, incógnita2, incógnita3,  ) {\displaystyle \ {\bigl (}\ x_{0},\ x_{1},\ x_{2},\ x_{3},\ \ldots \ {\bigr )}\ }definido por esta ecuación converge a límitenorteincógnitanorte=S  .{\displaystyle \ \lim _{n\to \infty }x_{n}={\sqrt {S~}}~.}

Esto es equivalente a utilizar el método de Newton para resolverincógnita2S=0{\displaystyle x^{2}-S=0}Este algoritmo converge cuadráticamente : el número de dígitos correctos deincógnitanorte{\displaystyle x_{n}}aproximadamente se duplica con cada iteración. [ 5 ]

Derivación

La idea básica es que si incógnita {\displaystyle \ x\ }es una sobreestimación de la raíz cuadrada de un número real positivo S {\displaystyle \ S\ }entonces  S incógnita {\displaystyle \ {\tfrac {\ S\ }{x}}\ }será una subestimación, y viceversa, por lo que cabe esperar que el promedio de estos dos números proporcione una mejor aproximación. (La demostración formal de esta afirmación se basa en la desigualdad de las medias aritmética y geométrica , que muestra que este promedio siempre sobreestima la raíz cuadrada, como se indica en el artículo sobre raíces cuadradas , lo que garantiza la convergencia).

Más precisamente, si incógnita {\displaystyle \ x\ }es nuestra suposición inicial de S  {\displaystyle \ {\sqrt {S~}}\ }y ε {\displaystyle \ \varepsilon \ }es el error en nuestra estimación tal que S=(incógnita+ε)2 ,{\displaystyle \ S=\left(x+\varepsilon \right)^{2}\ ,}Entonces podemos expandir el binomio como:  ( incógnita+ε )2=incógnita2+2incógnitaε+ε2{\displaystyle \ {\bigl (}\ x+\varepsilon \ {\bigr )}^{2}=x^{2}+2x\varepsilon +\varepsilon ^{2}} y resolver para el término de error

ε= Sincógnita2  2incógnita+ε  Sincógnita2 2incógnita ,{\displaystyle \varepsilon ={\frac {\ Sx^{2}\ }{\ 2x+\varepsilon \ }}\approx {\frac {\ Sx^{2}\ }{2x}}\ ,}si suponemos que εincógnita {\displaystyle \ \varepsilon \ll x~}

Por lo tanto, podemos compensar el error y actualizar nuestra antigua estimación como  incógnita+ε  incógnita+ Sincógnita2 2incógnita =  S+incógnita2 2incógnita =  S incógnita +incógnita 2  incógnitarmivismid .{\displaystyle \ x+\varepsilon \ \approx \ x+{\frac {\ Sx^{2}\ }{2x}}\ =\ {\frac {\ S+x^{2}\ }{2x}}\ =\ {\frac {\ {\frac {S}{\ x\ }}+x\ }{2}}\ \equiv \ x_{\mathsf {revisado}}~.} Dado que el error calculado no fue exacto, esta no es la respuesta correcta, sino que se convierte en nuestra nueva estimación para la siguiente ronda de corrección. El proceso de actualización se repite hasta obtener la precisión deseada.

Este algoritmo funciona igual de bien en los números p -ádicos , pero no se puede utilizar para identificar raíces cuadradas reales con raíces cuadradas p -ádicas; por ejemplo, se puede construir una secuencia de números racionales mediante este método que converge a +3 en los reales, pero a −3 en los 2-ádicos.

Implementación en Python

from decimal import Decimal , localcontext , getcontextTipoNúmero = entero | flotante | decimaldef sqrt_Heron (s : TipoNúmero ,precisión : int | Ninguno = Ninguno ,adivinar : NumberType | Ninguno = Ninguno) -> Decimal :""" Calcula sqrt(s) utilizando el método de Heron-Newton con precisión arbitraria. :param s: Número no negativo cuya raíz cuadrada se va a calcular. :param precision: Número de dígitos significativos. Utiliza por defecto el contexto decimal actual. La precisión mínima admitida es 2. (No se permite una precisión de 1 para evitar anomalías de redondeo). :param guess: Estimación inicial. Por defecto es s / 2. :return: Aproximación de sqrt(s) redondeada a la precisión especificada. """si s == 0 :devolver Decimal ( 0 )s = Decimal ( s )si s < 0 :Generar ValueError ( "sqrt(s) no está definido para números negativos." )Si precision es None :precisión = obtenercontexto () . prec # usar el contexto global actual si no se especifica# Imponer silenciosamente precisión mínimasi la precisión < 2 :precisión = 2Si guess es None :Suposición = Decimal ( s / 2 )guardia = 25 # dígitos adicionales temporales para estabilidad internamax_iter = 10_000# Contexto local: aislar los cambios de precisióncon localcontext () como ctx :ctx . prec = precisión + guardiaSuposición = ( suposición + s / suposición ) / 2para _ en rango ( iterador_máximo ):siguiente_suposición = ( suposición + s / suposición ) / 2# Detente cuando la mejora sea lo suficientemente pequeñaSi guess - next_guess < Decimal ( f "1e- { precisión } " ):romperadivinanza = siguiente_adivinanzademás :generar ArithmeticError ( f "El método Heron no convergió en { max_iter } iteraciones" )# Redondeo a precisión de objetivo (eliminando la protección)ctx.prec = precisióndevolver + siguiente_suposición

Ejemplo de cálculo

El siguiente ejemplo demuestra la ejecución de la sqrt_Heronfunción con diferentes entradas.

print ( f "1) { sqrt_Heron ( 125348 , precision = 7 , guess = 600 ) } " ) print ( f "2) { sqrt_Heron ( Decimal ( '3.1415926535897932384626433832795028841971693993' )) } " ) print ( f "3) { sqrt_Heron ( 2 , 1_000_157 ) } " ) print ( f "4) { sqrt_Heron ( 2 , 10_000_005 , 1.414 ) } " ) print ( f "5) { sqrt_Heron ( 2 , 100_000_000 , 1 ) } " )

Esto produce el siguiente resultado:

1) 354.0452 2) 1.772453850905516027298167483 3) 1.4142135623730950488016887242 ... 269732025731849141493880004856742892 4) 1.4142135623730950488016887242 ... ... 872480508054123572727872131589714262 5) 1.4142135623730950488016887242 ... ... ... 023678977744844723443287604232894971 

Anuncio 1)

El cálculo des{\displaystyle {\sqrt {s\,}}}paras=125348{\displaystyle s=125348}Para siete cifras significativas se sigue el siguiente curso: incógnita0=6102=600incógnita1=12(incógnita0+Sincógnita0)=12(600456666666+1253486004566666)404.456666666incógnita2=12(incógnita1+Sincógnita1)=12(404.456666666+125348404.456666666)357.186837334incógnita3=12(incógnita2+Sincógnita2)=12(357.186837334+125348357.186837334)354.059011038incógnita4=12(incógnita3+Sincógnita3)=12(354.059011038+125348354.059011038)354.045195124incógnita5=12(incógnita4+Sincógnita4)=12(354.045195124+125348354.045195124)354.045194855{\displaystyle {\begin{alignedat}{5}x_{0}&=6\cdot 10^{2}=600\\x_{1}&={\frac {1}{2}}\left(x_{0}+{\frac {S}{x_{0}}}\right)&&={\frac {1}{2}}\left(600{\phantom {456666666}}+{\frac {125348}{600}}{\phantom {4566666}}\right)&&\approx 404.456666666\\x_{2}&={\frac {1}{2}}\left(x_{1}+{\frac {S}{x_{1}}}\right)&&={\frac {1}{2}}\left(404.456666666+{\frac {125348}{404.456666666}}\right)&&\approx 357.186837334\\x_{3}&={\frac {1}{2}}\left(x_{2}+{\frac {S}{x_{2}}}\right)&&={\frac {1}{2}}\left(357.186837334+{\frac {125348}{357.186837334}}\right)&&\approx 354.059011038\\x_{4}&={\frac {1}{2}}\left(x_{3}+{\frac {S}{x_{3}}}\right)&&={\frac {1}{2}}\left(354.059011038+{\frac {125348}{354.059011038}}\right)&&\approx 354.045195124\\x_{5}&={\frac {1}{2}}\left(x_{4}+{\frac {S}{x_{4}}}\right)&&={\frac {1}{2}}\left(354.045195124+{\frac {125348}{354.045195124}}\right)&&\approx 354.045194855\end{alignedat}}}

Por lo tanto125348354.0452{\displaystyle {\sqrt {\,125348\,}}\approx 354.0452}a siete cifras significativas (redondeadas).

Ad 2) Cálculo (en 6 pasos de iteración) deπ×1046×1046{\displaystyle {\sqrt {\left\lfloor \pi \times 10^{46}\right\rfloor \times 10^{-46}}}}a precisión predeterminada. [ Nota 5 ]

Ad 3) Cálculo (en 22 pasos de iteración) de2{\displaystyle {\sqrt {2}}}hasta 1.000.157 dígitos. [ 6 ]

Anuncio 4) Cálculo (en 23 pasos de iteración) de2{\displaystyle {\sqrt {2}}}hasta 10.000.005 dígitos. [ 7 ]

Ad 5) Cálculo (en 28 pasos de iteración) de2{\displaystyle {\sqrt {2}}}hasta 100 millones de dígitos.

Parece ser que, para obtener estimaciones iniciales razonables, no se necesitan muchas iteraciones.

Notas

Explicación de las líneas 41 y 46

El método de Herón tiene la siguiente propiedad:

incógnitanorte<Sincógnitanorte+1>Sincógnitanorte=Sincógnitanorte+1=incógnitanorteincógnitanorte>Sincógnitanorte+1<incógnitanorte{\displaystyle {\begin{array}{rcl}x_{n}<{\sqrt {S}}&\implies &x_{n+1}>{\sqrt {S}}\\x_{n}={\sqrt {S}}&\implies &x_{n+1}=x_{n}\\x_{n}>{\sqrt {S}}&\implies &x_{n+1}<x_{n}\end{array}}}

En palabras sencillas: Una vez que la iteración produce un valor mayor queS{\displaystyle {\sqrt {S}}}(lo cual sucede inmediatamente siincógnita0>S{\displaystyle x_{0}>{\sqrt {S}}}, o después de un paso siincógnita0<S{\displaystyle x_{0}<{\sqrt {S}}}), cada estimación siguiente permanece por encimaS{\displaystyle {\sqrt {S}}}pero se hace más pequeño cada vez, por lo que la secuencia “se desliza hacia abajo” haciaS{\displaystyle {\sqrt {S}}}y converge.

En la línea 41 del programa, guessse establece un valor.S{\displaystyle \geq {\sqrt {S}}}. Luego, en la línea 46 del código,δnorte=incógnitanorteincógnitanorte+1{\displaystyle \delta _{n}=x_{n}-x_{n+1}}no puede ser negativo.

Justificación del criterio de parada

Utilizando la diferencia entre estimaciones sucesivas,

δnorte=incógnitanorteincógnitanorte+1{\displaystyle \delta _{n}=x_{n}-x_{n+1}},

como criterio de parada, el método asegura que la secuencia de aproximacionesincógnitanorte{\displaystyle x_{n}}está convergiendo hacia el valor verdaderoS{\displaystyle {\sqrt {S}}}Cuando las diferencias sucesivasδnorte{\displaystyle \delta _{n}}Cuando se vuelven suficientemente pequeños, se alcanza el objetivo dado. La idea clave es que el error absoluto

εnorte=incógnitanorteS{\displaystyle \varepsilon _{n}=x_{n}-{\sqrt {S}}},

está directamente relacionado con el tamaño de la mejora sucesivaδnorte{\displaystyle \delta _{n}}Específicamente, para los métodos iterativos que convergen lineal o cuadráticamente, existe una constantedo<1{\displaystyle C<1}de tal manera que

εnorte+1=doδnorte{\displaystyle \varepsilon _{n+1}=C\cdot \delta _{n}}.

Esta relación implica que comoδnorte{\displaystyle \delta _{n}}disminuye, el error absolutoεnorte{\displaystyle \varepsilon _{n}}también se vuelve más pequeño. Por lo tanto, detener la iteración cuandoδnorte{\displaystyle \delta _{n}}Si el valor cae por debajo de un umbral determinado, se garantiza que el error real se encuentre como máximo dentro de ese umbral.

Convergencia

Gráficos semilogarítmicos que comparan la velocidad de convergencia del método de Herón para hallar la raíz cuadrada de 100 con diferentes valores iniciales. Los valores negativos convergen a la raíz negativa, y los positivos a la raíz positiva. Nótese que los valores más cercanos a la raíz convergen más rápido, y todas las aproximaciones son sobreestimaciones. En el archivo SVG, coloque el cursor sobre un gráfico para visualizar sus puntos.

Supongamos que incógnita0>0  anorted  S>0 .{\displaystyle \ x_{0}>0~~{\mathsf {and}}~~S>0~.}Entonces, para cualquier número natural norte:incógnitanorte>0 .{\displaystyle \ n:x_{n}>0~.}Dejemos el error relativo en incógnitanorte {\displaystyle \ x_{n}\ }ser definido por  εnorte= incógnitanorte  S  1>1 {\displaystyle \ \varepsilon _{n}={\frac {~x_{n}\ }{\ {\sqrt {S~}}\ }}-1>-1\ } y por lo tanto  incógnitanorte=S (1+εnorte) .{\displaystyle \ x_{n}={\sqrt {S~}}\cdot \left(1+\varepsilon _{n}\right)~.}

Entonces se puede demostrar que  εnorte+1=εnorte22(1+εnorte)0 .{\displaystyle \ \varepsilon _{n+1}={\frac {\varepsilon _{n}^{2}}{2(1+\varepsilon _{n})}}\geq 0~.}

Y así que  εnorte+2min{  εnorte+12 2, εnorte+1 2 } {\displaystyle \ \varepsilon _{n+2}\leq \min \left\{\ {\frac {\ \varepsilon _{n+1}^{2}\ }{2}},{\frac {\ \varepsilon _{n+1}\ }{2}}\ \right\}\ } y, en consecuencia, esa convergencia está asegurada y es cuadrática .

Peor escenario para la convergencia

Si se utiliza la estimación aproximada anterior con el método babilónico, los casos menos precisos en orden ascendente son los siguientes: S= 1 ;incógnita0= 2 ;incógnita1= 1.250 ;ε1= 0,250 .S= 10 ;incógnita0= 2 ;incógnita1= 3.500 ;ε1< 0,107 .S= 10 ;incógnita0= 6 ;incógnita1= 3.833 ;ε1< 0,213 .S= 100 ;incógnita0= 6 ;incógnita1= 11.333 ;ε1< 0,134 .{\displaystyle {\begin{aligned}S&=\ 1\ ;&x_{0}&=\ 2\  ;&x_{1}&=\ 1.250\  ;&\varepsilon _{1}&=\ 0.250~.\\S&=\ 10\  ;&x_{0}&=\ 2\  ;&x_{1}&=\ 3.500\  ;&\varepsilon _{1}&<\ 0.107~.\\S&=\ 10\  ;&x_{0}&=\ 6\  ;&x_{1}&=\ 3.833\  ;&\varepsilon _{1}&<\ 0.213~.\\S&=\ 100\  ;&x_{0}&=\ 6\  ;&x_{1}&=\ 11.333\  ;&\varepsilon _{1}&<\ 0,134~.\end{aligned}}}

Por lo tanto, en cualquier caso, ε122.ε2<25<101 .ε3<211<103 .ε4<223<106 .ε5<247<1014 .ε6<295<1028 .ε7<2191<1057 .ε8<2383<10115 .{\displaystyle {\begin{aligned}\varepsilon _{1}&\leq 2^{-2}.\\\varepsilon _{2}&<2^{-5}<10^{-1}~.\\\varepsilon _{3}&<2^{-11}<10^{-3}~.\\\varepsilon _{4}&<2^{-23}<10^{-6}~.\\\varepsilon _{5}&<2^{-47}<10^{-14}~.\\\varepsilon _{6}&<2^{-95}<10^{-28}~.\\\varepsilon _{7}&<2^{-191}<10^{-57}~.\\\varepsilon _{8}&<2^{-383}<10^{-115}~.\end{aligned}}}

Los errores de redondeo ralentizarán la convergencia. Se recomienda mantener al menos un dígito adicional más allá de la precisión deseada. incógnitanorte {\displaystyle \ x_{n}\ }se calcula para evitar errores de redondeo significativos .

El método de Halley

Cuando en el programa anterior las estimaciones las líneas resaltadas 41 y 43

Suposición = ( suposición + s / suposición ) / 2# ...siguiente_suposición = ( suposición + s / suposición ) / 2

son reemplazados por

adivinanza *= ( adivinanza * adivinanza + 3 * s ) / ( 3 * adivinanza * adivinanza + s )# ...siguiente_suposición = suposición * ( suposición * suposición + 3 * s ) / ( 3 * suposición * suposición + s )

la función sqrt_Heronse transforma en una implementación del método de Halley , donde

incógnitanorte+1=incógnitanorteincógnitanorte2+3S3incógnitanorte2+S{\displaystyle x_{n+1}=x_{n}\cdot {\frac {x_{n}^{2}+3S}{3x_{n}^{2}+S}}}

Como en el método de Halley las estimaciones no siempre se mueven en una dirección, la línea 46 se convertirá en

si abs ( guess - next_guess ) < Decimal ( f "1e- { precision } " ):

El método de Halley converge más rápido —la tasa de convergencia a la raíz es cúbica, mejor que la cuadrática— iteración por iteración, pero implica cinco multiplicaciones por iteración (considerando la división como tres multiplicaciones). Los cinco cálculos de ejemplo se completan en 4, 4, 14, 15 y 19 iteraciones, respectivamente. En cambio, el método de Herón solo requiere una división, es decir, tres multiplicaciones, por lo que resulta ligeramente mejor a largo plazo.

Método Bakhshali

Este método para hallar una aproximación a una raíz cuadrada se describió en un antiguo manuscrito indio , llamado manuscrito Bakhshali . Es algebraicamente equivalente a dos iteraciones del método de Herón y, por lo tanto, converge cuárticamente, lo que significa que el número de dígitos correctos de la aproximación se cuadruplica aproximadamente con cada iteración. [ 8 ] La presentación original, utilizando notación moderna, es la siguiente: Para calcularS{\displaystyle {\sqrt {S}}}, dejarincógnita02{\displaystyle x_{0}^{2}}ser la aproximación inicial aS{\displaystyle S}. Luego, itere sucesivamente de la siguiente manera: anorte=Sincógnitanorte22incógnitanorte,incógnitanorte+1=incógnitanorte+anorte,incógnitanorte+2=incógnitanorte+1anorte22incógnitanorte+1.{\displaystyle {\begin{aligned}a_{n}&={\frac {S-x_{n}^{2}}{2x_{n}}},\\x_{n+1}&=x_{n}+a_{n},\\x_{n+2}&=x_{n+1}-{\frac {a_{n}^{2}}{2x_{n+1}}}.\end{aligned}}}

Los valoresincógnitanorte+1{\displaystyle x_{n+1}}yincógnitanorte+2{\displaystyle x_{n+2}}son exactamente iguales a los calculados por el método de Herón. Para ver esto, el segundo paso del método de Herón calcularía incógnitanorte+2=incógnitanorte+12+S2incógnitanorte+1=incógnitanorte+1+Sincógnitanorte+122incógnitanorte+1{\displaystyle x_{n+2}={\frac {x_{n+1}^{2}+S}{2x_{n+1}}}=x_{n+1}+{\frac {S-x_{n+1}^{2}}{2x_{n+1}}}} y podemos utilizar las definiciones deincógnitanorte+1{\displaystyle x_{n+1}}yanorte{\displaystyle a_{n}}para reorganizar el numerador en: Sincógnitanorte+12=S(incógnitanorte+anorte)2=Sincógnitanorte22incógnitanorteanorteanorte2=Sincógnitanorte2(Sincógnitanorte2)anorte2=anorte2.{\displaystyle {\begin{aligned}S-x_{n+1}^{2}&=S-(x_{n}+a_{n})^{2}\\&=S-x_{n}^{2}-2x_{n}a_{n}-a_{n}^{2}\\&=S-x_{n}^{2}-(S-x_{n}^{2})-a_{n}^{2}\\&=-a_{n}^{2}.\end{aligned}}}

Esto se puede utilizar para construir una aproximación racional a la raíz cuadrada comenzando con un número entero. Siincógnita0=norte{\displaystyle x_{0}=N}es un número entero elegido de tal manera quenorte2{\displaystyle N^{2}}está cerca deS{\displaystyle S}, yd=Snorte2{\displaystyle d=SN^{2}}es la diferencia cuyo valor absoluto se minimiza, entonces la primera iteración se puede escribir como: Snorte+d2norted28norte3+4norted=8norte4+8norte2d+d28norte3+4norted=norte4+6norte2S+S24norte3+4norteS=norte2(norte2+6S)+S24norte(norte2+S).{\displaystyle {\sqrt {S}}\approx N+{\frac {d}{2N}}-{\frac {d^{2}}{8N^{3}+4Nd}}={\frac {8N^{4}+8N^{2}d+d^{2}}{8N^{3}+4Nd}}={\frac {N^{4}+6N^{2}S+S^{2}}{4N^{3}+4NS}}={\frac {N^{2}(N^{2}+6S)+S^{2}}{4N(N^{2}+S)}}.}

El método de Bakhshali se puede generalizar al cálculo de una raíz arbitraria, incluidas las raíces fraccionarias. [ 9 ]

Se podría pensar que la segunda mitad del método Bakhshali podría usarse como una forma más simple de la iteración de Heron y usarse repetidamente, por ejemplo anorte+1=anorte22incógnitanorte+1,incógnitanorte+2=incógnitanorte+1+anorte+1,anorte+2=anorte+122incógnitanorte+2,incógnitanorte+3=incógnitanorte+2+anorte+2, etc.{\displaystyle {\begin{aligned}a_{n+1}&={\frac {-a_{n}^{2}}{2x_{n+1}}},&x_{n+2}&=x_{n+1}+a_{n+1},\\a_{n+2}&={\frac {-a_{n+1}^{2}}{2x_{n+2}}},&x_{n+3}&=x_{n+2}+a_{n+2},{\text{ etc.}}\end{aligned}}} Sin embargo, esto es numéricamente inestable . Sin ninguna referencia al valor de entrada original.S{\displaystyle S}, la precisión está limitada por la del cálculo original deanorte{\displaystyle a_{n}}y eso rápidamente se vuelve insuficiente.

Ejemplo

Usando el mismo ejemploS=125348{\displaystyle S=125348}como en el ejemplo del método de Herón , la primera iteración da incógnita0=600a0=12534860022×600=195.5433200incógnita1=600+(200)=400incógnita2=400(200)22×400=350{\displaystyle {\begin{alignedat}{3}x_{0}&=600\\[1ex]a_{0}&={\frac {125348-600^{2}}{2\times 600}}&&=-195.5433\approx -200\\[1ex]x_{1}&=600+(-200)&&={\phantom {-}}400\\[1ex]x_{2}&=400-{\frac {(-200)^{2}}{2\times 400}}&&={\phantom {-}}350\end{alignedat}}}

De igual modo, la segunda iteración da como resultado a2=12534835022×350=004.06857incógnita3=350+4.06857=354.06857incógnita4=354.068574.0685722×354.06857=354.045194{\displaystyle {\begin{alignedat}{3}a_{2}&={\frac {125348-350^{2}}{2\times 350}}&&={\phantom {00}}4.06857\\[1ex]x_{3}&=350+4.06857&&=354.06857\\[1ex]x_{4}&=354.06857-{\frac {4.06857^{2}}{2\times 354.06857}}&&=354.045194\end{alignedat}}} A diferencia del método de Herón,incógnita3{\displaystyle x_{3}}debe calcularse hasta 8 dígitos porque la fórmula paraincógnita4{\displaystyle x_{4}}no corrige ningún error enincógnita3{\displaystyle x_{3}}.

Cálculo dígito por dígito

Esta técnica proviene del trabajo de François Viète , publicado alrededor de 1600. [ 10 ] , y se basa en el teorema del binomio y es esencialmente un algoritmo inverso para resolver(incógnita+y)2=incógnita2+2incógnitay+y2{\displaystyle (x+y)^{2}=x^{2}+2xy+y^{2}}Es más lento que el método babilónico, pero tiene varias ventajas:

  • Puede resultar más sencillo para los cálculos manuales.
  • Se sabe que cada dígito de la raíz encontrada es correcto, es decir, no es necesario cambiarlo posteriormente.
  • Si la raíz cuadrada tiene una expansión que termina, el algoritmo finaliza después de encontrar el último dígito. Por lo tanto, se puede utilizar para comprobar si un número entero dado es un cuadrado perfecto .
  • El algoritmo funciona para cualquier base y, naturalmente, la forma en que procede depende de la base elegida.

Las desventajas son:

  • Se vuelve inmanejable para las raíces más altas.
  • No tolera estimaciones inexactas ni subcálculos; tales errores provocan que todos los dígitos siguientes del resultado sean incorrectos, a diferencia del método de Newton , que autocorrige cualquier error de aproximación.
  • Si bien el cálculo dígito por dígito es suficientemente eficiente en teoría, resulta demasiado costoso para su implementación en software. Cada iteración implica números mayores, lo que requiere más memoria, pero solo avanza la respuesta en un dígito correcto. Por lo tanto, el algoritmo tarda más tiempo por cada dígito adicional.

Los huesos de Napier incluyen una ayuda para la ejecución de este algoritmo. El algoritmo de raíz enésima desplazada es una generalización de este método.

Principio básico

Primero, consideremos el caso de hallar la raíz cuadrada de un número S , es decir, el cuadrado de un número de dos dígitos en base 10, XY , donde X es el dígito de las decenas e Y es el dígito de las unidades. Específicamente: S=(10incógnita+Y)2=100incógnita2+20incógnitaY+Y2.{\displaystyle S=\left(10X+Y\right)^{2}=100X^{2}+20XY+Y^{2}.}La letra S constará de 3 o 4 dígitos decimales.

Ahora, para comenzar el algoritmo dígito por dígito, dividimos los dígitos de S en dos grupos de dos dígitos, comenzando desde la derecha. Esto significa que el primer grupo será de 1 o 2 dígitos. Luego determinamos el valor de X como el dígito más grande tal que X 2 sea menor o igual que el primer grupo. Luego calculamos la diferencia entre el primer grupo y X 2 y comenzamos la segunda iteración concatenando el segundo grupo a él. Esto es equivalente a restar100incógnita2{\displaystyle 100X^{2}}de S , y nos quedamos conS=20incógnitaY+Y2{\displaystyle S'=20XY+Y^{2}}Dividimos S' entre 10, luego lo dividimos entre 2X y conservamos la parte entera para intentar adivinar Y. Concatenamos 2X con la Y tentativa y la multiplicamos por Y. Si nuestra suposición es correcta, esto equivale a calcular:(10(2incógnita)+Y)Y=20incógnitaY+Y2=S,{\displaystyle (10(2X)+Y)Y=20XY+Y^{2}=S',}Por lo tanto, el resto, es decir, la diferencia entre S' y el resultado, es cero; si el resultado es mayor que S' , disminuimos nuestra estimación en 1 y volvemos a intentarlo hasta que el resto sea 0. Dado que este es un caso simple donde la respuesta es una raíz cuadrada perfecta XY , el algoritmo se detiene aquí.

La misma idea se puede extender a cualquier cálculo de raíz cuadrada arbitrario a continuación. Supongamos que podemos encontrar la raíz cuadrada de S expresándola como una suma de n números positivos tales que S=(a1+a2+a3++anorte)2.{\displaystyle S=\left(a_{1}+a_{2}+a_{3}+\dots +a_{n}\right)^{2}.}

Al aplicar repetidamente la identidad básica (incógnita+y)2=incógnita2+2incógnitay+y2,{\displaystyle (x+y)^{2}=x^{2}+2xy+y^{2},} El término del lado derecho se puede expandir como (a1+a2+a3++anorte)2=a12+2a1a2+a22+2(a1+a2)a3+a32++anorte12+2(i=1norte1ai)anorte+anorte2=a12+[2a1+a2]a2+[2(a1+a2)+a3]a3++[2(i=1norte1ai)+anorte]anorte.{\displaystyle {\begin{aligned}&(a_{1}+a_{2}+a_{3}+\dotsb +a_{n})^{2}\\=&\,a_{1}^{2}+2a_{1}a_{2}+a_{2}^{2}+2(a_{1}+a_{2})a_{3}+a_{3}^{2}+\dots +a_{n-1}^{2}+2\left(\sum _{i=1}^{n-1}a_{i}\right)a_{n}+a_{n}^{2}\\=&\,a_{1}^{2}+[2a_{1}+a_{2}]a_{2}+[2(a_{1}+a_{2})+a_{3}]a_{3}+\dots +\left[2\left(\sum _{i=1}^{n-1}a_{i}\right)+a_{n}\right]a_{n}.\end{aligned}}}

Esta expresión nos permite encontrar la raíz cuadrada adivinando secuencialmente los valores deai{\displaystyle a_{i}}s. Supongamos que los númerosa1,,ametro1{\displaystyle a_{1},\ldots ,a_{m-1}}Como ya se ha adivinado, el m -ésimo término del lado derecho de la suma anterior viene dado porYmetro=[2PAGmetro1+ametro]ametro,{\displaystyle Y_{m}=\left[2P_{m-1}+a_{m}\right]a_{m},}dóndePAGmetro1=i=1metro1ai{\textstyle P_{m-1}=\sum _{i=1}^{m-1}a_{i}}es la raíz cuadrada aproximada encontrada hasta ahora. Ahora cada nueva suposiciónametro{\displaystyle a_{m}}debe satisfacer la recursión incógnitametro=incógnitametro1Ymetro,{\displaystyle X_{m}=X_{m-1}-Y_{m},} dóndeincógnitametro{\displaystyle X_{m}}es la suma de todos los términos despuésYmetro{\displaystyle Y_{m}}, es decir, el resto, de tal manera queincógnitametro0{\displaystyle X_{m}\geq 0}a pesar de1metronorte,{\displaystyle 1\leq m\leq n,}con inicializaciónincógnita0=S.{\displaystyle X_{0}=S.}Cuandoincógnitanorte=0,{\displaystyle X_{n}=0,}Se ha encontrado la raíz cuadrada exacta; si no, entonces la suma de laai{\displaystyle a_{i}}s proporciona una aproximación adecuada de la raíz cuadrada, conincógnitanorte{\displaystyle X_{n}}siendo el error de aproximación.

Por ejemplo, en el sistema numérico decimal tenemos S=(a110norte1+a210norte2++anorte110+anorte)2,{\displaystyle S=\left(a_{1}\cdot 10^{n-1}+a_{2}\cdot 10^{n-2}+\cdots +a_{n-1}\cdot 10+a_{n}\right)^{2},} dónde10nortei{\displaystyle 10^{ni}}son marcadores de posición y los coeficientesai{0,1,2,,9}{\displaystyle a_{i}\in \{0,1,2,\ldots ,9\}}. En cualquier etapa m del cálculo de la raíz cuadrada, la raíz aproximada encontrada hasta el momento,PAGmetro1{\displaystyle P_{m-1}}y el término de sumatoriaYmetro{\displaystyle Y_{m}}son dados por PAGmetro1=i=1metro1ai10nortei=10nortemetro+1i=1metro1ai10metroi1,{\displaystyle P_{m-1}=\sum _{i=1}^{m-1}a_{i}\cdot 10^{ni}=10^{n-m+1}\sum _{i=1}^{m-1}a_{i}\cdot 10^{mi-1},}Ymetro=[2PAGmetro1+ametro10nortemetro]ametro10nortemetro=[20i=1metro1ai10metroi1+ametro]ametro102(nortemetro).{\displaystyle Y_{m}=\left[2P_{m-1}+a_{m}\cdot 10^{n-m}\right]a_{m}\cdot 10^{n-m}=\left[20\sum _{i=1}^{m-1}a_{i}\cdot 10^{m-i-1}+a_{m}\right]a_{m}\cdot 10^{2(n-m)}.}

Aquí, dado que el valor posicional deYmetro{\displaystyle Y_{m}}es una potencia par de 10, solo necesitamos trabajar con el par de dígitos más significativos del restoincógnitametro1{\displaystyle X_{m-1}}, cuyo primer mandato esYmetro{\displaystyle Y_{m}}, en cualquier etapa m. La sección siguiente codifica este procedimiento.

Es obvio que se puede utilizar un método similar para calcular la raíz cuadrada en sistemas numéricos distintos del sistema decimal. Por ejemplo, encontrar la raíz cuadrada dígito por dígito en el sistema binario es bastante eficiente ya que el valor deai{\displaystyle a_{i}}se busca en un conjunto más pequeño de dígitos binarios {0,1}. Esto hace que el cálculo sea más rápido ya que en cada etapa el valor deYmetro{\displaystyle Y_{m}}es oYmetro=0{\displaystyle Y_{m}=0}paraametro=0{\displaystyle a_{m}=0}oYmetro=2PAGmetro1+1{\displaystyle Y_{m}=2P_{m-1}+1}paraametro=1{\displaystyle a_{m}=1}. El hecho de que solo tengamos dos opciones posibles paraametro{\displaystyle a_{m}}también hace que el proceso de decidir el valor deametro{\displaystyle a_{m}}en la m -ésima etapa del cálculo es más fácil. Esto se debe a que solo necesitamos comprobar siYmetroincógnitametro1{\displaystyle Y_{m}\leq X_{m-1}}paraametro=1.{\displaystyle a_{m}=1.}Si se cumple esta condición, entonces tomamosametro=1{\displaystyle a_{m}=1}; si no entoncesametro=0.{\displaystyle a_{m}=0.}Además, el hecho de que la multiplicación por 2 se realice mediante desplazamientos de bits a la izquierda facilita el cálculo.

Decimal (base 10)

Escribe el número original en forma decimal. Los números se escriben de forma similar al algoritmo de la división larga , y, como en la división larga, la raíz se escribirá en la línea superior. Ahora separa los dígitos en pares, comenzando desde el punto decimal y avanzando hacia la izquierda y hacia la derecha. El punto decimal de la raíz estará encima del punto decimal del cuadrado. Un dígito de la raíz aparecerá encima de cada par de dígitos del cuadrado.

Comenzando con el par de dígitos más a la izquierda, realice el siguiente procedimiento para cada par:

  1. Comenzando por la izquierda, baje el par de dígitos más significativos (el de más a la izquierda) que aún no se hayan usado (si ya se han usado todos los dígitos, escriba "00") y escríbalos a la derecha del resto del paso anterior (en el primer paso, no habrá resto). En otras palabras, multiplique el resto por 100 y sume los dos dígitos. Este será el valor actual c .
  2. Halla p , y y x , de la siguiente manera:
    • Sea p la parte de la raíz encontrada hasta ahora , sin tener en cuenta la coma decimal. (Para el primer paso, p = 0).
    • Determina el mayor dígito x tal queincógnita(20pag+incógnita)do{\displaystyle x(20p+x)\leq c}. Usaremos una nueva variable y = x (20 p + x ).
      • Nota: 20p + x es simplemente el doble de p , con el dígito x añadido a la derecha.
      • Nota: x se puede encontrar adivinando qué es c /(20· p ) y haciendo un cálculo de prueba de y , luego ajustando x hacia arriba o hacia abajo según sea necesario.
    • Coloca el dígitoincógnita{\displaystyle x}como el siguiente dígito de la raíz, es decir, encima de los dos dígitos del cuadrado que acabas de bajar. Por lo tanto, la siguiente p será la antigua p multiplicada por 10 más x .
  3. Resta y de c para formar un nuevo resto.
  4. Si el resto es cero y no hay más dígitos que bajar, el algoritmo ha terminado. De lo contrario, vuelva al paso 1 para una nueva iteración.

Ejemplos

Encuentra la raíz cuadrada de 152.2756.

 1 2. 3 4 / \/ 01 52.27 56 01 1·1 ≤ 1 < 2·2 x = 1 01 y = x·x = 1·1 = 1 00 52 22·2 ≤ 52 < 23·3 x = 2 00 44 y = (20+x)·x = 22·2 = 44 08 27 243·3 ≤ 827 < 244·4 x = 3 07 29 y = (240+x)·x = 243·3 = 729 98 56 2464·4 ≤ 9856 < 2465·5 x = 4 98 56 y = (2460 + x)·x = 2464·4 = 9856 00 00 El algoritmo finaliza: Respuesta = 12.34

Sistema numérico binario (base 2)

Esta sección utiliza el formalismo de la sección de cálculo dígito por dígito anterior , con la ligera variación que dejamosnorte2=(anorte++a0)2{\displaystyle N^{2}=(a_{n}+\dotsb +a_{0})^{2}}, con cada unoametro=2metro{\displaystyle a_{m}=2^{m}}oametro=0{\displaystyle a_{m}=0}Iteramos todo2metro{\displaystyle 2^{m}}, de2norte{\displaystyle 2^{n}}hasta20{\displaystyle 2^{0}}y construir una solución aproximadaPAGmetro=anorte+anorte1++ametro{\displaystyle P_{m}=a_{n}+a_{n-1}+\ldots +a_{m}}, la suma de todosai{\displaystyle a_{i}}para el cual hemos determinado el valor. Para determinar siametro{\displaystyle a_{m}}igual2metro{\displaystyle 2^{m}}o0{\displaystyle 0}, dejamosPAGmetro=PAGmetro+1+2metro{\displaystyle P_{m}=P_{m+1}+2^{m}}. SiPAGmetro2norte2{\displaystyle P_{m}^{2}\leq N^{2}}(es decir, el cuadrado de nuestra solución aproximada incluyendo2metro{\displaystyle 2^{m}}si no excede el cuadrado objetivo) entoncesametro=2metro{\displaystyle a_{m}=2^{m}}, de lo contrarioametro=0{\displaystyle a_{m}=0}yPAGmetro=PAGmetro+1{\displaystyle P_{m}=P_{m+1}}Para evitar la cuadraturaPAGmetro{\displaystyle P_{m}}En cada paso, almacenamos la diferencia.incógnitametro=norte2PAGmetro2{\displaystyle X_{m}=N^{2}-P_{m}^{2}}y actualizarlo incrementalmente mediante la configuraciónincógnitametro=incógnitametro+1Ymetro{\displaystyle X_{m}=X_{m+1}-Y_{m}}conYmetro=PAGmetro2PAGmetro+12=2PAGmetro+1ametro+ametro2{\displaystyle Y_{m}=P_{m}^{2}-P_{m+1}^{2}=2P_{m+1}a_{m}+a_{m}^{2}}. Inicialmente, establecimosanorte=PAGnorte=2norte{\displaystyle a_{n}=P_{n}=2^{n}}para el más grandenorte{\displaystyle n}con(2norte)2=4nortenorte2{\displaystyle (2^{n})^{2}=4^{n}\leq N^{2}}.

Como optimización adicional, almacenamosPAGmetro+12metro+1{\displaystyle P_{m+1}2^{m+1}}y(2metro)2{\displaystyle (2^{m})^{2}}, los dos términos deYmetro{\displaystyle Y_{m}}en caso de esoametro{\displaystyle a_{m}}es distinto de cero, en variables separadasdometro{\displaystyle c_{m}},dmetro{\displaystyle d_{m}}: dometro=PAGmetro+12metro+1{\displaystyle c_{m}=P_{m+1}2^{m+1}}dmetro=(2metro)2{\displaystyle d_{m}=(2^{m})^{2}}Ymetro={dometro+dmetrosi ametro=2metro0si ametro=0{\displaystyle Y_{m}={\begin{cases}c_{m}+d_{m}&{\text{if }}a_{m}=2^{m}\\0&{\text{if }}a_{m}=0\end{cases}}}

dometro{\displaystyle c_{m}}ydmetro{\displaystyle d_{m}}se puede actualizar de manera eficiente en cada paso: dometro1=PAGmetro2metro=(PAGmetro+1+ametro)2metro=PAGmetro+12metro+ametro2metro={dometro/2+dmetrosi ametro=2metrodometro/2si ametro=0{\displaystyle c_{m-1}=P_{m}2^{m}=(P_{m+1}+a_{m})2^{m}=P_{m+1}2^{m}+a_{m}2^{m}={\begin{cases}c_{m}/2+d_{m}&{\text{if }}a_{m}=2^{m}\\c_{m}/2&{\text{if }}a_{m}=0\end{cases}}}dmetro1=dmetro4{\displaystyle d_{m-1}={\frac {d_{m}}{4}}}

Tenga en cuenta que: do1=PAG020=PAG0=norte,{\displaystyle c_{-1}=P_{0}2^{0}=P_{0}=N,}que es el resultado final que devuelve la función siguiente.

Implementación

El programa Python calculaesqrt(norte)=norte.{\displaystyle \operatorname {isqrt} (n)=\lfloor {\sqrt {n}}\rfloor .}El algoritmo es un método dígito por dígito (bit por bit) para raíces cuadradas enteras . [ 11 ]

def isqrt ( x : int ) -> int : assert x >= 0 , "La entrada para la raíz cuadrada debe ser no negativa"op : int = x # X_(n+1) res : int = 0 # c_n# d_n que comienza en la mayor potencia de cuatro <= n uno : int = 1 mientras uno <= op : uno <<= 2 # Ahora 'uno' es la mayor potencia de cuatro <= x uno >>= 2# para dₙ … d₀ mientras uno != 0 : si op >= res + uno : # si X_(m+1) ≥ Y_m entonces a_m = 2^m op -= res + uno # X_m = X_(m+1) - Y_m res += 2 * uno # c_m = c_m + 2*d_m res //= 2 # c_(m-1) = c_m / 2 uno //= 4 # d_(m-1) = d_m / 4# c_(-1) devuelve res

Se pueden lograr algoritmos más rápidos, en binario y decimal o en cualquier otra base, utilizando tablas de búsqueda; en efecto, se intercambia más espacio de almacenamiento por un menor tiempo de ejecución . [ 12 ]

Identidad exponencial

Las calculadoras de bolsillo suelen implementar buenas rutinas para calcular la función exponencial y el logaritmo natural , y luego calcular la raíz cuadrada de S usando la identidad encontrada usando las propiedades de los logaritmos (lnincógnitanorte=nortelnincógnita{\displaystyle \ln x^{n}=n\ln x}) y exponenciales (milnincógnita=incógnita{\displaystyle e^{\ln x}=x}) : S=mi12lnS.{\displaystyle {\sqrt {S}}=e^{{\frac {1}{2}}\ln S}.} El denominador de la fracción corresponde a la raíz enésima . En el caso anterior, el denominador es 2, por lo que la ecuación indica que se debe hallar la raíz cuadrada. Esta misma identidad se utiliza al calcular raíces cuadradas con tablas de logaritmos o reglas de cálculo .

Un método iterativo de dos variables

Este método es aplicable para hallar la raíz cuadrada de0<S<3{\displaystyle 0<S<3\,\!}y converge mejor paraS1{\displaystyle S\approx 1}Sin embargo, esto no representa una limitación real para un cálculo informático, ya que en las representaciones de punto flotante y punto fijo en base 2, es trivial multiplicar.S{\displaystyle S\,\!}por una potencia entera de 4, y por lo tantoS{\displaystyle {\sqrt {S}}}mediante la potencia correspondiente de 2, cambiando el exponente o desplazando, respectivamente. Por lo tanto,S{\displaystyle S\,\!}se puede trasladar al rango12S<2{\textstyle {\tfrac {1}{2}}\leq S<2}Además, el método que se describe a continuación no emplea divisiones generales, sino únicamente sumas, restas, multiplicaciones y divisiones por potencias de dos, operaciones que, de nuevo, son muy fáciles de implementar. Una desventaja del método es que se acumulan errores numéricos, a diferencia de los métodos iterativos de una sola variable, como el babilónico.

El paso de inicialización de este método es a0=Sdo0=S1{\displaystyle {\begin{aligned}a_{0}&=S\\c_{0}&=S-1\end{aligned}}} mientras se leen los pasos iterativos anorte+1=anorteanortedonorte/2donorte+1=donorte2(donorte3)/4{\displaystyle {\begin{aligned}a_{n+1}&=a_{n}-a_{n}c_{n}/2\\c_{n+1}&=c_{n}^{2}(c_{n}-3)/4\end{aligned}}} Entonces,anorteS{\displaystyle a_{n}\to {\sqrt {S}}}(mientrasdonorte0{\displaystyle c_{n}\to 0}).

La convergencia dedonorte{\displaystyle c_{n}\,\!}y por lo tanto también deanorte{\displaystyle a_{n}\,\!}, es cuadrática.

La demostración del método es bastante sencilla. Primero, reescribimos la definición iterativa dedonorte{\displaystyle c_{n}}como 1+donorte+1=(1+donorte)(112donorte)2.{\displaystyle 1+c_{n+1}=(1+c_{n})(1-{\tfrac {1}{2}}c_{n})^{2}.} Entonces es sencillo demostrar por inducción que S(1+donorte)=anorte2{\displaystyle S(1+c_{n})=a_{n}^{2}} y por lo tanto la convergencia deanorte{\displaystyle a_{n}\,\!}al resultado deseadoS{\displaystyle {\sqrt {S}}}está garantizado por la convergencia dedonorte{\displaystyle c_{n}\,\!}a 0, lo cual a su vez se deduce de1<do0<2{\displaystyle -1<c_{0}<2\,\!}.

Este método fue desarrollado alrededor de 1950 por MV Wilkes , DJ Wheeler y S. Gill [ 13 ] para su uso en EDSAC , una de las primeras computadoras electrónicas. [ 14 ] Posteriormente, el método se generalizó, permitiendo el cálculo de raíces no cuadradas. [ 15 ]

Métodos iterativos para raíces cuadradas recíprocas

Los siguientes son métodos iterativos para encontrar la raíz cuadrada recíproca de S , que es1/S{\displaystyle 1/{\sqrt {S}}}Una vez encontrado, encontrarS{\displaystyle {\sqrt {S}}}mediante una simple multiplicación:S=S(1/S){\displaystyle {\sqrt {S}}=S\cdot (1/{\sqrt {S}})}Estas iteraciones solo implican multiplicación, no división. Por lo tanto, son más rápidas que el método babilónico . Sin embargo, no son estables. Si el valor inicial no está cerca de la raíz cuadrada recíproca, las iteraciones divergirán en lugar de converger. Por consiguiente, puede ser ventajoso realizar una iteración del método babilónico con una estimación aproximada antes de comenzar a aplicar estos métodos.

  • Aplicando el método de Newton a la ecuación(1/incógnita2)S=0{\displaystyle (1/x^{2})-S=0}produce un método que converge cuadráticamente utilizando tres multiplicaciones por paso:incógnitanorte+1=incógnitanorte2(3Sincógnitanorte2)=incógnitanorte(32S2incógnitanorte2).{\displaystyle x_{n+1}={\frac {x_{n}}{2}}\cdot (3-S\cdot x_{n}^{2})=x_{n}\cdot \left({\frac {3}{2}}-{\frac {S}{2}}\cdot x_{n}^{2}\right).}
  • Otra iteración se obtiene mediante el método de Halley , que es el método de Householder de segundo orden. Este converge cúbicamente , pero implica cinco multiplicaciones por iteración: [ Nota 6 ]ynorte=Sincógnitanorte2,{\displaystyle y_{n}=S\cdot x_{n}^{2},}yincógnitanorte+1=incógnitanorte8(15ynorte(103ynorte))=incógnitanorte(158ynorte(10838ynorte)).{\displaystyle x_{n+1}={\frac {x_{n}}{8}}\cdot (15-y_{n}\cdot (10-3\cdot y_{n}))=x_{n}\cdot \left({\frac {15}{8}}-y_{n}\cdot \left({\frac {10}{8}}-{\frac {3}{8}}\cdot y_{n}\right)\right).}
  • Si se realiza aritmética de punto fijo , la multiplicación por 3 y la división por 8 se pueden implementar usando desplazamientos y sumas. Si se utiliza punto flotante, el método de Halley se puede reducir a cuatro multiplicaciones por iteración mediante el preprocesamiento.3/8S{\textstyle {\sqrt {3/8}}S}y ajustando todas las demás constantes para compensar:ynorte=38Sincógnitanorte2,{\displaystyle y_{n}={\sqrt {\frac {3}{8}}}S\cdot x_{n}^{2},}yincógnitanorte+1=incógnitanorte(158ynorte(256ynorte)).{\displaystyle x_{n+1}=x_{n}\cdot \left({\frac {15}{8}}-y_{n}\cdot \left({\sqrt {\frac {25}{6}}}-y_{n}\right)\right).}

El algoritmo de Goldschmidt

El algoritmo de Goldschmidt es una extensión de la división de Goldschmidt , que recibe su nombre de Robert Elliot Goldschmidt, [ 16 ] [ 17 ] que se puede utilizar para calcular raíces cuadradas. Algunas computadoras utilizan el algoritmo de Goldschmidt para calcular simultáneamenteS{\displaystyle {\sqrt {S}}}y1/S{\displaystyle 1/{\sqrt {S}}}El algoritmo de Goldschmidt encuentraS{\displaystyle {\sqrt {S}}}más rápido que la iteración de Newton-Raphson en una computadora con una instrucción de multiplicación-suma fusionada y una unidad de punto flotante segmentada o dos unidades de punto flotante independientes. [ 18 ]

La primera forma de escribir el algoritmo de Goldschmidt comienza

b0=S{\displaystyle b_{0}=S}
Y01/S{\displaystyle Y_{0}\approx 1/{\sqrt {S}}}(normalmente mediante una búsqueda en tabla)
y0=Y0{\displaystyle y_{0}=Y_{0}}
incógnita0=Sy0{\displaystyle x_{0}=Sy_{0}}

y itera bnorte+1=bnorteYnorte2Ynorte+1=12(3bnorte+1)incógnitanorte+1=incógnitanorteYnorte+1ynorte+1=ynorteYnorte+1{\displaystyle {\begin{aligned}b_{n+1}&=b_{n}Y_{n}^{2}\\Y_{n+1}&={\tfrac {1}{2}}(3-b_{n+1})\\x_{n+1}&=x_{n}Y_{n+1}\\y_{n+1}&=y_{n}Y_{n+1}\end{aligned}}} hastabi{\displaystyle b_{i}}es suficientemente cercano a 1, o a un número fijo de iteraciones. Las iteraciones convergen a límitenorteincógnitanorte=S,{\displaystyle \lim _{n\to \infty }x_{n}={\sqrt {S}},}y límitenorteynorte=1/S.{\displaystyle \lim _{n\to \infty }y_{n}=1/{\sqrt {S}}.} Tenga en cuenta que es posible omitir cualquiera de los dos.incógnitanorte{\displaystyle x_{n}}yynorte{\displaystyle y_{n}}del cálculo, y si ambos son deseados entoncesincógnitanorte=Synorte{\displaystyle x_{n}=Sy_{n}}puede utilizarse al final en lugar de calcularlo en cada iteración.

Una segunda forma, que utiliza operaciones de multiplicación y suma fusionadas , comienza

y01/S{\displaystyle y_{0}\approx 1/{\sqrt {S}}}(normalmente mediante una búsqueda en tabla)
incógnita0=Sy0{\displaystyle x_{0}=Sy_{0}}
h0=12y0{\displaystyle h_{0}={\tfrac {1}{2}}y_{0}}

y itera rnorte=0,5incógnitanortehnorteincógnitanorte+1=incógnitanorte+incógnitanorternortehnorte+1=hnorte+hnorternorte{\displaystyle {\begin{aligned}r_{n}&=0.5-x_{n}h_{n}\\x_{n+1}&=x_{n}+x_{n}r_{n}\\h_{n+1}&=h_{n}+h_{n}r_{n}\end{aligned}}} hastari{\displaystyle r_{i}}está suficientemente cerca de 0, o de un número fijo de iteraciones. Esto converge a límitenorteincógnitanorte=S,{\displaystyle \lim _{n\to \infty }x_{n}={\sqrt {S}},}y límitenorte2hnorte=1/S.{\displaystyle \lim _{n\to \infty }2h_{n}=1/{\sqrt {S}}.}

Serie Taylor

Si N es una aproximación aS{\displaystyle {\sqrt {S}}}Se puede encontrar una mejor aproximación utilizando la serie de Taylor de la función raíz cuadrada :norte2+d=nortenorte=0(1)norte(2norte)¡(12norte)norte¡24nortednortenorte2norte=norte(1+d2norte2d28norte4+d316norte65d4128norte8+){\displaystyle {\sqrt {N^{2}+d}}=N\sum _{n=0}^{\infty }{\frac {(-1)^{n}(2n)!}{(1-2n)n!^{2}4^{n}}}{\frac {d^{n}}{N^{2n}}}=N\left(1+{\frac {d}{2N^{2}}}-{\frac {d^{2}}{8N^{4}}}+{\frac {d^{3}}{16N^{6}}}-{\frac {5d^{4}}{128N^{8}}}+\cdots \right)}

Como método iterativo , el orden de convergencia es igual al número de términos utilizados. Con dos términos, es idéntico al método babilónico . Con tres términos, cada iteración requiere casi tantas operaciones como la aproximación de Bakhshali , pero converge más lentamente. Por lo tanto, no es una forma de cálculo particularmente eficiente. Para maximizar la tasa de convergencia, elija N de modo que|d|norte2{\displaystyle {\frac {|d|}{N^{2}}}\,}es lo más pequeño posible.

Expansión continua de fracciones

La representación en fracción continua de un número real puede utilizarse en lugar de su expansión decimal o binaria, y esta representación tiene la propiedad de que la raíz cuadrada de cualquier número racional (que no sea ya un cuadrado perfecto) tiene una expansión periódica y repetitiva, similar a como los números racionales tienen expansiones repetitivas en el sistema de notación decimal.

Irracionales cuadráticos (números de la formaa+bdo{\displaystyle {\frac {a+{\sqrt {b}}}{c}}}(donde a , b y c son enteros), y en particular, las raíces cuadradas de enteros, tienen fracciones continuas periódicas . A veces, lo que se desea es encontrar no el valor numérico de una raíz cuadrada, sino su expansión en fracción continua y, por lo tanto, su aproximación racional. Sea S el número positivo para el cual se nos pide encontrar la raíz cuadrada. Entonces, suponiendo que a sea un número que sirve como estimación inicial y r el término restante, podemos escribir:S=a2+r.{\displaystyle S=a^{2}+r.}Dado que tenemosSa2=(S+a)(Sa)=r{\displaystyle S-a^{2}=({\sqrt {S}}+a)({\sqrt {S}}-a)=r}, podemos expresar la raíz cuadrada de S como S=a+ra+S.{\displaystyle {\sqrt {S}}=a+{\frac {r}{a+{\sqrt {S}}}}.}

Al aplicar esta expresión paraS{\displaystyle {\sqrt {S}}}Al término denominador de la fracción, tenemos: S=a+ra+(a+ra+S)=a+r2a+ra+S.{\displaystyle {\sqrt {S}}=a+{\frac {r}{a+(a+{\frac {r}{a+{\sqrt {S}}}})}}=a+{\frac {r}{2a+{\frac {r}{a+{\sqrt {S}}}}}}.}

Notación compacta : la expansión numerador/denominador para fracciones continuas (arriba) es engorrosa de escribir e integrar en sistemas de formato de texto. Por ello, los matemáticos han ideado varias notaciones alternativas, tales como: S=a+r2a+r2a+r2a+{\displaystyle {\sqrt {S}}=a+{\frac {r}{2a+}}\,{\frac {r}{2a+}}\,{\frac {r}{2a+}}\cdots }

Cuandor=1{\displaystyle r=1}En todo el texto, una notación aún más compacta es: [ Nota 7 ][a;2a,2a,2a,]{\displaystyle [a;2a,2a,2a,\cdots ]} Para fracciones continuas periódicas (que son todas las raíces cuadradas de cuadrados no perfectos), el período se representa solo una vez, con una línea superior para indicar una repetición no terminante de la parte subrayada: [ Nota 8 ][a;2a¯]{\displaystyle [a;{\overline {2a}}]}

Para √2 , el valor de a = 1 , por lo que su representación es: [1;2¯]{\displaystyle [1;{\overline {2}}]}

Procediendo de esta manera, obtenemos una fracción continua generalizada para la raíz cuadrada como S=a+r2a+r2a+r2a+{\displaystyle {\sqrt {S}}=a+{\cfrac {r}{2a+{\cfrac {r}{2a+{\cfrac {r}{2a+\ddots }}}}}}}

El primer paso para evaluar dicha fracción [ 19 ] y obtener una raíz es realizar sustituciones numéricas para la raíz del número deseado y el número de denominadores seleccionados. Por ejemplo, en forma canónica, r = 1 y para √2 , a = 1 , por lo que la fracción continua numérica para 3 denominadores es: 21+12+12+12{\displaystyle {\sqrt {2}}\approx 1+{\cfrac {1}{2+{\cfrac {1}{2+{\cfrac {1}{2}}}}}}}

El paso 2 consiste en simplificar la fracción continua de abajo hacia arriba, un denominador a la vez, para obtener una fracción racional cuyo numerador y denominador sean números enteros. La simplificación se realiza de la siguiente manera (tomando los tres primeros denominadores): 1+12+12+12=1+12+152=1+12+25=1+1125=1+512=1712{\displaystyle {\begin{aligned}1+{\cfrac {1}{2+{\cfrac {1}{2+{\cfrac {1}{2}}}}}}&=1+{\cfrac {1}{2+{\cfrac {1}{\frac {5}{2}}}}}\\&=1+{\cfrac {1}{2+{\cfrac {2}{5}}}}=1+{\cfrac {1}{\frac {12}{5}}}\\&=1+{\cfrac {5}{12}}={\frac {17}{12}}\end{aligned}}}

Finalmente (paso 3), divide el numerador por el denominador de la fracción racional para obtener el valor aproximado de la raíz: 17÷12=1.42{\displaystyle 17\div 12=1.42}redondeado a tres cifras decimales de precisión.

El valor real de √2 es 1,41 con tres cifras significativas. El error relativo es del 0,17%, por lo que la fracción racional es buena con casi tres cifras de precisión. Tomar más denominadores da aproximaciones sucesivamente mejores: cuatro denominadores dan como resultado la fracción4129=1.4137{\displaystyle {\frac {41}{29}}=1.4137}, con una precisión de casi 4 dígitos, etc.

Los siguientes son ejemplos de raíces cuadradas , sus fracciones continuas simples y sus primeros términos —llamados convergentes— hasta el denominador 99 inclusive :

En general, cuanto mayor sea el denominador de una fracción racional, mejor será la aproximación. También se puede demostrar que truncar una fracción continua produce una fracción racional que es la mejor aproximación a la raíz cuadrada de cualquier fracción con denominador menor o igual al denominador de dicha fracción ; por ejemplo, ninguna fracción con denominador menor o igual a 70 es una aproximación tan buena a √2 como 99/70 .

Aproximaciones que dependen de la representación de punto flotante

Un número se representa en formato de punto flotante comometro×bpag{\displaystyle m\times b^{p}}que también se llama notación científica . Su raíz cuadrada esmetro×bpag/2{\displaystyle {\sqrt {m}}\times b^{p/2}}y fórmulas similares se aplicarían a raíces cúbicas y logaritmos. A primera vista, esto no supone una mejora en la simplicidad, pero supongamos que solo se requiere una aproximación: entonces simplementebpag/2{\displaystyle b^{p/2}}es bueno hasta un orden de magnitud. A continuación, reconozca que algunas potencias, p , serán impares, por lo tanto, para 3141.59 = 3.14159 × 103 En lugar de trabajar con potencias fraccionarias de la base, multiplique la mantisa por la base y reste uno a la potencia para que sea par. La representación ajustada será el equivalente a 31.4159 × 102 de modo que la raíz cuadrada será31.4159 × 101 .

Si se toma la parte entera de la mantisa ajustada, solo pueden existir los valores del 1 al 99, que podrían usarse como índice en una tabla de 99 raíces cuadradas precalculadas para completar la estimación. Una computadora que use base dieciséis requeriría una tabla más grande, pero una que use base dos requeriría solo tres entradas: los bits posibles de la parte entera de la mantisa ajustada son 01 (la potencia es par, por lo que no hubo desplazamiento, recordando que un número de punto flotante normalizado siempre tiene un dígito de orden superior distinto de cero) o, si la potencia es impar, 10 u 11, que son los dos primeros bits de la mantisa original. Así, 6.25 = 110.01 en binario, normalizado a 1.1001 × 2 2 una potencia par por lo que los bits emparejados de la mantisa son 01, mientras que .625 = 0.101 en binario se normaliza a 1.01 × 2 −1 una potencia impar por lo que el ajuste es a 10.1 × 2 −2 y los bits emparejados son 10. Nótese que el bit de orden inferior de la potencia se refleja en el bit de orden superior de la mantisa de pares. Una potencia par tiene su bit de orden inferior cero y la mantisa ajustada comenzará con 0, mientras que para una potencia impar ese bit es uno y la mantisa ajustada comenzará con 1. Por lo tanto, cuando la potencia se divide a la mitad, es como si su bit de orden inferior se desplazara para convertirse en el primer bit de la mantisa de pares.

Una tabla con solo tres entradas podría ampliarse incorporando bits adicionales de la mantisa. Sin embargo, con las computadoras, en lugar de calcular una interpolación en una tabla, suele ser mejor encontrar un cálculo más simple que dé resultados equivalentes. Todo depende ahora de los detalles exactos del formato de la representación, además de las operaciones disponibles para acceder y manipular las partes del número. Por ejemplo, Fortran ofrece una EXPONENT(x)función para obtener la potencia. El esfuerzo invertido en idear una buena aproximación inicial se recupera evitando así las iteraciones adicionales del proceso de refinamiento que habrían sido necesarias para una aproximación deficiente. Dado que estas son pocas (una iteración requiere una división, una suma y una división por la mitad), la restricción es severa.

Muchas computadoras siguen la representación IEEE (o suficientemente similar), y se puede obtener una aproximación muy rápida a la raíz cuadrada para comenzar el método de Newton. La técnica que sigue se basa en el hecho de que el formato de punto flotante (en base dos) aproxima el logaritmo en base 2. Es decir,registro2(metro×2pag)=pag+registro2(metro){\displaystyle \log _{2}(m\times 2^{p})=p+\log _{2}(m)}

Así, para un número de punto flotante de precisión simple de 32 bits en formato IEEE (donde cabe destacar que la potencia tiene un sesgo de 127 añadido para la forma representada), se puede obtener el logaritmo aproximado interpretando su representación binaria como un entero de 32 bits y escalándolo por223{\displaystyle 2^{-23}}y eliminando un sesgo de 127, es decir incógnitaentero223127registro2(incógnita).{\displaystyle x_{\text{int}}\cdot 2^{-23}-127\approx \log _{2}(x).}

Por ejemplo, 1.0 está representado por un número hexadecimal 0x3F800000 , que representaría1065353216=127×223{\displaystyle 1065353216=127\times 2^{23}}si se toma como un número entero. Usando la fórmula anterior se obtiene1065353216×223127=0{\displaystyle 1065353216\times 2^{-23}-127=0}, como se esperaba deregistro2(1.0){\displaystyle \log _{2}(1.0)}. De manera similar se obtiene 0,5 a partir de 1,5 ( 0x3FC00000 ).

Para obtener la raíz cuadrada, divida el logaritmo por 2 y convierta el valor de nuevo. El siguiente programa demuestra la idea. El bit menos significativo del exponente se permite intencionalmente que se propague a la mantisa. Una forma de justificar los pasos de este programa es asumir que b es el sesgo del exponente y n es el número de bits almacenados explícitamente en la mantisa y luego demostrar que ((12(incógnitaentero/2norteb))+b)2norte=12(incógnitaentero2norte)+(12(b+1))2norte.{\displaystyle \left(\left({\tfrac {1}{2}}\left(x_{\text{int}}/2^{n}-b\right)\right)+b\right)\cdot 2^{n}={\tfrac {1}{2}}\left(x_{\text{int}}-2^{n}\right)+\left({\tfrac {1}{2}}\left(b+1\right)\right)\cdot 2^{n}.}

/* Se asume que float está en el formato de punto flotante de precisión simple IEEE 754 */ #include <stdint.h>unión FloatUInt { float f ; uint32_t i ; }float sqrtApprox ( float z ) { union FloatUInt val = { z }; // Convertir tipo, conservando el patrón de bits /*  * Para justificar el siguiente código, demuestre que  *  * ((((val.i / 2^m) - b) / 2) + b) * 2^m = ((val.i - 2^m) / 2) + ((b + 1) / 2) * 2^m)  *  * donde  *  * b = sesgo del exponente  * m = número de bits de la mantisa  */ val . i -= 1 << 23 ; // Restar 2^m. val . i >>= 1 ; // Dividir por 2. val . i += 1 << 29 ; // Sumar ((b + 1) / 2) * 2^m.// Interpretar de nuevo como float return val . f ; }

Las tres operaciones matemáticas que forman el núcleo de la función anterior se pueden expresar en una sola línea. Se puede agregar un ajuste adicional para reducir el error relativo máximo. Por lo tanto, las tres operaciones, sin incluir la conversión de tipo, se pueden reescribir como

val.i = ( 1 << 29 ) + ( val.i >> 1 ) - ( 1 << 22 ) + a ;

donde a es un sesgo para ajustar los errores de aproximación. Por ejemplo, con a = 0 los resultados son precisos para potencias pares de 2 (por ejemplo, 1.0), pero para otros números los resultados serán ligeramente demasiado grandes (por ejemplo, 1.5 para 2.0 en lugar de 1.414... con un error del 6%). Con a = 0x4B0D2 , el error relativo máximo se minimiza a ±3.5%. Para a = 0, la aproximación es mayor o igual que la raíz cuadrada de val para cualquier valor de val .

Si la aproximación se va a utilizar como una estimación inicial para el método de Newton a la ecuación(1/incógnita2)S=0{\displaystyle (1/x^{2})-S=0}, entonces se prefiere la forma recíproca que se muestra en la siguiente sección.

La aproximación se puede optimizar aún más combinando la suma de a , 1 << 29 y -1 << 22 en una sola operación. Esto da como resultado (val.i >> 1) + 0x1FBB4F2Epara a = 0x4B0D2 que minimiza la corrección de error y para a = 0.(val.i >> 1) + 1FC00000

Recíproco de la raíz cuadrada

A continuación se incluye una variante de la rutina anterior, que se puede utilizar para calcular el recíproco de la raíz cuadrada, es decir,incógnita1/2{\displaystyle x^{-1/2}}En cambio, fue escrito por Greg Walsh. La aproximación de desplazamiento entero produjo un error relativo de menos del 4%, y el error se redujo aún más a 0,15% con una iteración del método de Newton en la siguiente línea. [ 20 ] En gráficos por computadora es una forma muy eficiente de normalizar un vector.

unión FloatInt { float x ; int i ; };float inv_sqrt ( float x ) { float xhalf = 0.5f * x ; union FloatInt u ; u.x = x ; u.i = 0x5f375a86 - ( u.i >> 1 ) ; // La siguiente línea se puede repetir cualquier número de veces para aumentar la precisión u.x = u.x * ( 1.5f - xhalf * u.x * u.x ) ; return u.x ; }

Algunos componentes de hardware VLSI implementan la raíz cuadrada inversa utilizando una estimación polinómica de segundo grado seguida de una iteración de Goldschmidt . [ 21 ]

Cuadrado negativo o complejo

Si S  <  0, entonces su raíz cuadrada principal es S=|S|i.{\displaystyle {\sqrt {S}}={\sqrt {\vert S\vert }}\,\,i\,.}

Si S  = a + bi donde a y b son reales y b ≠ 0, entonces su raíz cuadrada principal es    S=|S|+a2+sgn(b)|S|a2i.{\displaystyle {\sqrt {S}}={\sqrt {\frac {\vert S\vert +a}{2}}}\,+\,\operatorname {sgn}(b){\sqrt {\frac {\vert S\vert -a}{2}}}\,\,i\,.}

Esto se puede verificar elevando al cuadrado la raíz. [ 22 ] [ 23 ] Aquí |S|=a2+b2{\displaystyle \vert S\vert ={\sqrt {a^{2}+b^{2}}}}

es el módulo de S. La raíz cuadrada principal de un número complejo se define como la raíz con la parte real no negativa.

Véase también

Notas

  1. Los factores dos y seis se utilizan porque se aproximan a las medias geométricas de los valores mínimo y máximo posibles con el número de dígitos dado:110=1041,78{\displaystyle {\sqrt {{\sqrt {1}}\cdot {\sqrt {10}}}}={\sqrt[{4}]{10}}\approx 1.78\,}y10100=100045.62{\displaystyle {\sqrt {{\sqrt {10}}\cdot {\sqrt {100}}}}={\sqrt[{4}]{1000}}\approx 5.62\,}.
  2. La estimación sin redondear tiene un error absoluto máximo de 2,65 en 100 y un error relativo máximo de 26,5% en y=1, 10 y 100.
  3. Si el número está exactamente a la mitad entre dos cuadrados, como 30.5, adivina el número mayor, que en este caso es 6.
  4. Esta es, incidentalmente, la ecuación de la recta tangente a y = x 2 en y = 1.
  5. La precisión predeterminada puede variar de 1 adecimal.MAX_PREC
  6. Ver arriba .
  7. ver: Fracción continua#Notaciones
  8. ver: Fracción continua periódica

Referencias

  1. Jackson 2011 .
  2. Fowler y Robson 1998 .
  3. 1 2 Heath 1921 .
  4. Baumann 2024 .
  5. Johnson 2015 .
  6. Nemiroff y Bonnell 1994 .
  7. Nemiroff y Bonnell 1994a .
  8. Bailey y Borwein 2012 .
  9. Simplemente Curioso 2018 .
  10. Herrero Piñeyro, PJ; Linero Bas, A.; Massa Esteve, MR; Mellado Romero, A. (2023). "Un problema sobre la aproximación de n-raíces basado en el trabajo de Viète" (PDF) . MATerials MATemàtics . 5 : 1–27 . Recuperado el 15 de noviembre de 2025 .Véase la página 8.
  11. Guy 1985 .
  12. ^ Steinarson, Corbit y Hendry 2003 .
  13. Wilkes, Wheeler y Gill 1951 , págs. 146 . 
  14. Campbell-Kelly 2009 .
  15. Gower 1958 .
  16. Goldschmidt, Robert E. (1964). Aplicaciones de la división por convergencia (PDF) (Tesis). Disertación de maestría. MIT OCLC 34136725. Archivado (PDF) del original el 10 de diciembre de 2015. Recuperado el 15 de septiembre de 2015 . 
  17. "Autores" . IBM Journal of Research and Development . 11 : 125–127 . 1967. doi : 10.1147/rd.111.0125 . Archivado del original el 18 de julio de 2018.
  18. Markstein 2004 .
  19. Sardina 2007 , pág. 10, 2.3j.
  20. Lomont 2003 .
  21. Piñeiro & Díaz Bruguera 2002 .
  22. Abramowitz y Stegun 1970 , pág. 17, Sección 3.7.26.
  23. Cooke 2008 , pág. 59.

Bibliografía

  • Abramowitz, Milton ; Stegun, Irene A. , eds. (1970) [1964]. Manual de funciones matemáticas con fórmulas, gráficas y tablas matemáticas . Serie de Matemáticas Aplicadas 55 (Novena edición  ). Washington D.C.: Dover Publications . ISBN 978-0-486-61272-0. LCCN 64-60036 . MR 0167642 .  
  • Bailey, David; Borwein, Jonathan (2012). " Raíces cuadradas de la antigua India: un ejercicio de paleomatemáticas forenses" (PDF) . American Mathematical Monthly . Vol.  119, n.º  8, págs. 646–657 . Consultado el 14 de septiembre de 2017 . 
  • Campbell-Kelly, Martin (septiembre de 2009). "Origen de la computación". Scientific American . 301 (3): 62– 69. Bibcode : 2009SciAm.301c..62C . doi : 10.1038/scientificamerican0909-62 . JSTOR 26001527. PMID 19708529 .  
  • Cooke, Roger (2008). Álgebra clásica: su naturaleza, orígenes y usos . John Wiley and Sons. ISBN 978-0-470-25952-8.
  • Fowler, David; Robson, Eleanor (1998). "Aproximaciones de raíz cuadrada en matemáticas babilónicas antiguas: YBC 7289 en contexto" (PDF) . Historia Mathematica . 25 (4): 376. doi : 10.1006/hmat.1998.2209 .
  • Gower, John C. (1958). "Una nota sobre un método iterativo para la extracción de raíces" . The Computer Journal . 1 (3): 142– 143. doi : 10.1093/comjnl/1.3.142 .
  • Guy, Martin (1985). "Raíz cuadrada de enteros rápida mediante el algoritmo del ábaco del Sr. Woo" . Universidad de Kent en Canterbury (UKC). Archivado del original el 6 de marzo de 2012. Recuperado el 5 de octubre de 2025 .
  • Heath, Thomas (1921). Historia de las matemáticas griegas, vol. 2. Oxford: Clarendon Press. págs. 323-324 . 
  • Jackson, Terence (1 de julio de 2011). "95.42 Raíces cuadradas irracionales de números naturales: un enfoque geométrico" . The Mathematical Gazette . 95 (533): 327–330 . doi : 10.1017/S0025557200003193 . ISSN 0025-5572 . S2CID 123995083 .  
  • Johnson, SG (4 de febrero de 2015). "Raíces cuadradas mediante el método de Newton" (PDF) . Curso 18.335 del MIT: Introducción a los métodos numéricos . Recuperado el 12 de octubre de 2025 .
  • Lomont, Chris (2003). "Raíz cuadrada inversa rápida" (PDF) .
  • Markstein, Peter (noviembre de 2004). División de software y raíz cuadrada utilizando algoritmos de Goldschmidt (PDF) . 6.ª Conferencia sobre números reales y computadoras . Dagstuhl , Alemania. CiteSeerX 10.1.1.85.9648 . 
  • Nemiroff, Robert; Bonnell, Jerry (1994). "La raíz cuadrada de 2 a 1 millón de dígitos" . Imagen astronómica del día de la NASA . Recuperado el 21 de octubre de 2025 .
  • Nemiroff, Robert; Bonnell, Jerry (1994a). "La raíz cuadrada de 2 a 10 millones de dígitos" . Imagen astronómica del día de la NASA . Recuperado el 21 de octubre de 2025 .
  • Piñeiro, José-Alejandro; Díaz Bruguera, Javier (diciembre de 2002). "Cálculo de alta velocidad y doble precisión de recíprocos, divisiones, raíces cuadradas e inversas de raíces cuadradas" . IEEE Transactions on Computers . 51 (12): 1377– 1388. Bibcode : 2002ITCmp..51.1377P . doi : 10.1109/TC.2002.1146704 .
  • Sardina, Manny (2007). "Método general para la extracción de raíces mediante fracciones continuas (plegadas)" . Surrey (Reino Unido).
  • Simply Curious (5 de junio de 2018). "Profundizando en el manuscrito de Bakhshali" . Blog de Simply Curious . Consultado el 21 de diciembre de 2020 .
  • Steinarson, Arne; Corbit, Dann; Hendry, Mateo (2003). "Función de raíz cuadrada entera" .
  • Wilkes, MV ; Wheeler, DJ ; Gill, S. (1951). La preparación de programas para una computadora digital electrónica . Oxford: Addison-Wesley. pág.  262. OCLC 475783493 . 
  • Baumann, Claude (2024). "Jugando con la raíz cuadrada" (PDF) . computarium.lcd . Consultado el 28 de enero de 2026 .
  • Weisstein, Eric W. "Algoritmos de raíz cuadrada" . MathWorld .
  • Jarvis, AF "Raíces cuadradas por sustracción" (PDF) . afjarvis.staff.shef.ac.uk . Archivado del original (PDF) el 10 de abril de 2022.
  • Radovic, Andrija. "Algoritmo de raíz cuadrada entera" . andrijar.com . Consultado el 7 de noviembre de 2025 .
  • Egbert, William E. (mayo de 1977). "Algoritmos de calculadora personal I: raíces cuadradas" (PDF) . Hewlett-Packard Journal : 22–24 . Recuperado el 7 de noviembre de 2025 .
  • "Calculadora para aprender la raíz cuadrada" . calculatorsquareroot.com . Consultado el 7 de noviembre de 2025 .