将稀疏数组从 matlab 传递到 Eigen (C++），然后再传递回 matlab

passing sparse arrays from matlab to Eigen (C++) and back to matlab?

本文关键字：matlab C++ 然后 数组 Eigen      更新时间：2023-10-16

下面是一个 mex 代码，它使用 Eigen 将 matlab 中的密集数组 g 和 G 相乘。当 g 稀疏时如何执行此操作？

#include <iostream>
#include <Eigen/Dense>
#include "mex.h"
using Eigen::MatrixXd;
using namespace Eigen;
/*gateway function*/
void mexFunction( int nlhs, mxArray *plhs[],
        int nrhs, const mxArray *prhs[]) {
    int nRows=(int)mxGetM(prhs[0]);
    int nCols=nRows;
    double* g=mxGetPr(prhs[0]);
    double* Gr=mxGetPr(prhs[1]);
    Map<MatrixXd> gmap (g, nRows, nCols );
    Map<MatrixXd> Grmap (Gr, nRows, nCols );
    plhs[0] = mxCreateDoubleMatrix(nRows, nCols, mxREAL);
    Map<MatrixXd> resultmap (mxGetPr(plhs[0]), nRows, nCols); 
    resultmap = gmap*Grmap; 
}

您可以使用这些函数在 MATLAB 和 Eigen* 之间传递稀疏(压缩(双矩阵：

#include "mex.h"
#include <Eigen/Sparse>
#include <type_traits>
#include <limits>
using namespace Eigen;
typedef SparseMatrix<double,ColMajor,std::make_signed<mwIndex>::type> MatlabSparse;

Map<MatlabSparse >
matlab_to_eigen_sparse(const mxArray * mat)
{
    mxAssert(mxGetClassID(mat) == mxDOUBLE_CLASS,
    "Type of the input matrix isn't double");
    mwSize     m = mxGetM (mat);
    mwSize     n = mxGetN (mat);
    mwSize    nz = mxGetNzmax (mat);
    /*Theoretically fails in very very large matrices*/
    mxAssert(nz <= std::numeric_limits< std::make_signed<mwIndex>::type>::max(),
    "Unsupported Data size."
    );
    double  * pr = mxGetPr (mat);
    MatlabSparse::StorageIndex* ir = reinterpret_cast<MatlabSparse::StorageIndex*>(mxGetIr (mat));
    MatlabSparse::StorageIndex* jc = reinterpret_cast<MatlabSparse::StorageIndex*>(mxGetJc (mat));
    Map<MatlabSparse> result (m, n, nz, jc, ir, pr);
    return result;
}
mxArray* 
eigen_to_matlab_sparse(const Ref<const MatlabSparse,StandardCompressedFormat>& mat)
{
    mxArray * result = mxCreateSparse (mat.rows(), mat.cols(), mat.nonZeros(), mxREAL);
    const MatlabSparse::StorageIndex* ir = mat.innerIndexPtr();
    const MatlabSparse::StorageIndex* jc = mat.outerIndexPtr();
    const double* pr = mat.valuePtr();
    mwIndex * ir2 = mxGetIr (result);
    mwIndex * jc2 = mxGetJc (result);
    double  * pr2 = mxGetPr (result);
    for (mwIndex i = 0; i < mat.nonZeros(); i++) {
        pr2[i] = pr[i];
        ir2[i] = ir[i];
    }
    for (mwIndex i = 0; i < mat.cols() + 1; i++) {
        jc2[i] = jc[i];
    }
    return result;
}

  • 读写从这里采用的MATLAB/Octave稀疏矩阵。

  • 感谢@chtz的建议。

相关文章：