在 C++ 中优化模算术计算答案

【问题标题】：optimizing modular arithmetic computations in C++在 C++ 中优化模算术计算
【发布时间】：2010-12-13 21:57:30
【问题描述】：

我正在开发一些模板化的线性代数代码矩阵系数类型。一种可能的类型是要做的类模运算，天真地实现如下：

template<typename val_t> // `val_t` is an integer type
class Modular 
{
  val_t val_;
  static val_t modulus_;
public:
  Modular(const val_t& value) : val_(value) { };
  static void global_set_modulus(const val_t& modulus) { modulus_ = modulus; };

  Modular<val_t>& operator=(const Modular<val_t>& other) { val_ = other.val_; return *this; }

  Modular<val_t>& operator+=(const Modular<val_t>& other) { val_ += other.val_; val_ %= modulus_; return *this; }
  Modular<val_t>& operator-=(const Modular<val_t>& other) { val_ -= other.val_; val_ %= modulus_; return *this; }
  Modular<val_t>& operator*=(const Modular<val_t>& other) { val_ *= other.val_; val_ %= modulus_; return *this; }
  Modular<val_t>& operator/=(const Modular<val_t>& other) { val_ *= other.inverse().val_; val_ %= modulus_; return *this; }

  friend Modular<val_t> operator+(const Modular<val_t>& a, const Modular<val_t>& b) { return Modular<val_t>((a.val_ + b.val_) % Modular<val_t>::modulus_); };
  friend Modular<val_t> operator-(const Modular<val_t>& a, const Modular<val_t>& b) { return Modular<val_t>((a.val_ - b.val_) % Modular<val_t>::modulus_); };
  // ...etc.
};

但是，当程序以Modular<int> 系数运行时，它是几个比使用int 系数运行时慢几倍。

我应该在“模块化”类中更改哪些内容为了获得最佳性能？

例如，是否可以优化 a*b + c*d 之类的表达式到(a.val_*b.val_ + c.val_*d.val_) % modulus，而不是obvious：

(((a.val_*b.val_) % modulus) + ((c.val_*d.val_ % modulus) % modulus) % modulus)

【问题讨论】：

它在什么编程环境中运行？你打开编译器优化了吗？
您可能希望在乘法之前执行归约，以防止溢出。
如果模数是 2 的幂，您可以使用“x & (modulus-1)”代替“x % 模数”。请注意，负 x 的结果不同，尽管 (-10 % 8 是 -2，但 -10 & (8-1) 是 6)
@Tom Sirgedas：如果在编译时知道模数，任何理智的编译器都会这样做。它被称为“强度降低”，因为它用较弱的运算符替换了一个较强的运算符。
@Jørgen Fogh：你是对的。由于负数的不同行为，我持怀疑态度，但诀窍是 AND 80000007h 然后 OR 与 FFFFFFF8h 用于负数（但这需要一个分支）。

标签： c++ performance math optimization modulo

【解决方案1】：

是的。有可能的。您要查找的是“表达式模板”并从那里开始。从那时起，您将不得不构建一些元程序逻辑来优化/简化表达式。远非一项微不足道的任务，但这不是你所要求的。

NVM - 这很简单：

int count = 0;
int modulus() { count++; return 10; }

template < typename T >
struct modular
{
  modular(T v) : value_(v) {}

  T value() const { return value_; }
  void value(T v) { value_ = v; }

  typedef modular<T> modular_type;
  typedef T raw_type;
private:
  T value_;
};

template < typename LH, typename RH >
struct multiply
{
  multiply(LH l, RH r) : lh(l), rh(r) {}

  typedef typename LH::modular_type modular_type;
  typedef typename LH::raw_type raw_type;

  raw_type value() const { return lh.value() * rh.value(); }

  operator modular_type () const { return modular_type(value() % modulus()); }

private:
  LH lh; RH rh;
};

template < typename LH, typename RH >
struct add
{
  add(LH l, RH r) : lh(l), rh(r) {}

  typedef typename LH::modular_type modular_type;
  typedef typename LH::raw_type raw_type;

  raw_type value() const { return lh.value() + rh.value(); }
  operator modular_type () const { return modular_type(value() % modulus()); }

private:
  LH lh; RH rh;
};

template < typename LH, typename RH >
add<LH,RH> operator+(LH const& lh, RH const& rh)
{
  return add<LH,RH>(lh,rh);
}

template < typename LH, typename RH >
multiply<LH,RH> operator*(LH const& lh, RH const& rh)
{
  return multiply<LH,RH>(lh,rh);
}

#include <iostream>

int main()
{
  modular<int> a = 5;
  modular<int> b = 7;
  modular<int> c = 3;
  modular<int> d = 8;

  std::cout << (5*7+3*8) % 10 << std::endl;

  modular<int> result = a * b + c * d;
  std::cout << result.value() << std::endl;

  std::cout << count << std::endl;

  std::cin.get();
}

如果你很聪明，你会在模块化的构造函数中使用 %，所以它总是模块化的；您还需要检查以确保 LH 和 RH 与 SFINAE 废话兼容，以防止操作员在任何时候杀死它。您还可以将模数作为模板参数并提供元函数来访问它。无论如何......你去吧。

编辑：顺便说一句，您可以使用相同的技术使您的矩阵计算得更快。不是为一系列操作中的每个操作创建一个新矩阵，而是制作这些东西，然后在分配结果时最终逐个元素地进行数学运算。互联网上有关于它的论文和所有东西，将它与 FORTRAN 等进行比较。是元编程的第一个用途之一，例如 C++ 中的模板使用。同样在书中http://www.amazon.com/Scientific-Engineering-Introduction-Advanced-Techniques/dp/0201533936

【讨论】：

或者只是将模数移到转换运算符中。
+1 用于解决根本问题，而不是建议编译器可能已经或可能尚未执行的个别优化。
非常感谢！我不敢指望一个如此详细和有用的答案。 :-)

【解决方案2】：

模数在加法上是可分配的。因此 A % N + B % N == (A + B) % N

关于负操作数或模数，我上次检查 C++ 标准时将结果留给供应商自行决定。所以以上可能不适用于否定。

【讨论】：

对于整数，是的。对于 n 位整数......不是那么多。
如果 A+B 溢出 A % N + B % N 可能不等于 (A + B) % N，因为 A % N + B % N 可能不会溢出。
另外，等式的右边应该是 (A%N+B%N)%N。
我的意思不是严格相等，我的意思是全等 mod N。所以如果没有溢出，我们可以说 (A%N+B%N)%N == (A+B) %N

【解决方案3】：

如果不知道该库的用途，我无法确定，但在我看来，这可能太低级了。

由于您担心性能，我假设您的矩阵非常大。这意味着通过使用更快的算法，您可能会看到比尝试优化这样的东西更大的速度提升。不管你做什么，int 系数可能会更快。

即使您节省了一些 mod 操作，加速也只会是一个常数因子，可能不到 10 倍。对于大多数矩阵运算，优化缓存可能会给您带来更多好处。

我的建议是分析一下哪些操作太慢，然后用谷歌搜索该操作并查看存在哪些算法（例如，Strassen algorithm 用于乘法）。您应该知道矩阵有多大以及它们是稀疏的还是密集的。

无论如何，如果您必须询问这些内容，您最好还是使用现有的库。

【讨论】：