实现 Eigen 库伪逆函数的 MEX 文件崩溃

Question

我正在尝试在 Matlab MEX 文件中实现 Eigen 库伪逆函数。它编译成功但是当我运行它时崩溃。

我正在尝试关注 the FAQ on how to implement a pseudo-inverse function using the Eigen library。

常见问题解答建议将其作为方法添加到 JacobiSVD class，但由于您不能在 C++ 中这样做，我将其添加到子 class。它编译成功，但随后崩溃且没有错误消息。如果我用 .pinv 调用注释掉该行，它会成功输出 "hi" 而不会崩溃，这就是问题所在。对于运行，我只是编译它（如 test.cpp），然后在命令行输入 test。我在 MacOS 10.14.5 和 Eigen 3.3.7 下使用 Matlab R2019a。在我的完整代码中，我还收到了很多关于 pinv 代码的奇怪错误消息，但在我进行故障排除之前，我需要这个简单的测试用例才能工作。这完全超出了我对 C++ 理解的极限。任何帮助表示赞赏。

#include "mex.h"
#include <Eigen/Dense>
#include <Eigen/SVD>
#include <cstdlib>
#include <cmath>
#include <iostream>

using namespace Eigen;
using namespace std;

//
class JacobiSVDext : public JacobiSVD<MatrixXf> {
    typedef SVDBase<JacobiSVD<MatrixXf>> Base;
    public:
    using JacobiSVD::JacobiSVD; //inherit constructors //
    MatrixXf pinv() //http://eigen.tuxfamily.org/index.php?title=FAQ
    {
        eigen_assert(m_isInitialized && "SVD is not initialized.");
        double  pinvtoler=1.e-6; // choose your tolerance wisely!
        JacobiSVDext::SingularValuesType singularValues_inv=m_singularValues;
        for ( long i=0; i<m_workMatrix.cols(); ++i) {
            if ( m_singularValues(i) > pinvtoler )
                singularValues_inv(i)=1.0/m_singularValues(i);
            else singularValues_inv(i)=0;
        }
        return m_matrixV*singularValues_inv.asDiagonal()*m_matrixU.transpose();
    };
};

/* The gateway function */
void mexFunction( int nlhs, mxArray *plhs[],
             int nrhs, const mxArray *prhs[])
{
    MatrixXf X = MatrixXf::Random(5, 5);
    JacobiSVDext svd(X);
    MatrixXf Y=svd.pinv();
    cout << Y << endl;
    cout << "hi" << endl;
}

预期的结果是输出随机矩阵的伪逆，也是"hi"。相反，它会在没有错误消息的情况下崩溃。

Answer 1

构造Eigen::JacobiSVD对象时，您请求计算矩阵U和V失败。默认情况下，不计算这些。显然，如果未计算这些矩阵，则访问这些矩阵将导致分段违规。

参见 the documentation to the constructor。第二个输入参数必须指定 ComputeFullU | ComputeFullV 或 ComputeThinU | ComputeThinV。计算 pseudo-inverse 时最好使用薄矩阵，因为不需要其他矩阵。

我不会从 JacobiSVD class 派生只是为了添加一个方法。相反，我会简单地编写一个免费函数。这既更容易，又允许您仅使用 Eigen API.

的文档部分

我写了下面的 MEX-file，它按预期工作（使用我已有的代码进行此计算）。它的作用相同，但方式略有不同，可避免编写显式循环。不确定这种写法是不是很清楚，但是可以用。

// Compile with:
//    mex -v test.cpp -I/usr/local/include/eigen3

#include "mex.h"
#include <Eigen/Dense>
#include <Eigen/SVD>
#include <cstdlib>
#include <cmath>
#include <iostream>

Eigen::MatrixXf PseudoInverse(Eigen::MatrixXf matrix) {
   Eigen::JacobiSVD< Eigen::MatrixXf > svd( matrix, Eigen::ComputeThinU | Eigen::ComputeThinV );
   float tolerance = 1.0e-6f * float(std::max(matrix.rows(), matrix.cols())) * svd.singularValues().array().abs()(0);
   return svd.matrixV()
         * (svd.singularValues().array().abs() > tolerance).select(svd.singularValues().array().inverse(), 0).matrix().asDiagonal()
         * svd.matrixU().adjoint();
}

void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    Eigen::MatrixXf X = Eigen::MatrixXf::Random(5, 5);
    Eigen::MatrixXf Y = PseudoInverse(X);
    std::cout << Y << '\n';
    std::cout << "hi\n";
}

实现 Eigen 库伪逆函数的 MEX 文件崩溃

MEX-file implementing Eigen library pseudo-inverse function crashes

c++

matlab

mex

eigen