mat_mult.ml

(** Run:
      corebuild -pkg core_bench mat_mult.native && sleep 1 && ./mat_mult.native
*)

open Core.Std
open Core_bench.Std

module Array_mat = struct
  type mat3 = float array array
  type vec3 = float array

  let rand_vec3 () =
    let open Random in
    [| float 10.0; float 10.0; float 10.0; |]
  ;;

  let rand_mat3 () =
    [| rand_vec3 (); rand_vec3 (); rand_vec3 (); |]
  ;;

  let zero () =
    [| [| 0.0; 0.0; 0.0; |]
     ; [| 0.0; 0.0; 0.0; |]
     ; [| 0.0; 0.0; 0.0; |]
    |]
  ;;

  let ( * ) m1 m2 =
    let m3 = zero () in
    for i = 0 to 2 do
      for j = 0 to 2 do
        for k = 0 to 2 do
          m3.(i).(j) <- m3.(i).(j) +. m1.(i).(k) *. m2.(k).(j)
        done
      done
    done;
    m3
  ;;

  let ( *| ) mat vec =
    let v3 = [| 0.0; 0.0; 0.0; |] in
    for i = 0 to 2 do
      for j = 0 to 2 do
        v3.(i) <- v3.(i) +. mat.(i).(j) *. vec.(j)
      done
    done
  ;;

  module Expanded = struct
    let ( * ) m1 m2 =
      [| [| m1.(0).(0) *. m2.(0).(0) +. m1.(0).(1) *. m2.(1).(0) +. m1.(0).(2) *. m2.(2).(0)
          ; m1.(0).(0) *. m2.(0).(1) +. m1.(0).(1) *. m2.(1).(1) +. m1.(0).(2) *. m2.(2).(1)
          ; m1.(0).(0) *. m2.(0).(2) +. m1.(0).(1) *. m2.(1).(2) +. m1.(0).(2) *. m2.(2).(2) |]
       ; [| m1.(1).(0) *. m2.(0).(0) +. m1.(1).(1) *. m2.(1).(0) +. m1.(1).(2) *. m2.(2).(0)
          ; m1.(1).(0) *. m2.(0).(1) +. m1.(1).(1) *. m2.(1).(1) +. m1.(1).(2) *. m2.(2).(1)
          ; m1.(1).(0) *. m2.(0).(2) +. m1.(1).(1) *. m2.(1).(2) +. m1.(1).(2) *. m2.(2).(2) |]
       ; [| m1.(2).(0) *. m2.(0).(0) +. m1.(2).(1) *. m2.(1).(0) +. m1.(2).(2) *. m2.(2).(0)
          ; m1.(2).(0) *. m2.(0).(1) +. m1.(2).(1) *. m2.(1).(1) +. m1.(2).(2) *. m2.(2).(1)
          ; m1.(2).(0) *. m2.(0).(2) +. m1.(2).(1) *. m2.(1).(2) +. m1.(2).(2) *. m2.(2).(2) |]
      |]
    ;;

    let ( *| ) mat vec =
      [| mat.(0).(0) *. vec.(0) +. mat.(0).(1) *. vec.(1) +. mat.(0).(2) *. vec.(2)
       ; mat.(1).(0) *. vec.(0) +. mat.(1).(1) *. vec.(1) +. mat.(1).(2) *. vec.(2)
       ; mat.(2).(0) *. vec.(0) +. mat.(2).(1) *. vec.(1) +. mat.(2).(2) *. vec.(2)
      |]
    ;;
  end
end

