MatrixMultiplyTest

using System.Diagnostics;
using BenchmarkDotNet.Attributes;
using BenchmarkDotNet.Running;

namespace MatrixMultiplyTest;

public class Program
{
    private const int MatrixSize = 1000;

    private static readonly double[,] TestMatrixA = new double[MatrixSize, MatrixSize];
    private static readonly double[,] TestMatrixB = new double[MatrixSize, MatrixSize];

    static Program()
    {
        var stopwatch = Stopwatch.StartNew();
        var rand = new Random(0);
        for (var i = 0; i < MatrixSize; i++)
        {
            for (var j = 0; j < MatrixSize; j++)
            {
                TestMatrixA[i, j] = rand.NextDouble();
                TestMatrixB[i, j] = rand.NextDouble();
            }
        }
        stopwatch.Stop();
        Console.WriteLine($"Initialized in {stopwatch.ElapsedMilliseconds} ms.");
    }

    private static double[,] ParallelMultiplicationTest(double[,] matrixA, double[,] matrixB)
    {
        var rowsA = matrixA.GetLength(0);
        var colsA = matrixA.GetLength(1);
        var rowsB = matrixB.GetLength(0);
        var colsB = matrixB.GetLength(1);
        Debug.Assert(colsA == rowsB);

        var result = new double[rowsA, colsB];
        Parallel.For(0, rowsA, i =>
        {
            for (var j = 0; j < colsB; j++)
            {
                result[i, j] = 0.0;
                for (var k = 0; k < colsA; k++)
                {
                    result[i, j] += matrixA[i, k] * matrixB[k, j];
                }
            }
        });
        return result;
    }

    private static double[,] SerialMultiplicationTest(double[,] matrixA, double[,] matrixB)
    {
        var rowsA = matrixA.GetLength(0);
        var colsA = matrixA.GetLength(1);
        var rowsB = matrixB.GetLength(0);
        var colsB = matrixB.GetLength(1);
        Debug.Assert(colsA == rowsB);

        var result = new double[rowsA, colsB];
        for (var i = 0; i < rowsA; i++)
        {
            for (var j = 0; j < colsB; j++)
            {
                result[i, j] = 0.0;
                for (var k = 0; k < colsA; k++)
                {
                    result[i, j] += matrixA[i, k] * matrixB[k, j];
                }
            }
        }
        return result;
    }

    [Benchmark]
    public void ParallelMultiplication() => ParallelMultiplicationTest(TestMatrixA, TestMatrixB);

    [Benchmark]
    public void SerialMultiplication() => SerialMultiplicationTest(TestMatrixA, TestMatrixB);

    static void Main()
    {
        Validate();
#if RELEASE
        _ = BenchmarkRunner.Run<Program>();
#endif
    }

    private static void Validate()
    {
        var sStopwatch = Stopwatch.StartNew();
        var sResult = SerialMultiplicationTest(TestMatrixA, TestMatrixB);
        sStopwatch.Stop();

        var pStopwatch = Stopwatch.StartNew();
        var pResult = ParallelMultiplicationTest(TestMatrixA, TestMatrixB);
        pStopwatch.Stop();

        Console.WriteLine($"Parallel time: {pStopwatch.ElapsedMilliseconds} ms. " +
                          $"Serial time: {sStopwatch.ElapsedMilliseconds} ms.");

        if (pResult.GetLength(0) != sResult.GetLength(0))
        {
            throw new Exception($"pResult.GetLength(0) {pResult.GetLength(0)} != sResult.GetLength(0) {sResult.GetLength(0)}");
        }

        if (pResult.GetLength(1) != sResult.GetLength(1))
        {
            throw new Exception($"pResult.GetLength(1) {pResult.GetLength(1)} != sResult.GetLength(1) {sResult.GetLength(1)}");
        }

        var maxDiff = double.MinValue;
        for (var i = 0; i < pResult.GetLength(0); i++)
        {
            for (var j = 0; j < pResult.GetLength(1); j++)
            {
                var diff = Math.Abs(pResult[i, j] - sResult[i, j]);
                if (diff > maxDiff) maxDiff = diff;
            }
        }
        Console.WriteLine($"maxDiff = {maxDiff}");

        for (var i = 0; i < pResult.GetLength(0); i++)
        {
            for (var j = 0; j < pResult.GetLength(1); j++)
            {
                // ReSharper disable once CompareOfFloatsByEqualityOperator
                if (pResult[i, j] != sResult[i, j])
                {
                    throw new Exception($"pResult[{i}, {j}] {pResult[i, j]} != sResult[{i}, {j}] {sResult[i, j]}");
                }
            }
        }

        Console.WriteLine("Validated");
    }
}