Actualmente estoy escribiendo un código donde tengo algo como:

double a = SomeCalculation1();
double b = SomeCalculation2();

if (a < b)
    DoSomething2();
else if (a > b)
    DoSomething3();

Y luego, en otros lugares, puede que necesite hacer igualdad:

double a = SomeCalculation3();
double b = SomeCalculation4();

if (a == 0.0)
   DoSomethingUseful(1 / a);
if (b == 0.0)
   return 0; // or something else here

En resumen, tengo muchas matemáticas de punto flotante y necesito hacer varias comparaciones para las condiciones. No puedo convertirlo en matemáticas enteras porque tal cosa no tiene sentido en este contexto.

He leído antes que las comparaciones de punto flotante pueden ser poco confiables, ya que puede tener cosas como esta sucediendo:

double a = 1.0 / 3.0;
double b = a + a + a;
if ((3 * a) != b)
    Console.WriteLine("Oh no!");

En resumen, me gustaría saber: ¿Cómo puedo comparar de manera confiable números de punto flotante (menor que, mayor que, igualdad)?

El rango de números que estoy usando es aproximadamente de 10E-14 a 10E6, por lo que necesito trabajar tanto con números pequeños como grandes.

Lo he etiquetado como independiente del idioma porque estoy interesado en cómo puedo lograrlo sin importar el idioma que esté usando.

La comparación de mayor / menor no es realmente un problema a menos que esté trabajando justo al borde del límite de flotación / doble precisión.

Para una comparación "difusa es igual", esto (código Java, debería ser fácil de adaptar) es lo que se me ocurrió para The Floating-Point Guide después de mucho trabajo y teniendo en cuenta muchas críticas:

public static boolean nearlyEqual(float a, float b, float epsilon) {
    final float absA = Math.abs(a);
    final float absB = Math.abs(b);
    final float diff = Math.abs(a - b);

    if (a == b) { // shortcut, handles infinities
        return true;
    } else if (a == 0 || b == 0 || diff < Float.MIN_NORMAL) {
        // a or b is zero or both are extremely close to it
        // relative error is less meaningful here
        return diff < (epsilon * Float.MIN_NORMAL);
    } else { // use relative error
        return diff / (absA + absB) < epsilon;
    }
}

Viene con un conjunto de pruebas. Debe descartar inmediatamente cualquier solución que no lo haga, porque está virtualmente garantizado que fallará en algunos casos extremos, como tener un valor 0, dos valores muy pequeños opuestos a cero o infinitos.

Una alternativa (consulte el enlace anterior para obtener más detalles) es convertir los patrones de bits de los flotantes a números enteros y aceptar todo dentro de una distancia entera fija.

En cualquier caso, probablemente no exista una solución perfecta para todas las aplicaciones. Idealmente, desarrollaría / adaptaría el suyo propio con un conjunto de pruebas que cubra sus casos de uso reales.

TL; DR

• Utilice la siguiente función en lugar de la solución actualmente aceptada para evitar algunos resultados no deseados en ciertos casos límite, al tiempo que es potencialmente más eficiente.
• Conozca la imprecisión esperada que tiene en sus números y aliméntelos en consecuencia en la función de comparación.

bool nearly_equal(
  float a, float b,
  float epsilon = 128 * FLT_EPSILON, float relth = FLT_MIN)
  // those defaults are arbitrary and could be removed
{
  assert(std::numeric_limits<float>::epsilon() <= epsilon);
  assert(epsilon < 1.f);

  if (a == b) return true;

  auto diff = std::abs(a-b);
  auto norm = std::min((std::abs(a) + std::abs(b)), std::numeric_limits<float>::max());
  // or even faster: std::min(std::abs(a + b), std::numeric_limits<float>::max());
  // keeping this commented out until I update figures below
  return diff < std::max(relth, epsilon * norm);
}

Gráficos, ¿por favor?

Al comparar números de coma flotante, hay dos "modos".

El primero es el modo relativo , donde la diferencia entre xy yse considera relativa a su amplitud |x| + |y|. Cuando se traza en 2D, da el siguiente perfil, donde verde significa igualdad de xy y. (Tomé una epsilonde 0.5 con fines ilustrativos).

El modo relativo es lo que se utiliza para valores de puntos flotantes "normales" o "suficientemente grandes". (Más sobre eso más adelante).

El segundo es un modo absoluto , cuando simplemente comparamos su diferencia con un número fijo. Da el siguiente perfil (nuevamente con un epsilon0,5 y un relth1 para la ilustración).

Este modo absoluto de comparación es el que se utiliza para valores de coma flotante "diminutos".

Ahora la pregunta es, ¿cómo unimos esos dos patrones de respuesta?

