0

I am working to encapsulate a C library using f2py (to later use them in python/numpy). After getting help here to get the first encapsulation working I am now trying to encapsulate more complex types, specifically C functions returning arrays allocated via malloc.

This is the C library mockup (c_lib.c):

#include <stdlib.h>

double multiply_numbers(double a, double b);
double *get_array(double a);

double multiply_numbers(double a, double b) {
        return a*b;
}

double *get_array(double a) {
        double *retval;
        int i;

        retval = malloc(100*sizeof(a));

        for(i=0; i<100; i++) {
                retval[i] = i*a;
        }

        return retval;
}

I can succesfully encapsulate the basic multiply_numbers function (f_mod.f90):

module test_c_lib

        use iso_c_binding
        implicit none
        
contains

        subroutine multiply(x, y, z)

                use iso_c_binding

                real(8), intent(in)        :: x
                real(8), intent(in)        :: y
                real(8), intent(out)       :: z
        
                ! Interface to C function
                interface                        
                        real(c_double) function c_multiply_numbers(a, b) bind(C, name="multiply_numbers")
                                import
                                real(c_double), value :: a,b
                        end function
                end interface

                ! Call C function               
                z = c_multiply_numbers(x,y)                

        end subroutine
        
end module

I compile everything using this makefile:

f_mod.so:       f_mod.f90 c_lib.o
                f2py -c f_mod.f90 c_lib.o -m f_mod

c_lib.o:        c_lib.c
                gcc -c -fpic c_lib.c -o c_lib.o

But I am struggling to understand how would I code the subroutine that wraps the C function get_array.

Specifically,

  • how do I receive a 100 items array?
  • which would be the interface definition?

Note: I have not tried anything because I do not where to start.

EDIT

With francescalous comments I have been able to encapsulate the function and make it work. It is worth to note that f2py might require the extra step of assigning the Fortran pointer to a Fortran array. I am not sure on this point but I can not get f2py compiled when trying to output the Fortran pointer directly.

This is the final wrapper subroutine in Fortran for f2py:

    subroutine get_array(x, z)

            use iso_c_binding

            real(8), intent(in)        :: x
            real(8), intent(out)       :: z(100)
            
            ! Auxiliary variables to convert C pointer to Fortran Array
            type(c_ptr)                :: ret_c_ptr
            real(8), pointer           :: f_ptr(:)
    
            ! Interface to C function
            interface                        
                    type(c_ptr) function c_get_array(a) bind(C, name="get_array")
                            import
                            real(c_double), value :: a
                    end function
            end interface

            ! Call C function               
            ret_c_ptr = c_get_array(x)                
            call c_f_pointer(ret_c_ptr, f_ptr, [100])
            
            ! Assign Fortran pointer to Fortran output array       
            z = f_ptr

    end subroutine
M.E.
  • 4,955
  • 4
  • 49
  • 128
  • 1
    In summary, the Fortran interface has function result of `type(c_ptr)` which then gives you a Fortran array pointer with `c_f_pointer` using shape `[100]`. – francescalus Jul 12 '21 at 16:35
  • Can I use the returned pointer as a regular Fortran array (like I do in C)? The intention is that the output variable from the wrapping Fortran subroutine returns a 100 element Fortran array so `f2py` can then encapsulate the array into a numpy ndarray. If so, how do I convert a Fortran pointer to a 100 elements `real(8)/double` into a 100 elements `real(8)` Fortran array? – M.E. Jul 12 '21 at 16:36
  • 1
    `call c_f_pointer(c_ptr, f_array, [100]) ! f_array becomes a (pointer) array with 100 elements` – francescalus Jul 12 '21 at 16:44
  • Ok but if I do: `real(8), intent(out) :: z(100) ... call c_f_pointer(c_ptr, z, [100]) ....` then I get after compiling with f2py the following error `48 | call c_f_pointer(ptr, z, [100]) Error: Argument FPTR at (1) to C_F_POINTER must be a pointer`. Is it possible to get the data into an actual array? If I return that via f2py will I get a numpy vector? – M.E. Jul 13 '21 at 06:09
  • 1
    You'll need to declare `z` as a pointer: `real(c_double), pointer :: z(:)`. That is an "actual array". – francescalus Jul 13 '21 at 06:14
  • Thanks, that helps a lot, I managed to get it compiled and it works. However, I get an error while trying to output directly the Fortran pointer. I am not 100% sure about the cause but I have used an assignment between the Fortran pointer and the array (seems to work). I guess this is due to the limitations of `f2py` not being able to handle Fortran pointers directly into numpy arrays. See the edit for details in case you are curious about it. – M.E. Jul 13 '21 at 06:27

0 Answers0