(bind-func mmmul
(lambda (m1:float* r1:i64 c1:i64 m2:float* r2:i64 c2:i64 result:float*)
(if (= c1 r2)
(let ((r 0) (c 0) (k 0))
(dotimes (r r1)
(dotimes (c c2)
(pset! result (+ c (* r c2)) 0.0)
(dotimes (k c1)
(pset! result (+ c (* r c2))
(+ (pref result (+ c (* r c2)))
(* (pref m1 (+ k (* r c1)))
(pref m2 (+ (* k c2) c))))))))
result)
(begin (printf "bad matrix multiplication [%lld,%lld]*[%lld,%lld]\n"
r1 c1 r2 c2)
(cast null float*)))))