Есть ли способ улучшить взаимное (деление 1 по X) относительно скорости, если точность не имеет решающего значения?
Итак, мне нужно вычислить 1/X. Есть ли способ обхода, поэтому я теряю точность, но делаю это быстрее?
Быстрое 1/X-деление (обратное)
Ответ 1
Я считаю, что то, что он искал, - это более эффективный способ приближения 1.0/x вместо некоторого технического определения приближения, в котором говорится, что вы можете использовать 1 как очень неточный ответ. Я также считаю, что это удовлетворяет это.
#ifdef __cplusplus
#include <cstdint>
#else
#include <stdint.h>
#endif
__inline__ double __attribute__((const)) reciprocal( double x ) {
union {
double dbl;
#ifdef __cplusplus
std::uint_least64_t ull;
#else
uint_least64_t ull;
#endif
} u;
u.dbl = x;
u.ull = ( 0xbfcdd6a18f6a6f52ULL - u.ull ) >> 1;
// pow( x, -0.5 )
u.dbl *= u.dbl; // pow( pow(x,-0.5), 2 ) = pow( x, -1 ) = 1.0 / x
return u.dbl;
}
__inline__ float __attribute__((const)) reciprocal( float x ) {
union {
float single;
#ifdef __cplusplus
std::uint_least32_t uint;
#else
uint_least32_t uint;
#endif
} u;
u.single = x;
u.uint = ( 0xbe6eb3beU - u.uint ) >> 1;
// pow( x, -0.5 )
u.single *= u.single; // pow( pow(x,-0.5), 2 ) = pow( x, -1 ) = 1.0 / x
return u.single;
}
Хм....... Я был бы уверен, если бы производители ЦП знали, что вы могли бы приблизить обратное значение только с помощью одного умножения, вычитания и сдвига битов, когда они проектировали ЦП.... хмм........,
Что касается тестирования производительности, инструкции по аппаратному обеспечению x 2 в сочетании с инструкциями по аппаратному вычитанию выполняются так же быстро, как и инструкции по оборудованию 1.0/x на современных компьютерах (мои тесты проводились на Intel i7, но я бы предположил, что аналогичные результаты для других процессоров), Однако если бы этот алгоритм был внедрен в аппаратное обеспечение как новая инструкция по сборке, то увеличение скорости, вероятно, было бы достаточно хорошим, чтобы эта инструкция была довольно практичной.
Для получения дополнительной информации об этом методе, эта реализация основана на замечательном "быстром" алгоритме обратного квадратного корня.
Как Pharap обратил мое внимание, чтение неактивного свойства из объединения является неопределенным поведением, поэтому есть два возможных решения, которые я разработал из его полезных комментариев, чтобы избежать неопределенного поведения. Первое решение больше похоже на неприятный трюк для обхода языковой семантики, которая практически не лучше оригинального решения.
#ifdef __cplusplus
#include <cstdint>
#else
#include <stdint.h>
#endif
__inline__ double __attribute__((const)) reciprocal( double x ) {
union {
double dbl[2];
#ifdef __cplusplus
std::uint_least64_t ull[2];
#else
uint_least64_t ull[2];
#endif
} u;
u.dbl[0] = x; // dbl is now the active property, so only dbl can be read now
u.ull[1] = 0;//trick to set ull to the active property so that ull can be read
u.ull][0] = ( 0xbfcdd6a18f6a6f52ULL - u.ull[0] ) >> 1;
u.dbl[1] = 0; // now set dbl to the active property so that it can be read
u.dbl[0] *= u.dbl[0];
return u.dbl[0];
}
__inline__ float __attribute__((const)) reciprocal( float x ) {
union {
float flt[2];
#ifdef __cplusplus
std::uint_least32_t ull[2];
#else
uint_least32_t ull[2];
#endif
} u;
u.flt[0] = x; // now flt is active
u.uint[1] = 0; // set uint to be active for reading and writing
u.uint[0] = ( 0xbe6eb3beU - u.uint[0] ) >> 1;
u.flt[1] = 0; // set flt to be active for reading and writing
u.flt[0] *= u.flt[0];
return u.flt[0];
}
Второе возможное решение гораздо более приемлемо, потому что оно полностью избавляет от профсоюзов. Однако это решение будет намного медленнее, если оно не будет должным образом оптимизировано компилятором. Но, с другой стороны, решение ниже будет полностью независимым от порядка следования байтов:
- что байты имеют ширину 8 бит
- эти байты являются наименьшей атомарной единицей на целевой машине.
- двойники имеют ширину 8 байт, а числа с плавающей запятой - 4 байта.
#ifdef __cplusplus
#include <cstdint>
#include <cstring>
#define stdIntWithEightBits std::uint8_t
#define stdIntSizeOfFloat std::uint32_t
#define stdIntSizeOfDouble std::uint64_t
#else
#include <stdint.h>
#include <string.h>
#define stdIntWithEightBits uint8_t
#define stdIntSizeOfFloat uint32_t
#define stdIntSizeOfDouble uint64_t
#endif
__inline__ double __attribute__((const)) reciprocal( double x ) {
double byteIndexFloat = 1.1212798184631136e-308;//00 08 10 18 20 28 30 38 bits
stdIntWithEightBits* byteIndexs = reinterpret_cast<stdIntWithEightBits*>(&byteIndexFloat);
stdIntWithEightBits* inputBytes = reinterpret_cast<stdIntWithEightBits*>(&x);
stdIntSizeOfDouble inputAsUll = (
(inputBytes[0] << byteIndexs[0]) |
(inputBytes[1] << byteIndexs[1]) |
(inputBytes[2] << byteIndexs[2]) |
(inputBytes[3] << byteIndexs[3]) |
(inputBytes[4] << byteIndexs[4]) |
(inputBytes[5] << byteIndexs[5]) |
(inputBytes[6] << byteIndexs[6]) |
(inputBytes[7] << byteIndexs[7])
);
inputAsUll = ( 0xbfcdd6a18f6a6f52ULL - inputAsUll ) >> 1;
double outputDouble;
const stdIntWithEightBits outputBytes[] = {
inputAsUll >> byteIndexs[0],
inputAsUll >> byteIndexs[1],
inputAsUll >> byteIndexs[2],
inputAsUll >> byteIndexs[3],
inputAsUll >> byteIndexs[4],
inputAsUll >> byteIndexs[5],
inputAsUll >> byteIndexs[6],
inputAsUll >> byteIndexs[7]
};
memcpy(&outputDouble, &outputBytes, 8);
return outputDouble * outputDouble;
}
__inline__ float __attribute__((const)) reciprocal( float x ) {
float byteIndexFloat = 7.40457e-40; // 0x00 08 10 18 bits
stdIntWithEightBits* byteIndexs = reinterpret_cast<stdIntWithEightBits*>(&byteIndexFloat);
stdIntWithEightBits* inputBytes = reinterpret_cast<stdIntWithEightBits*>(&x);
stdIntSizeOfFloat inputAsInt = (
(inputBytes[0] << byteIndexs[0]) |
(inputBytes[1] << byteIndexs[1]) |
(inputBytes[2] << byteIndexs[2]) |
(inputBytes[3] << byteIndexs[3])
);
inputAsInt = ( 0xbe6eb3beU - inputAsInt ) >> 1;
float outputFloat;
const stdIntWithEightBits outputBytes[] = {
inputAsInt >> byteIndexs[0],
inputAsInt >> byteIndexs[1],
inputAsInt >> byteIndexs[2],
inputAsInt >> byteIndexs[3]
};
memcpy(&outputFloat, &outputBytes, 4);
return outputFloat * outputFloat;
}
Отказ от ответственности: Наконец, обратите внимание, что я больше новичок в C++. В связи с этим я приветствую с распростертыми объятиями любые передовые практики, правильное форматирование или внесение ясности в смысл до конца как улучшения качества этого ответа для всех, кто его читает, так и расширения моих знаний о C++ для всех моих долгие годы (если, конечно, я не попаду в автомобильную аварию завтра и умру).
Ответ 2
Во-первых, убедитесь, что это не случай преждевременной оптимизации. Вы знаете, что это ваше узкое место?
Как говорит Мистик, 1/x можно рассчитать очень быстро. Убедитесь, что вы не используете тип данных double для 1 или делителя. Поплавки намного быстрее.
Тем не менее, бенчмарк, бенчмарк, бенчмарк. Не тратьте свое время на часовые часы на численную теорию, чтобы обнаружить, что источником низкой производительности является доступ к IO.
Ответ 3
Прежде всего, если вы включите оптимизацию компилятора, компилятор, скорее всего, оптимизирует вычисление (например, вытащить его из цикла). Чтобы увидеть эту оптимизацию, вам нужно создать и запустить в режиме Release.
Подразделение может быть тяжелее, чем умножение (но комментатор указал, что ответчики так же быстро, как умножение на современные процессоры, и в этом случае это неверно для вашего случая), поэтому, если у вас есть 1/X где-то внутри цикла (и более одного раза), вы можете помочь, кэшируя результат внутри цикла (float Y = 1.0f/X;), а затем используя Y. (Оптимизация компилятора может сделать это в любом случае.)
Кроме того, некоторые формулы могут быть переработаны для удаления деления или других неэффективных вычислений. Для этого вы можете опубликовать более крупное вычисление. Даже там сама программа или алгоритм иногда могут быть реструктурированы, чтобы предотвратить частое попадание вредных циклов.
Какую точность можно пожертвовать? Если по возможности вам нужен только порядок, вы можете легко получить это с помощью оператора модуля или побитовых операций.
Однако, в общем, нет возможности ускорить деление. Если бы это было, компиляторы уже делали бы это.
Ответ 4
Это должно сделать это с помощью нескольких предварительно развернутых итераций нонтона, оцененных как многочлен Horner, который использует операции накопления с накоплением многократного накопления, выполняемые большинством современных CPU в одном цикле Clk (каждый раз):
float inv_fast(float x) {
union { float f; int i; } v;
float w, sx;
int m;
sx = (x < 0) ? -1:1;
x = sx * x;
v.i = (int)(0x7EF127EA - *(uint32_t *)&x);
w = x * v.f;
// Efficient Iterative Approximation Improvement in horner polynomial form.
v.f = v.f * (2 - w); // Single iteration, Err = -3.36e-3 * 2^(-flr(log2(x)))
// v.f = v.f * ( 4 + w * (-6 + w * (4 - w))); // Second iteration, Err = -1.13e-5 * 2^(-flr(log2(x)))
// v.f = v.f * (8 + w * (-28 + w * (56 + w * (-70 + w *(56 + w * (-28 + w * (8 - w))))))); // Third Iteration, Err = +-6.8e-8 * 2^(-flr(log2(x)))
return v.f * sx;
}
Тонкая печать: ближе к 0, приближение не так хорошо, так что либо вы, программист, должны проверить производительность, либо ограничить ввод данных до минимума, прежде чем приступать к аппаратной части. то есть быть ответственным!
Ответ 5
Я тестировал эти методы на Arduino NANO на скорость и "точность".
Основной расчет должен был установить переменные, Y = 15 000 000 и Z = 65 535
(в моем реальном случае Y - это константа, а Z может варьироваться от 65353 до 3000, поэтому полезный тест)
Время вычислений на Arduino было установлено путем установки минимального значения пин-кода, затем высокого значения, полученного при вычислении, а затем снова низкого уровня и сравнения времени с логическим анализатором. ЗА 100 ЦИКЛОВ. С переменными в виде целых чисел без знака: -
Y * Z takes 0.231 msec
Y / Z takes 3.867 msec.
With variables as floats:-
Y * Z takes 1.066 msec
Y / Z takes 4.113 msec.
Basic Bench Mark and ( 15,000,000/65535 = 228.885 via calculator.)
Используя обратный алгоритм {Джек Гиффинс}:
Y * reciprocal(Z) takes 1.937msec which is a good improvement, but accuracy less so 213.68.
Используя {nimig18s} float inv_fast:
Y* inv_fast(Z) takes 5.501 msec accuracy 228.116 with single iteration
Y* inv_fast(Z) takes 7.895 msec accuracy 228.883 with second iteration
Использование Википедии Q_rsqrt (на которую указывает {Джек Гиффин})
Y * Q*rsqrt(Z) takes 6.104 msec accuracy 228.116 with single iteration
All entertaining but ultimately disappointing!
Ответ 6
Самый быстрый способ, которым я знаю, - использовать операции SIMD. http://msdn.microsoft.com/en-us/library/796k1tty(v=vs.90).aspx