1 ;;; calc-math.el --- mathematical functions for Calc
3 ;; Copyright (C) 1990, 1991, 1992, 1993, 2001, 2002, 2003, 2004,
4 ;; 2005, 2006, 2007 Free Software Foundation, Inc.
6 ;; Author: David Gillespie <daveg@synaptics.com>
7 ;; Maintainer: Jay Belanger <jay.p.belanger@gmail.com>
9 ;; This file is part of GNU Emacs.
11 ;; GNU Emacs is free software; you can redistribute it and/or modify
12 ;; it under the terms of the GNU General Public License as published by
13 ;; the Free Software Foundation; either version 3, or (at your option)
16 ;; GNU Emacs is distributed in the hope that it will be useful,
17 ;; but WITHOUT ANY WARRANTY; without even the implied warranty of
18 ;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 ;; GNU General Public License for more details.
21 ;; You should have received a copy of the GNU General Public License
22 ;; along with GNU Emacs; see the file COPYING. If not, write to the
23 ;; Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
24 ;; Boston, MA 02110-1301, USA.
30 ;; This file is autoloaded from calc-ext.el.
36 ;;; Find out how many 9s in 9.9999... will give distinct Emacs floats,
37 ;;; then back off by one.
39 (defvar math-emacs-precision
42 (xx (+ x (* 9 (expt 10 (- n))))))
47 (setq xx (+ x (* 9 (expt 10 (- n)))))))
49 "The number of digits in an Emacs float.")
51 ;;; Find the largest power of 10 which is an Emacs float,
52 ;;; then back off by one so that any float d.dddd...eN
53 ;;; is an Emacs float, for acceptable d.dddd....
55 (defvar math-largest-emacs-expt
57 (while (condition-case nil
62 (while (condition-case nil
67 "The largest exponent which Calc will convert to an Emacs float.")
69 (defvar math-smallest-emacs-expt
71 (while (condition-case nil
76 (while (condition-case nil
81 "The smallest exponent which Calc will convert to an Emacs float.")
83 (defun math-use-emacs-fn (fn x)
84 "Use the native Emacs function FN to evaluate the Calc number X.
85 If this can't be done, return NIL."
87 (<= calc-internal-prec math-emacs-precision)
89 (let* ((fx (math-float x))
90 (xpon (+ (nth 2 x) (1- (math-numdigs (nth 1 x))))))
91 (and (<= math-smallest-emacs-expt xpon)
92 (<= xpon math-largest-emacs-expt)
97 (string-to-number (math-format-number (math-float x))))))
100 (defun calc-sqrt (arg)
103 (if (calc-is-inverse)
104 (calc-unary-op "^2" 'calcFunc-sqr arg)
105 (calc-unary-op "sqrt" 'calcFunc-sqrt arg))))
107 (defun calc-isqrt (arg)
110 (if (calc-is-inverse)
111 (calc-unary-op "^2" 'calcFunc-sqr arg)
112 (calc-unary-op "isqt" 'calcFunc-isqrt arg))))
115 (defun calc-hypot (arg)
118 (calc-binary-op "hypt" 'calcFunc-hypot arg)))
125 (defun calc-log10 (arg)
127 (calc-hyperbolic-func)
130 (defun calc-log (arg)
133 (if (calc-is-inverse)
134 (calc-binary-op "alog" 'calcFunc-alog arg)
135 (calc-binary-op "log" 'calcFunc-log arg))))
137 (defun calc-ilog (arg)
140 (if (calc-is-inverse)
141 (calc-binary-op "alog" 'calcFunc-alog arg)
142 (calc-binary-op "ilog" 'calcFunc-ilog arg))))
144 (defun calc-lnp1 (arg)
149 (defun calc-exp (arg)
152 (if (calc-is-hyperbolic)
153 (if (calc-is-inverse)
154 (calc-unary-op "lg10" 'calcFunc-log10 arg)
155 (calc-unary-op "10^" 'calcFunc-exp10 arg))
156 (if (calc-is-inverse)
157 (calc-unary-op "ln" 'calcFunc-ln arg)
158 (calc-unary-op "exp" 'calcFunc-exp arg)))))
160 (defun calc-expm1 (arg)
163 (if (calc-is-inverse)
164 (calc-unary-op "ln+1" 'calcFunc-lnp1 arg)
165 (calc-unary-op "ex-1" 'calcFunc-expm1 arg))))
170 (if (calc-is-inverse)
171 (if (calc-is-hyperbolic)
172 (if calc-symbolic-mode
173 (calc-pop-push-record 0 "phi" '(var phi var-phi))
174 (calc-pop-push-record 0 "phi" (math-phi)))
175 (if calc-symbolic-mode
176 (calc-pop-push-record 0 "gmma" '(var gamma var-gamma))
177 (calc-pop-push-record 0 "gmma" (math-gamma-const))))
178 (if (calc-is-hyperbolic)
179 (if calc-symbolic-mode
180 (calc-pop-push-record 0 "e" '(var e var-e))
181 (calc-pop-push-record 0 "e" (math-e)))
182 (if calc-symbolic-mode
183 (calc-pop-push-record 0 "pi" '(var pi var-pi))
184 (calc-pop-push-record 0 "pi" (math-pi)))))))
186 (defun calc-sin (arg)
189 (if (calc-is-hyperbolic)
190 (if (calc-is-inverse)
191 (calc-unary-op "asnh" 'calcFunc-arcsinh arg)
192 (calc-unary-op "sinh" 'calcFunc-sinh arg))
193 (if (calc-is-inverse)
194 (calc-unary-op "asin" 'calcFunc-arcsin arg)
195 (calc-unary-op "sin" 'calcFunc-sin arg)))))
197 (defun calc-arcsin (arg)
202 (defun calc-sinh (arg)
204 (calc-hyperbolic-func)
207 (defun calc-arcsinh (arg)
210 (calc-hyperbolic-func)
213 (defun calc-sec (arg)
216 (if (calc-is-hyperbolic)
217 (calc-unary-op "sech" 'calcFunc-sech arg)
218 (calc-unary-op "sec" 'calcFunc-sec arg))))
220 (defun calc-sech (arg)
222 (calc-hyperbolic-func)
225 (defun calc-cos (arg)
228 (if (calc-is-hyperbolic)
229 (if (calc-is-inverse)
230 (calc-unary-op "acsh" 'calcFunc-arccosh arg)
231 (calc-unary-op "cosh" 'calcFunc-cosh arg))
232 (if (calc-is-inverse)
233 (calc-unary-op "acos" 'calcFunc-arccos arg)
234 (calc-unary-op "cos" 'calcFunc-cos arg)))))
236 (defun calc-arccos (arg)
241 (defun calc-cosh (arg)
243 (calc-hyperbolic-func)
246 (defun calc-arccosh (arg)
249 (calc-hyperbolic-func)
252 (defun calc-csc (arg)
255 (if (calc-is-hyperbolic)
256 (calc-unary-op "csch" 'calcFunc-csch arg)
257 (calc-unary-op "csc" 'calcFunc-csc arg))))
259 (defun calc-csch (arg)
261 (calc-hyperbolic-func)
264 (defun calc-sincos ()
267 (if (calc-is-inverse)
268 (calc-enter-result 1 "asnc" (list 'calcFunc-arcsincos (calc-top-n 1)))
269 (calc-enter-result 1 "sncs" (list 'calcFunc-sincos (calc-top-n 1))))))
271 (defun calc-tan (arg)
274 (if (calc-is-hyperbolic)
275 (if (calc-is-inverse)
276 (calc-unary-op "atnh" 'calcFunc-arctanh arg)
277 (calc-unary-op "tanh" 'calcFunc-tanh arg))
278 (if (calc-is-inverse)
279 (calc-unary-op "atan" 'calcFunc-arctan arg)
280 (calc-unary-op "tan" 'calcFunc-tan arg)))))
282 (defun calc-arctan (arg)
287 (defun calc-tanh (arg)
289 (calc-hyperbolic-func)
292 (defun calc-arctanh (arg)
295 (calc-hyperbolic-func)
298 (defun calc-cot (arg)
301 (if (calc-is-hyperbolic)
302 (calc-unary-op "coth" 'calcFunc-coth arg)
303 (calc-unary-op "cot" 'calcFunc-cot arg))))
305 (defun calc-coth (arg)
307 (calc-hyperbolic-func)
310 (defun calc-arctan2 ()
313 (calc-enter-result 2 "atn2" (cons 'calcFunc-arctan2 (calc-top-list-n 2)))))
315 (defun calc-conj (arg)
318 (calc-unary-op "conj" 'calcFunc-conj arg)))
320 (defun calc-imaginary ()
323 (calc-pop-push-record 1 "i*" (math-imaginary (calc-top-n 1)))))
325 (defun calc-to-degrees (arg)
328 (calc-unary-op ">deg" 'calcFunc-deg arg)))
330 (defun calc-to-radians (arg)
333 (calc-unary-op ">rad" 'calcFunc-rad arg)))
336 (defun calc-degrees-mode (arg)
340 (calc-change-mode 'calc-angle-mode 'deg)
341 (message "Angles measured in degrees")))
342 ((= arg 2) (calc-radians-mode))
343 ((= arg 3) (calc-hms-mode))
344 (t (error "Prefix argument out of range"))))
346 (defun calc-radians-mode ()
349 (calc-change-mode 'calc-angle-mode 'rad)
350 (message "Angles measured in radians")))
353 ;;; Compute the integer square-root floor(sqrt(A)). A > 0. [I I] [Public]
354 ;;; This method takes advantage of the fact that Newton's method starting
355 ;;; with an overestimate always works, even using truncating integer division!
356 (defun math-isqrt (a)
357 (cond ((Math-zerop a) a)
358 ((not (math-natnump a))
359 (math-reject-arg a 'natnump))
361 (math-isqrt-small a))
363 (math-normalize (cons 'bigpos (cdr (math-isqrt-bignum (cdr a))))))))
365 (defun calcFunc-isqrt (a)
367 (math-isqrt (math-floor a))
368 (math-floor (math-sqrt a))))
371 ;;; This returns (flag . result) where the flag is t if A is a perfect square.
372 (defun math-isqrt-bignum (a) ; [P.l L]
373 (let ((len (length a)))
375 (let* ((top (nthcdr (- len 2) a)))
376 (math-isqrt-bignum-iter
378 (math-scale-bignum-digit-size
380 (1+ (math-isqrt-small
381 (+ (* (nth 1 top) math-bignum-digit-size) (car top)))))
383 (let* ((top (nth (1- len) a)))
384 (math-isqrt-bignum-iter
386 (math-scale-bignum-digit-size
387 (list (1+ (math-isqrt-small top)))
390 (defun math-isqrt-bignum-iter (a guess) ; [l L l]
391 (math-working "isqrt" (cons 'bigpos guess))
392 (let* ((q (math-div-bignum a guess))
393 (s (math-add-bignum (car q) guess))
394 (g2 (math-div2-bignum s))
395 (comp (math-compare-bignum g2 guess)))
397 (math-isqrt-bignum-iter a g2)
398 (cons (and (= comp 0)
399 (math-zerop-bignum (cdr q))
403 (defun math-zerop-bignum (a)
406 (while (eq (car (setq a (cdr a))) 0))
409 (defun math-scale-bignum-digit-size (a n) ; [L L S]
415 (defun math-isqrt-small (a) ; A > 0. [S S]
416 (let ((g (cond ((>= a 1000000) 10000)
421 (while (< (setq g2 (/ (+ g (/ a g)) 2)) g)
428 ;;; Compute the square root of a number.
429 ;;; [T N] if possible, else [F N] if possible, else [C N]. [Public]
432 (and (Math-zerop a) a)
433 (and (math-known-nonposp a)
434 (math-imaginary (math-sqrt (math-neg a))))
436 (let ((sqrt (math-isqrt-small a)))
437 (if (= (* sqrt sqrt) a)
439 (if calc-symbolic-mode
440 (list 'calcFunc-sqrt a)
441 (math-sqrt-float (math-float a) (math-float sqrt))))))
442 (and (eq (car-safe a) 'bigpos)
443 (let* ((res (math-isqrt-bignum (cdr a)))
444 (sqrt (math-normalize (cons 'bigpos (cdr res)))))
447 (if calc-symbolic-mode
448 (list 'calcFunc-sqrt a)
449 (math-sqrt-float (math-float a) (math-float sqrt))))))
450 (and (eq (car-safe a) 'frac)
451 (let* ((num-res (math-isqrt-bignum (cdr (Math-bignum-test (nth 1 a)))))
452 (num-sqrt (math-normalize (cons 'bigpos (cdr num-res))))
453 (den-res (math-isqrt-bignum (cdr (Math-bignum-test (nth 2 a)))))
454 (den-sqrt (math-normalize (cons 'bigpos (cdr den-res)))))
455 (if (and (car num-res) (car den-res))
456 (list 'frac num-sqrt den-sqrt)
457 (if calc-symbolic-mode
458 (if (or (car num-res) (car den-res))
459 (math-div (if (car num-res)
460 num-sqrt (list 'calcFunc-sqrt (nth 1 a)))
462 den-sqrt (list 'calcFunc-sqrt (nth 2 a))))
463 (list 'calcFunc-sqrt a))
464 (math-sqrt-float (math-float a)
465 (math-div (math-float num-sqrt) den-sqrt))))))
466 (and (eq (car-safe a) 'float)
467 (if calc-symbolic-mode
468 (if (= (% (nth 2 a) 2) 0)
469 (let ((res (math-isqrt-bignum
470 (cdr (Math-bignum-test (nth 1 a))))))
472 (math-make-float (math-normalize
473 (cons 'bigpos (cdr res)))
475 (signal 'inexact-result nil)))
476 (signal 'inexact-result nil))
477 (math-sqrt-float a)))
478 (and (eq (car-safe a) 'cplx)
479 (math-with-extra-prec 2
480 (let* ((d (math-abs a))
481 (imag (math-sqrt (math-mul (math-sub d (nth 1 a))
484 (math-sqrt (math-mul (math-add d (nth 1 a)) '(float 5 -1)))
485 (if (math-negp (nth 2 a)) (math-neg imag) imag)))))
486 (and (eq (car-safe a) 'polar)
488 (math-sqrt (nth 1 a))
489 (math-mul (nth 2 a) '(float 5 -1))))
490 (and (eq (car-safe a) 'sdev)
491 (let ((sqrt (math-sqrt (nth 1 a))))
493 (math-div (nth 2 a) (math-mul sqrt 2)))))
494 (and (eq (car-safe a) 'intv)
495 (not (math-negp (nth 2 a)))
496 (math-make-intv (nth 1 a) (math-sqrt (nth 2 a)) (math-sqrt (nth 3 a))))
497 (and (eq (car-safe a) '*)
498 (or (math-known-nonnegp (nth 1 a))
499 (math-known-nonnegp (nth 2 a)))
500 (math-mul (math-sqrt (nth 1 a)) (math-sqrt (nth 2 a))))
501 (and (eq (car-safe a) '/)
502 (or (and (math-known-nonnegp (nth 2 a))
503 (math-div (math-sqrt (nth 1 a)) (math-sqrt (nth 2 a))))
504 (and (math-known-nonnegp (nth 1 a))
505 (not (math-equal-int (nth 1 a) 1))
506 (math-mul (math-sqrt (nth 1 a))
507 (math-sqrt (math-div 1 (nth 2 a)))))))
508 (and (eq (car-safe a) '^)
509 (math-known-evenp (nth 2 a))
510 (math-known-realp (nth 1 a))
511 (math-abs (math-pow (nth 1 a) (math-div (nth 2 a) 2))))
512 (let ((inf (math-infinitep a)))
514 (math-mul (math-sqrt (math-infinite-dir a inf)) inf)))
516 (calc-record-why 'numberp a)
517 (list 'calcFunc-sqrt a))))
518 (defalias 'calcFunc-sqrt 'math-sqrt)
520 (defun math-infinite-dir (a &optional inf)
521 (or inf (setq inf (math-infinitep a)))
522 (math-normalize (math-expr-subst a inf 1)))
524 (defun math-sqrt-float (a &optional guess) ; [F F F]
525 (if calc-symbolic-mode
526 (signal 'inexact-result nil)
527 (math-with-extra-prec 1 (math-sqrt-raw a guess))))
529 (defun math-sqrt-raw (a &optional guess) ; [F F F]
530 (if (not (Math-posp a))
533 (let ((ldiff (- (math-numdigs (nth 1 a)) 6)))
534 (or (= (% (+ (nth 2 a) ldiff) 2) 0) (setq ldiff (1+ ldiff)))
535 (setq guess (math-make-float (math-isqrt-small
536 (math-scale-int (nth 1 a) (- ldiff)))
537 (/ (+ (nth 2 a) ldiff) 2)))))
538 (math-sqrt-float-iter a guess)))
540 (defun math-sqrt-float-iter (a guess) ; [F F F]
541 (math-working "sqrt" guess)
542 (let ((g2 (math-mul-float (math-add-float guess (math-div-float a guess))
544 (if (math-nearly-equal-float g2 guess)
546 (math-sqrt-float-iter a g2))))
548 ;;; True if A and B differ only in the last digit of precision. [P F F]
549 (defun math-nearly-equal-float (a b)
550 (let ((ediff (- (nth 2 a) (nth 2 b))))
551 (cond ((= ediff 0) ;; Expanded out for speed
552 (setq ediff (math-add (Math-integer-neg (nth 1 a)) (nth 1 b)))
554 (and (not (consp ediff))
557 (= (math-numdigs (nth 1 a)) calc-internal-prec))))
559 (setq ediff (math-add (Math-integer-neg (nth 1 b))
560 (math-scale-int (nth 1 a) 1)))
561 (and (not (consp ediff))
564 (= (math-numdigs (nth 1 b)) calc-internal-prec)))
566 (setq ediff (math-add (Math-integer-neg (nth 1 a))
567 (math-scale-int (nth 1 b) 1)))
568 (and (not (consp ediff))
571 (= (math-numdigs (nth 1 a)) calc-internal-prec))))))
573 (defun math-nearly-equal (a b) ; [P N N] [Public]
574 (setq a (math-float a))
575 (setq b (math-float b))
576 (if (eq (car a) 'polar) (setq a (math-complex a)))
577 (if (eq (car b) 'polar) (setq b (math-complex b)))
578 (if (eq (car a) 'cplx)
579 (if (eq (car b) 'cplx)
580 (and (or (math-nearly-equal-float (nth 1 a) (nth 1 b))
581 (and (math-nearly-zerop-float (nth 1 a) (nth 2 a))
582 (math-nearly-zerop-float (nth 1 b) (nth 2 b))))
583 (or (math-nearly-equal-float (nth 2 a) (nth 2 b))
584 (and (math-nearly-zerop-float (nth 2 a) (nth 1 a))
585 (math-nearly-zerop-float (nth 2 b) (nth 1 b)))))
586 (and (math-nearly-equal-float (nth 1 a) b)
587 (math-nearly-zerop-float (nth 2 a) b)))
588 (if (eq (car b) 'cplx)
589 (and (math-nearly-equal-float a (nth 1 b))
590 (math-nearly-zerop-float a (nth 2 b)))
591 (math-nearly-equal-float a b))))
593 ;;; True if A is nearly zero compared to B. [P F F]
594 (defun math-nearly-zerop-float (a b)
596 (<= (+ (math-numdigs (nth 1 a)) (nth 2 a))
597 (1+ (- (+ (math-numdigs (nth 1 b)) (nth 2 b)) calc-internal-prec)))))
599 (defun math-nearly-zerop (a b) ; [P N R] [Public]
600 (setq a (math-float a))
601 (setq b (math-float b))
602 (if (eq (car a) 'cplx)
603 (and (math-nearly-zerop-float (nth 1 a) b)
604 (math-nearly-zerop-float (nth 2 a) b))
605 (if (eq (car a) 'polar)
606 (math-nearly-zerop-float (nth 1 a) b)
607 (math-nearly-zerop-float a b))))
609 ;;; This implementation could be improved, accuracy-wise.
610 (defun math-hypot (a b)
611 (cond ((Math-zerop a) (math-abs b))
612 ((Math-zerop b) (math-abs a))
613 ((not (Math-scalarp a))
614 (if (math-infinitep a)
615 (if (math-infinitep b)
620 (calc-record-why 'scalarp a)
621 (list 'calcFunc-hypot a b)))
622 ((not (Math-scalarp b))
623 (if (math-infinitep b)
625 (calc-record-why 'scalarp b)
626 (list 'calcFunc-hypot a b)))
627 ((and (Math-numberp a) (Math-numberp b))
628 (math-with-extra-prec 1
629 (math-sqrt (math-add (calcFunc-abssqr a) (calcFunc-abssqr b)))))
630 ((eq (car-safe a) 'hms)
631 (if (eq (car-safe b) 'hms) ; this helps sdev's of hms forms
632 (math-to-hms (math-hypot (math-from-hms a 'deg)
633 (math-from-hms b 'deg)))
634 (math-to-hms (math-hypot (math-from-hms a 'deg) b))))
635 ((eq (car-safe b) 'hms)
636 (math-to-hms (math-hypot a (math-from-hms b 'deg))))
638 (defalias 'calcFunc-hypot 'math-hypot)
640 (defun calcFunc-sqr (x)
645 (defun math-nth-root (a n)
646 (cond ((= n 2) (math-sqrt a))
650 (let ((root (math-nth-root-integer a n)))
653 (and (not calc-symbolic-mode)
654 (math-nth-root-float (math-float a) n
655 (math-float (cdr root)))))))
656 ((eq (car-safe a) 'frac)
657 (let* ((num-root (math-nth-root-integer (nth 1 a) n))
658 (den-root (math-nth-root-integer (nth 2 a) n)))
659 (if (and (car num-root) (car den-root))
660 (list 'frac (cdr num-root) (cdr den-root))
661 (and (not calc-symbolic-mode)
664 (math-div-float (math-float (cdr num-root))
665 (math-float (cdr den-root))))))))
666 ((eq (car-safe a) 'float)
667 (and (not calc-symbolic-mode)
668 (math-nth-root-float a n)))
669 ((eq (car-safe a) 'polar)
670 (let ((root (math-nth-root (nth 1 a) n)))
671 (and root (list 'polar root (math-div (nth 2 a) n)))))
674 ;; The variables math-nrf-n, math-nrf-nf and math-nrf-nfm1 are local
675 ;; to math-nth-root-float, but are used by math-nth-root-float-iter,
676 ;; which is called by math-nth-root-float.
679 (defvar math-nrf-nfm1)
681 (defun math-nth-root-float (a math-nrf-n &optional guess)
682 (math-inexact-result)
683 (math-with-extra-prec 1
684 (let ((math-nrf-nf (math-float math-nrf-n))
685 (math-nrf-nfm1 (math-float (1- math-nrf-n))))
686 (math-nth-root-float-iter a (or guess
688 1 (/ (+ (math-numdigs (nth 1 a))
693 (defun math-nth-root-float-iter (a guess)
694 (math-working "root" guess)
695 (let ((g2 (math-div-float (math-add-float (math-mul math-nrf-nfm1 guess)
697 a (math-ipow guess (1- math-nrf-n))))
699 (if (math-nearly-equal-float g2 guess)
701 (math-nth-root-float-iter a g2))))
703 ;; The variable math-nri-n is local to math-nth-root-integer, but
704 ;; is used by math-nth-root-int-iter, which is called by
705 ;; math-nth-root-int.
708 (defun math-nth-root-integer (a math-nri-n &optional guess) ; [I I S]
709 (math-nth-root-int-iter a (or guess
710 (math-scale-int 1 (/ (+ (math-numdigs a)
714 (defun math-nth-root-int-iter (a guess)
715 (math-working "root" guess)
716 (let* ((q (math-idivmod a (math-ipow guess (1- math-nri-n))))
717 (s (math-add (car q) (math-mul (1- math-nri-n) guess)))
718 (g2 (math-idivmod s math-nri-n)))
719 (if (Math-natnum-lessp (car g2) guess)
720 (math-nth-root-int-iter a (car g2))
721 (cons (and (equal (car g2) guess)
726 (defun calcFunc-nroot (x n)
727 (calcFunc-pow x (if (integerp n)
734 ;;;; Transcendental functions.
736 ;;; All of these functions are defined on the complex plane.
737 ;;; (Branch cuts, etc. follow Steele's Common Lisp book.)
739 ;;; Most functions increase calc-internal-prec by 2 digits, then round
740 ;;; down afterward. "-raw" functions use the current precision, require
741 ;;; their arguments to be in float (or complex float) format, and always
742 ;;; work in radians (where applicable).
744 (defun math-to-radians (a) ; [N N]
745 (cond ((eq (car-safe a) 'hms)
746 (math-from-hms a 'rad))
747 ((memq calc-angle-mode '(deg hms))
748 (math-mul a (math-pi-over-180)))
751 (defun math-from-radians (a) ; [N N]
752 (cond ((eq calc-angle-mode 'deg)
754 (math-div a (math-pi-over-180))
755 (list 'calcFunc-deg a)))
756 ((eq calc-angle-mode 'hms)
757 (math-to-hms a 'rad))
760 (defun math-to-radians-2 (a) ; [N N]
761 (cond ((eq (car-safe a) 'hms)
762 (math-from-hms a 'rad))
763 ((memq calc-angle-mode '(deg hms))
764 (if calc-symbolic-mode
765 (math-div (math-mul a '(var pi var-pi)) 180)
766 (math-mul a (math-pi-over-180))))
769 (defun math-from-radians-2 (a) ; [N N]
770 (cond ((memq calc-angle-mode '(deg hms))
771 (if calc-symbolic-mode
772 (math-div (math-mul 180 a) '(var pi var-pi))
773 (math-div a (math-pi-over-180))))
778 ;;; Sine, cosine, and tangent.
780 (defun calcFunc-sin (x) ; [N N] [Public]
781 (cond ((and (integerp x)
782 (if (eq calc-angle-mode 'deg)
785 (aref [0 1 0 -1] (math-mod (/ x 90) 4)))
787 (math-with-extra-prec 2
788 (math-sin-raw (math-to-radians (math-float x)))))
791 (math-with-extra-prec 2
792 (let* ((xx (math-to-radians (math-float (nth 1 x))))
793 (xs (math-to-radians (math-float (nth 2 x))))
794 (sc (math-sin-cos-raw xx)))
795 (math-make-sdev (car sc) (math-mul xs (cdr sc)))))
796 (math-make-sdev (calcFunc-sin (nth 1 x))
797 (math-mul (nth 2 x) (calcFunc-cos (nth 1 x))))))
798 ((and (eq (car x) 'intv) (math-intv-constp x))
799 (calcFunc-cos (math-sub x (math-quarter-circle nil))))
800 ((equal x '(var nan var-nan))
802 (t (calc-record-why 'scalarp x)
803 (list 'calcFunc-sin x))))
805 (defun calcFunc-cos (x) ; [N N] [Public]
806 (cond ((and (integerp x)
807 (if (eq calc-angle-mode 'deg)
810 (aref [1 0 -1 0] (math-mod (/ x 90) 4)))
812 (math-with-extra-prec 2
813 (math-cos-raw (math-to-radians (math-float x)))))
816 (math-with-extra-prec 2
817 (let* ((xx (math-to-radians (math-float (nth 1 x))))
818 (xs (math-to-radians (math-float (nth 2 x))))
819 (sc (math-sin-cos-raw xx)))
820 (math-make-sdev (cdr sc) (math-mul xs (car sc)))))
821 (math-make-sdev (calcFunc-cos (nth 1 x))
822 (math-mul (nth 2 x) (calcFunc-sin (nth 1 x))))))
823 ((and (eq (car x) 'intv) (math-intv-constp x))
824 (math-with-extra-prec 2
825 (let* ((xx (math-to-radians (math-float x)))
826 (na (math-floor (math-div (nth 2 xx) (math-pi))))
827 (nb (math-floor (math-div (nth 3 xx) (math-pi))))
828 (span (math-sub nb na)))
829 (if (memq span '(0 1))
830 (let ((int (math-sort-intv (nth 1 x)
831 (math-cos-raw (nth 2 xx))
832 (math-cos-raw (nth 3 xx)))))
835 (math-make-intv (logior (nth 1 x) 2)
838 (math-make-intv (logior (nth 1 x) 1)
842 (list 'intv 3 -1 1)))))
843 ((equal x '(var nan var-nan))
845 (t (calc-record-why 'scalarp x)
846 (list 'calcFunc-cos x))))
848 (defun calcFunc-sincos (x) ; [V N] [Public]
850 (math-with-extra-prec 2
851 (let ((sc (math-sin-cos-raw (math-to-radians (math-float x)))))
852 (list 'vec (cdr sc) (car sc)))) ; the vector [cos, sin]
853 (list 'vec (calcFunc-sin x) (calcFunc-cos x))))
855 (defun calcFunc-tan (x) ; [N N] [Public]
856 (cond ((and (integerp x)
857 (if (eq calc-angle-mode 'deg)
862 (math-with-extra-prec 2
863 (math-tan-raw (math-to-radians (math-float x)))))
866 (math-with-extra-prec 2
867 (let* ((xx (math-to-radians (math-float (nth 1 x))))
868 (xs (math-to-radians (math-float (nth 2 x))))
869 (sc (math-sin-cos-raw xx)))
870 (if (and (math-zerop (cdr sc)) (not calc-infinite-mode))
872 (calc-record-why "*Division by zero")
873 (list 'calcFunc-tan x))
874 (math-make-sdev (math-div-float (car sc) (cdr sc))
875 (math-div-float xs (math-sqr (cdr sc)))))))
876 (math-make-sdev (calcFunc-tan (nth 1 x))
878 (math-sqr (calcFunc-cos (nth 1 x)))))))
879 ((and (eq (car x) 'intv) (math-intv-constp x))
880 (or (math-with-extra-prec 2
881 (let* ((xx (math-to-radians (math-float x)))
882 (na (math-floor (math-div (math-sub (nth 2 xx)
885 (nb (math-floor (math-div (math-sub (nth 3 xx)
889 (math-sort-intv (nth 1 x)
890 (math-tan-raw (nth 2 xx))
891 (math-tan-raw (nth 3 xx))))))
892 '(intv 3 (neg (var inf var-inf)) (var inf var-inf))))
893 ((equal x '(var nan var-nan))
895 (t (calc-record-why 'scalarp x)
896 (list 'calcFunc-tan x))))
898 (defun calcFunc-sec (x)
899 (cond ((and (integerp x)
900 (eq calc-angle-mode 'deg)
906 (eq calc-angle-mode 'rad)
910 (math-with-extra-prec 2
911 (math-sec-raw (math-to-radians (math-float x)))))
914 (math-with-extra-prec 2
915 (let* ((xx (math-to-radians (math-float (nth 1 x))))
916 (xs (math-to-radians (math-float (nth 2 x))))
917 (sc (math-sin-cos-raw xx)))
918 (if (and (math-zerop (cdr sc))
919 (not calc-infinite-mode))
921 (calc-record-why "*Division by zero")
922 (list 'calcFunc-sec x))
923 (math-make-sdev (math-div-float '(float 1 0) (cdr sc))
925 (math-mul xs (car sc))
926 (math-sqr (cdr sc)))))))
927 (math-make-sdev (calcFunc-sec (nth 1 x))
930 (calcFunc-sin (nth 1 x)))
931 (math-sqr (calcFunc-cos (nth 1 x)))))))
932 ((and (eq (car x) 'intv)
933 (math-intv-constp x))
934 (math-with-extra-prec 2
935 (let* ((xx (math-to-radians (math-float x)))
936 (na (math-floor (math-div (math-sub (nth 2 xx)
939 (nb (math-floor (math-div (math-sub (nth 3 xx)
942 (naa (math-floor (math-div (nth 2 xx) (math-pi-over-2))))
943 (nbb (math-floor (math-div (nth 3 xx) (math-pi-over-2))))
944 (span (math-sub nbb naa)))
945 (if (not (equal na nb))
946 '(intv 3 (neg (var inf var-inf)) (var inf var-inf))
947 (let ((int (math-sort-intv (nth 1 x)
948 (math-sec-raw (nth 2 xx))
949 (math-sec-raw (nth 3 xx)))))
951 (if (math-evenp (math-div (math-add naa 1) 2))
952 (math-make-intv (logior (nth 1 int) 2)
955 (math-make-intv (logior (nth 1 int) 1)
959 ((equal x '(var nan var-nan))
961 (t (calc-record-why 'scalarp x)
962 (list 'calcFunc-sec x))))
964 (defun calcFunc-csc (x)
965 (cond ((and (integerp x)
966 (eq calc-angle-mode 'deg)
967 (= (% (- x 90) 180) 0))
968 (if (= (% (- x 90) 360) 0)
972 (math-with-extra-prec 2
973 (math-csc-raw (math-to-radians (math-float x)))))
976 (math-with-extra-prec 2
977 (let* ((xx (math-to-radians (math-float (nth 1 x))))
978 (xs (math-to-radians (math-float (nth 2 x))))
979 (sc (math-sin-cos-raw xx)))
980 (if (and (math-zerop (car sc))
981 (not calc-infinite-mode))
983 (calc-record-why "*Division by zero")
984 (list 'calcFunc-csc x))
985 (math-make-sdev (math-div-float '(float 1 0) (car sc))
987 (math-mul xs (cdr sc))
988 (math-sqr (car sc)))))))
989 (math-make-sdev (calcFunc-csc (nth 1 x))
992 (calcFunc-cos (nth 1 x)))
993 (math-sqr (calcFunc-sin (nth 1 x)))))))
994 ((and (eq (car x) 'intv)
995 (math-intv-constp x))
996 (math-with-extra-prec 2
997 (let* ((xx (math-to-radians (math-float x)))
998 (na (math-floor (math-div (nth 2 xx) (math-pi))))
999 (nb (math-floor (math-div (nth 3 xx) (math-pi))))
1000 (naa (math-floor (math-div (nth 2 xx) (math-pi-over-2))))
1001 (nbb (math-floor (math-div (nth 3 xx) (math-pi-over-2))))
1002 (span (math-sub nbb naa)))
1003 (if (not (equal na nb))
1004 '(intv 3 (neg (var inf var-inf)) (var inf var-inf))
1005 (let ((int (math-sort-intv (nth 1 x)
1006 (math-csc-raw (nth 2 xx))
1007 (math-csc-raw (nth 3 xx)))))
1009 (if (math-evenp (math-div naa 2))
1010 (math-make-intv (logior (nth 1 int) 2)
1013 (math-make-intv (logior (nth 1 int) 1)
1017 ((equal x '(var nan var-nan))
1019 (t (calc-record-why 'scalarp x)
1020 (list 'calcFunc-csc x))))
1022 (defun calcFunc-cot (x) ; [N N] [Public]
1023 (cond ((and (integerp x)
1024 (if (eq calc-angle-mode 'deg)
1025 (= (% (- x 90) 180) 0)
1029 (math-with-extra-prec 2
1030 (math-cot-raw (math-to-radians (math-float x)))))
1033 (math-with-extra-prec 2
1034 (let* ((xx (math-to-radians (math-float (nth 1 x))))
1035 (xs (math-to-radians (math-float (nth 2 x))))
1036 (sc (math-sin-cos-raw xx)))
1037 (if (and (math-zerop (car sc)) (not calc-infinite-mode))
1039 (calc-record-why "*Division by zero")
1040 (list 'calcFunc-cot x))
1041 (math-make-sdev (math-div-float (cdr sc) (car sc))
1042 (math-div-float xs (math-sqr (car sc)))))))
1043 (math-make-sdev (calcFunc-cot (nth 1 x))
1045 (math-sqr (calcFunc-sin (nth 1 x)))))))
1046 ((and (eq (car x) 'intv) (math-intv-constp x))
1047 (or (math-with-extra-prec 2
1048 (let* ((xx (math-to-radians (math-float x)))
1049 (na (math-floor (math-div (nth 2 xx) (math-pi))))
1050 (nb (math-floor (math-div (nth 3 xx) (math-pi)))))
1052 (math-sort-intv (nth 1 x)
1053 (math-cot-raw (nth 2 xx))
1054 (math-cot-raw (nth 3 xx))))))
1055 '(intv 3 (neg (var inf var-inf)) (var inf var-inf))))
1056 ((equal x '(var nan var-nan))
1058 (t (calc-record-why 'scalarp x)
1059 (list 'calcFunc-cot x))))
1061 (defun math-sin-raw (x) ; [N N]
1062 (cond ((eq (car x) 'cplx)
1063 (let* ((expx (math-exp-raw (nth 2 x)))
1064 (expmx (math-div-float '(float 1 0) expx))
1065 (sc (math-sin-cos-raw (nth 1 x))))
1067 (math-mul-float (car sc)
1068 (math-mul-float (math-add-float expx expmx)
1070 (math-mul-float (cdr sc)
1071 (math-mul-float (math-sub-float expx expmx)
1073 ((eq (car x) 'polar)
1074 (math-polar (math-sin-raw (math-complex x))))
1075 ((Math-integer-negp (nth 1 x))
1076 (math-neg-float (math-sin-raw (math-neg-float x))))
1077 ((math-lessp-float '(float 7 0) x) ; avoid inf loops due to roundoff
1078 (math-sin-raw (math-mod x (math-two-pi))))
1079 (t (math-sin-raw-2 x x))))
1081 (defun math-cos-raw (x) ; [N N]
1082 (if (eq (car-safe x) 'polar)
1083 (math-polar (math-cos-raw (math-complex x)))
1084 (math-sin-raw (math-sub (math-pi-over-2) x))))
1086 (defun math-sec-raw (x) ; [N N]
1087 (cond ((eq (car x) 'cplx)
1088 (let* ((x (math-mul x '(float 1 0)))
1089 (expx (math-exp-raw (nth 2 x)))
1090 (expmx (math-div-float '(float 1 0) expx))
1091 (sh (math-mul-float (math-sub-float expx expmx) '(float 5 -1)))
1092 (ch (math-mul-float (math-add-float expx expmx) '(float 5 -1)))
1093 (sc (math-sin-cos-raw (nth 1 x)))
1095 (math-mul-float (math-sqr (car sc))
1097 (math-mul-float (math-sqr (cdr sc))
1099 (and (not (eq (nth 1 d) 0))
1101 (math-div-float (math-mul-float (cdr sc) ch) d)
1102 (math-div-float (math-mul-float (car sc) sh) d)))))
1103 ((eq (car x) 'polar)
1104 (math-polar (math-sec-raw (math-complex x))))
1106 (let ((cs (math-cos-raw x)))
1109 (math-div-float '(float 1 0) cs))))))
1111 (defun math-csc-raw (x) ; [N N]
1112 (cond ((eq (car x) 'cplx)
1113 (let* ((x (math-mul x '(float 1 0)))
1114 (expx (math-exp-raw (nth 2 x)))
1115 (expmx (math-div-float '(float 1 0) expx))
1116 (sh (math-mul-float (math-sub-float expx expmx) '(float 5 -1)))
1117 (ch (math-mul-float (math-add-float expx expmx) '(float 5 -1)))
1118 (sc (math-sin-cos-raw (nth 1 x)))
1120 (math-mul-float (math-sqr (car sc))
1122 (math-mul-float (math-sqr (cdr sc))
1124 (and (not (eq (nth 1 d) 0))
1126 (math-div-float (math-mul-float (car sc) ch) d)
1127 (math-div-float (math-mul-float (cdr sc) sh) d)))))
1128 ((eq (car x) 'polar)
1129 (math-polar (math-csc-raw (math-complex x))))
1131 (let ((sn (math-sin-raw x)))
1134 (math-div-float '(float 1 0) sn))))))
1136 (defun math-cot-raw (x) ; [N N]
1137 (cond ((eq (car x) 'cplx)
1138 (let* ((x (math-mul x '(float 1 0)))
1139 (expx (math-exp-raw (nth 2 x)))
1140 (expmx (math-div-float '(float 1 0) expx))
1141 (sh (math-mul-float (math-sub-float expx expmx) '(float 5 -1)))
1142 (ch (math-mul-float (math-add-float expx expmx) '(float 5 -1)))
1143 (sc (math-sin-cos-raw (nth 1 x)))
1147 (and (not (eq (nth 1 d) 0))
1150 (math-mul-float (car sc) (cdr sc))
1154 (math-mul-float sh ch)
1156 ((eq (car x) 'polar)
1157 (math-polar (math-cot-raw (math-complex x))))
1159 (let ((sc (math-sin-cos-raw x)))
1160 (if (eq (nth 1 (car sc)) 0)
1161 (math-div (cdr sc) 0)
1162 (math-div-float (cdr sc) (car sc)))))))
1165 ;;; This could use a smarter method: Reduce x as in math-sin-raw, then
1166 ;;; compute either sin(x) or cos(x), whichever is smaller, and compute
1167 ;;; the other using the identity sin(x)^2 + cos(x)^2 = 1.
1168 (defun math-sin-cos-raw (x) ; [F.F F] (result is (sin x . cos x))
1169 (cons (math-sin-raw x) (math-cos-raw x)))
1171 (defun math-tan-raw (x) ; [N N]
1172 (cond ((eq (car x) 'cplx)
1173 (let* ((x (math-mul x '(float 2 0)))
1174 (expx (math-exp-raw (nth 2 x)))
1175 (expmx (math-div-float '(float 1 0) expx))
1176 (sc (math-sin-cos-raw (nth 1 x)))
1177 (d (math-add-float (cdr sc)
1178 (math-mul-float (math-add-float expx expmx)
1180 (and (not (eq (nth 1 d) 0))
1182 (math-div-float (car sc) d)
1183 (math-div-float (math-mul-float (math-sub-float expx
1185 '(float 5 -1)) d)))))
1186 ((eq (car x) 'polar)
1187 (math-polar (math-tan-raw (math-complex x))))
1189 (let ((sc (math-sin-cos-raw x)))
1190 (if (eq (nth 1 (cdr sc)) 0)
1191 (math-div (car sc) 0)
1192 (math-div-float (car sc) (cdr sc)))))))
1194 (defun math-sin-raw-2 (x orgx) ; This avoids poss of inf recursion. [F F]
1195 (let ((xmpo2 (math-sub-float (math-pi-over-2) x)))
1196 (cond ((Math-integer-negp (nth 1 xmpo2))
1197 (math-neg-float (math-sin-raw-2 (math-sub-float x (math-pi))
1199 ((math-lessp-float (math-pi-over-4) x)
1200 (math-cos-raw-2 xmpo2 orgx))
1201 ((math-lessp-float x (math-neg (math-pi-over-4)))
1202 (math-neg (math-cos-raw-2 (math-add (math-pi-over-2) x) orgx)))
1203 ((math-nearly-zerop-float x orgx) '(float 0 0))
1204 (calc-symbolic-mode (signal 'inexact-result nil))
1205 (t (math-sin-series x 6 4 x (math-neg-float (math-sqr-float x)))))))
1207 (defun math-cos-raw-2 (x orgx) ; [F F]
1208 (cond ((math-nearly-zerop-float x orgx) '(float 1 0))
1209 (calc-symbolic-mode (signal 'inexact-result nil))
1210 (t (let ((xnegsqr (math-neg-float (math-sqr-float x))))
1212 (math-add-float '(float 1 0)
1213 (math-mul-float xnegsqr '(float 5 -1)))
1214 24 5 xnegsqr xnegsqr)))))
1216 (defun math-sin-series (sum nfac n x xnegsqr)
1217 (math-working "sin" sum)
1218 (let* ((nextx (math-mul-float x xnegsqr))
1219 (nextsum (math-add-float sum (math-div-float nextx
1220 (math-float nfac)))))
1221 (if (math-nearly-equal-float sum nextsum)
1223 (math-sin-series nextsum (math-mul nfac (* n (1+ n)))
1224 (+ n 2) nextx xnegsqr))))
1227 ;;; Inverse sine, cosine, tangent.
1229 (defun calcFunc-arcsin (x) ; [N N] [Public]
1231 ((and (eq x 1) (eq calc-angle-mode 'deg)) 90)
1232 ((and (eq x -1) (eq calc-angle-mode 'deg)) -90)
1233 (calc-symbolic-mode (signal 'inexact-result nil))
1235 (math-with-extra-prec 2
1236 (math-from-radians (math-arcsin-raw (math-float x)))))
1238 (math-make-sdev (calcFunc-arcsin (nth 1 x))
1242 (math-sub 1 (math-sqr (nth 1 x))))))))
1244 (math-sort-intv (nth 1 x)
1245 (calcFunc-arcsin (nth 2 x))
1246 (calcFunc-arcsin (nth 3 x))))
1247 ((equal x '(var nan var-nan))
1249 (t (calc-record-why 'numberp x)
1250 (list 'calcFunc-arcsin x))))
1252 (defun calcFunc-arccos (x) ; [N N] [Public]
1254 ((and (eq x 0) (eq calc-angle-mode 'deg)) 90)
1255 ((and (eq x -1) (eq calc-angle-mode 'deg)) 180)
1256 (calc-symbolic-mode (signal 'inexact-result nil))
1258 (math-with-extra-prec 2
1259 (math-from-radians (math-arccos-raw (math-float x)))))
1261 (math-make-sdev (calcFunc-arccos (nth 1 x))
1265 (math-sub 1 (math-sqr (nth 1 x))))))))
1267 (math-sort-intv (nth 1 x)
1268 (calcFunc-arccos (nth 2 x))
1269 (calcFunc-arccos (nth 3 x))))
1270 ((equal x '(var nan var-nan))
1272 (t (calc-record-why 'numberp x)
1273 (list 'calcFunc-arccos x))))
1275 (defun calcFunc-arctan (x) ; [N N] [Public]
1277 ((and (eq x 1) (eq calc-angle-mode 'deg)) 45)
1278 ((and (eq x -1) (eq calc-angle-mode 'deg)) -45)
1280 (math-with-extra-prec 2
1281 (math-from-radians (math-arctan-raw (math-float x)))))
1283 (math-make-sdev (calcFunc-arctan (nth 1 x))
1286 (math-add 1 (math-sqr (nth 1 x)))))))
1288 (math-sort-intv (nth 1 x)
1289 (calcFunc-arctan (nth 2 x))
1290 (calcFunc-arctan (nth 3 x))))
1291 ((equal x '(var inf var-inf))
1292 (math-quarter-circle t))
1293 ((equal x '(neg (var inf var-inf)))
1294 (math-neg (math-quarter-circle t)))
1295 ((equal x '(var nan var-nan))
1297 (t (calc-record-why 'numberp x)
1298 (list 'calcFunc-arctan x))))
1300 (defun math-arcsin-raw (x) ; [N N]
1301 (let ((a (math-sqrt-raw (math-sub '(float 1 0) (math-sqr x)))))
1302 (if (or (memq (car x) '(cplx polar))
1303 (memq (car a) '(cplx polar)))
1304 (math-with-extra-prec 2 ; use extra precision for difficult case
1305 (math-mul '(cplx 0 -1)
1306 (math-ln-raw (math-add (math-mul '(cplx 0 1) x) a))))
1307 (math-arctan2-raw x a))))
1309 (defun math-arccos-raw (x) ; [N N]
1310 (math-sub (math-pi-over-2) (math-arcsin-raw x)))
1312 (defun math-arctan-raw (x) ; [N N]
1313 (cond ((memq (car x) '(cplx polar))
1314 (math-with-extra-prec 2 ; extra-extra
1316 (math-ln-raw (math-add 1 (math-mul '(cplx 0 1) x)))
1317 (math-ln-raw (math-add 1 (math-mul '(cplx 0 -1) x))))
1319 ((Math-integer-negp (nth 1 x))
1320 (math-neg-float (math-arctan-raw (math-neg-float x))))
1322 (calc-symbolic-mode (signal 'inexact-result nil))
1323 ((math-equal-int x 1) (math-pi-over-4))
1324 ((math-equal-int x -1) (math-neg (math-pi-over-4)))
1325 ((math-lessp-float '(float 414214 -6) x) ; if x > sqrt(2) - 1, reduce
1326 (if (math-lessp-float '(float 1 0) x)
1327 (math-sub-float (math-mul-float (math-pi) '(float 5 -1))
1328 (math-arctan-raw (math-div-float '(float 1 0) x)))
1329 (math-sub-float (math-mul-float (math-pi) '(float 25 -2))
1330 (math-arctan-raw (math-div-float
1331 (math-sub-float '(float 1 0) x)
1332 (math-add-float '(float 1 0)
1334 (t (math-arctan-series x 3 x (math-neg-float (math-sqr-float x))))))
1336 (defun math-arctan-series (sum n x xnegsqr)
1337 (math-working "arctan" sum)
1338 (let* ((nextx (math-mul-float x xnegsqr))
1339 (nextsum (math-add-float sum (math-div-float nextx (math-float n)))))
1340 (if (math-nearly-equal-float sum nextsum)
1342 (math-arctan-series nextsum (+ n 2) nextx xnegsqr))))
1344 (defun calcFunc-arctan2 (y x) ; [F R R] [Public]
1347 (math-with-extra-prec 2
1348 (math-from-radians (math-arctan2-raw (math-float y)
1350 (calc-record-why 'anglep x)
1351 (list 'calcFunc-arctan2 y x))
1352 (if (and (or (math-infinitep x) (math-anglep x))
1353 (or (math-infinitep y) (math-anglep y)))
1368 (calcFunc-arctan2 y x)
1369 '(var nan var-nan)))
1370 (calc-record-why 'anglep y)
1371 (list 'calcFunc-arctan2 y x))))
1373 (defun math-arctan2-raw (y x) ; [F R R]
1374 (cond ((math-zerop y)
1375 (if (math-negp x) (math-pi)
1376 (if (or (math-floatp x) (math-floatp y)) '(float 0 0) 0)))
1380 (math-neg (math-pi-over-2))))
1382 (math-arctan-raw (math-div-float y x)))
1384 (math-add-float (math-arctan-raw (math-div-float y x))
1387 (math-sub-float (math-arctan-raw (math-div-float y x))
1390 (defun calcFunc-arcsincos (x) ; [V N] [Public]
1391 (if (and (Math-vectorp x)
1393 (calcFunc-arctan2 (nth 2 x) (nth 1 x))
1394 (math-reject-arg x "*Two-element vector expected")))
1398 ;;; Exponential function.
1400 (defun calcFunc-exp (x) ; [N N] [Public]
1402 ((and (memq x '(1 -1)) calc-symbolic-mode)
1403 (if (eq x 1) '(var e var-e) (math-div 1 '(var e var-e))))
1405 (math-with-extra-prec 2 (math-exp-raw (math-float x))))
1406 ((eq (car-safe x) 'sdev)
1407 (let ((ex (calcFunc-exp (nth 1 x))))
1408 (math-make-sdev ex (math-mul (nth 2 x) ex))))
1409 ((eq (car-safe x) 'intv)
1410 (math-make-intv (nth 1 x) (calcFunc-exp (nth 2 x))
1411 (calcFunc-exp (nth 3 x))))
1412 ((equal x '(var inf var-inf))
1414 ((equal x '(neg (var inf var-inf)))
1416 ((equal x '(var nan var-nan))
1418 (t (calc-record-why 'numberp x)
1419 (list 'calcFunc-exp x))))
1421 (defun calcFunc-expm1 (x) ; [N N] [Public]
1423 ((math-zerop x) '(float 0 0))
1424 (calc-symbolic-mode (signal 'inexact-result nil))
1426 (math-with-extra-prec 2
1427 (let ((x (math-float x)))
1428 (if (and (eq (car x) 'float)
1429 (math-lessp-float x '(float 1 0))
1430 (math-lessp-float '(float -1 0) x))
1431 (math-exp-minus-1-raw x)
1432 (math-add (math-exp-raw x) -1)))))
1433 ((eq (car-safe x) 'sdev)
1435 (let ((ex (calcFunc-expm1 (nth 1 x))))
1436 (math-make-sdev ex (math-mul (nth 2 x) (math-add ex 1))))
1437 (math-make-sdev (calcFunc-expm1 (nth 1 x))
1438 (math-mul (nth 2 x) (calcFunc-exp (nth 1 x))))))
1439 ((eq (car-safe x) 'intv)
1440 (math-make-intv (nth 1 x)
1441 (calcFunc-expm1 (nth 2 x))
1442 (calcFunc-expm1 (nth 3 x))))
1443 ((equal x '(var inf var-inf))
1445 ((equal x '(neg (var inf var-inf)))
1447 ((equal x '(var nan var-nan))
1449 (t (calc-record-why 'numberp x)
1450 (list 'calcFunc-expm1 x))))
1452 (defun calcFunc-exp10 (x) ; [N N] [Public]
1455 (math-pow '(float 1 1) x)))
1457 (defun math-exp-raw (x) ; [N N]
1458 (cond ((math-zerop x) '(float 1 0))
1459 (calc-symbolic-mode (signal 'inexact-result nil))
1461 (let ((expx (math-exp-raw (nth 1 x)))
1462 (sc (math-sin-cos-raw (nth 2 x))))
1464 (math-mul-float expx (cdr sc))
1465 (math-mul-float expx (car sc)))))
1466 ((eq (car x) 'polar)
1467 (let ((xc (math-complex x)))
1469 (math-exp-raw (nth 1 xc))
1470 (math-from-radians (nth 2 xc)))))
1471 ((math-use-emacs-fn 'exp x))
1472 ((or (math-lessp-float '(float 5 -1) x)
1473 (math-lessp-float x '(float -5 -1)))
1474 (if (math-lessp-float '(float 921035 1) x)
1476 (if (math-lessp-float x '(float -921035 1))
1478 (let* ((two-x (math-mul-float x '(float 2 0)))
1479 (hint (math-scale-int (nth 1 two-x) (nth 2 two-x)))
1480 (hfrac (math-sub-float x (math-mul-float (math-float hint)
1482 (math-mul-float (math-ipow (math-sqrt-e) hint)
1483 (math-add-float '(float 1 0)
1484 (math-exp-minus-1-raw hfrac)))))
1485 (t (math-add-float '(float 1 0) (math-exp-minus-1-raw x)))))
1487 (defun math-exp-minus-1-raw (x) ; [F F]
1488 (math-exp-series x 2 3 x x))
1490 (defun math-exp-series (sum nfac n xpow x)
1491 (math-working "exp" sum)
1492 (let* ((nextx (math-mul-float xpow x))
1493 (nextsum (math-add-float sum (math-div-float nextx
1494 (math-float nfac)))))
1495 (if (math-nearly-equal-float sum nextsum)
1497 (math-exp-series nextsum (math-mul nfac n) (1+ n) nextx x))))
1503 (defun calcFunc-ln (x) ; [N N] [Public]
1504 (cond ((math-zerop x)
1505 (if calc-infinite-mode
1506 '(neg (var inf var-inf))
1507 (math-reject-arg x "*Logarithm of zero")))
1510 (math-with-extra-prec 2 (math-ln-raw (math-float x))))
1511 ((eq (car-safe x) 'sdev)
1512 (math-make-sdev (calcFunc-ln (nth 1 x))
1513 (math-div (nth 2 x) (nth 1 x))))
1514 ((and (eq (car-safe x) 'intv) (or (Math-posp (nth 2 x))
1515 (Math-zerop (nth 2 x))
1516 (not (math-intv-constp x))))
1517 (let ((calc-infinite-mode t))
1518 (math-make-intv (nth 1 x) (calcFunc-ln (nth 2 x))
1519 (calcFunc-ln (nth 3 x)))))
1520 ((equal x '(var e var-e))
1522 ((and (eq (car-safe x) '^)
1523 (equal (nth 1 x) '(var e var-e))
1524 (math-known-realp (nth 2 x)))
1527 (if (equal x '(var nan var-nan))
1529 '(var inf var-inf)))
1530 (t (calc-record-why 'numberp x)
1531 (list 'calcFunc-ln x))))
1533 (defun calcFunc-log10 (x) ; [N N] [Public]
1534 (cond ((math-equal-int x 1)
1535 (if (math-floatp x) '(float 0 0) 0))
1536 ((and (Math-integerp x)
1538 (let ((res (math-integer-log x 10)))
1540 (setq x (cdr res)))))
1542 ((and (eq (car-safe x) 'frac)
1544 (let ((res (math-integer-log (nth 2 x) 10)))
1546 (setq x (- (cdr res))))))
1549 (if calc-infinite-mode
1550 '(neg (var inf var-inf))
1551 (math-reject-arg x "*Logarithm of zero")))
1552 (calc-symbolic-mode (signal 'inexact-result nil))
1554 (math-with-extra-prec 2
1555 (let ((xf (math-float x)))
1556 (if (eq (nth 1 xf) 0)
1557 (math-reject-arg x "*Logarithm of zero"))
1558 (if (Math-integer-posp (nth 1 xf))
1559 (if (eq (nth 1 xf) 1) ; log10(1*10^n) = n
1560 (math-float (nth 2 xf))
1561 (let ((xdigs (1- (math-numdigs (nth 1 xf)))))
1563 (math-div-float (math-ln-raw-2
1564 (list 'float (nth 1 xf) (- xdigs)))
1566 (math-float (+ (nth 2 xf) xdigs)))))
1567 (math-div (calcFunc-ln xf) (math-ln-10))))))
1568 ((eq (car-safe x) 'sdev)
1569 (math-make-sdev (calcFunc-log10 (nth 1 x))
1571 (math-mul (nth 1 x) (math-ln-10)))))
1572 ((and (eq (car-safe x) 'intv) (or (Math-posp (nth 2 x))
1573 (not (math-intv-constp x))))
1574 (math-make-intv (nth 1 x)
1575 (calcFunc-log10 (nth 2 x))
1576 (calcFunc-log10 (nth 3 x))))
1578 (if (equal x '(var nan var-nan))
1580 '(var inf var-inf)))
1581 (t (calc-record-why 'numberp x)
1582 (list 'calcFunc-log10 x))))
1584 (defun calcFunc-log (x &optional b) ; [N N N] [Public]
1585 (cond ((or (null b) (equal b '(var e var-e)))
1586 (math-normalize (list 'calcFunc-ln x)))
1587 ((or (eq b 10) (equal b '(float 1 1)))
1588 (math-normalize (list 'calcFunc-log10 x)))
1590 (if calc-infinite-mode
1591 (math-div (calcFunc-ln x) (calcFunc-ln b))
1592 (math-reject-arg x "*Logarithm of zero")))
1594 (if calc-infinite-mode
1595 (math-div (calcFunc-ln x) (calcFunc-ln b))
1596 (math-reject-arg b "*Logarithm of zero")))
1597 ((math-equal-int b 1)
1598 (if calc-infinite-mode
1599 (math-div (calcFunc-ln x) 0)
1600 (math-reject-arg b "*Logarithm base one")))
1601 ((math-equal-int x 1)
1602 (if (math-floatp b) '(float 0 0) 0))
1603 ((and (Math-ratp x) (Math-ratp b)
1604 (math-posp x) (math-posp b)
1605 (let* ((sign 1) (inv nil)
1606 (xx (if (Math-lessp 1 x)
1610 (bb (if (Math-lessp 1 b)
1612 (setq sign (- sign))
1614 (res (if (Math-lessp xx bb)
1615 (setq inv (math-integer-log bb xx))
1616 (math-integer-log xx bb))))
1619 (math-div 1 (* sign (cdr res)))
1620 (* sign (cdr res)))))))
1622 (calc-symbolic-mode (signal 'inexact-result nil))
1623 ((and (Math-numberp x) (Math-numberp b))
1624 (math-with-extra-prec 2
1625 (math-div (math-ln-raw (math-float x))
1626 (math-log-base-raw b))))
1627 ((and (eq (car-safe x) 'sdev)
1629 (math-make-sdev (calcFunc-log (nth 1 x) b)
1632 (math-log-base-raw b)))))
1633 ((and (eq (car-safe x) 'intv) (or (Math-posp (nth 2 x))
1634 (not (math-intv-constp x)))
1636 (math-make-intv (nth 1 x)
1637 (calcFunc-log (nth 2 x) b)
1638 (calcFunc-log (nth 3 x) b)))
1639 ((or (eq (car-safe x) 'intv) (eq (car-safe b) 'intv))
1640 (math-div (calcFunc-ln x) (calcFunc-ln b)))
1641 ((or (math-infinitep x)
1643 (math-div (calcFunc-ln x) (calcFunc-ln b)))
1644 (t (if (Math-numberp b)
1645 (calc-record-why 'numberp x)
1646 (calc-record-why 'numberp b))
1647 (list 'calcFunc-log x b))))
1649 (defun calcFunc-alog (x &optional b)
1650 (cond ((or (null b) (equal b '(var e var-e)))
1651 (math-normalize (list 'calcFunc-exp x)))
1652 (t (math-pow b x))))
1654 (defun calcFunc-ilog (x b)
1655 (if (and (math-natnump x) (not (eq x 0))
1656 (math-natnump b) (not (eq b 0)))
1658 (math-reject-arg x "*Logarithm base one")
1659 (if (Math-natnum-lessp x b)
1661 (cdr (math-integer-log x b))))
1662 (math-floor (calcFunc-log x b))))
1664 (defun math-integer-log (x b)
1665 (let ((pows (list b))
1669 (while (not (Math-lessp x pow))
1670 (setq pows (cons pow pows)
1671 pow (math-sqr pow)))
1672 (setq n (lsh 1 (1- (length pows)))
1675 (while (and (setq pows (cdr pows))
1678 next (math-mul pow (car pows)))
1679 (or (Math-lessp x next)
1682 (cons (equal pow x) sum)))
1685 (defvar math-log-base-cache nil)
1686 (defun math-log-base-raw (b) ; [N N]
1687 (if (not (and (equal (car math-log-base-cache) b)
1688 (eq (nth 1 math-log-base-cache) calc-internal-prec)))
1689 (setq math-log-base-cache (list b calc-internal-prec
1690 (math-ln-raw (math-float b)))))
1691 (nth 2 math-log-base-cache))
1693 (defun calcFunc-lnp1 (x) ; [N N] [Public]
1694 (cond ((Math-equal-int x -1)
1695 (if calc-infinite-mode
1696 '(neg (var inf var-inf))
1697 (math-reject-arg x "*Logarithm of zero")))
1699 ((math-zerop x) '(float 0 0))
1700 (calc-symbolic-mode (signal 'inexact-result nil))
1702 (math-with-extra-prec 2
1703 (let ((x (math-float x)))
1704 (if (and (eq (car x) 'float)
1705 (math-lessp-float x '(float 5 -1))
1706 (math-lessp-float '(float -5 -1) x))
1707 (math-ln-plus-1-raw x)
1708 (math-ln-raw (math-add-float x '(float 1 0)))))))
1709 ((eq (car-safe x) 'sdev)
1710 (math-make-sdev (calcFunc-lnp1 (nth 1 x))
1711 (math-div (nth 2 x) (math-add (nth 1 x) 1))))
1712 ((and (eq (car-safe x) 'intv) (or (Math-posp (nth 2 x))
1713 (not (math-intv-constp x))))
1714 (math-make-intv (nth 1 x)
1715 (calcFunc-lnp1 (nth 2 x))
1716 (calcFunc-lnp1 (nth 3 x))))
1718 (if (equal x '(var nan var-nan))
1720 '(var inf var-inf)))
1721 (t (calc-record-why 'numberp x)
1722 (list 'calcFunc-lnp1 x))))
1724 (defun math-ln-raw (x) ; [N N] --- must be float format!
1725 (cond ((eq (car-safe x) 'cplx)
1727 (math-mul-float (math-ln-raw
1728 (math-add-float (math-sqr-float (nth 1 x))
1729 (math-sqr-float (nth 2 x))))
1731 (math-arctan2-raw (nth 2 x) (nth 1 x))))
1732 ((eq (car x) 'polar)
1733 (math-polar (list 'cplx
1734 (math-ln-raw (nth 1 x))
1735 (math-to-radians (nth 2 x)))))
1736 ((Math-equal-int x 1)
1738 (calc-symbolic-mode (signal 'inexact-result nil))
1739 ((math-posp (nth 1 x)) ; positive and real
1740 (let ((xdigs (1- (math-numdigs (nth 1 x)))))
1741 (math-add-float (math-ln-raw-2 (list 'float (nth 1 x) (- xdigs)))
1742 (math-mul-float (math-float (+ (nth 2 x) xdigs))
1745 (math-reject-arg x "*Logarithm of zero"))
1746 ((eq calc-complex-mode 'polar) ; negative and real
1748 (list 'cplx ; negative and real
1749 (math-ln-raw (math-neg-float x))
1751 (t (list 'cplx ; negative and real
1752 (math-ln-raw (math-neg-float x))
1755 (defun math-ln-raw-2 (x) ; [F F]
1756 (cond ((math-lessp-float '(float 14 -1) x)
1757 (math-add-float (math-ln-raw-2 (math-mul-float x '(float 5 -1)))
1759 (t ; now .7 < x <= 1.4
1760 (math-ln-raw-3 (math-div-float (math-sub-float x '(float 1 0))
1761 (math-add-float x '(float 1 0)))))))
1763 (defun math-ln-raw-3 (x) ; [F F]
1764 (math-mul-float (math-ln-raw-series x 3 x (math-sqr-float x))
1767 ;;; Compute ln((1+x)/(1-x))
1768 (defun math-ln-raw-series (sum n x xsqr)
1769 (math-working "log" sum)
1770 (let* ((nextx (math-mul-float x xsqr))
1771 (nextsum (math-add-float sum (math-div-float nextx (math-float n)))))
1772 (if (math-nearly-equal-float sum nextsum)
1774 (math-ln-raw-series nextsum (+ n 2) nextx xsqr))))
1776 (defun math-ln-plus-1-raw (x)
1777 (math-lnp1-series x 2 x (math-neg x)))
1779 (defun math-lnp1-series (sum n xpow x)
1780 (math-working "lnp1" sum)
1781 (let* ((nextx (math-mul-float xpow x))
1782 (nextsum (math-add-float sum (math-div-float nextx (math-float n)))))
1783 (if (math-nearly-equal-float sum nextsum)
1785 (math-lnp1-series nextsum (1+ n) nextx x))))
1787 (defconst math-approx-ln-10
1789 (math-read-number-simple "2.302585092994045684018"))
1790 "An approximation for ln(10).")
1792 (math-defcache math-ln-10 math-approx-ln-10
1793 (math-ln-raw-2 '(float 1 1)))
1795 (defconst math-approx-ln-2
1797 (math-read-number-simple "0.693147180559945309417"))
1798 "An approximation for ln(2).")
1800 (math-defcache math-ln-2 math-approx-ln-2
1801 (math-ln-raw-3 (math-float '(frac 1 3))))
1805 ;;; Hyperbolic functions.
1807 (defun calcFunc-sinh (x) ; [N N] [Public]
1809 (math-expand-formulas
1811 (list '/ (list '- (list 'calcFunc-exp x)
1812 (list 'calcFunc-exp (list 'neg x))) 2)))
1814 (if calc-symbolic-mode (signal 'inexact-result nil))
1815 (math-with-extra-prec 2
1816 (let ((expx (math-exp-raw (math-float x))))
1817 (math-mul (math-add expx (math-div -1 expx)) '(float 5 -1)))))
1818 ((eq (car-safe x) 'sdev)
1819 (math-make-sdev (calcFunc-sinh (nth 1 x))
1820 (math-mul (nth 2 x) (calcFunc-cosh (nth 1 x)))))
1822 (math-sort-intv (nth 1 x)
1823 (calcFunc-sinh (nth 2 x))
1824 (calcFunc-sinh (nth 3 x))))
1825 ((or (equal x '(var inf var-inf))
1826 (equal x '(neg (var inf var-inf)))
1827 (equal x '(var nan var-nan)))
1829 (t (calc-record-why 'numberp x)
1830 (list 'calcFunc-sinh x))))
1831 (put 'calcFunc-sinh 'math-expandable t)
1833 (defun calcFunc-cosh (x) ; [N N] [Public]
1835 (math-expand-formulas
1837 (list '/ (list '+ (list 'calcFunc-exp x)
1838 (list 'calcFunc-exp (list 'neg x))) 2)))
1840 (if calc-symbolic-mode (signal 'inexact-result nil))
1841 (math-with-extra-prec 2
1842 (let ((expx (math-exp-raw (math-float x))))
1843 (math-mul (math-add expx (math-div 1 expx)) '(float 5 -1)))))
1844 ((eq (car-safe x) 'sdev)
1845 (math-make-sdev (calcFunc-cosh (nth 1 x))
1847 (calcFunc-sinh (nth 1 x)))))
1848 ((and (eq (car x) 'intv) (math-intv-constp x))
1849 (setq x (math-abs x))
1850 (math-sort-intv (nth 1 x)
1851 (calcFunc-cosh (nth 2 x))
1852 (calcFunc-cosh (nth 3 x))))
1853 ((or (equal x '(var inf var-inf))
1854 (equal x '(neg (var inf var-inf)))
1855 (equal x '(var nan var-nan)))
1857 (t (calc-record-why 'numberp x)
1858 (list 'calcFunc-cosh x))))
1859 (put 'calcFunc-cosh 'math-expandable t)
1861 (defun calcFunc-tanh (x) ; [N N] [Public]
1863 (math-expand-formulas
1865 (let ((expx (list 'calcFunc-exp x))
1866 (expmx (list 'calcFunc-exp (list 'neg x))))
1868 (list '/ (list '- expx expmx) (list '+ expx expmx))))))
1870 (if calc-symbolic-mode (signal 'inexact-result nil))
1871 (math-with-extra-prec 2
1872 (let* ((expx (calcFunc-exp (math-float x)))
1873 (expmx (math-div 1 expx)))
1874 (math-div (math-sub expx expmx)
1875 (math-add expx expmx)))))
1876 ((eq (car-safe x) 'sdev)
1877 (math-make-sdev (calcFunc-tanh (nth 1 x))
1879 (math-sqr (calcFunc-cosh (nth 1 x))))))
1881 (math-sort-intv (nth 1 x)
1882 (calcFunc-tanh (nth 2 x))
1883 (calcFunc-tanh (nth 3 x))))
1884 ((equal x '(var inf var-inf))
1886 ((equal x '(neg (var inf var-inf)))
1888 ((equal x '(var nan var-nan))
1890 (t (calc-record-why 'numberp x)
1891 (list 'calcFunc-tanh x))))
1892 (put 'calcFunc-tanh 'math-expandable t)
1894 (defun calcFunc-sech (x) ; [N N] [Public]
1896 (math-expand-formulas
1898 (list '/ 2 (list '+ (list 'calcFunc-exp x)
1899 (list 'calcFunc-exp (list 'neg x))))))
1901 (if calc-symbolic-mode (signal 'inexact-result nil))
1902 (math-with-extra-prec 2
1903 (let ((expx (math-exp-raw (math-float x))))
1904 (math-div '(float 2 0) (math-add expx (math-div 1 expx))))))
1905 ((eq (car-safe x) 'sdev)
1906 (math-make-sdev (calcFunc-sech (nth 1 x))
1908 (math-mul (calcFunc-sech (nth 1 x))
1909 (calcFunc-tanh (nth 1 x))))))
1910 ((and (eq (car x) 'intv) (math-intv-constp x))
1911 (setq x (math-abs x))
1912 (math-sort-intv (nth 1 x)
1913 (calcFunc-sech (nth 2 x))
1914 (calcFunc-sech (nth 3 x))))
1915 ((or (equal x '(var inf var-inf))
1916 (equal x '(neg (var inf var-inf))))
1918 ((equal x '(var nan var-nan))
1920 (t (calc-record-why 'numberp x)
1921 (list 'calcFunc-sech x))))
1922 (put 'calcFunc-sech 'math-expandable t)
1924 (defun calcFunc-csch (x) ; [N N] [Public]
1925 (cond ((eq x 0) (math-div 1 0))
1926 (math-expand-formulas
1928 (list '/ 2 (list '- (list 'calcFunc-exp x)
1929 (list 'calcFunc-exp (list 'neg x))))))
1931 (if calc-symbolic-mode (signal 'inexact-result nil))
1932 (math-with-extra-prec 2
1933 (let ((expx (math-exp-raw (math-float x))))
1934 (math-div '(float 2 0) (math-add expx (math-div -1 expx))))))
1935 ((eq (car-safe x) 'sdev)
1936 (math-make-sdev (calcFunc-csch (nth 1 x))
1938 (math-mul (calcFunc-csch (nth 1 x))
1939 (calcFunc-coth (nth 1 x))))))
1941 (if (and (Math-negp (nth 2 x))
1942 (Math-posp (nth 3 x)))
1943 '(intv 3 (neg (var inf var-inf)) (var inf var-inf))
1944 (math-sort-intv (nth 1 x)
1945 (calcFunc-csch (nth 2 x))
1946 (calcFunc-csch (nth 3 x)))))
1947 ((or (equal x '(var inf var-inf))
1948 (equal x '(neg (var inf var-inf))))
1950 ((equal x '(var nan var-nan))
1952 (t (calc-record-why 'numberp x)
1953 (list 'calcFunc-csch x))))
1954 (put 'calcFunc-csch 'math-expandable t)
1956 (defun calcFunc-coth (x) ; [N N] [Public]
1957 (cond ((eq x 0) (math-div 1 0))
1958 (math-expand-formulas
1960 (let ((expx (list 'calcFunc-exp x))
1961 (expmx (list 'calcFunc-exp (list 'neg x))))
1963 (list '/ (list '+ expx expmx) (list '- expx expmx))))))
1965 (if calc-symbolic-mode (signal 'inexact-result nil))
1966 (math-with-extra-prec 2
1967 (let* ((expx (calcFunc-exp (math-float x)))
1968 (expmx (math-div 1 expx)))
1969 (math-div (math-add expx expmx)
1970 (math-sub expx expmx)))))
1971 ((eq (car-safe x) 'sdev)
1972 (math-make-sdev (calcFunc-coth (nth 1 x))
1974 (math-sqr (calcFunc-sinh (nth 1 x))))))
1976 (if (and (Math-negp (nth 2 x))
1977 (Math-posp (nth 3 x)))
1978 '(intv 3 (neg (var inf var-inf)) (var inf var-inf))
1979 (math-sort-intv (nth 1 x)
1980 (calcFunc-coth (nth 2 x))
1981 (calcFunc-coth (nth 3 x)))))
1982 ((equal x '(var inf var-inf))
1984 ((equal x '(neg (var inf var-inf)))
1986 ((equal x '(var nan var-nan))
1988 (t (calc-record-why 'numberp x)
1989 (list 'calcFunc-coth x))))
1990 (put 'calcFunc-coth 'math-expandable t)
1992 (defun calcFunc-arcsinh (x) ; [N N] [Public]
1994 (math-expand-formulas
1996 (list 'calcFunc-ln (list '+ x (list 'calcFunc-sqrt
1997 (list '+ (list '^ x 2) 1))))))
1999 (if calc-symbolic-mode (signal 'inexact-result nil))
2000 (math-with-extra-prec 2
2001 (math-ln-raw (math-add x (math-sqrt-raw (math-add (math-sqr x)
2003 ((eq (car-safe x) 'sdev)
2004 (math-make-sdev (calcFunc-arcsinh (nth 1 x))
2007 (math-add (math-sqr (nth 1 x)) 1)))))
2009 (math-sort-intv (nth 1 x)
2010 (calcFunc-arcsinh (nth 2 x))
2011 (calcFunc-arcsinh (nth 3 x))))
2012 ((or (equal x '(var inf var-inf))
2013 (equal x '(neg (var inf var-inf)))
2014 (equal x '(var nan var-nan)))
2016 (t (calc-record-why 'numberp x)
2017 (list 'calcFunc-arcsinh x))))
2018 (put 'calcFunc-arcsinh 'math-expandable t)
2020 (defun calcFunc-arccosh (x) ; [N N] [Public]
2022 ((and (eq x -1) calc-symbolic-mode)
2024 ((and (eq x 0) calc-symbolic-mode)
2025 (math-div (math-mul '(var pi var-pi) '(var i var-i)) 2))
2026 (math-expand-formulas
2028 (list 'calcFunc-ln (list '+ x (list 'calcFunc-sqrt
2029 (list '- (list '^ x 2) 1))))))
2031 (if calc-symbolic-mode (signal 'inexact-result nil))
2032 (if (Math-equal-int x -1)
2033 (math-imaginary (math-pi))
2034 (math-with-extra-prec 2
2035 (if (or t ; need to do this even in the real case!
2036 (memq (car-safe x) '(cplx polar)))
2037 (let ((xp1 (math-add 1 x))) ; this gets the branch cuts right
2039 (math-add x (math-mul xp1
2046 (math-add x (math-sqrt-raw (math-add (math-sqr x)
2047 '(float -1 0)))))))))
2048 ((eq (car-safe x) 'sdev)
2049 (math-make-sdev (calcFunc-arccosh (nth 1 x))
2052 (math-add (math-sqr (nth 1 x)) -1)))))
2054 (math-sort-intv (nth 1 x)
2055 (calcFunc-arccosh (nth 2 x))
2056 (calcFunc-arccosh (nth 3 x))))
2057 ((or (equal x '(var inf var-inf))
2058 (equal x '(neg (var inf var-inf)))
2059 (equal x '(var nan var-nan)))
2061 (t (calc-record-why 'numberp x)
2062 (list 'calcFunc-arccosh x))))
2063 (put 'calcFunc-arccosh 'math-expandable t)
2065 (defun calcFunc-arctanh (x) ; [N N] [Public]
2067 ((and (Math-equal-int x 1) calc-infinite-mode)
2069 ((and (Math-equal-int x -1) calc-infinite-mode)
2070 '(neg (var inf var-inf)))
2071 (math-expand-formulas
2073 (list 'calcFunc-ln (list '+ 1 x))
2074 (list 'calcFunc-ln (list '- 1 x))) 2))
2076 (if calc-symbolic-mode (signal 'inexact-result nil))
2077 (math-with-extra-prec 2
2078 (if (or (memq (car-safe x) '(cplx polar))
2080 (math-mul (math-sub (math-ln-raw (math-add '(float 1 0) x))
2081 (math-ln-raw (math-sub '(float 1 0) x)))
2083 (if (and (math-equal-int x 1) calc-infinite-mode)
2085 (if (and (math-equal-int x -1) calc-infinite-mode)
2086 '(neg (var inf var-inf))
2087 (math-mul (math-ln-raw (math-div (math-add '(float 1 0) x)
2090 ((eq (car-safe x) 'sdev)
2091 (math-make-sdev (calcFunc-arctanh (nth 1 x))
2093 (math-sub 1 (math-sqr (nth 1 x))))))
2095 (math-sort-intv (nth 1 x)
2096 (calcFunc-arctanh (nth 2 x))
2097 (calcFunc-arctanh (nth 3 x))))
2098 ((equal x '(var nan var-nan))
2100 (t (calc-record-why 'numberp x)
2101 (list 'calcFunc-arctanh x))))
2102 (put 'calcFunc-arctanh 'math-expandable t)
2105 ;;; Convert A from HMS or degrees to radians.
2106 (defun calcFunc-rad (a) ; [R R] [Public]
2107 (cond ((or (Math-numberp a)
2109 (math-with-extra-prec 2
2110 (math-mul a (math-pi-over-180))))
2112 (math-from-hms a 'rad))
2114 (math-make-sdev (calcFunc-rad (nth 1 a))
2115 (calcFunc-rad (nth 2 a))))
2116 (math-expand-formulas
2117 (math-div (math-mul a '(var pi var-pi)) 180))
2118 ((math-infinitep a) a)
2119 (t (list 'calcFunc-rad a))))
2120 (put 'calcFunc-rad 'math-expandable t)
2122 ;;; Convert A from HMS or radians to degrees.
2123 (defun calcFunc-deg (a) ; [R R] [Public]
2124 (cond ((or (Math-numberp a)
2126 (math-with-extra-prec 2
2127 (math-div a (math-pi-over-180))))
2129 (math-from-hms a 'deg))
2131 (math-make-sdev (calcFunc-deg (nth 1 a))
2132 (calcFunc-deg (nth 2 a))))
2133 (math-expand-formulas
2134 (math-div (math-mul 180 a) '(var pi var-pi)))
2135 ((math-infinitep a) a)
2136 (t (list 'calcFunc-deg a))))
2137 (put 'calcFunc-deg 'math-expandable t)
2139 (provide 'calc-math)
2141 ;;; arch-tag: c7367e8e-d0b8-4f70-8577-2fb3f31dbb4c
2142 ;;; calc-math.el ends here