[SciPy-Dev] f2py, the fortran integer type, and npy_intp

[Kurt Smith]

On Sat, Jul 10, 2010 at 6:33 PM, Sturla Molden <sturla@molden.no> wrote:
> Charles R Harris skrev:
>>
>> I don't think the bindings don't help, we need to have the default
>> integers in the existing pre-FORTRAN77 code in scipy compiled as
>> Py_ssize_d, and then f2py needs to be modified to generate appropriate
>> Python bindings. That's a lot of work. Even with the bindings I think
>> one would need to have a script to rewrite the FORTRAN code since the
>> c-type corresponding to Py_ssize_d isn't fixed.
> We should also beware that the problem applies to real numbers as well
> as integer. You cannot rely on standard mappings from C float to REAL
> and C double to DOUBLE PRECISION.
>
> With pre-Fortran 90, there is no way of controlling this portably. With
> Fortran 90 and later, we have the standard methods selected_real_kind
> and selected_int_kind that returns (compiler dependent) "kind" numbers,
> which can be used to declare real and integers with specific precitions.
>
>  integer, parameter :: single = selected_real_kind(p=6, r=37)
>  integer, parameter :: double = selected_real_kind(p=13)
>
>  integer, parameter :: npy_int8 = selected_int_kind(2)
>  integer, parameter :: npy_int16 = selected_int_kind(4)
>  integer, parameter :: npy_int32 = selected_int_kind(9)
>  integer, parameter :: npy_int64 = selected_int_kind(18)
>
> Now we can declare an npy_int32 like this:
>
>  integer(kind=npy_int32) :: i
>
> Still we have no ide what npy_intp would map to. We can do this in
> Fortran 2003:
>
>  use, intrinsic ::  iso_c_binding
>  integer, parameter :: npy_intp = c_intptr_t

Trouble is that 'c_intptr_t' isn't defined for some fortran compilers.
 So you first have to find which one of 'int', 'long int' or 'long
long int' corresponds to 'npy_intp' in C code, then declare a
parameter in fortran:

integer, parameter :: fort_npy_intp = c_long

for example.  Not a big deal.

>
>  integer(kind=npy_intp) ::  i
>
> Real numers and other integer sre also easier:
>
>  integer, parameter :: float = c_float
>  integer, parameter :: double = c_double
>  integer, parameter :: npy_int32 = c_int32_t
>
> Such declarations can be put in a module, and subsequently imported to
> Fortran 90.
>
> It might be that f2c is the only cure for old Fortran code. The other
> option is to write a Fortran 2003 wrapper to interface C. Then this
> wrapper will call the old Fortran code. We then need to declare the old
> Fortran routine (as an interface block) to Fortran 2003. The Fortran
> compiler is then smart enough to do the correct conversions, including
> making a local copy of an array if that is needed.

Umm, unless I misunderstand you, this isn't how many fortran compilers
behave (ifort, gfortran & g95).

Try the following:

subroutine outer()
    integer(kind=8) :: arr(10,10)

    call inner(arr)

    contains

    subroutine inner(arr)
        integer(kind=4), intent(inout) :: arr(:,:)
        arr = 1
    end subroutine
end subroutine

Compiling with gfortran 4.4.4:

    call inner(arr)
               1
Error: Type mismatch in argument 'arr' at (1); passed INTEGER(8) to INTEGER(4)

But perhaps I misunderstand you.

>
> Wasn't Kurt Smith working on this for a GSOC project?

Yes :-)

It is coming along, and being actively developed.  Many niggling
issues to work out.  We're at near-parity with f2py at the moment,
with a few enhancements, but no 'pyf' interface modules (yet).

http://fwrap.sourceforge.net/

http://fortrancython.wordpress.com/

Kurt