cgit.sxemacs.org Git - packages/blob - xemacs-packages/calc/calc-vec.el

   1 ;; Calculator for GNU Emacs, part II [calc-vec.el]
   2 ;; Copyright (C) 1990, 1991, 1992, 1993 Free Software Foundation, Inc.
   3 ;; Written by Dave Gillespie, daveg@synaptics.com.
   4
   5 ;; This file is part of GNU Emacs.
   6
   7 ;; GNU Emacs is distributed in the hope that it will be useful,
   8 ;; but WITHOUT ANY WARRANTY.  No author or distributor
   9 ;; accepts responsibility to anyone for the consequences of using it
  10 ;; or for whether it serves any particular purpose or works at all,
  11 ;; unless he says so in writing.  Refer to the GNU Emacs General Public
  12 ;; License for full details.
  13
  14 ;; Everyone is granted permission to copy, modify and redistribute
  15 ;; GNU Emacs, but only under the conditions described in the
  16 ;; GNU Emacs General Public License.   A copy of this license is
  17 ;; supposed to have been given to you along with GNU Emacs so you
  18 ;; can know your rights and responsibilities.  It should be in a
  19 ;; file named COPYING.  Among other things, the copyright notice
  20 ;; and this notice must be preserved on all copies.
  21
  22
  23
  24 ;; This file is autoloaded from calc-ext.el.
  25 (require 'calc-ext)
  26
  27 (require 'calc-macs)
  28
  29 (defun calc-Need-calc-vec () nil)
  30
  31
  32 (defun calc-display-strings (n)
  33   (interactive "P")
  34   (calc-wrapper
  35    (message (if (calc-change-mode 'calc-display-strings n t t)
  36                 "Displaying vectors of integers as quoted strings."
  37               "Displaying vectors of integers normally.")))
  38 )
  39
  40
  41 (defun calc-pack (n)
  42   (interactive "P")
  43   (calc-wrapper
  44    (let* ((nn (if n 1 2))
  45           (mode (if n (prefix-numeric-value n) (calc-top-n 1)))
  46           (mode (if (and (Math-vectorp mode) (cdr mode)) (cdr mode)
  47                   (if (integerp mode) mode
  48                     (error "Packing mode must be an integer or vector of integers"))))
  49           (num (calc-pack-size mode))
  50           (items (calc-top-list num nn)))
  51      (calc-enter-result (+ nn num -1) "pack" (calc-pack-items mode items))))
  52 )
  53
  54 (defun calc-pack-size (mode)
  55   (cond ((consp mode)
  56          (let ((size 1))
  57            (while mode
  58              (or (integerp (car mode)) (error "Vector of integers expected"))
  59              (setq size (* size (calc-pack-size (car mode)))
  60                    mode (cdr mode)))
  61            (if (= size 0)
  62                (error "Zero dimensions not allowed")
  63              size)))
  64         ((>= mode 0) mode)
  65         (t (or (cdr (assq mode '((-3 . 3) (-13 . 1) (-14 . 3) (-15 . 6))))
  66                2)))
  67 )
  68
  69 (defun calc-pack-items (mode items)
  70   (cond ((consp mode)
  71          (if (cdr mode)
  72              (let* ((size (calc-pack-size (cdr mode)))
  73                     (len (length items))
  74                     (new nil)
  75                     p row)
  76                (while (> len 0)
  77                  (setq p (nthcdr (1- size) items)
  78                        row items
  79                        items (cdr p)
  80                        len (- len size))
  81                  (setcdr p nil)
  82                  (setq new (cons (calc-pack-items (cdr mode) row) new)))
  83                (calc-pack-items (car mode) (nreverse new)))
  84            (calc-pack-items (car mode) items)))
  85         ((>= mode 0)
  86          (cons 'vec items))
  87         ((= mode -3)
  88          (if (and (math-objvecp (car items))
  89                   (math-objvecp (nth 1 items))
  90                   (math-objvecp (nth 2 items)))
  91              (if (and (math-num-integerp (car items))
  92                       (math-num-integerp (nth 1 items)))
  93                  (if (math-realp (nth 2 items))
  94                      (cons 'hms items)
  95                    (error "Seconds must be real"))
  96                (error "Hours and minutes must be integers"))
  97            (math-normalize (list '+
  98                                  (list '+
  99                                        (if (eq calc-angle-mode 'rad)
 100                                            (list '* (car items)
 101                                                  '(hms 1 0 0))
 102                                          (car items))
 103                                        (list '* (nth 1 items) '(hms 0 1 0)))
 104                                  (list '* (nth 2 items) '(hms 0 0 1))))))
 105         ((= mode -13)
 106          (if (math-realp (car items))
 107              (cons 'date items)
 108            (if (eq (car-safe (car items)) 'date)
 109                (car items)
 110              (if (math-objvecp (car items))
 111                  (error "Date value must be real")
 112                (cons 'calcFunc-date items)))))
 113         ((memq mode '(-14 -15))
 114          (let ((p items))
 115            (while (and p (math-objvecp (car p)))
 116              (or (math-integerp (car p))
 117                  (error "Components must be integers"))
 118              (setq p (cdr p)))
 119            (if p
 120                (cons 'calcFunc-date items)
 121              (list 'date (math-dt-to-date items)))))
 122         ((or (eq (car-safe (car items)) 'vec)
 123              (eq (car-safe (nth 1 items)) 'vec))
 124          (let* ((x (car items))
 125                 (vx (eq (car-safe x) 'vec))
 126                 (y (nth 1 items))
 127                 (vy (eq (car-safe y) 'vec))
 128                 (z nil)
 129                 (n (1- (length (if vx x y)))))
 130            (and vx vy
 131                 (/= n (1- (length y)))
 132                 (error "Vectors must be the same length"))
 133            (while (>= (setq n (1- n)) 0)
 134              (setq z (cons (calc-pack-items
 135                             mode
 136                             (list (if vx (car (setq x (cdr x))) x)
 137                                   (if vy (car (setq y (cdr y))) y)))
 138                            z)))
 139            (cons 'vec (nreverse z))))
 140         ((= mode -1)
 141          (if (and (math-realp (car items)) (math-realp (nth 1 items)))
 142              (cons 'cplx items)
 143            (if (and (math-objectp (car items)) (math-objectp (nth 1 items)))
 144                (error "Components must be real"))
 145            (math-normalize (list '+ (car items)
 146                                  (list '* (nth 1 items) '(cplx 0 1))))))
 147         ((= mode -2)
 148          (if (and (math-realp (car items)) (math-anglep (nth 1 items)))
 149              (cons 'polar items)
 150            (if (and (math-objectp (car items)) (math-objectp (nth 1 items)))
 151                (error "Components must be real"))
 152            (math-normalize (list '* (car items)
 153                                  (if (math-anglep (nth 1 items))
 154                                      (list 'polar 1 (nth 1 items))
 155                                    (list 'calcFunc-exp
 156                                          (list '*
 157                                                (math-to-radians-2
 158                                                 (nth 1 items))
 159                                                (list 'polar
 160                                                      1
 161                                                      (math-quarter-circle
 162                                                       nil)))))))))
 163         ((= mode -4)
 164          (let ((x (car items))
 165                (sigma (nth 1 items)))
 166            (if (or (math-scalarp x) (not (math-objvecp x)))
 167                (if (or (math-anglep sigma) (not (math-objvecp sigma)))
 168                    (math-make-sdev x sigma)
 169                  (error "Error component must be real"))
 170              (error "Mean component must be real or complex"))))
 171         ((= mode -5)
 172          (let ((a (car items))
 173                (m (nth 1 items)))
 174            (if (and (math-anglep a) (math-anglep m))
 175                (if (math-posp m)
 176                    (math-make-mod a m)
 177                  (error "Modulus must be positive"))
 178              (if (and (math-objectp a) (math-objectp m))
 179                  (error "Components must be real"))
 180              (list 'calcFunc-makemod a m))))
 181         ((memq mode '(-6 -7 -8 -9))
 182          (let ((lo (car items))
 183                (hi (nth 1 items)))
 184            (if (and (or (math-anglep lo) (eq (car lo) 'date)
 185                         (not (math-objvecp lo)))
 186                     (or (math-anglep hi) (eq (car hi) 'date)
 187                         (not (math-objvecp hi))))
 188                (math-make-intv (+ mode 9) lo hi)
 189              (error "Components must be real"))))
 190         ((eq mode -10)
 191          (if (math-zerop (nth 1 items))
 192              (error "Denominator must not be zero")
 193            (if (and (math-integerp (car items)) (math-integerp (nth 1 items)))
 194                (math-normalize (cons 'frac items))
 195              (if (and (math-objectp (car items)) (math-objectp (nth 1 items)))
 196                  (error "Components must be integers"))
 197              (cons 'calcFunc-fdiv items))))
 198         ((memq mode '(-11 -12))
 199          (if (and (math-realp (car items)) (math-integerp (nth 1 items)))
 200              (calcFunc-scf (math-float (car items)) (nth 1 items))
 201            (if (and (math-objectp (car items)) (math-objectp (nth 1 items)))
 202                (error "Components must be integers"))
 203            (math-normalize
 204             (list 'calcFunc-scf
 205                   (list 'calcFunc-float (car items))
 206                   (nth 1 items)))))
 207         (t
 208          (error "Invalid packing mode: %d" mode)))
 209 )
 210
 211 (defun calc-unpack (mode)
 212   (interactive "P")
 213   (calc-wrapper
 214    (let ((calc-unpack-with-type t))
 215      (calc-pop-push-record-list 1 "unpk" (calc-unpack-item
 216                                           (and mode
 217                                                (prefix-numeric-value mode))
 218                                           (calc-top)))))
 219 )
 220
 221 (defun calc-unpack-type (item)
 222   (cond ((eq (car-safe item) 'vec)
 223          (1- (length item)))
 224         ((eq (car-safe item) 'intv)
 225          (- (nth 1 item) 9))
 226         (t
 227          (or (cdr (assq (car-safe item) '( (cplx . -1) (polar . -2)
 228                                            (hms . -3) (sdev . -4) (mod . -5)
 229                                            (frac . -10) (float . -11)
 230                                            (date . -13) )))
 231              (error "Argument must be a composite object"))))
 232 )
 233
 234 (defun calc-unpack-item (mode item)
 235   (cond ((not mode)
 236          (if (or (and (not (memq (car-safe item) '(frac float cplx polar vec
 237                                                         hms date sdev mod
 238                                                         intv)))
 239                       (math-objvecp item))
 240                  (eq (car-safe item) 'var))
 241              (error "Argument must be a composite object or function call"))
 242          (if (eq (car item) 'intv)
 243              (cdr (cdr item))
 244            (cdr item)))
 245         ((> mode 0)
 246          (let ((dims nil)
 247                type new row)
 248            (setq item (list item))
 249            (while (> mode 0)
 250              (setq type (calc-unpack-type (car item))
 251                    dims (cons type dims)
 252                    new (calc-unpack-item nil (car item)))
 253              (while (setq item (cdr item))
 254                (or (= (calc-unpack-type (car item)) type)
 255                    (error "Inconsistent types or dimensions in vector elements"))
 256                (setq new (append new (calc-unpack-item nil (car item)))))
 257              (setq item new
 258                    mode (1- mode)))
 259            (if (cdr dims) (setq dims (list (cons 'vec (nreverse dims)))))
 260            (cond ((eq calc-unpack-with-type 'pair)
 261                   (list (car dims) (cons 'vec item)))
 262                  (calc-unpack-with-type
 263                   (append item dims))
 264                  (t item))))
 265         ((eq calc-unpack-with-type 'pair)
 266          (let ((calc-unpack-with-type nil))
 267            (list mode (cons 'vec (calc-unpack-item mode item)))))
 268         ((= mode -3)
 269          (if (eq (car-safe item) 'hms)
 270              (cdr item)
 271            (error "Argument must be an HMS form")))
 272         ((= mode -13)
 273          (if (eq (car-safe item) 'date)
 274              (cdr item)
 275            (error "Argument must be a date form")))
 276         ((= mode -14)
 277          (if (eq (car-safe item) 'date)
 278              (math-date-to-dt (math-floor (nth 1 item)))
 279            (error "Argument must be a date form")))
 280         ((= mode -15)
 281          (if (eq (car-safe item) 'date)
 282              (append (math-date-to-dt (nth 1 item))
 283                      (and (not (math-integerp (nth 1 item)))
 284                           (list 0 0 0)))
 285            (error "Argument must be a date form")))
 286         ((eq (car-safe item) 'vec)
 287          (let ((x nil)
 288                (y nil)
 289                res)
 290            (while (setq item (cdr item))
 291              (setq res (calc-unpack-item mode (car item))
 292                    x (cons (car res) x)
 293                    y (cons (nth 1 res) y)))
 294            (list (cons 'vec (nreverse x))
 295                  (cons 'vec (nreverse y)))))
 296         ((= mode -1)
 297          (if (eq (car-safe item) 'cplx)
 298              (cdr item)
 299            (if (eq (car-safe item) 'polar)
 300                (cdr (math-complex item))
 301              (if (Math-realp item)
 302                  (list item 0)
 303                (error "Argument must be a complex number")))))
 304         ((= mode -2)
 305          (if (or (memq (car-safe item) '(cplx polar))
 306                  (Math-realp item))
 307              (cdr (math-polar item))
 308            (error "Argument must be a complex number")))
 309         ((= mode -4)
 310          (if (eq (car-safe item) 'sdev)
 311              (cdr item)
 312            (list item 0)))
 313         ((= mode -5)
 314          (if (eq (car-safe item) 'mod)
 315              (cdr item)
 316            (error "Argument must be a modulo form")))
 317         ((memq mode '(-6 -7 -8 -9))
 318          (if (eq (car-safe item) 'intv)
 319              (cdr (cdr item))
 320            (list item item)))
 321         ((= mode -10)
 322          (if (eq (car-safe item) 'frac)
 323              (cdr item)
 324            (if (Math-integerp item)
 325                (list item 1)
 326              (error "Argument must be a rational number"))))
 327         ((= mode -11)
 328          (if (eq (car-safe item) 'float)
 329              (list (nth 1 item) (math-normalize (nth 2 item)))
 330            (error "Expected a floating-point number")))
 331         ((= mode -12)
 332          (if (eq (car-safe item) 'float)
 333              (list (calcFunc-mant item) (calcFunc-xpon item))
 334            (error "Expected a floating-point number")))
 335         (t
 336          (error "Invalid unpacking mode: %d" mode)))
 337 )
 338 (setq calc-unpack-with-type nil)
 339
 340 (defun calc-diag (n)
 341   (interactive "P")
 342   (calc-wrapper
 343    (calc-enter-result 1 "diag" (if n
 344                                    (list 'calcFunc-diag (calc-top-n 1)
 345                                          (prefix-numeric-value n))
 346                                  (list 'calcFunc-diag (calc-top-n 1)))))
 347 )
 348
 349 (defun calc-ident (n)
 350   (interactive "NDimension of identity matrix = ")
 351   (calc-wrapper
 352    (calc-enter-result 0 "idn" (if (eq n 0)
 353                                   '(calcFunc-idn 1)
 354                                 (list 'calcFunc-idn 1
 355                                       (prefix-numeric-value n)))))
 356 )
 357
 358 (defun calc-index (n &optional stack)
 359   (interactive "NSize of vector = \nP")
 360   (calc-wrapper
 361    (if (consp stack)
 362        (calc-enter-result 3 "indx" (cons 'calcFunc-index (calc-top-list-n 3)))
 363      (calc-enter-result 0 "indx" (list 'calcFunc-index
 364                                        (prefix-numeric-value n)))))
 365 )
 366
 367 (defun calc-build-vector (n)
 368   (interactive "NSize of vector = ")
 369   (calc-wrapper
 370    (calc-enter-result 1 "bldv" (list 'calcFunc-cvec
 371                                      (calc-top-n 1)
 372                                      (prefix-numeric-value n))))
 373 )
 374
 375 (defun calc-cons (arg)
 376   (interactive "P")
 377   (calc-wrapper
 378    (if (calc-is-hyperbolic)
 379        (calc-binary-op "rcns" 'calcFunc-rcons arg)
 380      (calc-binary-op "cons" 'calcFunc-cons arg)))
 381 )
 382
 383
 384 (defun calc-head (arg)
 385   (interactive "P")
 386   (calc-wrapper
 387    (if (calc-is-inverse)
 388        (if (calc-is-hyperbolic)
 389            (calc-unary-op "rtai" 'calcFunc-rtail arg)
 390          (calc-unary-op "tail" 'calcFunc-tail arg))
 391      (if (calc-is-hyperbolic)
 392          (calc-unary-op "rhed" 'calcFunc-rhead arg)
 393        (calc-unary-op "head" 'calcFunc-head arg))))
 394 )
 395
 396 (defun calc-tail (arg)
 397   (interactive "P")
 398   (calc-invert-func)
 399   (calc-head arg)
 400 )
 401
 402 (defun calc-vlength (arg)
 403   (interactive "P")
 404   (calc-wrapper
 405    (if (calc-is-hyperbolic)
 406        (calc-unary-op "dims" 'calcFunc-mdims arg)
 407      (calc-unary-op "len" 'calcFunc-vlen arg)))
 408 )
 409
 410 (defun calc-arrange-vector (n)
 411   (interactive "NNumber of columns = ")
 412   (calc-wrapper
 413    (calc-enter-result 1 "arng" (list 'calcFunc-arrange (calc-top-n 1)
 414                                      (prefix-numeric-value n))))
 415 )
 416
 417 (defun calc-vector-find (arg)
 418   (interactive "P")
 419   (calc-wrapper
 420    (let ((func (cons 'calcFunc-find (calc-top-list-n 2))))
 421      (calc-enter-result
 422       2 "find"
 423       (if arg (append func (list (prefix-numeric-value arg))) func))))
 424 )
 425
 426 (defun calc-subvector ()
 427   (interactive)
 428   (calc-wrapper
 429    (if (calc-is-inverse)
 430        (calc-enter-result 3 "rsvc" (cons 'calcFunc-rsubvec
 431                                          (calc-top-list-n 3)))
 432      (calc-enter-result 3 "svec" (cons 'calcFunc-subvec (calc-top-list-n 3)))))
 433 )
 434
 435 (defun calc-reverse-vector (arg)
 436   (interactive "P")
 437   (calc-wrapper
 438    (calc-unary-op "rev" 'calcFunc-rev arg))
 439 )
 440
 441 (defun calc-mask-vector (arg)
 442   (interactive "P")
 443   (calc-wrapper
 444    (calc-binary-op "vmsk" 'calcFunc-vmask arg))
 445 )
 446
 447 (defun calc-expand-vector (arg)
 448   (interactive "P")
 449   (calc-wrapper
 450    (if (calc-is-hyperbolic)
 451        (calc-enter-result 3 "vexp" (cons 'calcFunc-vexp (calc-top-list-n 3)))
 452      (calc-binary-op "vexp" 'calcFunc-vexp arg)))
 453 )
 454
 455 (defun calc-sort ()
 456   (interactive)
 457   (calc-slow-wrapper
 458    (if (calc-is-inverse)
 459        (calc-enter-result 1 "rsrt" (list 'calcFunc-rsort (calc-top-n 1)))
 460      (calc-enter-result 1 "sort" (list 'calcFunc-sort (calc-top-n 1)))))
 461 )
 462
 463 (defun calc-grade ()
 464   (interactive)
 465   (calc-slow-wrapper
 466    (if (calc-is-inverse)
 467        (calc-enter-result 1 "rgrd" (list 'calcFunc-rgrade (calc-top-n 1)))
 468      (calc-enter-result 1 "grad" (list 'calcFunc-grade (calc-top-n 1)))))
 469 )
 470
 471 (defun calc-histogram (n)
 472   (interactive "NNumber of bins: ")
 473   (calc-slow-wrapper
 474    (if calc-hyperbolic-flag
 475        (calc-enter-result 2 "hist" (list 'calcFunc-histogram
 476                                          (calc-top-n 2)
 477                                          (calc-top-n 1)
 478                                          (prefix-numeric-value n)))
 479      (calc-enter-result 1 "hist" (list 'calcFunc-histogram
 480                                        (calc-top-n 1)
 481                                        (prefix-numeric-value n)))))
 482 )
 483
 484 (defun calc-transpose (arg)
 485   (interactive "P")
 486   (calc-wrapper
 487    (calc-unary-op "trn" 'calcFunc-trn arg))
 488 )
 489
 490 (defun calc-conj-transpose (arg)
 491   (interactive "P")
 492   (calc-wrapper
 493    (calc-unary-op "ctrn" 'calcFunc-ctrn arg))
 494 )
 495
 496 (defun calc-cross (arg)
 497   (interactive "P")
 498   (calc-wrapper
 499    (calc-binary-op "cros" 'calcFunc-cross arg))
 500 )
 501
 502 (defun calc-remove-duplicates (arg)
 503   (interactive "P")
 504   (calc-wrapper
 505    (calc-unary-op "rdup" 'calcFunc-rdup arg))
 506 )
 507
 508 (defun calc-set-union (arg)
 509   (interactive "P")
 510   (calc-wrapper
 511    (calc-binary-op "unio" 'calcFunc-vunion arg '(vec) 'calcFunc-rdup))
 512 )
 513
 514 (defun calc-set-intersect (arg)
 515   (interactive "P")
 516   (calc-wrapper
 517    (calc-binary-op "intr" 'calcFunc-vint arg '(vec) 'calcFunc-rdup))
 518 )
 519
 520 (defun calc-set-difference (arg)
 521   (interactive "P")
 522   (calc-wrapper
 523    (calc-binary-op "diff" 'calcFunc-vdiff arg '(vec) 'calcFunc-rdup))
 524 )
 525
 526 (defun calc-set-xor (arg)
 527   (interactive "P")
 528   (calc-wrapper
 529    (calc-binary-op "xor" 'calcFunc-vxor arg '(vec) 'calcFunc-rdup))
 530 )
 531
 532 (defun calc-set-complement (arg)
 533   (interactive "P")
 534   (calc-wrapper
 535    (calc-unary-op "cmpl" 'calcFunc-vcompl arg))
 536 )
 537
 538 (defun calc-set-floor (arg)
 539   (interactive "P")
 540   (calc-wrapper
 541    (calc-unary-op "vflr" 'calcFunc-vfloor arg))
 542 )
 543
 544 (defun calc-set-enumerate (arg)
 545   (interactive "P")
 546   (calc-wrapper
 547    (calc-unary-op "enum" 'calcFunc-venum arg))
 548 )
 549
 550 (defun calc-set-span (arg)
 551   (interactive "P")
 552   (calc-wrapper
 553    (calc-unary-op "span" 'calcFunc-vspan arg))
 554 )
 555
 556 (defun calc-set-cardinality (arg)
 557   (interactive "P")
 558   (calc-wrapper
 559    (calc-unary-op "card" 'calcFunc-vcard arg))
 560 )
 561
 562 (defun calc-unpack-bits (arg)
 563   (interactive "P")
 564   (calc-wrapper
 565    (if (calc-is-inverse)
 566        (calc-unary-op "bpck" 'calcFunc-vpack arg)
 567      (calc-unary-op "bupk" 'calcFunc-vunpack arg)))
 568 )
 569
 570 (defun calc-pack-bits (arg)
 571   (interactive "P")
 572   (calc-invert-func)
 573   (calc-unpack-bits arg)
 574 )
 575
 576
 577 (defun calc-rnorm (arg)
 578   (interactive "P")
 579   (calc-wrapper
 580    (calc-unary-op "rnrm" 'calcFunc-rnorm arg))
 581 )
 582
 583 (defun calc-cnorm (arg)
 584   (interactive "P")
 585   (calc-wrapper
 586    (calc-unary-op "cnrm" 'calcFunc-cnorm arg))
 587 )
 588
 589 (defun calc-mrow (n &optional nn)
 590   (interactive "NRow number: \nP")
 591   (calc-wrapper
 592    (if (consp nn)
 593        (calc-enter-result 2 "mrow" (cons 'calcFunc-mrow (calc-top-list-n 2)))
 594      (setq n (prefix-numeric-value n))
 595      (if (= n 0)
 596          (calc-enter-result 1 "getd" (list 'calcFunc-getdiag (calc-top-n 1)))
 597        (if (< n 0)
 598            (calc-enter-result 1 "rrow" (list 'calcFunc-mrrow
 599                                              (calc-top-n 1) (- n)))
 600          (calc-enter-result 1 "mrow" (list 'calcFunc-mrow
 601                                            (calc-top-n 1) n))))))
 602 )
 603
 604 (defun calc-mcol (n &optional nn)
 605   (interactive "NColumn number: \nP")
 606   (calc-wrapper
 607    (if (consp nn)
 608        (calc-enter-result 2 "mcol" (cons 'calcFunc-mcol (calc-top-list-n 2)))
 609      (setq n (prefix-numeric-value n))
 610      (if (= n 0)
 611          (calc-enter-result 1 "getd" (list 'calcFunc-getdiag (calc-top-n 1)))
 612        (if (< n 0)
 613            (calc-enter-result 1 "rcol" (list 'calcFunc-mrcol
 614                                              (calc-top-n 1) (- n)))
 615          (calc-enter-result 1 "mcol" (list 'calcFunc-mcol
 616                                            (calc-top-n 1) n))))))
 617 )
 618
 619
 620 ;;;; Vectors.
 621
 622 (defun calcFunc-mdims (m)
 623   (or (math-vectorp m)
 624       (math-reject-arg m 'vectorp))
 625   (cons 'vec (math-mat-dimens m))
 626 )
 627
 628
 629 ;;; Apply a function elementwise to vector A.  [V X V; N X N] [Public]
 630 (defun math-map-vec (f a)
 631   (if (math-vectorp a)
 632       (cons 'vec (mapcar f (cdr a)))
 633     (funcall f a))
 634 )
 635
 636 (defun math-dimension-error ()
 637   (calc-record-why "*Dimension error")
 638   (signal 'wrong-type-argument nil)
 639 )
 640
 641
 642 ;;; Build a vector out of a list of objects.  [Public]
 643 (defun calcFunc-vec (&rest objs)
 644   (cons 'vec objs)
 645 )
 646
 647
 648 ;;; Build a constant vector or matrix.  [Public]
 649 (defun calcFunc-cvec (obj &rest dims)
 650   (math-make-vec-dimen obj dims)
 651 )
 652
 653 (defun math-make-vec-dimen (obj dims)
 654   (if dims
 655       (if (natnump (car dims))
 656           (if (or (cdr dims)
 657                   (not (math-numberp obj)))
 658               (cons 'vec (copy-sequence
 659                           (make-list (car dims)
 660                                      (math-make-vec-dimen obj (cdr dims)))))
 661             (cons 'vec (make-list (car dims) obj)))
 662         (math-reject-arg (car dims) 'fixnatnump))
 663     obj)
 664 )
 665
 666 (defun calcFunc-head (vec)
 667   (if (and (Math-vectorp vec)
 668            (cdr vec))
 669       (nth 1 vec)
 670     (calc-record-why 'vectorp vec)
 671     (list 'calcFunc-head vec))
 672 )
 673
 674 (defun calcFunc-tail (vec)
 675   (if (and (Math-vectorp vec)
 676            (cdr vec))
 677       (cons 'vec (cdr (cdr vec)))
 678     (calc-record-why 'vectorp vec)
 679     (list 'calcFunc-tail vec))
 680 )
 681
 682 (defun calcFunc-cons (head tail)
 683   (if (Math-vectorp tail)
 684       (cons 'vec (cons head (cdr tail)))
 685     (calc-record-why 'vectorp tail)
 686     (list 'calcFunc-cons head tail))
 687 )
 688
 689 (defun calcFunc-rhead (vec)
 690   (if (and (Math-vectorp vec)
 691            (cdr vec))
 692       (let ((vec (copy-sequence vec)))
 693         (setcdr (nthcdr (- (length vec) 2) vec) nil)
 694         vec)
 695     (calc-record-why 'vectorp vec)
 696     (list 'calcFunc-rhead vec))
 697 )
 698
 699 (defun calcFunc-rtail (vec)
 700   (if (and (Math-vectorp vec)
 701            (cdr vec))
 702       (nth (1- (length vec)) vec)
 703     (calc-record-why 'vectorp vec)
 704     (list 'calcFunc-rtail vec))
 705 )
 706
 707 (defun calcFunc-rcons (head tail)
 708   (if (Math-vectorp head)
 709       (append head (list tail))
 710     (calc-record-why 'vectorp head)
 711     (list 'calcFunc-rcons head tail))
 712 )
 713
 714
 715
 716 ;;; Apply a function elementwise to vectors A and B.  [O X O O] [Public]
 717 (defun math-map-vec-2 (f a b)
 718   (if (math-vectorp a)
 719       (if (math-vectorp b)
 720           (let ((v nil))
 721             (while (setq a (cdr a))
 722               (or (setq b (cdr b))
 723                   (math-dimension-error))
 724               (setq v (cons (funcall f (car a) (car b)) v)))
 725             (if a (math-dimension-error))
 726             (cons 'vec (nreverse v)))
 727         (let ((v nil))
 728           (while (setq a (cdr a))
 729             (setq v (cons (funcall f (car a) b) v)))
 730           (cons 'vec (nreverse v))))
 731     (if (math-vectorp b)
 732         (let ((v nil))
 733           (while (setq b (cdr b))
 734             (setq v (cons (funcall f a (car b)) v)))
 735           (cons 'vec (nreverse v)))
 736       (funcall f a b)))
 737 )
 738
 739
 740
 741 ;;; "Reduce" a function over a vector (left-associatively).  [O X V] [Public]
 742 (defun math-reduce-vec (f a)
 743   (if (math-vectorp a)
 744       (if (cdr a)
 745           (let ((accum (car (setq a (cdr a)))))
 746             (while (setq a (cdr a))
 747               (setq accum (funcall f accum (car a))))
 748             accum)
 749         0)
 750     a)
 751 )
 752
 753 ;;; Reduce a function over the columns of matrix A.  [V X V] [Public]
 754 (defun math-reduce-cols (f a)
 755   (if (math-matrixp a)
 756       (cons 'vec (math-reduce-cols-col-step f (cdr a) 1 (length (nth 1 a))))
 757     a)
 758 )
 759
 760 (defun math-reduce-cols-col-step (f a col cols)
 761   (and (< col cols)
 762        (cons (math-reduce-cols-row-step f (nth col (car a)) col (cdr a))
 763              (math-reduce-cols-col-step f a (1+ col) cols)))
 764 )
 765
 766 (defun math-reduce-cols-row-step (f tot col a)
 767   (if a
 768       (math-reduce-cols-row-step f
 769                                  (funcall f tot (nth col (car a)))
 770                                  col
 771                                  (cdr a))
 772     tot)
 773 )
 774
 775
 776
 777 (defun math-dot-product (a b)
 778   (if (setq a (cdr a) b (cdr b))
 779       (let ((accum (math-mul (car a) (car b))))
 780         (while (setq a (cdr a) b (cdr b))
 781           (setq accum (math-add accum (math-mul (car a) (car b)))))
 782         accum)
 783     0)
 784 )
 785
 786
 787 ;;; Return the number of elements in vector V.  [Public]
 788 (defun calcFunc-vlen (v)
 789   (if (math-vectorp v)
 790       (1- (length v))
 791     (if (math-objectp v)
 792         0
 793       (list 'calcFunc-vlen v)))
 794 )
 795
 796 ;;; Get the Nth row of a matrix.
 797 (defun calcFunc-mrow (mat n)   ; [Public]
 798   (if (Math-vectorp n)
 799       (math-map-vec (function (lambda (x) (calcFunc-mrow mat x))) n)
 800     (if (and (eq (car-safe n) 'intv) (math-constp n))
 801         (calcFunc-subvec mat
 802                          (math-add (nth 2 n) (if (memq (nth 1 n) '(2 3)) 0 1))
 803                          (math-add (nth 3 n) (if (memq (nth 1 n) '(1 3)) 1 0)))
 804       (or (and (integerp (setq n (math-check-integer n)))
 805                (> n 0))
 806           (math-reject-arg n 'fixposintp))
 807       (or (Math-vectorp mat)
 808           (math-reject-arg mat 'vectorp))
 809       (or (nth n mat)
 810           (math-reject-arg n "*Index out of range"))))
 811 )
 812
 813 (defun calcFunc-subscr (mat n &optional m)
 814   (setq mat (calcFunc-mrow mat n))
 815   (if m
 816       (if (math-num-integerp n)
 817           (calcFunc-mrow mat m)
 818         (calcFunc-mcol mat m))
 819     mat)
 820 )
 821
 822 ;;; Get the Nth column of a matrix.
 823 (defun math-mat-col (mat n)
 824   (cons 'vec (mapcar (function (lambda (x) (elt x n))) (cdr mat)))
 825 )
 826
 827 (defun calcFunc-mcol (mat n)   ; [Public]
 828   (if (Math-vectorp n)
 829       (calcFunc-trn
 830        (math-map-vec (function (lambda (x) (calcFunc-mcol mat x))) n))
 831     (if (and (eq (car-safe n) 'intv) (math-constp n))
 832         (if (math-matrixp mat)
 833             (math-map-vec (function (lambda (x) (calcFunc-mrow x n))) mat)
 834           (calcFunc-mrow mat n))
 835       (or (and (integerp (setq n (math-check-integer n)))
 836                (> n 0))
 837           (math-reject-arg n 'fixposintp))
 838       (or (Math-vectorp mat)
 839           (math-reject-arg mat 'vectorp))
 840       (or (if (math-matrixp mat)
 841               (and (< n (length (nth 1 mat)))
 842                    (math-mat-col mat n))
 843             (nth n mat))
 844           (math-reject-arg n "*Index out of range"))))
 845 )
 846
 847 ;;; Remove the Nth row from a matrix.
 848 (defun math-mat-less-row (mat n)
 849   (if (<= n 0)
 850       (cdr mat)
 851     (cons (car mat)
 852           (math-mat-less-row (cdr mat) (1- n))))
 853 )
 854
 855 (defun calcFunc-mrrow (mat n)   ; [Public]
 856   (and (integerp (setq n (math-check-integer n)))
 857        (> n 0)
 858        (< n (length mat))
 859        (math-mat-less-row mat n))
 860 )
 861
 862 ;;; Remove the Nth column from a matrix.
 863 (defun math-mat-less-col (mat n)
 864   (cons 'vec (mapcar (function (lambda (x) (math-mat-less-row x n)))
 865                      (cdr mat)))
 866 )
 867
 868 (defun calcFunc-mrcol (mat n)   ; [Public]
 869   (and (integerp (setq n (math-check-integer n)))
 870        (> n 0)
 871        (if (math-matrixp mat)
 872            (and (< n (length (nth 1 mat)))
 873                 (math-mat-less-col mat n))
 874          (math-mat-less-row mat n)))
 875 )
 876
 877 (defun calcFunc-getdiag (mat)   ; [Public]
 878   (if (math-square-matrixp mat)
 879       (cons 'vec (math-get-diag-step (cdr mat) 1))
 880     (calc-record-why 'square-matrixp mat)
 881     (list 'calcFunc-getdiag mat))
 882 )
 883
 884 (defun math-get-diag-step (row n)
 885   (and row
 886        (cons (nth n (car row))
 887              (math-get-diag-step (cdr row) (1+ n))))
 888 )
 889
 890 (defun math-transpose (mat)   ; [Public]
 891   (let ((m nil)
 892         (col (length (nth 1 mat))))
 893     (while (> (setq col (1- col)) 0)
 894       (setq m (cons (math-mat-col mat col) m)))
 895     (cons 'vec m))
 896 )
 897
 898 (defun calcFunc-trn (mat)
 899   (if (math-vectorp mat)
 900       (if (math-matrixp mat)
 901           (math-transpose mat)
 902         (math-col-matrix mat))
 903     (if (math-numberp mat)
 904         mat
 905       (math-reject-arg mat 'matrixp)))
 906 )
 907
 908 (defun calcFunc-ctrn (mat)
 909   (calcFunc-conj (calcFunc-trn mat))
 910 )
 911
 912 (defun calcFunc-pack (mode els)
 913   (or (Math-vectorp els) (math-reject-arg els 'vectorp))
 914   (if (and (Math-vectorp mode) (cdr mode))
 915       (setq mode (cdr mode))
 916     (or (integerp mode) (math-reject-arg mode 'fixnump)))
 917   (condition-case err
 918       (if (= (calc-pack-size mode) (1- (length els)))
 919           (calc-pack-items mode (cdr els))
 920         (math-reject-arg els "*Wrong number of elements"))
 921     (error (math-reject-arg els (nth 1 err))))
 922 )
 923
 924 (defun calcFunc-unpack (mode thing)
 925   (or (integerp mode) (math-reject-arg mode 'fixnump))
 926   (condition-case err
 927       (cons 'vec (calc-unpack-item mode thing))
 928     (error (math-reject-arg thing (nth 1 err))))
 929 )
 930
 931 (defun calcFunc-unpackt (mode thing)
 932   (let ((calc-unpack-with-type 'pair))
 933     (calcFunc-unpack mode thing))
 934 )
 935
 936 (defun calcFunc-arrange (vec cols)   ; [Public]
 937   (setq cols (math-check-fixnum cols t))
 938   (if (math-vectorp vec)
 939       (let* ((flat (math-flatten-vector vec))
 940              (mat (list 'vec))
 941              next)
 942         (if (<= cols 0)
 943             (nconc mat flat)
 944           (while (>= (length flat) cols)
 945             (setq next (nthcdr cols flat))
 946             (setcdr (nthcdr (1- cols) flat) nil)
 947             (setq mat (nconc mat (list (cons 'vec flat)))
 948                   flat next))
 949           (if flat
 950               (setq mat (nconc mat (list (cons 'vec flat)))))
 951           mat)))
 952 )
 953
 954 (defun math-flatten-vector (vec)   ; [L V]
 955   (if (math-vectorp vec)
 956       (apply 'append (mapcar 'math-flatten-vector (cdr vec)))
 957     (list vec))
 958 )
 959
 960 (defun calcFunc-vconcat (a b)
 961   (math-normalize (list '| a b))
 962 )
 963
 964 (defun calcFunc-vconcatrev (a b)
 965   (math-normalize (list '| b a))
 966 )
 967
 968 (defun calcFunc-append (v1 v2)
 969   (if (and (math-vectorp v1) (math-vectorp v2))
 970       (append v1 (cdr v2))
 971     (list 'calcFunc-append v1 v2))
 972 )
 973
 974 (defun calcFunc-appendrev (v1 v2)
 975   (calcFunc-append v2 v1)
 976 )
 977
 978
 979 ;;; Copy a matrix.  [Public]
 980 (defun math-copy-matrix (m)
 981   (if (math-vectorp (nth 1 m))
 982       (cons 'vec (mapcar 'copy-sequence (cdr m)))
 983     (copy-sequence m))
 984 )
 985
 986 ;;; Convert a scalar or vector into an NxN diagonal matrix.  [Public]
 987 (defun calcFunc-diag (a &optional n)
 988   (and n (not (integerp n))
 989        (setq n (math-check-fixnum n)))
 990   (if (math-vectorp a)
 991       (if (and n (/= (length a) (1+ n)))
 992           (list 'calcFunc-diag a n)
 993         (if (math-matrixp a)
 994             (if (and n (/= (length (elt a 1)) (1+ n)))
 995                 (list 'calcFunc-diag a n)
 996               a)
 997           (cons 'vec (math-diag-step (cdr a) 0 (1- (length a))))))
 998     (if n
 999         (cons 'vec (math-diag-step (make-list n a) 0 n))
1000       (list 'calcFunc-diag a)))
1001 )
1002
1003 (defun calcFunc-idn (a &optional n)
1004   (if n
1005       (if (math-vectorp a)
1006           (math-reject-arg a 'numberp)
1007         (calcFunc-diag a n))
1008     (if (integerp calc-matrix-mode)
1009         (calcFunc-idn a calc-matrix-mode)
1010       (list 'calcFunc-idn a)))
1011 )
1012
1013 (defun math-mimic-ident (a m)
1014   (if (math-square-matrixp m)
1015       (calcFunc-idn a (1- (length m)))
1016     (if (math-vectorp m)
1017         (if (math-zerop a)
1018             (cons 'vec (mapcar (function (lambda (x)
1019                                            (if (math-vectorp x)
1020                                                (math-mimic-ident a x)
1021                                              a)))
1022                                (cdr m)))
1023           (math-dimension-error))
1024       (calcFunc-idn a)))
1025 )
1026
1027 (defun math-diag-step (a n m)
1028   (if (< n m)
1029       (cons (cons 'vec
1030                   (nconc (make-list n 0)
1031                          (cons (car a)
1032                                (make-list (1- (- m n)) 0))))
1033             (math-diag-step (cdr a) (1+ n) m))
1034     nil)
1035 )
1036
1037 ;;; Create a vector of consecutive integers. [Public]
1038 (defun calcFunc-index (n &optional start incr)
1039   (if (math-messy-integerp n)
1040       (math-float (calcFunc-index (math-trunc n) start incr))
1041     (and (not (integerp n))
1042          (setq n (math-check-fixnum n)))
1043     (let ((vec nil))
1044       (if start
1045           (progn
1046             (if (>= n 0)
1047                 (while (>= (setq n (1- n)) 0)
1048                   (setq vec (cons start vec)
1049                         start (math-add start (or incr 1))))
1050               (while (<= (setq n (1+ n)) 0)
1051                 (setq vec (cons start vec)
1052                       start (math-mul start (or incr 2)))))
1053             (setq vec (nreverse vec)))
1054         (if (>= n 0)
1055             (while (> n 0)
1056               (setq vec (cons n vec)
1057                     n (1- n)))
1058           (let ((i -1))
1059             (while (>= i n)
1060               (setq vec (cons i vec)
1061                     i (1- i))))))
1062       (cons 'vec vec)))
1063 )
1064
1065 ;;; Find an element in a vector.
1066 (defun calcFunc-find (vec x &optional start)
1067   (setq start (if start (math-check-fixnum start t) 1))
1068   (if (< start 1) (math-reject-arg start 'posp))
1069   (setq vec (nthcdr start vec))
1070   (let ((n start))
1071     (while (and vec (not (Math-equal x (car vec))))
1072       (setq n (1+ n)
1073             vec (cdr vec)))
1074     (if vec n 0))
1075 )
1076
1077 ;;; Return a subvector of a vector.
1078 (defun calcFunc-subvec (vec start &optional end)
1079   (setq start (math-check-fixnum start t)
1080         end (math-check-fixnum (or end 0) t))
1081   (or (math-vectorp vec) (math-reject-arg vec 'vectorp))
1082   (let ((len (1- (length vec))))
1083     (if (<= start 0)
1084         (setq start (+ len start 1)))
1085     (if (<= end 0)
1086         (setq end (+ len end 1)))
1087     (if (or (> start len)
1088             (<= end start))
1089         '(vec)
1090       (setq vec (nthcdr start vec))
1091       (if (<= end len)
1092           (let ((chop (nthcdr (- end start 1) (setq vec (copy-sequence vec)))))
1093             (setcdr chop nil)))
1094       (cons 'vec vec)))
1095 )
1096
1097 ;;; Remove a subvector from a vector.
1098 (defun calcFunc-rsubvec (vec start &optional end)
1099   (setq start (math-check-fixnum start t)
1100         end (math-check-fixnum (or end 0) t))
1101   (or (math-vectorp vec) (math-reject-arg vec 'vectorp))
1102   (let ((len (1- (length vec))))
1103     (if (<= start 0)
1104         (setq start (+ len start 1)))
1105     (if (<= end 0)
1106         (setq end (+ len end 1)))
1107     (if (or (> start len)
1108             (<= end start))
1109         vec
1110       (let ((tail (nthcdr end vec))
1111             (chop (nthcdr (1- start) (setq vec (copy-sequence vec)))))
1112         (setcdr chop nil)
1113         (append vec tail))))
1114 )
1115
1116 ;;; Reverse the order of the elements of a vector.
1117 (defun calcFunc-rev (vec)
1118   (if (math-vectorp vec)
1119       (cons 'vec (reverse (cdr vec)))
1120     (math-reject-arg vec 'vectorp))
1121 )
1122
1123 ;;; Compress a vector according to a mask vector.
1124 (defun calcFunc-vmask (mask vec)
1125   (if (math-numberp mask)
1126       (if (math-zerop mask)
1127           '(vec)
1128         vec)
1129     (or (math-vectorp mask) (math-reject-arg mask 'vectorp))
1130     (or (math-constp mask) (math-reject-arg mask 'constp))
1131     (or (math-vectorp vec) (math-reject-arg vec 'vectorp))
1132     (or (= (length mask) (length vec)) (math-dimension-error))
1133     (let ((new nil))
1134       (while (setq mask (cdr mask) vec (cdr vec))
1135         (or (math-zerop (car mask))
1136             (setq new (cons (car vec) new))))
1137       (cons 'vec (nreverse new))))
1138 )
1139
1140 ;;; Expand a vector according to a mask vector.
1141 (defun calcFunc-vexp (mask vec &optional filler)
1142   (or (math-vectorp mask) (math-reject-arg mask 'vectorp))
1143   (or (math-constp mask) (math-reject-arg mask 'constp))
1144   (or (math-vectorp vec) (math-reject-arg vec 'vectorp))
1145   (let ((new nil)
1146         (fvec (and filler (math-vectorp filler))))
1147     (while (setq mask (cdr mask))
1148       (if (math-zerop (car mask))
1149           (setq new (cons (or (if fvec
1150                                   (car (setq filler (cdr filler)))
1151                                 filler)
1152                               (car mask)) new))
1153         (setq vec (cdr vec)
1154               new (cons (or (car vec) (car mask)) new))))
1155     (cons 'vec (nreverse new)))
1156 )
1157
1158
1159 ;;; Compute the row and column norms of a vector or matrix.  [Public]
1160 (defun calcFunc-rnorm (a)
1161   (if (and (Math-vectorp a)
1162            (math-constp a))
1163       (if (math-matrixp a)
1164           (math-reduce-vec 'math-max (math-map-vec 'calcFunc-cnorm a))
1165         (math-reduce-vec 'math-max (math-map-vec 'math-abs a)))
1166     (calc-record-why 'vectorp a)
1167     (list 'calcFunc-rnorm a))
1168 )
1169
1170 (defun calcFunc-cnorm (a)
1171   (if (and (Math-vectorp a)
1172            (math-constp a))
1173       (if (math-matrixp a)
1174           (math-reduce-vec 'math-max
1175                            (math-reduce-cols 'math-add-abs a))
1176         (math-reduce-vec 'math-add-abs a))
1177     (calc-record-why 'vectorp a)
1178     (list 'calcFunc-cnorm a))
1179 )
1180
1181 (defun math-add-abs (a b)
1182   (math-add (math-abs a) (math-abs b))
1183 )
1184
1185
1186 ;;; Sort the elements of a vector into increasing order.
1187 (defun calcFunc-sort (vec)   ; [Public]
1188   (if (math-vectorp vec)
1189       (cons 'vec (sort (copy-sequence (cdr vec)) 'math-beforep))
1190     (math-reject-arg vec 'vectorp))
1191 )
1192
1193 (defun calcFunc-rsort (vec)   ; [Public]
1194   (if (math-vectorp vec)
1195       (cons 'vec (nreverse (sort (copy-sequence (cdr vec)) 'math-beforep)))
1196     (math-reject-arg vec 'vectorp))
1197 )
1198
1199 (defun calcFunc-grade (grade-vec)
1200   (if (math-vectorp grade-vec)
1201       (let* ((len (1- (length grade-vec))))
1202         (cons 'vec (sort (cdr (calcFunc-index len)) 'math-grade-beforep)))
1203     (math-reject-arg grade-vec 'vectorp))
1204 )
1205
1206 (defun calcFunc-rgrade (grade-vec)
1207   (if (math-vectorp grade-vec)
1208       (let* ((len (1- (length grade-vec))))
1209         (cons 'vec (nreverse (sort (cdr (calcFunc-index len))
1210                                    'math-grade-beforep))))
1211     (math-reject-arg grade-vec 'vectorp))
1212 )
1213
1214 (defun math-grade-beforep (i j)
1215   (math-beforep (nth i grade-vec) (nth j grade-vec))
1216 )
1217
1218
1219 ;;; Compile a histogram of data from a vector.
1220 (defun calcFunc-histogram (vec wts &optional n)
1221   (or n (setq n wts wts 1))
1222   (or (Math-vectorp vec)
1223       (math-reject-arg vec 'vectorp))
1224   (if (Math-vectorp wts)
1225       (or (= (length vec) (length wts))
1226           (math-dimension-error)))
1227   (or (natnump n)
1228       (math-reject-arg n 'fixnatnump))
1229   (let ((res (make-vector n 0))
1230         (vp vec)
1231         (wvec (Math-vectorp wts))
1232         (wp wts)
1233         bin)
1234     (while (setq vp (cdr vp))
1235       (setq bin (car vp))
1236       (or (natnump bin)
1237           (setq bin (math-floor bin)))
1238       (and (natnump bin)
1239            (< bin n)
1240            (aset res bin (math-add (aref res bin)
1241                                    (if wvec (car (setq wp (cdr wp))) wts)))))
1242     (cons 'vec (append res nil)))
1243 )
1244
1245
1246 ;;; Set operations.
1247
1248 (defun calcFunc-vunion (a b)
1249   (if (Math-objectp a)
1250       (setq a (list 'vec a))
1251     (or (math-vectorp a) (math-reject-arg a 'vectorp)))
1252   (if (Math-objectp b)
1253       (setq b (list b))
1254     (or (math-vectorp b) (math-reject-arg b 'vectorp))
1255     (setq b (cdr b)))
1256   (calcFunc-rdup (append a b))
1257 )
1258
1259 (defun calcFunc-vint (a b)
1260   (if (and (math-simple-set a) (math-simple-set b))
1261       (progn
1262         (setq a (cdr (calcFunc-rdup a)))
1263         (setq b (cdr (calcFunc-rdup b)))
1264         (let ((vec (list 'vec)))
1265           (while (and a b)
1266             (if (math-beforep (car a) (car b))
1267                 (setq a (cdr a))
1268               (if (Math-equal (car a) (car b))
1269                   (setq vec (cons (car a) vec)
1270                         a (cdr a)))
1271               (setq b (cdr b))))
1272           (nreverse vec)))
1273     (calcFunc-vcompl (calcFunc-vunion (calcFunc-vcompl a)
1274                                       (calcFunc-vcompl b))))
1275 )
1276
1277 (defun calcFunc-vdiff (a b)
1278   (if (and (math-simple-set a) (math-simple-set b))
1279       (progn
1280         (setq a (cdr (calcFunc-rdup a)))
1281         (setq b (cdr (calcFunc-rdup b)))
1282         (let ((vec (list 'vec)))
1283           (while a
1284             (while (and b (math-beforep (car b) (car a)))
1285               (setq b (cdr b)))
1286             (if (and b (Math-equal (car a) (car b)))
1287                 (setq a (cdr a)
1288                       b (cdr b))
1289               (setq vec (cons (car a) vec)
1290                     a (cdr a))))
1291           (nreverse vec)))
1292     (calcFunc-vcompl (calcFunc-vunion (calcFunc-vcompl a) b)))
1293 )
1294
1295 (defun calcFunc-vxor (a b)
1296   (if (and (math-simple-set a) (math-simple-set b))
1297       (progn
1298         (setq a (cdr (calcFunc-rdup a)))
1299         (setq b (cdr (calcFunc-rdup b)))
1300         (let ((vec (list 'vec)))
1301           (while (or a b)
1302             (if (and a
1303                      (or (not b)
1304                          (math-beforep (car a) (car b))))
1305                 (setq vec (cons (car a) vec)
1306                       a (cdr a))
1307               (if (and a (Math-equal (car a) (car b)))
1308                   (setq a (cdr a))
1309                 (setq vec (cons (car b) vec)))
1310               (setq b (cdr b))))
1311           (nreverse vec)))
1312     (let ((ca (calcFunc-vcompl a))
1313           (cb (calcFunc-vcompl b)))
1314       (calcFunc-vunion (calcFunc-vcompl (calcFunc-vunion ca b))
1315                        (calcFunc-vcompl (calcFunc-vunion a cb)))))
1316 )
1317
1318 (defun calcFunc-vcompl (a)
1319   (setq a (math-prepare-set a))
1320   (let ((vec (list 'vec))
1321         (prev '(neg (var inf var-inf)))
1322         (closed 2))
1323     (while (setq a (cdr a))
1324       (or (and (equal (nth 2 (car a)) '(neg (var inf var-inf)))
1325                (memq (nth 1 (car a)) '(2 3)))
1326           (setq vec (cons (list 'intv
1327                                 (+ closed
1328                                    (if (memq (nth 1 (car a)) '(0 1)) 1 0))
1329                                 prev
1330                                 (nth 2 (car a)))
1331                           vec)))
1332       (setq prev (nth 3 (car a))
1333             closed (if (memq (nth 1 (car a)) '(0 2)) 2 0)))
1334     (or (and (equal prev '(var inf var-inf))
1335              (= closed 0))
1336         (setq vec (cons (list 'intv (+ closed 1)
1337                               prev '(var inf var-inf))
1338                         vec)))
1339     (math-clean-set (nreverse vec)))
1340 )
1341
1342 (defun calcFunc-vspan (a)
1343   (setq a (math-prepare-set a))
1344   (if (cdr a)
1345       (let ((last (nth (1- (length a)) a)))
1346         (math-make-intv (+ (logand (nth 1 (nth 1 a)) 2)
1347                            (logand (nth 1 last) 1))
1348                         (nth 2 (nth 1 a))
1349                         (nth 3 last)))
1350     '(intv 2 0 0))
1351 )
1352
1353 (defun calcFunc-vfloor (a &optional always-vec)
1354   (setq a (math-prepare-set a))
1355   (let ((vec (list 'vec)) (p a) (prev nil) b mask)
1356     (while (setq p (cdr p))
1357       (setq mask (nth 1 (car p))
1358             a (nth 2 (car p))
1359             b (nth 3 (car p)))
1360       (and (memq mask '(0 1))
1361            (not (math-infinitep a))
1362            (setq mask (logior mask 2))
1363            (math-num-integerp a)
1364            (setq a (math-add a 1)))
1365       (setq a (math-ceiling a))
1366       (and (memq mask '(0 2))
1367            (not (math-infinitep b))
1368            (setq mask (logior mask 1))
1369            (math-num-integerp b)
1370            (setq b (math-sub b 1)))
1371       (setq b (math-floor b))
1372       (if (and prev (Math-equal (math-sub a 1) (nth 3 prev)))
1373           (setcar (nthcdr 3 prev) b)
1374         (or (Math-lessp b a)
1375             (setq vec (cons (setq prev (list 'intv mask a b)) vec)))))
1376     (setq vec (nreverse vec))
1377     (math-clean-set vec always-vec))
1378 )
1379
1380 (defun calcFunc-vcard (a)
1381   (setq a (calcFunc-vfloor a t))
1382   (or (math-constp a) (math-reject-arg a "*Set must be finite"))
1383   (let ((count 0))
1384     (while (setq a (cdr a))
1385       (if (eq (car-safe (car a)) 'intv)
1386           (setq count (math-add count (math-sub (nth 3 (car a))
1387                                                 (nth 2 (car a))))))
1388       (setq count (math-add count 1)))
1389     count)
1390 )
1391
1392 (defun calcFunc-venum (a)
1393   (setq a (calcFunc-vfloor a t))
1394   (or (math-constp a) (math-reject-arg a "*Set must be finite"))
1395   (let ((p a) next)
1396     (while (cdr p)
1397       (setq next (cdr p))
1398       (if (eq (car-safe (nth 1 p)) 'intv)
1399           (setcdr p (nconc (cdr (calcFunc-index (math-add
1400                                                  (math-sub (nth 3 (nth 1 p))
1401                                                            (nth 2 (nth 1 p)))
1402                                                  1)
1403                                                 (nth 2 (nth 1 p))))
1404                            (cdr (cdr p)))))
1405       (setq p next))
1406     a)
1407 )
1408
1409 (defun calcFunc-vpack (a)
1410   (setq a (calcFunc-vfloor a t))
1411   (if (and (cdr a)
1412            (math-negp (if (eq (car-safe (nth 1 a)) 'intv)
1413                           (nth 2 (nth 1 a))
1414                         (nth 1 a))))
1415       (math-reject-arg (nth 1 a) 'posp))
1416   (let ((accum 0))
1417     (while (setq a (cdr a))
1418       (if (eq (car-safe (car a)) 'intv)
1419           (if (equal (nth 3 (car a)) '(var inf var-inf))
1420               (setq accum (math-sub accum
1421                                     (math-power-of-2 (nth 2 (car a)))))
1422             (setq accum (math-add accum
1423                                   (math-sub
1424                                    (math-power-of-2 (1+ (nth 3 (car a))))
1425                                    (math-power-of-2 (nth 2 (car a)))))))
1426         (setq accum (math-add accum (math-power-of-2 (car a))))))
1427     accum)
1428 )
1429
1430 (defun calcFunc-vunpack (a &optional w)
1431   (or (math-num-integerp a) (math-reject-arg a 'integerp))
1432   (if w (setq a (math-clip a w)))
1433   (if (math-messy-integerp a) (setq a (math-trunc a)))
1434   (let* ((calc-number-radix 2)
1435          (neg (math-negp a))
1436          (aa (if neg (math-sub -1 a) a))
1437          (str (if (eq aa 0)
1438                   ""
1439                 (if (consp aa)
1440                     (math-format-bignum-binary (cdr aa))
1441                   (math-format-binary aa))))
1442          (zero (if neg ?1 ?0))
1443          (one (if neg ?0 ?1))
1444          (len (length str))
1445          (vec (list 'vec))
1446          (pos (1- len)) pos2)
1447     (while (>= pos 0)
1448       (if (eq (aref str pos) zero)
1449           (setq pos (1- pos))
1450         (setq pos2 pos)
1451         (while (and (>= pos 0) (eq (aref str pos) one))
1452           (setq pos (1- pos)))
1453         (setq vec (cons (if (= pos (1- pos2))
1454                             (- len pos2 1)
1455                           (list 'intv 3 (- len pos2 1) (- len pos 2)))
1456                         vec))))
1457     (if neg
1458         (setq vec (cons (list 'intv 2 len '(var inf var-inf)) vec)))
1459     (math-clean-set (nreverse vec)))
1460 )
1461
1462 (defun calcFunc-rdup (a)
1463   (if (math-simple-set a)
1464       (progn
1465         (and (Math-objectp a) (setq a (list 'vec a)))
1466         (or (math-vectorp a) (math-reject-arg a 'vectorp))
1467         (setq a (sort (copy-sequence (cdr a)) 'math-beforep))
1468         (let ((p a))
1469           (while (cdr p)
1470             (if (Math-equal (car p) (nth 1 p))
1471                 (setcdr p (cdr (cdr p)))
1472               (setq p (cdr p)))))
1473         (cons 'vec a))
1474     (math-clean-set (math-prepare-set a)))
1475 )
1476
1477 (defun math-prepare-set (a)
1478   (if (Math-objectp a)
1479       (setq a (list 'vec a))
1480     (or (math-vectorp a) (math-reject-arg a 'vectorp))
1481     (setq a (cons 'vec (sort (copy-sequence (cdr a)) 'math-beforep))))
1482   (let ((p a) res)
1483
1484     ;; Convert all elements to non-empty intervals.
1485     (while (cdr p)
1486       (if (eq (car-safe (nth 1 p)) 'intv)
1487           (if (math-intv-constp (nth 1 p))
1488               (if (and (memq (nth 1 (nth 1 p)) '(0 1 2))
1489                        (Math-equal (nth 2 (nth 1 p)) (nth 3 (nth 1 p))))
1490                   (setcdr p (cdr (cdr p)))
1491                 (setq p (cdr p)))
1492             (math-reject-arg (nth 1 p) 'constp))
1493         (or (Math-anglep (nth 1 p))
1494             (eq (car (nth 1 p)) 'date)
1495             (equal (nth 1 p) '(var inf var-inf))
1496             (equal (nth 1 p) '(neg (var inf var-inf)))
1497             (math-reject-arg (nth 1 p) 'realp))
1498         (setcar (cdr p) (list 'intv 3 (nth 1 p) (nth 1 p)))
1499         (setq p (cdr p))))
1500
1501     ;; Combine redundant intervals.
1502     (setq p a)
1503     (while (cdr (cdr p))
1504       (if (or (memq (setq res (math-compare (nth 3 (nth 1 p))
1505                                             (nth 2 (nth 2 p))))
1506                     '(-1 2))
1507               (and (eq res 0)
1508                    (memq (nth 1 (nth 1 p)) '(0 2))
1509                    (memq (nth 1 (nth 2 p)) '(0 1))))
1510           (setq p (cdr p))
1511         (setq res (math-compare (nth 3 (nth 1 p)) (nth 3 (nth 2 p))))
1512         (setcdr p (cons (list 'intv
1513                               (+ (logand (logior (nth 1 (nth 1 p))
1514                                                  (if (Math-equal
1515                                                       (nth 2 (nth 1 p))
1516                                                       (nth 2 (nth 2 p)))
1517                                                      (nth 1 (nth 2 p))
1518                                                    0))
1519                                          2)
1520                                  (logand (logior (if (memq res '(1 0 2))
1521                                                      (nth 1 (nth 1 p)) 0)
1522                                                  (if (memq res '(-1 0 2))
1523                                                      (nth 1 (nth 2 p)) 0))
1524                                          1))
1525                               (nth 2 (nth 1 p))
1526                               (if (eq res 1)
1527                                   (nth 3 (nth 1 p))
1528                                 (nth 3 (nth 2 p))))
1529                         (cdr (cdr (cdr p))))))))
1530   a
1531 )
1532
1533 (defun math-clean-set (a &optional always-vec)
1534   (let ((p a) res)
1535     (while (cdr p)
1536       (if (and (eq (car-safe (nth 1 p)) 'intv)
1537                (Math-equal (nth 2 (nth 1 p)) (nth 3 (nth 1 p))))
1538           (setcar (cdr p) (nth 2 (nth 1 p))))
1539       (setq p (cdr p)))
1540     (if (and (not (cdr (cdr a)))
1541              (eq (car-safe (nth 1 a)) 'intv)
1542              (not always-vec))
1543         (nth 1 a)
1544       a))
1545 )
1546
1547 (defun math-simple-set (a)
1548   (or (and (Math-objectp a)
1549            (not (eq (car-safe a) 'intv)))
1550       (and (Math-vectorp a)
1551            (progn
1552              (while (and (setq a (cdr a))
1553                          (not (eq (car-safe (car a)) 'intv))))
1554              (null a))))
1555 )
1556
1557
1558
1559
1560 ;;; Compute a right-handed vector cross product.  [O O O] [Public]
1561 (defun calcFunc-cross (a b)
1562   (if (and (eq (car-safe a) 'vec)
1563            (= (length a) 4))
1564       (if (and (eq (car-safe b) 'vec)
1565                (= (length b) 4))
1566           (list 'vec
1567                 (math-sub (math-mul (nth 2 a) (nth 3 b))
1568                           (math-mul (nth 3 a) (nth 2 b)))
1569                 (math-sub (math-mul (nth 3 a) (nth 1 b))
1570                           (math-mul (nth 1 a) (nth 3 b)))
1571                 (math-sub (math-mul (nth 1 a) (nth 2 b))
1572                           (math-mul (nth 2 a) (nth 1 b))))
1573         (math-reject-arg b "*Three-vector expected"))
1574     (math-reject-arg a "*Three-vector expected"))
1575 )
1576
1577
1578
1579
1580
1581 (defun math-read-brackets (space-sep close)
1582   (and space-sep (setq space-sep (not (math-check-for-commas))))
1583   (math-read-token)
1584   (while (eq exp-token 'space)
1585     (math-read-token))
1586   (if (or (equal exp-data close)
1587           (eq exp-token 'end))
1588       (progn
1589         (math-read-token)
1590         '(vec))
1591     (let ((save-exp-pos exp-pos)
1592           (save-exp-old-pos exp-old-pos)
1593           (save-exp-token exp-token)
1594           (save-exp-data exp-data)
1595           (vals (let ((exp-keep-spaces space-sep))
1596                   (if (or (equal exp-data "\\dots")
1597                           (equal exp-data "\\ldots"))
1598                       '(vec (neg (var inf var-inf)))
1599                     (catch 'syntax (math-read-vector))))))
1600       (if (stringp vals)
1601           (if space-sep
1602               (let ((error-exp-pos exp-pos)
1603                     (error-exp-old-pos exp-old-pos)
1604                     vals2)
1605                 (setq exp-pos save-exp-pos
1606                       exp-old-pos save-exp-old-pos
1607                       exp-token save-exp-token
1608                       exp-data save-exp-data)
1609                 (let ((exp-keep-spaces nil))
1610                   (setq vals2 (catch 'syntax (math-read-vector))))
1611                 (if (and (not (stringp vals2))
1612                          (or (assoc exp-data '(("\\ldots") ("\\dots") (";")))
1613                              (equal exp-data close)
1614                              (eq exp-token 'end)))
1615                     (setq space-sep nil
1616                           vals vals2)
1617                   (setq exp-pos error-exp-pos
1618                         exp-old-pos error-exp-old-pos)
1619                   (throw 'syntax vals)))
1620             (throw 'syntax vals)))
1621       (if (or (equal exp-data "\\dots")
1622               (equal exp-data "\\ldots"))
1623           (progn
1624             (math-read-token)
1625             (setq vals (if (> (length vals) 2)
1626                            (cons 'calcFunc-mul (cdr vals)) (nth 1 vals)))
1627             (let ((exp2 (if (or (equal exp-data close)
1628                                 (equal exp-data ")")
1629                                 (eq exp-token 'end))
1630                             '(var inf var-inf)
1631                           (math-read-expr-level 0))))
1632               (setq vals
1633                     (list 'intv
1634                           (if (equal exp-data ")") 2 3)
1635                           vals
1636                           exp2)))
1637             (if (not (or (equal exp-data close)
1638                          (equal exp-data ")")
1639                          (eq exp-token 'end)))
1640                 (throw 'syntax "Expected `]'")))
1641         (if (equal exp-data ";")
1642             (let ((exp-keep-spaces space-sep))
1643               (setq vals (cons 'vec (math-read-matrix (list vals))))))
1644         (if (not (or (equal exp-data close)
1645                      (eq exp-token 'end)))
1646             (throw 'syntax "Expected `]'")))
1647       (or (eq exp-token 'end)
1648           (math-read-token))
1649       vals))
1650 )
1651
1652 (defun math-check-for-commas (&optional balancing)
1653   (let ((count 0)
1654         (pos (1- exp-pos)))
1655     (while (and (>= count 0)
1656                 (setq pos (string-match
1657                            (if balancing "[],[{}()<>]" "[],[{}()]")
1658                            exp-str (1+ pos)))
1659                 (or (/= (aref exp-str pos) ?,) (> count 0) balancing))
1660       (cond ((memq (aref exp-str pos) '(?\[ ?\{ ?\( ?\<))
1661              (setq count (1+ count)))
1662             ((memq (aref exp-str pos) '(?\] ?\} ?\) ?\>))
1663              (setq count (1- count)))))
1664     (if balancing
1665         pos
1666       (and pos (= (aref exp-str pos) ?,))))
1667 )
1668
1669 (defun math-read-vector ()
1670   (let* ((val (list (math-read-expr-level 0)))
1671          (last val))
1672     (while (progn
1673              (while (eq exp-token 'space)
1674                (math-read-token))
1675              (and (not (eq exp-token 'end))
1676                   (not (equal exp-data ";"))
1677                   (not (equal exp-data close))
1678                   (not (equal exp-data "\\dots"))
1679                   (not (equal exp-data "\\ldots"))))
1680       (if (equal exp-data ",")
1681           (math-read-token))
1682       (while (eq exp-token 'space)
1683         (math-read-token))
1684       (let ((rest (list (math-read-expr-level 0))))
1685         (setcdr last rest)
1686         (setq last rest)))
1687     (cons 'vec val))
1688 )
1689
1690 (defun math-read-matrix (mat)
1691   (while (equal exp-data ";")
1692     (math-read-token)
1693     (while (eq exp-token 'space)
1694       (math-read-token))
1695     (setq mat (nconc mat (list (math-read-vector)))))
1696   mat
1697 )
1698