Lecture de tableaux à partir de fichiers .npy dans Fortran 90

Question

J'utilise Python pour générer des données initiales, sous la forme de tableaux 2D tels que 'X', puis en utilisant Fortran pour faire des calculs avec eux. Initialement, lorsque les tailles de tableau étaient d'environ 10 000 x 10 000, np.savetxt fonctionnait bien en termes de vitesse. Mais une fois que je commence à augmenter les dimensions des tableaux, la vitesse ralentit significativement pour le savetxt. J'ai donc essayé np.save et cela entraîne des économies beaucoup plus rapides, mais le fichier est enregistré au format .npy. Comment pourrais-je lire un tel fichier dans Fortran pour reconstruire le tableau original? D'après ce que j'ai lu, binaire se traduit généralement par la plus faible consommation d'espace et les vitesses les plus rapides.

Dans Fortran 90,

open(10,file='/home/X.npy') 

read(10,*)X 

Je reçois l'erreur suivante: Erreur d'exécution Fortran: nombre réel Bad au point 1 de l'entrée de la liste

modifier: En python je fais cela,

import numpy as np 
x = np.linspace(-50,50,10) 
y = np.linspace(-50,50,10) 
X,Y = np.meshgrid(x,y) 
np.save('X',X) 

Puis en Fortran je fais cela,

program read_unformatted_binary_from_python 
    use iso_c_binding 
    implicit none 

    integer, parameter :: DPR = selected_real_kind(p=15) 
    character(len=*), parameter :: filename = 'X.npy' 

    integer, parameter :: N = 10 
    real(DPR), allocatable, dimension(:, :) :: X 

    allocate(X(N, N)) 


    open(40, file=filename, status='old', access='stream', form='unformatted') 
    read(40) X 
    close(40) 

    write(*,*) X 

end program read_unformatted_binary_from_python 

La sortie commence par 1.8758506894003703E-309 + 171 1.1711999948023422E 5.2274167985502976E-037 8.4474009688929314E + 252 + 180 2.345660E 9.9215260506210473E + 247 + 233 2.1620996769994603E 7.5805790251297605E-096 3.4756671925988047E-152 6.5549091408576423E-260 -50,000000000000000 - 38.888888888888886 -27.777777777777779 et ainsi de suite

Cette erreur ne se produit pas lors de l'utilisation du format .bin, uniquement lors de l'utilisation de np.save, ce qui entraîne la présence de npy.

Answer 1

Vladimir F est correct que vous voulez l'accès «en flux» d'un fichier «binaire brut» dans Fortran. Voici un MWE:

Python

import numpy as np 
A = np.random.rand(10000, 10000) 
print(A.sum()) 
A.tofile('data.bin') 

Fortran

program read_unformatted_binary_from_python 
    use iso_c_binding 
    implicit none 

    integer, parameter :: DPR = selected_real_kind(p=15) 
    character(len=*), parameter :: filename = 'data.bin' 

    integer, parameter :: N = 10000 
    real(DPR), allocatable, dimension(:, :) :: dat 

    allocate(dat(N, N)) 

    open(40, file=filename, status='old', access='stream', form='unformatted') 
    read(40) dat 
    close(40) 

    write(*,*) sum(dat) 

end program read_unformatted_binary_from_python 

Mon exemple Fortran est peut-être un peu plus longtemps que nécessaire, comme je l'utilise de nombreux systèmes différents et suites du compilateur, et la haine aussi grande statique tableaux (je suis un utilisateur Fortran après tout). J'ai codé cela très rapidement avec Python 2.7.x, Numpy 13.x, et gfortran de Homebrew GCC 6.3.0_1 sur un MacBook Pro, mais cela devrait fonctionner sur tous les systèmes.

MISE À JOUR: Un soin particulier doit être pris avec les formes et les tailles de réseau ici. Si dat est alloué pour être plus grand que ce qui est dans le fichier, que le streaming read doit essayer de remplir le tableau entier, appuyez sur le symbole EOF et émettre une erreur. En Python, la méthode np.fromfile() lira jusqu'à EOF puis retournera un tableau 1D avec la longueur appropriée, même si A était à l'origine multidimensionnel. En effet, un binaire brut n'a pas de métadonnées et est simplement une chaîne d'octets contiguës depuis la RAM.

En conséquence, les lignes suivantes en Python produisent le même fichier:

A = np.random.rand(10000, 10000) 
A.tofile('file.whatever') 
A.ravel().tofile('file.whatever') 
A.reshape((100, 1000, 1000)).tofile('file.whatever') 

Ce fichier peut être lu et remodelé comme:

B = np.fromfile('file.whatever').reshape(A.shape) 
B = np.fromfile('file.whatever').reshape((100, 1000, 100, 10)) 
# or something like 
B = np.fromfile('file.whatever') # just a 1D array 
B.resize(A.shape) # resized in-place 

En Fortran, il est facile à lire un fichier brut entier utilisant l'accès au flux sans connaître sa taille à l'avance, bien que vous ayez évidemment besoin d'une sorte d'entrée de l'utilisateur pour remodeler les données:

program read_unformatted_binary_from_python 
    use iso_c_binding 
    implicit none 

    integer, parameter :: DPR = selected_real_kind(p=15) 
    character(len=*), parameter :: filename = 'data.bin' 
    integer :: N = 10000, bytes, reals, M 
    real(DPR), allocatable :: A(:,:), D(:, :), zeros(:) 
    real(DPR), allocatable, target :: B(:) 
    real(DPR), pointer :: C(:, :) 

    allocate(A(N, N)) 

    open(40, file=filename, status='old', access='stream', form='unformatted') 

    read(40) A 
    write(*,*) 'sum of A', sum(A) 

    inquire(unit=40, size=bytes) 
    reals = bytes/8 
    allocate(B(reals)) 

    read(40, pos=1) B 
    write(*,*) 'sum of B', sum(B) 

    ! now reshape B in-place assuming the user wants the first dimension 
    ! (which would be the outer dimension in Python) to be length 100 
    N = 100 
    if(mod(reals, N) == 0) then 
     M = reals/N 
     call C_F_POINTER (C_LOC(B), C, [N, M]) 
     write(*, *) 'sum of C', sum(C) 
     write(*, *) 'shape of C', shape(C) 
    else 
     write(*,*) 'file size is not divisible by N!, did not point C to B' 
    end if 

    ! now reshape B out-of-place with first dimension length 9900, and 
    ! pad the result so that there is no size mismatch error 
    N = 9900 
    M = reals/N 
    if(mod(reals, N) > 0) M=M+1 

    allocate(D(N, M)) 
    allocate(zeros(N), source=real(0.0, DPR)) 
    D = reshape(B, [N, M], pad=zeros) 

    write(*,*) 'sum of D', sum(D) 
    write(*,*) 'shape of D', shape(D) 

    ! obviously you can also work with shape(A) in fortran the same way you 
    ! would use A.shape() in Python, if you already knew that A was the 
    ! correct shape of the data 
    deallocate(D) 
    allocate(D, mold=A) 
    D = reshape(B, shape(A)) 
    write(*,*) 'sum of D', sum(D) 
    write(*,*) 'shape of D', shape(D) 

    ! or, just directly read it in, skipping the whole reshape B part 
    read(40, pos=1) D 
    write(*,*) 'sum of D', sum(D) 

    close(40) 

end program read_unformatted_binary_from_python