;;; -*- mode: lisp -*-
;;; $Id: matrix.poplisp,v 1.0 2002/05/03 14:16:00 dada Exp $

(proclaim '(optimize (speed 3) (space 0) (compilation-speed 0) (debug 0) (safety 0)))

(defun matmul (a b c n m k)
  (declare (optimize (speed 3) (safety 0) (debug 0))
           (type (simple-array (unsigned-byte 32) (*)) a b c)
           (fixnum n m k))
  (let ((sum 0)
        (i1 (- m))
        (k2 0))
    (declare (type (unsigned-byte 32) sum) (type fixnum i1 k2))
    (dotimes (i n c)
      (declare (fixnum i))
      (setf i1 (+ i1 m)) ;; i1=i*m
      (dotimes (j k)
        (declare (fixnum j))
        (setf sum 0)
        (setf k2 (- k))
        (dotimes (l m)
          (declare (fixnum l))
          (setf k2 (+ k2 k)) ;; k2= l*k
          (setf sum (the (unsigned-byte 32) (+ (the (unsigned-byte 32) sum) 
                                               (the (unsigned-byte 32) (* (aref a (+ i1 l))
                                                                          (aref b (+ k2 j))))))))
        (setf (aref c (+ i1 j)) sum)))))

(defun make-matrix (rows cols)
  (declare (type (unsigned-byte 32) rows cols)
           (optimize (speed 3) (safety 0))); (hcl:fixnum-safety 0)))
  (let* ((space (* rows cols))
         (matrix (make-array space
                             :element-type '(unsigned-byte 32))))
    (declare (type (simple-array (unsigned-byte 32) (*)) matrix)
             (fixnum space))
    (loop :for i :of-type fixnum :from 0 :below space
          :do (setf (aref matrix i) (1+ i)))
    matrix))

(let ((n (parse-integer (or (car pop11::poparglist) "1"))))
(declare (fixnum n)    
     (optimize (speed 3) (debug 0) (safety 0)))
(let* ((m1 (make-matrix 30 30))
   (m2 (make-matrix 30 30))
   (m3 (make-matrix 30 30))
   (mm (make-array '(30 30) :element-type '(unsigned-byte 32) :displaced-to m3)))
  (loop repeat n do (matmul m1 m2 m3 30 30 30))
  (format t "~A ~A ~A ~A~%"
      (aref mm 0 0) (aref mm 2 3) (aref mm 3 2) (aref mm 4 4))))