Race the matrix multiplication methods!

MatrixRace.cpp

#include <ctime>
#include <cstring>
#include <iostream>
#include <vector>
#include <random>
using namespace std;

template<typename T> struct Matrix4x4
{
    T v[16];
};

template<typename T> constexpr Matrix4x4<T> Identity4x4()
{
    return {1,0,0,0,0,1,0,0,0,0,1,0,0,0,0,1};
}

template<typename T> bool operator==(
    const Matrix4x4<T>& a,
    const Matrix4x4<T>& b)
{
    return !memcmp(a.v, b.v, sizeof(a.v));
}

template<typename T> bool operator!=(
    const Matrix4x4<T>& a,
    const Matrix4x4<T>& b)
{
    return memcmp(a.v, b.v, sizeof(a.v));
}

template<typename T> Matrix4x4<T> operator*(Matrix4x4<T> a, Matrix4x4<T> b)
{
    Matrix4x4<T> result;
    
    for (int i = 0; i < 4; ++i)
    {
        for (int j = 0; j < 16; j += 4)
        {
            int index = j + i;
            result.v[index] = 0;
            
            for (int k = 0; k < 4; ++k)
            {
                result.v[index] += a.v[k * 4 + i] * b.v[j + k];
            }
        }
    }
    
    return result;
}

template<typename T> void Multiply4x4(T* result, const T* a, const T* b)
{
    for (int i = 0; i < 4; ++i)
    {
        for (int j = 0; j < 16; j += 4)
        {
            int index = j + i;
            result[index] = 0;
            
            for (int k = 0; k < 4; ++k)
            {
                result[index] += a[k * 4 + i] * b[j + k];
            }
        }
    }
}

template<typename T> ostream& operator<<(ostream& stream, Matrix4x4<T> matrix)
{
    for (int i = 0; i < 4; ++i)
    {
        for (int j = i; j < 16; j += 4)
        {
            stream << '[' << matrix.v[j] << ']';
        }
        
        stream << '\n';
    }

    return stream;
}

Matrix4x4<float> Random4x4(mt19937& mt)
{
    Matrix4x4<float> result;
    uniform_real_distribution<float> distribution(-128.0f, 128.0f);
    
    for (auto& value : result.v) value = distribution(mt);
    return result;
}

int64_t GetNanoseconds()
{
    timespec ts;
    clock_gettime(CLOCK_MONOTONIC, &ts);
    return (int64_t)ts.tv_sec * 1000000000 + ts.tv_nsec;
}

struct Timer
{
    int64_t nanoseconds;
    bool running;
};

Timer NewTimer()
{
    return {GetNanoseconds(), true};
}

Timer Stopped(Timer timer)
{
    return timer.running ?
        Timer{GetNanoseconds() - timer.nanoseconds, false} :
        timer;
}

ostream& operator<<(ostream& stream, Timer timer)
{
    int64_t milliseconds = timer.running ?
        (GetNanoseconds() - timer.nanoseconds) / 1000000 :
        timer.nanoseconds / 1000000;
    return stream << milliseconds << "ms";
}

int main(int argc, char** argv)
{
    vector<Matrix4x4<float>> stacks[4];
    
    constexpr int Count = 1 << 22;
    for (auto& s : stacks) s.resize(Count);
    
    auto timer = NewTimer();
    cout << "Generating random values...";
    cout.flush();
    mt19937 mt;
    for (int i = 0; i < Count; ++i)
    {
        stacks[0][i] = Random4x4(mt);
        stacks[1][i] = Random4x4(mt);
    }
    
    cout << timer << "\nMultiplying by value... ";
    cout.flush();
    timer = NewTimer();
    for (int i = 0; i < Count; ++i)
    {
        stacks[2][i] = stacks[0][i] * stacks[1][i];
    }
    
    cout << timer << "\nMultiplying by pointer... ";
    cout.flush();
    timer = NewTimer();
    for (int i = 0; i < Count; ++i)
    {
        Multiply4x4(stacks[3][i].v, stacks[0][i].v, stacks[1][i].v);
    }
    
    cout << timer << "\nresult: " << (stacks[2] == stacks[3]) << endl;
    
    return 0;
}