I currently have a Fortran function that I want to optimize using SciPy wrapping it using Ctypes. Is this possible? Perhaps I've done something wrong in my implementation. For example, assume I have:

cost.f90

module cost_fn
  use iso_c_binding, only: c_float
  implicit none

contains

  function sin_2_cos(x,y) bind(c)
      real(c_float) :: x, y, sin_2_cos
      sin_2_cos = sin(x)**2 * cos(y)
  end function sin_2_cos

end module cost_fn

that I compile with:

gfortran -fPIC -shared -g -o cost.so cost.f90

and then try to find a (local) minimum with:

cost.py

#!/usr/bin/env python

from ctypes import *
import numpy as np
import scipy.optimize as sopt

cost = cdll.LoadLibrary('./cost.so')

cost.sin_2_cos.argtypes = [POINTER(c_float), POINTER(c_float)]
cost.sin_2_cos.restype = c_float

def f2(x):
    return cost.sin_2_cos(c_float(x[0]), c_float(x[1]))
    # return np.sin(x[0])**2 * np.cos(x[1])

# print(f2( [1, 1] ))
# print(f2( [0.5 * np.pi, np.pi] ))

print( sopt.minimize( f2, (1.0, 1.0), options={'disp': True}, tol=1e-8) )

I expect a local minimum f2(pi / 2, pi) = -1. When I call f2 with the cost.sin_2_cos return value, the "minimimum" is just given at the initial guess of (1, 1). If I call f2 with the "Python" return value, optimize finds the correct minimum.

I've tried redefining sin_2_cos to take dimension(2) array input, but was seeing similar behavior. Perhaps I need to call sin_2_cos directly with minimize (but then how would I specify c_float for the arguments)? Any thoughts are appreciated!

Edit: To a comment below, note that the two commented print(f2(...)) lines produce the expected values. Thus, I believe the Fortran function is being properly called through the Python f2 function.