;;; -*-  Mode: Lisp; Package: Maxima; Syntax: Common-Lisp; Base: 10 -*- ;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;     The data in this file contains enhancments.                    ;;;;;
;;;                                                                    ;;;;;
;;;  Copyright (c) 1984,1987 by William Schelter,University of Texas   ;;;;;
;;;     All rights reserved                                            ;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;     (c) Copyright 1982 Massachusetts Institute of Technology         ;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

(in-package "MAXIMA")
(macsyma-module series)

(declare-top (genprefix ps)
	 (*lexpr factor $gcd $rat)
	 (fixnum %n %nn)
	 (special var *n *a *m *c *index $cauchysum *gcd*
		  nn* dn* $ratsimpexpons *infsumsimp
		  *ratexp splist *var usexp $verbose ans *trigred
		  *form indl *noexpand $ratexpand))

(load-macsyma-macros rzmac)

;
;******************************************************************************
;				driver 	stage
;******************************************************************************
;
;		the following throw labels are used
;
;psex -- for throws on failure to expand
;

(defmfun $powerseries (l var *pt)
	(cond ((signp e *pt)
	       (cond ((symbolp var) (seriesexpand* l))
		     ((numberp var) (improper-arg-err var '$powerseries))
		     (t (sbstpt l 'x (gensym) var var))))
	      ((eq *pt '$inf)
	       (sbstpt l (m// 1 'x) (gensym) var (div* 1 var)))
	      (t (sbstpt l (m+ 'x *pt) (gensym)
			 var (simplifya (m- var *pt) nil)))))

(defun sbstpt (exp sexp var var1 usexp)
	(setq sexp (subst var 'x sexp))
	(setq exp (MAXIMA-SUBSTITUTE sexp var1 exp))
	(MAXIMA-SUBSTITUTE usexp var (seriesexpand* exp)))

(defun seriesexpand* (x)
  (let (*n *a *m *c
	   (*index (gensumindex))
	   ($cauchysum t)
	   ($ratsimpexpons t) $ratexpand
	   *infsumsimp *ratexp *trigred *noexpand)
       (mkind *index '$integer)
       (setq x (catch 'psex (sp2expand (seriespass1 x))))
       (cond ((and x (atom x)) x)
	     ((and x (not (eq (car x) 'err))) x)
	     ($verbose
	      (mtell "Unable to expand for the following reason:")
	      (cond ((null x) (mtell "~%no reason given") '|&Unable to expand|)
		    (t (cdr x))))
	     (t '|&Unable to expand|))))

(defun out-of (e)
  (let  ((e      (cond ((and (boundp '*var) *var)
			(subst (list '(mexpt) *var *gcd*) var e))
		       (t e)))
	 (var      (cond ((and (boundp '*var) *var)) (t var))))
    (cond ((and (boundp 'usexp) usexp)
	   (subst usexp var e))
	  (t e))))

(defun show-exp (e)
  (mtell "~%~%~M~%~%"
	 (list '(mlable) nil (out-of e))))

(defun seriespass1 (e)
  (let ((w (sratsimp (sp1 e))))
       (when $verbose
	     (terpri)
	     (mtell "In the first simplification we have returned:")
	     (show-exp w))
       w))

;
;*****************************************************************************
;	pass two		expansion phase
;*****************************************************************************
;

(defun sp2expand (exp)
    (cond ((or (free exp var) (atom exp)) exp)
	  ((mbagp exp) (cons (car exp) (mapcar #'sp2expand (cdr exp))))
	  ((sratp exp var) (ratexp exp))
	  ((eq (caar exp) 'mexpt) (sp2expt exp))
	  ((oldget (caar exp) 'sp2) (sp2sub (sp2trig exp) (cadr exp)))
	  ((poly? exp var) exp)
	  ((eq (caar exp) 'mtimes) (m*l (mapcar #'sp2expand (cdr exp))))
	  ((eq (caar exp) 'mplus) (m+l (mapcar #'sp2expand (cdr exp))))
	  ((eq (caar exp) '%log) (sp2log (cadr exp)))
	  ((eq (caar exp) '%derivative) (sp2diff (cadr exp) (cddr exp)))
	  ((eq (caar exp) '%integrate)
	   (sp2integ (cadr exp) (caddr exp) (cdddr exp)))
	  ((memq (caar exp) '(%product %sum))
	   (list* (car exp) (sp2expand (cadr exp)) (cddr exp)))
	  (t (list '(%sum)
		   (m* (m^ var *index)
		       (m^ (list '(mfactorial) *index) -1)
		       (list '(%at) (list '(%derivative) exp var *index)
			     (list '(mequal) var 0)))
		   *index 0 '$inf))))

(defun sp2sub (s exp)
  (cond ((smono exp var) (MAXIMA-SUBSTITUTE exp 'sp2var (simplify s)))
	(t (throw 'psex (list 'err '(mtext)
			       "Tried to MAXIMA-SUBSTITUTE the expansion of  "
			       (out-of exp)
			       " into an expansion")))))

(defun ratexp (exp)
  (let (nn* dn* *gcd*)
       (if $verbose
	   (mtell "trying to do a rational function expansion of~%~M"
		  (list '(mlable) nil exp)))
       (numden exp)
       (sratexpnd nn* dn*)))

(defun sratexpnd (n d)
    (let ((ans (list nil))
	  (splist)
	  (linpat
	   '((mtimes) ((coefftt) (cc not-zero-free var))
		      ((mexpt) ((mplus) ((coeffpt) 
					 (w m1 ((mexpt) (x equal var)
							(m not-zero-free var)))
					 (c freevar))
					((coeffpp) (a freevar)))
			       (n not-zero-free var)))))
	(cond ((and (not (equal n 1)) (smono n var))
	       (m* n (sratexpnd 1 d)))
	      ((free d var)
	       (cond ((poly? n var)
		      (m// n d))
		     ((m1 n linpat)
		      (m* (srbinexpnd (cdr ans)) (div* 1 d)))
		     ((throw 'psex nil))))
	      ((smonop d var)
	       (cond ((mplusp n)
		      (m+l (mapcar #'(lambda (q) (div* q d)) (cdr n))))
		     (t (m// n d))))
	      ((not (equal 1 (setq *gcd* (ggcd (nconc (exlist n) (exlist d))))))
	       (sratsubst *gcd* n d))
	      ((and (equal n 1)
		    (prog2 (setq d (let (($ratfac t))
					(ratdisrep ($rat (factor d) var))))
			   (m1 d linpat)))
	       (m// (srbinexpnd (cdr ans)) (cdr (assq 'cc (cdr ans)))))
	      (t
	       (and *ratexp (throw 'psex nil))
	       (if (not (eq (caar d) 'mtimes)) (ratexand1 n d))
	       (do ((p (cdr d) (cdr p)))
		   ((null p) (ratexand1 n d))
		   (cond ((or (eq (car p) var)
			      (and (mexptp (car p)) (eq (cadaar p) var)))
		 	  (return (m* (sratexpnd n (meval (div* d (car p))))
				      (list '(mexpt) (car p) -1))))))))))

(defun ratexand1 (n d)
    (and $verbose
	 (prog2 (mtell "Trying partial fraction expansion of ")
		(show-exp (list '(mquotient) n d))
		(terpri)))
  (funcall #'(lambda (*ratexp) 
    (m+l (mapcar #'ratexp
		 (funcall #'(lambda (l)
			  (cond ((eq (caar l) 'mplus)
				 (and $verbose
				      (mtell "which is ")
				      (show-exp l))
				 (cdr l))
				(t (throw 'psex
					   '(err (mtext)
						 "Partial fraction expansion failed")))))
		  ($partfrac (div* n d) var)))))
   t))

(defun sratsubst (gcd num den)
    (and $verbose
	 (prog2 (mtell "We are substituting for the occurrences of")
		(show-exp (list '(mexpt) var gcd))
		(mtell "in")
	 	(show-exp (list '(mquotient) num den))
		(terpri)))
    (funcall #'(lambda (var* *var)
	(setq num (MAXIMA-SUBSTITUTE (m^ var* (m^ gcd -1)) *var num)
	      den (MAXIMA-SUBSTITUTE (m^ var* (m^ gcd -1)) *var den))
	(MAXIMA-SUBSTITUTE (m^ *var gcd) var*
		    (funcall #'(lambda (var) (sratexpnd num den)) var*)))
     (gensym) var))

(defun ggcd (l)
	(cond ((null l) 1)
	      ((null (cdr l)) (car l))
	      ((equal 1 (car l)) 1)
	      (t ($gcd (ggcd (cdr l)) (car l)))))

(defun exlist (exp)
	(cond ((null exp) nil)
	      ((atom exp)
	       (and (eq exp var) (ncons 1)))
	      ((and (not (atom (car exp))) (eq (caar exp) 'mplus))
	       (exlist (cdr exp)))
	      ((smono (car exp) var)
	       (cond ((equal *n 0) (exlist (cdr exp)))
		     (t (cons *n (exlist (cdr exp))))))
	      (t (exlist (cdr exp)))))

(defun srbinexpnd (ans)
    (let ((n (cdr (assq 'n ans)))
	  (a (cdr (assq 'a ans)))
	  (m (cdr (assq 'm ans)))
	  (c (cdr (assq 'c ans))))
	(cond ((integerp n) (srintegexpd n a m c))
	      (t (list '(%sum)
		       (m// (m* (m^ (m* c var) (m* m *index))
				(m^ a (m- n *index)))
			    (m* (list '($beta) (m- n (m1- *index)) (m1+ *index))
				(m1+ n)))
		       *index 0 '$inf)))))

(defun psp2form (coeff exp bas)
	(list '(%sum) (m* coeff (m^ var exp)) *index bas '$inf))

(defun srintegexpd (n a m c)
    (and $verbose
	 (prog2 (mtell "Using a special rule for expressions of form ")
		(show-exp '((mexpt) ((mplus) $a ((mtimes) $c ((mexpt) $var $m)))
				    ((mminus) $n)))
		(mtell "Here we have")
		(show-exp (list '(mlist) (list '(mequal) '$n n) (list '(mequal) '$a a)
					 (list '(mequal) '$c c) (list '(mequal) '$m m)))))
    (cond ((= n 1)
	   (psp2form
	    (m* (m^ a (m* -1 (m+ 1 *index)))
		(m^ (m* -1 c) *index))
	    (if (equal m 1) *index (m* *index m))
	    0))
	  ((= 2 n)
	   (psp2form (m* (m+ 1 *index)
			 (m^ c (m* -1 (m+ 2 *index)))
			 (m^ (m* -1 a) *index))
		     (if (equal m 1) *index (m* *index m))
		     0))
	  (t (psp2form (m* (do ((nn (f1- n) (f1- nn))
				(l nil (cons (list '(mplus) *index nn) l)))
			       ((= nn 0)
				(m*l (cons (m// 1 (factorial (f1- n))) l))))
			   (m^ (m* -1 c (m^ a -1)) *index)
			   (m^ a (f- n)))
		       (if (equal m 1) *index (m* *index m))
		       0))))

(defun sratp (a var)
	(cond ((atom a) t)
	      ((memq (caar a) '(mplus mtimes)) (sandmap (cdr a)))
	      ((eq (caar a) 'mexpt)
	       (cond ((free (cadr a) var) (free (caddr a) var))
		     ((smono a var) t)
		     ((and (free (caddr a) var) (sratp (cadr a) var)))))
	      (t (free (cadr a) var))))

(defun sandmap (l) (or (null l) (and (sratp (car l) var) (sandmap (cdr l)))))

(defun sp2trig (exp) (subst *index '*index (oldget (caar exp) 'sp2)))

(defun sp2log (e)
    (funcall #'(lambda (exp *a *n)
	(cond ((or (atom e) (free e var)) (list '(%log) e))
	      ((null (smono exp var)) (throw 'psex nil))
	      ((or (and (numberp *a)
			(minusp *a)
			(setq *a (minus *a)))
		   (and (mtimesp *a)
			(numberp (cadr *a))
			(minusp (cadr *a))
			(setq *a (simptimes
				   (list* (car *a) (minus (cadr *a)) (cddr *a))
				   1 t))))
	       (list '(%sum)
		     (m* -1
			 (m^ (m* (m^ var *n) *a) *index)
			 (m^ *index -1))
		     *index 1 '$inf))
	      (t (list '(%sum)
		       (m* -1
			   (m^ (m* -1 *a (m^ var *n)) *index)
			   (m^ *index -1))
		       *index 1 '$inf))))
     (m- e 1) nil nil))

(defun sp2expt (exp)
       (cond ((and (numberp (caddr exp)) (mexptp (cadr exp)))
	      (sp2expt (m^ (cadr (cadr exp))
			   (m* (caddr exp) (caddr (cadr exp))))))
	     ((and (free (caddr exp) var)
		   (signp g (caddr exp))
		   (lessp (caddr exp) $maxposex))
	      (m*l (dup (sp2expand (cadr exp)) (caddr exp))))
	     ((free (cadr exp) var)
	      (sp2sub (subst *index
			     '*index
			     (subst (cond ((eq (cadr exp) '$%e) 'sp2var)
					  (t (list '(mtimes)
						   (list '(mlog) (cadr exp))
						   'sp2var)))
				    'sp2var
				    (get 'mexpt 'sp2)))
		      (caddr exp)))
	     (t (throw 'psex nil))))

(defun dup (x %n) (if (= %n 1) (ncons x) (cons x (dup x (f1- %n)))))

(defun sp2diff (exp l)
  (let (indl)
    (cond ((free exp var)
	   (sp3form (sp2expand exp) (cons '(%derivative) l)))
	  (t (do ((l l (cddr l)) (ll))
		 ((null l)
		  (if ll (sp3form exp (cons '(%derivative) (nreverse ll)))
			 exp))
		 (cond ((eq (car l) var)
			(do ((%i (cadr l) (f1- %i)))
			    ((= %i 0) exp)
			    (setq indl nil
				  exp (sp2diff1 (sp2expand exp) nil nil))))
		       (t (setq ll (list* (cadr l) (car l) ll)))))))))

(defun sp2diff1 (exp ind lol)	;ind is a list of the indices
				;lol is a list of the lower limits
    (cond ((atom exp) (sdiff exp var))
	  ((eq (caar exp) 'mplus)
	   (cons '(mplus)
		 (mapcar #'(lambda (q) (sp2diff1 q ind lol))
			 (cdr exp))))
	  ((eq (caar exp) '%sum)
	   (setq indl (cons (copy-list (cddr exp)) indl))
	   (sp2diff1 (cadr exp)
		     (cons (caddr exp) ind)
		     (cons (cadddr exp) lol)))
	  (t (sp2diff2 exp ind lol))))

(defun sp2diff2 (exp ind lol)
    (let (e fr)
	(setq e (m2 exp '((mtimes) ((coefftt) (fr freevar))
				   ((coefftt) (e true)))
		    nil)
	      fr (cdr (assq 'fr e))
	      e  (cdr (assq 'e e)))
	(sp3reconst
	 (cond ((and (mexptp e) (eq (cadr e) var))
		(cond ((equal 0 (mbinding (ind lol)
					  (meval (m* fr (caddr e)))))
		       (m* (sp3substp1 ind ind (m* fr (caddr e))) e))
		      ((mgrp 1 (mbinding (ind lol)
					 (simplify (mevalatoms (caddr e)))))
		       (m* fr (caddr e) (m^ (cadr e) (m- (caddr e) 1))))
		      (t (sdiff exp var))))
	       (t (sdiff exp var))))))

(defun sp2integ (exp v l)
    (if (null l)
	(if (eq var v)
	    (sp2integ1 ($expand (sp2expand exp)))
	    (sp3form (sp2expand exp) (list '(%integrate) v)))
	(sp2integ2 exp v (car l) (cadr l))))

(defun sp2integ1 (exp)
  (let (pair)
    (cond ((ratp exp var) (ratint exp var))
	  ((eq (caar exp) 'mplus)
	   (cons '(mplus) (mapcar #'sp2integ1 (cdr exp))))
	  ((and (eq (caar exp) 'mtimes)
		(prog2 (setq pair (partition exp var 1))
		       (not (equal (car pair) 1))))
	   (mul2* (car pair) (sp2integ1 (cdr pair))))
	  ((and (eq (caar exp) 'mtimes)
		(prog2 (setq exp ($intosum exp)) nil)))
	  ((or (not (eq (caar exp) '%sum)) (not (isinop (cadr exp) '%sum)))
	   (sinint exp var))
	  (t (let ((indl (ncons (cddr exp))))
	       (sp2integ12 (cadr exp) (ncons (caddr exp)) (ncons (cadddr exp))))))))

(defun sp2integ12 (exp ind lol)
    (cond ((atom exp)
	   (sp3reconst (ratint exp var)))
	  ((eq (caar exp) 'mplus)
	   (sp3reconst
		(m+l (mapcar #'(lambda (q) (sp2integ13 q ind lol))
			     (cdr exp)))))
	  ((eq (caar exp) '%sum)
	   (setq indl (cons (cddr exp) indl))
	   (sp2integ12 (cadr exp)
		       (cons (caddr exp) ind)
		       (cons (cadddr exp) lol)))
	  (t (sp3reconst (sp2integ13 exp ind lol)))))

(defun sp2integ13 (exp ind lol)
    (let (e fr)
	(setq e (m2 exp '((mtimes) ((coefftt) (fr freevar))
				   ((coefftt) (e true)))
		    nil)
	      fr (cdr (assq 'fr e))
	      e  (cdr (assq 'e e)))
	(cond ((and (mexptp e) (eq (cadr e) var))
	       (cond ((mgrp -1 (mbinding (ind lol)
					 (meval (caddr e))))
		      (m* (sp3substpn ind ind (m* fr (caddr e)) -1) e))
		     (t (sinint exp var))))
	      (t (sinint exp var)))))

(defun sp2integ2 (exp v lo hi)
    (if (eq v var) (setq v (gensym) exp (subst v var exp)))
    (cond ((and (free lo var) (free hi var))
	   (cond ((free exp var)
		  (list '(%integrate) (subst var v exp) var lo hi))
		 (t (sp3form (sp2expand exp)
			     (list '(%integrate) v lo hi)))))
	  (t (m+ (sp2sub (setq exp (sp2expand (subst var v exp))) hi)
		 (m* -1 (sp2sub exp lo))))))

;
;************************************************************************************
;	phase three		miscellaneous garbage and final simplification
;************************************************************************************
;

(defun sp3reconst (e)
  (do ((l indl (cdr l)) (e e (list* '(%sum) e (car l))))
      ((null l) e)))

(defun sp3substpn (vars vals exp n)
  (sp3subst vars (mapcar #'(lambda (q) (add2* q n)) vals) exp))

(defun sp3substp1 (vars vals exp) (sp3substpn vars vals exp 1))

(defun sp3subst (vars vals exp)
    (simplify (sublis (mapcar #'cons (cdr vars) (cdr vals))
		      (subst (car vals) (car vars) exp))))

(defun sp3form (e *form) (sp3form1 e))

(defun sp3form1 (e)
    (cond ((atom e) (list* (car *form) e (cdr *form)))
	  ((eq (caar e) 'mplus)
	   (cons '(mplus) (mapcar #'sp3form1 (cdr e))))
	  ((eq (caar e) '%sum)
	   (list* '(%sum) (sp3form1 (cadr e)) (cddr e)))
	  (t (list* (car *form) e (cdr *form)))))

; These are the series expansions for circular functions

(defprop %sin
	 ((%sum) ((mtimes)
		  ((mexpt) -1 *index)
		  ((mexpt) ((mfactorial) ((mplus) ((mtimes) 2 *index) 1)) -1)
		  ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) 1)))
		 *index 0 $inf)
	 sp2)

(defprop %cos
	 ((%sum) ((mtimes) ((mexpt) -1 *index)
			   ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
			   ((mexpt) sp2var ((mtimes) 2 *index)))
		 *index 0 $inf)
	 sp2)

(defprop %tan
	 ((%sum) ((mtimes) ((mexpt) -1 ((mplus) *index -1))
			   ((mexpt) 2 ((mtimes) 2 *index))
			   ((mplus) ((mexpt) 2 ((mtimes) 2 *index)) -1)
			   ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
			   (($bern) ((mtimes) 2 *index))
			   ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) -1)))
		 *index 0 $inf)
	 sp2)

(defprop %csc
	 ((%sum) ((mtimes) 2 
		  ((mexpt) -1 ((mplus) *index -1))
		  ((mplus) ((mexpt) 2 ((mplus) ((mtimes) 2 *index) -1)) -1)
		  ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
		  (($bern) ((mtimes) 2 *index))
		  ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) -1)))
		 *index 0 $inf)
	 sp2)

(defprop %cot
	 ((%sum) ((mtimes)
		  ((mexpt) -1 *index)
		  ((mexpt) 2 ((mtimes) 2 *index))
		  ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
		  (($bern) ((mtimes) 2 *index))
		  ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) -1)))
		 *index 0 $inf)
	 sp2)

(defprop %sec
	 ((%sum) ((mtimes) ((mexpt) -1 *index)
			   ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
			   (($euler) ((mtimes) 2 *index))
			   ((mexpt) sp2var ((mtimes) 2 *index)))
		 *index 0 $inf)
	 sp2)

; These are the series definitions of exponential functions.

(defprop mexpt
	 ((%sum)
	  ((mtimes) ((mexpt) ((mfactorial) *index) -1) ((mexpt) sp2var *index))
	  *index 0 $inf)
	 sp2)

(defprop %sinh
	 ((%sum) ((mtimes)
		  ((mexpt) ((mfactorial) ((mplus) ((mtimes) 2 *index) 1)) -1)
		  ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) 1)))
		 *index 0 $inf)
	 sp2)

(defprop %cosh
	 ((%sum) ((mtimes)
		  ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
		  ((mexpt) sp2var ((mtimes) 2 *index)))
		 *index 0 $inf)
	 sp2)

(defprop %tanh
	((%sum)
	 ((mtimes) ((mexpt) 4 *index)
		   ((mplus) ((mexpt) 4 *index) -1)
		   (($bern) ((mtimes) 2 *index))
		   ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) -1))
		   ((mexpt)
		    ((mfactorial) ((mtimes) 2 *index))
		    -1))
	 *index 0 $inf)
	sp2)

(defprop %coth
	((%sum)
	 ((mtimes) ((mexpt) 4 *index)
		   (($bern) ((mtimes) 2 *index))
		   ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
		   ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) -1)))
	 *index 0 $inf)
	sp2)

(defprop %sech 
	((%sum)
	 ((mtimes) (($euler) ((mtimes) 2 *index))
		   ((mexpt) ((mfactorial) ((mtimes) 2 *index)) -1)
		   ((mexpt) sp2var ((mtimes) 2 *index)))
	 *index 0 $inf)
	sp2)

(defprop %csch
	((%sum)
	 ((mtimes) -2 ((mplus) ((mexpt) 2 ((mplus) ((mtimes) 2 *index) -1)) -1)
		   ((mexpt) ((mfactorial) ((mtimes) *index 2)) -1)
		   (($bern) ((mtimes) 2 *index))
		   ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) -1)))
	 *index 0 $inf)
	sp2)

;arc trigonometric function expansions

(defprop %asin
	((%sum)
	 ((mtimes) ((%genfact) ((mplus) ((mtimes) 2 *index) -1) *index 2)
		   ((mexpt) ((%genfact) ((mtimes) 2 *index) *index 2) -1)
	   	   ((mexpt) ((mplus) ((mtimes) 2 *index) 1) -1)
		   ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) 1)))
	 *index 0 $inf)
	sp2)

(defprop %atan
	((%sum)
	 ((mtimes) ((mexpt) -1 *index)
		   ((mexpt) ((mplus) ((mtimes) 2 *index) 1) -1)
		   ((mexpt) sp2var ((mplus) ((mtimes) 2 *index) 1)))
	 *index 0 $inf)
	sp2)