src/zero.lisp

   1 ;;; -*-  Mode: Lisp; Package: Maxima; Syntax: Common-Lisp; Base: 10 -*- ;;;;
   2 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
   3 ;;;     The data in this file contains enhancements.                   ;;;;;
   4 ;;;                                                                    ;;;;;
   5 ;;;  Copyright (c) 1984,1987 by William Schelter,University of Texas   ;;;;;
   6 ;;;     All rights reserved                                            ;;;;;
   7 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
   8 ;;;     (c) Copyright 1982 Massachusetts Institute of Technology         ;;;
   9 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  10
  11 (in-package :maxima)
  12
  13 (macsyma-module zero)
  14
  15 (defmfun $zeroequiv (exp var)
  16   (declare (special var ))
  17   (prog (r s v varlist genvar)
  18      (declare (special s v))
  19      (setq exp (specrepcheck exp))
  20      (setq r (let ($listconstvars) ($listofvars exp)))
  21      (if (and (cdr r) (or (cddr r) (not (alike1 (cadr r) var))))
  22          (return '$dontknow))
  23      (setq exp ($exponentialize exp))
  24      (setq r (sdiff exp var))
  25      (if (isinop r '%derivative) (return '$dontknow))
  26      ($rat r)
  27      (setq r ($rat exp))
  28      (setq s (car r))
  29      (setq v (ratnumerator (cdr r)))
  30      (return (zeroequiv1 v))))
  31
  32 (defun zeroequiv1 (v)
  33   (declare (special var v s))
  34   (prog (v1 v2 coeff deg)
  35      (declare (special v1 v2))
  36      (if (atom v) (return (equal v 0)))
  37      coeffloop (if (null (cdr v)) (return t))
  38      (setq deg (cadr v))
  39      (if (equal deg 0) (return (zeroequiv1 (caddr v))))
  40      (setq coeff (caddr v))
  41      (when (zeroequiv1 coeff)
  42        (setq v (cons (car v) (cdddr v)))
  43        (go coeffloop))
  44      (setq v1 ($rat (sdiff (ratdisrep (cons s (cons v (caddr v)))) var)))
  45      (setq v2 (cadr ($rat (ratdisrep v1))))
  46      (if (equal (pdegree v2 (car v)) (cadr v))
  47          (return (zeroequiv2 v)))
  48      (if (< (pdegree v2 (car v)) (cadr v))
  49          (return (if (zeroequiv1 v2) (zeroequiv2 v))))
  50      (return '$dontknow)))
  51
  52 (defun zeroequiv2 (v)
  53   (declare (special var v s))
  54   (prog (r r1 r2)
  55      (declare (special r1 r2))
  56      (setq r (sin (* 1e-3 (random 1000.))))
  57      (setq v (maxima-substitute r var (ratdisrep (cons s (cons v 1)))))
  58      (setq v (meval '(($ev) v $numer)))
  59      (cond ((and (numberp v) (< (abs v) (* r 1e-2)))
  60             (return t))
  61            ((numberp v) (return nil)))
  62      (if (and (free v '$%i) (not (isinop v '%log)))
  63          (return '$dontknow))
  64      (setq r1 ($realpart v))
  65      (setq r1 (meval '(($ev) r1 $numer)))
  66      (if (not (numberp r1)) (return '$dontknow))
  67      (setq r2 ($imagpart v))
  68      (setq r2 (meval '(($ev) r2 $numer)))
  69      (if (not (numberp r2)) (return '$dontknow))
  70      (cond ((and (< (abs r1) (* r 1e-2))
  71                  (< (abs r2) (* r 1e-2)))
  72             (return t))
  73            (t (return nil)))))