Update copyright year.
[bpt/emacs.git] / lisp / calc / calc-arith.el
index 3a436cb..2372b0e 100644 (file)
@@ -1,6 +1,7 @@
 ;;; calc-arith.el --- arithmetic functions for Calc
 
-;; Copyright (C) 1990, 1991, 1992, 1993, 2001 Free Software Foundation, Inc.
+;; Copyright (C) 1990, 1991, 1992, 1993, 2001, 2002, 2003, 2004,
+;;   2005, 2006 Free Software Foundation, Inc.
 
 ;; Author: David Gillespie <daveg@synaptics.com>
 ;; Maintainer: Jay Belanger <belanger@truman.edu>
                                       calcFunc-max calcFunc-min))
 
 (defvar math-real-if-arg-functions '(calcFunc-sin calcFunc-cos
-                                    calcFunc-tan calcFunc-arctan
+                                    calcFunc-tan calcFunc-sec
+                                     calcFunc-csc calcFunc-cot
+                                     calcFunc-arctan
                                     calcFunc-sinh calcFunc-cosh
-                                    calcFunc-tanh calcFunc-exp
+                                    calcFunc-tanh calcFunc-sech
+                                     calcFunc-csch calcFunc-coth
+                                     calcFunc-exp
                                     calcFunc-gamma calcFunc-fact))
 
 (defvar math-integer-functions '(calcFunc-idiv
     (real number)
     (number)
     (scalar)
+    (sqmatrix matrix vector)
     (matrix vector)
     (vector)
     (const)))
   (and (not (Math-scalarp a))
        (not (math-known-scalarp a t))))
 
