anchor.diff:         do (setf (row-major-aref result i) (* A (row-major-aref B i))))
       result))
   
  +(defmethod .* ((A number) (B simple-array))
  +  (declare (optimize speed)
  +	   (type double-float A)
  +	   (type (simple-array double-float) B))
  +  (with-result (result (array-dimensions B))
  +    (loop for i from 0 below (array-total-size B)
  +       do (setf (row-major-aref result i) (* A (row-major-aref B i))))
  +    result))
  +
  +(defmethod .* ((A number) (B array))
  +  (declare (optimize speed)
  +	   (type double-float A)
  +	   (type (array double-float) B))
  +  (with-result (result (array-dimensions B))
  +    (loop for i from 0 below (array-total-size B)
  +       do (setf (row-major-aref result i) (* A (row-major-aref B i))))
  +    result))
  +
   (defmethod .* ((A array) (B number))
     (.* B A))
   
   (defmethod .* ((A number) (B vector))
  -  (with-result (result (list (length B)))
  -    (do-matrix (B (i))
  -      (setf (aref result i) (* A (aref B i))))))
  +  (declare (optimize speed)
  +	   (type double-float A)
  +	   (type (vector double-float) B))
  +  (do-matrix-with-result (B (i) result (list (length B)))
  +    (setf (aref result i) (* A (aref B i)))))
  +
  +(defmethod .* ((A number) (B simple-vector))
  +  (declare (optimize speed)
  +	   (type double-float A)
  +	   (type simple-vector B))
  +  (do-matrix-with-result (B (i) result (list (length B)))
  +    (setf (aref result i) (* A (aref B i)))))
   
   (defmethod .* ((A vector) (B number))
     (.* B A))
   
  +(defgeneric .- (A B))
  +
   (defmethod .- (A B)
  -  (.+ A (.* -1 B)))
  +  (.+ A (.* -1d0 B)))
   
   (defmethod norm ((matrix array))
     (.* (1/ (loop for element across (as-vector matrix)
   	     maximize element))
         matrix))
  +
  +(defgeneric .^ (A B))
   
   (defmethod .^ ((A array) (p integer))
     (assert (squarep A))
  
positive.code: 
  (defun minor-matrix (matrix &rest subscripts)
    ;; Result matrix is of dimensions N-1xN-1 compared to the original matrix
    (with-result (result (mapcar #'1- (array-dimensions matrix)) (array-element-type matrix))
      (let ((row-major 0))
  	 (do-matrix (matrix indices)
  	   ;; Don't copy a number if any of the indices matches a sprcified subscript
  	   (unless (reduce (lambda (a b) (or a b)) (mapcar #'= indices subscripts))
  	     (setf (row-major-aref result row-major) (apply #'aref matrix indices))
  	     (incf row-major))))))
  
  (defun random-matrix (&rest dimensions)
    (with-result (result dimensions)
      (do-matrix (result subscripts)
        (setf (apply #'aref result subscripts) (random 1.0d0)))))
  
  (defun mapply (func data)
    (with-result (result (array-dimensions data))
  
negative.code: 	 sum (aref M i i))))
  
  (defun transpose (matrix)
    (do-matrix-with-result (result (i j) result (nreverse (array-dimensions matrix)))
      (setf (aref result i j) (aref matrix j i))))
  
  (defgeneric .+ (A B)
    (:documentation "Matrix addition function, will also add single numbers elementwise to the matrix."))
  