Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

Problems with mxSetPr in a mex file

53 views
Skip to first unread message

Sira

unread,
Dec 30, 2011, 8:40:09 AM12/30/11
to
Hello, I am a newbie in the use of mex. This is my problem, I would like to access to a complex matrix S as in matlab S(:,:,2), but without copying all the data. To access position 1 in C of the second parameter of my code, I wrote the following code:

double* SIPR = mxGetPr(prhs[1]);
double* SIPI = mxGetPi(prhs[1]);

mxArray *si_c = mxCreateDoubleMatrix(n, n, mxCOMPLEX);
mxSetPr(si_c,SIPR+ 1*n*n);
mxSetPi(si_c, SIPI+ 1*n*n);


note that the matrix has dimension n by n.
After executing this code, the values of the input matrix prhs[1] change and the code explodes. What am I doing wrong? Is there a easier way to do this?

Thank you for your response,
Sira.







James Tursa

unread,
Dec 30, 2011, 1:38:08 PM12/30/11
to
"Sira" wrote in message <jdkevp$1i6$1...@newscl01ah.mathworks.com>...
1) You can't attach the data pointers of one mxArray to another mxArray like that. It will confuse the MATLAB memory manager since it no longer knows where all of the shared data copies are.

2) You can't legally attach *any* data pointer to an mxArray that was not the result of some MATLAB memory manager API function. In your case, SIPR and SIPI are addresses that were the result of a MATLAB memory manager call, but SPIR + 1*nn and SIPI + 1*nn are NOT. On older versions of MATLAB the code you have will execute but downstream MATLAB will bomb as it tries to free memory starting at address SIPR + 1*nn and SIPI + 1*nn when the variable is cleared. On later versions of MATLAB your program will bomb at the mxSetPr statements themselves with an assertion type of fault since it checks the addresses passed to it in the mxSetXX calls.

3) Even if you could do what you are trying to do, your use of mxSetPr and mxSetPi will cause a memory leak since these functions do NOT free up the memory behind the current data pointers before they set the new ones.

You are apparently trying to create an mxArray that has data pointers that point to the interior data of another mxArray. You can't legally do that. There *are* ways to fake MATLAB out temporarily to get away with this, but they involve hacking into the mxArray structure itself and manually tinkering with the fields directly. And even this approach can bomb MATLAB unless you really know what you are doing and are very careful. E.g., see this FEX submission by Bruno Luong:

http://www.mathworks.com/matlabcentral/fileexchange/24576-inplacearray-a-semi-pointer-package-for-matlab

Q: Why are you doing this? Are you trying to gain some speed advantage for some calculations, or what? Maybe there is another way to go about your particular problem that accesses the 2D pages of your variable without actually forming a separate mxArray from them.

James Tursa

Sira

unread,
Dec 30, 2011, 3:56:08 PM12/30/11
to
"James Tursa" wrote in message <jdl0eg$1fn$1...@newscl01ah.mathworks.com>...
The reason to do this is that the matrix in prhs[1] has size nxnxM and I would like to access to these 2D matrices of size nxn without copying the whole matrix, and this is why the position of :
SIPR+ 1*n*n
should correspond to the position S(:,:,2). So, if there is any other way of doing this it would be perfect to me.

James Tursa

unread,
Dec 30, 2011, 4:36:08 PM12/30/11
to
"Sira" wrote in message <jdl8h8$rg8$1...@newscl01ah.mathworks.com>...
> "James Tursa" wrote in message <jdl0eg$1fn$1...@newscl01ah.mathworks.com>...
> >
> > Q: Why are you doing this? Are you trying to gain some speed advantage for some calculations, or what? Maybe there is another way to go about your particular problem that accesses the 2D pages of your variable without actually forming a separate mxArray from them.
>
> The reason to do this is that the matrix in prhs[1] has size nxnxM and I would like to access to these 2D matrices of size nxn without copying the whole matrix, and this is why the position of :
> SIPR+ 1*n*n
> should correspond to the position S(:,:,2). So, if there is any other way of doing this it would be perfect to me.

You didn't answer my question. You simply told me again *what* you were trying to do. I already knew that. I am asking *why* you want to do this. Since there is no way to legally do what you are attempting to do in a mex routine, I would like to know why you are trying to do it in the first place. Why can't you operate on the 2D pages directly in your mex routine without making mxArrays out of them? Why do you need separate mxArrays of the 2D pages to work with? Are you trying to pass them back into the MATLAB workspace? Are you trying to call a MATLAB function with them as input(s)? Or what?

James Tursa

Sira

unread,
Dec 31, 2011, 6:56:07 AM12/31/11
to
"James Tursa" wrote in message <jdlas8$5fm$1...@newscl01ah.mathworks.com>...
Sorry, I didn't understand your question. I need to implement the following formula:
B = sum_i (A*S(:,:,i)*A)^0.5

where all matrices A,B,S are complex. Because I will need to iterate it, I would like it to be as fast as possible, so I thought the best option was to use the BLAS method ZHEMM (the matrices are also Hermitian), but to do that I need to isolate the 2D pages....

James Tursa

