1

I am trying to wrap some fortran code into python using the ctypes library but have am having major issues with the data transfer. When I print the data from my python script, it looks substantially different from when I print it within the fortran code. Can someone please help me figure out what is going on here? I tried playing around with the datatypes which did not fix the solution and so far I have not found any other SO questions that have addressed my issue. Also please note that using f2py will not work for the end product that this example refers too.

Below is an example of my code:

temp.f90

       subroutine test(num_mod_sel,goodval,x,y,coeff,coeff_flag)
       implicit none


! Input/Output
       integer num_mod_sel, goodval
       real, dimension(goodval,num_mod_sel) :: x
       real, dimension(goodval)             :: y
       real, dimension(num_mod_sel)         :: coeff
       integer, dimension(num_mod_sel)      :: coeff_flag


      print*, num_mod_sel,goodval,x


      return
      end subroutine test

!===================================================================================================!

The above f90 code is compiled with:

gfortran -shared -fPIC -o temp.so temp.f90

test.py

from ctypes import CDLL,cdll, POINTER, c_int, c_float
import numpy as np


def test_fcode(num_mod_sel,goodval,x,y,coeff,coeff_flag):

    fortran = CDLL('./temp.so')

    fortran.test_.argtypes = [ POINTER(c_int),
                                     POINTER(c_int),
                                     POINTER(c_float),
                                     POINTER(c_float),
                                     POINTER(c_float),
                                     POINTER(c_int) ]

    fortran.test_.restype = None

    num_mod_sel_ = c_int(num_mod_sel)
    goodval_     = c_int(goodval)
    x_           = x.ctypes.data_as(POINTER(c_float))
    y_           = y.ctypes.data_as(POINTER(c_float))
    coeff_       = coeff.ctypes.data_as(POINTER(c_float))
    coeff_flag_  = coeff_flag.ctypes.data_as(POINTER(c_int))

    fortran.test_(num_mod_sel_,goodval_,x_,y_,coeff_,coeff_flag_)


    
 #Create some test data   
    
num_mod_sel = 4
goodval = 10
x = np.full((num_mod_sel,goodval),999.,dtype=float)
x[:] = np.random.rand(num_mod_sel,goodval)
y = np.full(goodval,999.,dtype=float)
y[:] = np.random.rand(goodval)

coeff = np.empty(num_mod_sel,dtype=float)
coeff_flag = np.empty(num_mod_sel,dtype=int)


#Run the fortran code
test_fcode(num_mod_sel,goodval,x,y,coeff,coeff_flag)

print(x) from the python code:

[[0.36677304 0.8734628 0.72076823 0.20234787 0.91754331 0.26591916 0.46325577 0.00334941 0.98890871 0.3284262 ] [0.15428096 0.24979671 0.97374747 0.83996786 0.59849493 0.55188578 0.9668523 0.98441142 0.50954678 0.22003844] [0.54362548 0.42636074 0.65118397 0.69455346 0.30531619 0.88668116 0.97278714 0.29046492 0.64851937 0.64885967] [0.31798739 0.37279389 0.88855305 0.38754276 0.94985151 0.56566525 0.99488508 0.13812829 0.0940132 0.07921261]]

print*, x from f90:

-2.91465824E-17 1.68338645 13.0443134 1.84336567 -7.44153724E-34 1.80519199 -2.87629426E+27 1.57734776 -1297264.38 1.85438573 -236487.531 1.63295949 -1.66118658E-33 1.73162782 -6.73423983E-09 0.919681191 -1.09687280E+21 1.87222707 5.50313165E+09 1.66421306 8.38275158E+34 1.52928090 -2.15154066E-13 1.62479663 3.88800366E+30 1.86843681 127759.977 1.83499193 -3.55062879E+15 1.77462363 2.43241945E+19 1.76297140 3.16150975E-03 1.86671305 1.35183692E+21 1.87110281 1.74403865E-31 1.75238669 9.85857248E-02 1.59503841 -2.33541620E+30 1.79045486 -1.86185171E+11 1.78229403 4.23132255E-20 1.81525886 2.96771497E-04 1.82888138 -4.55096013E-26 1.86097753 0.00000000 3.68934881E+19 -7.37626273E+15 1.58494916E+29 0 -1064355840 -646470284 -536868869

2
  • 2
    For "dtype=float" in numpy the Python "float" is interpreted as "numpy.float_" which is usually same as "numpy.float64" (C double) while Fortran's "real" is more likely "numpy. float32". Commented Sep 28, 2021 at 21:11
  • @MichaelButscher is correct: in the Python code, the dtype for x, y and coeff must be np.float32. Michael, make your comment an answer! Commented Sep 28, 2021 at 23:24

1 Answer 1

2

The problem is a mismatch of datatypes.

The Fortran real is usually a 32 bit float (C float) while numpy interprets the Python datatype float as numpy.float_ which is an alias of numpy.float64, the C double with 64 bits.

Solution: In Python use numpy.float32 as dtype for numpy array creation.

Sign up to request clarification or add additional context in comments.

Comments

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.