Управление кэшем для умножения разреженных матриц с использованием OpenMP - PullRequest
Новая
       41

Управление кэшем для умножения разреженных матриц с использованием OpenMP

0 голосов
/ 29 октября 2018

У меня проблемы с тем, что я считаю ложным кешированием, я получаю лишь небольшое ускорение при использовании следующего кода по сравнению с несравненной версией.

matrix1 и matrix2 - это разреженные матрицы в структурев формате (row, col, val). 

void pMultiply(struct SparseRow *matrix1, struct SparseRow *matrix2, int m1Rows, int m2Rows, struct SparseRow **result) {

*result = malloc(1 * sizeof(struct SparseRow));

int resultNonZeroEntries = 0;

#pragma omp parallel for atomic
for(int i = 0; i < m1Rows; i++)
{
    int curM1Row = matrix1[i].row;
    int curM1Col = matrix1[i].col;
    float curM1Value = matrix1[i].val;

    for(int j = 0; j < m2Rows; j++)
    {

        int curM2Row = matrix2[j].row;
        int curM2Col = matrix2[j].col;
        float curM2Value = matrix2[j].val;

        if(curM1Col == curM2Row)
        {
            *result = realloc(*result, 
            (sizeof(struct SparseRow)*(resultNonZeroEntries+1)));

            (*result)[resultNonZeroEntries].row = curM1Row;
            (*result)[resultNonZeroEntries].col = curM2Col;
            (*result)[resultNonZeroEntries].val += curM1Value*curM2Value;
            resultNonZeroEntries++;
            break;
        }

    }
}

1 Ответ

0 голосов
/ 29 октября 2018

Несколько вопросов там:

  • Как упомянул Брайан Брохерс, пункт #pragma omp atomic должен быть помещен непосредственно перед строкой, которая должна быть защищена от состояния гонки.

  • Перераспределение памяти на каждом шаге, вероятно, снижает производительность. Если память нельзя перераспределить на месте и ее нужно скопировать в другое место, это будет медленным. Это также источник ошибок, так как значение указателя result изменено. Другие потоки продолжают работать, пока происходит перераспределение, и могут пытаться получить доступ к памяти по «старому» адресу, или несколько потоков могут попытаться перераспределить results одновременно. Размещение всей части сложения realloc + в критической секции было бы более безопасным, но по существу сериализовало бы функцию для чего угодно, кроме проверки равенства индексов строк / столбцов за счет значительных накладных расходов. Потоки должны работать с локальным буфером, а затем объединить их результаты на более позднем этапе. Перераспределение должно выполняться блоками достаточного размера. 

    // Make sure this will compile even without openmp + include memcpy
#include <string.h>
#ifdef _OPENMP
   #define thisThread omp_thread_num()
   #define nThreads omp_num_threads()
#else
   #define  thisThread 0
   #define  nThreads 1
#endif

// shared variables
int totalNonZero,*copyIndex,*threadNonZero;
#pragma omp parallel
{
// each thread now initialize a local buffer and local variables 
int localNonZero = 0;
int allocatedSize = 1024;
SparseRow *localResult = malloc(allocatedSize  * sizeof(*SparseRow));

// one thread initialize an array
#pragma omp single
{
threadNonZero=malloc(nThreads*sizeof(int));copyIndex=malloc((nThreads+1)*sizeof(int));
}

#pragma omp for
for (int i = 0; i < m1Rows; i++){
    /* 
     * do the same as your initial code but:
     * realloc an extra 1024 lines each time localNonZeros exceeds allocatedSize
     * fill the local buffer and increment the localNonZeros counter
     * this is safe, no need to use critical / atomic clauses
     */
    }

copyIndex[thisThread]=localNonZero; //put number of non zero into a shared variable
#pragma omp barrier

// Wrap_up : check how many non zero values for each thread, allocate the output and check where each thread will copy its local buffer
#pragma omp single
{
    copyIndex[0]=0;
    for (int i=0; i<nThreads; ii++)
        copyIndex[i+1]=localNonZero[i]+copyIndex[i];
    result=malloc(copyIndex[nThreads+1]*sizeof(*SparseRow));
}

// Copy the results from local to global result   memcpy(&result[copyIndex[thisThread]],localResult,localNonZero*sizeof(*SparseRow);

// Free memory
free(localResult);
#pragma omp single
{
free(copyIndex);
free(localNonZero);
}
} // end parallel

  • Обратите внимание, что алгоритм будет генерировать дубликаты, например, если первая матрица содержит значения в позициях (1,10) и (1,20), а вторая (10,5) и (20,5), в результате будет две (1,5) строки. В какой-то момент понадобится функция сжатия, которая объединит повторяющиеся строки.

...