Extend basic operations

doc/high-level.md

@@ -64,14 +64,20 @@ Magicl provides some "block matrix" constructors: these construct matrices from
 
 ### Basic Operations
 
-| MAGICL     | MATLAB   | NumPy            | Description                 |
-|------------|----------|------------------|-----------------------------|
-| `(@ a b)`  | `a * b`  | `a @ b`          | Matrix multiplication       |
-| `(.+ a b)` | `a + b`  | `a + b`          | Element-wise add            |
-| `(.- a b)` | `a - b`  | `a - b`          | Element-wise subtract       |
-| `(.* a b)` | `a .* b` | `a * b`          | Element-wise multiply       |
-| `(./ a b)` | `a./b`   | `a/b`            | Element-wise divide         |
-| `(.^ a b)` | `a.^b`   | `np.power(a,b)`  | Element-wise exponentiation |
+| MAGICL       | MATLAB     | NumPy             | Description                    |
+|--------------|------------|-------------------|--------------------------------|
+| `(@ a b)`    | `a * b`    | `a @ b`           | Matrix multiplication          |
+| `(.+ a b)`   | `a + b`    | `a + b`           | Element-wise add               |
+| `(.- a b)`   | `a - b`    | `a - b`           | Element-wise subtract          |
+| `(.* a b)`   | `a .* b`   | `a * b`           | Element-wise multiply          |
+| `(./ a b)`   | `a./b`     | `a/b`             | Element-wise divide            |
+| `(.^ a b)`   | `a.^b`     | `np.power(a,b)`   | Element-wise exponentiation    |
+| `(.exp a)`   | `exp(a)`   | `np.exp(a)`       | Element-wise exponential       |
+| `(.log a)`   | `log(a)`   | `np.log(a)`       | Element-wise natural logarithm |
+| `(.max a b)` | `max(a,b)` | `np.maximum(a,b)` | Element-wise maximum           |
+| `(.min a b)` | `min(a,b)` | `np.minimum(a,b)` | Element-wise minimum           |
+
+Note: Elementwise operators with two arguments (e.g. `.+`, `.-`) also act as expected when one argument is a number; for example `(.- A 5.0)` subtracts 5.0 from each element of A.
 
 ### Linear Algebra
 

src/high-level/abstract-tensor.lisp

@@ -201,6 +201,40 @@ If TARGET is not specified then a new tensor is created with the same element ty
                   target dims)))
         target))))
 
