How to Access Allocatable Array Data from a Foreign Language Library

All,

I’m very new to Python, and C/C++ for that matter… I’ve been a Fortran programmer since the dawn of time, in an environment that didn’t really require expertise in interoperability with other languages. Times have changed!

I’m trying to build a Python script that accesses array data from a foreign language library (a Fortran DLL) that dynamically allocates memory for the array of interest. I gather from various online searches that this is not trivial to do, if even possible.

The Fortran function discussed here is a simplified example of the one I actually need to integrate with a Python script. But if I can get this example to work, I can make the real one work. The constraint is that the Fortran library itself is immutable… either it works as-is, or it doesn’t.

Although this is a Python forum, I’ll include the Fortran code here for completeness. It’s a simple function that fills an array whose size is initially undertermined:

      Integer Function SimpleDLL(OutInt,OutArray) Bind(C, name="SimpleDLL")
      USE, INTRINSIC :: ISO_C_BINDING
      IMPLICIT NONE
      !DIR$ ATTRIBUTES DLLEXPORT :: SimpleDLL
      Integer(C_INT), INTENT(OUT) :: OutInt
      Real(C_FLOAT), Allocatable, INTENT(OUT) :: OutArray(:)
      Real :: r,sum
      Integer :: i,j
      
      i = 0
      sum = 0
      Call RANDOM_SEED()
      Do While (sum <= 5)
        Call RANDOM_NUMBER(r)
        sum = sum + r
        i = i + 1
      EndDo
      Allocate(OutArray(i))
      Do j = 1,i
        OutArray(j) = FLOAT(j)**(i/5.)
      EndDo
      OutInt = i
      SimpleDLL = 5/i
      
      End Function SimpleDLL

The above is compiled into a DLL that is interoperable with C via the BIND(C) directive (and hence, in principle, is interoperable with Python). I have verified that this DLL can be succesfully called from a Fortran driver that declares the function argument corresponding to OutArray to be an allocatable C_FLOAT array.

This is my attempt to call this DLL function from a Python script (I’m also not at all fluent with ctypes!):

# Python script that calls SimpleDLL
import ctypes as ct
import numpy as np

# import the dll
flib = ct.CDLL('./SimpleDLL.dll')

# declare the DLL function's result type
flib.SimpleDLL.restype = ct.c_int

# define argument types
c_int_p = ct.POINTER(ct.c_int)
c_float_p = ct.POINTER(ct.c_float)

# initialize arguments
InInt = ct.c_int(0)
InArray = ct.c_float()

# call the DLL function
result = flib.SimpleDLL(ct.byref(InInt), ct.byref(InArray))

# show the results
print('Result = ',result)
print('Array Size = ' + str(InInt.value))
print('Array = ' + str(InArray.value))

The output for a typical execution of this script is:

Result = 0
Array Size = 11
Array = 5.925913876308824e-33

The script is evidently communicating with the DLL (I successfully get back the function result and the array size), but what I really need is the data in the returned array, which I clearly don’t know how to extract. Is it even possible from Python? Thanks for any info-

I’ve never done any FORTRAN/C interop myself, but I’m guessing that this translates to

int SimpleDLL(int *OutInt, float **OutFloat);

in C - float ** being a pointer to pointer to float, i.e., you’re passing the function a location where you would like it to store the location of the array it allocates.

As to how to express that in ctypes, I don’t have time to look it up right now, but maybe this helps as a hint (assuming my hunch is correct).

(Note there might not be any easy way to free the allocated array from Python)

I gave it a proper go, first from C and then from Python once I’d figured out what was going on.

I modified your test function to print some debug output so that I know what I should expect

Integer(C_INT) Function SimpleDLL(OutInt,OutArray) Bind(C, name="SimpleDLL")
    USE, INTRINSIC :: ISO_C_BINDING
    IMPLICIT NONE
    Integer(C_INT), INTENT(OUT) :: OutInt
    Real(C_FLOAT), Allocatable, INTENT(OUT) :: OutArray(:)
    Real :: r,sum
    Integer :: i,j

    i = 1
    sum = 0
    Call RANDOM_SEED()

    Do While (sum <= 5)
        Call RANDOM_NUMBER(r)
        sum = sum + r
        i = i + 1
    EndDo


    Allocate(OutArray(i))

    Do j = 1,i
        OutArray(j) = FLOAT(j)**(i/5.)
        print *, "value", j, "=", OutArray(j)
    EndDo
    
    OutInt = i 
    SimpleDLL = 5/i 
    print *, "OutInt=", i, "SimpleDLL", SimpleDLL

End Function SimpleDLL

At first the results I was getting calling this function made no sense at all. Not even the integer was being passed out correctly. Not until I commented out the allocation part of the Fortran code.

That’s how I discovered that a Fortran allocatable array uses more than one pointer’s worth of space (on my Linux machine, using GFortran, it appears to use at least 6, or 48 bytes rather than 8) — at first, Fortran was overwriting memory that I hadn’t intended for it to overwrite.

Anyway, this is the C code with which I got it working:

extern int SimpleDLL(int *OutInt, float **OutFloat);

#include <stdio.h>
#include <string.h>

#define NOTHING_SIZE 8

