CODEWITHSUNDEEP
c++matlabmex
Jame
Jame

Reputation: 3854

How to represent a 2D array of MATLAB in mex code

I have a MATLAB code which is

%% Inputs are theta and h (size NxM)
alpha=zeros(N,M);
h_tmp=zeros(N,M);
h_tmp(1:N-1,:)=h(2:N ,:);
for i=1:N
    alpha(i,:)=theta.*(h_tmp(i,:)+h(i,:));
end

By using vectorized method, the above code can be

alpha = theta .* [h(1:N-1,:) + h(2:N,:); h(N,:)];

To speed up the code, I want to rewrite it in MEX file using C++. The main different between MATLAB and C++ in 2D array is row-major order (MATLAB) and column-major order(C++)

double  *h, *alpha, *h_temp;
int N,M;
double theta;    
N      = (int) mxGetN(prhs[0]); //cols
M      = (int) mxGetM(prhs[0]); //rows
h      = (double *)mxGetData(prhs[0]);
theta  = (double)*mxGetPr(prhs[1]);
/* Initial zeros matrix*/
plhs[0]   = mxCreateDoubleMatrix(M, N, mxREAL);  alpha = mxGetPr(plhs[0]);
//////////////Compute alpha/////////    
for (int rows=0; rows < M; rows++) {
    //h[N*rows+cols] is h_tmp
    for (int cols=0; cols < N; cols++) {        
         alpha[N*rows+cols]=theta*(h[N*rows+cols+1]+h[N*rows+cols]);
    }
}

Are my Mex code and MATLAB code equivalent? If not, could you help me to fix it?

Upvotes: 1

Views: 1069

Answers (1)

user2457516
user2457516

Reputation:

Besides the corrections from the comments to your question, there is one minor difference. What is missing is that you skip h(N,:) in the Matlab code, where in the C code iteration over the code is done until cols < N, which (due to 0 indexing in C) will also process the last element in each column.

#include "mex.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
    double  *h, *alpha, *h_temp;
    int num_columns, num_rows;
    double theta;    
    num_columns = (int) mxGetN(prhs[0]); //cols
    num_rows    = (int) mxGetM(prhs[0]); //rows
    h           = (double *)mxGetData(prhs[0]);
    theta       = (double)*mxGetPr(prhs[1]);
    /* Initial zeros matrix*/
    plhs[0]   = mxCreateDoubleMatrix(num_rows, num_columns, mxREAL);  alpha = mxGetPr(plhs[0]);
    //////////////Compute alpha/////////
    // there are num_rows many elements in each column
    // and num_columns many rows. Matlab stores column first.
    // h[0] ... h[num_rows-1] == h(:,1)
    int idx; // to help make code cleaner
    for (int column_idx=0; column_idx < num_columns; column_idx++) {
        //iterate over each column
        for (int row_idx=0; row_idx < num_rows-1; row_idx++) {// exclude h(end,row_idx)
            //for each row in a column do
            idx = num_columns * column_idx + row_idx;
            alpha[idx]= theta * (h[idx+1] + h[idx]);
        }
    }
    //the last column wasn't modified and set to 0 upon initialization.
    //set it now
    for(int rows = 0; rows < num_rows; rows++) {
        alpha[num_columns*rows+(num_rows-1)] = theta * h[num_columns*rows+(num_rows-1)];
    }
}

Note that I decided to rename some of the variables, so that I think it becomes easier to read.

Edit: Removed the suggestion with prhs[0] = plhs[0], as suggested in the comments to this answer. One may get away with this under certain circumstances, but in general it is no good practice when coding matlab .mex functions and it may crash Matlab.

Upvotes: 1

Related Questions