matmul.cpp

#include <chrono>
#include <iostream>
#include <random>

const int N = 1000;
const int TRIALS = 10;

// [memo] When "N" equals 1000, # of flops is 2000000000
const double nFlops = 2.0 * 1000 * 1000 * 1000;

// [memo] MacBookPro-13-2020's peak flops
const double pFlops = 32.0 * 2 * 4 * 1000 * 1000 * 1000;

// NxN matrices A, B and C
double A[N][N] = {}, B[N][N] = {}, C[N][N] = {};

// Calculate C = AB
void run() {
    // Implement here

#pragma omp parallel for
    for (int i = 0; i < N; ++i) {
        for (int j = 0; j < N; ++j) {
            C[i][j] = 0;
        }
    }

#pragma omp parallel for
    for (int i = 0; i < N; i++) {
        for (int k = 0; k < N; k++) {
            for (int j = 0; j < N; j++) {
                C[i][j] += A[i][k] * B[k][j];
            }
        }
    }
}

// Check if result is correct
void check_result() {
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            double sum = 0.0;
            for (int k = 0; k < N; k++) {
                sum += A[i][k] * B[k][j];
            }

            if (std::abs(C[i][j] - sum) > 1e-9) {
                std::cout << "Result is incorrect!" << std::endl;
                std::cout << "C[i][j] is " << C[i][j] << std::endl;
                std::cout << "sum is " << sum << std::endl;
                return;
            }
        }
    }

    std::cout << "Result is correct." << std::endl;
}

int main() {
    // Initialize A and B with random values
    std::random_device seed_gen;
    std::mt19937 engine(seed_gen());
    std::uniform_real_distribution<> dist(-1.0, 1.0);

    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            A[i][j] = dist(engine);
            B[i][j] = dist(engine);
        }
    }

    // Run benchmark for TRIALS
    for (int i = 0; i < TRIALS; i++) {
        // Start timer
        const std::chrono::steady_clock::time_point start =
                std::chrono::steady_clock::now();

        run();

        // Stop timer
        const std::chrono::steady_clock::time_point end =
                std::chrono::steady_clock::now();

        // Calculate and print runtime
        std::chrono::duration<double, std::milli> runtime = end - start;
        std::cout << "runtime: " << runtime.count() << " [ms]" << std::endl;

        double aflops = 0;
        aflops = nFlops / runtime.count() * 1000;
        std::cout << "Attained FLOP/s " << aflops << std::endl;

        double e;
        e = aflops / pFlops * 100;
        std::cout << "Efficiency: " << e << " [%]" << std::endl;

        check_result();
    }
}