int main()
{
    printf("Calling SimpleDLL\n");
    int i = 0;
    struct {
        long nothing1[NOTHING_SIZE];
        float *f;
        long nothing2[NOTHING_SIZE];
    } s;
    memset(&s, 0, sizeof(s));

    printf("Have %d and %p\n", i, s.f);

    int rv = SimpleDLL(&i, &s.f);

    printf("Returned %d, passed out %d and %p\n", rv, i, s.f);

    for (int j = 0; j < i; ++j) {
        printf(" array[%d] = %f\n", j, (double)s.f[j]);
    }

    printf("The nothing around the float is:\n");
    for (int j = 0; j < NOTHING_SIZE; ++j) {
        printf("before[%d] = %#lx\n", j, s.nothing1[j]);
    }
    for (int j = 0; j < NOTHING_SIZE; ++j) {
        printf("after[%d] = %#lx\n", j, s.nothing2[j]);
    }
    printf("Address of s: %p, and of f: %p\n", &s, s.f);
    printf("all good\n");
    return 0;
}

As you can see, in the end I added some extra space around the pointer I’m passing to Fortran. Really, the C signature of the function should be something like

struct fortran_allocatable_f32 {
    float *address;
   /* unknown other fields */
};
extern int SimpleDLL(int *OutInt, struct fortran_allocatable_f32 *OutFloat);

The C program above prints (e.g.)

Calling SimpleDLL
Have 0 and (nil)
 value           1 =   1.00000000    
 value           2 =   4.59479380    
 value           3 =   11.2115793    
 value           4 =   21.1121273    
 value           5 =   34.4932442    
 value           6 =   51.5148926    
 value           7 =   72.3128891    
 value           8 =   97.0058670    
 value           9 =   125.699501    
 value          10 =   158.489334    
 value          11 =   195.462723    
 OutInt=          11 SimpleDLL           0
Returned 0, passed out 11 and 0x1c14d20
 array[0] = 1.000000
 array[1] = 4.594794
 array[2] = 11.211579
 array[3] = 21.112127
 array[4] = 34.493244
 array[5] = 51.514893
 array[6] = 72.312889
 array[7] = 97.005867
 array[8] = 125.699501
 array[9] = 158.489334
 array[10] = 195.462723
The nothing around the float is:
before[0] = 0
before[1] = 0
before[2] = 0
before[3] = 0
before[4] = 0
before[5] = 0
before[6] = 0
before[7] = 0
after[0] = 0x4
after[1] = 0x403000100000000
after[2] = 0x1
after[3] = 0xb
after[4] = 0x4
after[5] = 0
after[6] = 0
after[7] = 0
Address of s: 0x7fff6b0ad9c0, and of f: 0x1c14d20
all good

I don’t know what the extra information Fortran writes out is, and, crucially, I don’t know how much space it will take on other compilers and architectures, but some of it (after[3]) appears to be the size of the array.

Armed with this knowledge, it’s not too hard to get it to work with ctypes:

import ctypes as ct

flib = ct.CDLL('./alloctest.so')

flib.SimpleDLL.restype = ct.c_int

# Define the arguments
length = ct.c_int(0)
# Allocate some extra memory for FORTRAN to write to
floatp_16arr_t = ct.POINTER(ct.c_float) * 16
outarray_mem = floatp_16arr_t()

result = flib.SimpleDLL(ct.byref(length), ct.byref(outarray_mem[0]))

print("result =", result)
print("length =", length)
# Get the data out of the C pointer
outarray_ptr = outarray_mem[0]
res_list = [ outarray_ptr[i] for i in range(length.value) ]
print("result =", res_list)

(Note that you will have to go back into Fortran to free the memory allocated)

tjol,

Thanks so much for your feedback! My ignorance of C is such that I have only the most feeble grasp of what you’ve done here, but clearly you’ve identified a workable approach. I was trying various ways of using C pointers without success (hence the anomalous do-nothing definitions of c_int_p and c_float_p in my original posted Python script).

I’ll need a while to digest your methodology, and see if I can reproduce it on my end. I gather from your closing comment that the DLL itself must deallocate the memory it allocated for OutArray. That would require modification of the DLL (and technically violate
my immutable DLL constraint), but it would be fairly trivial… e.g., add an optional input argument that, if present, would simply trigger deallocation of OutArray and return.

Thanks again for your investigation and insightful replies. I’ll post again with results of my next attempt, successful or not.

That would be ideal, and may be necessary. (The traditional way to do this in C APIs is to have a separate subroutine to deallocate/free/clean up pointers and arrays returned by the library).

The issue is that the memory is allocated using the Fortran allocator used by the DLL, so it must be deallocated by whatever allocator that is. It may be possible to use a separate Fortran DLL with a simple deallocation subroutine, but I think at least on Windows it’s possible separate DLLs get separate pools of memory, so that might not work.

tjol,

My apologies for the extended delay, too many other fires to put out!

First let me report the very good news that I was able to reproduce your results using the python script you provided. It works flawlessly! Moreover, I rewrote the SimpleDLL function to create a more complex multidimensional array, which was also accessible from the same python script, as is. The resulting list in python (res_list) even retains the Fortran convention of the first index being the most rapidly varying, which makes sense since it is simply listing the contents of addresses that are contiguous in memory.

Most of my subsequent effort has been to try to reproduce your C code and see about the extra space required around the pointer for the allocatable array on my specific system (Windows 7 64-bit). I was utterly unable to do this. I created a C++ project in Visual Studio 2017, pasted in your code, and linked (or tried linking to) my SimpleDLL fortran library. No matter what I do, it insists that the reference in the C code to SimpleDLL is an unresolved external. I tried manually editing Additional Dependencies, added the SimpleDLL.lib file as a resource, etc, etc, no luck. I finally gave up.

In any case, I’m entirely delighted with how well this Fortran allocatable array extraction worked. I can assure you that, after a pretty extensive review of search results of other software development forums, there is widespread belief that what’s been demonstrated here is not possible in Python. Thanks again for your help!