home *** CD-ROM | disk | FTP | other *** search
/ Enter 2004 January / enter-2004-01.iso / files / maxima-5.9.0.exe / {app} / share / maxima / 5.9.0 / src / defint.lisp < prev    next >
Encoding:
Text File  |  2003-02-09  |  80.4 KB  |  2,690 lines
  1. ;;; -*-  Mode: Lisp; Package: Maxima; Syntax: Common-Lisp; Base: 10 -*- ;;;;
  2. ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  3. ;;;     The data in this file contains enhancments.                    ;;;;;
  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. (macsyma-module defint)
  13.  
  14. ;;;          this is the definite integration package. 
  15. ;    defint does definite integration by trying to find an
  16. ;appropriate method for the integral in question.  the first thing that
  17. ;is looked at is the endpoints of the problem. 
  18. ;
  19. ;    i(grand,var,a,b) will be used for integrate(grand,var,a,b)
  20.  
  21. ;;references are to evaluation of definite integrals by symbolic
  22. ;;manipulation by paul s. wang.
  23. ;;
  24. ;;    nointegrate is a macsyma level flag which inhibits indefinite
  25. ;integration.
  26. ;    abconv is a macsyma level flag which inhibits the absolute
  27. ;convergence test.
  28. ;
  29. ;    $defint is the top level function that takes the user input
  30. ;and does minor changes to make the integrand ready for the package.
  31. ;
  32. ;    next comes defint, which is the function that does the
  33. ;integration.  it is often called recursively from the bowels of the
  34. ;package.  defint does some of the easy cases and dispatches to:
  35. ;
  36. ;    dintegrate.  this program first sees if the limits of
  37. ;integration are 0,inf or minf,inf.  if so it sends the problem to
  38. ;ztoinf or mtoinf, respectivly.
  39. ;    else, dintegrate tries:
  40. ;
  41. ;    intsc1 - does integrals of sin's or cos's or exp(%i var)'s
  42. ;         when the interval is 0,2 %pi or 0,%pi.
  43. ;         method is conversion to rational function and find
  44. ;         residues in the unit circle. [wang, pp 107-109]
  45. ;
  46. ;    ratfnt - does rational functions over finite interval by
  47. ;         doing polynomial part directly, and converting
  48. ;         the rational part to an integral on 0,inf and finding
  49. ;         the answer by residues.
  50. ;
  51. ;    zto1   - i(x^(k-1)*(1-x)^(l-1),x,0,1) = beta(k,l)  or
  52. ;         i(log(x)*x^(x-1)*(1-x)^(l-1),x,0,1) = psi...
  53. ;         [wang, pp 116,117]
  54. ;
  55. ;    dintrad- i(x^m/(a*x^2+b*x+c)^(n+3/2),x,0,inf) [wang, p 74]
  56. ;
  57. ;    dintlog- i(log(g(x))*f(x),x,0,inf) = 0 (by symetry) or
  58. ;         tries an integration by parts.  (only routine to
  59. ;         try integration by parts) [wang, pp 93-95]
  60. ;
  61. ;    dintexp- i(f(exp(x)),x,a,b) = i(f(x+1),x,a',b')
  62. ;
  63. ;dintegrate also tries indefinite integration based on certain 
  64. ;predicates (such as abconv) and tries breaking up the integrand
  65. ;over a sum or tries a change of variable.
  66. ;
  67. ;    ztoinf is the routine for doing integrals over the range 0,inf.
  68. ;          it goes over a series of routines and sees if any will work:
  69. ;
  70. ;       scaxn  - sc(b*x^n) (sc stands for sin or cos) [wang, pp 81-83]
  71. ;
  72. ;       ssp    - a*sc^n(r*x)/x^m  [wang, pp 83,84]
  73. ;
  74. ;       zmtorat- rational function. done by multiplication by plog(-x)
  75. ;            and finding the residues over the keyhole contour
  76. ;            [wang, pp 59-61]
  77. ;
  78. ;       log*rat- r(x)*log^n(x) [wang, pp 89-92]
  79. ;
  80. ;       logquad0 log(x)/(a*x^2+b*x+c) uses formula
  81. ;            i(log(x)/(x^2+2*x*a*cos(t)+a^2),x,0,inf) =
  82. ;            t*log(a)/sin(t).  a better formula might be
  83. ;            i(log(x)/(x+b)/(x+c),x,0,inf) = 
  84. ;            (log^2(b)-log^2(c))/(2*(b-c))
  85. ;
  86. ;       batapp - x^(p-1)/(b*x^n+a)^m uses formula related to the beta
  87. ;            function [wang, p 71]
  88. ;            there is also a special case when m=1 and a*b<0
  89. ;            see [wang, p 65]
  90. ;
  91. ;          sinnu  - x^-a*n(x)/d(x) [wang, pp 69-70]
  92. ;
  93. ;       ggr    - x^r*exp(a*x^n+b) 
  94. ;
  95. ;       dintexp- see dintegrate
  96. ;
  97. ;     ztoinf also tries 1/2*mtoinf if the integrand is an even function
  98. ;
  99. ; mtoinf is the routine for doing integrals on minf,inf.  
  100. ;        it too tries a series of routines and sees if any succeed.
  101. ;
  102. ;     scaxn  - when the integrand is an even function, see ztoinf
  103. ;
  104. ;     mtosc  - exp(%i*m*x)*r(x) by residues on either the upper half
  105. ;          plane or the lower half plane, depending on whether
  106. ;          m is positive or negative.
  107. ;
  108. ;     zmtorat- does rational function by finding residues in upper
  109. ;              half plane
  110. ;
  111. ;     dintexp- see dintegrate
  112. ;
  113. ;     rectzto%pi2 - poly(x)*rat(exp(x)) by finding residues in
  114. ;               rectangle [wang, pp98-100]
  115. ;
  116. ;     ggrm   - x^r*exp((x+a)^n+b)
  117. ;
  118. ;   mtoinf also tries 2*ztoinf if the integrand is an even function.
  119.  
  120. (load-macsyma-macros rzmac)
  121.  
  122. (DECLARE-TOP(*lexpr $DIFF $LIMIT $SUBSTITUTE $EZGCD $RATSIMP context)
  123.      (*expr subfunmake $coeff $logcontract $radcan $makegamma
  124.         $constantp $subvarp MAXIMA-SUBSTITUTE freeof ith
  125.         $oddp $hipow $multthru $xthru $num $denom 
  126.         stripdollar MAXIMA-FIND sdiff partition
  127.         constant free mapatom
  128.  
  129.         $ratdisrep ratdisrep $ratp ratp ratnumerator 
  130.         sratsimp ratdenominator $ratsubst ratnump ratcoef
  131.         pterm rdis pdis ratrep newvar pdivide pointergp
  132.               
  133.         $factor factor $sqfr oddelm zerop1
  134.  
  135.         $asksign asksign $sign ask-integer assume forget
  136.               
  137.         $residue residue res res1 polelist partnum
  138.  
  139.         solve solvex sinint
  140.               
  141.         $rectform $realpart $imagpart trisplit cabs
  142.               
  143.         among involve notinvolve  
  144.         numden* ratgreaterp
  145.         subin polyinx genfind xor fmt polyp numden andmapcar
  146.         abless1 even1 rddeg tansc radicalp deg simplerd
  147.         no-err-sub oscip %einvolve sin-sq-cos-sq-sub)
  148.               
  149. ;;;rsn* is in comdenom. does a ratsimp of numerator.
  150.      (SPECIAL *DEF2* PCPRNTD MTOINF* RSN* SEMIRAT*
  151.           SN* SD* LEADCOEF CHECKFACTORS 
  152.           *NODIVERG RD* EXP1
  153.           UL1 LL1 *DFLAG BPTU BPTD PLM* ZN ZD
  154.           *UPDN UL LL EXP PE* PL* RL* PL*1 RL*1
  155.           LOOPSTOP* VAR NN* ND* DN* p*
  156.           IND* FACTORS RLM*
  157.           PLOGABS *ZEXPTSIMP? SCFLAG
  158.           sin-cos-recur rad-poly-recur dintlog-recur
  159.           dintexp-recur defintdebug defint-assumptions
  160.           current-assumptions
  161.           global-defint-assumptions)
  162.      
  163.      (ARRAY* (NOTYPE *I* 1 *J* 1))
  164.      (GENPREFIX DEF)
  165.      (muzzled t)
  166.      ;expvar
  167.      (special $intanalysis $abconvtest $noprincipal $nointegrate)
  168.      ;impvar
  169.      (special $solveradcan $solvetrigwarn *roots *failures 
  170.          $logabs $tlimswitch $maxposex $maxnegex
  171.          $trigsign $savefactors $radexpand $breakup $%emode
  172.          $float $exptsubst dosimp context rp-polylogp
  173.          %P%I HALF%PI %PI2 HALF%PI3 VARLIST genvar
  174.          $domain $m1pbranch errorsw errrjfflag raterr
  175.          limitp $algebraic
  176.          ;;LIMITP T Causes $ASKSIGN to do special things
  177.          ;;For DEFINT like eliminate epsilon look for prin-inf
  178.          ;;take realpart and imagpart.
  179.          integer-info
  180.          ;;If LIMITP is non-null ask-integer conses 
  181.          ;;its assumptions onto this list.
  182.          generate-atan2))
  183.                    ;If this switch is () then RPART returns ATAN's
  184.                            ;instead of ATAN2's
  185.  
  186. (declare-top(special infinities real-infinities infinitesimals))
  187. (cond ;These are really defined in LIMIT but DEFINT uses them also.
  188.   ((not (boundp 'infinities))
  189.   (setq INFINITIES '($INF $MINF $INFINITY))
  190.   (setq REAL-INFINITIES '($INF $MINF))
  191.   (setq INFINITESIMALS '($ZEROA $ZEROB))))
  192.  
  193. (defmvar defintdebug () "If true Defint prints out debugging information")
  194.  
  195. (DEFMVAR INTEGERL NIL
  196.   "An integer-list for non-atoms found out to be INTEGERs")
  197.  
  198. (DEFMVAR NONINTEGERL NIL
  199.   "A non-integer-list for non-atoms found out to be NONINTEGERs")
  200.  
  201. (DEFUN $DEFINT (EXP VAR LL UL)
  202.  (let ((global-defint-assumptions ())
  203.        (integer-info ()) (integerl integerl) (nonintegerl nonintegerl))
  204.  (with-new-context (context)
  205.  (unwind-protect
  206.        (let ((defint-assumptions ())  (*def2* ())  (rad-poly-recur ())
  207.          (sin-cos-recur ())  (dintexp-recur ())  (dintlog-recur 0.)
  208.              (ans nil)  (orig-exp exp)  (orig-var var)
  209.          (orig-ll ll)  (orig-ul ul) 
  210.          (pcprntd nil)  (*NODIVerg nil)  ($logabs t)  (limitp t)
  211.          (rp-polylogp ())
  212.          ($domain '$real) ($m1pbranch ())) ;Try this out.
  213.  
  214.         (FIND-FUNCTION '$LIMIT)
  215.         (make-global-assumptions) ;sets global-defint-assumptions
  216.         (FIND-FUNCTION '$RESIDUE)
  217.         (SETQ EXP (RATDISREP EXP))
  218.         (SETQ VAR (RATDISREP VAR))
  219.         (SETQ LL (RATDISREP LL))
  220.         (SETQ UL (RATDISREP UL))
  221.         (COND (($CONSTANTP VAR)
  222.            (merror "Variable of integration not a variable: ~M"
  223.                VAR))
  224.           (($SUBVARP VAR)  (SETQ VAR (STRIPDOLLAR (CAAR VAR)))
  225.                    (SETQ EXP ($SUBSTITUTE VAR orig-VAR EXP))))
  226.         (COND ((NOT (ATOM VAR))
  227.            (merror "Improper variable of integration: ~M" VAR))
  228.           ((OR (AMONG VAR UL)
  229.                (AMONG VAR LL)) 
  230.            (SETQ VAR (STRIPDOLLAR VAR))
  231.            (SETQ EXP ($SUBSTITUTE VAR orig-VAR EXP))))
  232.         (cond ((not (equal ($ratsimp ($imagpart ll)) 0))
  233.            (merror "Defint: Lower limit of integration must be real."))
  234.           ((not (equal ($ratsimp ($imagpart ul)) 0))
  235.            (merror
  236.             "Defint: Upper limit of integration must be real.")))
  237.  
  238.         (COND ((SETQ ANS (DEFINT EXP VAR LL UL))
  239.            (SETQ ANS (SUBST orig-VAR VAR ANS))
  240.            (COND ((atom ans)  ans)
  241.              ((and (free ans '%limit)
  242.                    (free ans '%integrate)
  243.                    (OR (not (free ans '$INF))
  244.                    (not (free ans '$MINF))
  245.                    (not (free ans '$INFINITY))))
  246.               (diverg))
  247.              ((not (free ans '$und))
  248.              `((%integrate) ,orig-exp ,orig-var ,orig-ll ,orig-ul))
  249.                   (t ANS)))
  250.           (t `((%integrate) ,orig-exp ,orig-var ,orig-ll ,orig-ul))))
  251.        (forget-global-assumptions)))))
  252.  
  253. (DEFUN EEZZ (EXP LL UL)
  254.        (COND ((OR (polyinx EXP VAR nil)
  255.           (CATCH 'PIN%EX (PIN%EX EXP)))
  256.           (SETQ EXP (ANTIDERIV EXP))
  257. ;;;If antideriv can't do it, returns nil
  258. ;;;use limit to evaluate every answer returned by antideriv.
  259.           (COND ((NULL EXP) NIL)
  260.             (t (INTSUBS EXP LL UL))))))
  261. ;;;Hack the expression up for exponentials.
  262.  
  263. (defun sinintp (expr var)
  264. ;;; Is this expr a candidate for SININT ?
  265.        (let ((expr (factor expr))
  266.          (numer nil)
  267.          (denom nil))
  268.         (setq numer ($num expr))
  269.         (setq denom ($denom expr))
  270.         (cond ((polyinx numer var nil)
  271.            (cond ((and (polyinx denom var nil)
  272.                    (deg-lessp denom var 2))
  273.               t)))
  274. ;;;ERF type things go here.
  275.           ((let ((exponent (%einvolve numer)))
  276.             (and (polyinx exponent var nil)
  277.                  (deg-lessp exponent var 2)))
  278.            (cond ((free denom var)
  279.               t))))))
  280.  
  281. (defun deg-lessp (expr var power)
  282.        (cond  ((or (atom expr) 
  283.            (mnump expr)) t)
  284.           ((or (mtimesp expr) 
  285.            (mplusp expr))
  286.            (do ((ops (cdr expr) (cdr ops)))
  287.            ((null ops) t)
  288.            (cond ((not (deg-lessp (car ops) var power))
  289.               (return ())))))
  290.           ((mexptp expr)
  291.            (and (or (not (alike1 (cadr expr) var))
  292.             (and (numberp (caddr expr))
  293.                  (not (eq (asksign (m+ power (m- (caddr expr))))
  294.                       '$negative))))
  295.             (deg-lessp (cadr expr) var power)))))
  296.  
  297. (DEFUN ANTIDERIV (A)
  298.        (let ((limitp ())
  299.          (Ans ())
  300.          (generate-atan2 ()))
  301.      (setq ans (SININT A VAR))
  302.      (COND ((AMONG '%INTEGRATE Ans)  NIL) 
  303.            (T (SIMPLIFY Ans)))))
  304.  
  305. ;;;This routine tries to take a limit a couple of ways.
  306. (defmfun get-limit nargs
  307.   (let ((ans (apply 'limit-no-err (listify nargs)))
  308.     (val ()) (var ()) (exp ()) (dir ()))
  309.     (cond ((and ans (not (among '%limit ans)))  ans)
  310.       (t (cond ((and (or (equal nargs 3) (equal nargs 4))
  311.              (memq (setq val (arg 3)) '($inf $minf)))
  312.             (setq var (arg 2))
  313.             (setq exp (MAXIMA-SUBSTITUTE (m^t var -1) var (arg 1)))
  314.             (cond ((eq val '$inf)  (setq dir '$plus))
  315.               (t (setq dir '$minus)))
  316.             ;(setq ans (apply 'limit-no-err `(,exp ,var 0 ,dir)))
  317.             (setq ans (limit-no-err exp var 0 dir))
  318.             (cond ((not (among '%limit ans))  ans)
  319.               (t ()))))))))
  320.  
  321. ;(defun limit-no-err nargs
  322. ;  (let ((errorsw t)  (ans ()))
  323. ;    (setq ans (catch 'errorsw (apply '$limit (listify nargs))))
  324. ;    (cond ((not (eq ans t))  ans)
  325. ;      (t nil))))
  326. (defun limit-no-err (&rest argvec)
  327.   (declare (special errorsw))
  328.   (let ((errorsw t) (ans nil))
  329.     (setq ans (catch 'errorsw (apply '$limit argvec)))
  330.     (if (eq ans t) nil ans)))
  331.  
  332. (DEFUN INTCV (NV IND FLAG)
  333.   (let ((D (bx**n+a nv))
  334.     (*ROOTS ())  (*FAILURES ())  ($BREAKUP ()))
  335.     (cond ((AND (EQ UL '$INF)
  336.         (EQUAL LL 0)
  337.         (EQUAL (CADR D) 1)) ())
  338.       (t (SOLVE (m+t 'YX (m*t -1. NV)) VAR 1.)
  339.          (COND (*ROOTS (SETQ D (SUBST VAR 'YX (CADDAR *ROOTS)))
  340.                (COND (FLAG (INTCV2 D IND NV))
  341.                  (T (INTCV1 D IND NV))))
  342.            (t ()))))))
  343.  
  344. (DEFUN INTCV1 (D IND NV) 
  345.        (COND ((AND (INTCV2 D IND NV)
  346.            (NOT (alike1 LL1 UL1)))
  347.           (let ((*DEF2* t))
  348.            (DEFINT EXP1 VAR LL1 UL1)))))
  349.  
  350. (DEFUN INTCV2 (D IND NV)
  351.        (INTCV3 D IND NV)
  352.        (AND (COND ((AND (ZEROP1 (m+ LL UL))
  353.             (EVENFN NV VAR))
  354.            (SETQ EXP1 (m* 2 EXP1)
  355.              LL1 (LIMCP NV VAR 0 '$PLUS)))
  356.           (T (SETQ LL1 (LIMCP NV VAR LL '$PLUS))))
  357.         (SETQ UL1 (LIMCP NV VAR UL '$MINUS))))
  358.  
  359. (DEFUN LIMCP (A B C D) 
  360.        (let ((Ans ($LIMIT A B C D)))
  361.         (COND ((NOT (OR (null ans)
  362.                 (among '%limit ans)
  363.                 (AMONG '$IND Ans)
  364.                 (AMONG '$UND Ans)))
  365.            Ans))))
  366.  
  367. (DEFUN INTCV3 (D IND NV)
  368.        (SETQ NN* ($RATSIMP (SDIFF D VAR)))
  369.        (SETQ EXP1 (SUBST 'YX NV EXP))
  370.        (SETQ EXP1 (m* NN* (COND (IND EXP)
  371.                 (T (SUBST D VAR EXP1)))))
  372.        (SETQ EXP1 (sRATSIMP (SUBST VAR 'YX EXP1))))
  373.  
  374. (DEFUN DEFINT (EXP VAR LL UL)
  375.  (let ((old-assumptions defint-assumptions)  (current-assumptions ()))
  376.   (unwind-protect
  377.    (prog ()
  378.       (setq current-assumptions (make-defint-assumptions 'noask))
  379.       (let ((exp (resimplify exp))            
  380.         (var (resimplify var))
  381.         ($exptsubst t)
  382.         (loopstop* 0)
  383.         ;; D (not used? -- cwh)
  384.         ANS NN* DN* ND* $NOPRINCIPAL)
  385.        (cond ((setq ans (defint-list exp var ll ul))
  386.           (return ans))
  387.          ((OR (ZEROP1 EXP)
  388.               (ALIKE1 UL LL))
  389.           (RETURN 0.))
  390.          ((NOT (AMONG VAR EXP))
  391.           (COND ((OR (MEMQ UL '($INF $MINF))
  392.                  (MEMQ LL '($INF $MINF)))
  393.              (Diverg))
  394.             (T (SETQ ANS (m* EXP (m+ UL (m- LL))))
  395.                (return ANS)))))
  396.        (let* ((EXP (RMCONST1 EXP))
  397.           (C (CAR EXP))
  398.           (EXP (%i-out-of-denom (CDR EXP))))
  399.          (COND ((AND (NOT $NOINTEGRATE)
  400.              (NOT (ATOM EXP))
  401.              (or (among 'mqapply exp)
  402.                  (NOT (MEMQ (CAAR EXP)
  403.                     '(MEXPT MPLUS MTIMES %SIN %COS
  404.                         %TAN %SINH %COSH %TANH
  405.                         %LOG %ASIN %ACOS %ATAN
  406.                         %COT %ACOT %SEC 
  407.                         %ASEC %CSC %ACSC 
  408.                         %DERIVATIVE)))))
  409.             (COND ((SETQ ANS (ANTIDERIV EXP))
  410.                (SETQ ANS (INTSUBS ANS LL UL))
  411.                (return (m* C ANS)))
  412.               (t (return nil)))))
  413.          (SETQ EXP (TANSC EXP))
  414.          (cond ((setq  ans (initial-analysis exp var ll ul))
  415.             (return (m* c ans))))
  416.          (return nil))))
  417.    (restore-defint-assumptions old-assumptions current-assumptions))))
  418.  
  419. (defun defint-list (exp var ll ul)
  420.        (COND ((AND (NOT (ATOM EXP)) 
  421.            (MEMQ (CAAR EXP)
  422.              '(MEQUAL MLIST $MATRIX)))
  423.           (let ((ANS (CONS (CAR EXP)
  424.                    (MAPCAR
  425.                 #'(LAMBDA (SUB-EXP)
  426.                     (DEFINT SUB-EXP VAR LL UL))
  427.                 (CDR EXP)))))
  428.            (COND (ANS (simplify ANS))
  429.              (T NIL))))
  430.          (t nil)))
  431.  
  432. (defun initial-analysis (exp var ll ul)
  433.        (let ((pole (cond ((not $intanalysis)
  434.               '$no)          ;don't do any checking.
  435.              (t (poles-in-interval exp var ll ul)))))
  436.         (cond ((eq pole '$no)
  437.            (cond ((and (oddfn exp var)
  438.                    (or (and (eq ll '$minf)
  439.                     (eq ul '$inf))
  440.                    (eq ($sign (m+ ll ul))
  441.                        '$zero)))  0)
  442.              (t (parse-integrand exp var ll ul))))
  443.           ((eq pole '$unknown)  ())
  444.           (t (principal-value-integral exp var ll ul pole)))))
  445.  
  446. (defun parse-integrand (exp var ll ul)
  447.   (let (ans)
  448.        (COND ((SETQ ANS (EEZZ EXP LL UL))  ans)
  449.          ((and (RATP EXP VAR)
  450.            (setq ans (method-by-limits EXP VAR LL UL)))  ans)
  451.          ((and (mplusp EXP)
  452.            (setq ans (INTBYTERM EXP T)))  ans)
  453.          ((setq ans (method-by-limits EXP VAR LL UL))  ans)
  454.          (t ()))))
  455.  
  456. (DEFUN RMCONST1 (E)
  457.   (COND ((AMONG VAR E) 
  458.      (PARTITION E VAR 1))
  459.     (T (CONS E 1))))
  460.  
  461.  
  462. (defun method-by-limits (exp var ll ul)
  463.  (let ((old-assumptions defint-assumptions))
  464.    (setq current-assumptions (make-defint-assumptions 'noask))
  465. ;;Should be a PROG inside of unwind-protect, but Multics has a compiler
  466. ;;bug wrt. and I want to test this code now.
  467.    (unwind-protect
  468.        (COND ((and (and (EQ UL '$INF)
  469.             (eq ll '$minf))  (mtoinf exp var)))
  470.          ((and (and (eq ul '$inf)
  471.             (equal ll 0.))  (ztoinf exp var)))
  472. ;;;This seems((and (and (eq ul '$inf)
  473. ;;;fairly losing    (setq exp (subin (m+ ll var) exp))
  474. ;;;            (setq ll 0.))
  475. ;;;           (ztoinf exp var)))
  476.          ((and (equal ll 0.)
  477.             (freeof var ul) (eq ($asksign ul) '$pos) (zto1 exp)))
  478. ;         ((and (and (equal ul 1.)
  479. ;            (equal ll 0.))  (zto1 exp)))
  480.          (t (dintegrate exp var ll ul)))
  481.        (restore-defint-assumptions old-assumptions defint-assumptions))))
  482.        
  483.  
  484. (DEFUN DINTEGRATE (EXP VAR LL UL)
  485.    (let ((ans nil) (arg nil) (scflag nil) 
  486.      (*dflag nil) ($%emode t))
  487. ;;;NOT COMPLETE for sin's and cos's.
  488.     (cond ((and (not sin-cos-recur)
  489.             (oscip exp)
  490.             (SETQ SCFLAG T)
  491.             (INTSC1 ll ul exp)))
  492.           ((and (not rad-poly-recur)
  493.             (notinvolve exp '(%log))
  494.             (not (%einvolve exp))
  495.             (method-radical-poly exp var ll ul)))
  496.           ((and (not (equal dintlog-recur 2.))
  497.             (SETQ arg (INVOLVE EXP '(%log)))
  498.             (DINTLOG EXP arg)))
  499.           ((and (not dintexp-recur)
  500.             (SETQ arg (%EINVOLVE exp))
  501.             (DINTEXP EXP var)))
  502.           ((and (NOT (RATP EXP VAR)) 
  503.             (SETQ ANS ($EXPAND EXP))
  504.             (NOT (ALIKE1 ANS EXP))
  505.             (INTBYTERM ANS T)))
  506.           ((setq ans (antideriv exp))
  507.            (intsubs ans ll ul))
  508.           (t nil))))
  509.  
  510. (defun method-radical-poly (exp var ll ul)
  511. ;;;Recursion stopper
  512.   (let ((rad-poly-recur t)   ;recursion stopper
  513.     (result ()))
  514.     (cond ((and (sinintp EXP VAR) 
  515.         (setq result (antideriv exp))
  516.         (intsubs result ll ul)))
  517.       ((and (RATP EXP VAR)
  518.         (setq result (RATFNT EXP))))
  519.       ((and (setq result (antideriv exp))
  520.         (intsubs result ll ul)))
  521.       ((AND (NOT SCFLAG)
  522.         (NOT (EQ UL '$INF))
  523.         (radic exp var)
  524.         (KINDP34)
  525.         (setq result (CV EXP))))
  526.       (t ()))))
  527.  
  528. ;;; LIMIT loss can't set logabs to true for these cases.
  529. (defun principal-value-integral (exp var ll ul poles)
  530.   (let (($logabs ())  (anti-deriv ()))
  531.     (cond ((not (null (setq anti-deriv (antideriv exp))))
  532.        (cond ((not (null poles))
  533.           (order-limits 'ask)
  534.           (cond ((take-principal anti-deriv ll ul poles))
  535.             (t ()))))))))
  536.  
  537. (defun take-principal (anti-deriv ll ul poles &aux ans merged-list)
  538.   (setq anti-deriv (cond ((involve anti-deriv '(%log))
  539.               ($logcontract anti-deriv))
  540.              (t anti-deriv)))
  541.   (setq ans 0.)
  542.   (setq merged-list (interval-list poles ll ul))
  543.   (do ((current-pole (cdr merged-list) (cdr current-pole))
  544.        (previous-pole merged-list (cdr previous-pole)))
  545.       ((null current-pole)  t)
  546.     (setq ans (m+ ans        
  547.           (intsubs anti-deriv (m+ (caar previous-pole) 'epsilon)
  548.                (m+ (caar current-pole) (m- 'epsilon))))))
  549.                
  550. ;;;Hack answer to simplify "Correctly".
  551.     (cond ((not (freeof '%log ans)) 
  552.        (setq ans ($logcontract ans))))
  553.     (setq ans (get-limit (get-limit ans 'epsilon 0 '$plus) 'prin-inf '$inf))
  554. ;;;Return setion.
  555.     (cond ((or (null ans)
  556.            (not (free ans '$infinity)) 
  557.            (not (free ans '$ind)))  ())
  558.       ((or (among '$minf ans)
  559.            (among '$inf ans)
  560.            (among '$und ans))
  561.        (diverg))
  562.       (t (principal) ans)))
  563.  
  564. (defun interval-list (pole-list ll ul)
  565.   (let ((first (car (first pole-list)))
  566.     (last (caar (last pole-list))))
  567.     (cond ((eq ul last)  
  568.        (if (eq ul '$inf)
  569.            (setq pole-list (subst 'prin-inf '$inf pole-list))))
  570.       (t (if (eq ul '$inf) 
  571.          (setq ul 'prin-inf))
  572.          (setq pole-list (append pole-list (list (cons ul 'ignored))))))
  573.     (cond ((eq ll first) 
  574.        (if (eq ll '$minf)
  575.            (setq pole-list (subst (m- 'prin-inf) '$minf pole-list))))
  576.       (t (if (eq ll '$minf)
  577.          (setq ll (m- 'prin-inf)))
  578.          (setq pole-list (append (list (cons ll 'ignored)) pole-list)))))
  579.   pole-list)
  580.  
  581. (DEFUN CV (EXP)
  582.   (if (not (or (real-infinityp ll) (real-infinityp ul)))
  583.       (method-by-limits (INTCV3 (M// (m+t LL (m*t UL VAR))
  584.                      (m+t 1. VAR)) NIL 'YX)
  585.             VAR 0. '$INF)
  586.       ()))
  587.  
  588. (DEFUN RATFNT (EXP)
  589.    (let ((e (pqr exp)))
  590.      (COND ((EQUAL 0. (CAR E))  (CV EXP))
  591.        ((EQUAL 0. (CDR E))  (EEZZ (CAR E) LL UL))
  592.        (T (m+t (EEZZ (CAR E) LL UL)
  593.            (CV (M// (CDR E) DN*)))))))
  594.  
  595. (DEFUN PQR (E)
  596.   (let ((VARLIST (list var)))
  597.        (NEWVAR E)
  598.        (SETQ E (CDR (RATREP* E)))
  599.        (SETQ DN* (PDIS (RATDENOMINATOR E)))
  600.        (SETQ E (PDIVIDE (RATNUMERATOR E) (RATDENOMINATOR E)))
  601.        (CONS (SIMPLIFY (RDIS (CAR E))) (SIMPLIFY (RDIS (CADR E))))))
  602.  
  603.  
  604. (DEFUN INTBYTERM (EXP *NODIVERG)
  605.   (let ((saved-exp exp))
  606.     (COND ((mplusp EXP)
  607.        (let ((ans (CATCH 'Divergent 
  608.                   (ANDMAPCAR #'(LAMBDA (new-exp) 
  609.                          (let ((*def2* t))
  610.                            (DEFINT new-exp VAR LL UL)))
  611.                      (CDR EXP)))))
  612.         (COND ((NULL ans) NIL)
  613.               ((EQ ans 'Divergent)
  614.                (let ((*NODIVerg nil))
  615.              (cond ((setq ans (ANTIDERIV saved-EXP))
  616.                 (intsubs ans ll ul))
  617.                    (t nil))))
  618.               (T (sRATSIMP (m+l ans))))))
  619. ;;;If leadop isn't plus don't do anything.
  620.       (t nil))))
  621.  
  622. (DEFUN KINDP34 NIL
  623.   (NUMDEN EXP)
  624.   (let* ((d dn*)
  625.      (a (COND ((and (ZEROP1 ($LIMIT D VAR LL '$PLUS))
  626.             (eq (limit-pole (m+ exp (m+ (m- LL) VAR)) var LL '$PLUS) '$yes))
  627.            t)
  628.           (t nil)))
  629.      (b (COND ((and (ZEROP1 ($LIMIT D VAR UL '$MINUS))
  630.             (eq (limit-pole (m+ exp (m+ UL (m- VAR))) var UL '$MINUS) '$yes))
  631.            t)
  632.           (t nil))))
  633.     (or a b)))
  634.  
  635. (DEFUN Diverg NIL
  636.   (COND (*NODIVERG (THROW 'Divergent 'Divergent))
  637.     (T (MERROR "Integral is divergent"))))
  638.  
  639. (defun make-defint-assumptions (ask-or-not)
  640.   (cond ((null (order-limits ask-or-not))  ())
  641.     (t (mapc 'forget defint-assumptions)
  642.        (setq defint-assumptions ())
  643.        (let ((sign-ll (cond ((eq ll '$inf)  '$pos)
  644.                 ((eq ll '$minf) '$neg)
  645.                 (t ($sign ($limit ll)))))
  646.          (sign-ul (cond ((eq ul '$inf)  '$pos)
  647.                 ((eq ul '$minf)  '$neg)
  648.                 (t ($sign ($limit ul)))))
  649.          (sign-ul-ll (cond ((and (eq ul '$inf)
  650.                      (not (eq ll '$inf)))  '$pos)
  651.                    ((and (eq ul '$minf)
  652.                      (not (eq ll '$minf)))  '$neg)
  653.                    (t ($sign ($limit (m+ ul (m- ll))))))))
  654.         (cond ((eq sign-ul-ll '$pos)
  655.                (setq defint-assumptions
  656.                  `(,(assume `((mgreaterp) ,var ,ll))
  657.                    ,(assume `((mgreaterp) ,ul ,var)))))
  658.               ((eq sign-ul-ll '$neg)
  659.                (setq defint-assumptions
  660.                  `(,(assume `((mgreaterp) ,var ,ul))
  661.                    ,(assume `((mgreaterp) ,ll ,var))))))
  662.         (cond ((and (eq sign-ll '$pos)
  663.                 (eq sign-ul '$pos))
  664.                (setq defint-assumptions
  665.                  `(,(assume `((mgreaterp) ,var 0))
  666.                    ,@defint-assumptions)))
  667.               ((and (eq sign-ll '$neg)
  668.                 (eq sign-ul '$neg))
  669.                (setq defint-assumptions
  670.                  `(,(assume `((mgreaterp) 0 ,var))
  671.                    ,@defint-assumptions)))
  672.               (t defint-assumptions))))))
  673.  
  674. (defun restore-defint-assumptions (old-assumptions assumptions)
  675.   (do ((llist assumptions (cdr llist)))
  676.       ((null llist) t)
  677.       (forget (car llist)))
  678.   (do ((llist old-assumptions (cdr llist)))
  679.       ((null llist) t)
  680.       (assume (car llist))))
  681.  
  682. (defun make-global-assumptions ()
  683.         (setq global-defint-assumptions
  684.           (cons (assume '((mgreaterp) *z* 0.))
  685.             global-defint-assumptions))
  686. ;;; *Z* is a "zero parameter" for this package.
  687. ;;; Its also used to transform.
  688. ;;  limit(exp,var,val,dir) -- limit(exp,tvar,0,dir)
  689.         (setq global-defint-assumptions
  690.           (cons (assume '((mgreaterp) epsilon 0.))
  691.             global-defint-assumptions))       
  692.         (setq global-defint-assumptions
  693.           (cons (assume '((mlessp) epsilon 1.0e-8))
  694.             global-defint-assumptions)) 
  695. ;;; EPSILON is used in principal vaule code to denote the familiar
  696. ;;; mathematical entity.
  697.         (setq global-defint-assumptions
  698.           (cons (assume '((mgreaterp) prin-inf 1.0e+8))
  699.             global-defint-assumptions)))
  700. ;;; PRIN-INF Is a special symbol in the principal value code used to
  701. ;;; denote an end-point which is proceeding to infinity.
  702.  
  703. (defun forget-global-assumptions ()
  704.        (do ((llist global-defint-assumptions (cdr llist)))
  705.        ((null llist) t)
  706.        (forget (car llist)))
  707.        (cond ((not (null integer-info))
  708.           (do ((llist integer-info (cdr llist)))
  709.           ((null llist) t)
  710.           (I-$remove `(,(cadar llist) ,(caddar llist)))))))
  711.  
  712. (DEFUN order-limits (ask-or-not)
  713.        (cond ((or (not (equal ($imagpart ll) 0))
  714.           (not (equal ($imagpart ul) 0)))  ())
  715.          (t (COND ((ALIKE1 LL (m*t -1 '$INF))
  716.                (SETQ LL '$MINF)))
  717.         (COND ((ALIKE1 UL (m*t -1 '$INF))
  718.                (SETQ UL '$MINF)))
  719.         (cond ((alike1 ll (m*t -1 '$minf))
  720.                (setq ll '$inf)))
  721.         (cond ((alike1 ul (m*t -1 '$minf))
  722.                (setq ul '$inf)))
  723.         (COND ((EQ UL '$INF) NIL)
  724.               ((EQ LL '$MINF)
  725.                (SETQ EXP (SUBIN (m- VAR) EXP))
  726.                (SETQ LL (m- ul))
  727.                (SETQ UL '$INF))
  728.               ((eq ll '$inf)
  729.                (setq ll ul)
  730.                (setq ul '$inf)
  731.                (setq exp (m- exp))))
  732. ;;;Fix limits so that ll < ul. 
  733.         (let ((D (COMPLM ask-or-not)))
  734.              (COND ((EQUAL D -1.)
  735.                 (SETQ exp (m- exp))
  736.                 (SETQ D LL)
  737.                 (SETQ LL UL)
  738.                 (SETQ UL D))
  739.                (T T))))))
  740.  
  741. (DEFUN COMPLM (ask-or-not)
  742.   (let ((askflag (cond ((eq ask-or-not 'ask)  t)
  743.                (t nil)))
  744.     (A ()))
  745.     (COND ((ALIKE1 UL LL)  0.)
  746.       ((EQ (SETQ A (cond (askflag ($asksign ($limit (m+t UL (m- LL)))))
  747.                  (t ($sign ($limit (m+t UL (m- LL)))))))
  748.            '$pos)  1.)
  749.       ((EQ A '$neg)  -1.)
  750.       (T 1.))))
  751.  
  752.  
  753. (DEFUN INTSUBS (E A B)
  754.   (cond ((easy-subs e a b))
  755.     (t (setq current-assumptions
  756.          (make-defint-assumptions 'ask)) ;get forceful!
  757.        (let ((generate-atan2 ())  ($algebraic t)
  758.          (rpart ())  (ipart ()))
  759.          (desetq (rpart . ipart) (cond ((not (free e '$%i))
  760.                         (trisplit e))
  761.                        (t (cons e 0))))
  762.          (cond ((not (equal (sratsimp ipart) 0))  
  763.             (let ((rans (cond ((limit-subs rpart a b))
  764.                       (t (m-
  765.                       `((%limit) ,rpart ,var ,b $minus)
  766.                       `((%limit) ,rpart ,var ,a $plus)))))
  767.               (ians (cond ((limit-subs ipart a b))
  768.                       (t (m-
  769.                       `((%limit) ,ipart ,var ,b $minus)
  770.                       `((%limit) ,ipart ,var ,a $plus))))))
  771.               (m+ rans (m* '$%i ians))))
  772.            (t (setq rpart (sratsimp rpart))
  773.               (cond ((limit-subs rpart a b))
  774.                 (t (same-sheet-subs rpart a b)))))))))
  775.  
  776. (defun easy-subs (e ll ul)
  777.     (cond ((or (infinityp ll) (infinityp ul)) ())
  778.       (t (cond ((polyinx e var ())
  779.             (let ((ll-val (no-err-sub ll e))
  780.               (ul-val (no-err-sub ul e)))
  781.               (cond ((and ll-val ul-val)  (m- ul-val ll-val))
  782.                 (t ()))))
  783.            (t ())))))
  784.  
  785. (defun limit-subs (e ll ul)
  786.   (cond ((not (free e '%atan))  ())
  787.     (t (setq e ($multthru e))
  788.        (let ((a1 ($limit e var ll '$plus))    
  789.          (a2 ($limit e var ul '$minus)))
  790.          (cond ((MEMQ A1 '($INF $MINF $INFINITY ))
  791.             (COND ((MEMQ A2 '($INF $MINF $INFINITY))
  792.                (COND ((EQ A2 A1)  (Diverg))
  793.                  (T ())))
  794.               (T (Diverg))))
  795.            ((MEMQ A2 '($INF $MINF $INFINITY))  (Diverg))
  796.            ((OR (MEMQ A1 '($UND $IND))
  797.             (MEMQ A2 '($UND $IND)))  ())
  798.            (T (m- A2 A1)))))))
  799.  
  800. ;;;This function works only on things with ATAN's in them now.
  801. (defun same-sheet-subs (exp ll ul &aux ans)
  802.   (let ((poles (atan-poles exp ll ul)))
  803. ;POLES -> ((mlist) ((mequal) ((%atan) foo) replacement) ......) 
  804. ;We can then use $SUBSTITUTE
  805.     (setq ans ($limit exp var ll '$plus))
  806.     (setq exp (sratsimp ($substitute poles exp)))
  807.     (m- ($limit exp var ul '$minus) ans)))
  808.  
  809. (defun atan-poles (exp ll ul)
  810.   `((mlist) ,@(atan-pole1 exp ll ul)))
  811.  
  812. (defun atan-pole1 (exp ll ul &aux ipart)
  813.   (cond 
  814.    ((mapatom exp)  ())
  815.    ((matanp exp) ;neglect multiplicity and '$unknowns for now.
  816.     (desetq (exp . ipart) (trisplit exp))
  817.     (cond 
  818.      ((not (equal (sratsimp ipart) 0))  ())
  819.      (t (let ((pole (poles-in-interval (let (($algebraic t))
  820.                      (sratsimp (cadr exp)))
  821.                        var ll ul)))
  822.       (cond ((and pole (not (or (eq pole '$unknown)
  823.                     (eq pole '$no))))
  824.          (do ((l pole (cdr l)) (llist ()))
  825.              ((null l)  llist)
  826.            (cond 
  827.             ((eq (caar l) ll)  t) ;Skip this one by definition.
  828.             (t (let ((low-lim ($limit (cadr exp) var (caar l) '$minus))
  829.                  (up-lim ($limit (cadr exp) var (caar l) '$plus)))
  830.              (cond ((and (not (eq low-lim up-lim))
  831.                      (real-infinityp low-lim)
  832.                      (real-infinityp up-lim))
  833.                 (let ((change (if (eq low-lim '$minf)
  834.                           (m- '$%pi)
  835.                            '$%pi)))
  836.                   (setq llist (cons `((mequal simp) ,exp  ,(m+ exp change))
  837.                            llist)))))))))))))))
  838.    (t (do ((l (cdr exp) (cdr l))
  839.        (llist ()))
  840.       ((null l)  llist)
  841.     (setq llist (append llist (atan-pole1 (car l) ll ul)))))))
  842.  
  843. (DEFUN DIFAPPLY (N D S FN1)
  844.   (PROG (K M R $NOPRINCIPAL) 
  845.     (COND ((eq ($asksign (m+ (DEG D) (m- S) (m- 2.)))  '$neg)
  846.            (RETURN NIL)))
  847.     (SETQ $NOPRINCIPAL T)
  848.     (COND ((OR (NOT (mexptp D))
  849.            (NOT (NUMBERP (SETQ R (CADDR D)))))
  850.            (RETURN NIL))
  851.           ((AND (EQUAL N 1.)
  852.             (EQ FN1 'MTORAT)
  853.             (EQUAL 1. (DEG (CADR D))))
  854.            (RETURN 0.)))
  855.     (SETQ M (DEG (SETQ D (CADR D))))
  856.     (SETQ K (m// (m+ S 2.) M))
  857.     (COND ((eq (ask-integer (m// (m+ S 2.) M) '$any)  '$yes)
  858.            NIL)
  859.           (T (SETQ K (M+ 1 K))))
  860.     (COND ((eq ($sign (m+ r (m- k))) '$pos)
  861.            (RETURN (DIFFHK FN1 N D K (m+ R (m- K)))))))) 
  862.  
  863. (DEFUN DIFFHK (FN1 N D R M)
  864.   (PROG (D1 *DFLAG) 
  865.     (SETQ *DFLAG T)
  866.     (SETQ D1 (FUNCALL FN1 N
  867.               (m^ (m+t '*Z* D) R)
  868.               (m* R (DEG D))))
  869.     (COND (D1 (RETURN (DIFAP1 D1 R '*Z* M 0.)))))) 
  870.  
  871. (DEFUN PRINCIPAL NIL
  872.        (COND ($NOPRINCIPAL (Diverg))
  873.          ((NOT PCPRNTD)
  874.           (PRINC "Principal Value")
  875.           (SETQ PCPRNTD T))))
  876.  
  877. (DEFUN RIB (E S)
  878.   (let (*UPDN C) 
  879.     (COND ((OR (MNUMP E) (CONSTANT E))
  880.        (SETQ BPTU (CONS E BPTU)))
  881.       (t (SETQ E (RMCONST1 E))
  882.          (SETQ C (CAR E))
  883.          (SETQ NN* (CDR E))
  884.          (SETQ ND* S)
  885.          (SETQ E (CATCH 'PTIMES%E (PTIMES%E NN* ND*)))
  886.          (COND ((NULL E) NIL)
  887.            (T (SETQ E (m* C E))
  888.               (COND (*UPDN (SETQ BPTU (CONS E BPTU)))
  889.                 (T (SETQ BPTD (CONS E BPTD))))))))))
  890.  
  891. (DEFUN PTIMES%E (TERM N)
  892.        (COND ((POLYINX TERM VAR NIL) TERM)
  893.          ((AND (mexptp TERM)
  894.            (EQ (CADR TERM) '$%E)
  895.            (POLYINX (CADDR TERM) VAR NIL)
  896.            (eq ($sign (m+ (DEG ($realpart (CADDR TERM))) -1.))
  897.                '$neg)
  898.            (eq ($sign (m+ (DEG (SETQ NN* ($imagpart (CADDR TERM)))) 
  899.                      -2.))
  900.                '$neg))
  901.           (COND ((EQ ($askSIGN (RATCOEF NN* VAR)) '$pos) 
  902.              (SETQ *UPDN T))
  903.             (T (SETQ *UPDN NIL)))
  904.           TERM)
  905.          ((AND (mtimesp TERM)
  906.            (SETQ NN* (POLFACTORS TERM))
  907.            (OR (NULL (CAR NN*))
  908.                (eq ($sign (m+ n (m- (deg (car nn*)))))
  909.                '$pos))
  910.            (PTIMES%E (CADR NN*) N)
  911.            TERM))
  912.          (T (THROW 'PTIMES%E NIL))))
  913.  
  914. (DEFUN CSEMIDOWN (N D VAR)
  915.   (let ((pcprntd t)) ;Not sure what to do about PRINCIPAL values here.
  916.     (PRINCIP (RES N D #'LOWERHALF #'(lambda (X)
  917.                           (cond ((equal ($imagpart x) 0)  t)
  918.                         (t ())))))))
  919.  
  920. (DEFUN LOWERHALF (J) (EQ ($askSIGN ($imagpart J)) '$neg)) 
  921.  
  922. (DEFUN UPPERHALF (J) (EQ ($askSIGN ($imagpart J)) '$pos)) 
  923.  
  924.  
  925. (DEFUN CSEMIUP (N D VAR)
  926.   (let ((pcprntd t)) ;I'm not sure what to do about PRINCIPAL values here.
  927.     (PRINCIP (RES N D #'UPPERHALF #'(lambda (X)
  928.                       (cond ((equal ($imagpart x) 0)  t)
  929.                         (t ())))))))
  930.  
  931. (DEFUN PRINCIP (N)
  932.        (cond ((null n) nil)
  933.          (t (m*t '$%I ($RECTFORM (m+ (COND ((CAR N)
  934.                         (m*t 2. (CAR N)))
  935.                            (T 0.))
  936.                      (COND ((CADR N)
  937.                         (PRINCIPAL)
  938.                         (CADR N))
  939.                            (T 0.))))))))
  940.  
  941.  
  942. (DEFUN SCONVERT (E)
  943.        (COND ((ATOM E) E)
  944.          ((POLYINX E VAR NIL) E)
  945.          ((EQ (CAAR E) '%SIN)
  946.           (m* '((RAT) -1. 2.)
  947.           '$%I
  948.           (m+t (m^t '$%E (m*t '$%I (CADR E)))
  949.                (m- (m^t '$%E (m*t (m- '$%I) (CADR E)))))))
  950.          ((EQ (CAAR E) '%COS)
  951.           (mul* '((RAT) 1. 2.)
  952.             (m+t (m^t '$%E (m*t '$%I (CADR E)))
  953.              (m^t '$%E (m*t (m- '$%I) (CADR E))))))
  954.          (T (simplify
  955.          (CONS (LIST (CAAR E)) (MAPCAR #'SCONVERT (CDR E)))))))
  956.  
  957. (DEFUN POLFACTORS (EXP)
  958.   (let (poly rest)
  959.        (COND ((mplusp exp)  nil)
  960.          (t (cond ((mtimesp EXP)
  961.                (SETQ EXP (REVERSE (CDR EXP))))
  962.               (T (SETQ EXP (LIST EXP))))
  963.         (mapc #'(lambda (term)
  964.              (cond ((polyinx term var nil)
  965.                 (push term poly))
  966.                    (t (push term rest))))
  967.               exp)
  968.         (list (m*l poly) (m*l rest))))))
  969.  
  970. (DEFUN ESAP (E)
  971.        (PROG (D) 
  972.          (COND ((ATOM E) (RETURN E))
  973.            ((NOT (AMONG '$%E E)) (RETURN E))
  974.            ((AND (mexptp E)
  975.              (EQ (CADR E) '$%E))
  976.             (SETQ D ($imagpart (CADDR E)))
  977.             (RETURN (m* (m^t '$%E ($realpART (CADDR E)))
  978.                    (m+ `((%COS) ,D)
  979.                       (m*t '$%I `((%SIN) ,D))))))
  980.            (T (RETURN (simplify (CONS (LIST (CAAR E))
  981.                           (MAPCAR #'ESAP (CDR E)))))))))
  982.  
  983. (DEFUN MTOSC (GRAND)
  984.   (NUMDEN GRAND)
  985.   (let ((n nn*)
  986.     (d dn*)
  987.     plf bptu bptd s upans downans)
  988.     (COND ((not (or (POLYINX D VAR NIL)
  989.             (AND (SETQ GRAND (%EINVOLVE D))
  990.              (AMONG '$%I GRAND)
  991.              (POLYINX (SETQ D ($RATSIMP (M// D (m^t '$%E GRAND))))
  992.                   VAR
  993.                   NIL)
  994.              (SETQ N (M// N (m^t '$%E GRAND))))))  nil)
  995.       ((EQUAL (SETQ S (DEG D)) 0)  NIL)
  996. ;;;Above tests for applicability of this method.
  997.       ((and (or (SETQ PLF (POLFACTORS N))  t)
  998.         (SETQ N ($EXPAND (COND ((CAR PLF)
  999.                     (m*t 'X* (SCONVERT (CADR PLF))))
  1000.                        (T (SCONVERT N)))))
  1001.         (COND ((mplusp N)  (SETQ N (CDR N)))
  1002.               (T (SETQ N (LIST N))))
  1003.         (do ((do-var n (cdr do-var)))
  1004.             ((null do-var) t)
  1005.           (cond ((rib (car do-var) s))
  1006.             (t (return nil))))
  1007. ;;;Function RIB sets up the values of BPTU and BPTD
  1008.         (COND ((CAR PLF)
  1009.                (SETQ BPTU (SUBST (CAR PLF) 'X* BPTU))
  1010.                (SETQ BPTD (SUBST (CAR PLF) 'X* BPTD))
  1011.                t) ;CROCK, CROCK. This is TERRIBLE code.
  1012.               (t t))
  1013. ;;;If there is BPTU then CSEMIUP must succeed.
  1014. ;;;Likewise for BPTD.
  1015.         (cond (bptu (cond ((setq upans (csemiup (m+l bptu) d var)))
  1016.                   (t nil)))
  1017.               (t (setq upans 0)))
  1018.         (cond (bptd (cond ((setq downans (csemidown (m+l bptd) d var)))
  1019.                   (t nil)))
  1020.               (t (setq downans 0))))
  1021.        (sratsimp (m* '$%pi (m+ upans (m- downans))))))))
  1022.  
  1023.  
  1024. (DEFUN EVENFN (E var)
  1025.        (let ((temp (m+ (m- E) (cond ((atom var)
  1026.                      ($substitute (m- var) var e))
  1027.                     (t ($ratsubst (m- var) var E))))))
  1028.         (cond ((zerop1 temp)
  1029.            t)
  1030.           ((zerop1 ($ratsimp temp))
  1031.            t)
  1032.           (t nil))))
  1033.         
  1034. (DEFUN ODDFN (E VAR)       
  1035.        (let ((temp (m+ e (cond ((atom var)
  1036.                 ($SUBSTITUTE (m- VAR) var E))
  1037.                    (t ($ratsubst (m- var) var e))))))
  1038.         (cond ((zerop1 temp)
  1039.            t)
  1040.           ((zerop1 ($ratsimp temp))
  1041.            t)
  1042.           (t nil))))
  1043.  
  1044. (DEFUN ZTOINF (GRAND VAR)
  1045.        (PROG (N D SN* SD* VARLIST
  1046.         S NC DC
  1047.         ANS R $SAVEFACTORS CHECKFACTORS temp test-var)
  1048.          (SETQ $SAVEFACTORS T SN* (SETQ SD* (LIST 1.)))
  1049.          (COND ((eq ($sign (m+ LOOPSTOP* -1))
  1050.             '$pos)
  1051.             (RETURN NIL))
  1052.            ((SETQ temp (OR (SCAXN GRAND)
  1053.                    (SSP GRAND))) 
  1054.             (RETURN temp))
  1055.            ((INVOLVE GRAND '(%SIN %COS %TAN))
  1056.             (SETQ GRAND (SCONVERT GRAND))
  1057.             (GO ON)))
  1058.  
  1059.          (COND ((POLYINX GRAND VAR NIL)
  1060.             (Diverg))
  1061.            ((AND (RATP GRAND VAR)
  1062.              (mtimesp GRAND)
  1063.              (ANDMAPCAR #'SNUMDEN (CDR GRAND)))
  1064.             (SETQ NN* (M*L SN*)
  1065.               SN* NIL)
  1066.             (SETQ DN* (M*L SD*)
  1067.               SD* NIL))
  1068.            (t (NUMDEN GRAND)))
  1069. ;;;
  1070. ;;;New section.
  1071.     (SETQ N (RMCONST1 NN*))
  1072.     (SETQ D (RMCONST1 DN*))
  1073.     (SETQ NC (CAR N))
  1074.     (SETQ N (CDR N))
  1075.     (SETQ DC (CAR D))
  1076.     (SETQ D (CDR D))
  1077.     (COND ((POLYINX D VAR NIL) 
  1078.        (SETQ S (DEG D)))
  1079.       (T (GO FINDOUT)))
  1080.     (COND ((AND (SETQ R (FINDP N))
  1081.         (eq (ask-integer R '$integer) '$yes)
  1082.         (SETQ test-var (BXM D S))
  1083.         (SETQ ans (APPLY 'FAN (CONS (m+ 1. R) test-var))))
  1084.        (RETURN (m* (M// NC DC) ($RATSIMP ans))))
  1085.       ((and (RATP GRAND VAR)
  1086.         (SETQ ANS (ZMTORAT N (COND ((mtimesp d) d)
  1087.                        (T ($SQFR d)))
  1088.                    S #'ZTORAT)))
  1089.        (RETURN (m* (M// NC DC) ANS)))
  1090.       ((AND (EVENFN D VAR) 
  1091.         (SETQ NN* (P*LOGNXP N S)))
  1092.        (SETQ ANS (LOG*RAT (CAR NN*) D (CADR NN*)))
  1093.        (RETURN (m* (M// NC DC) ANS)))
  1094.       ((INVOLVE GRAND '(%LOG))
  1095.        (COND ((SETQ ANS (LOGQUAD0 GRAND))
  1096.           (RETURN (m* (M// NC DC) ANS)))
  1097.          (T (RETURN NIL)))))
  1098.  FINDOUT
  1099.     (COND ((SETQ temp (BATAPP GRAND)) 
  1100.        (RETURN temp))
  1101.       (T nil))
  1102.   ON
  1103.     (COND ((let ((MTOINF* nil))
  1104.         (SETQ temp (GGR GRAND T)))
  1105.        (RETURN temp))
  1106.       ((mplusp GRAND)
  1107.        (COND ((let ((*NODIVERG t))
  1108.                (SETQ ANS (CATCH 'Divergent
  1109.                      (ANDMAPCAR #'(LAMBDA (G)
  1110.                             (ZTOINF G VAR))
  1111.                             (CDR GRAND)))))
  1112.           (COND ((EQ ANS 'Divergent) NIL)
  1113.             (T (RETURN (sRATSIMP (m+l ANS)))))))))
  1114.  
  1115.     (COND ((AND (EVENFN GRAND VAR)
  1116.         (SETQ LOOPSTOP* (M+ 1 LOOPSTOP*))
  1117.         (SETQ ANS (method-by-limits grand var '$minf '$inf)))
  1118.        (return (m*t '((RAT) 1. 2.) ANS)))
  1119.       (T (RETURN NIL)))))
  1120.    
  1121. (DEFUN ZTORAT (N D S)
  1122.        (COND ((AND (NULL *DFLAG)
  1123.            (SETQ S (DIFAPPLY N D NN* #'ZTORAT)))
  1124.           S)
  1125.          ((SETQ N (let ((PLOGABS ()))
  1126.             (KEYHOLE (m* `((%PLOG) ,(m- VAR)) N) D VAR)))
  1127.           (m- N))
  1128.          (T (MERROR "Keyhole failed"))))
  1129.  
  1130. (SETQ *DFLAG NIL) 
  1131.  
  1132. (DEFUN LOGQUAD0 (EXP)
  1133.   (let ((a ()) (b ())  (c ()))
  1134.     (COND ((SETQ EXP (LOGQUAD EXP))
  1135.        (SETQ A (CAR EXP) B (CADR EXP) C (CADDR EXP))
  1136.        ($asksign b) ;let the data base know about the sign of B.
  1137.        (COND ((EQ ($askSIGN C) '$pos)
  1138.           (SETQ C (m^ (M// C A) '((RAT) 1. 2.)))
  1139.           (setq b (simplify 
  1140.                `((%ACOS) ,(add* 'epsilon (M// B (mul* 2. A C))))))
  1141.           (setq a (M// (m* B `((%LOG) ,C))
  1142.                    (mul* A (Simplify `((%SIN) ,B)) C)))
  1143.           (get-limit a 'epsilon 0 '$plus))))
  1144.       (t ()))))
  1145.     
  1146. (DEFUN LOGQUAD (EXP)
  1147.   (let ((VARLIST (list var)))
  1148.     (NEWVAR EXP)
  1149.     (SETQ EXP (CDR (RATREP* EXP)))
  1150.     (COND ((AND (ALIKE1 (PDIS (CAR EXP))
  1151.             `((%LOG) ,VAR))
  1152.         (NOT (ATOM (CDR EXP)))
  1153.         (EQUAL (CADR (CDR EXP)) 2.)
  1154.         (not (equal (pterm (cddr exp) 0.) 0.)))
  1155.        (SETQ EXP (MAPCAR 'PDIS (CDR (ODDELM (CDR EXP)))))))))
  1156.  
  1157. (DEFUN MTOINF (GRAND VAR) 
  1158.   (PROG (ANS SD* SN* P* PE* N D S NC DC $SAVEFACTORS CHECKFACTORS temp)
  1159.     (SETQ $SAVEFACTORS T)
  1160.     (SETQ SN* (SETQ SD* (LIST 1.)))
  1161.     (COND ((eq ($sign (m+ LOOPSTOP* -1)) '$pos)
  1162.        (RETURN NIL))
  1163.       ((INVOLVE GRAND '(%SIN %COS))
  1164.        (COND ((AND (EVENFN GRAND VAR)
  1165.                (or (SETQ temp (SCAXN GRAND))
  1166.                (setq temp (ssp grand))))
  1167.           (RETURN (m*t 2. temp)))
  1168.          ((SETQ temp (MTOSC GRAND))
  1169.           (RETURN temp))
  1170.          (T (GO EN))))
  1171.       ((among '$%i (%EINVOLVE GRAND))
  1172.        (COND ((SETQ temp (MTOSC GRAND))
  1173.           (RETURN temp))
  1174.          (T (GO EN)))))
  1175.     (COND ((POLYINX GRAND VAR NIL)
  1176.        (Diverg))
  1177.       ((AND (RATP GRAND VAR)
  1178.         (mtimesp GRAND)
  1179.         (ANDMAPCAR #'SNUMDEN (CDR GRAND)))
  1180.        (SETQ NN* (M*L SN*) SN* NIL)
  1181.        (SETQ DN* (M*L SD*) SD* NIL))
  1182.       (t (numden grand)))
  1183.     (SETQ N (RMCONST1 NN*))
  1184.     (SETQ D (RMCONST1 DN*))
  1185.     (SETQ NC (CAR N))
  1186.     (SETQ N (CDR N))
  1187.     (SETQ DC (CAR D))
  1188.     (SETQ D (CDR D))
  1189.     (COND ((POLYINX D VAR NIL)
  1190.        (SETQ S (DEG D))))
  1191.     (COND ((AND (NOT (%einvolve grand))
  1192.         (NOTINVOLVE EXP '(%SINH %COSH %TANH))
  1193.         (SETQ P* (FINDP N))
  1194.         (eq (ask-integer P* '$integer) '$yes)
  1195.         (SETQ PE* (BXM D S)))
  1196.        (COND ((AND (eq (ask-integer (CADDR PE*) '$even) '$yes)
  1197.                (eq (ask-integer P* '$even) '$yes))
  1198.           (COND ((SETQ ANS (APPLY 'FAN (CONS (m+ 1. P*) PE*)))
  1199.              (SETQ ANS (m*t 2. ANS))
  1200.              (RETURN (m* (M// NC DC) ANS)))))
  1201.          ((EQUAL (CAR PE*) 1.)
  1202.           (COND ((AND (SETQ ANS (APPLY 'FAN (CONS (m+ 1. P*) PE*)))
  1203.                   (SETQ NN* (FAN (m+ 1. P*)
  1204.                          (CAR PE*)
  1205.                          (m* -1.(CADR PE*))
  1206.                          (CADDR PE*)
  1207.                          (CADDDR PE*))))
  1208.              (SETQ ANS (m+ ANS (m*t (m^ -1. P*) NN*)))
  1209.              (RETURN (m* (M// NC DC) ANS))))))))
  1210.     (COND ((RATP GRAND VAR)
  1211.        (SETQ ANS (m*t '$%PI (ZMTORAT N (COND ((mtimesp d) d)
  1212.                          (T ($SQFR d)))
  1213.                      S
  1214.                      #'MTORAT)))
  1215.        (RETURN (m* (M// NC DC) ANS)))
  1216.       ((AND (OR (%einvolve grand)
  1217.             (INVOLVE GRAND '(%SINH %COSH %TANH)))
  1218.         (P*PIN%EX N) ;setq's P* and PE*...Barf again.
  1219.         (SETQ ANS (CATCH 'PIN%EX (PIN%EX D))))
  1220.        (COND ((NULL P*)
  1221.           (RETURN (DINTEXP GRAND VAR)))
  1222.          ((AND (ZEROP1 (get-LIMIT GRAND VAR '$INF))
  1223.                (ZEROP1 (get-LIMIT GRAND VAR '$MINF))
  1224.                (SETQ ANS (RECTZTO%PI2 (M*L P*) (M*L PE*) D)))
  1225.           (RETURN (m* (M// NC DC) ANS)))
  1226.          (T (Diverg)))))
  1227.     EN
  1228.     (COND ((SETQ ANS (GGRM GRAND)) 
  1229.        (RETURN ANS))
  1230.       ((AND (EVENFN GRAND VAR)
  1231.         (SETQ LOOPSTOP* (M+ 1 LOOPSTOP*))
  1232.         (SETQ ANS (method-by-limits grand var 0 '$inf)))
  1233.        (RETURN (m*t 2. ANS)))
  1234.       (T (RETURN NIL)))))
  1235.  
  1236. (DEFUN LINPOWER0 (EXP VAR)
  1237.        (COND ((AND (SETQ EXP (LINPOWER EXP VAR))
  1238.            (eq (ask-integer (CADDR EXP) '$even)
  1239.                '$yes)
  1240.            (RATGREATERP 0. (CAR EXP)))
  1241.           EXP))) 
  1242.  
  1243. ;;; given (b*x+a)^n+c returns  (a b n c)
  1244. (DEFUN LINPOWER (EXP VAR)
  1245.  (let (LINPART DEG LC C VARLIST) 
  1246.       (COND ((NOT (POLYP EXP))   NIL)
  1247.         (t (let ((varlist (list var)))
  1248.          (NEWVAR EXP)
  1249.          (SETQ LINPART (CADR (RATREP* EXP)))
  1250.          (COND ((ATOM LINPART)
  1251.             NIL)
  1252.                (t (SETQ DEG (CADR LINPART)) 
  1253. ;;;get high degree of poly
  1254.               (SETQ LINPART ($DIFF EXP VAR (m+ deg -1))) 
  1255. ;;;diff down to linear.
  1256.               (SETQ LC (SDIFF LINPART VAR))    
  1257. ;;;all the way to constant.
  1258.               (SETQ LINPART ($RATSIMP (M// LINPART lc))) 
  1259.               (SETQ LC ($RATSIMP (M// LC `((MFACTORIAL) ,DEG))))
  1260. ;;;get rid of factorial from differentiation.
  1261.               (SETQ C ($RATSIMP (m+ EXP (m* (m- LC)
  1262.                             (m^ LINPART DEG)))))))
  1263. ;;;Sees if can be expressed as (a*x+b)^n + part freeof x.
  1264.          (COND ((NOT (AMONG VAR c))
  1265.             `(,LC ,LINPART ,DEG ,C))
  1266.                (t nil)))))))
  1267.  
  1268. (DEFUN MTORAT (N D S)
  1269.   (let ((SEMIRAT* t)) 
  1270.        (COND ((AND (NULL *DFLAG)
  1271.            (SETQ S (DIFAPPLY N D s #'MTORAT)))
  1272.           S)
  1273.          (T (CSEMIUP N D VAR)))))
  1274.  
  1275. (DEFUN ZMTORAT (N D S FN1)
  1276.  (PROG (C) 
  1277.    (COND ((eq ($sign (m+ S (M+ 1 (SETQ NN* (DEG N))))) 
  1278.           '$neg)
  1279.       (Diverg))
  1280.      ((eq ($sign (m+ s -4))
  1281.           '$neg)
  1282.       (GO ON)))
  1283.    (SETQ D ($FACTOR D))
  1284.    (SETQ C (RMCONST1 D))
  1285.    (SETQ D (CDR C))
  1286.    (SETQ C (CAR C))
  1287.    (COND
  1288.     ((mtimesp D)
  1289.      (SETQ D (CDR D))
  1290.      (SETQ N (PARTNUM N D))
  1291.      (let ((RSN* t))
  1292.       (SETQ N ($XTHRU (M+L
  1293.                (MAPCAR #'(LAMBDA (A B) 
  1294.                     (M// (FUNCALL FN1 (CAR A) B (DEG B)) 
  1295.                          (CADR A)))
  1296.                    N
  1297.                    D)))))
  1298.      (RETURN (COND (C (M// N C)) 
  1299.            (T N)))))
  1300.    ON
  1301.  
  1302.    (SETQ N (FUNCALL FN1 N D S))
  1303.    (RETURN  (sRATSIMP (COND (C  (M// N C))
  1304.                    (T N))))))
  1305.  
  1306. (DEFUN PFRNUM (F G N N2 VAR)
  1307.   (let ((VARLIST (list var))  GENVAR)
  1308.        (SETQ F (POLYFORM F)
  1309.          G (POLYFORM G)
  1310.          N (POLYFORM N)
  1311.          N2 (POLYFORM N2))
  1312.        (SETQ VAR (CAADR (RATREP* VAR)))
  1313.        (SETQ F (RESPROG0 F G N N2))
  1314.        (LIST (LIST (PDIS (CADR F)) (PDIS (CDDR F)))
  1315.          (LIST (PDIS (CAAR F)) (PDIS (CDAR F))))))
  1316.  
  1317. (DEFUN POLYFORM (E)
  1318.        (PROG (F D)
  1319.          (NEWVAR E)
  1320.          (SETQ F (RATREP* E))
  1321.          (AND (EQUAL (CDDR F) 1)
  1322.           (RETURN (CADR F)))
  1323.          (AND (EQUAL (LENGTH (SETQ D (CDDR F))) 3)
  1324.           (NOT (AMONG (CAR D)
  1325.                   (CADR F)))
  1326.           (RETURN (LIST (CAR D)
  1327.                 (f* -1 (CADR D))
  1328.                 (PTIMES (CADR F) (CADDR D)))))
  1329.          (MERROR "Bug from PFRNUM in RESIDU")))
  1330.  
  1331. (DEFUN PARTNUM (N DL)
  1332.   (let ((n2 1)  ANS NL)
  1333.        (do ((dl dl (cdr dl)))
  1334.        ((null (cdr dl))
  1335.         (nconc ans (ncons (list n n2))))
  1336.        (SETQ NL (PFRNUM (CAR DL) (M*L (CDR DL)) N N2 VAR))
  1337.        (SETQ ANS (NCONC ANS (NCONS (CAR NL))))
  1338.        (SETQ N2 (CADADR NL) N (CAADR NL) NL NIL))))
  1339.  
  1340. (DEFUN GGRM (E)
  1341.   (PROG (POLY EXPO MTOINF* MB  VARLIST  GENVAR L C GVAR) 
  1342.     (SETQ VARLIST (LIST VAR))
  1343.     (SETQ MTOINF* T)
  1344.     (COND ((AND (SETQ EXPO (%EINVOLVE E))
  1345.         (POLYP (SETQ POLY ($RATSIMP (M// E (m^t '$%E EXPO)))))
  1346.         (SETQ L (CATCH 'GGRM (GGR (m^t '$%E EXPO) NIL))))
  1347.        (SETQ MTOINF* NIL)
  1348.        (SETQ MB (m- (SUBIN 0. (CADR L))))
  1349.        (SETQ POLY (m+ (SUBIN (m+t MB VAR) POLY)
  1350.               (SUBIN (m+t MB (m*t -1. VAR)) POLY))))
  1351.       (T (RETURN NIL)))
  1352.     (SETQ EXPO (CADDR L)
  1353.       C (CADDDR L)
  1354.       L (m* -1. (CAR L))
  1355.       E NIL)
  1356.     (NEWVAR POLY)
  1357.     (SETQ POLY (CDR (RATREP* POLY)))
  1358.     (SETQ MB (m^ (PDIS (CDR POLY)) -1.) 
  1359.       POLY (CAR POLY))
  1360.     (SETQ GVAR (CAADR (RATREP* VAR)))
  1361.     (COND ((OR (ATOM POLY)
  1362.            (POINTERGP GVAR (CAR POLY))) 
  1363.        (SETQ POLY (LIST 0. POLY)))
  1364.       (T (SETQ POLY (CDR POLY))))
  1365.     (return (do ((poly poly (cddr poly)))
  1366.         ((null poly)
  1367.          (mul* (m^t '$%E C) (m^t EXPO -1.) MB (M+L E)))
  1368.         (SETQ E (CONS (GGRM1 (CAR POLY) (PDIS (CADR POLY)) L EXPO)
  1369.                   E))))))
  1370.  
  1371. (DEFUN GGRM1 (D K A B)
  1372.        (SETQ B (M// (m+t 1. D) B))
  1373.        (m* K `((%GAMMA) ,B) (m^ A (m- B)))) 
  1374.  
  1375. (DEFUN RADIC (E V) 
  1376. ;;;If rd* is t the m^ts must just be free of var.
  1377. ;;;If rd* is () the m^ts must be mnump's.
  1378.        (let ((rd* ())) 
  1379.         (RADICALP E V)))
  1380.  
  1381. (DEFUN KEYHOLE (N D VAR)
  1382.        (let ((SEMIRAT* ()))
  1383.      (SETQ N (RES N D #'(LAMBDA (J) 
  1384.                   (OR (NOT (equal ($imagpart j) 0))
  1385.                   (EQ ($askSIGN J) '$neg)))
  1386.                   #'(LAMBDA (J)
  1387.                   (COND ((EQ ($askSIGN J) '$pos)
  1388.                      T)
  1389.                     (T (Diverg))))))
  1390.      (let ((RSN* t))
  1391.        ($rectform ($multthru (m+ (COND ((CAR N) (CAR N))
  1392.                        (T 0.))
  1393.                      (COND ((CADR N) (CADR N))
  1394.                        (T 0.))))))))
  1395.  
  1396. (DEFUN SKR (E)
  1397.        (PROG (M R K) 
  1398.          (COND ((ATOM E) (RETURN NIL)))
  1399.          (SETQ E (PARTITION E VAR 1))
  1400.          (SETQ M (CAR E))
  1401.          (SETQ E (CDR E))
  1402.          (COND ((SETQ R (SINRX E)) (RETURN (LIST M R 1)))
  1403.            ((AND (mexptp E)
  1404.              (eq (ask-integer (SETQ K (CADDR E)) '$integer) '$yes)
  1405.              (SETQ R (SINRX (CADR E))))
  1406.             (RETURN (LIST M R K)))))) 
  1407.  
  1408. (DEFUN SINRX (E)
  1409.        (COND ((EQ (CAAR E) '%SIN)
  1410.           (COND ((EQ (CADR E) VAR) 1.)
  1411.             ((AND (SETQ E (PARTITION (CADR E) VAR 1)) (EQ (CDR E) VAR))
  1412.              (CAR E)))))) 
  1413.  
  1414. (DECLARE-TOP(SPECIAL N)) 
  1415.  
  1416.  
  1417. (DEFUN SSP (EXP)
  1418.  (PROG (U N C) 
  1419.        (setq exp ($substitute (m^t `((%sin) ,var) 2.)
  1420.                   (m+t 1. (m- (m^t `((%cos) ,var) 2.)))
  1421.                   exp))
  1422.        (numden exp)
  1423.        (setq u nn*)
  1424.        (COND ((AND (SETQ N (FINDP DN*))
  1425.            (eq (ask-integer N '$integer) '$yes))
  1426.           (COND ((SETQ C (SKR U)) 
  1427.              (RETURN (SCMP C N)))
  1428.             ((AND (mplusp U)
  1429.               (SETQ C (ANDMAPCAR #'SKR (CDR U))))
  1430.              (RETURN (m+l (MAPCAR #'(LAMBDA (J) (SCMP J N))
  1431.                       C)))))))))
  1432.  
  1433. (DECLARE-TOP(UNSPECIAL N)) 
  1434.  
  1435. (DEFUN SCMP (C N)
  1436.        (m* (CAR C) (m^ (CADR C) (m+ N -1)) `((%SIGNUM) ,(CADR C))
  1437.        (SINSP (CADDR C) N)))
  1438.  
  1439. (DEFUN SEVN (N)
  1440.        (m* HALF%PI ($MAKEGAMMA `((%BINOMIAL) ,(m+t (m+ N -1) '((RAT) -1. 2.))
  1441.                          ,(m+ N -1)))))
  1442.  
  1443.  
  1444. (DEFUN SFORX (N)
  1445.   (COND ((EQUAL N 1.) 
  1446.      HALF%PI) 
  1447.     (T (BYGAMMA (m+ N -1) 0.)))) 
  1448.  
  1449. (DEFUN SINSP (L K)
  1450.   (let ((I ())
  1451.     (J ()))
  1452.     (COND ((eq ($sign (m+ L (m- (m+ K -1))))
  1453.            '$neg)
  1454.        (Diverg))
  1455.       ((NOT (EVEN1 (m+ L K)))
  1456.        NIL)
  1457.       ((EQUAL K 2.)
  1458.        (SEVN (m// L 2.)))
  1459.       ((EQUAL K 1.) 
  1460.        (SFORX L))
  1461.       ((eq ($sign  (m+ K -2.))
  1462.            '$pos)
  1463.        (SETQ I (m* (m+ K -1) (SETQ J (m+ K -2.))))
  1464.        (m+ (m* L (m+ L -1) (m^t I -1.) (SINSP (m+ L -2.) J))
  1465.            (m* (m- (m^ L 2)) (m^t I -1.)
  1466.            (SINSP L J)))))))
  1467.  
  1468. (DEFUN FPART (A)
  1469.   (COND ((NULL A) 0.)
  1470.     ((NUMBERP A) 0.)
  1471.     ((MNUMP A)
  1472.      (LIST (CAR A) (REMAINDER (CADR A) (CADDR A)) (CADDR A)))
  1473.     ((AND (ATOM A) (ABLESS1 A)) A)
  1474.     ((AND (mplusp A)
  1475.           (NULL (CDDDR A))
  1476.           (ABLESS1 (CADDR A)))
  1477.      (CADDR A)))) 
  1478.  
  1479. (DEFUN THRAD (E)
  1480.        (COND ((POLYINX E VAR NIL) 0.)
  1481.          ((AND (mexptp E) 
  1482.            (EQ (CADR E) VAR)
  1483.            (MNUMP (CADDR E)))
  1484.           (FPART (CADDR E)))
  1485.          ((mtimesp E)
  1486.           (m+l (MAPCAR #'THRAD E)))))
  1487.  
  1488.  
  1489. ;;; THE FOLLOWING FUNCTION IS FOR TRIG FUNCTIONS OF THE FOLLOWING TYPE:
  1490. ;;; LOWER LIMIT=0 B A MULTIPLE OF %PI SCA FUNCTION OF SIN (X) COS (X)
  1491. ;;; B<=%PI2
  1492.  
  1493. (DEFUN PERIOD (P E VAR)
  1494.   (and (ALIKE1 (no-err-sub var E) (setq e (NO-ERR-SUB (m+ P VAR) E)))
  1495.        ;; means there was no error
  1496.        (not (eq e t))))
  1497.  
  1498. (DEFUN INFR (A)
  1499.   (let ((var '$%i)
  1500.     (r (subin 0. a))
  1501.     c)
  1502.        (SETQ C (SUBIN 1. (m+ A (m*t -1. R))))
  1503.        (SETQ A (IGPRT (m* '((RAT) 1. 2.) C)))
  1504.        (CONS A (m+ R (m*t (M+ C (m* -2. A)) '$%PI)))))
  1505.  
  1506. (DEFUN IGPRT (R) 
  1507.        (M+ R (m* -1. (FPART R)))) 
  1508.  
  1509.  
  1510. ;;;Try making exp(%i*var) --> yy, if result is rational then do integral
  1511. ;;;around unit circle. Make corrections for limits of integration if possible.
  1512. (DEFUN SCRAT (SC B)
  1513.   (let* ((exp-form (sconvert sc)) ;Exponentialize
  1514.      (rat-form (MAXIMA-SUBSTITUTE 'yy (m^t '$%e (m*t '$%i var))
  1515.                    exp-form))) ;Try to make Rational fun.
  1516.     (COND ((AND (RATP rat-form 'YY)
  1517.             (NOT (AMONG VAR rat-form)))
  1518.            (COND ((ALIKE1 B %PI2) 
  1519.               (let ((ans (ZTO%PI2 rat-form 'YY)))
  1520.                (cond (ans ans)
  1521.                  (t nil))))
  1522.              ((AND (EQ B '$%PI)
  1523.                (EVENFN exp-form VAR))
  1524.               (let ((ans (ZTO%PI2 rat-form 'YY)))
  1525.                (cond (ans (m*t '((RAT) 1. 2.) ans))
  1526.                  (t nil))))
  1527.              ((AND (ALIKE1 B HALF%PI)
  1528.                (EVENFN exp-form VAR)
  1529.                (ALIKE1 rat-form 
  1530.                    (NO-ERR-SUB (m+t '$%PI (m*t -1. VAR))
  1531.                             rat-form)))
  1532.                (let ((ans (ZTO%PI2 rat-form 'yy)))
  1533.                 (cond (ans (m*t '((RAT) 1. 4.) ans))
  1534.                       (t nil)))))))))
  1535.  
  1536. ;;; Do integrals of sin and cos. this routine makes sure lower limit
  1537. ;;; is zero.
  1538. (defun INTSC1 (A B E)
  1539.   (let ((limit-diff (m+ b (m* -1 a)))
  1540.     ($%emode t)
  1541.     ($trigsign t)
  1542.     (sin-cos-recur t)) ;recursion stopper
  1543.      (PROG (ANS D NZP2 L) 
  1544.        (COND ((OR (NOT (MNUMP (M// limit-diff '$%PI)))
  1545.               (NOT (PERIOD %PI2 E VAR)))
  1546.           (RETURN NIL))
  1547.          ((not (equal a 0.))
  1548.           (setq e (MAXIMA-SUBSTITUTE (m+ a var) var e))
  1549.           (setq a 0.)
  1550.           (setq b limit-diff)))
  1551. ;;;Multiples of 2*%pi in limits.
  1552.        (COND ((eq (ask-integer (SETQ D (let (($float nil))
  1553.                          (m// limit-diff %pi2))) '$integer)
  1554.               '$yes)
  1555.           (SETQ ANS (m* D (COND ((SETQ ANS (INTSC E %PI2 VAR))
  1556.                      (RETURN ans))
  1557.                     (T (RETURN NIL)))))))
  1558.           (COND ((RATGREATERP %PI2 B)
  1559.              (RETURN (INTSC E B VAR)))
  1560.             (T (SETQ L A) 
  1561.                (SETQ A 0.)))
  1562.           (SETQ B (INFR B))
  1563.           (COND ((NULL L) 
  1564.              (SETQ NZP2 (CAR B))
  1565.              (SETQ limit-diff 0.)
  1566.              (GO OUT)))
  1567.           (SETQ L (INFR L))
  1568.           (SETQ limit-diff
  1569.             (m*t -1. (COND ((SETQ ANS (INTSC E (CDR L) VAR)) 
  1570.                     ANS)
  1571.                        (T (RETURN NIL)))))
  1572.           (SETQ NZP2 (m+ (CAR B) (m- (CAR L))))
  1573.              OUT  (SETQ ANS (add* (COND ((ZEROP1 NZP2) 0.)
  1574.                     ((SETQ ANS (INTSC E %PI2 VAR))
  1575.                      (m*t NZP2 ANS))
  1576.                     (T (RETURN NIL)))
  1577.                   (COND ((ZEROP1 (CDR B)) 0.)
  1578.                     ((SETQ ANS (INTSC E (CDR B) VAR))
  1579.                      ANS)
  1580.                     (T (RETURN NIL)))
  1581.                   limit-diff))
  1582.          (RETURN ANS))))
  1583.  
  1584. (DEFUN INTSC (SC B VAR)
  1585.        (COND ((EQ ($sign B) '$neg)
  1586.           (SETQ B (m*t -1. B))
  1587.           (SETQ SC (m* -1. (SUBIN (m*t -1. VAR) SC)))))
  1588.        (SETQ SC (PARTITION SC VAR 1))
  1589.        (COND ((SETQ B (INTSC0 (CDR SC) B VAR))
  1590.           (m* (resimplify (CAR SC)) B))))
  1591.  
  1592. (DEFUN INTSC0 (SC B VAR)
  1593.   (let ((nn* (scprod sc))
  1594.     (dn* ()))
  1595.        (COND (NN* (COND ((ALIKE1 B HALF%PI)
  1596.              (BYGAMMA (CAR NN*) (CADR NN*)))
  1597.             ((EQ B '$%PI)
  1598.              (cond ((eq (real-branch (cadr nn*) -1.) '$yes)
  1599.                  (m* (m+ 1. (m^ -1 (CADR NN*)))
  1600.                      (BYGAMMA (car nn*) (cadr nn*))))))
  1601.             ((ALIKE1 B %PI2)
  1602.              (COND ((or (AND (eq (ask-integer (CAR NN*) '$even)
  1603.                          '$yes)
  1604.                      (eq (ask-integer (CADR NN*) '$even)
  1605.                          '$yes))
  1606.                      (and (ratnump (car nn*))
  1607.                       (eq (real-branch (car nn*) -1.)
  1608.                           '$yes)
  1609.                       (ratnump (cadr nn*))
  1610.                       (eq (real-branch (cadr nn*) -1.)
  1611.                           '$yes)))
  1612.                  (m* 4.    (BYGAMMA (car nn*) (cadr nn*))))
  1613.                 ((or (eq (ask-integer (car nn*) '$odd) '$yes)
  1614.                      (eq (ask-integer (cadr nn*) '$odd) '$yes))
  1615.                  0.)
  1616.                 (t nil)))
  1617.              ((ALIKE1 B HALF%PI3)
  1618.               (m* (m+ 1. (m^ -1 (CADR NN*)) (m^ -1 (m+l NN*)))
  1619.                   (bygamma (car nn*) (cadr nn*))))))
  1620.          (t (cond ((AND (OR (EQ B '$%PI)
  1621.                 (ALIKE1 B %PI2)
  1622.                 (ALIKE1 B HALF%PI))
  1623.                 (SETQ DN* (SCRAT SC B)))
  1624.                DN*)
  1625.               ((SETQ NN* (ANTIDERIV SC))
  1626.                (sin-cos-intsubs nn* var 0. b))
  1627.               (t ()))))))
  1628.  
  1629. ;;;Is careful about substitution of limits where the denominator may be zero
  1630. ;;;because of various assumptions made.
  1631. (defun sin-cos-intsubs (exp var ll ul)
  1632.   (cond ((mplusp exp)
  1633.      (m+l (mapcar #'sin-cos-intsubs1 (cdr exp))))
  1634.     (t (sin-cos-intsubs1 exp))))
  1635.  
  1636. (defun sin-cos-intsubs1 (exp)     
  1637.   (let* ((rat-exp ($rat exp))
  1638.      (num (pdis (cadr rat-exp)))
  1639.      (denom (pdis (cddr rat-exp))))
  1640.     (cond ((not (equal (intsubs num ll ul) 0.))
  1641.        (intsubs exp ll ul))
  1642.       ((not (equal ($asksign denom) '$zero))
  1643.        0.)
  1644.       (t (let (($%piargs ()))
  1645.            (intsubs exp ll ul))))))
  1646.  
  1647. (DEFUN SCPROD (E)
  1648.  (let ((great-minus-1 #'(lambda (temp)
  1649.               (ratgreaterp temp -1)))
  1650.        m n)
  1651.    (COND
  1652.     ((SETQ M (Powerofx E `((%SIN) ,VAR) great-minus-1 VAR))
  1653.      (list m 0.))
  1654.     ((SETQ N (Powerofx E `((%COS) ,VAR) great-minus-1 VAR))
  1655.      (SETQ M 0.)
  1656.      (list 0. n))
  1657.     ((AND (mtimesp E)
  1658.       (OR (SETQ M (Powerofx (CADR E) `((%SIN) ,VAR) great-minus-1 VAR))
  1659.           (SETQ N (Powerofx (CADR E) `((%COS) ,VAR) great-minus-1 VAR)))
  1660.       (COND
  1661.        ((NULL M)
  1662.         (SETQ M (Powerofx (CADDR E) `((%SIN) ,VAR) great-minus-1 VAR)))
  1663.        (T (SETQ N (Powerofx (CADDR E) `((%COS) ,VAR) great-minus-1 VAR))))
  1664.       (NULL (CDDDR E)))
  1665.      (list m n))
  1666.     (T ()))))
  1667.  
  1668. (defun real-branch (exponent value)
  1669. ;;;Says wether (m^t value exponent) has at least one real branch.
  1670. ;;;Only works for values of 1 and -1 now.
  1671. ;;;Returns $yes $no $unknown.
  1672.        (cond ((equal value 1.)
  1673.           '$yes)
  1674.          ((eq (ask-integer exponent '$integer) '$yes)
  1675.           '$yes)
  1676.          ((ratnump exponent)
  1677.           (cond ((eq ($oddp (caddr exponent)) t)
  1678.              '$yes)
  1679.             (t '$no)))
  1680.          (t '$unknown)))
  1681.  
  1682. (DEFUN BYGAMMA (M N)
  1683.        (let ((one-half (m//t 1. 2.)))
  1684.        (m* one-half `(($BETA) ,(m* one-half (m+t 1. M))
  1685.                   ,(m* one-half (m+t 1. N))))))
  1686.  
  1687. ;Seems like Guys who call this don't agree on what it should return.
  1688. (DEFUN Powerofx (E X P VAR)
  1689.        (SETQ E (COND ((NOT (AMONG VAR E)) NIL)
  1690.              ((ALIKE1 E X) 1.)
  1691.              ((ATOM E) NIL)
  1692.              ((AND (mexptp E)
  1693.                (ALIKE1 (CADR E) X)
  1694.                (NOT (AMONG VAR (CADDR E))))
  1695.               (CADDR E))))
  1696.        (COND ((NULL E) NIL)
  1697.          ((FUNCALL P E) E))) 
  1698.  
  1699. (DECLARE-TOP(SPECIAL L C K)) 
  1700.  
  1701. (COMMENT (THE FOLLOWING FUNC IS NOT COMPLETE)) 
  1702.  
  1703. ;(DEFUN ZTO1 (E)
  1704. ;  (prog (ans k l)
  1705. ;    (COND ((NOTINVOLVE E '(%SIN %COS %TAN %LOG))
  1706. ;       (cond ((SETQ ANS (BATAP E))
  1707. ;          (return ans)))))
  1708. ;    (cond ((AND (NOTINVOLVE E '(%SIN %COS %TAN))
  1709. ;        (AMONG '%LOG E))
  1710. ;       (COND ((SETQ ANS (BATAP (M// E `((%LOG) ,VAR))))
  1711. ;          (SETQ K NN* L DN*)
  1712. ;          (SETQ ANS (m* ANS
  1713. ;                (m+ (subfunmake '$PSI '(0) (list K))
  1714. ;                    (m* -1. (subfunmake '$PSI
  1715. ;                            '(0)
  1716. ;                            (ncons (m+ K
  1717. ;                                   L)))))))
  1718. ;          (return ans)))))))
  1719.  
  1720.  
  1721.  
  1722. (DEFUN BATA0 (E)
  1723.   (COND ((ATOM E) NIL)
  1724.     ((AND (mtimesp E)
  1725.           (NULL (CDDDR E))
  1726.           (OR (AND (SETQ K (FINDP (CADR E)))
  1727.                (SETQ C (BXM (CADDR E) (POLYINX (CADDR E) VAR NIL))))
  1728.           (AND (SETQ K (FINDP (CADDR E)))
  1729.                (SETQ C (BXM (CADR E) (POLYINX (CADR E) VAR NIL))))))
  1730.      T)
  1731.     ((SETQ C (BXM E (POLYINX E VAR NIL)))
  1732.      (SETQ K 0.)))) 
  1733.  
  1734. ;(DEFUN BATAP (E) 
  1735. ;  (PROG (K C L) 
  1736. ;    (COND ((NOT (BATA0 E)) (RETURN NIL))
  1737. ;      ((AND (EQUAL -1. (CADDDR C))
  1738. ;        (EQ ($askSIGN (SETQ K (m+ 1. K)))
  1739. ;            '$pos)
  1740. ;        (EQ ($askSIGN (SETQ L (m+ 1. (CAR C))))
  1741. ;            '$pos)
  1742. ;        (ALIKE1 (CADR C)
  1743. ;            (m^ UL (CADDR C)))
  1744. ;        (SETQ E (CADR C))
  1745. ;        (EQ ($askSIGN (SETQ C (CADDR C))) '$pos))
  1746. ;       (RETURN (M// (m* (m^ UL (m+t K (m* C (m+t -1. L))))
  1747. ;                `(($BETA) ,(SETQ NN* (M// K C))
  1748. ;                      ,(SETQ DN* L)))
  1749. ;            C))))))
  1750.  
  1751.  
  1752. ;;; Integrals of the form i(log(x)^m*x^k*(a+b*x^n)^l,x,0,ul) with
  1753. ;;; ul>0, b<0, a=-b*ul^n, k>-1, l>-1, n>0, m a nonnegative integer. 
  1754. ;;; These integrals are essentially partial derivatives of the 
  1755. ;;; Beta function (i.e. the Eulerian integral of the first kind).
  1756. ;;; Note that, currently, with the default setting intanalysis:true,
  1757. ;;; this function might not even be called for some of these integrals.
  1758. ;;; However, this can be palliated by setting intanalysis:false. 
  1759.  
  1760. (defun zto1 (e)                
  1761.   (when (or (mtimesp e) (mexptp e))
  1762.     (let ((m 0) (log (list '(%log) var)))
  1763.       (flet ((set-m (p) (setq m p)))
  1764.         (find-if #'(lambda (fac) (powerofx fac log #'set-m var)) (cdr e)))
  1765.       (when (and (freeof var m) 
  1766.          (eq (ask-integer m '$integer) '$yes)
  1767.          (not (eq ($asksign m) '$neg))) 
  1768.     (setq e (m//t e (list '(mexpt) log m)))
  1769.     (multiple-value-bind
  1770.      (k/n l n b) (batap-new e)
  1771.      (when k/n
  1772.        (let ((beta (simplify (list '($beta) k/n l)))
  1773.          (m (if (eq ($asksign m) '$zero) 0 m)))
  1774.          ;; The result looks like B(k/n,l) ( ... ).
  1775.          ;; Perhaps, we should just $factor, instead of
  1776.          ;; pulling out beta like this.
  1777.          (m*t
  1778.           beta
  1779.           ($fullratsimp
  1780.            (m//t
  1781.         (m*t
  1782.          (m^t (m-t b) (m1-t l))
  1783.          (m^t ul (m*t n (m1-t l)))
  1784.          (m^t n (m-t (m1+t m)))
  1785.          ($at ($diff (m*t (m^t ul (m*t n var)) (list '($beta) var l))
  1786.                  var m) (list '(mequal) var k/n))) 
  1787.         beta))))))))))
  1788.  
  1789.  
  1790. ;;; If e is of the form given below, make the obvious change
  1791. ;;; of variables (substituting ul*x^(1/n) for x) in order to reduce
  1792. ;;; e to the usual form of the integrand in the Eulerian
  1793. ;;; integral of the first kind.
  1794. ;;; N. B: The old version of ZTO1 completely ignored this  
  1795. ;;; substitution; the log(x)s were just thrown in, which,
  1796. ;;; of course would give wrong results.
  1797.  
  1798. (defun batap-new (e) 
  1799.   (let (k c) 
  1800.     (declare (special k c))
  1801.     ;; Parse e
  1802.     (when (bata0 e)
  1803.       (multiple-value-bind
  1804.        ;; e=x^k*(a+b*x^n)^l
  1805.        (l a n b)  (values-list c)
  1806.        (when (and (freeof var k) (freeof var n) (freeof var l)
  1807.           (alike1 a (m-t (m*t b (m^t ul n))))
  1808.           (eq ($asksign b) '$neg)
  1809.           (eq ($asksign (setq k (m1+t k))) '$pos)
  1810.           (eq ($asksign (setq l (m1+t l))) '$pos)
  1811.           (eq ($asksign n) '$pos))
  1812.      (values (m//t k n) l n b))))))
  1813.  
  1814.  
  1815.  
  1816. (DEFUN BATAPP (E)
  1817.        (PROG (K C D L AL) 
  1818.          (COND ((NOT (OR (EQUAL LL 0) (EQ LL '$MINF)))
  1819.             (SETQ E (SUBIN (m+ LL VAR) E))))
  1820.          (COND ((NOT (BATA0 E)) (RETURN NIL))
  1821.            ((AND (RATGREATERP (SETQ AL (CADDR C)) 0.)
  1822.              (EQ ($askSIGN (SETQ K (M// (m+ 1. K)
  1823.                                AL)))
  1824.                  '$pos)
  1825.              (RATGREATERP (SETQ L (m* -1. (CAR C)))
  1826.                       K)
  1827.              (EQ ($askSIGN (m* (SETQ D (CADR C))
  1828.                        (SETQ C (CADDDR C))))
  1829.                  '$pos))
  1830.             (SETQ L (m+ L (m*t -1. K)))
  1831.             (RETURN (M// `(($BETA) ,K ,L)
  1832.                   (mul* AL (m^ C K) (m^ D L)))))))) 
  1833.  
  1834. (DECLARE-TOP(UNSPECIAL L C K)) 
  1835.  
  1836. (DEFUN GAMMA1 (C A B D)
  1837.        (m* (m^t '$%E D)
  1838.        (m^ (m* B (m^ A (SETQ C (M// (m+t C 1.) B))))
  1839.            -1.)
  1840.        `((%GAMMA) ,C)))
  1841.        
  1842. (DEFUN ZTO%PI2 (GRAND VAR)
  1843.        (let ((result (UNITCIR ($RATSIMP (M// GRAND VAR)) VAR)))
  1844.         (cond (result (sratsimp (m* (m- '$%I) result)))
  1845.           (t nil))))
  1846.  
  1847. (DEFUN UNITCIR (GRAND VAR)
  1848.   (NUMDEN GRAND)
  1849.   (let ((result (PRINCIP (RES NN* DN* #'(LAMBDA (PT)
  1850.                       (RATGREATERP 1 (CABS PT)))
  1851.                           #'(LAMBDA (PT)
  1852.                       (ALIKE1 1 (CABS PT)))))))
  1853.        (cond (result (m* '$%PI result))
  1854.          (t nil))))
  1855.  
  1856.  
  1857. (DEFUN LOGX1 (EXP LL UL)
  1858.   (let ((arg nil))
  1859.     (COND
  1860.      ((AND (NOTINVOLVE EXP '(%SIN %COS %TAN %ATAN %ASIN %ACOS))
  1861.        (SETQ ARG (INVOLVE EXP '(%LOG))))
  1862.       (COND ((EQ ARG VAR)
  1863.          (COND ((RATGREATERP 1. LL)
  1864.             (COND ((NOT (EQ UL '$INF))
  1865.                (INTCV1 (m^t '$%E (m- VAR)) () (m- `((%LOG) ,VAR))))
  1866.               (T (INTCV1 (m^t '$%E VAR) () `((%LOG) ,VAR)))))))
  1867.         (t (INTCV ARG NIL nil)))))))
  1868.  
  1869.  
  1870. (DEFUN SCAXN (E)
  1871.   (let (IND S G) 
  1872.        (COND ((ATOM E)  NIL)
  1873.          ((AND (OR (EQ (CAAR E) '%SIN)
  1874.                (EQ (CAAR E) '%COS))
  1875.            (SETQ IND (CAAR E))
  1876.            (SETQ E (BX**N (CADR E))))
  1877.           (COND ((EQUAL (CAR E) 1.)  '$IND)
  1878.             ((ZEROP (SETQ S (let ((sign ($askSIGN (CADR E))))
  1879.                      (cond ((eq sign '$pos) 1)
  1880.                            ((eq sign '$neg) -1)
  1881.                            ((eq sign '$zero) 0)))))
  1882.              NIL)
  1883.             ((not (EQ ($askSIGN (m+ -1 (CAR E)))  '$pos))
  1884.              nil)
  1885.             (t (SETQ G (GAMMA1 0. (m* S (CADR E)) (CAR E) 0.))
  1886.                (SETQ E (m* G `((,IND) ,(M// HALF%PI (CAR E))))) 
  1887.                (m* (COND ((AND (EQ IND '%SIN)
  1888.                        (EQUAL S -1.))
  1889.                   -1.)
  1890.                  (T 1.))
  1891.                E)))))))
  1892.               
  1893.  
  1894. (COMMENT THIS IS THE SECOND PART OF THE DEFINITE INTEGRAL PACKAGE) 
  1895.  
  1896. (DECLARE-TOP(SPECIAL VAR PLM* PL* RL* PL*1 RL*1)) 
  1897.  
  1898. (DEFUN P*LOGNXP (A S)
  1899.   (let (B) 
  1900.     (COND ((NOT (AMONG '%LOG A)) 
  1901.        ())
  1902.       ((AND (POLYINX (SETQ B (MAXIMA-SUBSTITUTE 1. `((%LOG) ,VAR) A))
  1903.              VAR T)
  1904.         (eq ($sign (m+ S (M+ 1 (DEG B))))
  1905.             '$pos)
  1906.         (EVENFN B VAR)
  1907.         (SETQ A (LOGNXP ($RATSIMP (M// A B)))))
  1908.        (LIST B A)))))
  1909.  
  1910. (DEFUN LOGNXP (A)
  1911.    (COND ((ATOM A) NIL)
  1912.      ((AND (EQ (CAAR A) '%LOG) 
  1913.            (EQ (CADR A) VAR)) 1.)
  1914.      ((AND (mexptp A)
  1915.            (NUMBERP (CADDR A))
  1916.            (LOGNXP (CADR A)))
  1917.       (CADDR A)))) 
  1918.  
  1919. (COMMENT CHECK THE FOLLOWING FUNCTION FOR UNUSED PROG VAR A) 
  1920.  
  1921. (DEFUN LOGCPI0 (N D)
  1922.   (PROG (PL DP) 
  1923.     (SETQ PL (POLELIST D #'UPPERHALF #'(LAMBDA (J)
  1924.                          (COND ((ZEROP1 J) NIL)
  1925.                            ((equal ($imagpart j) 0)
  1926.                             T)))))
  1927.     (cond ((null pl)
  1928.            (return nil)))
  1929.     (SETQ FACTORS (CAR PL) 
  1930.           PL (CDR PL))
  1931.     (COND ((OR (CADR PL)
  1932.            (CADDR PL))
  1933.            (SETQ DP (SDIFF D VAR))))
  1934.     (COND ((SETQ PLM* (CAR PL))
  1935.            (SETQ RLM* (RESIDUE N (COND (LEADCOEF FACTORS)
  1936.                        (T D))
  1937.                    PLM*))))
  1938.     (COND ((SETQ PL* (CADR PL)) 
  1939.            (SETQ RL* (RES1 N DP PL*))))
  1940.     (COND ((SETQ PL*1 (CADDR PL))
  1941.            (SETQ RL*1 (RES1 N DP PL*1))))
  1942.     (RETURN (m*t (m//t 1. 2.)
  1943.              (m*t '$%PI 
  1944.               (PRINCIP 
  1945.                (LIST (COND ((SETQ NN* (APPEND RL* RLM*))
  1946.                     (M+L NN*)))
  1947.                  (COND (RL*1 (M+L RL*1))))))))))
  1948.  
  1949. (DEFUN LOGNX2 (NN DN PL RL)
  1950.   (do ((pl pl (cdr pl))
  1951.        (rl rl (cdr rl))
  1952.        (ans ()))
  1953.       ((or (null pl)
  1954.        (null rl))  ans)
  1955.     (SETQ ANS (CONS (m* DN (car rl) (m^ `((%PLOG) ,(car pl)) NN))
  1956.             ANS))))
  1957.  
  1958. (DEFUN LOGCPJ (N D I)
  1959.  (SETQ N (APPEND
  1960.       (COND (PLM* (LIST (mul* (m*t '$%i %pi2)
  1961.                   (M+L
  1962.                    (RESIDUE (m* (m^ `((%PLOG) ,VAR) I)
  1963.                         N)
  1964.                         D
  1965.                         PLM*))))))
  1966.       (LOGNX2 I (m*t '$%i %pi2) PL* RL*)
  1967.       (LOGNX2 I %P%I PL*1 RL*1)))
  1968.  (COND ((NULL N) 0)
  1969.        (T (SIMPLIFY (M+L N)))))
  1970.  
  1971.  
  1972. ;;should replace these references to *i* and *j* to symbol-value arrays.
  1973. ;;here and in SUMI, and LOGCPI.  These are the only references in this file.
  1974. ;;I did change I to *I* 
  1975.  
  1976. #-cl ;in case other lisps don't understand internal declares.
  1977. (declare-top(special *i* *j*))
  1978.  
  1979. (DEFUN LOG*RAT (N D M)
  1980.   (PROG (LEADCOEF FACTORS C PLM* PL* RL* PL*1 RL*1 RLM*)
  1981.     (declare (special *i* *j*))
  1982. ;    (ARRAY *I* T (M+ 1 M))
  1983. ;    (ARRAY *J* T (M+ 1 M))
  1984.     (setq *i* (*ARRAY nil T (M+ 1 M)))
  1985.     (setq *j* (*ARRAY nil T (M+ 1 M)))
  1986.     (SETQ C 0.)
  1987.     (STORE (aref *J* C) 0.)
  1988.     (do ((c 0. (m+ 1 c)))
  1989.         ((equal c m)
  1990.          (return (logcpi n d m)))
  1991.      (STORE (aref *I* C) (LOGCPI N D C))
  1992.      (STORE (aref *J* C) (LOGCPJ N FACTORS C)))))
  1993.  
  1994. (DEFUN LOGCPI (N D C)
  1995.   (declare (special *j*))
  1996.   (COND ((EQUAL C 0.)
  1997.      (LOGCPI0 N D))
  1998.     (T (m* '((RAT) 1. 2.)
  1999.            (m+ (aref *J* C) (m* -1. (SUMI C)))))))
  2000.  
  2001. (DEFUN SUMI (C)
  2002.   (declare (special *i*))
  2003.   (do ((k 1. (m+ 1 k))
  2004.        (ans ()))
  2005.       ((equal k c)
  2006.        (m+l ans))
  2007.      (SETQ ANS (CONS (mul* ($MAKEGAMMA `((%BINOMIAL) ,C ,K))
  2008.                (m^t '$%PI K)
  2009.                (m^t '$%I K)
  2010.                (aref *I* (m+ C (m- K))))
  2011.              ANS))))
  2012.  
  2013. (DEFUN FAN (P M A N B)
  2014.   (let ((povern (m// p n))
  2015.     (ab (m// a b)))
  2016.     (COND
  2017.      ((OR (eq (ask-integer POVERN '$integer) '$yes)
  2018.       (NOT (equal ($imagpart ab) 0)))  ())
  2019.      (t (let ((IND ($askSIGN AB)))
  2020.       (COND ((eq IND '$zero) NIL)
  2021.         ((eq ind '$neg) nil)
  2022.         ((NOT (RATGREATERP M POVERN)) NIL)
  2023.         (t (M// (m* '$%pi
  2024.                 ($MAKEGAMMA `((%BINOMIAL) ,(m+ -1. M (m- POVERN))
  2025.                               ,(m+t -1. M)))
  2026.                 `((mabs) ,(m^ A (m+ POVERN (m- m)))))
  2027.             (m* (m^ b povern)
  2028.                 N
  2029.                 `((%SIN) ,(m*t '$%PI POVERN)))))))))))
  2030.  
  2031.  
  2032. ;;Makes a new poly such that np(x)-np(x+2*%i*%pi)=p(x).
  2033. ;;Constructs general POLY of degree one higher than P with
  2034. ;;arbitrary coeff. and then solves for coeffs by equating like powers
  2035. ;;of the varibale of integration.
  2036. ;;Can probably be made simpler now.
  2037.  
  2038. (DEFUN MAKPOLY (P)
  2039.   (let ((n (deg p))  (ans ())  (varlist ())  (gp ())  (cl ())  (zz ()))
  2040.     (SETQ ANS (GENPOLY (M+ 1 N))) ;Make poly with gensyms of 1 higher deg.
  2041.     (setq cl (cdr ans)) ;Coefficient list
  2042.     (SETQ VARLIST (APPEND CL (LIST VAR))) ;Make VAR most important.
  2043.     (SETQ GP (CAR ANS)) ;This is the poly with gensym coeffs.
  2044. ;;;Now, poly(x)-poly(x+2*%i*%pi)=p(x), P is the original poly.
  2045.     (SETQ ANS (m+ GP (SUBIN (m+t (m*t '$%i %pi2) VAR) (m- GP)) (m- P)))
  2046.     (NEWVAR ANS)
  2047.     (SETQ ANS (RATREP* ANS)) ;Rational rep with VAR leading.
  2048.     (SETQ ZZ (COEFSOLVE N CL (COND ((NOT (EQ (CAADR ANS) ;What is Lead Var.
  2049.                          (GENFIND (CAR ANS) VAR)))
  2050.                     (LIST 0 (CADR ANS)));No VAR in ans.
  2051.                    ((CDADR ANS)))));The real Poly.
  2052.     (if (or (null zz) (null gp)) 
  2053.       -1
  2054.       ($SUBSTITUTE ZZ GP))));Substitute Values for gensyms.
  2055.  
  2056. (DEFUN COEFSOLVE (N CL E)    
  2057.   (DO (($BREAKUP)
  2058.        (EQL (NCONS (PDIS (PTERM E N))) (CONS (PDIS (PTERM E M)) EQL))
  2059.        (M (m+ N -1) (m+ M -1)))
  2060.       ((SIGNP L M) (SOLVEX EQL CL NIL NIL))))
  2061.  
  2062. (DEFUN RECTZTO%PI2 (P PE D)
  2063. (PROG (DP N PL A B C denom-exponential)
  2064. (if (not (and (SETQ denom-exponential (CATCH 'PIN%EX (PIN%EX D)))
  2065.           (%e-integer-coeff pe)
  2066.           (%e-integer-coeff d)))
  2067.   (return ()))
  2068. (SETQ N (m* (COND ((NULL P) -1.)
  2069.           (T ($EXPAND (m*t '$%i %pi2 (MAKPOLY P)))))
  2070.         PE))
  2071.   (let ((var 'z*)
  2072.     (leadcoef ()))
  2073.     (SETQ PL (CDR (POLELIST DENOM-EXPONENTIAL #'(LAMBDA (J) 
  2074.                     (OR (NOT (equal ($imagpart J) 0))
  2075.                     (EQ ($askSIGN ($realpart J)) '$neg)))
  2076.                 #'(LAMBDA (J)
  2077.                     (NOT (EQ ($askSIGN ($realpart J)) '$zero)))))))
  2078.   (COND ((NULL PL)  (RETURN NIL))
  2079.     ((OR (CADR PL)
  2080.          (CADDR PL))  (SETQ DP (SDIFF D VAR))))
  2081.   (COND ((CADR PL)  (SETQ B (MAPCAR #'log-imag-0-2%pi (CADR PL)))
  2082.             (SETQ B (RES1 N DP B))
  2083.             (setq b (m+l b)))
  2084.     (t (setq b 0.)))
  2085.    (COND ((CADDR PL)
  2086.       (let ((temp (MAPCAR #'log-imag-0-2%pi (CADDR PL))))
  2087.         (setq c (append temp (mapcar #'(lambda (j) 
  2088.                          (m+ (m*t '$%i %pi2) j))
  2089.                      temp)))
  2090.         (SETQ C (RES1 N DP C))
  2091.         (setq c (m+l c))))
  2092.      (t (setq c 0.)))
  2093.    (COND ((CAR PL)
  2094.       (let ((poles (mapcar #'log-imag-0-2%pi (caar pl)))
  2095.         (exp (m// n (subst (m^t '$%e var) 'z* denom-exponential))))
  2096.            (SETQ A (mapcar #'(lambda (j) 
  2097.                    ($RESIDUE exp var j))
  2098.                    poles))
  2099.            (setq a (m+l a))))
  2100.      (t (setq a 0.)))
  2101.    (RETURN (sRATSIMP (m+ a b (m* '((rat) 1. 2.) c))))))
  2102.  
  2103. (DEFUN GENPOLY (I)
  2104.   (do ((i i (m+ i -1))
  2105.        (c (gensym) (gensym))
  2106.        (cl ())
  2107.        (ans ()))
  2108.       ((zerop i)
  2109.        (cons (m+l ans) cl))
  2110.     (SETQ ANS (CONS (m* C (m^t VAR I)) ANS))
  2111.     (SETQ CL (CONS C CL))))
  2112.  
  2113. (DECLARE-TOP(SPECIAL *FAILFLAG *LHFLAG LHV *INDICATOR CNT *DISCONFLAG)) 
  2114.  
  2115. (defun %e-integer-coeff (exp)
  2116.   (cond ((mapatom exp) t)
  2117.     ((and (mexptp exp)
  2118.           (eq (cadr exp) '$%e)
  2119.           (eq (ask-integer ($coeff (caddr exp) var) '$integer)
  2120.           '$yes))  t)
  2121.     (t (andmapc '%e-integer-coeff (cdr exp)))))
  2122.  
  2123. (DEFUN WLINEARPOLY (E VAR)
  2124.        (COND ((AND (SETQ E (POLYINX E VAR T))
  2125.            (EQUAL (DEG E) 1.))
  2126.           (SUBIN 1. E)))) 
  2127.  
  2128. (DECLARE-TOP(SPECIAL E $EXPONENTIALIZE))
  2129.  
  2130. (DEFUN PIN%EX (EXP)
  2131.   (PIN%EX0 (COND ((NOTINVOLVE EXP '(%SINH %COSH %TANH)) EXP)
  2132.          (T (SETQ EXP (let (($EXPONENTIALIZE t))
  2133.                 ($expand EXP)))))))
  2134.  
  2135. (DEFUN PIN%EX0 (E)
  2136.   (COND ((NOT (AMONG VAR E))  E)
  2137.     ((ATOM E)  (THROW 'PIN%EX NIL))
  2138.     ((AND (mexptp E)
  2139.           (EQ (CADR E)  '$%E))
  2140.      (COND ((EQ (CADDR E) VAR)  'Z*)
  2141.            ((let ((linterm (wlinearpoly (caddr e) var)))
  2142.           (and linterm
  2143.                (m* (subin 0 e) (m^t 'z* linterm)))))
  2144.            (T (THROW 'PIN%EX NIL))))
  2145.     ((mtimesp E)  (m*l (MAPCAR #'PIN%EX0 (CDR E))))
  2146.     ((mplusp E)  (m+l (MAPCAR #'PIN%EX0 (CDR E))))
  2147.     (T (THROW 'PIN%EX NIL))))
  2148.  
  2149. (DECLARE-TOP (UNSPECIAL E)) 
  2150.  
  2151. (DEFUN P*PIN%EX (ND*)
  2152.        (SETQ ND* ($FACTOR ND*))
  2153.        (COND ((POLYINX ND* VAR NIL) (SETQ P* (CONS ND* P*)) T)
  2154.          ((CATCH 'PIN%EX (PIN%EX ND*)) (SETQ PE* (CONS ND* PE*)) T)
  2155.          ((mtimesp ND*)
  2156.           (andMAPCAR #'P*PIN%EX (CDR ND*)))))
  2157.  
  2158. (DEFUN FINDSUB (P)
  2159.        (COND ((FINDP P) NIL)
  2160.          ((SETQ ND* (BX**N P)) 
  2161.           (m^t VAR (CAR ND*)))
  2162.          ((SETQ P (BX**N+A P))
  2163.           (m* (CADDR P) (m^t VAR (CADR P)))))) 
  2164.  
  2165. (DEFUN FUNCLOGOR%E (E)
  2166.   (PROG (ANS ARG NVAR R) 
  2167.     (COND ((OR (RATP E VAR)
  2168.            (INVOLVE E '(%SIN %COS %TAN))
  2169.            (NOT (SETQ ARG (XOR (AND (SETQ ARG (INVOLVE E '(%LOG)))
  2170.                     (SETQ R '%LOG))
  2171.                    (%EINVOLVE E)))))
  2172.        (RETURN NIL)))
  2173. AG (SETQ NVAR (COND ((EQ R '%LOG) `((%LOG) ,ARG))
  2174.             (T (m^t '$%E ARG))))
  2175.    (SETQ ANS (MAXIMA-SUBSTITUTE (m^t 'YX -1.) (m^t NVAR -1.) (MAXIMA-SUBSTITUTE 'YX nvar E)))
  2176.    (COND ((NOT (AMONG VAR ANS))  (RETURN (LIST (SUBST VAR 'YX ANS) NVAR)))
  2177.      ((AND (NULL R) 
  2178.            (SETQ ARG (FINDSUB ARG)))
  2179.       (GO AG)))))
  2180.  
  2181. (DEFUN DINTBYPART (U V A B)
  2182. ;;;SINCE ONLY CALLED FROM DINTLOG TO get RID OF LOGS - IF LOG REMAINS, QUIT
  2183.        (let ((ad (antideriv v)))
  2184.         (COND ((or (NULL AD)
  2185.                (INVOLVE AD '(%LOG)))
  2186.            NIL)
  2187.           (t (let ((P1 (m* U AD))
  2188.                (P2 (m* AD (SDIFF U VAR))))
  2189.               (let ((P1-part1 (get-LIMit P1 VAR B '$MINUS))
  2190.                 (p1-part2 (get-LIMIT P1 VAR A '$PLUS)))
  2191.                    (cond ((or (null p1-part1)
  2192.                       (null p1-part2))
  2193.                       nil)
  2194.                      (t (let ((P2 (Let ((*DEF2* t))
  2195.                                (DEFINT P2 VAR A B))))
  2196.                          (COND (P2 (add* p1-part1 
  2197.                                  (m- p1-part2)
  2198.                                  (m- P2)))
  2199.                            (t nil)))))))))))
  2200.  
  2201. (DEFUN DINTEXP (EXP arg &aux ans)
  2202.   (let ((dintexp-recur t))    ;recursion stopper
  2203.     (COND ((and (sinintp exp var) ;To be moved higher in the code.
  2204.         (setq ans (antideriv exp))
  2205.         (setq ans (intsubs ans ll ul))))
  2206.       ((setq ANS (FUNCLOGOR%E EXP))
  2207.        (COND ((AND (EQUAL LL 0.) (EQ UL '$INF))
  2208.           (SETQ EXP (SUBIN (m+t 1. arg) (CAR ANS)))
  2209.           (SETQ ANS (m+t -1. (CADR ANS))))
  2210.          (T (SETQ EXP (CAR ANS))
  2211.             (SETQ ANS (CADR ANS))))
  2212.          (INTCV ANS T NIL)))))
  2213.  
  2214. (DEFUN DINTLOG (EXP ARG)
  2215. (let ((dintlog-recur (f1+ dintlog-recur))) ;recursion stopper
  2216.  (PROG (ANS D) 
  2217.        (COND ((AND (EQ UL '$INF)
  2218.            (EQUAL LL 0.)
  2219.            (EQ ARG VAR)
  2220.            (EQUAL 1. ($RATSIMP (M// EXP (m* (m- (SUBIN (m^t VAR -1.)
  2221.                                    EXP))
  2222.                             (m^t VAR -2.))))))
  2223.           (RETURN 0.))
  2224.          ((setq ans (antideriv exp))
  2225.           (return (intsubs ans ll ul)))
  2226.          ((SETQ ANS (LOGX1 EXP LL UL))
  2227.           (RETURN ANS)))
  2228.        (SETQ ANS (M// EXP `((%LOG) ,ARG)))
  2229.        (COND ((INVOLVE ANS '(%LOG)) (RETURN NIL))
  2230.          ((AND (EQ ARG VAR)
  2231.            (EQUAL 0. (NO-ERR-SUB 0. ANS))
  2232.            (SETQ D (let ((*DEF2* t))
  2233.                 (DEFINT (m* ANS (m^t VAR '*Z*))
  2234.                     VAR LL UL))))
  2235.           (RETURN (DERIVAT '*Z* 1. D 0.)))
  2236.          ((SETQ ANS (DINTBYPART `((%LOG) ,ARG) ANS LL UL))
  2237.           (RETURN ANS))))))
  2238.  
  2239. (DEFUN DERIVAT (VAR N E PT) (SUBIN PT (APPLY '$DIFF (LIST E VAR N)))) 
  2240.  
  2241. ;;;MAYBPC RETURNS (COEF EXPO CONST)
  2242. (DEFUN MAYBPC (E VAR)
  2243.        (COND (MTOINF* (THROW 'GGRM (LINPOWER0 E VAR)))
  2244.          ((AND (NOT MTOINF*)
  2245.            (NULL (SETQ E (BX**N+A E))))    ;bx**n+a --> (a b n) or nil.
  2246.           NIL)                ;with var being x.
  2247.          ((AND (AMONG '$%I (CADDR E))
  2248.            (ZEROP1 ($realpart (CADDR E)))
  2249.            (SETQ ZN ($imagPART (CADDR E)))
  2250.            (EQ ($askSIGN (CADR E)) '$pos))
  2251.           (COND ((EQ ($askSIGN ZN) '$neg)
  2252.              (SETQ VAR -1.)
  2253.              (SETQ ZN (m- ZN)))
  2254.             (T (SETQ VAR 1.)))
  2255.           (SETQ ZD (m^t '$%E (M// (mul* VAR '$%I '$%PI (m+t 1. ND*))
  2256.                       (m*t 2. (CADR E)))))
  2257.           `(,ZN ,(CADR E) ,(CAR E)))
  2258.          ((AND (OR (EQ (SETQ VAR ($askSIGN ($realPART (CADDR E)))) '$neg)
  2259.                (EQUAL VAR '$zero))
  2260.            (equal ($imagpart (CADR E)) 0)
  2261.            (RATGREATERP (CADR E) 0.))
  2262.           `(,(m- (CADDR E)) ,(CADR E) ,(CAR E)))))
  2263.  
  2264. (DEFUN GGR (E IND)
  2265.  (PROG (C ZD ZN NN* DN* ND* DOSIMP $%EMODE) 
  2266.    (SETQ ND* 0.)
  2267.    (COND (IND (SETQ E ($EXPAND E))
  2268.           (COND ((AND (mplusp E)
  2269.               (let ((*NODIVERG t))
  2270.                    (SETQ E (CATCH 'Divergent
  2271.                            (ANDMAPCAR
  2272.                         #'(LAMBDA (J) 
  2273.                             (GGR J NIL))
  2274.                         (CDR E))))))
  2275.              (COND ((EQ E 'Divergent) NIL)
  2276.                (T (RETURN (sRATSIMP (CONS '(MPLUS) E)))))))))
  2277.    (SETQ E (RMCONST1 E))
  2278.    (SETQ C (CAR E))
  2279.    (SETQ E (CDR E))
  2280.    (COND ((SETQ E (GGR1 E VAR))
  2281.       (SETQ E (APPLY 'GAMMA1 E))
  2282.       (COND (ZD (SETQ $%EMODE T)
  2283.             (SETQ DOSIMP T)
  2284.             (SETQ E (m* ZD E))))))
  2285.    (COND (E (RETURN (m* C E))))))
  2286.  
  2287.  
  2288.  
  2289. (DEFUN GGR1 (E VAR) 
  2290.        (COND ((ATOM E) NIL)
  2291.          ((AND (mexptp E)
  2292.            (EQ (CADR E) '$%E))
  2293.           (COND ((SETQ E (MAYBPC (CADDR E) VAR)) (CONS 0. E))))
  2294.          ((AND (mtimesp E)
  2295.            (NULL (CDDDR E))
  2296.            (OR (AND (SETQ DN* (XTORTERM (CADR E) VAR))
  2297.                 (RATGREATERP (SETQ ND* ($realPART DN*))
  2298.                      -1.)
  2299.                 (SETQ NN* (GGR1 (CADDR E) VAR)))
  2300.                (AND (SETQ DN* (XTORTERM (CADDR E) VAR))
  2301.                 (RATGREATERP (SETQ ND* ($realPART DN*))
  2302.                      -1.)
  2303.                 (SETQ NN* (GGR1 (CADR E) VAR)))))
  2304.           (RPLACA NN* DN*))))
  2305.  
  2306.  
  2307. (DEFUN BX**N+A (E)
  2308. ;;; returns list of (a b n) or nil.
  2309.        (COND ((EQ E VAR) 
  2310.           (LIST 0. 1. 1.))
  2311.          ((OR (ATOM E) 
  2312.           (MNUMP E)) ())
  2313.          (t (let ((a (NO-ERR-SUB 0. E)))
  2314.              (cond ((null a)  ())
  2315.                (t (SETQ E (m+ E (m*t -1. A)))
  2316.                   (COND ((SETQ E (BX**N E))
  2317.                      (CONS A E))
  2318.                     (t ()))))))))
  2319.  
  2320. (DEFUN BX**N (E)
  2321. ;;;returns a list (n e) or nil.
  2322.        (let ((N ()))
  2323.         (AND (SETQ N (XEXPONGET E VAR))
  2324.          (NOT (AMONG VAR
  2325.                  (SETQ E (let (($maxposex 1)
  2326.                        ($maxnegex 1))
  2327.                        ($EXPAND (M// E (m^t VAR N)))))))
  2328.          (LIST N E))))
  2329.  
  2330. (DEFUN XEXPONGET (E NN*)
  2331.        (COND ((ATOM E) (COND ((EQ E VAR) 1.)))
  2332.          ((MNUMP E) NIL)
  2333.          ((AND (mexptp E)
  2334.            (EQ (CADR E) NN*)
  2335.            (NOT (AMONG NN* (CADDR E))))
  2336.           (CADDR E))
  2337.          (T (ORMAPC #'(LAMBDA (J)
  2338.                 (XEXPONGET J NN*))
  2339.             (CDR E)))))
  2340.  
  2341.  
  2342. ;;; given (b*x^n+a)^m returns (m a n b)
  2343. (DEFUN BXM (E IND)
  2344. (let (m r)
  2345.   (COND ((OR (ATOM E)
  2346.          (MNUMP E)
  2347.          (INVOLVE E '(%LOG %SIN %COS %TAN))
  2348.          (%EINVOLVE E))  NIL)
  2349.     ((mtimesp E)  NIL)
  2350.     ((mexptp E)  (COND ((AMONG VAR (CADDR E))  NIL)
  2351.                ((SETQ R (BX**N+A (CADR E))) 
  2352.                 (CONS (caddr e) R))))
  2353.     ((SETQ R (BX**N+A E))  (CONS 1. R))
  2354.     ((not (null IND))
  2355. ;;;Catches Unfactored forms.
  2356.      (SETQ M (M// (SDIFF E VAR) E))
  2357.      (NUMDEN M)
  2358.      (SETQ M NN*)
  2359.      (setq r dn*)
  2360.      (COND 
  2361.       ((AND (SETQ R (BX**N+A ($ratsimp r)))
  2362.         (NOT (AMONG VAR (SETQ M (M// M (m* (CADR R) (CADDR R)
  2363.                            (m^t VAR (m+t -1.
  2364.                                 (CADR R))))))))
  2365.         (SETQ E (M// (SUBIN 0. E) (m^t (CAR R) M))))
  2366.        (COND ((equal E 1.)
  2367.           (CONS M R))
  2368.          (T (SETQ E (m^ E (M// 1. M)))
  2369.             (LIST M (m* E (CAR R)) (CADR R) 
  2370.               (m* E (CADDR R))))))))
  2371.     (t ()))))
  2372.  
  2373. ;;;Is E = VAR raised to some power? If so return power or 0.
  2374. (DEFUN FINDP (E) 
  2375.        (COND ((NOT (AMONG VAR E)) 0.)
  2376.          (T (XTORTERM E VAR))))
  2377.  
  2378. (DEFUN XTORTERM (E VAR1)
  2379. ;;;Is E = VAR1 raised to some power? If so return power.
  2380.        (COND ((ALIKE1 E VAR1) 1.)
  2381.          ((ATOM E) NIL)
  2382.          ((AND (mexptp E)
  2383.            (ALIKE1 (CADR E) VAR1)
  2384.            (NOT (AMONG VAR (CADDR E))))
  2385.           (CADDR E)))) 
  2386.  
  2387. (DEFUN TBF (L)
  2388.        (m^ (m* (m^ (CADDR L) '((RAT) 1. 2.))
  2389.            (m+ (CADR L) (m^ (m* (CAR L) (CADDR L))
  2390.                 '((RAT) 1. 2.))))
  2391.        -1.))
  2392.  
  2393. (defun radbyterm (d l)
  2394.   (do ((l l (cdr l))
  2395.        (ans ()))
  2396.       ((null l)
  2397.        (m+l ans))
  2398.     (let (((const . integrand) (rmconst1 (car l))))
  2399.      (setq ans (cons (m* const (dintrad0 integrand d))
  2400.              ans)))))
  2401.  
  2402. (DEFUN SQDTC (E IND)
  2403.        (PROG (A B C VARLIST) 
  2404.          (SETQ VARLIST (LIST VAR))
  2405.          (NEWVAR E)
  2406.          (SETQ E (CDADR (RATREP* E)))
  2407.          (SETQ C (PDIS (PTERM E 0.)))
  2408.          (SETQ B (m*t (m//t 1. 2.) (PDIS (PTERM E 1.))))
  2409.          (SETQ A (PDIS (PTERM E 2.)))
  2410.          (COND ((AND (EQ ($askSIGN (m+ B (m^ (m* A C)
  2411.                          '((RAT) 1. 2.))))
  2412.                  '$pos)
  2413.              (OR (AND IND
  2414.                   (NOT (EQ ($askSIGN A) '$neg))
  2415.                   (EQ ($askSIGN C) '$pos))
  2416.                  (AND (EQ ($askSIGN A) '$pos)
  2417.                   (NOT (EQ ($askSIGN C) '$neg)))))
  2418.             (RETURN (LIST A B C)))))) 
  2419.  
  2420. (DEFUN DIFAP1 (E PWR VAR M PT)
  2421.        (M// (mul* (COND ((eq (ask-integer M '$even) '$yes)
  2422.               1.)
  2423.              (T -1.))
  2424.            `((%GAMMA) ,PWR)
  2425.            (DERIVAT VAR M E PT))
  2426.          `((%GAMMA) ,(m+ PWR M)))) 
  2427.  
  2428. (DEFUN SQRTINVOLVE (E)
  2429.        (COND ((ATOM E) NIL)
  2430.          ((MNUMP E) NIL)
  2431.          ((AND (mexptp E) 
  2432.            (AND (MNUMP (CADDR E))
  2433.             (NOT (NUMBERP (CADDR E)))
  2434.             (EQUAL (CADDR (CADDR E)) 2.))
  2435.            (AMONG VAR (CADR E)))
  2436.           (CADR E))
  2437.          (T (ORMAPC #'SQRTINVOLVE (CDR E))))) 
  2438.  
  2439. (DEFUN BYDIF (R S D)
  2440.        (let ((B 1)  p)
  2441.         (SETQ D (m+ (m*t '*Z* VAR) D))
  2442.         (COND ((OR (ZEROP1 (SETQ P (m+ S (m*t -1. R))))
  2443.                (AND (ZEROP1 (m+ 1. P))
  2444.                 (SETQ B VAR)))
  2445.            (DIFAP1 (DINTRAD0 B (m^ D '((RAT) 3. 2.)))
  2446.                     '((RAT) 3. 2.) '*Z* R 0.))
  2447.            ((EQ ($askSIGN P) '$pos)
  2448.             (DIFAP1 (DIFAP1 (DINTRAD0 1. (m^ (m+t 'Z** D)
  2449.                              '((RAT) 3. 2.)))
  2450.                     '((RAT) 3. 2.) '*Z* R 0.)
  2451.                 '((RAT) 3. 2.) 'Z** P 0.)))))
  2452.  
  2453. (DEFUN DINTRAD0 (N D)
  2454.        (let (L R S) 
  2455.      (COND ((AND (mexptp D) 
  2456.              (EQUAL (DEG (CADR D)) 2.))
  2457.         (COND ((ALIKE1 (CADDR D) '((RAT) 3. 2.))
  2458.                (COND ((AND (EQUAL N 1.)
  2459.                    (SETQ L (SQDTC (CADR D) T)))
  2460.                   (TBF L))
  2461.                  ((AND (EQ N VAR)
  2462.                    (SETQ L (SQDTC (CADR D) NIL)))
  2463.                   (TBF (REVERSE L)))))
  2464.               ((AND (SETQ R (FINDP N))
  2465.                 (OR (EQ ($askSIGN (m+ -1. (m-  R) (m*t 2.
  2466.                                    (CADDR D))))
  2467.                     '$pos)
  2468.                 (Diverg))
  2469.                 (SETQ S (m+ '((RAT) -3. 2.) (CADDR D)))
  2470.                 (EQ ($askSIGN S) '$pos)
  2471.                 (eq (ask-integer S '$integer) '$yes))
  2472.                (BYDIF R S (CADR D)))
  2473.               ((POLYINX N VAR NIL)
  2474.                (RADBYTERM D (CDR N))))))))
  2475.  
  2476.  
  2477. ;;;Looks at the IMAGINARY part of a log and puts it in the interval 0 2*%pi.
  2478. (defun log-imag-0-2%pi (x)
  2479.   (let ((plog (simplify `((%plog) ,x))))
  2480.     (cond ((not (free plog '%plog))
  2481.        (subst '%log '%plog plog))
  2482.       (t (let (((real . imag) (trisplit plog)))
  2483.            (cond ((eq ($asksign imag) '$neg)
  2484.               (setq imag (m+ imag %pi2)))
  2485.              ((eq ($asksign (m- imag %pi2)) '$pos)
  2486.               (setq imag (m- imag %pi2)))
  2487.              (t t))
  2488.            (m+ real (m* '$%i imag)))))))
  2489.  
  2490.         
  2491. ;;; Temporary fix for a lacking in taylor, which loses with %i in denom.
  2492. ;;; Besides doesn't seem like a bad thing to do in general.
  2493. (defun %i-out-of-denom (exp)
  2494.        (let ((denom ($denom exp))
  2495.          (den-conj nil))
  2496.         (cond ((among '$%i denom)
  2497.            (setq den-conj (MAXIMA-SUBSTITUTE (m- '$%i) '$%i denom))
  2498.            (setq exp (m* den-conj ($ratsimp (m// exp den-conj))))
  2499.            (setq exp (simplify ($multthru  (sratsimp exp)))))
  2500.           (t exp))))
  2501.  
  2502. ;;; LL and UL must be real otherwise this routine return $UNKNOWN.
  2503. ;;; Returns $no $unknown or a list of poles in the interval (ll ul)
  2504. ;;; for exp w.r.t. var.
  2505. ;;; Form of list ((pole . multiplicity) (pole1 . multiplicity) ....)
  2506. (defun poles-in-interval (exp var ll ul)
  2507.  (let* ((denom (cond ((mplusp exp)
  2508.               ($denom (sratsimp exp)))
  2509.              ((and (mexptp exp)
  2510.                (free (caddr exp) var)
  2511.                (eq ($asksign (caddr exp)) '$neg))
  2512.               (m^ (cadr exp) (m- (caddr exp))))
  2513.              (t ($denom exp))))
  2514.     (roots (real-roots denom var))
  2515.     (ll-pole (limit-pole exp var ll '$plus))
  2516.     (ul-pole (limit-pole exp var ul '$minus)))
  2517.        (cond ((or (eq roots '$failure)
  2518.           (null ll-pole)
  2519.           (null ul-pole))   '$unknown)
  2520.          ((and (eq roots '$no)
  2521.            (eq ll-pole '$no)
  2522.            (eq ul-pole '$no))  '$no)
  2523.          (t (cond ((equal roots '$no)
  2524.                (setq roots ())))
  2525.         (do ((dummy roots (cdr dummy))
  2526.              (pole-list (cond ((not (eq ll-pole '$no))
  2527.                        `((,ll . 1)))
  2528.                       (t nil))))
  2529.             ((null dummy)
  2530.              (cond ((not (eq ul-pole '$no))
  2531.                 (sort-poles (push `(,ul . 1) pole-list)))
  2532.                ((not (null pole-list))
  2533.                 (sort-poles pole-list))
  2534.                (t '$no)))
  2535.             (let* ((soltn (caar dummy))
  2536.                ;; (multiplicity (cdar dummy)) (not used? -- cwh)
  2537.                (root-in-ll-ul (in-interval soltn ll ul)))
  2538.               (cond ((eq root-in-ll-ul '$no) '$no)
  2539.                 ((eq root-in-ll-ul '$yes)
  2540.                  (let ((lim-ans (is-a-pole exp soltn)))
  2541.                       (cond ((null lim-ans)
  2542.                          (return '$unknown))
  2543.                         ((equal lim-ans 0)
  2544.                          '$no)
  2545.                         (t (push (car dummy)
  2546.                              pole-list))))))))))))
  2547.  
  2548.  
  2549. ;;;Returns $YES if there is no pole and $NO if there is one.
  2550. (defun limit-pole (exp var limit direction)
  2551.   (let ((ans (cond ((memq limit '($minf $inf))
  2552.             (cond ((eq (special-convergent-formp exp limit) '$yes)
  2553.                '$no)
  2554.               (t (get-limit (m* exp var) var limit direction))))
  2555.            (t '$no))))
  2556.        (cond ((eq ans '$no)   '$no)
  2557.           ((null ans)   nil)
  2558.           ((equal ans 0.)   '$no)
  2559.           (t '$yes))))
  2560.  
  2561. ;;;Takes care of forms that the ratio test fails on.
  2562. (defun special-convergent-formp (exp limit)
  2563.   (cond ((not (oscip exp))  '$no)
  2564.     ((or (eq (sc-converg-form exp limit) '$yes)
  2565.          (eq (exp-converg-form exp limit) '$yes))
  2566.      '$yes)
  2567.     (t  '$no)))
  2568.  
  2569. (defun exp-converg-form (exp limit)
  2570.   (let (exparg)
  2571.     (setq exparg (%einvolve exp))
  2572.     (cond ((or (null exparg)
  2573.            (freeof '$%i exparg))
  2574.        '$no)
  2575.       (t (cond
  2576.           ((and (freeof '$%i 
  2577.                 (%einvolve 
  2578.                  (setq exp 
  2579.                    (sratsimp (m// exp (m^t '$%e exparg))))))
  2580.             (equal (get-limit exp var limit)  0))
  2581.             '$yes)
  2582.           (t '$no))))))
  2583.  
  2584. (defun sc-converg-form (exp limit)       
  2585.  (prog (scarg trigpow)
  2586.   (setq exp ($expand exp))
  2587.   (setq scarg (involve (sin-sq-cos-sq-sub exp) '(%sin %cos)))
  2588.   (cond ((null scarg) (return '$no))
  2589.     ((and (polyinx scarg var ())
  2590.           (eq ($asksign (m- ($hipow scarg var) 1)) '$pos)) (return '$yes))
  2591.     ((not (freeof var (sdiff scarg var)))
  2592.      (return '$no))
  2593.     ((and (setq trigpow ($hipow exp `((%sin) ,scarg)))
  2594.           (eq (ask-integer trigpow '$odd) '$yes)
  2595.           (equal (get-limit (m// exp `((%sin) ,scarg)) var limit)
  2596.              0))
  2597.      (return '$yes))
  2598.     ((and (setq trigpow ($hipow exp `((%cos) ,scarg)))
  2599.           (eq (ask-integer trigpow '$odd) '$yes)
  2600.           (equal (get-limit (m// exp `((%cos) ,scarg)) var limit)
  2601.              0))
  2602.      (return '$yes))
  2603.     (t (return '$no)))))
  2604.  
  2605. (defun is-a-pole (exp soltn)
  2606.        (get-limit ($radcan 
  2607.            (m* (MAXIMA-SUBSTITUTE (m+ 'epsilon soltn) var exp)
  2608.              'epsilon))
  2609.           'epsilon 0 '$plus))
  2610.  
  2611. (defun in-interval (place ll ul)
  2612.     ;; real values for ll and ul; place can be imaginary.
  2613.   (let ((order (ask-greateq ul ll)))
  2614.     (cond ((eq order '$yes))
  2615.       ((eq order '$no) (let ((temp ul)) (setq ul ll ll temp)))
  2616.       (t (merror "Incorrect limits given to DEFINT:~%~M"
  2617.              (list '(mlist simp) ll ul)))))
  2618.   (if (not (equal ($imagpart place) 0))
  2619.       '$no
  2620.       (let ((lesseq-ul (ask-greateq ul place))
  2621.         (greateq-ll (ask-greateq place ll)))
  2622.     (if (and (eq lesseq-ul '$yes) (eq greateq-ll '$yes)) '$yes '$no))))
  2623.  
  2624. (defun real-roots (exp var)
  2625.   (let (($solvetrigwarn (cond (defintdebug t) ;Rest of the code for
  2626.                   (t ())))           ;TRIGS in denom needed.
  2627.     ($solveradcan (cond ((or (among '$%i exp)
  2628.                  (among '$%e exp)) t)
  2629.                 (t nil)))
  2630.     *roots *failures) ;special vars for solve.
  2631.        (cond ((not (among var exp))   '$no)
  2632.          (t (solve exp var 1)
  2633.         (cond (*failures '$failure)
  2634.               (t (do ((dummy *roots (cddr dummy))
  2635.                   (rootlist))
  2636.                  ((null dummy)
  2637.                   (cond ((not (null rootlist))
  2638.                      rootlist)
  2639.                     (t '$no)))
  2640.                (cond ((equal ($imagpart (caddar dummy)) 0)
  2641.                   (setq rootlist 
  2642.                     (cons (cons 
  2643.                            ($rectform (caddar dummy))
  2644.                            (cadr dummy))
  2645.                           rootlist)))))))))))
  2646.  
  2647. (defun ask-greateq (x y)
  2648. ;;; Is x > y. X or Y can be $MINF or $INF, zeroA or zeroB.
  2649.        (let ((x (cond ((among 'zeroa x)
  2650.                (subst 0 'zeroa x))
  2651.               ((among 'zerob x)
  2652.                (subst 0 'zerob x))
  2653.               ((among 'epsilon x)
  2654.                (subst 0 'epsilon x))
  2655.               ((or (among '$inf x)
  2656.                (among '$minf x))
  2657.                ($limit x))
  2658.               (t x)))
  2659.          (y (cond ((among 'zeroa y)
  2660.                (subst 0 'zeroa y))
  2661.               ((among 'zerob y)
  2662.                (subst 0 'zerob y))
  2663.               ((among 'epsilon y)
  2664.                (subst 0 'epsilon y))
  2665.               ((or (among '$inf y)
  2666.                (among '$minf y))
  2667.                ($limit y))
  2668.               (t y))))
  2669.         (cond ((eq x '$inf)
  2670.            '$yes)
  2671.           ((eq x '$minf)
  2672.            '$no)
  2673.           ((eq y '$inf)
  2674.            '$no)
  2675.           ((eq y '$minf)
  2676.            '$yes)
  2677.           (t (let ((ans ($asksign (m+ x (m- y)))))
  2678.               (cond ((memq ans '($zero $pos))
  2679.                  '$yes)
  2680.                 ((eq ans '$neg)
  2681.                  '$no)
  2682.                 (t '$unknown)))))))
  2683.  
  2684. (defun sort-poles (pole-list)
  2685.        (sort pole-list #'(lambda (x y)
  2686.                 (cond ((eq (ask-greateq (car x) (car y))
  2687.                        '$yes)
  2688.                        nil)
  2689.                       (t t)))))
  2690.