No fixed arrays in Fortran / Mex possible?
Roland Neumann
I cannot access t then using t(2,1) or other indices (Compiler error). If I define it double precision, allocatable,dimension(:,:) :: t then I get an segmentation fault with the mxCopyPtr... command. Is there a mxCopyPtrToMtx command?
Seemingly at the end of my skills here, I just want to get a matrix like 10000x30 from Matlab command in a mex Fortran file. I could dimension it with real*8 t(30000) but then it's too large for 20x10 and too small for 20000x30. All online help examples just define it fix and never access like t(2,1), Why?
NZTideMan
Not sure what you're doing wrong, but this is how I do it:
real*8,allocatable :: y(:,:)
-
-
! Bring array into Fortran
n=mxGetM(prhs(1)) ! No of data
ncol=mxGetN(prhs(1)) ! No of columns of data
allocate (y(n,ncol))
ptr=mxGetPr(prhs(1))
call mxCopyPtrtoReal8(ptr,y,n*ncol)
James Tursa
> I understand that you have to call mxCopyPtrToReal8(t, tr, size) to copy a given matrix from a matlab command to the fortran mex file. Though it seems not allowed to define t as real*8 t(*).
>
> I cannot access t then using t(2,1) or other indices (Compiler error). If I define it double precision, allocatable,dimension(:,:) :: t then I get an segmentation fault with the mxCopyPtr... command. Is there a mxCopyPtrToMtx command?
>
> Seemingly at the end of my skills here, I just want to get a matrix like 10000x30 from Matlab command in a mex Fortran file. I could dimension it with real*8 t(30000) but then it's too large for 20x10 and too small for 20000x30. All online help examples just define it fix and never access like t(2,1), Why?
There are basically two ways to do it, depending on whether you want to use the matrix as a read-only matrix, or be able to modify it inside the mex routine. You had the right idea with allocatable, but you must allocate it manually before copying into it. That was why you got a seg fault. Here are the two methods (caveat, I am typing this off the top of my head and don't have access to Fortran on the computer I am using right now ... so code is not tested yet. Also, this is bare bones and I have not included any argument checking etc.):
1) Copying the input to an allocatable matrix, which can then be changed if desired.
#include "fintrf.h"
subroutine mexFunction(nlhs, plhs, nrhs, prhs)
integer*4 nlhs, nrhs
mwPointer plhs(*), prhs(*)
mwPointer mxGetPr
integer*4 mxGetM, mxGetN
real*8, allocatable :: t(:,:)
integer*4 m, n
m = mxGetM(prhs(1))
n = mxGetN(prhs(1))
allocate(t(m,n))
call mxCopyPtrToReal8(mxGetPr(prhs(1)), t, m*n)
:
! use t here
:
deallocate(t)
:
end subroutine mexFunction
2) Using the input as a read-only matrix. No copying required (faster).
#include "fintrf.h"
subroutine mexFunction(nlhs, plhs, nrhs, prhs)
integer*4 nlhs, nrhs
mwPointer plhs(*), prhs(*)
mwPointer mxGetPr, pr
integer*4 m, n
m = mxGetM(prhs(1))
n = mxGetN(prhs(1))
pr = mxGetPr(prhs(1))
call sub(%VAL(pr), m, n)
:
end subroutine mexFunction
subroutine sub(t,m,n)
integer*4 m, n
real*8, intent(in) :: t(m,n)
:
! use t here
:
return
end subroutine sub
James Tursa
Roland Neumann
real*4, allocatable :: weight(:)
The call didnt get the right result in weight:
allocate(weight(zeil))
tempzeiger = mxGetPr(prhs(5))
call mxCopyPtrToReal4(tempzeiger,weight, zeil)
Strangely, this code works:
allocate(weight(zeil)) !is real*4
allocate(tempvekt(zeil)) !is real*8
tempzeiger = mxGetPr(prhs(5))
call mxCopyPtrToReal8(tempzeiger, tempvekt, zeil)
weight = tempvekt
deallocate(tempvekt)
Can't I use the mxCopy Functions with other than real8? Maybe the answer to that problem points me to my current problem, after some fortran subroutine calls the correct values of weight disappear. So any hint would be appreciated
NZTideMan
My version of Matlab (2006a) does not have a mxCopyPtrToReal4
subroutine available, so I always use the second version that you have
given. I don't think that is "strange" at all.
James Tursa
> Thank you both. It works like a charm. Unfortunately only until the fortran routine (from netlib.org) is called. It expects real*4, so I defined:
> real*4, allocatable :: weight(:)
> The call didnt get the right result in weight:
> allocate(weight(zeil))
> tempzeiger = mxGetPr(prhs(5))
> call mxCopyPtrToReal4(tempzeiger,weight, zeil)
All of the mxCopy___ functions are simple memory copy functions. *No* type conversion of any kind is performed. So it is not at all surprising that the above code, where you are doing a bit-for-bit memory copy from a real*8 (i.e., the double from MATLAB) into a real*4, gives a garbage result.
> Strangely, this code works:
> allocate(weight(zeil)) !is real*4
> allocate(tempvekt(zeil)) !is real*8
> tempzeiger = mxGetPr(prhs(5))
> call mxCopyPtrToReal8(tempzeiger, tempvekt, zeil)
> weight = tempvekt
> deallocate(tempvekt)
Not strange at all. Here your are copying the real*8 from MATLAB (the double) into a real*8 Fortran variable. The bit-for-bit copy works because they are both exactly the same type. Then you assign tempvekt, your real*8 Fortran variable, to weight, a real*4 Fortran variable. This works fine because the Fortran language itself is doing the type conversion for you as part of the assignment statement. But you are doing too much work here. The first copy is completely wasteful and not necessary ... you don't need this copy at all. Better to do just one copy/conversion as follows:
#include "fintrf.h"
subroutine mexFunction(nlhs, plhs, nrhs, prhs)
integer*4 nlhs, nrhs
mwPointer plhs(*), prhs(*)
mwPointer mxGetPr, pr
integer*4 mxGetM, mxGetN
integer*4 m, n
real*4, allocatable :: weight(:,:)
m = mxGetM(prhs(1))
n = mxGetN(prhs(1))
pr = mxGetPr(prhs(1))
allocate(weight(m,n))
call r4eqr8(weight, %VAL(pr), m, n)
:
! use weight here
:
deallocate(weight)
end subroutine mexFunction
!-----
subroutine r4eqr8(r4, r8, m, n)
integer*4 m, n
real*4 r4(m,n)
real*8, intent(in) :: r8(m,n)
r4 = r8
return
end subroutine r4eqr8
> Can't I use the mxCopy Functions with other than real8? Maybe the answer to that problem points me to my current problem, after some fortran subroutine calls the correct values of weight disappear. So any hint would be appreciated
You are going to have to be more explicit here. What do you mean "disappear"? Can you show some code that demonstrates this? Like I said before, all of the mxCopy___ functions are simple memory copy functions, nothing else. Are you trying to copy the results of your Fortran code back into a plhs(1) that is returned to MATLAB? If so, you will have to do the reverse of the conversion I have shown above ... a simple memory copy will not work. You will need a conversion as well if you want to return a double type.
James Tursa
James Tursa
>
> My version of Matlab (2006a) does not have a mxCopyPtrToReal4
> subroutine available
Since all the mxCopy___ functions are just memory copy functions, you can use mxCopyPtrToInteger4 instead (if it is available). This will copy single variables into real*4 variables just fine. e.g., I routinely use mxCopyPtrToReal8 (and its complement) for integer*8 variables since there is no mxCopyPtrToInteger8 function available.
James Tursa
Roland Neumann
> ...All of the mxCopy___ functions are simple memory copy functions. *No* type conversion of any kind is performed...
This was the most important part, which is completely missing in the Matlab documentation. Thank you very much for it.
Besides, I found the solution for my disappearing parameters problem... I made a very stupid error. The subroutine call was 2lines, I forgot the , after the last parameter in the 1st line. So it assumed it was a new variable (I know now, implicit none...) with just a strange name and one missing...
Anyway, two things (am I allowed to discuss it here or each in a new thread to better be found for others?): The %VAL() construct doesn't compile in Ifort 9.1.
Second, the statistics code I call contains lots of nice outputs write(6, 11). Is there a way to display them in matlab or in a file at least? I tried to define
>character*100 ::str
in each subroutine. Then I replaced the 6 above with str and put a mexWarnMsgTxt(str). It leads to Matlab collapse with some kind of overflow though. But till this point I'm impressed how fast fortran runs. Thank you very much, Roland
James Tursa
>
> Anyway, two things (am I allowed to discuss it here or each in a new thread to better be found for others?): The %VAL() construct doesn't compile in Ifort 9.1.
The %VAL() construct compiles just fine in Ifort 9.1 ... I have that exact compiler. %VAL() has been an unofficial part of most (all?) Fortran compilers for probably decades. My guess is you are misusing it.
Remember for most practical purposes Fortran is basically a pass-by-reference language (as opposed to C which is pass-by-value, e.g.), so when you make any type of function or subroutine call it is the addresses of the arguments that get passed to the routine (in C, copies of the arguments themselves are passed). But sometimes you want the actual value of the argument to be passed, not the address of the argument. That is the case with the example I posted earlier. The result of a mxGetPr call is an address which gets stored in the pr variable (equivalent to a C pointer). I want that actual address value contained in pr to be passed to the subroutine, not the address of the pr variable itself. So I use the %VAL() construct. The restriction is that %VAL() can *only* be used in the argument list itself. i.e., you cannot use it in any other context like an assignment (e.g., A = %VAL(pr)
would be invalid), because it is strictly a compiler instruction on how to pass arguments. If you can't get it to work feel free to post the offending code and I can look at it.
>
> Second, the statistics code I call contains lots of nice outputs write(6, 11). Is there a way to display them in matlab or in a file at least? I tried to define
> >character*100 ::str
> in each subroutine. Then I replaced the 6 above with str and put a mexWarnMsgTxt(str). It leads to Matlab collapse with some kind of overflow though. But till this point I'm impressed how fast fortran runs. Thank you very much, Roland
You're stuck. It is a pain in the neck to try to modify a large chunk of Fortran code and line-by-line change every write(6,*) ___ into a write(str,*) ___ and mexPrintf() combination, isn't it? My guess is you don't need this output interactively, so I would suggest the following:
Put this statement at the front of your mexFunction routine:
integer*4 mexPrintf, k
character*200 ::str
:
open(unit=6,file='temp.txt',status='unknown')
and leave all of your write(6,__) statements as-is. Then just before mexFunction returns back to MATLAB, do this:
rewind 6
do
read(6,'(a)',end=500) line
k = mexPrintf(line//achar(13))
enddo
500 close(unit=6)
James Tursa
James Tursa
> do
> read(6,'(a)',end=500) line
> k = mexPrintf(line//achar(13))
> enddo
> 500 close(unit=6)
>
Typo. That last part should have str, not line:
rewind 6
do
read(6,'(a)',end=500) str
k = mexPrintf(str//achar(13))
Roland Neumann
> Typo. That last part should have str, not line:
Until now, I learnt so much from you, I even found that before. Thank you so much. Especially the output redirection was a good idea. Everything runs fine (hopefully in the future aswell).
Roland Neumann
My Fortran routine returns a variable a(x, y, z). y is an positive integer. How can I return its value to Matlab? mxCreateArray only does 2 dims. I wouldn't like to return just the part a(x,z) for y=1. With this, I also get a nice crash where whole Matlab just disappears. Segmentation fault was yesterday...
Steven Lord
Try mxCreateNumericArray:
http://www.mathworks.com/access/helpdesk/help/techdoc/apiref/mxcreatenumericarray.html
mxCreateNumericMatrix is limited to 2D. mxCreateNumericArray is not, as
shown in the example on that reference page.
--
Steve Lord
sl...@mathworks.com
James Tursa
> Now I came to a interesting point. Three-dimensional arrays....
>
> My Fortran routine returns a variable a(x, y, z). y is an positive integer. How can I return its value to Matlab? mxCreateArray only does 2 dims. I wouldn't like to return just the part a(x,z) for y=1. With this, I also get a nice crash where whole Matlab just disappears. Segmentation fault was yesterday...
In addition to Steven Lords's comments, I would point out that you again have the choice of copying your Fortran variable into the mxArray plhs variable, OR you can create your mxArray plhs variable up front and use the %VAL() construct to pass the mxGetPr(plhs(1)) result to your routine. That way the results automatically get stored in the return plhs variable and you will avoid the extra unnecessary copy.
James Tursa
Roland Neumann
Thomas Clark
Thanks!
Tom Clark