How to loop through matrix elements in mex c++ function for MATLAB?

Question

I am trying to index write an external c++ function for MATLAB using mex to manipulate matrices, and am not able to use multidimensional indexing. There are examples provided here, but I have not found how to fix the problem I describe below. I have a sample matrix:

>> mat
mat =
 1    10
 2    20
 3    30
 4    40
 5    50

Currently I use a linear index through the matrix which works:

#include <mex.h>
#include <iostream>
using namespace std;
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) 
{ 
 //1.get pointer to input graph_list and allocate it    
    double *graph_list = mxGetPr(prhs[0]);
    mwSize mrows = mxGetM(prhs[0]);
    mwSize ncols = mxGetN(prhs[0]);  
    cout<< mrows<<" rows\n";
    cout<< ncols<<" cols\n";
    int mm, nn;
    for (nn=0;nn<ncols;nn++) {
        for (mm=0;mm<mrows;mm++){
            cout << graph_list[nn*(mrows) +mm]  <<"\n";            
        }
    }     
}

This produces:

>> mexTryAlex(mat)    
5 rows
2 cols
1
2
3
4
5
10
20
30
40
50

When I change the definition of graph_list and try 2D indexing to graph_list there is a compilation error with mex:

double **graph_list = mxGetPr(prhs[0]);
cout << graph_list[nn][mm];

EDIT: here is the error message received

>> mex mexTryAlex.cpp
Warning: You are using gcc version "4.4.3-4ubuntu5)".  The version
     currently supported with MEX is "4.3.4".
     For a list of currently supported compilers see: 
     http://www.mathworks.com/support/compilers/current_release/
mexTryAlex.cpp: In function ‘void mexFunction(int, mxArray**, int, const mxArray**)’:
mexTryAlex.cpp:16: error: cannot convert ‘double*’ to ‘double**’ in initialization
mex: compile of ' "mexTryAlex.cpp"' failed.
??? Error using ==> mex at 208
Unable to complete successfully.

Answer 1

The compiler says it all.

In C, a 2D array is like an array of arrays. Therefore, a 2D array is fundamentally different from a 1D array; it's an array of pointers, in which each element contains a pointer to an array (hence a double-pointer, double**).

You're asking mxGetPr() to return a double**, but it returns a double*, e.g., a pointer to the first element of a 1D array. This 1D array can only be indexed linearly.

My guess is that MATLAB does this it this way to keep indexing arrays consistent -- would you really expect/want a double**** for a 4-D array?

Moreover, mxGetPr() cannot be overloaded by return type (it's C after all).

In order to be able to double-index a 1D array, you could sneak in a little macro:

#define A(i,j) A[(i) + (j)*numrows]

and use it like so

double *A = mxGetPr(...);
int numrows = 4;   /* or get with mxGetM() or so) */

double blah = A(3,2); /* call to MACRO */

Obviously, as with all macros, there's a few things to look out for:

1. there's no bounds checking
2. C is 0-based, and Matlab 1-based, making all indices different
3. All arrays will have to be called 'A'

You could write a function to mitigate these drawbacks:

double getValue(double** array, int row, int* dims);

(or use mxCalcSingleSubscript as pointed out by Shai), but that doesn't really improve expressive power IMHO:

double blah = getValue(array, 3,4, dims);
/* or the ugliness from mxCalcSingleSubscript(); */

You could also write in C++, make a Matrix-type class with an operator(), construct it with the pointer and dimensions from mxGetPr() and mxGetDims() etc., compile in Matlab using g++ or equivalent, but that introduces a whole host of other problems and adds way much more complexity than needed for most cases.

Therefore, to avoid all of this mess, I just always compute the index in-place :)

Answer 2

Having a matrix class is the by far easiest way to deal with this kind of problems. There are plenty to choose from so don't bother writing your own. Armadillo is fairly good and also integrates with LAPACK if you use that. http://arma.sourceforge.net/docs.html

See example below

#include <mex.h>
#include <iostream>
#include <armadillo>
using namespace std;
using namespace arma;

//creates an armadillo matrix from a matlab matrix
mat armaMatrix(const mxArray *matlabMatrix[]){
    mwSize mrows = mxGetM(matlabMatrix[0]);
    mwSize ncols = mxGetN(matlabMatrix[0]);
    double *values = mxGetPr(matlabMatrix[0]);

    return mat(values, nrows, ncols);
}

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) 
{
    mat graph_list = armaMatrix(prhs);

    //print the matrix
    cout << graph_list<<"\n";
    //print the first column
    cout << graph_list(span::all,0) <<"\n";
}