4

我有以下代码(我为稍大的代码片段道歉,这是我能够减少问题的最小示例):

#include <Eigen/Dense>
#include <complex>
#include <iostream>
#include <typeinfo>

// Dynamic Matrix over Scalar field
template <typename Scalar> 
using DynMat = Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>;

// Dynamic column vector over Scalar field
template <typename Scalar>
using DynVect = Eigen::Matrix<Scalar, Eigen::Dynamic, 1>;

// Returns the D x D Identity matrix over the field Derived::Scalar
// deduced from the expression Eigen::MatrixBase<Derived>& A
template<typename Derived>
DynMat<typename Derived::Scalar> Id(const Eigen::MatrixBase<Derived>& A, std::size_t D)
{   
    DynMat<typename Derived::Scalar> result =
            DynMat<typename Derived::Scalar>::Identity(D, D);

    return result;
}

int main()
{
    //using ScalarField = std::complex<double>; // same issue even if I use complex numbers
    using ScalarField = double; // we use doubles in this example

    // A double dynamic matrix (i.e. MatrixXd)
    DynMat<ScalarField> Foo; // used to deduce the type in Id<>()

    // A double dynamic column vector (i.e. VectorXd)
    DynVect<ScalarField> v(4);
    v << 1., 0. , 0. ,0.; // plug in some values into it

    // Make sure that Id(Foo, 4) correctly deduces the template parameters
    std::cout << "Id(Foo, 4) is indeed the 4 x 4 identiy matrix over the ScalarField of "
              << "typeid().name(): " << typeid(ScalarField).name() << std::endl;
    std::cout << Id(Foo, 4) << std::endl; // Indeed the 4 x 4 complex Identity matrix

    // Use auto type deduction for GenMatProduct, junk is displayed. Why?!
    std::cout << std::endl << "Use auto type deduction for GenMatProduct,\
                 sometimes junk is displayed. Why?!" << std::endl;
    auto autoresult = Id(Foo, 4) * v; // evaluated result must be identically equal to v
    for(int i=0; i<10; i++)
    {
            std::cout << autoresult.transpose(); // thought 1 0 0 0 is the result, but NO, junk
            std::cout << " has norm: " << autoresult.norm() << std::endl; // junk
    }

    // Use implicit cast to Dynamic Matrix, works fine
    std::cout << std::endl << "Use implicit cast to Dynamic Matrix, works fine" << std::endl;
    DynMat<ScalarField> castresult = Id(Foo, 4) * v; // evaluated result must be identically equal to v
    for(int i=0; i<10; i++)
    {
            std::cout << castresult.transpose(); // 1 0 0 0, works ok
            std::cout << " has norm: " << castresult.norm() << std::endl; // ok
    }
}

主要思想是模板函数Id<>()将 Eigen 表达式A 作为参数,连同 size D,并在表达式的标量字段上生成单位矩阵A。此功能本身可以正常工作。但是,当我在具有推导类型的 Eigen 乘积中使用它时auto,例如在 line auto autoresult = Id(Foo, 4) * v中,我希望将向量乘以v 单位矩阵,因此最终结果应该是一个表达式,在评估时应该完全等于v. 但事实并非如此,请参阅第一个for循环,每当我显示结果并计算其范数时,我大部分时间都是垃圾。另一方面,如果我隐式转换 Eigen 积Id(Foo, 4) * v对于动态矩阵,一切正常,结果得到正确评估。

我在 OS X Yosemite 上使用 Eigen 3.2.2,并在 g++4.9.1 和 Apple LLVM 6.0 版(clang-600.0.54)(基于 LLVM 3.5svn)中得到相同的奇怪行为

问题:

  • 我不明白第一个for循环中发生了什么,为什么在我使用时std::cout,甚至在我使用该norm方法时都没有评估产品?我错过了什么吗?这里没有涉及混叠,我真的对发生的事情感到困惑。我知道 Eigen 使用惰性求值,并在需要时对表达式求值,但这里似乎并非如此。这个问题对我来说非常重要,因为我有很多与 类似的函数,Id<>()auto 推导表达式中使用时可能会失败。

问题经常发生,但并非总是如此。但是,如果您将程序运行 3-4 次,您肯定会看到它。

我用来编译和运行它的命令是:

clang++ (g++) -std=c++11 -isystem ./eigen_3.2.2/ testeigen.cpp -otesteigen; ./testeigen

我在实际运行中得到的典型输出是:

Id(Foo, 4) is indeed the 4 x 4 identiy matrix over the ScalarField of typeid().name(): d
1 0 0 0
0 1 0 0
0 0 1 0
0 0 0 1

Use GenMatProduct, sometimes junk is displayed. Why?!
1 0 0 0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf
3.10504e+231 3.10504e+231 3.95253e-323            0 has norm: inf

Use implicit cast to Dynamic Matrix, works fine
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1
1 0 0 0 has norm: 1

即使我eval()使用

  std::cout << autoresult.eval().transpose(); // thought 1 0 0 0 is the result, but NO, junk
  std::cout << " has norm: " << autoresult.eval().norm() << std::endl; // junk

我得到了同样奇怪的行为。

4

2 回答 2

3

问题是它Id()返回一个临时值,它通过引用存储在表示表达式的对象中Id(Foo, 4) * v。因此,在 auto 语句之后,autoresult存储对死对象的引用。如果您不想要抽象表达式而是实际结果,请不要使用auto或调用eval来强制评估:

auto autoresult = (Id(Foo, 4) * v).eval();

第三种选择是使返回的对象Id()可用于进一步计算:

auto id4 = Id(Foo,4);
auto autoresult = id4 * v;

但在这种情况下,无论您何时使用autoresult,都会重新评估产品,以下将输出不同的结果:

cout << autoresult;
v.setRandom();
cout << autoresult;
于 2014-11-14T07:50:09.703 回答
1

它可能有一个惰性评估类型,只能安全评估一次。您可以使用以下方法捕获它:

auto autoresultmatrix = autoresult.eval()
于 2014-11-02T23:02:44.153 回答