Posted By

benjamin on 10/29/09


Tagged

c matlab IT numerical computing armadillo itpp mex


Versions (?)

Fast Scientific Computation in Mex Functions Using Armadillo


 / Published in: C++
 

URL: http://www.myoutsourcedbrain.com/2009/08/fast-scientific-computation-in-mex.html

This code is a demonstration of how to do fast scientific computation in mex functions using the armadillo library. See my blog post for details and leave comments there.

  1. #include "mex.h"
  2. #include "math.h"
  3. #include<armadillo>
  4.  
  5. using namespace arma;
  6.  
  7. void matlab2arma(mat& A, const mxArray *mxdata){
  8. // delete [] A.mem; // don't do this!
  9. access::rw(A.mem)=mxGetPr(mxdata);
  10. access::rw(A.n_rows)=mxGetM(mxdata); // transposed!
  11. access::rw(A.n_cols)=mxGetN(mxdata);
  12. access::rw(A.n_elem)=A.n_rows*A.n_cols;
  13. };
  14.  
  15. void freeVar(mat& A, const double *ptr){
  16. access::rw(A.mem)=ptr;
  17. access::rw(A.n_rows)=1; // transposed!
  18. access::rw(A.n_cols)=1;
  19. access::rw(A.n_elem)=1;
  20. };
  21.  
  22. void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
  23. {
  24. if (nrhs != 2)
  25. mexErrMsgTxt("Incorrect number of input arguments");
  26. if (nlhs != 1)
  27. mexErrMsgTxt("Incorrect number of output arguments");
  28.  
  29. mat D1(1,1);
  30. const double* D1mem=access::rw(D1.mem);
  31. matlab2arma(D1,prhs[0]); // First create the matrix, then change it to point to the matlab data.
  32.  
  33. mat D2(1,1);
  34. const double* D2mem=access::rw(D2.mem);
  35. matlab2arma(D2,prhs[1]);
  36.  
  37. // check if the input corresponds to what you are expecting
  38. if( D1.n_rows != D2.n_rows )
  39. mexErrMsgTxt("Columns of D1 and D2 must be of equal length!");
  40.  
  41. if( D1.n_cols != D2.n_cols )
  42. mexErrMsgTxt("Rows of D1 and D2 must be of equal length!");
  43.  
  44. plhs[0] = mxCreateDoubleMatrix(D1.n_rows, D1.n_cols, mxREAL);
  45. mat output(1,1);
  46. const double* outputmem=access::rw(output.mem);
  47. matlab2arma(output,plhs[0]);
  48.  
  49. output=D1+D2;
  50. // output.print();
  51.  
  52. freeVar(D1,D1mem); // Change back the pointers!!
  53. freeVar(D2,D2mem);
  54. freeVar(output,outputmem);
  55. return;
  56. }

Report this snippet  

You need to login to post a comment.