我想知道将大型 fortran 数组（5000 x 5000 实单精度数）写入文件的最佳方法。我正在尝试保存数值计算的结果以供以后使用，因此它们不需要重复。根据计算，每个数字 5000 x 5000 x 4bytes 是 100 Mb，是否可以将其保存为只有 100Mb 的形式？有没有办法将 fortran 数组保存为二进制文件并将其读回以供以后使用？
我注意到将数字保存到文本文件会产生比所保存数据类型的大小大得多的文件。这是因为数字被保存为字符吗？
我熟悉的唯一写入文件的方法是
open (unit=41, file='outfile.txt')
do i=1,len
do j=1,len
write(41,*) Array(i,j)
end do
end do
虽然我想有更好的方法来做到这一点。如果有人可以向我指出一些资源或示例来证明我能够有效地（就内存而言）编写和读取更大的文件，那就太好了。 谢谢！
【问题讨论】：
以二进制形式编写数据文件，除非您实际上要读取输出 - 而且您不会读取 250 万个元素的数组。
使用二进制的原因有三个，重要性在降低：
准确性问题可能是最明显的。当您将（二进制）浮点数转换为十进制数的字符串表示时，您不可避免地会在某些时候截断。如果您确定当您将文本值读回浮点值时，您肯定会得到相同的值，那没关系；但这实际上是一个微妙的问题，需要仔细选择格式。使用默认格式，各种编译器以不同程度的质量执行此任务。 This blog post，从游戏程序员的角度写的，很好地涵盖了这些问题。
让我们考虑一个小程序，它针对各种格式，将单精度实数写入字符串，然后再次将其读回，跟踪它遇到的最大错误。我们只是从 0 到 1，以机器 epsilon 为单位。代码如下：
program testaccuracy
character(len=128) :: teststring
integer, parameter :: nformats=4
character(len=20), parameter :: formats(nformats) = &
[ '( E11.4)', '( E13.6)', '( E15.8)', '(E17.10)' ]
real, dimension(nformats) :: errors
real :: output, back
real, parameter :: delta=epsilon(output)
integer :: i
errors = 0
output = 0
do while (output < 1)
do i=1,nformats
write(teststring,FMT=formats(i)) output
read(teststring,*) back
if (abs(back-output) > errors(i)) errors(i) = abs(back-output)
enddo
output = output + delta
end do
print *, 'Maximum errors: '
print *, formats
print *, errors
print *, 'Trying with default format: '
errors = 0
output = 0
do while (output < 1)
write(teststring,*) output
read(teststring,*) back
if (abs(back-output) > errors(1)) errors(1) = abs(back-output)
output = output + delta
end do
print *, 'Error = ', errors(1)
end program testaccuracy
当我们运行它时，我们得到：
$ ./accuracy
Maximum errors:
( E11.4) ( E13.6) ( E15.8) (E17.10)
5.00082970E-05 5.06639481E-07 7.45058060E-09 0.0000000
Trying with default format:
Error = 7.45058060E-09
请注意，即使使用小数点后 8 位的格式（考虑到single precision reals are only accurate to 6-7 decimal places，我们可能认为这已经足够了）我们也不会得到精确的副本，相差大约 1e-8。而且这个编译器的默认格式没有给我们准确的往返浮点值；引入了一些错误！如果您是视频游戏程序员，那么这种准确度可能就足够了。但是，如果您正在对湍流进行随时间变化的模拟，那可能绝对不行，特别是如果在引入误差的位置存在一些偏差，或者如果误差发生在应该是守恒的量中。
请注意，如果您尝试运行此代码，您会发现它需要很长时间才能完成。这是因为，也许令人惊讶的是，性能是浮点数文本输出的另一个真正问题。考虑以下简单程序，它只是将 5000 × 5000 实数数组的示例写成文本和未格式化的二进制文件：
program testarray
implicit none
integer, parameter :: asize=5000
real, dimension(asize,asize) :: array
integer :: i, j
integer :: time, u
forall (i=1:asize, j=1:asize) array(i,j)=i*asize+j
call tick(time)
open(newunit=u,file='test.txt')
do i=1,asize
write(u,*) (array(i,j), j=1,asize)
enddo
close(u)
print *, 'ASCII: time = ', tock(time)
call tick(time)
open(newunit=u,file='test.dat',form='unformatted')
write(u) array
close(u)
print *, 'Binary: time = ', tock(time)
contains
subroutine tick(t)
integer, intent(OUT) :: t
call system_clock(t)
end subroutine tick
! returns time in seconds from now to time described by t
real function tock(t)
integer, intent(in) :: t
integer :: now, clock_rate
call system_clock(now,clock_rate)
tock = real(now - t)/real(clock_rate)
end function tock
end program testarray
这里是时序输出，用于写入磁盘或 ramdisk：
Disk:
ASCII: time = 41.193001
Binary: time = 0.11700000
Ramdisk
ASCII: time = 40.789001
Binary: time = 5.70000000E-02
请注意，当写入磁盘时，二进制输出的速度是 ASCII 的 352 倍，而对于 ramdisk，它接近 700 倍。这有两个原因 - 一是您可以一次写出所有数据，而不必循环；另一个是生成浮点数的字符串十进制表示是一个非常微妙的操作，需要对每个值进行大量计算。
最后是数据大小；上面示例中的文本文件（在我的系统上）大约是二进制文件大小的 4 倍。
现在，二进制输出确实存在问题。特别是，原始 Fortran（或者，就此而言，C）二进制输出非常脆弱。如果您更改平台，或者您的数据大小发生变化，您的输出可能不再有任何好处。向输出中添加新变量会破坏文件格式，除非您总是在文件末尾添加新数据，并且您无法提前知道从您的协作者（谁可能是你，三个月前）。使用像NetCDF 这样的库可以避免二进制输出的大部分缺点，这些库编写的自描述二进制文件比原始二进制文件更“面向未来”。更好的是，由于它是一个标准，许多工具都可以读取 NetCDF 文件。
网上有很多NetCDF教程；我们的是here。一个使用 NetCDF 的简单示例给出了与原始二进制文件相似的时间：
$ ./array
ASCII: time = 40.676998
Binary: time = 4.30000015E-02
NetCDF: time = 0.16000000
但给你一个很好的自我描述文件：
$ ncdump -h test.nc
netcdf test {
dimensions:
X = 5000 ;
Y = 5000 ;
variables:
float Array(Y, X) ;
Array:units = "ergs" ;
}
文件大小与原始二进制文件大致相同：
$ du -sh test.*
96M test.dat
96M test.nc
382M test.txt
代码如下：
program testarray
implicit none
integer, parameter :: asize=5000
real, dimension(asize,asize) :: array
integer :: i, j
integer :: time, u
forall (i=1:asize, j=1:asize) array(i,j)=i*asize+j
call tick(time)
open(newunit=u,file='test.txt')
do i=1,asize
write(u,*) (array(i,j), j=1,asize)
enddo
close(u)
print *, 'ASCII: time = ', tock(time)
call tick(time)
open(newunit=u,file='test.dat',form='unformatted')
write(u) array
close(u)
print *, 'Binary: time = ', tock(time)
call tick(time)
call writenetcdffile(array)
print *, 'NetCDF: time = ', tock(time)
contains
subroutine tick(t)
integer, intent(OUT) :: t
call system_clock(t)
end subroutine tick
! returns time in seconds from now to time described by t
real function tock(t)
integer, intent(in) :: t
integer :: now, clock_rate
call system_clock(now,clock_rate)
tock = real(now - t)/real(clock_rate)
end function tock
subroutine writenetcdffile(array)
use netcdf
implicit none
real, intent(IN), dimension(:,:) :: array
integer :: file_id, xdim_id, ydim_id
integer :: array_id
integer, dimension(2) :: arrdims
character(len=*), parameter :: arrunit = 'ergs'
integer :: i, j
integer :: ierr
i = size(array,1)
j = size(array,2)
! create the file
ierr = nf90_create(path='test.nc', cmode=NF90_CLOBBER, ncid=file_id)
! define the dimensions
ierr = nf90_def_dim(file_id, 'X', i, xdim_id)
ierr = nf90_def_dim(file_id, 'Y', j, ydim_id)
! now that the dimensions are defined, we can define variables on them,...
arrdims = (/ xdim_id, ydim_id /)
ierr = nf90_def_var(file_id, 'Array', NF90_REAL, arrdims, array_id)
! ...and assign units to them as an attribute
ierr = nf90_put_att(file_id, array_id, "units", arrunit)
! done defining
ierr = nf90_enddef(file_id)
! Write out the values
ierr = nf90_put_var(file_id, array_id, array)
! close; done
ierr = nf90_close(file_id)
return
end subroutine writenetcdffile
end program testarray
【讨论】：
以“未格式化”的形式打开文件进行读写，并在不提供格式的情况下读写数据，如下面的程序所示。
program xunformatted
integer, parameter :: n = 5000, inu = 20, outu = 21
real :: x(n,n)
integer :: i
character (len=*), parameter :: out_file = "temp_num"
call random_seed()
call random_number(x)
open (unit=outu,form="unformatted",file=out_file,action="write")
do i=1,n
write (outu) x(i,:) ! write one row at a time
end do
print*,"sum(x) =",sum(x)
close (outu)
open (unit=inu,form="unformatted",file=out_file,action="read")
x = 0.0
do i=1,n
read (inu) x(i,:) ! read one row at a time
end do
print*,"sum(x) =",sum(x)
end program xunformatted
【讨论】：