Most Efficient Way of Using mexCallMATLAB in Converting Double* to mxArray*

I am writing a MEX code in which I need to use pinv function. I am trying to find a way to pass the array of type double to pinv using mexCallMATLAB in the most efficient way. Let's for the sake of example say the array is named G and its size is 100.

    double *G = (double*) mxMalloc( 100 * sizeof(double) );

where

    G[0] = G11; G[1] = G12;
    G[2] = G21; G[3] = G22;

Which means every four consecutive elements of G is a 2×2 matrix. G stores 25 different values of this 2×2 matrix.

I should note that these 2×2 matrices are not well-conditioned and they may contain all zero in their element. How can I use pinv function to calculate the pseudoinverse in the elements of G? For example, how can I pass the array to mexCallMATLAB in order to calculate the pseudoinverse of the first 2×2 matrix in G?

I thought of the following approach:

    mxArray *G_PINV_input     = mxCreateDoubleMatrix(2, 2, mxREAL);
    mxArray *G_PINV_output    = mxCreateDoubleMatrix(2, 2, mxREAL);
    double  *G_PINV_input_ptr = mxGetPr(G_PINV_input);

    memcpy( G_PINV_input_ptr, &G[0], 4 * sizeof(double));
    mexCallMATLAB(1, G_PINV_output, 1, G_PINV_input, "pinv");

I am not sure how good this approach is. Copying the values is not economical at all because the total number of elements in G in my actual application is large. Is there anyway to skip this copying?

Answer 1

[James Tursa]

Preliminary comment: You should not pre-allocate the G_PINV_output, since this will be automatically created by the mexCallMATLAB call itself. Just pass the uninitialized pointer in. They way you have it coded you actually have a memory leak since the result of your mxCreateDoubleMatrix call will get wiped out by the mexCallMATLAB call.

Now on to the real question ...

The answer from an official API function standpoint is no, you cannot do what you are trying to do without copying the data. The only memory address you can officially attach to an mxArray data area (e.g., pr or pi) is one that comes from an an API allocation routine. For example, you could attach G directly to an mxArray data area (e.g. using mxSetPr) since G came directly from an API routine, mxMalloc. That would allow you to pass the first 2x2 matrix to pinv. All fine and good. It doesn't matter that G actually points to a much larger data set ... simply using mxSetM and mxSetN to set a 2x2 size will work.

The problem comes when you try to do the same thing for the next 2x2 block, i.e. starting with G+4. If you pass this address to mxSetPr MATLAB will bomb with an assertion fault. The API routines check the incoming address to make sure they are on an internal list of allocated memory blocks. Since G+4 is not on the list (it is in the middle of an allocated block, but the G+4 address itself is not on the list) a real-time fault will be generated and MATLAB will bomb. So you can't do this neatly in a loop using simple pointer arithmetic and mxSetPr. You need to copy the data. There is a way to get around this by hacking into the mxArray structure itself and setting the pr field directly, bypassing the mxSetPr call entirely (and bypassing the assertion check). However, those details are beyond the scope of this comment and come with their own set of problems to avoid MATLAB crashes, so I will not post them here.

My advice is to simply live with the copying. This is not the bottleneck anyway. The mexCallMATLAB call overhead, plus the creation of the output mxArray variable, plus the pinv function, all will swamp the relatively itty-bitty time it takes to copy 4 doubles. If you really want to speed things up and it is important enough for you, then hand-code the pinv algorithm for the 2x2 yourself inside your mex function and avoid the mexCallMATLAB call completely.

I should also point out this this call has errors:

mexCallMATLAB(1, G_PINV_output, 1, G_PINV_input, "pinv");

it should be this instead:

mexCallMATLAB(1, &G_PINV_output, 1, &G_PINV_input, "pinv");

Since you declared the 2nd and 4th argument variables as mxArray * and not mxArray *[ ], you will need to add the & operator to get the correct level of indirection.