En la respuesta de Michael Borgwardt, el cambio se basa en el valor de diff, que debería estar por debajo relth( Float.MIN_NORMALen su respuesta). Esta zona de cambio se muestra sombreada en el siguiente gráfico.

Debido a que relth * epsilones más pequeño que relth, los parches verdes no se pegan, lo que a su vez le da a la solución una mala propiedad: podemos encontrar tripletes de números como ese x < y_1 < y_2y aún x == y2pero x != y1.

Tome este ejemplo sorprendente:

x  = 4.9303807e-32
y1 = 4.930381e-32
y2 = 4.9309825e-32

Tenemos x < y1 < y2, y de hecho y2 - xes más de 2000 veces más grande que y1 - x. Y sin embargo, con la solución actual,

nearlyEqual(x, y1, 1e-4) == False
nearlyEqual(x, y2, 1e-4) == True

Por el contrario, en la solución propuesta anteriormente, la zona de cambio se basa en el valor de |x| + |y|, que está representado por el cuadro sombreado a continuación. Asegura que ambas zonas se conecten con elegancia.

Además, el código anterior no tiene ramificaciones, lo que podría ser más eficiente. Tenga en cuenta que operaciones como maxy abs, que a priori necesitan ramificación, suelen tener instrucciones de montaje dedicadas. Por esta razón, creo que este enfoque es superior a otra solución que sería arreglar el de Michael nearlyEqualcambiando el interruptor de diff < reltha diff < eps * relth, lo que produciría esencialmente el mismo patrón de respuesta.

¿Dónde cambiar entre comparación relativa y absoluta?

El cambio entre esos modos se realiza alrededor relth, lo que se toma como FLT_MINen la respuesta aceptada. Esta elección significa que la representación de float32es lo que limita la precisión de nuestros números de coma flotante.

Esto no siempre tiene sentido. Por ejemplo, si los números que compara son el resultado de una resta, quizás algo en el rango de FLT_EPSILONtenga más sentido. Si son raíces cuadradas de números restados, la imprecisión numérica podría ser aún mayor.

Es bastante obvio cuando consideras comparar un punto flotante con 0. Aquí, cualquier comparación relativa fallará, porque |x - 0| / (|x| + 0) = 1. Por lo tanto, la comparación debe cambiar al modo absoluto cuando xestá en el orden de la imprecisión de su cálculo, y rara vez es tan bajo como FLT_MIN.

Este es el motivo de la introducción del relthparámetro anterior.

Además, al no multiplicar relthcon epsilon, la interpretación de este parámetro es simple y corresponde al nivel de precisión numérica que esperamos de esos números.

Ruido matemático

(guardado aquí principalmente para mi propio placer)

De manera más general, supongo que un operador de comparación de punto flotante con buen comportamiento =~debería tener algunas propiedades básicas.

Los siguientes son bastante obvios:

• auto-igualdad: a =~ a
• simetría: a =~ bimplicab =~ a
• invariancia por oposición: a =~ bimplica-a =~ -b

(No tenemos a =~ be b =~ cimplica a =~ c, =~no es una relación de equivalencia).

Agregaría las siguientes propiedades que son más específicas para las comparaciones de punto flotante

• si a < b < c, entonces a =~ cimplica a =~ b(los valores más cercanos también deben ser iguales)
• si a, b, m >= 0entonces a =~ bimplica a + m =~ b + m(valores más grandes con la misma diferencia también deberían ser iguales)
• si 0 <= λ < 1entonces a =~ bimplica λa =~ λb(quizás menos obvio para el argumento).

Esas propiedades ya dan fuertes restricciones sobre posibles funciones de casi igualdad. La función propuesta anteriormente los verifica. Quizás falten una o varias propiedades obvias.

Cuando se piensa =~en una relación de familia de igualdad =~[Ɛ,t]parametrizada por Ɛy relth, también se podría agregar

• si Ɛ1 < Ɛ2entonces a =~[Ɛ1,t] bimplica a =~[Ɛ2,t] b(igualdad para una tolerancia dada implica igualdad para una tolerancia más alta)
• si t1 < t2entonces a =~[Ɛ,t1] bimplica a =~[Ɛ,t2] b(igualdad para una imprecisión dada implica igualdad en una imprecisión mayor)

La solución propuesta también verifica estos.

Tuve el problema de comparar números de punto flotante A < By esto A > B es lo que parece funcionar:

if(A - B < Epsilon) && (fabs(A-B) > Epsilon)
{
    printf("A is less than B");
}

if (A - B > Epsilon) && (fabs(A-B) > Epsilon)
{
    printf("A is greater than B");
}

Las fábricas - valor absoluto - se encarga de si son esencialmente iguales.

Tenemos que elegir un nivel de tolerancia para comparar los números flotantes. Por ejemplo,

