C# Test matrix mul

using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Threading.Tasks;
using System.IO;



namespace MatrixTest1
{
    class Program
    {
        static void Main(string[] args)
        {
            const int TRIALS = 10;
            var A = Matrix.RandomMatrix(2000, 2000, 1.0);
            var B = Matrix.RandomMatrix(2000, 2000, 1.0);
            var mList = new List<Matrix>(TRIALS * 2);

            Console.WriteLine("Warming up, wait...");
            //give VM chance to warm-up
            for (int i = 0; i < TRIALS; i++)
                mList.Add(A * B);

            Console.WriteLine("Starting trials!");
            //actual trials
            var sw = new Stopwatch();
            sw.Restart();
            for (int i = 0; i < TRIALS; i++)
                mList.Add(A * B);
            sw.Stop();
            var sec = sw.Elapsed.TotalSeconds;
            Console.WriteLine($"{TRIALS} trials in {sec:F} seconds, {sec / TRIALS:F4}s per trial");
            Console.ReadLine();
        }
    }


    public class Matrix
    {
        public readonly int rows;
        public readonly int cols;
        private readonly double[][] mat;

        public Matrix(int iRows, int iCols)
        {
            rows = iRows;
            cols = iCols;
            mat = new double[rows][];
            for (int i = 0; i < rows; i++)
                mat[i] = new double[cols];
        }




        public double this[int iRow, int iCol]
        {
            get { return mat[iRow][iCol]; }
            set { mat[iRow][iCol] = value; }
        }

        private double[] this[int iRow]
        {
            get { return mat[iRow]; }
            set { mat[iRow] = value; }
        }

        public static Matrix RandomMatrix(int iRows, int iCols, double dispersion)
        {
            var random = new ThreadRandom();
            var matrix = new Matrix(iRows, iCols);
            for (int i = 0; i < iRows; i++)
                for (int j = 0; j < iCols; j++)
                    matrix[i, j] = random.NextDouble(-dispersion, dispersion);
            return matrix;
        }


        private static void CalcMulRow(Matrix result, Matrix m1, Matrix m2, int index)
        {
            var iRowA = m1[index];
            var iRowC = result[index];
            for (int k = 0; k < m1.cols; k++)
            {
                var kRowB = m2[k];
                var ikA = iRowA[k];
                for (int j = 0; j < result.cols; j++)
                    iRowC[j] += ikA * kRowB[j];
            }
        }


        private static Matrix Multiply(Matrix m1, Matrix m2)
        {
            if (m1.cols != m2.rows)
                throw new Exception("Wrong dimensions of matrix, m1.cols (" + m1.cols + ") != m2.rows (" + m2.rows + ")!");
            var result = new Matrix(m1.rows, m2.cols);
            Parallel.For(0, result.rows, (i) => CalcMulRow(result, m1, m2, i));
            return result;
        }


        public static Matrix operator *(Matrix m1, Matrix m2) => Matrix.Multiply(m1, m2);
    }


    public class ThreadRandom
    {
        private static readonly Random seedRNG = new Random();
        private readonly Random localRNG = new Random(seedRNG.Next());

        public int Next()
        {
            return localRNG.Next();
        }

        public int Next(int maxValue)
        {
            return localRNG.Next(maxValue);
        }

        public int Next(int minValue, int maxValue)
        {
            return localRNG.Next(minValue, maxValue);
        }

        public double NextDouble()
        {
            return localRNG.NextDouble();
        }

        public double NextDouble(double minValue, double maxValue)
        {
            var r = localRNG.NextDouble() * (maxValue - minValue);
            return minValue + r;
        }
    }
}