Быстрая итеративная запись в файл hdf5 - PullRequest
Новая
       115

Быстрая итеративная запись в файл hdf5

0 голосов
/ 26 октября 2019

Некоторое время назад я написал параллельный код FORTRAN, который диагонализирует очень большие плотные матрицы на суперкомпьютере. Эти матрицы считываются из плотных фрагментированных наборов данных hdf5. Теперь я хочу использовать этот код для очень разреженных матриц, которые построены с использованием Python.

Однако, когда я пытаюсь записать свои данные в плотный файл hdf5, это занимает очень много времени. Разреженная матрица состоит из 3x3 ненулевых блоков и сохраняется с использованием трех массивов: rows, cols, data. Я попытался написать каждый блок итеративно: 

fl = h5py.File(filepath, 'w')
dataset = fl.create_dataset("matrix", shape, dtype='d',
                            chunks=(60, 60), compression='szip',
                            fillvalue=0)

for row, col, val in zip(rows, cols, data):
    dataset[row*3: row*3 + 3, col*3: col*3 + 3] = val

fl.close()

Для маленькой матрицы, состоящей из 14848 ненулевых блоков (плотная форма (1536, 1536)), для записи требуется 2,6 секунды. И мне нужно написать матрицы, которые более чем в 100 раз больше (с гораздо большей разреженностью).

1 Ответ

1 голос
/ 26 октября 2019

Я не знаю, поможет ли это, как по скорости, так и по удобству, но:

scipy.sparse имеет формат сжатия блоков, который напоминает мне ваши данные. Это не совсем то же самое.

Из документов sparse.bsr_matrix: 

In [375]: >>> indptr = np.array([0, 2, 3, 6]) 
     ...: >>> indices = np.array([0, 2, 2, 0, 1, 2]) 
     ...: >>> data = np.array([1, 2, 3, 4, 5, 6]).repeat(4).reshape(6, 2, 2) 
     ...: M = sparse.bsr_matrix((data,indices,indptr), shape=(6, 6)) 
     ...:  
In [377]: M                                                                     
Out[377]: 
<6x6 sparse matrix of type '<class 'numpy.int64'>'
    with 24 stored elements (blocksize = 2x2) in Block Sparse Row format>
In [378]: M.data                                                                
Out[378]: 
array([[[1, 1],
        [1, 1]],

       [[2, 2],
        [2, 2]],

       [[3, 3],
        [3, 3]],

       [[4, 4],
        [4, 4]],

       [[5, 5],
        [5, 5]],

       [[6, 6],
        [6, 6]]])
In [379]: M.data.shape                                                          
Out[379]: (6, 2, 2)
In [380]: M.indptr                                                              
Out[380]: array([0, 2, 3, 6], dtype=int32)
In [381]: M.indices                                                             
Out[381]: array([0, 2, 2, 0, 1, 2], dtype=int32)

Это сжатый формат с indptr и indices вместо col иrow массивы. sparse не имеет блочной версии формата coo.

В любом случае, sparse имеет (относительно) быстрые методы преобразования между форматами. 

In [382]: Mo = M.tocoo()                                                        

In [384]: (Mo.row, Mo.col, Mo.data)                                             
Out[384]: 
(array([0, 0, 1, 1, 0, 0, 1, 1, 2, 2, 3, 3, 4, 4, 5, 5, 4, 4, 5, 5, 4, 4,
        5, 5], dtype=int32),
 array([0, 1, 0, 1, 4, 5, 4, 5, 4, 5, 4, 5, 0, 1, 0, 1, 2, 3, 2, 3, 4, 5,
        4, 5], dtype=int32),
 array([1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5, 6, 6,
        6, 6]))

Эти данныеможет использоваться для заполнения массива zeros одним выражением: 

In [385]: A = np.zeros((6,6),int)                                               
In [386]: A[Mo.row, Mo.col] = Mo.data                                           
In [387]: A                                                                     
Out[387]: 
array([[1, 1, 0, 0, 2, 2],
       [1, 1, 0, 0, 2, 2],
       [0, 0, 0, 0, 3, 3],
       [0, 0, 0, 0, 3, 3],
       [4, 4, 5, 5, 6, 6],
       [4, 4, 5, 5, 6, 6]])
In [388]: M.A                                                                   
Out[388]: 
array([[1, 1, 0, 0, 2, 2],
       [1, 1, 0, 0, 2, 2],
       [0, 0, 0, 0, 3, 3],
       [0, 0, 0, 0, 3, 3],
       [4, 4, 5, 5, 6, 6],
       [4, 4, 5, 5, 6, 6]])

https://docs.h5py.org/en/stable/high/dataset.html#fancy-indexing предупреждает, что необычное индексирование h5py может быть медленным, особенно если оно охватывает куски. Тем не менее, это может быть быстрее, чем итеративная запись срезов 3x3.

Таким образом, неизвестны:

  • как преобразовать формат вашего блока в bsr
  • скоростьbsr.tocoo() шаг
  • относительная скорость фантазии h5py запись
...