¿Cómo puedo encontrar la raíz cúbica de un número de manera eficiente? Creo que se puede usar el método de Newton-Raphson, pero no sé cómo adivinar la solución inicial mediante programación para minimizar el número de iteraciones.Sembrando la iteración de Newton para raíz cúbica eficientemente
Esta es una pregunta engañosamente compleja. Here es una buena encuesta de algunos posibles enfoques.
Podría sonar grosero, pero creo que en caso de que el enlace se acabe, debe agregar * algún * contenido a su respuesta. (Copia y pega si quieres) – Astrobleme
En vista de la "rotura del enlace" que superó a la respuesta aceptada, daré una respuesta más independiente centrándose en el tema de obtener rápidamente una conjetura inicial adecuada para la iteración superlineal.
La "encuesta" por metamerist (Wayback link) proporcionó algunas comparaciones de tiempo para varias combinaciones de valor inicial/iteración (se incluyen los métodos Newton y Halley). Sus referencias son a trabajos por W. Kahan, "Computing a Real Cube Root", y por K. Turkowski, "Computing the Cube Root".
metamarist actualiza la técnica de bits de tocar el violín era DEC-VAX de W. Kahan con este fragmento, que "asume enteros de 32 bits" y se basa en formato IEEE 754 en dobles "para generar las estimaciones iniciales con 5 bits de precisión" :
inline double cbrt_5d(double d)
{
const unsigned int B1 = 715094163;
double t = 0.0;
unsigned int* pt = (unsigned int*) &t;
unsigned int* px = (unsigned int*) &d;
pt[1]=px[1]/3+B1;
return t;
}
El código por K. Turkowski proporciona un poco más de precisión ("aproximadamente 6 bits") por una escala convencional poderes-de-dos en float fr, seguido de una aproximación cuadrática a su raíz cubo sobre intervalo [0,125 , 1.0):
/* Compute seed with a quadratic qpproximation */
fr = (-0.46946116F * fr + 1.072302F) * fr + 0.3812513F;/* 0.5<=fr<1 */
y una posterior restauración del exponente de dos (ajustado a un tercio). La extracción y restauración del exponente/mantisa utilizan math library calls a frexp y ldexp.
Comparación con otras aproximaciones "semilla" raíz cúbica
Para apreciar esas aproximaciones raíz cúbica que necesitamos para compararlas con otras formas posibles. Primero, los criterios para juzgar: consideramos la aproximación en el intervalo [1/8,1], y usamos el mejor error relativo (minimizando el máximo).
Es decir, si f(x) es una aproximación propuesta para x^{1/3}, que encuentran su error relativo:
error_rel = max | f(x)/x^(1/3) - 1 | on [1/8,1]
La aproximación más simple sería, por supuesto utilizar un único constante en el intervalo, y el mejor error relativo en ese caso se logra seleccionando f_0(x) = sqrt(2)/2, la media geométrica de los valores en los puntos finales. Esto le da 1.27 bits de precisión relativa, un punto de inicio rápido pero sucio para una iteración de Newton.
una mejor aproximación sería la mejor polinomio de primer grado:
f_1(x) = 0.6042181313*x + 0.4531635984
Esto da 4,12 bits de precisión relativa, una gran mejora pero por debajo de los 5-6 bits de precisión relativa prometidos por los respectivos métodos de Kahan y Turkowski. Pero está en el estadio y usa solo una multiplicación (y una adición).
Por último, ¿y si nos dejamos una división en lugar de una multiplicación?Resulta que con una división y dos "adiciones" podemos tener la mejor función lineal-fraccionada:
f_M(x) = 1.4774329094 - 0.8414323527/(x+0.7387320679)
que da 7.265 bits de precisión relativa.
De un vistazo esto parece ser un enfoque atractivo, pero una antigua regla general era tratar el costo de una división FP como tres multiplicaciones FP (y para ignorar las adiciones y sustracciones). Sin embargo, con los diseños actuales de FPU esto no es realista. Si bien el costo relativo de las multiplicaciones para agregar/restar ha bajado, en la mayoría de los casos a un factor de dos o incluso a la igualdad, el costo de la división no ha disminuido, pero a menudo ha aumentado hasta 7-10 veces el costo de la multiplicación. Por lo tanto, debemos ser mezquinos con nuestras operaciones de división.
static double cubeRoot(double num) {
double x = num;
if(num >= 0) {
for(int i = 0; i < 10 ; i++) {
x = ((2 * x * x * x) + num)/(3 * x * x);
}
}
return x;
}
Parece que ya se ha abordado la cuestión de optimización, pero me gustaría añadir una mejora de la función cubeRoot() publicado aquí, para otras personas de tropiezo en esta página en busca de un algoritmo de raíz cúbica rápida .
El algoritmo existente funciona bien, pero fuera del rango de 0-100 da resultados incorrectos.
Aquí hay una versión revisada que funciona con números entre -/+ 1 cuatrillón (1E15). Si necesita trabajar con números más grandes, simplemente use más iteraciones.
static double cubeRoot(double num){
boolean neg = (num < 0);
double x = Math.abs(num);
for(int i = 0, iterations = 60; i < iterations; i++){
x = ((2 * x * x * x) + num)/(3 * x * x);
}
if(neg){ return 0 - x; }
return x;
}
En cuanto a optimización, supongo que el cartel original estaba preguntando cómo predecir el número mínimo de iteraciones para un resultado exacto, dado un tamaño de entrada arbitraria. Pero parece que, en la mayoría de los casos generales, la ganancia de la optimización no vale la pena la complejidad añadida. Incluso con la función anterior, 100 iteraciones toman menos de 0.2 ms en el hardware de consumo promedio. Si la velocidad era de suma importancia, consideraría usar tablas de búsqueda pre calculadas. Pero esto proviene de un desarrollador de escritorio, no de un ingeniero de sistemas integrado.
¿Qué quiere decir con "número"? – starblue