+(defun math-known-square-matrixp (a)
+  (and (math-known-matrixp a)
+       (math-check-known-square-matrixp a)))
+
 ;;; Try to prove that A is a scalar (i.e., a non-vector).
 (defun math-check-known-scalarp (a)
   (cond ((Math-objectp a) t)
         (let ((decl (if (eq (car a) 'var)
                         (or (assq (nth 2 a) math-decls-cache)
                             math-decls-all)
-                      (assq (car a) math-decls-cache))))
-          (memq 'scalar (nth 1 decl))))))
+                      (assq (car a) math-decls-cache)))
+               val)
+           (cond
+            ((memq 'scalar (nth 1 decl))
+             t)
+            ((and (eq (car a) 'var)
+                  (boundp (nth 2 a))
+                  (setq val (symbol-value (nth 2 a))))
+             (math-check-known-scalarp val))
+            (t
+             nil))))))
 
 ;;; Try to prove that A is *not* a scalar.
 (defun math-check-known-matrixp (a)
         (let ((decl (if (eq (car a) 'var)
                         (or (assq (nth 2 a) math-decls-cache)
                             math-decls-all)
-                      (assq (car a) math-decls-cache))))
-          (memq 'vector (nth 1 decl))))))
-
+                      (assq (car a) math-decls-cache)))
+               val)
+           (cond
+            ((memq 'matrix (nth 1 decl))
+             t)
+            ((and (eq (car a) 'var)
+                  (boundp (nth 2 a))
+                  (setq val (symbol-value (nth 2 a))))
+             (math-check-known-matrixp val))
+            (t
+             nil))))))
+
+;;; Given that A is a matrix, try to prove that it is a square matrix.
+(defun math-check-known-square-matrixp (a)
+  (cond ((math-square-matrixp a)
+         t)
+        ((eq (car-safe a) '^)
+         (math-check-known-square-matrixp (nth 1 a)))
+        ((or
+          (eq (car-safe a) '*)
+          (eq (car-safe a) '+)
+          (eq (car-safe a) '-))
+         (and
+          (math-check-known-square-matrixp (nth 1 a))
+          (math-check-known-square-matrixp (nth 2 a))))
+        (t
+         (let ((decl (if (eq (car a) 'var)
+                         (or (assq (nth 2 a) math-decls-cache)
+                             math-decls-all)
+                       (assq (car a) math-decls-cache)))
+               val)
+           (cond
+            ((memq 'sqmatrix (nth 1 decl))
+             t)
+            ((and (eq (car a) 'var)
+                  (boundp (nth 2 a))
+                  (setq val (symbol-value (nth 2 a))))
+             (math-check-known-square-matrixp val))
+            ((and (or
+                   (integerp calc-matrix-mode)
+                   (eq calc-matrix-mode 'sqmatrix))
+                  (eq (car-safe a) 'var))
+             t)
+            ((memq 'matrix (nth 1 decl))
+             nil)
+            (t
+             nil))))))
 
 ;;; Try to prove that A is a real (i.e., not complex).
 (defun math-known-realp (a)
               (and (math-known-scalarp b)
                    (math-add (nth 1 a) b))))
       (and (eq (car-safe b) 'calcFunc-idn)
-          (= (length a) 2)
+          (= (length b) 2)
           (or (and (math-square-matrixp a)
                    (math-add a (math-mimic-ident (nth 1 b) a)))
               (and (math-known-scalarp a)
       (and (eq (car-safe b) '^)
           (Math-looks-negp (nth 2 b))
           (not (and (eq (car-safe a) '^) (Math-looks-negp (nth 2 a))))
+           (not (math-known-matrixp (nth 1 b)))
           (math-div a (math-normalize
                        (list '^ (nth 1 b) (math-neg (nth 2 b))))))
       (and (eq (car-safe a) '/)
                    (list 'calcFunc-idn (math-mul a (nth 1 b))))
               (and (math-known-matrixp a)
                    (math-mul a (nth 1 b)))))
+      (and (math-identity-matrix-p a t)
+           (or (and (eq (car-safe b) 'calcFunc-idn)
+                    (= (length b) 2)
+                    (list 'calcFunc-idn (math-mul 
+                                         (nth 1 (nth 1 a))
+                                         (nth 1 b))
+                          (1- (length a))))
+               (and (math-known-scalarp b)
+                    (list 'calcFunc-idn (math-mul 
+                                         (nth 1 (nth 1 a)) b)
+                          (1- (length a))))
+               (and (math-known-matrixp b)
+                    (math-mul (nth 1 (nth 1 a)) b))))
+      (and (math-identity-matrix-p b t)
+           (or (and (eq (car-safe a) 'calcFunc-idn)
+                    (= (length a) 2)
+                    (list 'calcFunc-idn (math-mul (nth 1 a) 
+                                                  (nth 1 (nth 1 b)))
+                          (1- (length b))))
+               (and (math-known-scalarp a)
+                    (list 'calcFunc-idn (math-mul a (nth 1 (nth 1 b))) 
+                          (1- (length b))))
+               (and (math-known-matrixp a)
+                    (math-mul a (nth 1 (nth 1 b))))))
       (and (math-looks-negp b)
           (math-mul (math-neg a) (math-neg b)))
       (and (eq (car-safe b) '-)
            (math-reject-arg b "*Division by zero"))
        a))))
 
+;; For math-div-symb-fancy
+(defvar math-trig-inverses
+  '((calcFunc-sin . calcFunc-csc)
+    (calcFunc-cos . calcFunc-sec)
+    (calcFunc-tan . calcFunc-cot)
+    (calcFunc-sec . calcFunc-cos)
+    (calcFunc-csc . calcFunc-sin)
+    (calcFunc-cot . calcFunc-tan)
+    (calcFunc-sinh . calcFunc-csch)
+    (calcFunc-cosh . calcFunc-sech)
+    (calcFunc-tanh . calcFunc-coth)
+    (calcFunc-sech . calcFunc-cosh)
+    (calcFunc-csch . calcFunc-sinh)
+    (calcFunc-coth . calcFunc-tanh)))
+
+(defvar math-div-trig)
+(defvar math-div-non-trig)
+
+(defun math-div-new-trig (tr)
+  (if math-div-trig
+      (setq math-div-trig
+            (list '* tr math-div-trig))
+    (setq math-div-trig tr)))
+
+(defun math-div-new-non-trig (ntr)
+  (if math-div-non-trig
+      (setq math-div-non-trig 
+            (list '* ntr math-div-non-trig))
+    (setq math-div-non-trig ntr)))
+
+(defun math-div-isolate-trig (expr)
+  (if (eq (car-safe expr) '*)
+      (progn
+        (math-div-isolate-trig-term (nth 1 expr))
+        (math-div-isolate-trig (nth 2 expr)))
+    (math-div-isolate-trig-term expr)))
+
+(defun math-div-isolate-trig-term (term)
+  (let ((fn (assoc (car-safe term) math-trig-inverses)))
+    (if fn
+        (math-div-new-trig
+         (cons (cdr fn) (cdr term)))
+      (math-div-new-non-trig term))))
+
 (defun math-div-symb-fancy (a b)
-  (or (and math-simplify-only
+  (or (and (math-known-matrixp b)
+           (math-mul a (math-pow b -1)))
+      (and math-simplify-only
           (not (equal a math-simplify-only))
           (list '/ a b))
       (and (Math-equal-int b 1) a)
                    (list 'calcFunc-idn (math-div a (nth 1 b))))
               (and (math-known-matrixp a)
                    (math-div a (nth 1 b)))))
+      (and math-simplifying
+           (let ((math-div-trig nil)
+                 (math-div-non-trig nil))
+             (math-div-isolate-trig b)
+             (if math-div-trig
+                 (if math-div-non-trig
+                     (math-div (math-mul a math-div-trig) math-div-non-trig)
+                   (math-mul a math-div-trig))
+               nil)))
       (if (and calc-matrix-mode
               (or (math-known-matrixp a) (math-known-matrixp b)))
          (math-combine-prod a b nil t nil)
               (math-mul-zero b a))))
       (list '/ a b)))
 
+;;; Division from the left.
+(defun calcFunc-ldiv (a b)
+  (if (math-known-scalarp a)
+      (math-div b a)
+    (math-mul (math-pow a -1) b)))
 
 (defun calcFunc-mod (a b)
   (math-normalize (list '% a b)))
           (cond ((and math-simplify-only
                       (not (equal a math-simplify-only)))
                  (list '^ a b))
+                 ((and (eq (car-safe a) '*)
+                       (or 
+                        (and
+                         (math-known-matrixp (nth 1 a))
+                         (math-known-matrixp (nth 2 a)))
+                        (and
+                         calc-matrix-mode
+                         (not (eq calc-matrix-mode 'scalar))
+                         (and (not (math-known-scalarp (nth 1 a)))
+                              (not (math-known-scalarp (nth 2 a)))))))
+                  (if (and (= b -1)
+                           (math-known-square-matrixp (nth 1 a))
+                           (math-known-square-matrixp (nth 2 a)))
+                      (math-mul (math-pow-fancy (nth 2 a) -1) 
+                                (math-pow-fancy (nth 1 a) -1))
+                    (list '^ a b)))
                 ((and (eq (car-safe a) '*)
                       (or (math-known-num-integerp b)
                           (math-known-nonnegp (nth 1 a))
         invb
         (math-looks-negp (nth 2 b)))
     (math-mul a (math-pow (nth 1 b) (math-neg (nth 2 b)))))
+   ((and math-simplifying
+         (math-combine-prod-trig a b)))
    (t (let ((apow 1) (bpow 1))
        (and (consp a)
             (cond ((and (eq (car a) '^)
                            (math-pow a apow)
                          (inexact-result (list '^ a apow)))))))))))
 
+(defun math-combine-prod-trig (a b)
+  (cond
+   ((and (eq (car-safe a) 'calcFunc-sin)
+         (eq (car-safe b) 'calcFunc-csc)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    1)
+   ((and (eq (car-safe a) 'calcFunc-sin)
+         (eq (car-safe b) 'calcFunc-sec)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-tan (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-sin)
+         (eq (car-safe b) 'calcFunc-cot)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-cos (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-cos)
+         (eq (car-safe b) 'calcFunc-sec)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    1)
+   ((and (eq (car-safe a) 'calcFunc-cos)
+         (eq (car-safe b) 'calcFunc-csc)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-cot (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-cos)
+         (eq (car-safe b) 'calcFunc-tan)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-sin (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-tan)
+         (eq (car-safe b) 'calcFunc-cot)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    1)
+   ((and (eq (car-safe a) 'calcFunc-tan)
+         (eq (car-safe b) 'calcFunc-csc)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-sec (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-sec)
+         (eq (car-safe b) 'calcFunc-cot)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-csc (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-sinh)
+         (eq (car-safe b) 'calcFunc-csch)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    1)
+   ((and (eq (car-safe a) 'calcFunc-sinh)
+         (eq (car-safe b) 'calcFunc-sech)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-tanh (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-sinh)
+         (eq (car-safe b) 'calcFunc-coth)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-cosh (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-cosh)
+         (eq (car-safe b) 'calcFunc-sech)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    1)
+   ((and (eq (car-safe a) 'calcFunc-cosh)
+         (eq (car-safe b) 'calcFunc-csch)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-coth (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-cosh)
+         (eq (car-safe b) 'calcFunc-tanh)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-sinh (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-tanh)
+         (eq (car-safe b) 'calcFunc-coth)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    1)
+   ((and (eq (car-safe a) 'calcFunc-tanh)
+         (eq (car-safe b) 'calcFunc-csch)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-sech (cdr a)))
+   ((and (eq (car-safe a) 'calcFunc-sech)
+         (eq (car-safe b) 'calcFunc-coth)
+         (= 0 (math-simplify (math-sub (cdr a) (cdr b)))))
+    (cons 'calcFunc-csch (cdr a)))
+   (t
+    nil)))
+
 (defun math-mul-or-div (a b ainv binv)
   (if (or (Math-vectorp a) (Math-vectorp b))
       (math-normalize