Как я могу дополнительно оптимизировать кеширование моего алгоритма умножения матриц? - PullRequest
Новая
       8

Как я могу дополнительно оптимизировать кеширование моего алгоритма умножения матриц?

0 голосов
/ 16 февраля 2019

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

Матрицы в основной столбец порядок.

Любые предложения приветствуются! 

const int blockSize = 16;

int min(int a, int b) {
        if (a < b)
                return a;
        else
                return b;
}

void readInMatrix(double* inputMatrix, double* outMatrix, int beginX, int beginY, int endX, int endY, int blockSize, int n) {

        int elementCount = blockSize*blockSize;
        for (int i = beginX; i < endX; i++) {

                for (int j = beginY; j < endY; j++) {

                        *outMatrix = *(inputMatrix + i * n + j);
                        outMatrix++;
                }


                for (int l = endY; l < beginY + blockSize; l++) {
                        *outMatrix = 0;
                        outMatrix++;
                }


        }

    for (int n = endX; n < beginX + blockSize; n++) {
                for (int m = 0; m < blockSize; m++) {
                        *outMatrix = 0;
                        outMatrix++;
                }

        }
    outMatrix -= elementCount;

}



void writeMatrix(double* inputMatrix, double* outMatrix, int beginX, int beginY, int endX, int endY, int blockSize, int n) {

        for (int i = beginX; i < endX; i++) {

                for (int j = beginY; j < endY; j++) {

                        *(outMatrix + n * i + j) = *inputMatrix;
                        inputMatrix++;
                }

                inputMatrix = inputMatrix - (endY - beginY) + blockSize;

        }

    inputMatrix -= (endY - beginY) * (endX - beginX);


}

void square_dgemm(int n, double* A, double* B, double* C) {

        if (n > blockSize) {

        double ac[blockSize*blockSize];
        double bc[blockSize*blockSize];
        double cc[blockSize*blockSize];



                for (int x = 0; x < n; x += blockSize) {
                        for (int y = 0; y < n; y += blockSize) {
                                readInMatrix(C, cc, x, y, min(x + blockSize, n), min(y + blockSize, n), blockSize, n);






                                for (int z = 0; z < n; z += blockSize) {

                                        readInMatrix(A, ac, z, y, min(z + blockSize, n), min(y + blockSize, n), blockSize, n);

                                        readInMatrix(B, bc, x, z, min(x + blockSize, n), min(z + blockSize, n), blockSize, n);



                                                //for x
        for (int i = 0; i < blockSize; i++) {
                // for y
                for (int j = 0; j < blockSize; j++) {
                        double cij = *(cc + i * blockSize + j);
                        for (int k = 0; k < blockSize; k+=4) {
                                cij += *(ac + j +  blockSize * k) * *(bc + i * blockSize + k) + *(ac + j +  blockSize * (k+1)) * *(bc + i * blockSize + k+1) + *(ac + j +  blockSize * (k+2)) * *(bc + i * blockSize + k+2) + *(ac + j +  b$

                        }
                        *(cc + i * blockSize + j) = cij;

                }

        }


                                }

                                writeMatrix(cc, C, x, y, min(x + blockSize, n), min(y + blockSize, n), blockSize, n);

                        }

                }

        }
    else {

                for (int i = 0; i < n; ++i)
                        /* For each column j of B */
                        for (int j = 0; j < n; ++j)
                        {
                                /* Compute C(i,j) */
                                double cij = C[i + j * n];
                                for (int k = 0; k < n; k++)
                                        cij += A[i + k * n] * B[k + j * n];
                                C[i + j * n] = cij;
                        }


        }


}

Результаты: Размер: 31 Mflop / s: 5612,97Процент: 13.08 Размер: 32 Mflop / s: 6313.5 Процент: 14.72 Размер: 96 Mflop / s: 6666.09 Процент: 15.54 Размер: 97 Mflop / s: 4480.09 Процент: 10.44 Размер: 127 Mflop / s: 6520.47 Процент: 15.20 Размер: 128Mflop / s: 6674,19%: 15,56 Размер: 129 Mflop / s: 4731,75%: 11,03 Размер: 191 Mflop / s: 6470,42%: 15,08 Размер: 192 Mflop / s: 6524,03%: 15,21 Размер: 229 Mflop / s: 5791,24%: 13.50 Размер: 255 Мфлоп / с: 6469.74 Процент: 15.08 Размер: 256 Мфлоп / с: 6216.15 Процент: 14.49 Размер: 257 Мфлоп / с: 5616.37 Процент: 13.09 Размер: 319 Мфлоп / с: 6522.59 Процент: 15.20 Размер: 320 Мфлоп/ s: 6281.06 Процент: 14.64 Размер: 321 Mflop / s: 5796.06 Процент: 13.51 Размер: 417 Mflop / s: 5982.09 Процент: 13.94 Размер: 479 Mflop / s: 6422.19 Процент: 14.97 Размер: 480 Mflop / s: 6502.51 Процент:15.16 Размер: 511 MФлоп / с: 6407,95 Процент: 14,94 Размер: 512 Мфлоп / с: 6407,58 Процент: 14,94 Размер: 639 Мфлоп / с: 6460,55 Процент: 15,06 Размер: 640 Мфлоп / с: 6173,96 Процент: 14,39 Размер: 767 Мфлоп / с: 6346: 14,79 Размер: 768 Мфлоп / с: 6490,66 Процент: 15,13 Размер: 769 Мфлоп / с: 6039,4 Процент: 14,08 Средний процент Пика = 14,3374

...