Как оптимизировать и ускорить умножение матрицы в C ++? - PullRequest
1 голос
/ 19 марта 2019

это оптимизированная реализация умножения матриц, и эта подпрограмма выполняет операцию умножения матриц. C: = C + A * B (где A, B и C - матрицы размером n на n, сохраненные в основном формате столбца) На выходе A и B сохраняют свои входные значения.

void matmul_optimized(int n, int *A, int *B, int *C)
{
    // to the effective bitwise calculation
    // save the matrix as the different type
    int i, j, k;
    int cij;
    for (i = 0; i < n; ++i) {
        for (j = 0; j < n; ++j) {
            cij = C[i + j * n]; // the initialization into C also, add separate additions to the product and sum operations and then record as a separate variable so there is no multiplication
            for (k = 0; k < n; ++k) {
                cij ^= A[i + k * n] & B[k + j * n]; // the multiplication of each terms is expressed by using & operator the addition is done by ^ operator.
            }
            C[i + j * n] = cij; // allocate the final result into C         }
    }
}

как мне ускорить умножение матрицы на основе вышеуказанной функции / метода?

эта функция проверяется до 2048 на матрице 2048.

функция matmul_optimized выполняется с помощью matmul.

#include <stdio.h>
#include <stdlib.h>

#include "cpucycles.c"
#include "helper_functions.c"
#include "matmul_reference.c"
#include "matmul_optimized.c"

int main()
{
    int i, j;
    int n = 1024; // Number of rows or columns in the square matrices
    int *A, *B;   // Input matrices
    int *C1, *C2; // Output matrices from the reference and optimized implementations

    // Performance and correctness measurement declarations
    long int CLOCK_start, CLOCK_end, CLOCK_total, CLOCK_ref, CLOCK_opt;
    long int COUNTER, REPEAT = 5;
    int difference;
    float speedup;

    // Allocate memory for the matrices
    A = malloc(n * n * sizeof(int));
    B = malloc(n * n * sizeof(int));
    C1 = malloc(n * n * sizeof(int));
    C2 = malloc(n * n * sizeof(int));

    // Fill bits in A, B, C1
    fill(A, n * n);
    fill(B, n * n);
    fill(C1, n * n);

    // Initialize C2 = C1
    for (i = 0; i < n; i++)
        for (j = 0; j < n; j++)
            C2[i * n + j] = C1[i * n + j];

    // Measure performance of the reference implementation
    CLOCK_total = 0;
    for (COUNTER = 0; COUNTER < REPEAT; COUNTER++)
    {
        CLOCK_start = cpucycles();
        matmul_reference(n, A, B, C1);
        CLOCK_end = cpucycles();
        CLOCK_total = CLOCK_total + CLOCK_end - CLOCK_start;
    }
    CLOCK_ref = CLOCK_total / REPEAT;
    printf("n=%d Avg cycle count for reference implementation = %ld\n", n, CLOCK_ref);

    // Measure performance of the optimized implementation
    CLOCK_total = 0;
    for (COUNTER = 0; COUNTER < REPEAT; COUNTER++)
    {
        CLOCK_start = cpucycles();
        matmul_optimized(n, A, B, C2);
        CLOCK_end = cpucycles();
        CLOCK_total = CLOCK_total + CLOCK_end - CLOCK_start;
    }
    CLOCK_opt = CLOCK_total / REPEAT;
    printf("n=%d Avg cycle count for optimized implementation = %ld\n", n, CLOCK_opt);

    speedup = (float)CLOCK_ref / (float)CLOCK_opt;

    // Check correctness by comparing C1 and C2
    difference = 0;
    for (i = 0; i < n; i++)
        for (j = 0; j < n; j++)
            difference = difference + C1[i * n + j] - C2[i * n + j];

    if (difference == 0)
        printf("Speedup factor = %.2f\n", speedup);
    if (difference != 0)
        printf("Reference and optimized implementations do not match\n");

    //print(C2, n);

    free(A);
    free(B);
    free(C1);
    free(C2);
    return 0;
}

Ответы [ 2 ]

0 голосов
/ 02 мая 2019

Оптимизация умножения матрицы на матрицу требует пристального внимания к ряду вопросов:

  • Во-первых, вы должны уметь использовать векторные инструкции. Только векторные инструкции могут получить доступ к параллелизму, присущему архитектуре. Таким образом, либо ваш компилятор должен иметь возможность автоматически отображать на векторные инструкции, либо вы должны делать это вручную, например, вызывая встроенную векторную библиотеку для инструкций AVX-2 (для архитектур x86).

  • Далее необходимо обратить особое внимание на иерархию памяти. Если вы этого не сделаете, ваша производительность может легко упасть до пика ниже 5%.

  • Как только вы сделаете это правильно, вы, надеюсь, разбите вычисление на достаточно маленькие вычислительные блоки, которые вы также можете распараллелить через OpenMP или pthreads.

Документ, который тщательно описывает все, что требуется, можно найти по адресу http://www.cs.utexas.edu/users/flame/laff/pfhp/LAFF-On-PfHP.html. (Это очень большая работа в процессе.) В конце всего этого у вас будет реализация, которая будет близка к производительность достигается с помощью высокопроизводительных библиотек, таких как Math Kernel Library (MKL) Intel или программного обеспечения BLAS-подобных библиотек (BLIS).

(И на самом деле вы МОЖЕТЕ также эффективно включить алгоритм Штрассена. Но это другая история, рассказанная в Разделе 3.5.3 этих заметок.) ​​

Вы можете найти следующую тему релевантной: Как BLAS достигает такой экстремальной производительности?

0 голосов
/ 19 марта 2019

Вы можете попробовать алгоритм типа Штрассен или Копперсмит-Виноград , а также здесь есть хороший пример .Или, может быть, попробуйте параллельные вычисления, например future :: task или std :: thread

...