final float TOLERANCE = 0.00001;
if (Math.abs(f1 - f2) < TOLERANCE)
    Console.WriteLine("Oh yes!");

Una nota. Tu ejemplo es bastante divertido.

double a = 1.0 / 3.0;
double b = a + a + a;
if (a != b)
    Console.WriteLine("Oh no!");

Algunas matemáticas aquí

a = 1/3
b = 1/3 + 1/3 + 1/3 = 1.

1/3 != 1

Oh si..

Quieres decir

if (b != 1)
    Console.WriteLine("Oh no!")

Idea que tuve para la comparación de punto flotante en Swift.

infix operator ~= {}

func ~= (a: Float, b: Float) -> Bool {
    return fabsf(a - b) < Float(FLT_EPSILON)
}

func ~= (a: CGFloat, b: CGFloat) -> Bool {
    return fabs(a - b) < CGFloat(FLT_EPSILON)
}

func ~= (a: Double, b: Double) -> Bool {
    return fabs(a - b) < Double(FLT_EPSILON)
}

Adaptación a PHP de la respuesta de Michael Borgwardt y bosonix:

class Comparison
{
    const MIN_NORMAL = 1.17549435E-38;  //from Java Specs

    // from http://floating-point-gui.de/errors/comparison/
    public function nearlyEqual($a, $b, $epsilon = 0.000001)
    {
        $absA = abs($a);
        $absB = abs($b);
        $diff = abs($a - $b);

        if ($a == $b) {
            return true;
        } else {
            if ($a == 0 || $b == 0 || $diff < self::MIN_NORMAL) {
                return $diff < ($epsilon * self::MIN_NORMAL);
            } else {
                return $diff / ($absA + $absB) < $epsilon;
            }
        }
    }
}

Debería preguntarse por qué está comparando los números. Si conoce el propósito de la comparación, también debe conocer la precisión requerida de sus números. Eso es diferente en cada situación y en cada contexto de aplicación. Pero en casi todos los casos prácticos hay un requisito absoluto precisión . Es muy raro que se aplique una precisión relativa.

Para dar un ejemplo: si su objetivo es dibujar un gráfico en la pantalla, entonces probablemente desee que los valores de punto flotante se comparen iguales si se asignan al mismo píxel en la pantalla. Si el tamaño de su pantalla es de 1000 píxeles y sus números están en el rango de 1e6, entonces probablemente querrá que 100 se compare con 200.

Dada la precisión absoluta requerida, el algoritmo se convierte en:

public static ComparisonResult compare(float a, float b, float accuracy) 
{
    if (isnan(a) || isnan(b))   // if NaN needs to be supported
        return UNORDERED;    
    if (a == b)                 // short-cut and takes care of infinities
        return EQUAL;           
    if (abs(a-b) < accuracy)    // comparison wrt. the accuracy
        return EQUAL;
    if (a < b)                  // larger / smaller
        return SMALLER;
    else
        return LARGER;
}

El consejo estándar es usar un pequeño valor "épsilon" (elegido probablemente según su aplicación), y considerar que los flotantes que están dentro de épsilon entre sí son iguales. por ejemplo, algo como

#define EPSILON 0.00000001

if ((a - b) < EPSILON && (b - a) < EPSILON) {
  printf("a and b are about equal\n");
}

Una respuesta más completa es complicada, porque el error de punto flotante es extremadamente sutil y confuso para razonar. Si realmente le importa la igualdad en un sentido preciso, probablemente esté buscando una solución que no implique el punto flotante.

Intenté escribir una función de igualdad con los comentarios anteriores en mente. Esto es lo que se me ocurrió:

Editar: Cambiar de Math.Max ​​(a, b) a Math.Max ​​(Math.Abs ​​(a), Math.Abs ​​(b))

static bool fpEqual(double a, double b)
{
    double diff = Math.Abs(a - b);
    double epsilon = Math.Max(Math.Abs(a), Math.Abs(b)) * Double.Epsilon;
    return (diff < epsilon);
}

Pensamientos Todavía necesito calcular un mayor que y un menor que también.

Debe tener en cuenta que el error de truncamiento es relativo. Dos números son aproximadamente iguales si su diferencia es tan grande como su ulp (Unidad en el último lugar).

Sin embargo, si realiza cálculos de punto flotante, su potencial de error aumenta con cada operación (especialmente, ¡cuidado con las restas!), Por lo que su tolerancia de error debe aumentar en consecuencia.

La mejor manera de comparar dobles para la igualdad / desigualdad es tomando el valor absoluto de su diferencia y comparándolo con un valor lo suficientemente pequeño (dependiendo de su contexto).

double eps = 0.000000001; //for instance

double a = someCalc1();
double b = someCalc2();

double diff = Math.abs(a - b);
if (diff < eps) {
    //equal
}