pérdida inesperada de precisión cuando se divide dobles

Question

tengo un getSlope función que toma como parámetros 4 dobles y devuelve otro doble calculada usando estos parámetros dados de la siguiente manera:

double QSweep::getSlope(double a, double b, double c, double d){ 
double slope; 
slope=(d-b)/(c-a); 
return slope; 
} 

El problema es que cuando se llama a esta función con argumentos por ejemplo:

getSlope(2.71156, -1.64161, 2.70413, -1.72219); 

el resultado devuelto es:

10.8557 

y este no es un buen resultado para mis cálculos. He calculado la pendiente usando Mathematica y el resultado de la pendiente para los mismos parámetros es:

10.8452 

o más dígitos de precisión:

10.845222072678331. 

El resultado devuelto por mi programa no es bueno en mis cálculos adicionales. Además, no entiendo cómo devuelve el programa 10.8557 a partir de 10.845222072678331 (suponiendo que este es el resultado aproximado de la división)? ¿Cómo puedo obtener el buen resultado para mi división?

gracias de antemano, Madalina

---

puedo imprimir el resultado utilizando la línea de comandos:

std::cout<<slope<<endl; 

Puede ser que mis parámetros no son tal vez bueno, ya que los leo de otro programa (que calcula un gráfico; después de leer estos parámetros desde su gráfica, los he visualizado para ver su valor, pero tal vez los vectores mostrados no tienen la misma precisión interna para el valor calculado ... No sé, es realmente extraño . Algunos numéricos aparecen los errores ..)

Cuando se calcula el gráfico del que estoy leyendo mis parámetros, se usan algunas librerías numéricas escritas en C++ (con plantillas). No se usa OpenGL para este cálculo.

gracias, Madalina

Answer 1

Podría ser que utiliza DirectX u OpenGL en su proyecto? Si es así, pueden desactivar la doble precisión y obtendrás resultados extraños.

Usted puede comprobar los ajustes de precisión con

std::sqrt(x) * std::sqrt(x) 

El resultado tiene que ser bastante cerca de x. Conocí este problema hace mucho tiempo y paso un mes revisando todas las fórmulas. Pero luego me he encontrado

D3DCREATE_FPU_PRESERVE 

Answer 2

El problema aquí es que (c-a) es pequeña, por lo que los errores de redondeo inherentes a las operaciones de punto flotante se magnifica en este ejemplo. Una solución general es volver a trabajar su ecuación para que no se divida por un número pequeño, aunque no estoy seguro de cómo lo haría aquí.

EDIT:

Neil tiene razón en su comentario a esta pregunta, calcula la respuesta en VB usando Dobles y tiene la misma respuesta que Mathematica.

Answer 3

el siguiente código:

#include <iostream> 
using namespace std; 

double getSlope(double a, double b, double c, double d){ 
    double slope; 
    slope=(d-b)/(c-a); 
    return slope; 
} 

int main() { 
    double s = getSlope(2.71156, -1.64161, 2.70413, -1.72219); 
    cout << s << endl; 
} 

da un resultado de 10,8452 con g ++. ¿Cómo está imprimiendo el resultado en su código?

Answer 4

Mejor Imprima también los argumentos. Cuando esté, como supongo, transfiriendo parámetros en notación decimal, perderá precisión para todos y cada uno de ellos. El problema es que 1/5 es una serie infinita en binario, por ejemplo, 0.2 se convierte en .001001001 .... Además, los decimales se cortan al convertir un flotador binario en una representación textual en decimal.

Además de eso, a veces el compilador elige la velocidad sobre la precisión. Esto debería ser un interruptor de compilador documentado.

Answer 5

Patrick parece ser correcto sobre (ca) siendo la causa principal:

db = -1,72219 - (-1,64161) = -0,08058

ca = 2, 70413 - 2,71156 = -0,00743

S = (db)/(ca) = -0,08058/-0,00743 = 10,845222

Usted comienza con seis dígitos de precisión a través de la resta obtienes una reducción a 3 y cuatro dígitos. Mi mejor suposición es que pierdes precisión adicional porque el número -0,00743 no se puede representar exactamente en un doble. Trate de usar las variables intermedias con una precisión más grande, como este:

double QSweep::getSlope(double a, double b, double c, double d) 
{ 
    double slope; 
    long double temp1, temp2; 

    temp1 = (d-b); 
    temp2 = (c-a); 
    slope = temp1/temp2; 

    return slope; 
} 

