Matrix ops comparison

core.clj

(ns matrix-comp.core)

(require '[clojure.math :refer (sqrt)])
(require '[clojure.core.matrix :refer (matrix mget add mul inverse det set-current-implementation dot mmul eseq)])
(require '[clojure.core.matrix.linear :refer (norm)])
(require '[criterium.core :refer :all])
(require '[fastmath.matrix :as fm])
(require '[fastmath.vector :as fv])
(import '[mikera.matrixx Matrix Matrixx])
(import '[mikera.vectorz Vector])
(import '[org.ejml.simple SimpleMatrix SimpleBase])
(import '[org.ejml.data DMatrixRMaj])

(set! *warn-on-reflection* true)
(set! *unchecked-math* true)
(set-current-implementation :vectorz)

(defmacro mybench [expr] `(do (println (quote ~expr)) (quick-bench ~expr) (println)))

(def a (matrix [[1.0 1.0 1.0 1.0] [1.0 2.0 2.0 2.0] [2.0 2.0 4.0 4.0] [3.0 3.0 3.0 8.0]]))
(def b (SimpleMatrix. 4 4 true (double-array [1.0 1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 2.0 4.0 4.0 3.0 3.0 3.0 8.0])))
(def c (Matrixx/create [[1.0 1.0 1.0 1.0] [1.0 2.0 2.0 2.0] [2.0 2.0 4.0 4.0] [3.0 3.0 3.0 8.0]]))
(def d (fm/mat4x4 1.0 1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 2.0 4.0 4.0 3.0 3.0 3.0 8.0))
(def u (matrix [2.0 3.0 5.0 7.0]))
(def v (SimpleMatrix. 4 1 true (double-array [2.0 3.0 5.0 7.0])))
(def w (Vector/create [2.0 3.0 5.0 7.0]))
(def z (fv/vec4 2.0 3.0 5.0 7.0))

(println "--------------------------------------------------------------------------------")
(mybench (matrix [[1.0 1.0 1.0 1.0] [0.0 2.0 2.0 2.0] [0.0 0.0 4.0 4.0] [0.0 0.0 0.0 8.0]]))
(mybench (SimpleMatrix. 4 4 true (double-array [1.0 1.0 1.0 1.0 0.0 2.0 2.0 2.0 0.0 0.0 4.0 4.0 0.0 0.0 0.0 8.0])))
(mybench (Matrixx/create [[1.0 1.0 1.0 1.0] [0.0 2.0 2.0 2.0] [0.0 0.0 4.0 4.0] [0.0 0.0 0.0 8.0]]))
(mybench (fm/mat4x4 1.0 1.0 1.0 1.0 0.0 2.0 2.0 2.0 0.0 0.0 4.0 4.0 0.0 0.0 0.0 8.0))

(println "--------------------------------------------------------------------------------")
(mybench (matrix [1.0 2.0 3.0 4.0]))
(mybench (SimpleMatrix. 4 1 true (double-array [1.0 2.0 3.0 4.0])))
(mybench (Vector/create [1.0 2.0 3.0 4.0]))
(mybench (fv/vec4 2.0 3.0 5.0 7.0))

(println "--------------------------------------------------------------------------------")
(mybench (add ^Vector u ^Vector u))
(mybench (.plus ^SimpleMatrix v ^SimpleMatrix v))
(mybench (let [result (.clone ^Matrix c)] (.add ^Matrix result ^Matrix c) result))
(mybench (fv/add z z))

(println "--------------------------------------------------------------------------------")
(mybench (inverse ^Matrix a))
(mybench (.invert ^SimpleMatrix b))
(mybench (.inverse ^Matrix c))
(mybench (fm/inverse d))

(println "--------------------------------------------------------------------------------")
(mybench (mul ^Matrix a ^Matrix a))
(mybench (.elementMult ^SimpleMatrix b ^SimpleMatrix b))
(mybench (let [result (.clone ^Matrix c)] (.multiply ^Matrix result ^Matrix c) result))
(mybench (fm/emulm d d))

(println "--------------------------------------------------------------------------------")
(mybench (mmul ^Matrix a ^Matrix a))
(mybench (.mult ^SimpleMatrix b ^SimpleMatrix b))
(mybench (.innerProduct ^Matrix c ^Matrix c))
(mybench (fm/mulm d d))

(println "--------------------------------------------------------------------------------")
(mybench (mmul ^Matrix a ^Vector u))
(mybench (.mult ^SimpleMatrix b ^SimpleMatrix v))
(mybench (.innerProduct ^Matrix c ^Vector w))
(mybench (fm/mulv d z))

(println "--------------------------------------------------------------------------------")
(mybench (dot ^Vector u ^Vector u))
(mybench (.dot ^SimpleMatrix v ^SimpleMatrix v))
(mybench (.innerProduct ^Vector w ^Vector w))
(mybench (fv/dot z z))

(println "--------------------------------------------------------------------------------")
(mybench (norm ^Matrix u))
(mybench (.normF ^SimpleMatrix v))
(mybench (fv/mag z))

(println "--------------------------------------------------------------------------------")
(mybench (det ^Matrix a))
(mybench (.determinant ^SimpleMatrix b))
(mybench (.determinant ^Matrix c))
(mybench (fm/det d))

(println "--------------------------------------------------------------------------------")
(mybench (mget ^Matrix a 1 1))
(mybench (.get ^SimpleMatrix b 1 1))
(mybench (.get ^Matrix a 1 1))
(mybench (d 1 1))

(println "--------------------------------------------------------------------------------")
(mybench (.getElements ^Matrix a))
(mybench (.data ^DMatrixRMaj (.getMatrix ^SimpleMatrix b)))
(mybench (.getElements ^Matrix c))
(mybench (fm/mat->array2d d))

(println "--------------------------------------------------------------------------------")
(def a3 (matrix [[1.0 1.0 1.0] [1.0 2.0 2.0] [2.0 2.0 4.0]]))
(def b3 (SimpleMatrix. 3 3 true (double-array [1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 4.0])))
(def c3 (Matrixx/create [[1.0 1.0 1.0] [1.0 2.0 2.0] [2.0 2.0 4.0]]))
(def d3 (fm/mat3x3 1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 4.0))
(mybench (inverse ^Matrix a3))
(mybench (.invert ^SimpleMatrix b3))
(mybench (.inverse ^Matrix c3))
(mybench (fm/inverse d3))

(println "--------------------------------------------------------------------------------")
(def a4 (matrix [[1.0 1.0 1.0 1.0] [1.0 2.0 2.0 2.0] [2.0 2.0 4.0 4.0] [3.0 3.0 3.0 8.0]]))
(def b4 (SimpleMatrix. 4 4 true (double-array [1.0 1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 2.0 4.0 4.0 3.0 3.0 3.0 8.0])))
(def c4 (Matrixx/create [[1.0 1.0 1.0 1.0] [1.0 2.0 2.0 2.0] [2.0 2.0 4.0 4.0] [3.0 3.0 3.0 8.0]]))
(def d4 (fm/mat4x4 1.0 1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 2.0 4.0 4.0 3.0 3.0 3.0 8.0))
(mybench (inverse ^Matrix a4))
(mybench (.invert ^SimpleMatrix b4))
(mybench (.inverse ^Matrix c4))
(mybench (fm/inverse d4))

(println "--------------------------------------------------------------------------------")
(def a5 (matrix [[1.0 1.0 1.0 1.0 1.0] [1.0 2.0 2.0 2.0 2.0] [2.0 2.0 4.0 4.0 4.0] [3.0 3.0 3.0 8.0 8.0] [4.0 4.0 4.0 4.0 16.0]]))
(def b5 (SimpleMatrix. 5 5 true (double-array [1.0 1.0 1.0 1.0 1.0 1.0 2.0 2.0 2.0 2.0 2.0 2.0 4.0 4.0 4.0 3.0 3.0 3.0 8.0 8.0 4.0 4.0 4.0 4.0 16.0])))
(def c5 (Matrixx/create [[1.0 1.0 1.0 1.0 1.0] [1.0 2.0 2.0 2.0 2.0] [2.0 2.0 4.0 4.0 4.0] [3.0 3.0 3.0 8.0 8.0] [4.0 4.0 4.0 4.0 16.0]]))
(mybench (inverse ^Matrix a5))
(mybench (.invert ^SimpleMatrix b5))
(mybench (.inverse ^Matrix c5))

project.clj

(defproject matrix-comp "0.1.0-SNAPSHOT"
  :description "FIXME: write description"
  :url "http://example.com/FIXME"
  :license {:name "EPL-2.0 OR GPL-2.0-or-later WITH Classpath-exception-2.0"
            :url "https://www.eclipse.org/legal/epl-2.0/"}
  :dependencies [[org.clojure/clojure "1.11.1"]
                 [net.mikera/core.matrix "0.63.0"]
                 [net.mikera/vectorz-clj "0.48.0"]
                 [org.ejml/ejml-all "0.43"]
                 [generateme/fastmath "2.2.1"]]
  :repl-options {:init-ns matrix-comp.core})

res.txt

--------------------------------------------------------------------------------
(matrix [[1.0 1.0 1.0 1.0] [0.0 2.0 2.0 2.0] [0.0 0.0 4.0 4.0] [0.0 0.0 0.0 8.0]])
             Execution time mean : 839.953820 ns
(SimpleMatrix. 4 4 true (double-array [1.0 1.0 1.0 1.0 0.0 2.0 2.0 2.0 0.0 0.0 4.0 4.0 0.0 0.0 0.0 8.0]))
             Execution time mean : 126.571723 ns
(Matrixx/create [[1.0 1.0 1.0 1.0] [0.0 2.0 2.0 2.0] [0.0 0.0 4.0 4.0] [0.0 0.0 0.0 8.0]])
             Execution time mean : 55.682185 ns
(fm/mat4x4 1.0 1.0 1.0 1.0 0.0 2.0 2.0 2.0 0.0 0.0 4.0 4.0 0.0 0.0 0.0 8.0)
             Execution time mean : 13.038411 ns
--------------------------------------------------------------------------------
(matrix [1.0 2.0 3.0 4.0])
             Execution time mean : 345.015969 ns
(SimpleMatrix. 4 1 true (double-array [1.0 2.0 3.0 4.0]))
             Execution time mean : 50.980032 ns
(Vector/create [1.0 2.0 3.0 4.0])
             Execution time mean : 10.422687 ns
(fv/vec4 2.0 3.0 5.0 7.0)
             Execution time mean : 3.364942 ns
--------------------------------------------------------------------------------
(add u u)
             Execution time mean : 13.722448 ns
(.plus v v)
             Execution time mean : 19.900253 ns
(let [result (.clone c)] (.add result c) result) ;; wrong test
             Execution time mean : 17.693223 ns
(fv/add z z)
             Execution time mean : 4.082166 ns
--------------------------------------------------------------------------------
(inverse a)
             Execution time mean : 368.827894 ns
(.invert b)
             Execution time mean : 58.952959 ns
(.inverse c)
             Execution time mean : 351.509533 ns
(fm/inverse d)
             Execution time mean : 28.930919 ns
--------------------------------------------------------------------------------
(mul a a)
             Execution time mean : 47.156526 ns
(.elementMult b b)
             Execution time mean : 24.918859 ns
(let [result (.clone c)] (.multiply result c) result)
             Execution time mean : 28.382067 ns
(fm/emulm d d)
             Execution time mean : 9.951789 ns
--------------------------------------------------------------------------------
(mmul a a)
             Execution time mean : 104.274636 ns
(.mult b b)
             Execution time mean : 81.108903 ns
(.innerProduct c c)
             Execution time mean : 106.270240 ns
(fm/mulm d d)
             Execution time mean : 19.046061 ns
--------------------------------------------------------------------------------
(mmul a u)
             Execution time mean : 22.328177 ns
(.mult b v)
             Execution time mean : 35.142648 ns
(.innerProduct c w)
             Execution time mean : 19.308981 ns
(fm/mulv d z)
             Execution time mean : 8.300744 ns
--------------------------------------------------------------------------------
(dot u u)
             Execution time mean : 6.465094 ns
(.dot v v)
             Execution time mean : 9.481777 ns
(.innerProduct w w)
             Execution time mean : 3.654975 ns
(fv/dot z z)
             Execution time mean : 2.979366 ns
--------------------------------------------------------------------------------
(norm u)
             Execution time mean : 14.393764 ns
(.normF v)
             Execution time mean : 6.301779 ns
(fv/mag z)
             Execution time mean : 3.320475 ns
--------------------------------------------------------------------------------
(det a)
             Execution time mean : 164.376744 ns
(.determinant b)
             Execution time mean : 5.752944 ns
(.determinant c)
             Execution time mean : 165.898599 ns
(fm/det d)
             Execution time mean : 8.452224 ns
--------------------------------------------------------------------------------
(mget a 1 1)
             Execution time mean : 3.925001 ns
(.get b 1 1)
             Execution time mean : 2.692314 ns
(.get a 1 1)
             Execution time mean : 3.208399 ns
(d 1 1) ;; through functional protocol, not optimized
             Execution time mean : 50.762566 ns
--------------------------------------------------------------------------------
(.getElements a)
             Execution time mean : 9.800249 ns
(.data (.getMatrix b))
             Execution time mean : 1.065920 ns
(.getElements c)
             Execution time mean : 13.100863 ns
(fm/mat->array2d d) ;; oops
             Execution time mean : 858.552802 ns
--------------------------------------------------------------------------------
(inverse a3)
             Execution time mean : 12.658557 ns
(.invert b3)
             Execution time mean : 37.919045 ns
(.inverse c3)
             Execution time mean : 10.052509 ns
(fm/inverse d3)
             Execution time mean : 9.854505 ns
--------------------------------------------------------------------------------
(inverse a4)
             Execution time mean : 355.921660 ns
(.invert b4)
             Execution time mean : 57.804481 ns
(.inverse c4)
             Execution time mean : 350.942934 ns
(fm/inverse d4)
             Execution time mean : 27.216275 ns
--------------------------------------------------------------------------------
(inverse a5)
             Execution time mean : 535.311499 ns
(.invert b5)
             Execution time mean : 121.917799 ns
(.inverse c5)
             Execution time mean : 517.591227 ns