Optimization Without Using Derivatives: the PRIMA Package, its Fortran Implementation, and Its Inclusion in SciPy

With pleasure!

If you want to start getting your hands dirty you could start with a very simple example. I reworked a minimum working example previously discussed here

Since I work on windows but also have linux through WSL I prepared it such that you can test the following combination A) windows with ifort B) Linux with gfortran. Also, let’s assume that for simplicity you have everything in the same folder:

fortran2python.f90
module fortran2python
  use iso_c_binding
  implicit none
  contains
    subroutine say_hello(X,dim,num_points) bind(c)
#if __INTEL_COMPILER
      !DEC$ ATTRIBUTES DLLEXPORT :: say_hello
#elif __GNUC__
      !GCC$ ATTRIBUTES DLLEXPORT :: say_hello
#endif      
    ! Variables
      integer(c_int), intent(in) :: num_points
      integer(c_int), intent(in) :: dim 
      real(c_double), intent(inout) :: X(dim,num_points) !! array of coordinates
      
      integer :: i
      !-------------------------------
      print *, 'hello from fortran'
      
      print *, 'Whole array: ',X(:,:)
      
      do i = 1, num_points
          print *, 'by data point', x(:,i)
      end do
      
      x(1,2) = 2*x(1,2)
      
  end subroutine
  
  end module

test.py
import os
import numpy as np
import ctypes

folderpath = '' 
if os.name == 'nt': #windows
    DLLname = 'fortran2python.dll'
elif os.name == 'posix': #Linux/Mac
    DLLname = 'libfortran2python.so'
api = np.ctypeslib.load_library(DLLname,folderpath)

X = np.array([
    [1.0, 2.0],
    [3.0, 4.0],
    [5.0, 6.0],
])

print('X shape:',X.shape)
print(X[0,0],X[0,1])
print(X[1,0],X[1,1])
print(X[2,0],X[2,1])

api.say_hello(  ctypes.byref(np.ctypeslib.as_ctypes( X )) ,
                ctypes.byref(ctypes.c_int(X.shape[1])) ,
                ctypes.byref(ctypes.c_int(X.shape[0])) 
             )

print('back in python')
print(X[0,0],X[0,1])
print(X[1,0],X[1,1])
print(X[2,0],X[2,1])

print('hello')

With those two files you can now compile the dynamic link library and run the test from python as follows

A) Windows + ifort

ifort -fpp -O3 -c fortran2python.f90
ifort -dll -exe:fortran2python.dll fortran2python.obj
python test.py

B) Linux + gfortran

gfortran -cpp -O3 -fpic -c fortran2python.f90
gfortran -shared -o libfortran2python.so fortran2python.o
python test.py

You can learn more here: Managing libraries (static and dynamic libraries) — Fortran Programming Language

few details to notice:

  • Had to use the -fpp/-cpp pre-processors because of these key lines:
    subroutine say_hello(X,dim,num_points) bind(c)
#if __INTEL_COMPILER
      !DEC$ ATTRIBUTES DLLEXPORT :: say_hello
#elif __GNUC__
      !GCC$ ATTRIBUTES DLLEXPORT :: say_hello
#endif      

intel and gnu compilers use different decorators, so in order to be cross-compiler compatible you will have to use this in your interface functions.

As mentioned in the page I linked before, the ‘-fPIC’ flag is necessary with GNU compiler in order to enable Position Independent Code, which is needed make the collection of functions in the library callable in any order.

Once you get comfortable with this interfacing with a very simple example you’ll have to put in the mix the “build system” part, you can go for CMake or meson (I’m most familiar with the first but the latter seems to be even easier)… As was said previously the SciPy team can help with that decision. Anyhow, at the end your work will be to provide the binary compiled library that would be as least intrusive as possible in their architecture. Since calling N dynamic (or static) libraries is basically how most applications work now a days it wont (shouldn’t) be a problem.

2 Likes