Regarding the input / output of python <-> fortran, it is difficult to find a cohesive description on the Internet, so I decided to summarize it as far as I know.
As an aside, np.savez and np.load are recommended for data input / output between pythons.
Output data according to the following procedure
reshapetofileAn example of python code that outputs a 3D double precision real number with little endian
import numpy as np
ix, jx, kx = 3, 4, 5
endian='<'
a = np.ones((ix,jx,kx))
a.reshape([ix*jx*kx],order='F').astype(endian+'d').tofile('test.dat')
To read this with fortran, do as follows. Suppose the calculator you are using is little endian. ʻAccess ='stream' is required for ʻopen in fortran
program test
implicit none
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.d0)), dimension(ix,jx,kx) :: a
open(newunit=idf, file='test.dat',form='unformatted',access='stream')
read(idf) a
close(idf)
stop
end program test
If you want to exchange big endian data, use python code
endian='>'
You can change it to. This ʻendian variable controls how to convert to binary with ʻas type
In fortran code
open(newunit=idf, file='test.dat',form='unformatted',access='stream',convert='big_endian')
And. If you want to output / input in single precision, rewrite the python code as follows
import numpy as np
ix, jx, kx = 3, 4, 5
endian='<'
a = np.ones((ix,jx,kx),dtype=np.float32)
a.reshape([ix*jx*kx],order='F').astype(endian+'f').tofile('test.dat')
Changed the parts of np.ones and ʻas type` to handle single precision.
The corresponding fortran code looks like this:
program test
implicit none
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.e0)), dimension(ix,jx,kx) :: a
open(newunit=idf, file='test.dat',form='unformatted',access='stream')
read(idf) a
close(idf)
stop
end program test
Only the part of real (KIND (0.e0)) was changed.
If you want to output multiple arrays to one file, follow the steps below.
np.arrayreshape to make a one-dimensional array.T.import numpy as np
ix, jx, kx = 3, 4, 5
endian='<'
a = np.ones((ix,jx,kx),dtype=np.float32)
b = np.ones((ix,jx,kx),dtype=np.float32)*2
c = np.ones((ix,jx,kx),dtype=np.float32)*3
np.array([a,b,c]).T.reshape([3*ix*jx*kx],order='F').astype(endian+'f').tofile('test.dat')
The fortran code to read is
program test
implicit none
integer :: i,j,k
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.e0)), dimension(ix,jx,kx) :: a,b,c
open(newunit=idf, file='test.dat',form='unformatted',access='stream')
read(idf) a,b,c
close(idf)
stop
end program test
Should be
On the other hand, when outputting arrays of different sizes
reshape to make a one-dimensional arraynp.concatenate to organize the arrays
You just have to follow the procedureThe python code looks like this
import numpy as np
ix, jx, kx = 3, 4, 5
lx, mx = 2, 4
endian='<'
a = np.ones((ix,jx,kx))
b = np.ones((lx,mx))*2
a1 = a.reshape([ix*jx*kx],order='F')
b1 = b.reshape([lx*mx],order='F')
np.concatenate([a1,b1]).astype(endian+'d').tofile('test.dat')
When reading from fortran, it is the same as before, but it is as follows.
program test
implicit none
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
integer, parameter :: lx = 2, mx = 4
real(KIND(0.d0)), dimension(ix,jx,kx) :: a
real(KIND(0.d0)), dimension(lx,mx) :: b
open(newunit=idf, file='test.dat',form='unformatted',access='stream')
read(idf) a,b
close(idf)
stop
end program test
The above method is easier for python code, but fortran has to edit the data headers by yourself unless ʻaccess ='stream'is used. It isscipy.io.FortranFile` that adjusts this part well.
Python code looks like this
import numpy as np
from scipy.io import FortranFile
ix, jx, kx = 3, 4, 5
endian='<'
a = np.ones((ix,jx,kx))
f = FortranFile('test.dat','w')
f.write_record(a)
f.close()
The Fortran code for reading is
program test
implicit none
integer :: i,j,k
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.d0)), dimension(ix,jx,kx) :: a
open(newunit=idf, file='test.dat',form='unformatted')
read(idf) a
close(idf)
stop
end program test
It should be done.
Outputting multiple arrays is easier than using Numpy
a = np.ones((ix,jx,kx))
b = np.ones((ix,jx,kx))*2
f = FortranFile('test.dat','w')
f.write_record(a,b)
f.close()
And just add it to the argument of write_record.
Fortran reading code
program test
implicit none
integer :: i,j,k
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.d0)), dimension(ix,jx,kx) :: a,b
open(newunit=idf, file='test.dat',form='unformatted')
read(idf) a,b
close(idf)
stop
end program test
All you have to do is arrange the arrays to be read side by side.
However, if you look at the scipy website
Since this is a non-standard file format, whose contents depend on the compiler and the endianness of the machine, caution is advised. Files from gfortran 4.8.0 and gfortran 4.1.2 on x86_64 are known to work.
Consider using Fortran direct-access files or files from the newer Stream I/O, which can be easily read by numpy.fromfile.
it is written like this. The fortran format seems to be highly compiler-dependent, and the method using numpy seems to be recommended.
In order to read fortran data with numpy, it is recommended to output with ʻaccess ='stream'` (it is not essential because you can adjust the header etc. yourself).
Define an array with fortran and output
program test
implicit none
integer :: i,j,k
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.d0)), dimension(ix,jx,kx) :: a
a = 1.d0
open(newunit=idf, file='test.dat',form='unformatted',access='stream')
write(idf) a
close(idf)
stop
end program test
The simplest way to load it in python is to use np.fromfile.
import numpy as np
ix, jx, kx = 3, 4, 5
endian='<'
f = open('test.dat','rb')
a = np.fromfile(f,endian+'d',ix*jx*kx)
a = a.reshape((ix,jx,kx),order='F')
f.close()
Even if you use np.fromfile, you can use np.dtype for more versatile reading. For example, do as follows
import numpy as np
ix, jx, kx = 3, 4, 5
endian='<'
f = open('test.dat','rb')
dtype = np.dtype([('a',endian+str(ix*jx*kx)+'d')])
a = np.fromfile(f,dtype=dtype)['a']
a = a.reshape((ix,jx,kx),order='F')
f.close()
It is more convenient to use np.dtype when dealing with multiple arrays.
First, fortran outputs multiple arrays to one file as shown below.
program test
implicit none
integer :: i,j,k
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.d0)), dimension(ix,jx,kx) :: a,b
a = 1.d0
b = 2.d0
open(newunit=idf, file='test.dat',form='unformatted',access='stream')
write(idf) a,b
close(idf)
stop
end program test
To load this in python:
import numpy as np
ix, jx, kx = 3, 4, 5
endian='<'
f = open('test.dat','rb')
dtype = np.dtype([('a',endian+str(ix*jx*kx)+'d'),('b',endian+str(ix*jx*kx)+'d')])
data = np.fromfile(f,dtype=dtype)
a = data['a'].reshape((ix,jx,kx),order='F')
b = data['b'].reshape((ix,jx,kx),order='F')
f.close()
It is easy to use scipy.io.FortranFile to read fortran data that does not use ʻaccess ='stream'`.
The output of fortran is as follows
program test
implicit none
integer :: i,j,k
integer :: idf
integer, parameter :: ix = 3, jx = 4, kx = 5
real(KIND(0.d0)), dimension(ix,jx,kx) :: a
a = 1.d0
open(newunit=idf, file='test.dat',form='unformatted')
write(idf) a
close(idf)
stop
end program test
The python code to read this looks like this
import numpy as np
from scipy.io import FortranFile
ix, jx, kx = 3, 4, 5
endian='<'
a = np.ones((ix,jx,kx))
b = np.ones((ix,jx,kx))*2
f = FortranFile('test.dat','r')
a = f.read_record(endian+'('+str(ix)+','+str(jx)+','+str(kx)+')d')
f.close()
It can also be read using np.dtype as follows:
import numpy as np
from scipy.io import FortranFile
ix, jx, kx = 3, 4, 5
endian='<'
a = np.ones((ix,jx,kx))
b = np.ones((ix,jx,kx))*2
f = FortranFile('test.dat','r')
dtype = np.dtype([('a',endian+str(ix*jx*kx)+'d')])
a = f.read_record(dtype=dtype)
a = a['a'].reshape(ix,jx,kx,order='F')
f.close()
Recommended Posts