Сегодня в обычном вопросе в комментариях в очередной раз всплыл вопрос, использовать ли float или double; мотивация при этом была та же, что возникала уже не раз: float должен работать быстрее.
Насколько это и в самом деле так? Стоит ли повышение скорости снижения точности?
Ставлю вычислительный эксперимент. Вычисление интеграла методом Симпсона с заранее недостижимой точностью, просто чтоб выполнить одинаковое количество вычислений. Специально не матрицы, чтобы всякие эффекты вроде размещения в кэше не играли особой роли. Код просто пережевывает значения с плавающей точкой; приведен ниже.
Компиляция VC++2022, опции командной строки, связанные с оптимизацией — /O2 /Ot /Oi /Ox /Oy /Ob2.
При компиляции в 64-разрадное приложение среднее (по пяти выполнениям) время выполнения на моей машине 2.35 с для float, 2.84 с для double. Получается, для стандартных, так сказать, настроек double на 20% медленнее.
Если использовать ключ /fp:fast, о котором говорится дословно результаты менее предсказуемы, то получим 0.65 с и 1.55 с (медленнее в 2.4 раза).
В случае 32-разрядного приложения разница — 3.90 с для float и 5.93 с для double, т.е. разница 50% (в полтора раза) и для /fp:fast 0.85 с и 1.55 с (в 1.8 раза).
Вопрос точности выполняемых вычислений остался за рамками эксперимента.
Выводы не делаю, просто чисто информативное сообщение. Каждый для своей программы (что важнее — скорость или точность/надежность) делает их сам. (Лично я для своих в основном научных расчетов вопрос давно решил в пользу точности.)
Интересно было бы сравнить результаты для других компиляторов.
Конструктивная критика принимается с благодарностью.
Вот код (суть по сути :) неважна, специально старался как можно больше использовать типы с плавающей точкой).
#include <iostream>
using namespace std;
//typedef float fptype;
typedef double fptype;
template<typename T, typename Func>
requires std::is_floating_point_v<T>
T Simpson(T a, T b, T eps, Func&& f)
{
T N = 4.;
T sum1 = 0., sum2 = 0.;
T x, h, index1, index2, func;
unsigned int i;
do {
h = (b-a)/N;
sum1 = sum2 = -(f(a)+f(b));
index1 = index2 = 4.;
i = 0;
for(x = a; x <= b+h/2.; x+=h) {
sum1 += (index1 = T{6.}-index1)*(func=f(x));
if (i++%2==0) {
sum2 += (index2 = T{6.}-index2)*func;
};
};
if (fabs((sum1-T{2.}*sum2)/sum1)<=eps)
return (h*sum1)/T{3.}+h*(T{2.}*sum2-sum1)/T{45.};
else N*= T{2.};
} while(N < (1<<20));
return (h*sum1)/T{3.};
};
fptype f(fptype x)
{
fptype s = 0;
for(int i = 0; i < 100; ++i) s += sin(x*i)*(cos(x*i)-sin(x*i));
fptype p = 1;
for(int i = 0; i < 10; ++i) p *= pow(abs(log(abs(x)+fptype(0.1))),fptype(1.)/(fptype(1.)+i));
p += s;
for(int i = 0; i < 10; ++i) p /= x+fptype(0.2);
return p + s;
}
int main(int argc, char * argv[])
{
cout << Simpson(fptype(0),fptype(100),fptype(1e-18),f);
}