A sparse matrix vector multiplication in C mex programming

Question

JM 2017-11-28

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/369788-a-sparse-matrix-vector-multiplication-in-c-mex-programming

评论： Royi Avital 2021-8-26

Hi everyone,

I'm new in C mex programming and I'm really stuck. I would like to write a code that computes a sparse matrix vector multiplication (sparse_matrix*vector=result_vector). Given how sparse matrices are managed in C mex, I have difficulties getting the non zero elements of each row (to do the multiplication with the vector), because the pointers Ir=mxGetIr(sparse_matrix) and Jc=mxGetJc(sparse_matrix) give the information about the columns only and not the rows.

I would be really grateful if somebody helps me out because I can't move on without solving this problem.

Thanks in advance.

Sincerely,

J

3 个评论
显示 1更早的评论隐藏 1更早的评论

James Tursa 2017-11-28

I'm with Matt J. Why are you doing this? And if for some reason you really do need to do this, what have you done so far? Where specifically are you stuck? The mxGetIr( ) gives the row indexing, btw. Is the vector full or sparse?

JM 2017-11-28

Thanks to both of you for you answers. Matt J I'm not coding in matlab but I'm coding in C and then I compile from matlab using mex engine, and as far as I know, there's no routine already implemented that does sparse matrix vector multiplication in C.

As you know, in order to carry out the multiplication, each row of the sparse matrix is multiplied by the vector which is a full one, and then the sum of the multiplication result gives the result that is stored in the result vector. So how can I, for each row of the sparse matrix, access the non zeros elements in order to multiply them with the vector ? I can't see how to do this with Ir and Jc because they contain the information about the columns only and not the rows.

Please let me know if I'm not clear enough.

Thanks in advance.

J

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

James Tursa 2017-11-29

2
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/369788-a-sparse-matrix-vector-multiplication-in-c-mex-programming#answer_293672

编辑：James Tursa 2017-12-6

在 MATLAB Online 中打开

Here is the bare bones C code for this. Most of the code consists of argument checking etc. The actual computation part of the code is only a few lines and is fairly trivial. It simply accumulates the individual element multiplies into the appropriate resulting element. The j for-loop loops over the columns of M, and the while loop does the individual element multiplies within each column.

/* Bare bones sparse matrix * full vector                               */
/* Produces a full vector result                                        */
/* Assumes all input elements are finite                                */
/* E.g., does not account for NaN or 0*inf etc.                         */
/* Programmer:  James Tursa                                             */
/* Includes ----------------------------------------------------------- */
#include "mex.h"
/* Gateway ------------------------------------------------------------ */
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    mwSize m, n, j, nrow;
    double *Mpr, *Vpr, *Cpr;
    mwIndex *Mir, *Mjc;
/* Argument checks */
    if( nlhs > 1 ) {
        mexErrMsgTxt("Too many outputs");
    }
    if( nrhs != 2 ) {
        mexErrMsgTxt("Need exactly two inputs");
    }
    if( !mxIsDouble(prhs[0]) || !mxIsSparse(prhs[0]) || mxIsComplex(prhs[0]) ) {
        mexErrMsgTxt("1st argument must be real double sparse matrix");
    }
    if( !mxIsDouble(prhs[1]) || mxIsSparse(prhs[1]) || mxIsComplex(prhs[1]) ||
        mxGetNumberOfDimensions(prhs[1])!=2 || mxGetN(prhs[1])!=1 ) {
        mexErrMsgTxt("2nd argument must be real double full column vector");
    }
    m = mxGetM(prhs[0]);
    n = mxGetN(prhs[0]);
    if( mxGetM(prhs[1])!=n ) {
        mexErrMsgTxt("Inner matrix dimensions must agree.");
    }
/* Create output */
    plhs[0] = mxCreateDoubleMatrix( m, 1, mxREAL );
/* Get data pointers */
    Mpr = mxGetPr(prhs[0]);
    Vpr = mxGetPr(prhs[1]);
    Cpr = mxGetPr(plhs[0]);
    Mir = mxGetIr(prhs[0]);
    Mjc = mxGetJc(prhs[0]);
/* Calculate result */
    for( j=0; j<n; j++ ) {
        nrow = Mjc[j+1] - Mjc[j];  /* Number of row elements for this column */
        while( nrow-- ) {
            Cpr[*Mir++] += *Mpr++ * Vpr[j];  /* Accumulate contribution of Vpr[j] */
        }
    }
}

2 个评论
显示无隐藏无

JM 2017-12-1

Thank you so much James, now I understand many things related to the handling of sparse matrices under C mex.

Cheers.

J

Royi Avital 2021-8-26

Isn't there an internal function to do this efficiently (Multi Threaded)?

请先登录，再进行评论。

A sparse matrix vector multiplication in C mex programming

3 个评论
显示 1更早的评论隐藏 1更早的评论

采纳的回答

2 个评论
显示无隐藏无

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

A sparse matrix vector multiplication in C mex programming

3 个评论 显示 1更早的评论隐藏 1更早的评论

采纳的回答

2 个评论 显示 无隐藏 无

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

3 个评论
显示 1更早的评论隐藏 1更早的评论

2 个评论
显示无隐藏无