Хитрости для получения быстрых квадратных корней и тому подобного достигают своей производительности, жертвуя точностью.(Ну, большинство из них.)
Вы уверены, что вам нужна точность double?Вы можете пожертвовать точностью достаточно легко:
double d = somevalue();
float c = 1.0f / ((float) d * (float) d);
В этом случае 1.0f абсолютно обязателен, если вы используете 1.0, вместо этого вы получите double точность.
Вы пытались включить "небрежную" математику на вашем компиляторе?На GCC вы можете использовать -ffast-math, есть аналогичные опции для других компиляторов.Небрежная математика может быть более чем достаточно для вашего приложения.( Редактировать: Я не увидел никакой разницы в итоговой сборке.)
Если вы используете GCC, рассматривали ли вы вопрос об использовании -mrecip?Существует функция «обратной оценки», которая имеет точность около 12 бит, но она намного быстрее.Вы можете использовать метод Ньютона-Рафсона, чтобы повысить точность результата.Опция -mrecip заставит компилятор автоматически сгенерировать для вас обратную оценку и шаги Ньютона-Рафсона, хотя вы всегда можете написать сборку самостоятельно, если хотите более точно настроить компромисс между точностью и производительностью.(Ньютон-Рафсон быстро сходится очень .) ( Редактировать: Мне не удалось получить GCC для генерации RCPSS. См. Ниже.)
Я нашел сообщение в блоге ( source ), в котором обсуждается точная проблема, с которой вы сталкиваетесь, и автор приходит к выводу, что такие методы, как метод Carmack, не могут конкурировать с инструкцией RCPSS (на которой установлен флаг -mrecip)GCC использует).
Причина , почему деление может быть таким медленным, заключается в том, что процессоры, как правило, имеют только один блок деления, и это часто не конвейерно.Таким образом, вы можете иметь несколько умножений в конвейере, которые выполняются одновременно, но деление невозможно, пока не закончится предыдущее деление.
Трюки, которые не работают
Метод Кармака: он устарел на современных процессорах, которые имеют коды обратной оценки.Для взаимных ссылок лучшая версия, которую я видел, дает только один бит точности - ничто по сравнению с 12 битами RCPSS.Я думаю, что это совпадение, что уловка работает так хорошо для взаимных квадратных корней;совпадение, которое вряд ли повторится.
Перемаркировка переменных.Что касается компилятора, то разница между 1.0/(x*x) и double x2 = x*x; 1.0/x2 очень мала.Я был бы удивлен, если бы вы нашли компилятор, который генерирует различный код для двух версий с оптимизацией, включенной даже до самого низкого уровня.
Использование pow.Библиотека pow - это монстр.При отключенном GCC -ffast-math вызов библиотеки довольно дорогой.При включенном -ffast-math в GCC вы получите точно такой же код сборки для pow(x, -2), что и для 1.0/(x*x), поэтому никаких преимуществ.
Обновление
Вот пример приближения Ньютона-Рафсона для обратного квадрата значения с плавающей запятой двойной точности.
static double invsq(double x)
{
double y;
int i;
__asm__ (
"cvtpd2ps %1, %0\n\t"
"rcpss %0, %0\n\t"
"cvtps2pd %0, %0"
: "=x"(y)
: "x"(x));
for (i = 0; i < RECIP_ITER; ++i)
y *= 2 - x * y;
return y * y;
}
К сожалению, с помощью тестов RECIP_ITER=1 на моем компьютере он немного медленнее~ 5%), чем простая версия 1.0/(x*x).Это быстрее (в 2 раза быстрее) с нулевыми итерациями, но тогда вы получите только 12 бит точности.Я не знаю, достаточно ли для вас 12 битов.
Я думаю, что одной из проблем здесь является то, что это слишком мало для микрооптимизации;в этом масштабе авторы компиляторов почти на равных с хакерами сборок.Возможно, если бы мы имели более широкую картину, мы могли бы увидеть способ сделать это быстрее.
Например, вы сказали, что -ffast-math вызвало нежелательную потерю точности;это может указывать на проблему числовой стабильности в алгоритме, который вы используете.При правильном выборе алгоритма многие проблемы можно решить с помощью float вместо double.(Конечно, вам может потребоваться более 24 бит. Я не знаю.)
Я подозреваю, что метод RCPSS срабатывает, если вы хотите вычислить несколько из них параллельно.