如何使用Z_2中的系数求解稀疏线性系统

Question

我想使用eigen求解具有Z_2系数的大稀疏线性方程系统。首先，我们尝试了不起作用的布尔类型，因为布尔值1 1 = 1，但我希望1 1 = 0。因此，解决方案可能是一种新的自定义标量类型。但是它如何工作？也欢迎其他软件包或软件的建议。

Answer 1

作为基本类型的运营商不能超载。您需要定制的标量类型。这是告诉您如何执行此操作的基本文件。

http://eigen.tuxfamily.org/dox/topiccustomizingeigen.html#customscalartype

基本上您需要做3件事。

1. 确保type t
支持公共操作员（ ， - ，*，/等）
2. 添加struct eigen :: Numtraits的专业化
3. 定义对您类型有意义的数学功能。这包括标准的标准，例如SQRT，POW，SIN，TAN，CONJ，REAL，IMAG等，以及特定于特异性的ABS2。（请参阅文件eigen/src/core/Mathfunctions.h）

实际上，您只能定义方程求解所需的那些运算符和数学功能。

上面的链接为adtl::adouble类型提供了一个简单的示例。该代码非常短，因为大多数操作已经定义得很好。在EIGEN源DIR unsupported/中，还有另一个第三方类型的mpfr::mpreal支持。您可以从此链接开始。

https://eigen.tuxfamily.org/dox/unsupported/group__mprealsupport__module.html

在特征源DIR中，这些文件与mpfr::mpreal支持有关。它们可能很有用。

./无支持/eigen/mprealsupport./unsupported/test/mpreal/mpreal.h./unsupported/test/mpreal_support.cpp

---

这是一个最低的工作示例，在z_2中使用矩阵乘法支持：

#include <iostream>
#include <Eigen/Eigen>
namespace Z2 {
struct Scalar {
  Scalar() :
      data(0) {
  }
  Scalar(bool x) :
      data(x) {
  }
  bool data;
  inline Scalar operator+=(const Scalar& a) {
    return data ^= a.data;
  }
};
inline Scalar abs(const Scalar& a) {
  return a;
}
inline Scalar operator+(const Scalar& a, const Scalar& b) {
  return a.data ^ b.data;
}
inline Scalar operator*(const Scalar& a, const Scalar& b) {
  return a.data & b.data;
}
template<class E, class CT>
std::basic_ostream<E, CT> &operator <<(std::basic_ostream<E, CT> &os,
                                       const Z2::Scalar &m) {
  return os << m.data;
}
}
namespace Eigen {
// follow all other traits of bool
template<> struct NumTraits<Z2::Scalar> : NumTraits<bool> {
  typedef Z2::Scalar Real;
  typedef typename internal::conditional<sizeof(Z2::Scalar) <= 2, float, double>::type NonInteger;
  typedef Z2::Scalar Nested;
};
}
int main(void) {
  using namespace Eigen;
  Matrix<Z2::Scalar, Dynamic, Dynamic> x(2, 2), y(2, 2), i(2, 2), j(2, 2);
  x.row(0) << 1, 1;
  y.col(0) << 1, 1;
  i.setIdentity();
  j = i.array() + 1;
  std::cout << "x=n" << x << std::endl;
  std::cout << "y=n" << y << std::endl;
  std::cout << "i=n" << i << std::endl;
  std::cout << "j=n" << j << std::endl;
  std::cout << "x+y=n" << x + y << std::endl;
  std::cout << "x.*y=n" << x.array() * y.array() << std::endl;
  std::cout << "y*j=n" << y * j << std::endl;
  std::cout << "abs(x)=n" << x.array().abs() << std::endl;
  return 0;
}

结果：

x=
1 1
0 0
y=
1 0
1 0
i=
1 0
0 1
j=
0 1
1 0
x+y=
0 1
1 0
x.*y=
1 0
0 0
y*j=
0 1
0 1
abs(x)=
1 1
0 0