Exchanging multi-dimensional array between Python and Fortran

Exchange multi-dimensional array data in binary format

Maintaining a very old program which was written in Fortran (F77), we want a part of the calculation to be finished in Python.

One should take care of:

  • Array size
  • Array indexing order
  • Data type

Ordering

When using scipy.io.FortranFile, you will always need to transpose the matrix.

When using numpy.ndarray.tofile, you will not need to transpose the matrix if you specify “order="F"”, but you will need transposition otherwise.

Running the Python program

When use the Intel ifort compiler, we use SYSTEM() to call the system command, and which is part of ifport.

       USE ifport
       SYSTEM("python3 xxx.py")

Example

A working prototype using scipy.io.FortranFile and numpy.ndarray.tofile is here.

The Fortran part

       PROGRAM BAKER

       USE ifport

       PARAMETER (mxnx = 3, mxny = 4, mxnz = 7)
       DATA nx, ny, nz / 2, 3, 7 /

       REAL*8 mat_from_fortran
       REAL*8 mat_scipy
       REAL*8 mat_tofile
       DIMENSION mat_from_fortran(mxnx, mxny, mxnz)
       DIMENSION mat_scipy(mxnx, mxny, mxnz)
       DIMENSION mat_to_file(mxnx, mxny, mxnz)

       OPEN(1001, FILE="params.dat", ACTION="write")
       WRITE(1001, '(3I9)') mxnx, mxny, mxnz
       WRITE(1001, '(3I9)') nx, ny, nz
       CLOSE(1001)

       DATA res / 0 /
       res = SYSTEM("python3 baker.py")

       OPEN(1002, FILE="out.bin", FORM="unformatted",
     ,            ACTION="write")
       WRITE(*,*) "Opening..."
       WRITE(1002) mat_from_fortran
       WRITE(*,*) "Written finished..."
       CLOSE(1002)

       WRITE(*,*) "Finished writing"

       OPEN(1002, FILE="scipy.bin", FORM="unformatted", 
     ,            ACTION="read")
       WRITE(*,*) "Opening..."
       READ(1002) mat_scipy
       WRITE(*,*) "Read finished..."
       CLOSE(1002)

       do i = 1, mxnx
         do j = 1, mxny
           do k = 1, mxnz
c            WRITE(*, '(3I2, F12.4)') i, j, k, mat(i,j,k)
           enddo
         enddo
       enddo

       WRITE(*, *) mat_scipy(1,1,6)
       WRITE(*, *) mat_scipy(1,3,5)
       WRITE(*, *) mat_scipy(3,2,5)
       WRITE(*, *) mat_scipy(3,4,7)

       OPEN(1002, FILE="tofile.bin", FORM="unformatted",
     ,            ACCESS="direct", ACTION="read", recl=8)
       WRITE(*,*) "Opening..."
       READ(1002, rec=1) mat_tofile
       WRITE(*,*) "Read finished..."
       CLOSE(1002)
       
       WRITE(*, *) mat_scipy(1,1,6)
       WRITE(*, *) mat_scipy(1,3,5)
       WRITE(*, *) mat_scipy(3,2,5)
       WRITE(*, *) mat_scipy(3,4,7)

       END

The Python part

import numpy, scipy.io

with open("params.dat") as fp:
    mxnx, mxny, mxnz = (int(i) for i in fp.readline().strip().split())
    nx, ny, nz = (int(i) for i in fp.readline().strip().split())

print('Python :D')
print(mxnx, mxny, mxnz)
print(nx, ny, nz)
    
mat = numpy.zeros((mxnx, mxny, mxnz), order="F", dtype=numpy.float64)
mat = numpy.random.rand(mxnx, mxny, mxnz).astype(numpy.double, order="F")

print(mat[0,0,5])
print(mat[0,2,4])
print(mat[2,1,4])
print(mat[2,3,6])
print(mat.shape)

fp = scipy.io.FortranFile("scipy.bin", "w")
fp.write_record(mat.T)
fp.close()

mat.tofile("tofile.bin")

Get the fortran matrix with correct shape

import scipy.io

def read_fortran_unformatted_array(fname, max_shape: tuple, real_shape: tuple):
    fp = scipy.io.FortranFile(fname, "r")
    mat = fp.read_record(numpy.double).reshape(max_shape, order="F")
    fp.close()
    return mat[tuple(numpy.indices(real_shape))]

Resources