batch_inverse.hpp¶
只做一次除法和线性次乘法，批量求一组非零域元素的逆元。
\[
\displaystyle y_i = x_i^{-1}
\]
Complexity: Time: \(O(n)\) field operations, one division for nonempty input. Space: \(O(n)\) output and \(O(1)\) additional workspace.
Implementation¶
当前头文件，省略 include guard；依赖见 #include。
/// @complexity Time: O(n) field operations, one division for nonempty input.
/// Space: O(n) output and O(1) additional workspace.
#include <array>
#include <cassert>
#include <vector>
namespace noya {
/// @brief Invert nonzero field elements with four independent product chains.
/// Each chain contains indices with the same remainder modulo four.
template <class T> std::vector<T> batch_inverse(const std::vector<T> &vs) {
if (vs.empty()) return {};
int n = int(vs.size());
std::vector<T> res(n);
std::array<T, 4> prod;
prod.fill(T(1));
#pragma GCC unroll 4
for (int i = 0; i < n; ++i) {
assert(vs[i] != T{});
res[i] = prod[i & 3];
prod[i & 3] *= vs[i];
}
// Invert the four chain products together, using just one division.
auto last = prod;
for (int i = 1; i < 4; ++i) prod[i] *= prod[i - 1];
T inv = T(1) / prod[3];
for (int i = 3; i > 0; --i) {
prod[i] = prod[i - 1] * inv;
inv *= last[i];
}
prod[0] = inv;
#pragma GCC unroll 4
for (int i = n - 1; i >= 0; --i) {
res[i] *= prod[i & 3];
prod[i & 3] *= vs[i];
}
return res;
}
} // namespace noya
// NOYA_BATCH_INVERSE_HPP
#ifndef NOYA_BATCH_INVERSE_HPP
#define NOYA_BATCH_INVERSE_HPP 1
/// @complexity Time: O(n) field operations, one division for nonempty input.
/// Space: O(n) output and O(1) additional workspace.
#include <array>
#include <cassert>
#include <vector>
namespace noya {
/// @brief Invert nonzero field elements with four independent product chains.
/// Each chain contains indices with the same remainder modulo four.
template <class T> std::vector<T> batch_inverse(const std::vector<T> &vs) {
if (vs.empty()) return {};
int n = int(vs.size());
std::vector<T> res(n);
std::array<T, 4> prod;
prod.fill(T(1));
#pragma GCC unroll 4
for (int i = 0; i < n; ++i) {
assert(vs[i] != T{});
res[i] = prod[i & 3];
prod[i & 3] *= vs[i];
}
// Invert the four chain products together, using just one division.
auto last = prod;
for (int i = 1; i < 4; ++i) prod[i] *= prod[i - 1];
T inv = T(1) / prod[3];
for (int i = 3; i > 0; --i) {
prod[i] = prod[i - 1] * inv;
inv *= last[i];
}
prod[0] = inv;
#pragma GCC unroll 4
for (int i = n - 1; i >= 0; --i) {
res[i] *= prod[i & 3];
prod[i & 3] *= vs[i];
}
return res;
}
} // namespace noya
#endif // NOYA_BATCH_INVERSE_HPP
/// @complexity Time: O(n) field operations, one division for nonempty input.
/// Space: O(n) output and O(1) additional workspace.
#include <array>
#include <cassert>
#include <vector>
namespace noya {
/// @brief Invert nonzero field elements with four independent product chains.
/// Each chain contains indices with the same remainder modulo four.
template <class T> std::vector<T> batch_inverse(const std::vector<T> &vs) {
if (vs.empty()) return {};
int n = int(vs.size());
std::vector<T> res(n);
std::array<T, 4> prod;
prod.fill(T(1));
#pragma GCC unroll 4
for (int i = 0; i < n; ++i) {
assert(vs[i] != T{});
res[i] = prod[i & 3];
prod[i & 3] *= vs[i];
}
// Invert the four chain products together, using just one division.
auto last = prod;
for (int i = 1; i < 4; ++i) prod[i] *= prod[i - 1];
T inv = T(1) / prod[3];
for (int i = 3; i > 0; --i) {
prod[i] = prod[i - 1] * inv;
inv *= last[i];
}
prod[0] = inv;
#pragma GCC unroll 4
for (int i = n - 1; i >= 0; --i) {
res[i] *= prod[i & 3];
prod[i & 3] *= vs[i];
}
return res;
}
} // namespace noya
// NOYA_BATCH_INVERSE_HPP