+;;; Extend bianary-operator to handle (TENSOR, NUMBER) and (NUMBER, TENSOR)
+;;; arguments Recall that, e.g., TENSOR + NUMBER is commutative, but
+;;; TENSOR - NUMBER isn't, so need two DEFMETHODs here
+;;;
+;;; N.B. This implementation still suffers from the existing issue that if the
+;;; type of the arguments differ then we can get a TYPE-ERROR from the SETF
+;;; Try adding a MATRIX/INT32 to a MATRIX/SINGLE-FLOAT
+
+;;  (TENSOR, NUMBER)
+(defmethod binary-operator ((function function) (source1 abstract-tensor) (source2 number) &optional target)
+  (let ((target (or target (copy-tensor source1))))
+    (map-indexes
+     (shape source1)
+     (lambda (&rest dims)
+       (apply #'(setf tref)
+              (funcall function
+                       (apply #'tref source1 dims)
+                       source2)
+              target dims)))
+    target))
+
+;; (NUMBER, TENSOR)
+(defmethod binary-operator ((function function) (source1 number) (source2 abstract-tensor) &optional target)
+  (let ((target (or target (copy-tensor source2))))
+    (map-indexes
+     (shape source2)
+     (lambda (&rest dims)
+       (apply #'(setf tref)
+              (funcall function
+                       source1
+                       (apply #'tref source2 dims))
+              target dims)))
+    target))
+
 (define-backend-function .+ (source1 source2 &optional target)
   "Add tensors elementwise, optionally storing the result in TARGET.
 If TARGET is not specified then a new tensor is created with the same element type as the first source tensor")
@@ -245,6 +279,47 @@ If TARGET is not specified then a new tensor is created with the same element ty
     ;; choice, like integers.
     (binary-operator #'expt source1 source2 target)))
 
+(define-extensible-function (.max max-lisp) (source1 source2 &optional target)
+  (:documentation "Apply MAX function to tensors elementwise, optionally storing the result in TARGET.
+If TARGET is not specified then a new tensor is created with the same element type as the first source tensor.
+If one argument is a NUMBER then apply MAX function to tensor element and number")
+  (:method (source1 source2 &optional target)
+    (binary-operator #'max source1 source2 target)))
+
+(define-extensible-function (.min min-lisp) (source1 source2 &optional target)
+  (:documentation "Apply MIN function to tensors elementwise, optionally storing the result in TARGET.
+If TARGET is not specified then a new tensor is created with the same element type as the first source tensor.
+If one argument is a NUMBER then apply MIN function to tensor element and number")
+  (:method (source1 source2 &optional target)
+    (binary-operator #'min source1 source2 target)))
+
+
+(defgeneric unary-operator (function source &optional target)
+  (:documentation "Perform a unary operator on tensor elementwise, optionally storing the result in TARGET.
+If TARGET is not specified then a new tensor is created with the same element type as the source tensor")
+  (:method ((function function) (source abstract-tensor) &optional target)
+    (let ((target (or target (copy-tensor source))))
+      (map-indexes
+       (shape source)
+       (lambda (&rest dims)
+         (apply #'(setf tref)
+                (funcall function
+                         (apply #'tref source dims))
+                target dims)))
+      target)))
+
+(define-extensible-function (.exp exp-lisp) (source &optional target)
+  (:documentation "Applies exponential function to tensor elementwise, optionally storing the restult in TARGET.
+If TARGET is not specified then a new tensor is created with the same element type as the source tensor")
+  (:method ((source abstract-tensor) &optional target)
+    (unary-operator #'exp source target)))
+
+(define-extensible-function (.log log-lisp) (source &optional target)
+  (:documentation "Applies natural logarithm to tensor elementwise, optionally storing the restult in TARGET.
+If TARGET is not specified then a new tensor is created with the same element type as the source tensor")
+  (:method ((source abstract-tensor) &optional target)
+    (unary-operator #'log source target)))
+
 
 (define-extensible-function (= =-lisp) (source1 source2 &optional epsilon)
   (:documentation "Check the equality of tensors with an optional EPSILON")

src/packages.lisp

@@ -123,13 +123,18 @@
 
            ;; Operators
            #:binary-operator
+           #:unary-operator
            #:.+
            #:.-
            #:.*
            #:./
            #:.^
            #:=
            #:map
+           #:.exp
+           #:.log
+           #:.max
+           #:.min
 
            ;; Matrix operators
            #:square-matrix-p

tests/abstract-tensor-tests.lisp

@@ -32,3 +32,68 @@
                                   :type 'double-float)))
     (loop :for i :below 15
           :do (= i (apply #'magicl:tref tensor (magicl::from-row-major-index i '(3 5)))))))
+
+(deftest test-tensor-number-ops ()
+  "Test that basic operations on a tensor and number give the expected results"
+  (let* ((input '(-1.0 -0.5 0.5 1.0))
+         (tensor (magicl:from-list input
+                                   '(2 2)
+                                   :type 'single-float)))
+
+    ;; .+
+    (is (magicl:= (magicl:.+ tensor 3.14)
+                  (magicl:from-list (loop for i in input collect (+ i 3.14))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))
+
+    ;; .-
+    (is (magicl:= (magicl:.- 3.14 tensor)
+                  (magicl:from-list (loop for i in input collect (- 3.14 i))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))
+
+    ;; .*
+    (is (magicl:= (magicl:.* 3.14 tensor)
+                  (magicl:from-list (loop for i in input collect (* 3.14 i))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))
+
+    ;; ./
+    (is (magicl:= (magicl:./ 3.14 tensor)
+                  (magicl:from-list (loop for i in input collect (/ 3.14 i))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))
+
+    ;; .max
+    (is (magicl:= (magicl:.max tensor 0.0)
+                  (magicl:from-list (loop for i in input collect (max 0.0 i))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))
+    ;; .min
+    (is (magicl:= (magicl:.min 0.0 tensor)
+                  (magicl:from-list (loop for i in input collect (min i 0.0))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))))
+
+(deftest test-unary-ops ()
+  "Test that basic unary operations on a tensor give the expected results"
+  (let* ((input '(-1.1 -0.4 0.5 1.0))
+         (tensor (magicl:from-list input
+                                   '(2 2)
+                                   :type 'single-float)))
+
+    ;; .exp
+    (is (magicl:= (magicl:.exp tensor)
+                  (magicl:from-list (loop for i in input collect (exp i))
+                                    (magicl:shape tensor)
+                                    :type 'single-float)))
+
+    ;; .log - Recall natural log of negative numbers are undefined!
+    (let* ((input '(0.1 0.5 1.1 2.0))
+           (tensor (magicl:from-list input
+                                     '(2 2)
+                                     :type 'single-float)))
+      (is (magicl:= (magicl:.log tensor)
+                    (magicl:from-list (loop for i in input collect (log i))
+                                      (magicl:shape tensor)
+                                      :type 'single-float))))))