Answer 6

He tratado con el flotador en lugar de doble y llego 10.845110 como resultado. Todavía se ve mejor que el resultado de madalina.

EDIT:

creo que sé por qué se obtiene este resultado. Si obtiene parámetros a, b, c y d de otro lugar y lo imprime, le da valores redondeados. Entonces, si lo coloca en Mathemtacia (o calc;)) le dará un resultado diferente.

He intentado cambiar un poco uno de tus parámetros. Cuando lo hice:

double c = 2.7041304; 

obtengo 10.845806. Solo agrego 0.0000004 a c! Así que creo que sus "errores" no son errores. Imprima a, b, c y d con mayor precisión y luego colóquelos en Mathematica.

Answer 7

Los resultados que obtiene son consistentes con la aritmética de 32 bits. Sin saber más acerca de su entorno, no es posible aconsejar qué hacer.

Suponiendo que el código que se muestra es lo que se está ejecutando, es decir, no está convirtiendo nada a cadenas o flotantes, entonces no hay una solución dentro de C++. Está fuera del código que ha mostrado y depende del entorno.

Como Patrick McDonald y Treb trajeron tanto la precisión de sus entradas como el error en a-c, pensé en echarle un vistazo. Una técnica para observar los errores de redondeo es la aritmética de intervalos, que hace que los límites superior e inferior que el valor representa sean explícitos (están implícitos en los números de coma flotante y se fijan a la precisión de la representación). Tratando cada valor como un límite superior e inferior, y extendiendo los límites por el error en la representación (aproximadamente x * 2^-53 para un valor doble x), obtienes un resultado que da los límites inferior y superior en el precisión de un valor, teniendo en cuenta los errores de precisión en el peor de los casos.

Por ejemplo, si tiene un valor en el rango [1.0, 2.0] y resta un valor en el rango [0.0, 1.0], entonces el resultado debe estar en el rango [debajo (0.0), arriba (2.0)] ya que el resultado mínimo es 1.0-1.0 y el máximo es 2.0-0.0. below y above son equivalentes a piso y techo, pero para el siguiente valor representable en lugar de enteros.

El uso de intervalos que representan el peor caso de doble redondeo:

getSlope(
a = [2.7115599999999995262:2.7115600000000004144], 
b = [-1.6416099999999997916:-1.6416100000000002357], 
c = [2.7041299999999997006:2.7041300000000005888], 
d = [-1.7221899999999998876:-1.7221900000000003317]) 
(d-b) = [-0.080580000000000526206:-0.080579999999999665783] 
(c-a) = [-0.0074300000000007129439:-0.0074299999999989383218] 

to double precision [10.845222072677243474:10.845222072679954195] 

Así que, aunque c-a es pequeño comparado con c o a, todavía es grande en comparación con doble redondeo, por lo que si se utilizara el doble peor imaginable redondeo de precisión, entonces usted puede confiar en que el valor sea preciso a 12 cifras - 10.8452220727. Ha perdido algunas cifras con doble precisión, pero aún está trabajando para algo más que la importancia de su entrada.

Pero si las entradas eran sólo exacta a los números de cifras significativas, a continuación, en lugar de ser el doble valor de 2.71156 +/- EPS, entonces el rango de entrada sería [2.711555,2.711565], para que pueda obtener el resultado:

getSlope(
a = [2.711555:2.711565], 
b = [-1.641615:-1.641605], 
c = [2.704125:2.704135], 
d = [-1.722195:-1.722185]) 
(d-b) = [-0.08059:-0.08057] 
(c-a) = [-0.00744:-0.00742] 

to specified accuracy [10.82930108:10.86118598] 

que es un rango mucho más amplio.

Pero tendrías que salir de tu camino para rastrear la precisión en los cálculos, y los errores de redondeo inherentes al punto flotante no son significativos en este ejemplo - es preciso para 12 figuras con el peor caso de doble redondeo de precisión.

Por otro lado, si sus entradas solo son conocidas por 6 cifras, en realidad no importa si obtiene 10.8557 o 10.8452. Ambos están dentro de [10.82930108: 10.86118598].

Answer 8

Si bien la discusión académica es excelente para aprender sobre las limitaciones de los lenguajes de programación, puede encontrar la solución más simple para el problema es una estructura de datos para arbitrary precision arithmetic.

Esto tendrá algunos gastos generales, pero debería poder encontrar algo con una precisión bastante garantizada.