module Record_mat = struct
  type mat3 = { m00 : float; m01 : float; m02 : float;
                m10 : float; m11 : float; m12 : float;
                m20 : float; m21 : float; m22 : float;
              }
  type vec3 = { v0 : float; v1 : float; v2 : float; }

  let vec3_of_array vec =
    { v0 = vec.(0); v1 = vec.(1); v2 = vec.(2); }
  ;;

  let mat3_of_array mat =
    {
      m00 = mat.(0).(0); m01 = mat.(0).(1); m02 = mat.(0).(2);
      m10 = mat.(1).(0); m11 = mat.(1).(1); m12 = mat.(1).(2);
      m20 = mat.(2).(0); m21 = mat.(2).(1); m22 = mat.(2).(2);
    }
  ;;

  let ( * ) m1 m2 =
    {
      m00 = m1.m00 *. m2.m00 +. m1.m01 *. m2.m10 +. m1.m02 *. m2.m20;
      m01 = m1.m00 *. m2.m01 +. m1.m01 *. m2.m11 +. m1.m02 *. m2.m21;
      m02 = m1.m00 *. m2.m02 +. m1.m01 *. m2.m12 +. m1.m02 *. m2.m22;
      m10 = m1.m10 *. m2.m00 +. m1.m11 *. m2.m10 +. m1.m12 *. m2.m20;
      m11 = m1.m10 *. m2.m01 +. m1.m11 *. m2.m11 +. m1.m12 *. m2.m21;
      m12 = m1.m10 *. m2.m02 +. m1.m11 *. m2.m12 +. m1.m12 *. m2.m22;
      m20 = m1.m20 *. m2.m00 +. m1.m21 *. m2.m10 +. m1.m22 *. m2.m20;
      m21 = m1.m20 *. m2.m01 +. m1.m21 *. m2.m11 +. m1.m22 *. m2.m21;
      m22 = m1.m20 *. m2.m02 +. m1.m21 *. m2.m12 +. m1.m22 *. m2.m22;
    }
  ;;

  let ( *| ) mat vec =
    {
      v0 = mat.m00 *. vec.v0 +. mat.m01 *. vec.v1 +. mat.m02 *. vec.v2;
      v1 = mat.m10 *. vec.v0 +. mat.m11 *. vec.v1 +. mat.m12 *. vec.v2;
      v2 = mat.m20 *. vec.v0 +. mat.m21 *. vec.v1 +. mat.m22 *. vec.v2;
    }
  ;;
end

let main () =
  Random.self_init ();
  let m1 = Array_mat.rand_mat3 () in
  let m2 = Array_mat.rand_mat3 () in
  let vec = Array_mat.rand_vec3 () in
  let m1' = Record_mat.mat3_of_array m1 in
  let m2' = Record_mat.mat3_of_array m2 in
  let vec' = Record_mat.vec3_of_array vec in
  Command.run (Bench.make_command [
      Bench.Test.create ~name:"Array_mat.mat_mult" (fun () ->
          ignore Array_mat.(m1 * m2));
      Bench.Test.create ~name:"Array_mat.Expanded.mat_mult" (fun () ->
          ignore Array_mat.Expanded.(m1 * m2));
      Bench.Test.create ~name:"Record_mat.mat_mult" (fun () ->
          ignore Record_mat.(m1' * m2'));
      Bench.Test.create ~name:"Array_mat.vec_mult" (fun () ->
          ignore Array_mat.(m1 *| vec));
      Bench.Test.create ~name:"Array_mat.Expanded.vec_mult" (fun () ->
          ignore Array_mat.Expanded.(m1 *| vec));
      Bench.Test.create ~name:"Record_mat.vec_mult" (fun () ->
          ignore Record_mat.(m1' *| vec'));
    ])
;;

let () = main ();;

Estimated testing time 3m (6 benchmarks x 30s). Change using -quota SECS.
┌─────────────────────────────┬──────────┬────────────┐
│ Name                        │ Time/Run │ Percentage │
├─────────────────────────────┼──────────┼────────────┤
│ Array_mat.mat_mult          │ 107.74ns │    100.00% │
│ Array_mat.Expanded.mat_mult │  41.12ns │     38.17% │
│ Record_mat.mat_mult         │  13.75ns │     12.77% │
│ Array_mat.vec_mult          │  25.49ns │     23.66% │
│ Array_mat.Expanded.vec_mult │  13.13ns │     12.19% │
│ Record_mat.vec_mult         │   7.23ns │      6.71% │
└─────────────────────────────┴──────────┴────────────┘