uniform random 3000x3000 matrix inversion and multiply benchmark tests

bench.clj

(ns clojure-matrix-bench.core
  (:use clatrix.core))

(defn- my-main []
  (let [n 3000
        m (rand n n)]
    (trace (* m (i m)))))

(defn -main []
  (print (time my-main)))

bench.cpp

include <iostream>
#include <Eigen/Dense>

#include <time.h>

using Eigen::MatrixXd;
using namespace std; 

int main()
{
  const int N=3000;
  MatrixXd m = MatrixXd::Random(N,N);
  clock_t old = clock();
  cout << "tr(m * m^-1) =" << endl << (m*m.inverse()).trace() << endl;
  double elapsed = ((double)clock()-old)/CLOCKS_PER_SEC;
  cout << "elapled = " << elapsed << endl;
  getchar();
}

bench.fs

open MathNet.Numerics.LinearAlgebra
open MathNet.Numerics.LinearAlgebra.Double
open MathNet.Numerics.Distributions

open System
open System.Diagnostics

MathNet.Numerics.Control.LinearAlgebraProvider <-
   new MathNet.Numerics.Algorithms.LinearAlgebra.Mkl.MklLinearAlgebraProvider()

let test n =
    let m : Matrix = 
        upcast DenseMatrix.randomCreate n n (ContinuousUniform (0.0,1.0))
    (m * m.Inverse()).Trace()

[<EntryPoint>]
let main argv = 
    let sw = new Stopwatch()
    sw.Start()
    let ans = test 3000
    sw.Stop()
    printfn "%A" ans
    printfn "elapsed = %d ms" <| sw.ElapsedMilliseconds
    0

bench.R

N <- 3000

pt <- proc.time()
m <- matrix(runif(N*N), nrow=N, ncol=N)
ans <- sum(diag(m %*% solve(m)))
pt <- proc.time()-pt

print(paste("ans =",ans))
print(paste("time = ",pt[3]))