Блог

 #include <tgmath.h>
#include <iostream>
int main(int argc, char** argv) {

        #define NUM1 -0.031679909079365576
        #define NUM2 -0.11491794452567111

        std::cout << "double precision :"<< std::endl;
        typedef std::numeric_limits< double > dbl;
        std::cout.precision(dbl::max_digits10);
        std::cout << std::hypot((double)NUM1, (double)NUM2);
        std::cout << " VS sqrt :" << sqrt((double )NUM1*(double )NUM1 
                                    (double )NUM2*(double )NUM2) << std::endl;

        std::cout << "long double precision :"<< std::endl;
        typedef std::numeric_limits<long double > ldbl;
        std::cout.precision(ldbl::max_digits10);
        std::cout << std::hypot((long double)NUM1, (long double)NUM2);
        std::cout << " VS sqrt :" << sqrt((long double )NUM1*(long double )NUM1   (long double )NUM2*(long double )NUM2);
}
 

Возвращает под Linux (Ubuntu 18.04 clang или gcc, любая оптимизация, glic 2.25):

двойная точность : 0,1192046585217293 ПРОТИВ sqrt : 0,1192046585217293 2

длинная двойная точность : 0,119204658521729311251 ПРОТИВ sqrt : 0,119204658521729311251

В соответствии с cppreference :

Реализации обычно гарантируют точность менее 1 ulp (единицы измерения в последнюю очередь): GNU, BSD, Open64 std::hypot(x, y) эквивалентно std::abs(std::complex(x,y)) POSIX указывает, что недостаточный поток может возникать только тогда, когда оба аргументасубнормальный и правильный результат также субнормальный (это запрещает наивные реализации)

Итак, hypot((double)NUM1, (double)NUM2) должен возвращать 0.11920465852172932, я полагаю (как наивная реализация sqrt). В Windows, использующей 64-разрядную версию MSVC, это так.

Почему мы видим эту разницу , используя glibc ? Как можно решить эту несогласованность ?

• 0.1192046585217293 2 представлено 0x1.e84324de1b576p-4 (как двойное)
• 0.1192046585217293 0 представлено 0x1.e84324de1b575p-4 (как двойное)
• 0.1192046585217293 1 1251 — это результат с длинным двойным значением, который, как мы можем предположить, верен с точностью до пары десятичных знаков. т.е. Точный результат ближе к округленному результату.

Эти битовые шаблоны FP отличаются только младшим битом мантиссы (он же significand), и точный результат находится между ними. Таким образом, каждый из них имеет ошибку округления менее 1 ulp, достигая того, к чему стремятся типичные реализации (включая glibc).

В отличие от «базовых» операций IEEE-754 (add / sub / mul / div / sqrt), hypot не требуется «правильно округлять». Это означает <= 0,5 ulp ошибки. Достижение этого было бы намного медленнее для операций, которые HW не предоставляет напрямую. (например, выполняйте вычисления с повышенной точностью, по крайней мере, с парой дополнительных определенно правильных битов, чтобы вы могли округлить до ближайшего double, как это делает аппаратное обеспечение для базовых операций)

Случается, что в этом случае наивный метод вычисления выдал правильно округленный результат, в то время как «безопасная» реализация glibc std::hypot (которая должна избегать недостаточного потока при возведении в квадрат небольших чисел перед добавлением) выдала результат с > 0,5, но <1 ulp ошибки.

Вы не указали, используете ли вы MSVC в 32-разрядном режиме.

Предположительно, 32-разрядный режим будет использовать x87 для FP-математики, что дает дополнительную временную точность. Хотя CRT-код некоторых версий MSVC устанавливает внутреннюю точность FPU x87 для округления до 53-битной мантиссы после каждой операции, поэтому он ведет себя как SSE2, используя actual double , за исключением более широкого диапазона экспонент. См. Сообщение в блоге Брюса Доусона.

Так что я не знаю, есть ли какая-либо причина, помимо удачи, что MSVC std::hypot получил для этого правильно округленный результат.

Обратите внимание, что long double в MSVC используется тот же тип, что и 64-разрядный double ; эта реализация C не предоставляет 80-разрядный аппаратный тип расширенной точности x86 / x86-64. (64-разрядная мантисса).