unread,
Jan 1, 2012, 2:59:09 AM1/1/12
to
"Sira" wrote in message <jdmt8n$qho$1...@newscl01ah.mathworks.com>...
You will NOT be able to call the BLAS ZHEMM routine without making a copy of your data. That is because all of the complex BLAS (and LAPACK) routines expect the complex data to be interleaved with the real data ala Fortran. MATLAB does not store real and imaginary data in interleaved fashion, but has separate memory blocks for the real and imaginary parts. So either you will need to copy your data into new memory blocks and call the complex BLAS routines, or perhaps call the real BLAS routines with the various real and imaginary pieces and combine the results into the complex result manually. E.g., you could do this:

mtimesx(mtimesx(A,S),A)

which would compute the 2D slice results via BLAS routines (it calls the real BLAS routines multiple times with the real and imaginary pieces and combines the results into the complex result). The mtimesx routine can be found on the FEX here:

http://www.mathworks.com/matlabcentral/fileexchange/25977-mtimesx-fast-matrix-multiply-with-multi-dimensional-support

As to the rest of the calculation, may I ask what the dimensions of all of your variables are? And could you be explicit as to whether the sum is inside the sqrt or ourside of the sqrt?

In any event, there is absolutely no reason to construct mxArrays of the 2D slices in order to calculate your desired result. You can do everything inside the mex routine without making data copies. I can help you with that once I know the dimensions involved and the order of the calculations.

James Tursa

Sira

unread,
Jan 5, 2012, 5:20:08 AM1/5/12
to
"James Tursa" wrote in message <jdp3od$esi$1...@newscl01ah.mathworks.com>...
Ok, understood, I will use mtimesx.
Regarding the formula B = sum_i [ (A*S(:,:,i)*A)^0.5 ], the square root is inside the sum and the dimensions of A are nxn, where n can be 1 or 3 and S is nxnxM, where the value of M depends but it is normally between 2 and 15. So, the space is not as important as the time it takes to do the computations, given that I will need to iterate this formula 20000 times, more specifically, I want to implement:

for j=1:20000
B = sum_i [ (A*S(:,:,i)*A)^0.5 ];
A = B^0.5;
end

For the pow(XX,0.5), is there a fast way of executing it in mex? I was planning to call the mexcallMATLAB, something like:

mxArray * squareRoot( mxArray *Input){
//inputs
mxArray* params[2];
params[0] = Input;
params[1] = mxCreateDoubleScalar(0.5); //power=0.5

mwSize n = mxGetM(Input);
//outputs
mxArray *Output[1];
Output[0] = mxCreateDoubleMatrix(n, n + n - 1, mxCOMPLEX);

/* get the square root */
mexCallMATLAB(1, Output, 2, params, "mpower");
return(Output[0]);
}

and then sum along "i" in the formula and every position using a regular C loop.

James Tursa

unread,
Jan 5, 2012, 11:27:07 AM1/5/12
to
"Sira" wrote in message <je3tgo$okp$1...@newscl01ah.mathworks.com>...
>
> Regarding the formula B = sum_i [ (A*S(:,:,i)*A)^0.5 ], the square root is inside the sum and the dimensions of A are nxn, where n can be 1 or 3 and S is nxnxM, where the value of M depends but it is normally between 2 and 15. So, the space is not as important as the time it takes to do the computations, given that I will need to iterate this formula 20000 times, more specifically, I want to implement:
>
> for j=1:20000
> B = sum_i [ (A*S(:,:,i)*A)^0.5 ];
> A = B^0.5;
> end

Just to be clear, for the n=3 case you are doing a matrix square root of a 3x3 for each index i? If speed is your ultimate goal, then I would probably opt to explicitly code the matrix multiply A*S(:,:,i)*A element-by-element in the C routine as this will likely be faster than anything else such as BLAS calls for only a 3x3 size. E.g., for the n=1 case you don't even need to use mexCallMATLAB to get the result. For the n=3 case you can do the matrix multiply explicitly, and use mexCallMATLAB on those results to get the matrix square roots. Which variables are complex?

James Tursa

Sira

unread,
Jan 5, 2012, 11:45:08 AM1/5/12
to
"James Tursa" wrote in message <je4j0r$62b$1...@newscl01ah.mathworks.com>...
> "Sira" wrote in message <je3tgo$okp$1...@newscl01ah.mathworks.com>...
> >
> > Regarding the formula B = sum_i [ (A*S(:,:,i)*A)^0.5 ], the square root is inside the sum and the dimensions of A are nxn, where n can be 1 or 3 and S is nxnxM, where the value of M depends but it is normally between 2 and 15. So, the space is not as important as the time it takes to do the computations, given that I will need to iterate this formula 20000 times, more specifically, I want to implement:
> >
> > for j=1:20000
> > B = sum_i [ (A*S(:,:,i)*A)^0.5 ];
> > A = B^0.5;
> > end
>
> Just to be clear, for the n=3 case you are doing a matrix square root of a 3x3 for each index i?

yes, it is not element by element, but a power function, I believe in matlab for a 3x3 matrix M, it should be mpower(M,0.5).

>If speed is your ultimate goal, then I would probably opt to explicitly code the matrix multiply A*S(:,:,i)*A element-by-element in the C routine as this will likely be faster than anything else such as BLAS calls for only a 3x3 size. E.g., for the n=1 case you don't even need to use mexCallMATLAB to get the result. For the n=3 case you can do the matrix multiply explicitly, and use mexCallMATLAB on those results to get the matrix square roots. Which variables are complex?
>
> James Tursa
All of the matrices are complex, A,S and B. Thanks a lot for your help!
0 new messages