[SciPy-Dev] f2py: questions on array arguments

[Sturla Molden]

Ramon Crehuet skrev:
>
> 1. The follwing fortran function:
>
> function f2(x,y)
> implicit none
> real,intent(in):: x,y
> real, dimension(3):: f2
> f2(1)=x+y**2
> f2(2)=sin(x*y)
> f2(3)=2*x-y
> end function f2
>
> gives a segmentation fault when called from python if it is not in a
> fortran module. If it is contained in a fortran module, it works fine
> and returns an array. That makes sense because fortran modules
> automatically generate an interface. However, I don't see that reflected
> in the .pyf files generated by f2py. So, is there a way to "correct" the
> function outside the module to work with f2py?

Wrap it in a subroutine f2wrap(x,y,output) and call that. C and Fortran 
interpop (which is what f2py does) can only go so far. C has no concept 
of functions returning arrays. The Fortran function must know where to 
write the return values, but how can you provide that information from 
C? So Fortran and C has different calling conventions here.

And the reason for this working inside a module is that f2py must write 
a Fortran wrapper for the module. With modules there are name mangling 
and stuff going on that are not interopable with C. But then you have a 
situation where f2 is called from Fortran, so it will not segfault.

>
> 2. I have read in several posts that automatic arrays do not work with
> f2py. So that something like:
>
> real function trace(m)
> real, dimension(:,:), intent(in) :: m
>
> Has to be converted into:
>
> real function trace2(m,n)
> integer, intent(in) :: n
> !f2py integer,intent(hide),depend(m) :: n=shape(m,0)
> real, dimension(n,n), intent(in) :: m
>

Again, here you have a C vs. Fortran incompatibility again. An assume 
shaped array like "m" in "trace" is actually a dope vector, i.e. a C 
struct not different from NumPy's PyArrayObject (though the binary 
layout is compiler dependent, and incompatible with NumPy). So we cannot 
pass a C pointer to "trace" and hope that it will work. But if we do 
know the C defintion of the Fortran compiler's dope vector, we can fill  
in the fields from a PyArrayObject and pass a pointer to the struct 
(there is a C library that attempts to do that, but I don't recall the 
name).

But in "trace2", the shape of "m" is explicit. Here Fortran will assume 
m is just a C pointer to the first element. Fortran 77 compilers usually 
treat explicit and assumed-ashed arrays differently for interopability 
with C and compatilibity with other F77 compilers. So in "trace2", we 
can pass a C pointer to the first element of "m", whereas in "trace" we 
cannot.

This difference is not required by the standard though. A Fortran 90 
compiler could use only dope vectors if it wanted to. (There are other 
non-standard ways of interfacing with C such as Cray pointers.) So this 
is a bit murky territory.

With Fortran 2003 we finally have a portable Fortran to C interface. Not 
all Fortran compilers implement this; gfortran does, whereas my 
favourite compiler Absoft does not.

Sturla