Какой самый эффективный способ сравнения с плавающей запятой и двойного сравнения?
Что было бы наиболее эффективным способом сравнить два double или два float ценности?
Просто делать это не правильно
bool CompareDoubles1 (double A, double B)
{
return A == B;
}
Но что-то вроде:
bool CompareDoubles2 (double A, double B)
{
diff = A - B;
return (diff < EPSILON) && (-diff < EPSILON);
}
Кажется, переработка отходов.
Кто-нибудь знает умный поплавковый компаратор?
35 ответов
Будьте предельно осторожны, используя любые другие предложения. Все зависит от контекста.
Я провел много времени, отслеживая ошибки в системе, которая предполагала a==b если |a-b|<epsilon, Основными проблемами были:
Неявная презумпция в алгоритме, если
a==bа такжеb==cзатемa==c,Использование одного и того же эпсилона для линий, измеряемых в дюймах, и линий, измеряемых в милях (0,001 дюйма). То есть
a==bно1000a!=1000b, (Вот почему почти Equal2sComplement запрашивает epsilon или max ULPS).Использование одного и того же эпсилона как для косинуса углов, так и для длины линий!
Использование такой функции сравнения для сортировки элементов в коллекции. (В этом случае использование встроенного оператора C++ == для двойных значений дает правильные результаты.)
Как я уже сказал: все зависит от контекста и ожидаемого размера a а также b,
КСТАТИ, std::numeric_limits<double>::epsilon() это "машина эпсилон". Это разница между 1.0 и следующим значением, представленным двойным. Я предполагаю, что это можно использовать в функции сравнения, но только если ожидаемые значения меньше 1. (Это в ответ на ответ @ cdv...)
Кроме того, если у вас есть в основном int арифметика в doubles (здесь мы используем double для хранения значений int в некоторых случаях) Ваша арифметика будет правильной. Например, 4.0/2.0 будет таким же, как 1.0+1.0. Это до тех пор, пока вы не делаете вещи, которые приводят к дробным частям (4.0/3.0) или не выходят за пределы размера int.
Сравнение со значением эпсилона - это то, что делает большинство людей (даже в программировании игр).
Вы должны немного изменить свою реализацию:
bool AreSame(double a, double b)
{
return fabs(a - b) < EPSILON;
}
Изменить: Кристер добавил стопку полезной информации по этой теме в недавнем сообщении в блоге. Наслаждаться.
Сравнение чисел с плавающей запятой для зависит от контекста. Поскольку даже изменение порядка операций может привести к различным результатам, важно знать, насколько "равными" должны быть числа.
Сравнение чисел с плавающей запятой по Брюсу Доусону - это хорошее место для начала, если взглянуть на сравнение с плавающей запятой.
Следующие определения взяты из Искусства компьютерного программирования от Кнута:
bool approximatelyEqual(float a, float b, float epsilon)
{
return fabs(a - b) <= ( (fabs(a) < fabs(b) ? fabs(b) : fabs(a)) * epsilon);
}
bool essentiallyEqual(float a, float b, float epsilon)
{
return fabs(a - b) <= ( (fabs(a) > fabs(b) ? fabs(b) : fabs(a)) * epsilon);
}
bool definitelyGreaterThan(float a, float b, float epsilon)
{
return (a - b) > ( (fabs(a) < fabs(b) ? fabs(b) : fabs(a)) * epsilon);
}
bool definitelyLessThan(float a, float b, float epsilon)
{
return (b - a) > ( (fabs(a) < fabs(b) ? fabs(b) : fabs(a)) * epsilon);
}
Конечно, выбор эпсилона зависит от контекста и определяет, насколько равными должны быть числа.
Другой метод сравнения чисел с плавающей запятой - посмотреть на ULP (единицы на последнем месте) чисел. Несмотря на то, что статья не касается конкретно сравнений, статья " То, что должен знать каждый компьютерщик" о числах с плавающей запятой, является хорошим ресурсом для понимания того, как работает число с плавающей запятой и каковы подводные камни, включая ULP.
Я обнаружил, что Google C++ Testing Framework содержит красивую кроссплатформенную реализацию на основе шаблонов NearEqual2sComplement, которая работает как с двойными, так и с плавающими числами. Учитывая, что он выпущен по лицензии BSD, использование его в вашем собственном коде не должно быть проблемой, если вы сохраняете лицензию. Я извлек приведенный ниже код из http://code.google.com/p/googletest/source/browse/trunk/include/gtest/internal/gtest-internal.h https://github.com/google/googletest/blob/master/googletest/include/gtest/internal/gtest-internal.h и добавил лицензию сверху.
Обязательно #define GTEST_OS_WINDOWS определите какое-либо значение (или измените код, в котором он используется, на что-то, что соответствует вашей кодовой базе - в конце концов, это BSD-лицензия).
Пример использования:
double left = // something
double right = // something
const FloatingPoint<double> lhs(left), rhs(right);
if (lhs.AlmostEquals(rhs)) {
//they're equal!
}
Вот код:
// Copyright 2005, Google Inc.
// All rights reserved.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// * Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
// * Redistributions in binary form must reproduce the above
// copyright notice, this list of conditions and the following disclaimer
// in the documentation and/or other materials provided with the
// distribution.
// * Neither the name of Google Inc. nor the names of its
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
// OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
// SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
// LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
// THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
// OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Authors: wan@google.com (Zhanyong Wan), eefacm@gmail.com (Sean Mcafee)
//
// The Google C++ Testing Framework (Google Test)
// This template class serves as a compile-time function from size to
// type. It maps a size in bytes to a primitive type with that
// size. e.g.
//
// TypeWithSize<4>::UInt
//
// is typedef-ed to be unsigned int (unsigned integer made up of 4
// bytes).
//
// Such functionality should belong to STL, but I cannot find it
// there.
//
// Google Test uses this class in the implementation of floating-point
// comparison.
//
// For now it only handles UInt (unsigned int) as that's all Google Test
// needs. Other types can be easily added in the future if need
// arises.
template <size_t size>
class TypeWithSize {
public:
// This prevents the user from using TypeWithSize<N> with incorrect
// values of N.
typedef void UInt;
};
// The specialization for size 4.
template <>
class TypeWithSize<4> {
public:
// unsigned int has size 4 in both gcc and MSVC.
//
// As base/basictypes.h doesn't compile on Windows, we cannot use
// uint32, uint64, and etc here.
typedef int Int;
typedef unsigned int UInt;
};
// The specialization for size 8.
template <>
class TypeWithSize<8> {
public:
#if GTEST_OS_WINDOWS
typedef __int64 Int;
typedef unsigned __int64 UInt;
#else
typedef long long Int; // NOLINT
typedef unsigned long long UInt; // NOLINT
#endif // GTEST_OS_WINDOWS
};
// This template class represents an IEEE floating-point number
// (either single-precision or double-precision, depending on the
// template parameters).
//
// The purpose of this class is to do more sophisticated number
// comparison. (Due to round-off error, etc, it's very unlikely that
// two floating-points will be equal exactly. Hence a naive
// comparison by the == operation often doesn't work.)
//
// Format of IEEE floating-point:
//
// The most-significant bit being the leftmost, an IEEE
// floating-point looks like
//
// sign_bit exponent_bits fraction_bits
//
// Here, sign_bit is a single bit that designates the sign of the
// number.
//
// For float, there are 8 exponent bits and 23 fraction bits.
//
// For double, there are 11 exponent bits and 52 fraction bits.
//
// More details can be found at
// http://en.wikipedia.org/wiki/IEEE_floating-point_standard.
//
// Template parameter:
//
// RawType: the raw floating-point type (either float or double)
template <typename RawType>
class FloatingPoint {
public:
// Defines the unsigned integer type that has the same size as the
// floating point number.
typedef typename TypeWithSize<sizeof(RawType)>::UInt Bits;
// Constants.
// # of bits in a number.
static const size_t kBitCount = 8*sizeof(RawType);
// # of fraction bits in a number.
static const size_t kFractionBitCount =
std::numeric_limits<RawType>::digits - 1;
// # of exponent bits in a number.
static const size_t kExponentBitCount = kBitCount - 1 - kFractionBitCount;
// The mask for the sign bit.
static const Bits kSignBitMask = static_cast<Bits>(1) << (kBitCount - 1);
// The mask for the fraction bits.
static const Bits kFractionBitMask =
~static_cast<Bits>(0) >> (kExponentBitCount + 1);
// The mask for the exponent bits.
static const Bits kExponentBitMask = ~(kSignBitMask | kFractionBitMask);
// How many ULP's (Units in the Last Place) we want to tolerate when
// comparing two numbers. The larger the value, the more error we
// allow. A 0 value means that two numbers must be exactly the same
// to be considered equal.
//
// The maximum error of a single floating-point operation is 0.5
// units in the last place. On Intel CPU's, all floating-point
// calculations are done with 80-bit precision, while double has 64
// bits. Therefore, 4 should be enough for ordinary use.
//
// See the following article for more details on ULP:
// http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm.
static const size_t kMaxUlps = 4;
// Constructs a FloatingPoint from a raw floating-point number.
//
// On an Intel CPU, passing a non-normalized NAN (Not a Number)
// around may change its bits, although the new value is guaranteed
// to be also a NAN. Therefore, don't expect this constructor to
// preserve the bits in x when x is a NAN.
explicit FloatingPoint(const RawType& x) { u_.value_ = x; }
// Static methods
// Reinterprets a bit pattern as a floating-point number.
//
// This function is needed to test the AlmostEquals() method.
static RawType ReinterpretBits(const Bits bits) {
FloatingPoint fp(0);
fp.u_.bits_ = bits;
return fp.u_.value_;
}
// Returns the floating-point number that represent positive infinity.
static RawType Infinity() {
return ReinterpretBits(kExponentBitMask);
}
// Non-static methods
// Returns the bits that represents this number.
const Bits &bits() const { return u_.bits_; }
// Returns the exponent bits of this number.
Bits exponent_bits() const { return kExponentBitMask & u_.bits_; }
// Returns the fraction bits of this number.
Bits fraction_bits() const { return kFractionBitMask & u_.bits_; }
// Returns the sign bit of this number.
Bits sign_bit() const { return kSignBitMask & u_.bits_; }
// Returns true iff this is NAN (not a number).
bool is_nan() const {
// It's a NAN if the exponent bits are all ones and the fraction
// bits are not entirely zeros.
return (exponent_bits() == kExponentBitMask) && (fraction_bits() != 0);
}
// Returns true iff this number is at most kMaxUlps ULP's away from
// rhs. In particular, this function:
//
// - returns false if either number is (or both are) NAN.
// - treats really large numbers as almost equal to infinity.
// - thinks +0.0 and -0.0 are 0 DLP's apart.
bool AlmostEquals(const FloatingPoint& rhs) const {
// The IEEE standard says that any comparison operation involving
// a NAN must return false.
if (is_nan() || rhs.is_nan()) return false;
return DistanceBetweenSignAndMagnitudeNumbers(u_.bits_, rhs.u_.bits_)
<= kMaxUlps;
}
private:
// The data type used to store the actual floating-point number.
union FloatingPointUnion {
RawType value_; // The raw floating-point number.
Bits bits_; // The bits that represent the number.
};
// Converts an integer from the sign-and-magnitude representation to
// the biased representation. More precisely, let N be 2 to the
// power of (kBitCount - 1), an integer x is represented by the
// unsigned number x + N.
//
// For instance,
//
// -N + 1 (the most negative number representable using
// sign-and-magnitude) is represented by 1;
// 0 is represented by N; and
// N - 1 (the biggest number representable using
// sign-and-magnitude) is represented by 2N - 1.
//
// Read http://en.wikipedia.org/wiki/Signed_number_representations
// for more details on signed number representations.
static Bits SignAndMagnitudeToBiased(const Bits &sam) {
if (kSignBitMask & sam) {
// sam represents a negative number.
return ~sam + 1;
} else {
// sam represents a positive number.
return kSignBitMask | sam;
}
}
// Given two numbers in the sign-and-magnitude representation,
// returns the distance between them as an unsigned number.
static Bits DistanceBetweenSignAndMagnitudeNumbers(const Bits &sam1,
const Bits &sam2) {
const Bits biased1 = SignAndMagnitudeToBiased(sam1);
const Bits biased2 = SignAndMagnitudeToBiased(sam2);
return (biased1 >= biased2) ? (biased1 - biased2) : (biased2 - biased1);
}
FloatingPointUnion u_;
};
РЕДАКТИРОВАТЬ: Этому посту 4 года. Это, вероятно, все еще действует, и код хорош, но некоторые люди нашли улучшения. Лучше иди получить последнюю версию AlmostEquals прямо из исходного кода Google Test, а не того, который я вставил сюда.
Для более глубокого подхода прочитайте Сравнение чисел с плавающей запятой. Вот фрагмент кода по этой ссылке:
// Usable AlmostEqual function
bool AlmostEqual2sComplement(float A, float B, int maxUlps)
{
// Make sure maxUlps is non-negative and small enough that the
// default NAN won't compare as equal to anything.
assert(maxUlps > 0 && maxUlps < 4 * 1024 * 1024);
int aInt = *(int*)&A;
// Make aInt lexicographically ordered as a twos-complement int
if (aInt < 0)
aInt = 0x80000000 - aInt;
// Make bInt lexicographically ordered as a twos-complement int
int bInt = *(int*)&B;
if (bInt < 0)
bInt = 0x80000000 - bInt;
int intDiff = abs(aInt - bInt);
if (intDiff <= maxUlps)
return true;
return false;
}
Понимание того, что это старая ветка, но эта статья является одной из самых простых, которые я обнаружил при сравнении чисел с плавающей запятой, и, если вы хотите узнать больше, она также содержит более подробные ссылки, а основной сайт охватывает полный спектр вопросов. Работа с числами с плавающей точкой Руководство по плавающей точке: Сравнение.
Мы можем найти несколько более практичную статью в разделе "Допуски с плавающей точкой", где отмечается, что есть тест абсолютной толерантности, который сводится к этому в C++:
bool absoluteToleranceCompare(double x, double y)
{
return std::fabs(x - y) <= std::numeric_limits<double>::epsilon() ;
}
и тест относительной толерантности:
bool relativeToleranceCompare(double x, double y)
{
double maxXY = std::max( std::fabs(x) , std::fabs(y) ) ;
return std::fabs(x - y) <= std::numeric_limits<double>::epsilon()*maxXY ;
}
В статье отмечается, что абсолютный тест не проходит, когда x а также y большие и терпят неудачу в относительном случае, когда они маленькие. Если предположить, что абсолютная и относительная толерантность одинаковы, комбинированный тест будет выглядеть так:
bool combinedToleranceCompare(double x, double y)
{
double maxXYOne = std::max( { 1.0, std::fabs(x) , std::fabs(y) } ) ;
return std::fabs(x - y) <= std::numeric_limits<double>::epsilon()*maxXYOne ;
}
Портативный способ получить эпсилон в C++
#include <limits>
std::numeric_limits<double>::epsilon()
Тогда функция сравнения становится
#include <cmath>
#include <limits>
bool AreSame(double a, double b) {
return std::fabs(a - b) < std::numeric_limits<double>::epsilon();
}
Я потратил довольно много времени на изучение материала в этой замечательной теме. Я сомневаюсь, что все хотят тратить так много времени, поэтому я бы выделил краткое изложение того, что я узнал, и решения, которое я реализовал.
Краткое резюме
- Есть две проблемы с плавающим сравнением: у вас ограниченная точность, и значение "приблизительно ноль" зависит от контекста (см. Следующий пункт).
- 1E-8 примерно такой же, как 1E-16? Если вы смотрите на шумные данные датчиков, то, вероятно, да, но если вы делаете молекулярное моделирование, то, возможно, нет! Итог: вам всегда нужно думать о значении толерантности в контексте конкретного вызова функции, а не просто делать его общей жестко закодированной константой приложения.
- Для общих библиотечных функций по-прежнему приятно иметь параметр с допуском по умолчанию. Типичный выбор
numeric_limits::epsilon()который совпадает с FLT_EPSILON в float.h. Это, однако, проблематично, потому что epsilon для сравнения значений, таких как 1.0, если не совпадает с epsilon для значений, таких как 1E9. FLT_EPSILON определен для 1.0. - Очевидная реализация для проверки, находится ли число в пределах допуска
fabs(a-b) <= epsilonоднако это не работает, потому что эпсилон по умолчанию определен для 1.0. Нам нужно масштабировать эпсилон вверх или вниз с точки зрения a и b. - Есть два решения этой проблемы: либо вы устанавливаете эпсилон пропорционально
max(a,b)или вы можете получить следующие представимые числа вокруг a, а затем посмотреть, попадает ли b в этот диапазон. Первый называется "относительным" методом, а позже - методом ULP. - Оба метода на самом деле дают сбой при сравнении с 0. В этом случае приложение должно предоставить правильный допуск.
Реализация служебных функций (C++11)
//implements relative method - do not use for comparing with zero
//use this most of the time, tolerance needs to be meaningful in your context
template<typename TReal>
static bool isApproximatelyEqual(TReal a, TReal b, TReal tolerance = std::numeric_limits<TReal>::epsilon())
{
TReal diff = std::fabs(a - b);
if (diff <= tolerance)
return true;
if (diff < std::fmax(std::fabs(a), std::fabs(b)) * tolerance)
return true;
return false;
}
//supply tolerance that is meaningful in your context
//for example, default tolerance may not work if you are comparing double with float
template<typename TReal>
static bool isApproximatelyZero(TReal a, TReal tolerance = std::numeric_limits<TReal>::epsilon())
{
if (std::fabs(a) <= tolerance)
return true;
return false;
}
//use this when you want to be on safe side
//for example, don't start rover unless signal is above 1
template<typename TReal>
static bool isDefinitelyLessThan(TReal a, TReal b, TReal tolerance = std::numeric_limits<TReal>::epsilon())
{
TReal diff = a - b;
if (diff < tolerance)
return true;
if (diff < std::fmax(std::fabs(a), std::fabs(b)) * tolerance)
return true;
return false;
}
template<typename TReal>
static bool isDefinitelyGreaterThan(TReal a, TReal b, TReal tolerance = std::numeric_limits<TReal>::epsilon())
{
TReal diff = a - b;
if (diff > tolerance)
return true;
if (diff > std::fmax(std::fabs(a), std::fabs(b)) * tolerance)
return true;
return false;
}
//implements ULP method
//use this when you are only concerned about floating point precision issue
//for example, if you want to see if a is 1.0 by checking if its within
//10 closest representable floating point numbers around 1.0.
template<typename TReal>
static bool isWithinPrecisionInterval(TReal a, TReal b, unsigned int interval_size = 1)
{
TReal min_a = a - (a - std::nextafter(a, std::numeric_limits<TReal>::lowest())) * interval_size;
TReal max_a = a + (std::nextafter(a, std::numeric_limits<TReal>::max()) - a) * interval_size;
return min_a <= b && max_a >= b;
}
Код, который вы написали, содержит ошибки:
return (diff < EPSILON) && (-diff > EPSILON);
Правильный код будет:
return (diff < EPSILON) && (diff > -EPSILON);
(... и да, это другое)
Интересно, не заставят ли вас в какой-то степени похабы ленивые оценки? Я бы сказал, что это зависит от компилятора. Вы можете попробовать оба. Если они в среднем эквивалентны, возьмите реализацию с потрясающими.
Если у вас есть какая-то информация о том, какой из двух чисел с большей вероятностью будет больше, чем другие, вы можете сыграть в порядке сравнения, чтобы лучше использовать ленивую оценку.
Наконец, вы можете получить лучший результат, вставив эту функцию. Хотя вряд ли что-то улучшится...
Редактировать: OJ, спасибо за исправление вашего кода. Я стер свой комментарий соответственно
`вернуть потрясающие (а - б)
Это хорошо, если:
- порядок величины ваших входов не сильно меняется
- очень небольшое количество противоположных знаков можно рассматривать как равные
Но в противном случае это приведет к неприятностям. Числа двойной точности имеют разрешение около 16 десятичных знаков. Если сравниваемые два числа больше по величине, чем EPSILON*1.0E16, то вы также можете сказать:
return a==b;
Я рассмотрю другой подход, который предполагает, что вам нужно беспокоиться о первой проблеме, и предположим, что вторая подойдет для вашего приложения. Решение будет что-то вроде:
#define VERYSMALL (1.0E-150)
#define EPSILON (1.0E-8)
bool AreSame(double a, double b)
{
double absDiff = fabs(a - b);
if (absDiff < VERYSMALL)
{
return true;
}
double maxAbs = max(fabs(a) - fabs(b));
return (absDiff/maxAbs) < EPSILON;
}
Это дорого в вычислительном отношении, но иногда это то, что требуется. Это то, что мы должны сделать в моей компании, потому что мы имеем дело с инженерной библиотекой, и входные данные могут варьироваться на несколько десятков порядков.
В любом случае, суть заключается в следующем (и относится практически ко всем программным проблемам): оцените свои потребности, а затем найдите решение, отвечающее вашим потребностям - не думайте, что простой ответ удовлетворит ваши потребности. Если после вашей оценки вы обнаружите, что fabs(a-b) < EPSILON хватит, отлично - пользуйся! Но помните о его недостатках и других возможных решениях.
Как отмечали другие, использование эпсилона с фиксированной экспонентой (например, 0,0000001) будет бесполезным для значений, отличных от значения эпсилона. Например, если вашими двумя значениями являются 10000.000977 и 10000, то между этими двумя числами НЕТ 32-битных значений с плавающей запятой - 10000 и 10000.000977 настолько близки, насколько это возможно, не будучи одинаковыми по битам. Здесь эпсилон менее 0,0009 не имеет смысла; Вы также можете использовать оператор прямого равенства.
Аналогично, когда эти два значения приближаются к размеру эпсилона, относительная ошибка возрастает до 100%.
Таким образом, попытка смешать число с фиксированной точкой, например 0,00001, со значениями с плавающей точкой (где показатель степени произвольный) является бессмысленным упражнением. Это будет работать только в том случае, если вы можете быть уверены, что значения операндов находятся в узкой области (т. Е. Близко к некоторому конкретному показателю степени), и если вы правильно выберете значение эпсилона для этого конкретного теста. Если вы вытаскиваете число из воздуха ("Эй! 0,00001 - это мало, значит, это должно быть хорошо!"), Вы обречены на ошибки в цифрах. Я потратил много времени на отладку плохого числового кода, где некоторые плохие чмо бросают в случайные значения epsilon, чтобы заставить работать еще один тестовый пример.
Если вы занимаетесь каким-либо числовым программированием и считаете, что вам нужно получить эпсилоны с фиксированной запятой, ПРОЧИТАЙТЕ СТАТЬЮ БРЮСА ПО СРАВНЕНИЮ Плавающих чисел.
Qt реализует две функции, можете ли вы поучиться у них:
static inline bool qFuzzyCompare(double p1, double p2)
{
return (qAbs(p1 - p2) <= 0.000000000001 * qMin(qAbs(p1), qAbs(p2)));
}
static inline bool qFuzzyCompare(float p1, float p2)
{
return (qAbs(p1 - p2) <= 0.00001f * qMin(qAbs(p1), qAbs(p2)));
}
И вам могут понадобиться следующие функции, так как
Обратите внимание, что сравнение значений, где p1 или p2 равно 0.0, не будет работать, равно как и сравнение значений, где одно из значений - NaN или бесконечность. Если одно из значений всегда 0.0, используйте взамен qFuzzyIsNull. Если одно из значений может составлять 0,0, одним из решений является добавление 1,0 к обоим значениям.
static inline bool qFuzzyIsNull(double d)
{
return qAbs(d) <= 0.000000000001;
}
static inline bool qFuzzyIsNull(float f)
{
return qAbs(f) <= 0.00001f;
}
Вот доказательство того, что с помощью std::numeric_limits::epsilon() это не ответ - он терпит неудачу для значений больше единицы:
Доказательство моего комментария выше:
#include <stdio.h>
#include <limits>
double ItoD (__int64 x) {
// Return double from 64-bit hexadecimal representation.
return *(reinterpret_cast<double*>(&x));
}
void test (__int64 ai, __int64 bi) {
double a = ItoD(ai), b = ItoD(bi);
bool close = std::fabs(a-b) < std::numeric_limits<double>::epsilon();
printf ("%.16f and %.16f %s close.\n", a, b, close ? "are " : "are not");
}
int main()
{
test (0x3fe0000000000000L,
0x3fe0000000000001L);
test (0x3ff0000000000000L,
0x3ff0000000000001L);
}
Запуск дает этот вывод:
0.5000000000000000 and 0.5000000000000001 are close.
1.0000000000000000 and 1.0000000000000002 are not close.
Обратите внимание, что во втором случае (один и чуть больше единицы) два входных значения находятся настолько близко, насколько это возможно, и все же сравниваются как не близкие. Таким образом, для значений больше 1.0 вы можете просто использовать тест на равенство. Исправленные эпсилоны не спасут вас при сравнении значений с плавающей точкой.
К сожалению, даже ваш "расточительный" код неверен. EPSILON - это наименьшее значение, которое можно добавить к 1,0 и изменить его значение. Значение 1.0 очень важно - большие числа не меняются при добавлении в EPSILON. Теперь вы можете масштабировать это значение до числа, которое вы сравниваете, чтобы определить, отличаются они или нет. Правильное выражение для сравнения двух чисел:
if (fabs(a - b) <= DBL_EPSILON * fmax(fabs(a), fabs(b)))
{
// ...
}
Это как минимум. В целом, однако, вы бы хотели учесть шум в ваших вычислениях и игнорировать несколько наименее значимых битов, поэтому более реалистичное сравнение будет выглядеть так:
if (fabs(a - b) <= 16 * DBL_EPSILON * fmax(fabs(a), fabs(b)))
{
// ...
}
Если для вас очень важна эффективность сравнения и вы знаете диапазон своих значений, вам следует вместо этого использовать числа с фиксированной запятой.
Вы должны выполнить эту обработку для сравнения с плавающей запятой, так как float нельзя точно сравнить, как целочисленные типы. Вот функции для различных операторов сравнения.
С плавающей точкой равно (==)
Я также предпочитаю технику вычитания, а не полагаться на
fabs() или же
abs(), но мне пришлось бы ускорить его профилирование на различных архитектурах от 64-битного ПК до микроконтроллера ATMega328 (Arduino), чтобы действительно увидеть, сильно ли это влияет на производительность.
Итак, давайте забудем обо всем этом, что касается абсолютных значений, и просто проведем вычитание и сравнение!
Изменено из примера Microsoft здесь:
/// @brief See if two floating point numbers are approximately equal.
/// @param[in] a number 1
/// @param[in] b number 2
/// @param[in] epsilon A small value such that if the difference between the two numbers is
/// smaller than this they can safely be considered to be equal.
/// @return true if the two numbers are approximately equal, and false otherwise
bool is_float_eq(float a, float b, float epsilon) {
return ((a - b) < epsilon) && ((b - a) < epsilon);
}
bool is_double_eq(double a, double b, double epsilon) {
return ((a - b) < epsilon) && ((b - a) < epsilon);
}
Пример использования:
constexpr float EPSILON = 0.0001; // 1e-4
is_float_eq(1.0001, 0.99998, EPSILON);
Я не совсем уверен, но мне кажется, что некоторые из критических замечаний по поводу подхода на основе эпсилон, как описано в комментариях ниже этого высоко оцененного ответа, могут быть решены с помощью переменной эпсилон, масштабируемой в соответствии с плавающей запятой. сравниваются значения, например:
float a = 1.0001;
float b = 0.99998;
float epsilon = std::max(std::fabs(a), std::fabs(b)) * 1e-4;
is_float_eq(a, b, epsilon);
Таким образом, эпсилон-значение масштабируется вместе со значениями с плавающей запятой и поэтому никогда не бывает настолько маленьким, что становится незначительным.
Для полноты картины добавим остальное:
Лучше чем (>) и меньше (<):
/// @brief See if floating point number `a` is > `b`
/// @param[in] a number 1
/// @param[in] b number 2
/// @param[in] epsilon a small value such that if `a` is > `b` by this amount, `a` is considered
/// to be definitively > `b`
/// @return true if `a` is definitively > `b`, and false otherwise
bool is_float_gt(float a, float b, float epsilon) {
return a > b + epsilon;
}
bool is_double_gt(double a, double b, double epsilon) {
return a > b + epsilon;
}
/// @brief See if floating point number `a` is < `b`
/// @param[in] a number 1
/// @param[in] b number 2
/// @param[in] epsilon a small value such that if `a` is < `b` by this amount, `a` is considered
/// to be definitively < `b`
/// @return true if `a` is definitively < `b`, and false otherwise
bool is_float_lt(float a, float b, float epsilon) {
return a < b - epsilon;
}
bool is_double_lt(double a, double b, double epsilon) {
return a < b - epsilon;
}
Больше или равно (>=), и меньше или равно (<=)
/// @brief Returns true if `a` is definitively >= `b`, and false otherwise
bool is_float_ge(float a, float b, float epsilon) {
return a > b - epsilon;
}
bool is_double_ge(double a, double b, double epsilon) {
return a > b - epsilon;
}
/// @brief Returns true if `a` is definitively <= `b`, and false otherwise
bool is_float_le(float a, float b, float epsilon) {
return a < b + epsilon;
}
bool is_double_le(double a, double b, double epsilon) {
return a < b + epsilon;
}
Смотрите также:
- Макро-формы некоторых из вышеперечисленных функций в моем репо здесь: utilities.h.
- ОБНОВЛЕНИЕ 29 НОЯБРЯ 2020: это незавершенная работа, и я собираюсь дать ему отдельный ответ, когда будет готов, но я создал лучшую, масштабированную эпсилон-версию всех функций на C в этом файле здесь: utilities.c. Взглянуть.
- Дополнительная литература, которую мне
нужно сделатьсейчас: пересмотр допусков с плавающей запятой, Кристер Эриксон
Общее сравнение чисел с плавающей запятой, как правило, не имеет смысла. Как сравнивать, действительно зависит от поставленной задачи. Во многих задачах числа достаточно дискретизированы, чтобы их можно было сравнивать с заданным допуском. К сожалению, есть столько же проблем, когда такой трюк не работает. В качестве одного примера рассмотрим работу с функцией Хевисайда (шага) рассматриваемого числа (на ум приходят цифровые акции), когда ваши наблюдения находятся очень близко к барьеру. Выполнение сравнения на основе допусков не принесло бы много пользы, поскольку оно фактически сдвинет проблему с исходного барьера на два новых. Опять же, не существует универсального решения для таких проблем, и конкретное решение может потребовать даже изменения численного метода для достижения устойчивости.
Мой класс основан на ранее опубликованных ответах. Очень похоже на код Google, но я использую смещение, при котором все значения NaN превышают 0xFF000000. Это позволяет быстрее проверить NaN.
Этот код предназначен для демонстрации концепции, а не для общего решения. Код Google уже показывает, как вычислить все значения, специфичные для платформы, и я не хотел дублировать все это. Я провел ограниченное тестирование этого кода.
typedef unsigned int U32;
// Float Memory Bias (unsigned)
// ----- ------ ---------------
// NaN 0xFFFFFFFF 0xFF800001
// NaN 0xFF800001 0xFFFFFFFF
// -Infinity 0xFF800000 0x00000000 ---
// -3.40282e+038 0xFF7FFFFF 0x00000001 |
// -1.40130e-045 0x80000001 0x7F7FFFFF |
// -0.0 0x80000000 0x7F800000 |--- Valid <= 0xFF000000.
// 0.0 0x00000000 0x7F800000 | NaN > 0xFF000000
// 1.40130e-045 0x00000001 0x7F800001 |
// 3.40282e+038 0x7F7FFFFF 0xFEFFFFFF |
// Infinity 0x7F800000 0xFF000000 ---
// NaN 0x7F800001 0xFF000001
// NaN 0x7FFFFFFF 0xFF7FFFFF
//
// Either value of NaN returns false.
// -Infinity and +Infinity are not "close".
// -0 and +0 are equal.
//
class CompareFloat{
public:
union{
float m_f32;
U32 m_u32;
};
static bool CompareFloat::IsClose( float A, float B, U32 unitsDelta = 4 )
{
U32 a = CompareFloat::GetBiased( A );
U32 b = CompareFloat::GetBiased( B );
if ( (a > 0xFF000000) || (b > 0xFF000000) )
{
return( false );
}
return( (static_cast<U32>(abs( a - b ))) < unitsDelta );
}
protected:
static U32 CompareFloat::GetBiased( float f )
{
U32 r = ((CompareFloat*)&f)->m_u32;
if ( r & 0x80000000 )
{
return( ~r - 0x007FFFFF );
}
return( r + 0x7F800000 );
}
};
Нашел еще одну интересную реализацию: https://en.cppreference.com/w/cpp/types/numeric_limits/epsilon
#include <cmath>
#include <limits>
#include <iomanip>
#include <iostream>
#include <type_traits>
#include <algorithm>
template<class T>
typename std::enable_if<!std::numeric_limits<T>::is_integer, bool>::type
almost_equal(T x, T y, int ulp)
{
// the machine epsilon has to be scaled to the magnitude of the values used
// and multiplied by the desired precision in ULPs (units in the last place)
return std::fabs(x-y) <= std::numeric_limits<T>::epsilon() * std::fabs(x+y) * ulp
// unless the result is subnormal
|| std::fabs(x-y) < std::numeric_limits<T>::min();
}
int main()
{
double d1 = 0.2;
double d2 = 1 / std::sqrt(5) / std::sqrt(5);
std::cout << std::fixed << std::setprecision(20)
<< "d1=" << d1 << "\nd2=" << d2 << '\n';
if(d1 == d2)
std::cout << "d1 == d2\n";
else
std::cout << "d1 != d2\n";
if(almost_equal(d1, d2, 2))
std::cout << "d1 almost equals d2\n";
else
std::cout << "d1 does not almost equal d2\n";
}
Я бы очень настороженно отнесся к любому из этих ответов, который включает вычитание с плавающей запятой (например, fabs(ab)
Несмотря на то, что ответ не переносимый, я думаю, что ответ Грома лучше всего помогает избежать этих проблем.
Это зависит от того, насколько точно вы хотите, чтобы сравнение было. Если вы хотите сравнить для того же числа, просто перейдите с ==. (Вы почти никогда не захотите делать это, если вы действительно хотите точно такой же номер.) На любой приличной платформе вы также можете сделать следующее:
diff= a - b; return fabs(diff)<EPSILON;
как fabs имеет тенденцию быть довольно быстрым Под довольно быстрым я подразумеваю, что это в основном побитовое И, так что лучше быть быстрым.
Целочисленные трюки для сравнения значений типа double и float хороши, но, как правило, затрудняют эффективную обработку различными конвейерами ЦП. И в наши дни это определенно не быстрее на определенных архитектурах по порядку из-за использования стека в качестве временного хранилища для значений, которые часто используются. (Load-hit-store для тех, кому не все равно.)
Я использую следующую функцию для сравнения чисел с плавающей запятой:
bool approximatelyEqual(double a, double b)
{
return fabs(a - b) <= ((fabs(a) < fabs(b) ? fabs(b) : fabs(a)) * std::numeric_limits<double>::epsilon());
}
На самом деле в числовом программном обеспечении есть случаи, когда вы хотите проверить, равны ли два числа с плавающей запятой. Я отправил это на подобный вопрос
/questions/4914999/sravnenie-ieee-s-plavayuschej-i-dvojnoj-dlya-ravenstva/4915007#4915007
Таким образом, вы не можете сказать, что "CompareDoubles1" является неправильным в целом.
Я придумал простой подход к настройке размера эпсилона в соответствии с размером сравниваемых чисел. Итак, вместо использования:
iif(abs(a - b) < 1e-6, "equal", "not")
если
iif(abs(a - b) < (10 ^ -abs(7 - log(a))), "equal", "not")
Я полагаю, что это не удовлетворяет все теоретические вопросы, обсуждаемые в других ответах, но у него есть то преимущество, что он представляет собой одну строку кода, поэтому его можно использовать в формуле Excel или запросе доступа без необходимости использования функции VBA.
Я поискал, использовали ли другие этот метод, но ничего не нашел. Я протестировал его в своем приложении, и, похоже, он работает хорошо. Таким образом, кажется, что это метод, который подходит для контекстов, которые не требуют сложности других ответов. Но мне интересно, есть ли у него проблема, о которой я не задумывался, потому что, похоже, никто другой ею не пользуется.
Если есть причина, по которой тест с журналом не подходит для простого сравнения чисел разного размера, пожалуйста, укажите почему в комментарии.
Я использую этот код. В отличие от приведенных выше ответов, это позволяет дать объяснение в комментариях к коду.
В первой версии сравниваются комплексные числа, так что ошибку можно объяснить с точки зрения угла между двумя «векторами» одинаковой длины в комплексной плоскости (что дает некоторое представление). Тогда оттуда следует правильная формула для двух действительных чисел.
https://github.com/CarloWood/ai-utils/blob/master/almost_equal.h
Последний тогда
template<class T>
typename std::enable_if<std::is_floating_point<T>::value, bool>::type
almost_equal(T x, T y, T const abs_relative_error)
{
return 2 * std::abs(x - y) <= abs_relative_error * std::abs(x + y);
}
где в основном (дважды) абсолютное значение того, что ближе всего к определению в литературе: относительная ошибка. Но это только выбор имени.
Что это на самом деле, наиболее четко видно в комплексной плоскости, я думаю. Если |х| = 1, а y лежит в окружности вокруг x диаметром
abs_relative_error, то они считаются равными.
Я использую этот код:
bool AlmostEqual(double v1, double v2)
{
return (std::fabs(v1 - v2) < std::fabs(std::min(v1, v2)) * std::numeric_limits<double>::epsilon());
}
По шкале количеств:
Если epsilon является малой долей величины величины (т.е. относительной величины) в некотором определенном физическом смысле и A а также B типы сопоставимы в том же смысле, чем я думаю, что вполне корректно следующее:
#include <limits>
#include <iomanip>
#include <iostream>
#include <cmath>
#include <cstdlib>
#include <cassert>
template< typename A, typename B >
inline
bool close_enough(A const & a, B const & b,
typename std::common_type< A, B >::type const & epsilon)
{
using std::isless;
assert(isless(0, epsilon)); // epsilon is a part of the whole quantity
assert(isless(epsilon, 1));
using std::abs;
auto const delta = abs(a - b);
auto const x = abs(a);
auto const y = abs(b);
// comparable generally and |a - b| < eps * (|a| + |b|) / 2
return isless(epsilon * y, x) && isless(epsilon * x, y) && isless((delta + delta) / (x + y), epsilon);
}
int main()
{
std::cout << std::boolalpha << close_enough(0.9, 1.0, 0.1) << std::endl;
std::cout << std::boolalpha << close_enough(1.0, 1.1, 0.1) << std::endl;
std::cout << std::boolalpha << close_enough(1.1, 1.2, 0.01) << std::endl;
std::cout << std::boolalpha << close_enough(1.0001, 1.0002, 0.01) << std::endl;
std::cout << std::boolalpha << close_enough(1.0, 0.01, 0.1) << std::endl;
return EXIT_SUCCESS;
}
/// testing whether two doubles are almost equal. We consider two doubles
/// equal if the difference is within the range [0, epsilon).
///
/// epsilon: a positive number (supposed to be small)
///
/// if either x or y is 0, then we are comparing the absolute difference to
/// epsilon.
/// if both x and y are non-zero, then we are comparing the relative difference
/// to epsilon.
bool almost_equal(double x, double y, double epsilon)
{
double diff = x - y;
if (x != 0 && y != 0){
diff = diff/y;
}
if (diff < epsilon && -1.0*diff < epsilon){
return true;
}
return false;
}
Я использовал эту функцию для своего небольшого проекта, и она работает, но обратите внимание на следующее:
Ошибка двойной точности может вас удивить. Допустим, epsilon = 1.0e-6, тогда 1.0 и 1.000001 НЕ должны считаться равными в соответствии с приведенным выше кодом, но на моей машине функция считает их равными, потому что 1.000001 не может быть точно переведен в двоичный формат, это, вероятно, 1.0000009xxx. Я проверяю это с 1.0 и 1.0000011, и на этот раз я получаю ожидаемый результат.
Как насчет этого?
template<typename T>
bool FloatingPointEqual( T a, T b ) { return !(a < b) && !(b < a); }
Я видел различные подходы, но никогда не видел этого, поэтому мне тоже любопытно услышать какие-либо комментарии!
Это еще одно решение с лямбда:
#include <cmath>
#include <limits>
auto Compare = [](float a, float b, float epsilon = std::numeric_limits<float>::epsilon()){ return (std::fabs(a - b) <= epsilon); };
Я пишу это для Java, но, возможно, вы найдете это полезным. Он использует длинные вместо двойных, но заботится о NaN, субнормальных и т. Д.
public static boolean equal(double a, double b) {
final long fm = 0xFFFFFFFFFFFFFL; // fraction mask
final long sm = 0x8000000000000000L; // sign mask
final long cm = 0x8000000000000L; // most significant decimal bit mask
long c = Double.doubleToLongBits(a), d = Double.doubleToLongBits(b);
int ea = (int) (c >> 52 & 2047), eb = (int) (d >> 52 & 2047);
if (ea == 2047 && (c & fm) != 0 || eb == 2047 && (d & fm) != 0) return false; // NaN
if (c == d) return true; // identical - fast check
if (ea == 0 && eb == 0) return true; // ±0 or subnormals
if ((c & sm) != (d & sm)) return false; // different signs
if (abs(ea - eb) > 1) return false; // b > 2*a or a > 2*b
d <<= 12; c <<= 12;
if (ea < eb) c = c >> 1 | sm;
else if (ea > eb) d = d >> 1 | sm;
c -= d;
return c < 65536 && c > -65536; // don't use abs(), because:
// There is a posibility c=0x8000000000000000 which cannot be converted to positive
}
public static boolean zero(double a) { return (Double.doubleToLongBits(a) >> 52 & 2047) < 3; }
Имейте в виду, что после ряда операций с плавающей точкой число может сильно отличаться от того, что мы ожидаем. Нет кода, чтобы это исправить.