diff --git a/changelog b/changelog
index 7ab42f6..55bebdf 100644
--- a/changelog
+++ b/changelog
@@ -1,3 +1,6 @@
+20120420 tpd src/axiom-website/patches.html 20120420.01.tpd.patch
+20120420 tpd src/input/Makefile add cohen.input
+20120420 tpd src/input/cohen.input Joel Cohen algebra example
 20120413 tpd src/axiom-website/patches.html 20120413.02.tpd.patch
 20120413 tpd src/axiom-website/download.html update download list
 20120413 tpd src/axiom-website/patches.html 20120413.01.tpd.patch
diff --git a/src/axiom-website/patches.html b/src/axiom-website/patches.html
index a023da0..2432cc5 100644
--- a/src/axiom-website/patches.html
+++ b/src/axiom-website/patches.html
@@ -3882,5 +3882,7 @@ faq add FAQ 53: Axiom won't build on Fedora<br/>
 Makefile.pamphlet add <<GCLOPTS-CUSTRELOC>> to all stanzas<br/>
 <a href="patches/20120413.02.tpd.patch">20120413.02.tpd.patch</a>
 src/axiom-website/download.html update download list<br/>
+<a href="patches/20120420.01.tpd.patch">20120420.01.tpd.patch</a>
+src/input/cohen.input Joel Cohen algebra example<br/>
  </body>
 </html>
diff --git a/src/input/Makefile.pamphlet b/src/input/Makefile.pamphlet
index 5cf6fad..38f0d42 100644
--- a/src/input/Makefile.pamphlet
+++ b/src/input/Makefile.pamphlet
@@ -301,7 +301,7 @@ REGRESSTESTS= ackermann.regress \
     carten.regress    cclass.regress   char.regress     ch.regress \
     chtheorem.regress classtalk.regress clements.regress \
     clifford.regress  clif.regress     cmds.regress \
-    coercels.regress  collect.regress \
+    coercels.regress  collect.regress  cohen.regress \
     complex.regress   complexfactor.regress conformal.regress \
     constant.regress  contfrac.regress contfrc.regress \
     curl.regress      curry.regress    cwmmt.regress \
@@ -649,7 +649,7 @@ FILES= ${OUT}/ackermann.input \
        ${OUT}/clifford.input ${OUT}/clements.input \
        ${OUT}/clif.input       ${OUT}/cmds.input \
        ${OUT}/coercels.input ${OUT}/collect.input    ${OUT}/color.input \
-       ${OUT}/complex.input  ${OUT}/complexfactor.input \
+       ${OUT}/complex.input  ${OUT}/complexfactor.input ${OUT}/cohen.input \
        ${OUT}/cone.input       ${OUT}/conformal.input \
        ${OUT}/constant.input \
        ${OUT}/contfrac.input ${OUT}/contfrc.input    ${OUT}/coordsys.input \
@@ -938,7 +938,7 @@ DOCFILES= \
   ${DOC}/clif.input.dvi        ${DOC}/cmds.input.dvi       \
   ${DOC}/coercels.input.dvi    ${DOC}/collect.input.dvi    \
   ${DOC}/color.input.dvi       ${DOC}/complex.input.dvi    \
-  ${DOC}/complexfactor.input.dvi \
+  ${DOC}/complexfactor.input.dvi ${DOC}/cohen.input.dvi    \
   ${DOC}/cone.input.dvi        ${DOC}/conformal.input.dvi  \
   ${DOC}/constant.input.dvi    ${DOC}/contfrac.input.dvi   \
   ${DOC}/contfrc.input.dvi     ${DOC}/coordsys.input.dvi   \
diff --git a/src/input/cohen.input.pamphlet b/src/input/cohen.input.pamphlet
new file mode 100644
index 0000000..9bbd52d
--- /dev/null
+++ b/src/input/cohen.input.pamphlet
@@ -0,0 +1,3317 @@
+\documentclass{article}
+\usepackage{axiom}
+\begin{document}
+\title{\$SPAD/src/algebra cohen.spad}
+\author{Timothy Daly}
+\maketitle
+\begin{abstract}
+\end{abstract}
+\newpage
+\tableofcontents
+\newpage
+
+\section{category CCAT CohenCategory}
+\begin{verbatim}
+  kind:(CExpr)->Boolean
+    ++ \spad{kind}(u) This operator returns the type of expressions
+    ++ (e.g. symbol, integer, fraction, +, *, ^, not, equal, less,
+    ++ lessorequal,greater, greaterequal, notequal, and, or ,not, set,
+    ++ list, and function names. E.g.
+    ++ kind(m*x+b) -> +
+    ++ If u is a compound expression, this returns the operator at the root
+    ++ of the expression tree.
+  numberOfOperand:(CExpr)->Integer
+    ++ \spad{numberOfOperand}(u) The operator returns the
+    ++ number of operands of the main operator u, e.g.
+    ++ numberOfOperands(a*x+b*x+c) -> 3
+    ++ If u is not a compound expression this returns "failed"
+  operand:(CExpr,Integer)->CExpr
+    ++ \spad{operand}(u,i) This operator returns the ith operand
+    ++ of the main operator u. e.g.
+    ++ operand(m*x+b,2) -> b
+    ++ If u is not a compound expression or u does not have an ith operand
+    ++ this returns "failed"
+    ++
+    ++ Operand(a+b+c,2) -> b
+    ++ Operand({a,b,c},3) -> c
+    ++ Operand(a=b,2) -> b
+  construct:(CExpr,List(CExpr))->CExpr
+    ++ \spad{construct}(f,L)->CExpr Let f be an operator (+, *, =
+    ++ etc.) or a symbol and let L=[a,b,...,c] be a list of expressions. This
+    ++ operator returns an expression with main operator f and operands 
+    ++ a,b,...,c. E.g.
+    ++ construct("+",[a,b,c]) -> a+b+c
+    ++ construct("*",[a+b,c+d,e+f]) -> (a+b)*(c+d)*(e+f)
+    ++ construct(g,[a,b,c]) -> g(a,b,c)
+  freeOf:(CExpr,CExpr)->Boolean
+    ++ \spad{freeOf}(u,t) Let u and t (for target) be mathematical 
+    ++ expressions. This operator returns false when t is identical to some
+    ++ complete subexpression of u and otherwise returns true. E.g.
+    ++ freeOf((a+b)c,a+b) -> false
+    ++
+    ++ freeeOf(a+b,b) -> false
+    ++ freeeOf(a+b,c) -> true
+    ++ freeeOf((a+b)*c,a+b) -> false
+    ++ freeeOf(sin(x)+2*x,sin(x)) -> false
+    ++ freeeOf((a+b+c)*d,a+b) -> true
+    ++ freeeOf((y+2*x-y)/x,x) -> true
+    ++ freeeOf((x*y)^2,x*y) -> true
+  substitute:(CExpr,CEq)->CExpr
+    ++ \spad{substitute}(u,t=r) Let u, t, and r be mathematical expressions.
+    ++ This operator forms a new exspression with each occurence of the target
+    ++ expression t in u replaced by the replacement expression r. The
+    ++ substitution occurs whenever t is structurally identical to a complete
+    ++ sub-expression of u. E.g.
+    ++ substitute((a+b)c,a+b=x) -> xc
+    ++
+    ++ substitute(2*x+1,x=b+1) -> 2*(b+1)+1
+  sequentialSubstitute:(CExpr,CEq) -> CExpr
+    ++ \spad{sequentialSubstitute}(u,L) let u be an expression and let L be
+    ++ a list of equations L={t1=r1,t2=r2,...tn=rn} where the targets ti are
+    ++ distinct. This operator retuns the expression u_n that is defined by
+    ++ the sequence of structural substitutions
+    ++  u_1 := substitute(u.t1=r1)
+    ++  u_2 := substitute(u.t2=r2)
+    ++  ...
+    ++  u_n := substitute(u.tn=rn)
+    ++ 
+    ++ u:=x^2+x+2
+    ++ v:=x^2+3*x-7
+    ++ w:=x^2-5*x+4
+    ++ algebraicExpand(sequentialSubstitute(u,[x=v,x=w]))
+    ++ x^8-20x^7+172x^6-830x^5+2439x^4-4390x^3+4573x^2-2365x+464
+    ++ algebraicExpand(sequentialSubstitute(u,[x=w,x=v]))
+    ++ x^8+12x^7+16x^6-234x^5-407x^4+2202x^3+1479x^2-10089x+7834
+  concurrentSubstitute:(CExpr,CEq) -> CExpr
+    ++ \spad{concurrentSubstitute}(u,S) let u be an expression and S be the
+    ++ set of equations S={t1=t1,t2=r2,...,tn=rn} where the targets ti are
+    ++ distinct. This operator returns a new expression defined by a
+    ++ recursive search thru the expression tree of u and compare each
+    ++ complete sub-expression v to each of the targets ti. If v is identical
+    ++ to some ti, substitute the corresponding ri for v
+  map:(CExpr,CExpr)->CExpr
+    ++ \spad{map}(f,u) Let u be a mathematical expression with 
+    ++ n=numberOfOperands(u) >= 1 and let f(x) be an operator.
+    ++ The \spad{map}(f,u) returns a new expression with main operator
+    ++ kind(u) and operands f(operand(u,1))...f(operand(u,n))
+  map:(CExpr,CExpr,List CId)->CExpr
+    ++ \spad{map}(f,u,c) Let u be a mathematical expression with 
+    ++ n=numberOfOperands(u) >= 1 and let f(x) be an operator and
+    ++ c be a list of identifiers
+    ++ The \spad{map}(f,u,c) returns a new expression with main operator
+    ++ kind(u) and operands g(operand(u,1),c[i],...g(operand(u,n),c[i])
+  degreeGPE:(CExpr,CExpr)->Integer
+    ++ \spad{degreeGPE}(u,x) Let u=c_1..c_r*x_1^(n_1)..x_m^(n_m) be a
+    ++ monomial with non-zero coefficient part. The degree of u with respect
+    ++ to x_i is denoted by deg(u,x_i)=n_i. By mathematical convention, the
+    ++ degree of the 0 monomial is minusInfinity. If u is a polynomial in x,
+    ++ that is a sum of monomials, then deg(u,x_i) is the maximum of the
+    ++ degrees of the monomials. If the generalized variable x_i is
+    ++ understood from context, we use the simpler notation deg(u). E.g.
+    ++ degreeGPE(sin^2(x)+b*sin(x)+c,sin(x)) -> 2
+  leadingCoefficientGPE:(CExpr,CId)->CExpr
+    ++ \spad{leadingCoefficientGPE}(u,x) Let u be a polynomial in x. The 
+    ++ leading coefficient of a polynomial u != 0 with respect to x is the
+    ++ sum of the coefficient parts of all monomials with variable parts
+    ++ x^deg(u,x). The zero polynomial has, by definition, leading
+    ++ coefficient zero. The leading coefficient is represented by lc(u,x)
+    ++ and when x is understood from context by lc(u). This operator returns
+    ++ lc(u,x). E.g.
+    ++ leadingCoefficientGPE(a*x+b*x+y,x)->a+b
+  monomialGPE:(CExpr,CId) -> Boolean
+    ++ \spad{monomialGPE}(u,x) Let u be an algebraic expression and let v be
+    ++ either a variable or s set of variables. This operator returns true
+    ++ whenever u is a polymonial in {x} or in the set of variabesl, otherwise
+    ++ false
+    ++ 
+    ++ monomialGPE(a*x^2*y^2,{x,y}) -> true
+    ++ monomialGPE(x^2+y^2,{x,y}) -> false
+  polynomialGPE:(CExpr,CId) -> Boolean
+    ++ \spad{polynomialGPE}(u,x) Let u be an algebraic expression and let v be
+    ++ either a variable or s set of variables. This operator returns true
+    ++ whenever u is a polymonial in {x} or in the set of variabesl, otherwise
+    ++ false
+    ++ 
+    ++ polynomialGPE(x^2+y^2,{x,y}) -> true
+    ++ polynomialGPE(sin(x)^2+2*sin(x)+3,sin(x)) -> true
+    ++ polynomialGPE(x/y+2*y,{x,y}) -> false
+    ++ polynomialGPE((x+1)*(x+3),x) -> false
+  coefficientGPE:(CExpr,CId,NonNegativeInteger) -> Union(CExpr,"failed")
+    ++ \spad{coefficientGPE}(u,x,j) Let u be an algebraic expression. If u
+    ++ is a GPE in a variable x and j >=0 then this operator returns the sum
+    ++ of the coefficient parts of all monomials of u with variable part x^j.
+    ++ If there is no monomial with variable part x^j, this returns 0. If u
+    ++ is not a polynomial in x, this returns "failed"
+    ++ 
+    ++ coefficientGPE(a*x^2+b*x+c,x,2) -> a
+    ++ coefficientGPE(2*x*y^2+5*x^2*y+7*x+9,x,1) -> 3*y^2+7
+    ++ coefficientGPE(3*x*y^2+5*x^2*y+7*x+9,x,3) -> 0
+    ++ coefficientGPE((2*sin(x))*x^2+(2*ln(x))*x+4,x,2) -> failed
+    ++ coefficientGPE(coefficientGPE(3*x*y^2+5*x^2*y+7*x+9,x,1),y,2) -> 3
+    ++ coefficientGPE(coefficientGPE(3*sin(x)*x^2+2*ln(x)*x+4,ln(x),1),x,1)->2
+    ++ coefficientGPE(coefficientGPE(3*sin(x)*x^2+2*ln(x)*x+4,x,1),ln(x),1)
+    ++    -> failed
+  variables:CExpr->List(CExpr)
+    ++ \spad{variables}(u) The polynomial structure of an algebraic 
+    ++ expression u depends on which expressions are chosen as the generalized
+    ++ variables. This operator selects a set of generalized variables so
+    ++ that the coefficients of all monomials in u are rational numbers. E.g.
+    ++ 
+    ++ variables(4*x^3+3*x^2*sin(x))->[x,sin(x)]
+    ++ variables(x^3+3*x^2*y+3*x*y^2+y^3) -> {x,y}
+    ++ variables(3*x*(x+1)*y^2*z^n) -> {x,x+1,y,z^n}
+    ++ variables(a*sin(x)^2+2*b*sin(x)+3*c) -> {a,b,c,sin(x)}
+    ++ variables(1/2) -> {}
+    ++ varialbes(sqrt(2)*x^2+sqrt(3)*x+sqrt(5)) -> {x,sqrt(2),sqrt(3),sqrt(5)}
+  deg:(CExpr) -> Integer
+    ++ \spad{deg}(u) If u is a monomial in a single variable then let
+    ++ u=c1..cr*x1^n1..xm^nm be a monomial with non-zero coefficient part.
+    ++ The degree of u is the sum of the exponents
+    ++ of the variables: deg(u) = n. The degree of the 0 monomial
+    ++ is minusInfinity. If u is a sum of monomials then deg(u) is the
+    ++ maximum of the degrees of the monomials. 
+  deg:(CExpr,List(CId)) -> Integer
+    ++ \spad{deg}(u,S) Let S = {x1,..xn} be a set of variables. Let
+    ++ u=c1..cr*x1^n1..xm^nm be a monomial with non-zero coefficient part.
+    ++ The degree of u with respect to the set S is the sum of the exponents
+    ++ of the variables: deg(u,S) = n1+n2+..._nm. The degree of the 0 monomial
+    ++ is minusInfinity. If u is a sum of monomials then deg(u,S) is the
+    ++ maximum of the degrees of the monomials. 
+  operandList:CExpr->List CExpr
+    ++ \spad{operandList}(u)
+  absoluteValue:CExpr->CExpr
+    ++ \spad{absoluteValue}(u)
+    ++ 
+    ++ absoluteValue(-2) -> 2
+  coefficientSv:(CExpr,CId,NonNegativeInteger) -> Union(CExpr,"failed")
+    ++ \spad{coefficientSv}(u,x,j) Let u be an algebraic expression. If u is
+    ++ a polynomial in x then this operator returns the coefficient u_j of
+    ++ x^j in u. If j > deg(u,x) then the result is 0. If u is not a 
+    ++ polynomial in x, the operator returns "failed"
+    ++ 
+    ++ coefficientSv(x^2_3*x+5,x,1) -> 3
+    ++ coefficientSv(2*x^3+3*x,x,4) -> 0
+    ++ coefficientSv(3,x,0) -> 3
+    ++ coefficientSv((x+1)*(x+3),x,2) -> failed
+  leadingCoefficientSv:(CExpr,CId)->Union(CExpr,"failed")
+    ++ \spad{leadingCoefficientSv}(u,x) Let u be an algebraic expression. If
+    ++ u is a polynomial in x then this operator returns the 
+    ++ leadingCoefficient(u,x). If u is not a polynomial in x, this returns
+    ++ "failed"
+    ++ 
+    ++ leadingCoefficientSv(x^2+3*x+5,x) -> 1
+    ++ leadingCoefficientSv(3,x) -> 3
+    ++ leadingCoefficientSv == coefficientSv(u,x,degreeSv(u,x)) 
+  max:List(CExpr)->CExpr
+    ++ \spad{max}(u)
+  min:List(CExpr)->CExpr
+    ++ \spad{min}(u)
+  numerator:CExpr->CExpr
+    ++ \spad{numerator}(u)
+  denominator:CExpr->CExpr
+    ++ \spad{denominator}(u)
+  integerQuotient:(CInt,CInt)->CInt
+    ++ \spad{integerQuotient}(a,b)
+  integerRemainder:(CInt,CInt)->CInt
+    ++ \spad{integerRemainder}(a,b)
+  integerGCD:(CInt,CInt)->CInt
+    ++ \spad{integerGCD}(a,b)
+  derivative:(CExpr,CExpr)->CExpr
+    ++ \spad{derivative}(u,v)
+    ++ 
+    ++ derivative(sin(x),x) -> cos(x)
+    ++ if u=x then derivative(u,x) -> 1
+    ++ if u=v^w then derivative(u,x) ->
+    ++   w*v^(w-1)*derivative(v,x)+derivative(w,x)*v^w+ln(v)
+    ++ if u is a sum and v=Operand(u,1) and w=u-v then 
+    ++   derivative(u,x) -> derivative(v,x)+derivative(w,x)
+    ++ if u is a product and v=Operand(u,1) and w=u/v then
+    ++   derivative(u,x) -> derivative(v,x)*w+v*derivative(w,x)
+    ++ if u = sin(x) then
+    ++   derivative(u,x) -> cos(v)*derivative(v,x)
+    ++ if freeOf(u,x) = true then
+    ++   derivative(u,x) -> 0
+  decimal:CExpr -> Float
+    ++ \spad{decimal}(u) Evaluate rational numbers, arithmetic operations, and
+    ++ numerical functions in an expression u to a real value
+    ++ 
+    ++ decimal(1/4) -> .25
+    ++ decimal(x+1/4) -> x+.25
+    ++ decimal(sin(2)+1/2) -> 1.409297
+  degree:(CExpr,CId) -> Integer
+    ++ \spad{degree}(u,x) Degree in x of polynomial expression u
+    ++ 
+    ++ degree(x^2+5*x+7,x) -> 2
+  coefficient:(CExpr,CId,Integer) -> CExpr
+    ++ \spad{coefficient}(u,x,j) Coefficient of x^j in a polynomial 
+    ++ expression u
+    ++ 
+    ++ Coefficient(x^2+5*x+7,x) -> 5
+  factor:CExpr -> CExpr
+    ++ \spad{factor}(u) polynomial factorization
+    ++ 
+    ++ factor(x^2+5*x+6) -> (x+2)*(x+3)
+  solve:(CExpr,CId) -> List(CExpr)
+    ++ \spad{solve}(u,x) Solution of an equation u for x
+    ++ 
+    ++ solve(a*x=b,x) -> [x = b/a]
+  solve:(List(CExpr),List(CId)) -> List(CExpr)
+    ++ \spad{solve}(u,x) Solution of a set of equations u a set of variables
+    ++ 
+    ++ solve([2*x+4*y=3,3*x-y=7],[x,y])->[x=31/14, y=-5/14]
+  limit:(CExpr,CId,CExpr) -> CExpr
+    ++ \spad{limit}(u,x,v)
+    ++ 
+    ++ limit(1/x,x,plusInfinity) -> 0
+  rationalSimplify:CExpr -> CExpr
+    ++ \spad{rationalSimplify}(u)
+    ++ 
+    ++ rationalSimplify(1/a+1/b) -> (a+b)/(a*b)
+    ++ rationalSimplify((x^2-1)/(x-1)) -> x+1
+    ++ rationalSimplify(1/(1/a+c/(a*b))+(a*b*c+a*c^2)/(b+c)^2) -> a
+
+  algebraicExpand:CExpr->CExpr
+    ++ \spad{algebraicExpand}(u) 
+    ++ 
+    ++ algebraicExpand((x*(y+1)^(3/2)+1)*(x*(y+1)^(3/2)-1)) -> x^2*(y_1)^3-1
+    ++                    = x^2*y^3+3*x^2*y^2+3*x^2*y+x^2-1
+    ++ algebraicExpand((x*(y+1)^(1/2)+1)^4) ->
+    ++    x^4*(y+1)^2 + 4*x^3*(y+1)^(3/2) + 6*x^2*(y+1) + 4*x*(y+1)^(1/2)+1
+    ++  = x^4*y^2+2*x^4*y+x^4+4*x^3*(y+1)^(3/2)+6*x^2*y+6*x^2+4*x*(y+1)^(1/2)+1
+    ++ algebraicExpand((x+2)*(x+3)) -> x^2+5*x+6
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p253
+    ++ Procedure Algebraic_expand(u)
+    ++ Input
+    ++   u : an algebraic expression where all exponents of powers are
+    ++       integers
+    ++ Output
+    ++   the expanded form of u
+    ++ Local Variables
+    ++   v,base,exponent
+    ++ Begin
+    ++   if Kind(u)="+" then
+    ++     v:=Operand(u,1)
+    ++     Return(Algebraic_expand(v)+Algebraic_expand(u-v))
+    ++   elseif Kind(u)="*" then
+    ++     v:=Operand(u,1)
+    ++     Return(Expand_product(Algegraic_expand(v),Algebraic_expand(u/v)))
+    ++   elseif Kind(u)="^" then
+    ++     base:=Operand(u,1)
+    ++     exponent:=Operand(u,2)
+    ++     if Kind(exponent)=integer and exponent >= 2 then
+    ++       Return(Expand_power(Algebraic_expand(base),exponent))
+    ++   Return(u)
+    ++ End
+
+  expandProduct(CExpr,CExpr)->CExpr
+    ++ \spad{expandProduct}(r,s)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p253
+    ++ Procedure Expand_product(r,s)
+    ++ Input
+    ++   r,s : expanded algebraic expressions, where all exponents of powers
+    ++         are integers
+    ++ Output
+    ++   the expanded form of r*s
+    ++ Local Variables
+    ++   f
+    ++ Begin
+    ++   if Kind(r)="+" then
+    ++     f:=Operand(r,1)
+    ++     Return(Expand_product(f,s)+Expand_product(r-f,s)
+    ++   elseif Kind(s)="+" then
+    ++     Return(Expand_product(s,r))
+    ++   else
+    ++     Return(r*s)
+    ++ End
+
+  expandPower(CExpr,NonNegativeInteger)->CExpr
+    ++ \spad{expandPower}(u,n)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p253
+    ++ Procedure Expand_power(u,n)
+    ++ Input
+    ++   u : an expanded algebraic expression where all exponents of powers
+    ++       are integers
+    ++   n : a non-negative integer
+    ++ Output
+    ++   the expanded form of u^n
+    ++ Local Variables
+    ++   f,r,k,s,c
+    ++ Begin
+    ++   if Kind(u)="+" then
+    ++     f:=Operand(u,1)
+    ++     r:=u-f
+    ++     s:=0
+    ++     for k:=0 to n do
+    ++       c:=n!/(k!*(n-k)!)
+    ++       s:=s+Expand_product(c*f^(n-k),Expand_power(r,k))
+    ++     Return(s)
+    ++   else
+    ++     Return(u^n)
+    ++ End
+
+
+  collectTerms:(CExpr,List(CId))->Union(CExpr,"failed")
+    ++ \spad{collectTerms}(u,s) Let u be an algebraic expression and s be
+    ++ a non-empty set of variables.
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p249
+    ++ Procedure Collect_terms(u,S)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   S : a non-empty set of generalized variables
+    ++ Output
+    ++   the collected form of u or the global symbol Undefined if u is
+    ++   not a polynomial in S
+    ++ Local Variables
+    ++   f,combined,i,j,N,T,v
+    ++ Begin
+    ++   if Kind(u) != "+" then
+    ++     if Coeff_var_monomial(u,S) = Undefined then
+    ++       Return(Undefined)
+    ++     else Return(u)
+    ++   else
+    ++     if u in S then Return(u)
+    ++     N:=0
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       f:=Coeff_var_monomial(Operand(u,i),S)
+    ++       if f = Undefined then Return(Undefined)
+    ++       else
+    ++         j:=1
+    ++         combined:=false
+    ++         while not combined and j <= N do
+    ++           if Operand(f,2) = Operand(T[j],2) then
+    ++             T[j]:=[Operand(f,1)+Operand(T[j],1),Operand(f,2)]
+    ++             combined:=true
+    ++           j:=j+1
+    ++         if not combined then
+    ++           T[N+1]:=f
+    ++           N:=N+1
+    ++     v:=0
+    ++     for j = 1 to N do
+    ++       v:=v+Operand(T[j],1)*Operand(T[j],2)
+    ++     Return(v)
+    ++ End
+
+  coefficientGPE:(CExpr,CId,PositiveInteger)->CExpr
+    ++ \spad{coefficientGPE}(u,x,j) Let u be a polynomial in x, and let j be
+    ++ a non-negative integer. This operator returns the sum of the 
+    ++ coefficient parts of all monomials of u with variable part x^j. E.g.
+    ++ coefficientGPE(a*x+b*x+y,x,1) -> a+b
+    ++ 
+    ++ coefficient_gpe(a*(x^2+1)^2+(x^2+1),x^2+1,1) -> 0
+    ++ coefficient_gpe(a*(x^2+1)^2+(x^2+1),x^2+1,0) -> x^2+1
+    ++ coefficient_gpe(2*(x^2)^2+3*(x^2),x^2,2)->0
+    ++ coefficient_gpe(2*(x^2)^2+3*(x^2),x^2,1)->3
+    ++ coefficient_gpe(2*(x^2)^2+3*(x^2),x^2,0)->2*x^4
+    ++ coefficient_gpe(3*x*y^2+5*x^2+y+7*x+9,x,1) -> 3*y^2+7
+    ++ coefficient_gpe(3*x*y^2+5*x^2+y+7*x+9,x,3) -> 0
+    ++ coefficient_gpe((2*sin(x))*x^2+(2*ln(x))*x+4,x,2) -> Undefined
+    ++ u := 3*x*y^2+5*x^2*y+78x+9
+    ++ coefficient_gpe(coefficient_gpe(u,x,1),y,2) -> 3
+    ++ v := 3*sin(x)*x^2+2*ln(x)*x+4
+    ++ coefficient_gpe(coefficient_gpe(v,ln(x),1),x,1) -> 2
+    ++
+    ++ Cohen, Joel "Elementary Algorithms" p233
+    ++ Procedure Coefficient_gpe(u,x,j)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a generalized variable
+    ++   j : a non-negative integer
+    ++ Output
+    ++   The coefficient of x^i in the polynomial u or the global symbol
+    ++   Undefined
+    ++ Local Variables
+    ++   i,c,f
+    ++ Begin
+    ++   if Kind(u) != "+" then
+    ++     f:=Coefficient_monomial_gpe(u,x)
+    ++     if f = Undefined then Return(Undefined)
+    ++     else
+    ++       if j = Operand(f,2) then Return(Operand(f,1))
+    ++       else Return(0)
+    ++   else
+    ++     if u = x then
+    ++       if j = 1 then Return(1) else Return(0)
+    ++     c:=0
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       f:=Coefficient_monomial_gpe(Operand(u,i),x)
+    ++       if f = Undefined then Return(Undefined)
+    ++       elseif Operand(f,2) = j then c = c+Operand(f,1)
+    ++     Return(c)
+    ++ End
+
+  coefficientMonomialGPE:(CExpr,CId)->Union(List(Integer),"failed")
+    ++ \spad{coefficientMonomialGPE}(u,v) When u is a polynomial in x, this
+    ++ procedure returns a list [c,m] where m is the degree of u in x and
+    ++ c is the coefficient part of the monomial. If u is not a monomial in x,
+    ++ the global symbol "failed" is returned.
+    ++
+    ++ Cohen, Joel "Elementary Algorithms" p232
+    ++ Procedure Coefficient_monomial_gpe(u,x)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a generalized variable
+    ++ Output
+    ++   The list [c,m] where m is the degree of the monomial and
+    ++   c is the coefficient of x^m or the global symbol Undefined
+    ++ Local Variables
+    ++   base, exponent, i, c, m, f
+    ++ Begin
+    ++   if u = x then
+    ++     Return([1,1])
+    ++   elseif Kind(u) = "^" then
+    ++     base:=Operand(u,1)
+    ++     exponent:=Operand(u,2)
+    ++     if base = x and Kind(exponent) = integer and exponent > 1 then
+    ++       Return([1,exponent])
+    ++   elseif Kind(u) = "+" then
+    ++     m:=0
+    ++     c:=u
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       f:=Coefficient_monomial_gpe(Operand(u,i),x)
+    ++       if f = Undefined then
+    ++         Return(Undefined)
+    ++       elseif Operand(f,2) != 0 then
+    ++         m:=Operand(f,2)
+    ++         c:=u/x^m
+    ++     Return([c,m])
+    ++   if Free_of(u,x) then
+    ++     Return([u,0])
+    ++   else
+    ++     Return(Undefined)
+    ++ End
+
+
+  monomialGPE:(CExpr,CId)->CExpr
+    ++ \spad{monomialGRE}(u,v) Let u be an algebraic expression, and let v
+    ++ be either a generalized variable x or a set S of generalized variables.
+    ++ This operator returns true whenever u is a polynomial in {x} or in S,
+    ++ and otherwise returns false.
+    ++
+    ++ Cohen, Joel "Elementary Algorithms" p227
+    ++ Procedure Monomial_gpe(u,v)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   v : a generalized variable or a set of generalized variables
+    ++ Output
+    ++   true or false
+    ++ Local Variables
+    ++   i,S,base,exponent
+    ++ Begin
+    ++   if Kind(v) != set then S:={v} else S:=v
+    ++   if u in S then
+    ++     Return(true)
+    ++   elseif Kind(u) = "^" then
+    ++     base:=Operand(u,1)
+    ++     exponent:=Operand(u,2)
+    ++     if base in S and Kind(exponent) = integer and exponent > 1 then
+    ++       Return(true)
+    ++   elseif Kind(u) = "+" then
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       if Monomial_gpe(Operand(u,i),S)=false then
+    ++         Return(false)
+    ++     Return(true)
+    ++   Return(Set_free_of(u,S))
+    ++ End
+
+  polynomialGPE:(CExpr,List CId)->Boolean
+    ++ \spad{polynomialGPE}(u,v) Let u be an algebraic expression and let v
+    ++ be either a generalized variable x or a set S of generalized variables.
+    ++ The \spad{polynomialGPE}(u,v) returns true whenenver u is a polynomial
+    ++ in x or in S and otherwise it returns false. E.g.
+    ++ polynomialGPE(x^2+y^2,[x,y]) -> true
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p228
+    ++ Procedure Polynomial_gpe(u,v):
+    ++ Input
+    ++   u : an algebraic expression
+    ++   v : a generalized variable or a set of generalized variables
+    ++ Output
+    ++   true or false
+    ++ Local Variables
+    ++   i,Si
+    ++ Begin
+    ++   if Kind(v) != set then S:={v} else S:=v;
+    ++   if Kind(u) != "+" then
+    ++     Return(Monomial_gpe(u,x))
+    ++   else
+    ++     if u in S then Return(true);
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       if Monomial_gpe(Operand(u,i),S) = false then
+    ++         Return(false);
+    ++     Return(true)
+    ++ End
+    ++ 
+    ++ Polynomial_gpe(x^2+y^2,{x,y}) -> true
+    ++ Polynomial_gpe(sin^2(x)+2*sin(x)+3,sin(x)) -> true
+    ++ Polynomial_gpe(x/y+2*y,{x,y}) -> false
+    ++ Polynomial_gpe((x+1)*(x+3),x) -> false 
+
+  rationalizeExpression:CExpr->CExpr
+    ++ \spad{rationalizeExpression}(u)
+    ++ 
+    ++ Rationalize_Expression((1+1/x)^2) -> ((x+1)^2)/x^2
+    ++ Rationalize_Expression((1+1/x)^(1/2)) -> ((x+1)/x)^(1/2)
+    ++ Rationalize_Expression(x/z+y/(z^2)) -> (z^2*x+z*y)/(z^3)
+    ++ Cohen, Joel "Elementary Algorithms" p265
+    ++ Procedure Rationalize_expression(u)
+    ++ Input
+    ++  u : an algebraic expression
+    ++ Output
+    ++   a rationalized form of u
+    ++ Local Variables
+    ++   f,g,r
+    ++ Begin
+    ++   if Kind(u) = "^" then
+    ++     Return(Rationalize_expression(Operand(u,1))^(Operand(u,2)))
+    ++   elseif Kind(u)="*" then
+    ++     f:=Operand(u,1)
+    ++     Return(rationalize_expression(f)*Rationalize_expression(u/f))
+    ++   elseif Kind(u)="+" then
+    ++     f:=Operand(u,1)
+    ++     g:=Rationalize_expression(f)
+    ++     r:=Rationalize_expression(u-f)
+    ++     Return(rationalize_sum(g,r))
+    ++   else
+    ++     Return(u)
+    ++ End
+
+  rationalizeSum:(CExpr,CExpr)->CExpr
+    ++ \spad{rationalizeSum}(u,v) computes the sum
+    ++ m/r + n/s -> (m*s+n*r)/(r*s)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p265
+    ++ Procedure Rationalize_sum(u,v)
+    ++ Input
+    ++   u,v : algebraic expression in rationalized form
+    ++ Output
+    ++   an algebraic expression in rationalized form
+    ++ Local Variables
+    ++   m,n,r,s
+    ++ Begin
+    ++   m:=Numerator(u)
+    ++   r:=Denominator(u)
+    ++   n:=Numerator(v)
+    ++   s:=Denominator(v)
+    ++   if r=1 and s=1 then
+    ++     Return(u+v)
+    ++   else
+    ++     Return(Rationalize_sum(m*s,n_r)/(r*s))
+    ++ End
+
+  expandExponential:CExpr -> CExpr
+    ++ \spad{expandExponential}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p278
+    ++ Procedure Expd_exp(u)
+    ++ Input
+    ++   u: an algebraic expression
+    ++ Output
+    ++   an algebraic expression in exponential-expanded form
+    ++ Local Variables
+    ++   v,A,f
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else
+    ++     v:=Map(Expd_exp,u)
+    ++     if Kind(v)=exp then
+    ++       A:=Operand(v,1)
+    ++       if Kind(A) = "+" then
+    ++         f:=Operand(A,1)
+    ++         Return(Expd_exp(exp(f))*Expd_exp(exp(A-f)))
+    ++     elseif Kind(A)="*" then
+    ++       f:=Operand(A,1)
+    ++       if Kind(f) = integer then
+    ++         Return(Expd_exp(exp(A/f))^f)
+    ++   Return(v)
+    ++ End
+
+  expandExponential1:CExpr -> CExpr
+    ++ \spad{expandExponential1}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p280
+    ++ Procedure Expand_exp(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   an algebraic expression in exponential-expanded form
+    ++ Local Variables
+    ++   v
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else 
+    ++     v:=Map(Expand_exp,u)
+    ++     if Kind(v)=exp then
+    ++       Return(Expand_exp_rules(Operand(v,1)))
+    ++     else
+    ++       Return(v)
+    ++ End
+
+  expandExponentialRules:CExpr -> CExpr
+    ++ \spad{expandExponentialRules}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p280
+    ++ Procedure Expand_exp_rules(A)
+    ++ Input
+    ++   A : an algebraic expression that is the argument of an exponential
+    ++       function
+    ++ Output
+    ++   the exponential-expanded form of exp(A)
+    ++ Local Variables
+    ++   f
+    ++ Begin
+    ++   if Kind(A)="+" then
+    ++     f:=Operand(A,1)
+    ++     Return(Expand_exp_rules(f)*Expand_exp_rules(A-f))
+    ++   elseif
+    ++     f:=Operand(A,1)
+    ++     if Kind(f) = integer then
+    ++       Return(Expand_exp_rules(A/f)^f)
+    ++   Return(exp(A))
+    ++ End
+
+  expandTrig:CExpr -> CExpr
+    ++ \spad{expandTrig}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p283
+    ++ Procedure Expand_trig(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   an algebraic expression in trigonometric-expanded form
+    ++ Local Variables
+    ++   v
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else
+    ++     v:=Map(expand_trig,u)
+    ++     if Kind(v)=sin then
+    ++       Return(Operand(Expand_trig_rules(Operand(v,1)),1))
+    ++     elseif Kind(v)=cos then
+    ++       Return(Operand(Expand_trig_rules(Operand(v,1)),2))
+    ++     else
+    ++       Return(v)
+    ++ End
+
+  expandTrigRules:CExpr -> List(CExpr)
+    ++ \spad{expandTrigRules}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p283
+    ++ Procedure Expand_trig_rules(A)
+    ++ Input
+    ++   A : an algebraic expression that is the argument of a sin or cos
+    ++ Output
+    ++   a two element list [s,c] where s and c are the 
+    ++   trigonometric-expanded forms of sin(A) and cos(A)
+    ++ Local variables
+    ++   f,r,s,c
+    ++ Begin
+    ++   if Kind(A)="+" then
+    ++     f:=Expand_trig_rules(Operand(A,1))
+    ++     r:=Expand_trig_rules(A-Operand(A,1))
+    ++     s:=Operand(f,1)*Operand(r,2)+Operand(f,2)*Operand(r,1)
+    ++     c:=Operand(f,2)*Operand(r,2)-Operand(f,1)*Operand(r,1)
+    ++     Return([s,c])
+    ++   elseif Kind(A)="*" then
+    ++     f:=Operand(A,1)
+    ++     if Kind(f)=integer then
+    ++       Return([Multiple_angle_sin(f,A/f),Multiple_angle_cos(f,A/f)])
+    ++   Return([sin(A),cos(A)])
+    ++ End
+
+  contractExp:CExpr -> CExpr
+    ++ \spad{contractExp}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p290
+    ++ Procedure Contract_exp(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   an algebraic expression in exponential-contracted form
+    ++ Local Variables
+    ++   v
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else
+    ++     v:=Map(Contract_exp,u)
+    ++     if Kind(v) in {"*", "^"} then
+    ++       Return(contract_exp_rules(v))
+    ++     else
+    ++       Return(v)
+    ++ End
+
+  contractExpRules:CExpr -> CExpr
+    ++ \spad{contractExpRules}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p291
+    ++ Procedure Contract_exp_rules(u)
+    ++ Input
+    ++   u : an algebraic expresson that is sent by either Contract_exp
+    ++       or a recursive call to this procedure
+    ++ Output
+    ++   an algebraic expression in exponential-contracted form
+    ++ Local Variables
+    ++   v,b,s,p,i,y
+    ++ Begin
+    ++   v:=Expand_main_op(u)
+    ++   if Kind(v)="^" then
+    ++     b:=Operand(v,1)
+    ++     v:=Operand(v,2
+    ++     if Kind(b)=exp then
+    ++       p:=Operand(b,1)*s
+    ++       if Kind(p) in {"*","^"} then
+    ++         p:=Constract_exp_rules(p)
+    ++       Return(exp(p))
+    ++     else
+    ++       Return(v)
+    ++   elseif Kind(v)="*" then
+    ++     p:=1
+    ++     s:=0
+    ++     for i:=1 to Number_of_operands(v) do
+    ++       y:=Operand(v,i)
+    ++       if Kind(y)=exp then
+    ++         s:=s+Operand(y,1)
+    ++       else
+    ++         p:=p*y
+    ++     Return(exp(s)*p)
+    ++   elseif Kind(v) = "+"
+    ++     s:=0
+    ++     for i:=1 to Number_of_operands(v) do
+    ++       y:=Operand(v,i)
+    ++       if Kind(y) in {"*","^"} then
+    ++         s:=s+Contract_exp_rules(y)
+    ++       else
+    ++         s:=s+y
+    ++     Return(s)
+    ++   else
+    ++     Return(v)
+    ++ End
+
+  contractTrig:CExpr -> CExpr
+    ++ \spad{contractTrig}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p296
+    ++ Procedure Contract_trig(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   an algebraic expression in trigonometric-contracted form
+    ++ Local Variables
+    ++   v
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else
+    ++     v:=Map(Contract_trig,u)
+    ++     if Kind(v) in {"*","^"} then
+    ++       Return(Contract_trig_rules(v))
+    ++     else
+    ++       Return(v)
+    ++ End
+
+  contractTrigRules:CExpr -> CExpr
+    ++ \spad{contractTrigRules}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p297
+    ++ Procedure Contract_trig_rules(u)
+    ++ Input
+    ++   u : an algebraic expression (sum, product, power) that is sent by
+    ++       either Contract_trig, Contract_trig_power, or a recursive call
+    ++       of this procedure
+    ++ Output
+    ++   an algebraic expression in trigonometric-contracted form
+    ++ Local Variables
+    ++   v,s,c,d,i,y
+    ++ Begin
+    ++   v:=Expand_main_op(u)
+    ++   if Kind(v) = "^" then
+    ++     Return(Contract_trig_power(v))
+    ++   elseif Kind(v) = "*" then
+    ++     s:=Separate_sin_cos(v)
+    ++     c:=Operand(s,1)
+    ++     d:=Operand(s,2)
+    ++     if d = 1 then
+    ++       Return(v)
+    ++     if Kind(d) in {sin,cos} then
+    ++       Return(v)
+    ++     elseif Kind(d) = "^" then
+    ++       Return(Expand_main_op(c*Contract_trig_power(d)))
+    ++     else
+    ++       Return(Expand_main_op(c*Contract_trig_product(d)))
+    ++   elseif Kind(v)="+" then
+    ++     s:=0
+    ++     for i:=1 to Number_of_operands(v) do
+    ++       y:=Operand(v,i)
+    ++       if Kind(y) in {"*","^"} then
+    ++         s:=s+Contract_trig_rules(y)
+    ++       else
+    ++        s:=s+y
+    ++     Return(s)
+    ++   else
+    ++     Return(v)
+    ++ End
+
+  contractTrigProduct:CExpr -> CExpr
+    ++ \spad{contractTrigProduct}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p298
+    ++ Procedure Construct_trig_product(u)
+    ++ Input
+    ++   u : a product of sines, cosines, and positive integer powers of
+    ++       sines and cosines
+    ++ Output
+    ++   the trigonometric-contracted form of u
+    ++ Local Variables
+    ++   A,B,theta,phi
+    ++ Begin
+    ++   if Number_of_operands(u) = 2 then
+    ++     A:=Operand(u,1)
+    ++     B:=Operand(u,2)
+    ++     if Kind(A)="^" then
+    ++       A:=Construct_trig_power(A)
+    ++       Return(Contruct_trig_rules(A*b))
+    ++     elseif Kind(B)="^" then
+    ++       B:=Construct_trig_power(B)
+    ++       Return(Contruct_trig_rules(A*B))
+    ++     else
+    ++       theta:=Operand(A,1)
+    ++       phi:=Operand(B,1)
+    ++       if Kind(A) = sin and Kind(B) = sin then
+    ++         Return(cos(theta-phi)/2 - cos(theta_phi)/2)
+    ++       elseif Kind(A) = cos and Kind(B) = cos then
+    ++         Return(cos(theta+phi)/2+cos(theta-phi)/2)
+    ++       elseif Kind(A) = sin and Kind(B) = cos then
+    ++         Return(sin(theta+phi)/2+sin(theta-phi)/2)
+    ++       elseif Kind(A) = cos and Kind(B) = sin then
+    ++         Return(sin(theta+phi)/2+sin(phi-theta)/2)
+    ++   else
+    ++     A:=Operand(u,1)
+    ++     B:=Contract_trig_product(u/A)
+    ++     Return(Contract_trig_rules(A*B))
+    ++ End
+
+  simplifyTrig:CExpr -> Union(CExpr,"failed")
+    ++ \spad{simplifyTrig}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p303
+    ++ Procedure Simplify_trig(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   either an algebraic expression in trigonometric contracted form
+    ++   or the global symbol Undefined
+    ++ Local Variables
+    ++   v,w,n,d
+    ++ Begin
+    ++   v:=Trig_substitute(u)
+    ++   w:=Rationalize_expression(v)
+    ++   n:=Expand_trig(Numerator(w))
+    ++   n:=Contract_trig(n)
+    ++   d:=Expand_trig(Denominator(w))
+    ++   d:=Contract_trig(d)
+    ++   if d = 0 then
+    ++     Return(Undefined)
+    ++   else
+    ++     Return(n/d)
+    ++ End
+
+  degreeSv:(CExpr,CId) -> Union(Integer,"failed")
+    ++ \spad{degreeSv}(u,x) Let u be an algebraic expresson. If u is a 
+    ++ polynomial in x, this operator returns deg(u,x). If u is not a
+    ++ polynomial in x, the operator returns the symbol "failed"
+    ++ 
+    ++ Degree_sv(3*x^2+4*x+5,x) -> 2
+    ++ Degree_sv(2*x^3,x) -> 3
+    ++ Degree_sv((x+1)*(x+3),x) -> Undefined
+    ++ Degree_sv(3,x) -> 0
+    ++ Cohen, Joel "Elementary Algorithms" p220
+    ++ Procedure degree_sv(u,x)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   deg(u,x) or the global symbol Undefined
+    ++ Local Variables
+    ++   d,i,f
+    ++ Begin
+    ++   d:=Degree_monomial_sv(u,x)
+    ++   if d != Undefined then
+    ++     Return(d)
+    ++   elseif Kind(u)="+" then
+    ++     d:=0
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       f:=Degree_monomial_sv(Operand(u,i),x)
+    ++       if f = Undefined then
+    ++         Return(Undefined)
+    ++       else
+    ++         d:=Max({d,f})
+    ++     Return(d)
+    ++   Return(Undefined)
+    ++ End
+
+  degreeMonomialSv:(CExpr,Cid) -> Union(Integer,"failed")
+    ++ \spad{degreeMonomialSv}(u,x)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p219
+    ++ Procedure Degree_monomial_sv(u,x)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   deg(u,x) or the global symbol Undefined
+    ++ Local Variables
+    ++   base, exponent, s,t
+    ++ Begin
+    ++   if u = 0 then
+    ++     Return(-oo)
+    ++   elseif Kind(u) in {integer,fraction} then
+    ++     Return(0)
+    ++   elseif u = x then
+    ++     Return(1)
+    ++   elseif Kind(u)="^" then
+    ++     base := Operand(u,1)
+    ++     exponent := Operand(u,2)
+    ++     if base = x and Kind(exponent) = integer and exponent > 1 then
+    ++       Return(exponent)
+    ++   elseif Kind(u)="+" then
+    ++     if Number_of_operands(y) = 2 then
+    ++       s := Degree_monomial_sv(Operand(u,1),x)
+    ++       t := Degree_monomial_sv(Operand(u,2),x)
+    ++       if s != Undefined and t != Undefined then
+    ++         Return(t)
+    ++   Return(Undefined)
+    ++ End
+
+  polynomialSv:(CExpr,Cid) -> Boolean
+    ++ \spad{polynomialSv}(u,x) Let u be an algebraic expression. This
+    ++ operator returns true when u is a polynomial in x and otherwise false.
+    ++ The suffix 'sv' stands for 'single variable'
+    ++ 
+    ++ Polynomial_sv(3*x^2+4*x+5,x) -> true
+    ++ Polynomial_sv(1/(x+1),x) -> false
+    ++ Polynomial_sv(a*x^2+b*x+c,x) -> false
+    ++ Cohen, Joel "Elementary Algorithms" p218
+    ++ Procedure Polynomial_sv(u,x)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   true or false
+    ++ Local Variables
+    ++   i
+    ++ Begin
+    ++   if Monomial_sv(u,x) then
+    ++     Return(true)
+    ++   elseif Kind(u)="+" then
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       if Monomial_sv(Operand(u,i),x) = false then
+    ++         Return(false)
+    ++     Return(true)
+    ++   Return(false)
+    ++ End
+
+  monomialSv:(CExpr,Cid) -> Boolean
+    ++ \spad{monomialSv}(u,x) Let u be an algebraic expression. This operator
+    ++ returns true when u is a monomial in x and otherwise returns false.
+    ++ The suffix 'Sv' stands for 'single variable'
+    ++ 
+    ++ monomialSv(2*x^3,x) -> true
+    ++ monomialSv(x+1,x) -> false
+    ++ Cohen, Joel "Elementary Algorithms" p217
+    ++ Procedure Monomial_sv(u,x)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   true or false
+    ++ Local Variables
+    ++   base,exponent
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction} then
+    ++     Return(true)
+    ++   elseif u = x then
+    ++     Return(true)
+    ++   elseif Kind(u)="^" then
+    ++     base := Operand(u,1)
+    ++     exponent := Operand(u,2)
+    ++     if base = x and Kind(exponent) = integer and exponent > 1 then
+    ++       Return(true)
+    ++   elseif Kind(u) = "+" then
+    ++     Return(Number_of_operand(u) = 2 and Monomial_sv(Operand(u,1),x)
+    ++            and Monomial_sv(Operand(u,2),x))
+    ++   Return(false)
+    ++ End
+
+  integral:(CExpr,Cid) -> Union(CExpr,"failed")
+    ++ \spad{integral}(u,x)
+    ++ 
+    ++ integral(cos(x),x) -> sin(x)
+    ++ Cohen, Joel "Elementary Algorithms" p206
+    ++ Procedure Integral(f,x)
+    ++ Input
+    ++   f : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   int f dx or the global symbol Fail
+    ++ Local Variables
+    ++   F, g
+    ++ Begin
+    ++   F:=Integral_table(f,x)
+    ++   if F:=Fail then
+    ++     F:=Linear_properties(f,x)
+    ++   if F:=Fail then
+    ++     F:=Substitution_method(f,x)
+    ++   if F:=Fail then
+    ++     g:=Algebraic_expand(f)
+    ++     if f != g then
+    ++       F:=Integral(g,x)
+    ++   Return(F)
+    ++ End
+
+  substitutionMethod:(CExpr,Cid) -> Union(CExpr,"failed")
+    ++ \spad{substitutionMethod}(f,x)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p206
+    ++ Procedure Substituion_method(f,x)
+    ++ Input
+    ++   f : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   int f dx or the global symbol Fail
+    ++ Local Variables
+    ++   P,F,i,u,g
+    ++ Global
+    ++  v
+    ++ Begin
+    ++   P:=Trial_substitution(f)
+    ++   F:=Fail
+    ++   i:=1
+    ++   while F = Fail and i <= Number_of_operands(P) do
+    ++     g:=Operand(P,i)
+    ++     if g != x and not Free_of(g,x) then
+    ++       u:=Substitute(f/Derivative(g,x),g=v)
+    ++       if Free_of(u,x) then
+    ++         F:=Substitute(Integral(u,v), v=g)
+    ++     i:=i+1
+    ++   Return(F)
+    ++ End
+
+  legendre:(NonNegativeInteger,CId) -> CExpr
+    ++ \spad{legendre}(n,x)
+    ++ 
+    ++ (Note that this recursive version is inefficient. See legendreIter)
+    ++ p_0(x) = 1
+    ++ p_1(x) = x
+    ++ p_2(x) = (1/2)*((2*2-1)*x*p_1(x)    -(2-1)*p_0(x)) = (3/2)*x^2-(1/2)
+    ++ p_n(x) = (1/n)*((2*n-1)*x*p_(n-1)(x)-(n-1)*p_(n-2)(x)) | n >= 2
+    ++ Cohen, Joel "Elementary Algorithms" p191
+    ++ Procedure Legendre(n,x)
+    ++ Input
+    ++   n : a non-negative integer
+    ++   x : a symbol
+    ++ Output
+    ++   p_n(x)
+    ++ Local Variables
+    ++   f,g
+    ++ Begin
+    ++   if n = 0 then
+    ++     Return(1)
+    ++   elseif n = 1 then
+    ++     Return(x)
+    ++   else
+    ++     f:=Legendre(n-1,x)
+    ++     g:=Legendre(n-2,x)
+    ++     Return(Algebraic_expand((1/n)*((2*n-1)*x_f-(n-1)*g)))
+    ++ End
+
+  legendreIter:(NonNegativeInteger,CId) -> CExpr
+    ++ \spad{legendre}(n,x)
+    ++ 
+    ++ (Note that this recursive version is inefficient. See legendreIter)
+    ++ p_0(x) = 1
+    ++ p_1(x) = x
+    ++ p_2(x) = (1/2)*((2*2-1)*x*p_1(x)    -(2-1)*p_0(x)) = (3/2)*x^2-(1/2)
+    ++ p_n(x) = (1/n)*((2*n-1)*x*p_(n-1)(x)-(n-1)*p_(n-2)(x)) | n >= 2
+    ++ Cohen, Joel "Elementary Algorithms" p191
+    ++ Procedure Legendre(n,x)
+    ++ Input
+    ++   n : a non-negative integer
+    ++   x : a symbol
+    ++ Output
+    ++   p_n(x)
+    ++ Local Variables
+    ++   f,g
+    ++ Begin
+    ++   if n = 0 then
+    ++     Return(1)
+    ++   elseif n = 1 then
+    ++     Return(x)
+    ++   else
+    ++     f:=1
+    ++     g:=x
+    ++     for i:=2 to n do
+    ++       p:=(1/i)*((2*i-1)*x*g-(i-1)*f
+    ++       f:=g
+    ++       g:=p
+    ++     Return(Algebraic_expand(p))
+    ++ End
+
+  trigSubstituteMap:CExpr -> CExpr
+    ++ \spad{trigSubstituteMap}(u)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p190
+    ++ Procedure Trig_substitute_map(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   a new expression where all instances of the functions 
+    ++   tan, cot, sec, and csc are replaced by the representations
+    ++   using sin and cos
+    ++ Local Variables
+    ++   U
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else
+    ++     U:=Map(Trig_substitute_map,u)
+    ++     if Kind(U) = tan then
+    ++       Return(sin(Operand(U,1))/cos(Operand(U,1)))
+    ++     if Kind(U) = cot then
+    ++       Return(cos(Operand(U,1))/sin(Operand(U,1)))
+    ++     if Kind(U) = sec then
+    ++       Return(1/cos(Operand(U,1)))
+    ++     if Kind(U) = csc then
+    ++       Return(1/sin(Operand(U,1)))
+    ++     else
+    ++       Return(U)
+    ++ End
+
+  trigSubstitute:CExpr -> CExpr
+    ++ \spad{trigSubstitute}(u) Replace trig operators with identities:
+    ++ tan(v) -> sin(v)/cos(v)
+    ++ cot(v) -> cos(v)/sin(v)
+    ++ sec(v) -> 1/cos(v)
+    ++ csc(v) -> 1/sin(v)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p187
+    ++ Procedure Trig_substitute(u)
+    ++ Input
+    ++   u : an algebraic expression
+    ++ Output
+    ++   a new expression with all instances of the functions
+    ++   tan, cot, sec, csc replaced by the representations using
+    ++   sin and cos
+    ++ Local Variables
+    ++   s,i,L
+    ++ Begin
+    ++   if Kind(u) in {integer, fraction, symbol} then
+    ++     Return(u)
+    ++   else
+    ++     L:=[]
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       L:=Join(L,[Trig_substitute(Operand(u,1))])
+    ++     if Kind(u) in {tan,cot,sec,csc} then
+    ++       s:=Operand(L,1)
+    ++       if Kind(u) = tan then
+    ++         Return(sin(s)/cos(s))
+    ++       if Kind(u) = cot then
+    ++         Return(cos(s)/sin(s))
+    ++       if Kind(u) = sec then
+    ++         Return(1/cos(s))
+    ++       if Kind(u) = csc then
+    ++         Return(1/sin(s))
+    ++     else
+    ++       Return(Construct(Kind(u),L))
+    ++ End
+
+  linearForm:(CExpr,Cid) -> Union(List(CExpr),"failed")
+    ++ \spad{linearForm}(u,x)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p181
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   the list [a,b] where a and b are algebraic expressions or the
+    ++   global symbol Fail
+    ++ Local Variables
+    ++   f,r
+    ++ Begin
+    ++   if u = x then
+    ++     Return([1,0])
+    ++   elseif Kind(u) in {symbol, integer, fraction} then
+    ++     Return([0,u])
+    ++   elseif Kind(u) = "*" then
+    ++     if Free_of(u,x) then
+    ++       Return([0,u])
+    ++     elseif Free_of(u/x,x) then
+    ++       Return([u/x,0])
+    ++     else
+    ++       Return(Fail)
+    ++   elseif Kind(u) = "+" then
+    ++     f:=Linear_form(Operand(u,1),x)
+    ++     if f = Fail then
+    ++       Return(Fail)
+    ++     else
+    ++       r:=Linear_form(u-Operand(u,1),x)
+    ++       if r = Fail then
+    ++         Return(Fail)
+    ++       else
+    ++         Return([Operand(f,1)+Operand(r,1),Operand(f,2)+Operand(r,2)])
+    ++   elseif Free_of(u,x) then
+    ++     Return([0,u])
+    ++   else
+    ++     Return(Fail)
+    ++ End
+
+  freeOf:(CExpr,CExpr) -> Boolean
+    ++ \spad{freeOf}(u,t)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p179
+    ++ Procedure Free_of(u,t)
+    ++ Input
+    ++   u,t : mathematical expressions
+    ++ Output
+    ++   true or false
+    ++ Local Variables
+    ++   i
+    ++ Begin
+    ++   if u = t then
+    ++     Return(false)
+    ++   elseif Kind(u) in {symbol, integer, real} then
+    ++     Return(true)
+    ++   else
+    ++     i:=1
+    ++     while i <= Number_of_operands(u) do
+    ++       if not Free_of(Operand(u,i),t) then
+    ++         Return(false)
+    ++       i:=i+1
+    ++     Return(true)
+    ++ End
+
+  completeSubExpressions:CExpr -> List(CExpr)
+    ++ \spad{completeSubExpressions}(u)
+    ++ 
+    ++ completeSubExpressions(a*(x+1)+3*cos(y)) ->
+    ++ [a*(x+1)+3*cos(y),a*(x+1),a,x+1,x,1,3*cos(y),3,cos(y),y]
+    ++ Cohen, Joel "Elementary Algorithms" p177
+    ++ Procedure Complete_sub_expressions(u)
+    ++ Input
+    ++   u : a mathematical expression
+    ++ Output
+    ++   the set of complete sub-expressions of u
+    ++ Local Variables
+    ++   s,i
+    ++ Begin
+    ++   if Kind(u) in {integer, symbol, read} then
+    ++     Return({u})
+    ++   else
+    ++     s:={u}
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       s:- s union Complete_sub_expressions(Operand(u,i))
+    ++     Return(s)
+    ++ End
+
+  recFact:NonNegativeInteger -> NonNegativeInteger
+    ++ \spad{recFact}(n)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p172
+    ++ Procedure Rec_fact(n)
+    ++ Input
+    ++   n : non-negative integer
+    ++ Output
+    ++   n!
+    ++ Local Variables
+    ++   f
+    ++ Begin
+    ++   if n = 0 then
+    ++     f:=1
+    ++   else
+    ++     f:=n*Rec_fact(n-1)
+    ++   Return(f)
+    ++ End
+
+  solveODE:(CExpr,CId,CId) -> Union(CExpr,"failed")
+    ++ \spad{solveODE}(w,x,y)
+    ++ 
+    ++ solveODE(derivative(y(x),x)=y(x),x,y) -> C*exp(x)
+    ++ Cohen, Joel "Elementary Algorithms" p162
+    ++ Procedure Solve_ode(w,x,y)
+    ++ Input
+    ++   w : a differential equation that can be transformed by rational
+    ++       simplification to the form M + N*(dy/dx) = 0
+    ++       where the derivative dy/dx is represented by the function
+    ++       form d(y,x)
+    ++   x,y : Symbols
+    ++ Output
+    ++   An implicit solution to the differential equation or the
+    ++   global symbol Fail
+    ++ Local Variables
+    ++   p,M,N,F
+    ++ Begin
+    ++   p:=Transform_ode(w,x,y)
+    ++   M:=Operand(p,1)
+    ++   N:=Operand(p,2)
+    ++   F:=Separable_ode(M,N,x,y)
+    ++   if F = Fail then
+    ++     F:-Solve_exaxt(M,N,x,y)
+    ++   Return(F)
+    ++ End
+
+  transformODE:(CExpr,CId,CId) -> List(CExpr)
+    ++ \spad{transformODE}(w,x,y)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p162
+    ++ Procedure Transform_ode(w,x,y)
+    ++ Input
+    ++   w : a differential equation that can be transformed by rational
+    ++       simplification to the form M + N*(dy/dx) = 0
+    ++       where the derivative dy/dx is represented by the function
+    ++       form d(y,x)
+    ++   x,y : Symbols
+    ++ Output
+    ++   the list [M,N]
+    ++ Local Variables
+    ++   v,n,M,N
+    ++ Begin
+    ++   v:=Rational_simplify(Operand(w,1)-Operand(w,2))
+    ++   n:=Numerator(v)
+    ++   M:=Coefficient(n,d(y,x),0)
+    ++   N:=Coefficient(n,d(y,x),1)
+    ++   Return([M,N])
+    ++ End
+
+  solveExact:(CExpr,CExpr,Cid,CId) -> Union(CExpr,"failed")
+    ++ \spad{solveExact}(m,n,x,y)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p163
+    ++ Procedure Solve_exact(M,N,x,y)
+    ++ Input
+    ++   M,N: algebraic expressions
+    ++   x,y: symbols
+    ++ Output
+    ++   An implicit solution to the differential equation or the global
+    ++   symbol Fail
+    ++ Local Variables
+    ++   My, Nx, d, u, F, G, g, h, hp
+    ++ Begin
+    ++   if N = 0 then
+    ++     Return(Fail)
+    ++   elseif M = 0 then
+    ++     Return(y=C)
+    ++   My:=Derivative(M,y)
+    ++   Nx:=Derivative(N,x)
+    ++   d:=My-Nx
+    ++   if d = 0 then
+    ++     u:=1
+    ++   else
+    ++     F:=Rational_simplify(d/N)
+    ++     if Free_of(F,y) then
+    ++       u:=exp(Integral(F,x))
+    ++       d:=0
+    ++     else
+    ++       G:=Rational_simplify(-d/M)
+    ++       if Free_of(G,x) then
+    ++         u:-exp(Integral(G,y))
+    ++         d:=0
+    ++   if d = 0 then
+    ++     g:=Integral(u*M,x)
+    ++     hp:=u*N-Derivative(g,y)
+    ++     h:=Integral(hp,y)
+    ++     Return(g+h=C)
+    ++   else
+    ++     Return(Fail)
+    ++ End
+
+  separateFactors:(CExpr,CExpr) -> List(CExpr)
+    ++ \spad{separateFactors}(u,x)
+    ++ 
+    ++ Separate_factors((c*x*sin(x))/2,x) -> [c/2,x*sin(x)]
+    ++ Cohen, Joel "Elementary Algorithms" p148
+    ++ Procedure Spearate_factors(u,x)
+    ++ Input
+    ++   u,x : algebraic expressions
+    ++ Output
+    ++   a list with two algebraic expressions
+    ++ Local Variables
+    ++   f, free_of_part,dependent_part,i
+    ++ Begin
+    ++   if Kind(u) = "*" then
+    ++     free_of_part:=1
+    ++     dependent_part:=1
+    ++     for i:=1 to Number_of_operands(u) do
+    ++       f:=Operand(u,i)
+    ++       if Free_of(f,x) then
+    ++         free_of_part:=f*free_of_part
+    ++       else
+    ++         dependent_part:=f*dependent_part
+    ++     Return([free_of_part,dependent_part])
+    ++   else
+    ++     if Free_of(u,x) then
+    ++       Return([u,1])
+    ++     else
+    ++       Return([1,u])
+    ++ End
+
+  evenOdd:(CExpr,CId) -> Symbol
+    ++ \spad{evenOdd}(u,x)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p144
+    ++ Procedure Even_odd(u,x)
+    ++ Input
+    ++   u : an algebraic expression
+    ++   x : a symbol
+    ++ Output
+    ++   one of the global symbols Even, Odd, or Unknown
+    ++ Local Variables
+    ++   v
+    ++ Begin
+    ++   v:=Substitute(u,x=-x)
+    ++   if u-v = 0 then
+    ++     Return(Even)
+    ++   elseif u+v = 0 then
+    ++     Return(Odd)
+    ++   else
+    ++     Return(Unknown)
+    ++ End
+
+  tangentLine:(CExpr,CId,CExpr) -> CExpr
+    ++ \spad{tangentLine}(f,x,a)
+    ++ 
+    ++ Cohen, Joel "Elementary Algorithms" p136
+    ++ Procedure Tangent_line(f,x,a)
+    ++ Input
+    ++   f : an algebraic expression (formula for a mathematical function)
+    ++   x : a symbol (independent variable)
+    ++   a : an algebraic expresion (point of tangency)
+    ++ Output
+    ++   an algebraic expression that is the formula for the tangent line
+    ++ Local Variables
+    ++   deriv,m,line
+    ++ Begin
+    ++   deriv:=Derivative(f,x)
+    ++   m:=Substitute(deriv,x=a)
+    ++   line:=Algebraic_expand(m*(x-a)+Substitute(f,x-a))
+    ++   Return(line)
+    ++ End
+\end{verbatim}
+
+
+\begin{chunk}{domain CINT CohenInteger}
+)abbrev domain CINT CohenInteger
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++ Elements of this domain are symbolic integers. Then do not act like
+++ ordinary integers since expressions involving them will not simplify
+++ automatically. Thus, if 1 and 2 are declared as CohenIntegers then\br
+++ \tab{5}1 + 2\br
+++ will not be automatically simplified to 3
+
+CohenInteger(): SetCategory with
+
+  coerce : Integer -> %
+ 
+  coerce : % -> Integer 
+
+ == add 
+
+  Rep := List Integer
+
+  coerce(i:Integer):% == [i]
+ 
+  coerce(cint:%):Integer == cint.1
+
+  coerce(cint:%):OutputForm == (cint.1)::OutputForm
+
+\end{chunk}
+
+\begin{chunk}{domain CID CohenIdentifier}
+)abbrev domain CID CohenIdentifier
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++
+++ Elements of this domain are symbolic identifiers. Then do not act like
+++ ordinary symbols since expressions involving them will not simplify
+++ automatically. Thus, if x is declared as a CohenIdentifier then\br
+++ \tab{5}x + x\br
+++ will not be automatically simplified to 2x. Note that identifiers
+++ have a "value cell" that has a setter and a getter. It defaults to 0.
+
+CohenIdentifier(): SetCategory with
+
+  coerce : Symbol -> %
+ 
+  coerce : % -> Symbol
+
+  set : (%,Integer) -> %
+
+  get : % -> Integer
+
+  name : % -> Symbol
+
+ == add 
+
+  Rep := Record(sym:Symbol,val:Integer)
+
+  coerce(s:Symbol):% == [s,0]
+ 
+  coerce(cid:%):Symbol == cid.sym
+
+  set(s:%,i:Integer):% == [s.sym,i]
+
+  get(cid:%):Integer == cid.val
+
+  name(cid:%):Symbol == cid.sym
+
+  coerce(cid:%):OutputForm == (cid.sym)::OutputForm
+
+\end{chunk}
+
+\begin{chunk}{domain COP CohenOperator}
+)abbrev domain COP CohenOperator
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++ CohenOperators are prefix, infix, or postfix. Prefix operators
+++ assume arity 1, infix assumes arity 2, and postfix assumes arity 0
+++ They are not applied automatically but are inert identifiers in an
+++ expression. You have to call a function to operate on them.
+++
+++ Note that functions can be associated with the operators.
+
+CohenOperator(): SetCategory with
+
+  copPrefix : Symbol -> %
+  copPrefix : (Symbol,NonNegativeInteger) -> %
+
+  copInfix : Symbol -> %
+  copInfix : (Symbol,NonNegativeInteger) -> %
+
+  copPostfix : Symbol -> %
+  copPostfix : (Symbol,NonNegativeInteger) -> %
+
+  name : % -> Symbol
+
+  arity : % -> NonNegativeInteger
+
+  prefix? : % -> Boolean
+
+  infix? : % -> Boolean
+
+  postfix? : % -> Boolean
+
+  coerce: % -> Symbol
+
+  coerce: Symbol -> %
+ 
+ == add 
+  
+  prefix   ==> 0
+  infix    ==> 1
+  postfix  ==> 2
+
+  Rep := Record(name:Symbol,arity:NonNegativeInteger,fix:NonNegativeInteger)
+
+  copPrefix(f:Symbol):% == [f,1,prefix]
+  copPrefix(f:Symbol,n:NonNegativeInteger):% == [f,n,prefix]
+
+  copInfix(f:Symbol):% == [f,2,infix]
+  copInfix(f:Symbol,n:NonNegativeInteger):% == [f,n,infix]
+
+  copPostfix(f:Symbol):% == [f,0,postfix]
+  copPostfix(f:Symbol,n:NonNegativeInteger):% == [f,n,postfix]
+
+  name(op:%):Symbol == op.name
+
+  arity(op:%):NonNegativeInteger == op.arity
+
+  prefix?(op:%):Boolean == op.fix = prefix
+ 
+  infix?(op:%):Boolean == op.fix = infix
+
+  postfix?(op:%):Boolean == op.fix = postfix
+
+  coerce(op:Symbol):% == copPrefix(op)
+
+  coerce(op:%):Symbol == op.name
+
+  coerce(cn:%):OutputForm == cn.name::OutputForm
+
+\end{chunk}
+
+\begin{chunk}{domain CFUNC CohenFunction}
+)abbrev domain CFUNC CohenFunction
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++ CohenOperators are prefix, infix, or postfix. Prefix operators
+++ assume arity 1, infix assumes arity 2, and postfix assumes arity 0
+++ They are not applied automatically but are inert identifiers in an
+++ expression. You have to call a function to operate on them.
+++
+++ Note that functions can be associated with the operators.
+
+CohenFunction(): SetCategory with
+
+  cFunction : Symbol -> %
+  cFunction : (Symbol,NonNegativeInteger) -> %
+
+  name : % -> Symbol
+
+  arity : % -> NonNegativeInteger
+
+  coerce: % -> Symbol
+
+  coerce: Symbol -> %
+ 
+ == add 
+  
+  function ==> 0
+
+  Rep := Record(name:Symbol,arity:NonNegativeInteger,fix:NonNegativeInteger)
+
+  cFunction(f:Symbol):% == [f,1,function]
+  cFunction(f:Symbol,n:NonNegativeInteger):% == [f,n,function]
+
+  name(op:%):Symbol == op.name
+
+  arity(op:%):NonNegativeInteger == op.arity
+
+  function?(op:%):Boolean == op.fix = function
+ 
+  coerce(op:%):Symbol == op.name
+
+  coerce(s:Symbol):% == cFunction(s)
+
+  coerce(cn:%):OutputForm == cn.name::OutputForm
+
+\end{chunk}
+
+\begin{chunk}{domain CNODE CohenNode}
+)abbrev domain CNODE CohenNode
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++ Elements of this domain are purely symbolic expressions represented
+++ in a tree structure. All operations on the elements are purely
+++ procedural so expressions are not simplified unless the simplification
+++ is requested. Thus, members of this domain can be used in program to
+++ explain the steps of a calculation, such as subtracting 1 from both
+++ sides of an equation and displaying the unsimplified result:\br
+++ \tab{5}x + 1 - 1 = 5 - 1
+
+CohenNode(): Export == Implementation where
+
+ CINT  ==> CohenInteger
+ CID   ==> CohenIdentifier
+ COP   ==> CohenOperator
+ CFUNC ==> CohenFunction
+ NNI   ==> NonNegativeInteger
+
+ Export ==> SetCategory with
+
+  intnode : CohenInteger -> %
+  coerce : CohenInteger -> %
+  retract : % -> CohenInteger
+  retractIfCan : % -> Union(CohenInteger,"failed")
+  cohenInteger? : % -> Boolean
+
+  idnode : CohenIdentifier -> %
+  coerce : CohenIdentifier -> %
+  retract : % -> CohenIdentifier
+  retractIfCan : % -> Union(CohenIdentifier,"failed")
+  cohenIdentifier? : % -> Boolean
+
+  opnode : CohenOperator -> %
+  coerce : CohenOperator -> %
+  retract : % -> CohenOperator
+  retractIfCan : % -> Union(CohenOperator,"failed")
+  cohenOperator? : % -> Boolean
+
+  fnnode : CohenFunction -> %
+  coerce : CohenFunction -> %
+  retract : % -> CohenFunction
+  retractIfCan : % -> Union(CohenFunction,"failed")
+  cohenFunction? : % -> Boolean
+
+ Implementation ==> add 
+
+  cinttag ==> 0
+  cidtag  ==> 1
+  coptag  ==> 2
+  cfntag  ==> 3
+
+  NODE ==> Record(int:CINT,id:CID,op:COP,fn:CFUNC)
+
+  Rep := Record(tag:NNI, node:NODE)
+
+  intnode(cint:CINT):% == 
+    rec:NODE :=
+     [cint,"failed"::Symbol::CID,"failed"::Symbol::COP,"failed"::Symbol::CFUNC]
+    [cinttag,rec]
+  coerce(cint:CINT):% == intnode(cint)
+  retract(cn:%):CINT == cn.node.int
+  retractIfCan(cn:%):Union(CINT,"failed") ==
+    cn.tag = cinttag =>  cn.node.int
+    "failed"
+  cohenInteger?(cn:%):Boolean == cn.tag = cinttag
+
+  idnode(cid:CID):% ==
+    rec:NODE :=
+     [0::Integer::CINT,cid,"failed"::Symbol::COP,"failed"::Symbol::CFUNC]
+    [cidtag,rec]
+  coerce(cid:CID):% == idnode(cid)
+  retract(cn:%):CID == cn.node.id
+  retractIfCan(cn:%):Union(CID,"failed") ==
+    cn.tag = cidtag =>  cn.node.id
+    "failed"
+  cohenIdentifier?(cn:%):Boolean == cn.tag = cidtag
+
+  opnode(cop:COP):% ==
+    rec:NODE := 
+     [0::Integer::CINT,"failed"::Symbol::CID,cop,"failed"::Symbol::CFUNC]
+    [coptag,rec]
+  coerce(cop:COP):% == opnode(cop)
+  retract(cn:%):COP == cn.node.op
+  retractIfCan(cn:%):Union(COP,"failed") ==
+    cn.tag = coptag =>  cn.node.op
+    "failed"
+  cohenOperator?(cn:%):Boolean == cn.tag = coptag
+
+  fnnode(cfn:CFUNC):% ==
+    rec:NODE :=
+     [0::Integer::CINT,"failed"::Symbol::CID,"failed"::Symbol::COP,cfn]
+    [cfntag,rec]
+  coerce(cfn:CFUNC):% == fnnode(cfn)
+  retract(cn:%):CFUNC == cn.node.fn
+  retractIfCan(cn:%):Union(CFUNC,"failed") ==
+    cn.tag = cfntag =>  cn.node.fn
+    "failed"
+  cohenFunction?(cn:%):Boolean == cn.tag = cfntag
+
+  coerce(cn:%):OutputForm ==
+    cn.tag = cinttag => (retract(cn)@CINT)::OutputForm
+    cn.tag = cidtag  => (retract(cn)@CID)::OutputForm
+    cn.tag = coptag  => (retract(cn)@COP)::OutputForm
+    cn.tag = cfntag  => (retract(cn)@CFUNC)::OutputForm
+    "failed"::OutputForm
+
+\end{chunk}
+
+\begin{chunk}{domain CREP CohenRepresentation}
+)abbrev domain CREP CohenRepresentation
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++ Elements of this domain are purely symbolic expressions represented
+++ in a tree structure. All operations on the elements are purely
+++ procedural so expressions are not simplified unless the simplification
+++ is requested. Thus, members of this domain can be used in program to
+++ explain the steps of a calculation, such as subtracting 1 from both
+++ sides of an equation and displaying the unsimplified result:\br
+++ \tab{5}x + 1 - 1 = 5 - 1
+
+CohenRepresentation(): SetCategory with
+  makeNode : CohenNode -> %
+  makeNode : (CohenNode, List %) -> %
+  empty?  : % -> Boolean
+  empty : () -> %
+  explain? : % -> Boolean
+  explain : (%,Boolean) -> Boolean
+  coerce: % -> OutputForm
+  getData : % -> CohenNode
+  getKids : % -> List %
+  parens? : % -> Boolean
+  parens : (%,Boolean) -> Boolean
+  
+ == add
+
+  NODE ==> Record(data: CohenNode, kids: List %)
+  Rep:= Union(node:NODE,empty:"empty")
+
+  t:%                      -- handle the empty case
+
+  getData(cr:%):CohenNode ==
+    cr.node.data
+
+  getKids(cr:%):List % ==
+    empty?(cr) => [empty()]
+    cr.node.kids
+
+  empty?(t:%):Boolean == t case empty
+  empty():% == ["empty"]
+
+  noisy:Boolean := false   -- default to being quiet
+
+  explain?(x:%):Boolean == noisy 
+  explain(x:%,onoff:Boolean):Boolean == noisy := onoff
+
+  parens:Boolean := false   -- default to being quiet
+
+  parens?(x:%):Boolean == noisy 
+  parens(x:%,onoff:Boolean):Boolean == parens := onoff
+
+  flatten(x:%):OutputForm ==
+    t1:=(x.node.data)@CohenNode
+    empty?(x.node.kids) => t1::OutputForm
+    empty?(first x.node.kids) => t1::OutputForm
+    cohenInteger?(t1)$CohenNode => retract(t1)@CohenInteger::OutputForm
+    cohenIdentifier?(t1)$CohenNode => retract(t1)@CohenIdentifier::OutputForm
+    cohenFunction?(t1)$CohenNode => 
+      t2:CohenFunction:=retract(t1)$CohenNode
+      t4:=[flatten(j) for j in x.node.kids]
+      hconcat [name(t2)::OutputForm,paren(commaSeparate(t4))]
+    cohenOperator?(t1)$CohenNode =>
+      t2:CohenOperator:=retract(t1)$CohenNode
+      prefix?(t2) => 
+        t4:=[flatten(j) for j in x.node.kids]
+        if parens
+        then paren(blankSeparate(cons(t2::OutputForm,t4)))::OutputForm
+        else blankSeparate(cons(t2::OutputForm,t4))::OutputForm
+      infix?(t2) => 
+        t4:=flatten(first x.node.kids)
+        t5:=[[name(t2)::OutputForm,flatten(j)] for j in rest(x.node.kids)]
+        t6:=cons(cons(t4,nil),t5)
+        t7:=concat t6
+        if parens
+        then paren(blankSeparate(t7))::OutputForm
+        else blankSeparate(t7)::OutputForm
+      postfix?(t2) => 
+        t4:=[flatten(j) for j in x.node.kids]
+        if parens
+        then paren(blankSeparate(append(t4,list(t2::OutputForm))))::OutputForm
+        else blankSeparate(append(t4,list(t2::OutputForm)))::OutputForm
+    "failed"::OutputForm
+
+  coerce(x:%):OutputForm == 
+    empty? x => "()"::OutputForm
+    flatten x
+    
+  makeNode(cnode:CohenNode):% ==  [[cnode,[empty()]]]
+
+  makeNode(cnode:CohenNode,babes:List %):% == [[cnode,babes]]
+
+\end{chunk}
+
+\begin{chunk}{domain CPOLY CohenPolynomial}
+)abbrev package CPOLY CohenPolynomial
+
+-- Author: Timothy Daly (based on Joel Cohen's books)
+-- Date Created: 7 April 2007
+-- Date Last Updated: 
+-- Basic Functions:
+-- Related Constructors: Fraction
+-- Also See:
+-- AMS Classifications:
+-- Keywords:
+-- References:
+-- "Computer Algebra and Symbolic Computation: Mathematical Methods"
+-- A.K. Peters, Ltd. ISBN 1-56881-159-4
+-- "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+-- A.K. Peters, Ltd. ISBN 1-56881-158-6
+-- Description:
+-- This type represents symbolic algebra as detailed in Joel Cohen's
+-- books on Computer Algebra and Symbolic Computation
+--
+-- Elements of this domain are purely symbolic expressions represented
+-- in a tree structure. All operations on the elements are purely
+-- procedural so expressions are not simplified unless the simplification
+-- is requested. Thus, members of this domain can be used in program to
+-- explain the steps of a calculation, such as subtracting 1 from both
+-- sides of an equation and displaying the unsimplified result:\br
+-- \tab{5}x + 1 - 1 = 5 - 1
+
+CohenPolynomial(): SetCategory with
+
+  poly2cohen : Polynomial(Fraction(Integer)) -> CohenRepresentation
+
+  cohen2poly : CohenRepresentation -> Polynomial(Fraction(Integer))
+
+  "+" : (CohenRepresentation,CohenRepresentation) -> CohenRepresentation
+
+  "-" : (CohenRepresentation,CohenRepresentation) -> CohenRepresentation
+
+  "*" : (CohenRepresentation,CohenRepresentation) -> CohenRepresentation
+
+ == add
+
+  ((a:CohenRepresentation) + (b:CohenRepresentation)):CohenRepresentation ==
+     makeNode(opnode(copInfix(_+)),[a,b])    
+
+  ((a:CohenRepresentation) - (b:CohenRepresentation)):CohenRepresentation ==
+     makeNode(opnode(copInfix(_-)),[a,b])    
+
+  ((a:CohenRepresentation) * (b:CohenRepresentation)):CohenRepresentation ==
+     makeNode(opnode(copInfix(_*)),[a,b])    
+
+  poly2cohen(poly:Polynomial(Fraction(Integer))):CohenRepresentation ==
+    hasvar:=variables poly
+    hasvar = [] => 
+      d1:Fraction(Integer):=retract(first coefficients poly)
+      numer:=numer d1
+      denom:=denom d1
+      denom = 1 => makeNode(intnode(numer::CohenInteger))
+      d2:CohenRepresentation:=makeNode(intnode(numer::CohenInteger))
+      d3:CohenRepresentation:=makeNode(intnode(denom::CohenInteger))
+      makeNode(opnode(copInfix(_/)),[d2,d3])
+    var:=first hasvar
+    deg:=degree(poly,var)
+    terms:CohenRepresentation:=empty()
+    monom:CohenRepresentation:=empty()
+    for c in 0..deg repeat
+      t1:Fraction(Integer):=retract(coefficient(poly,var,c))
+      numer:=numer t1
+      denom:=denom t1
+      if denom = 1 
+      then coef:=makeNode(intnode(numer::CohenInteger))
+      else
+        w1:CohenRepresentation:=makeNode(intnode(numer::CohenInteger))
+        w2:CohenRepresentation:=makeNode(intnode(denom::CohenInteger))
+        coef:=makeNode(opnode(copInfix(_/)),[w1,w2])
+      if t1 ~= 0 then
+        if c = 0
+        then 
+          monom:=coef
+        else
+          t2:=makeNode(idnode(var::CohenIdentifier))
+          t3:=makeNode(intnode(c::CohenInteger))
+          t4:=makeNode(opnode(copInfix(_^)),[t2,t3])
+          monom:=makeNode(opnode(copInfix(_*)),[coef,t4])
+        if empty?(terms)
+        then terms:=monom
+        else terms:=makeNode(opnode(copInfix(_+)),[monom,terms])
+    terms    
+
+  expon(x:Symbol,n:Integer):Polynomial(Fraction(Integer)) ==
+    monomial(1::Polynomial(Fraction(Integer)),x,n pretend NonNegativeInteger)
+
+  formPoly(cr:CohenRepresentation):Polynomial(Fraction(Integer)) ==
+    cn:CohenNode:=getData(cr)
+    ck:List(CohenRepresentation):=getKids(cr)
+    cohenInteger? cn => 
+      ((retract(cn)::CohenInteger)::Integer)::Polynomial(Fraction(Integer))
+    cohenIdentifier? cn => 
+      ((retract(cn)::CohenIdentifier)::Symbol)::Polynomial(Fraction(Integer))
+    cohenOperator? cn =>
+      t1:=retract(cn)::CohenOperator
+      t2:=arity(t1)
+      -- we have not yet defined prefix operators
+      prefix? t1 => 1::Polynomial(Fraction(Integer)) 
+      infix? t1 => 
+        name(t1) = _+ => formPoly(ck.1)+formPoly(ck.2)
+        name(t1) = _- => formPoly(ck.1)-formPoly(ck.2)
+        name(t1) = _* => formPoly(ck.1)*formPoly(ck.2)
+        name(t1) = _/ => 
+          numer:=getData(ck.1)
+          num:Union(CohenInteger,"failed"):=retractIfCan(numer)
+          num case "failed" => error "Only integer division supported"
+          numi:=(retract(numer)::CohenInteger)::Integer
+          denom:CohenNode:=getData(ck.2)
+          den:Union(CohenInteger,"failed"):=retractIfCan(denom)
+          den case "failed" => error "Only integer division supported"
+          deni:=(retract(denom)::CohenInteger)::Integer
+          (numi/deni)::Polynomial(Fraction(Integer)) 
+        name(t1) = _^ => 
+          base:=getData(ck.1)
+          power:CohenNode:=getData(ck.2)
+          pow:=(retract(power)::CohenInteger)::Integer
+          bas:=((retract(base)::CohenIdentifier)::Symbol)
+          expon(bas,pow)
+      -- we have not yet defined postfix operators
+      postfix? t1 => 64::Polynomial(Fraction(Integer)) 
+      128::Polynomial(Fraction(Integer)) 
+    cohenFunction? cn => 
+      t1:=retract(cn)::CohenFunction
+      t2:=arity(t1)
+      256::Polynomial(Fraction(Integer)) 
+    512::Polynomial(Fraction(Integer)) 
+
+  cohen2poly(cr:CohenRepresentation):Polynomial(Fraction(Integer)) ==
+    formPoly(cr)
+
+\end{chunk}
+
+\begin{verbatim}
+<<domain CEXPR CohenExpression>>=
+)abbrev domain CEXPR CohenExpression
+-- Author: Timothy Daly (based on Joel Cohen's books)
+-- Date Created: 7 April 2007
+-- Date Last Updated: 
+-- Basic Functions:
+-- Related Constructors: Fraction
+-- Also See:
+-- AMS Classifications:
+-- Keywords:
+-- References:
+-- "Computer Algebra and Symbolic Computation: Mathematical Methods"
+-- A.K. Peters, Ltd. ISBN 1-56881-159-4
+-- "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+-- A.K. Peters, Ltd. ISBN 1-56881-158-6
+-- Description:
+-- This type represents symbolic algebra as detailed in Joel Cohen's
+-- books on Computer Algebra and Symbolic Computation
+--
+-- Elements of this domain are purely symbolic expressions represented
+-- in a tree structure. All operations on the elements are purely
+-- procedural so expressions are not simplified unless the simplification
+-- is requested. Thus, members of this domain can be used in program to
+-- explain the steps of a calculation, such as subtracting 1 from both
+-- sides of an equation and displaying the unsimplified result:\br
+-- \tab{5}x + 1 - 1 = 5 - 1
+
+CohenExpression(CNODE:CohenNode): SetCategory with
+
+
+ == add 
+
+  Rep:= Tree CNODE
+
+  coerce(term:%):OutputForm ==
+    term.coef = 0 => 0::CINT()::OutputForm    
+    typeof(term) = inttag => term.int::OutputForm
+    typeof(term) = symtag => term.sym::OutputForm
+    "failed"::OutputForm
+
+\end{verbatim}
+
+\begin{chunk}{package CMANIP CohenManipulations}
+)abbrev package CMANIP CohenManipulations
+
+++ Author: Timothy Daly (based on Joel Cohen's books)
+++ Date Created: 7 April 2007
+++ Date Last Updated: 
+++ Basic Functions:
+++ Related Constructors: Fraction
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ "Computer Algebra and Symbolic Computation: Mathematical Methods"
+++ A.K. Peters, Ltd. ISBN 1-56881-159-4
+++ "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+++ A.K. Peters, Ltd. ISBN 1-56881-158-6
+++ Description:
+++ This type represents symbolic algebra as detailed in Joel Cohen's
+++ books on Computer Algebra and Symbolic Computation
+++
+++ Perform manipulations like constant folding
+
+CohenManipulations(): SetCategory with
+
+  constant? : CohenRepresentation -> Boolean
+
+  monomials : CohenRepresentation -> List CohenRepresentation
+
+  sum : List(CohenRepresentation) -> CohenRepresentation
+  sum : (List(CohenRepresentation),String) -> CohenRepresentation
+
+  constantFold : List(CohenRepresentation) -> List(CohenRepresentation)
+
+ == add 
+
+  constant?(cr:CohenRepresentation):Boolean ==
+    maybe:Union(CohenInteger,"failed"):=retractIfCan(getData(cr))
+    maybe case "failed" => false
+    true
+
+  monomials(cr:CohenRepresentation):List(CohenRepresentation) ==
+    cn:CohenNode:=getData(cr)
+    ck:List(CohenRepresentation):=getKids(cr)
+    cohenOperator? cn =>
+      op:CohenOperator:=retract(cn)
+      m:Symbol:=coerce(op)
+      m = '_+ => 
+        res:List(CohenRepresentation):=[]
+        for k in ck repeat res:=append(res,monomials k)
+        res
+      [cr]
+    [cr]
+
+  sum(cl:List(CohenRepresentation)):CohenRepresentation ==
+    sum(cl,"")
+
+  sum(cl:List(CohenRepresentation),msg:String):CohenRepresentation ==
+    messagePrint(msg)$OutputForm
+    folded:=constantFold(cl)
+    result:CohenRepresentation:=first folded
+    for monom in rest(folded) repeat                 -- remake the polynomial
+      result:=makeNode(opnode(copInfix('_+)),[result,monom])
+    result
+
+  constantFold(cl:List(CohenRepresentation)):List(CohenRepresentation) ==
+    result:List(CohenRepresentation):=[]
+    constantsum:Integer:=0
+    for term in cl repeat                            -- sum all the constants
+      if constant? term 
+      then
+        cohenint:CohenInteger:=retract(getData(term))
+        constantsum:=constantsum+coerce(cohenint)::Integer
+      else
+        result:=append(result,[term])
+    if constantsum > 0             -- make a nonzero constant into a monomial
+     then 
+       cint:CohenInteger:=coerce(constantsum)
+       result:=append(result,[makeNode(intnode(cint))])
+    #result = 0 =>                              -- everything collapsed to 0?
+      zero:Integer:=0
+      czero:CohenInteger:=coerce(zero)
+      [makeNode(intnode(czero))]
+    result
+
+\end{chunk}
+
+\begin{verbatim}
+)abbrev domain CPGM CohenProgram
+-- Author: Timothy Daly (based on Joel Cohen's books)
+-- Date Created: 7 April 2007
+-- Date Last Updated: 
+-- Basic Functions:
+-- Related Constructors: Fraction
+-- Also See:
+-- AMS Classifications:
+-- Keywords:
+-- References:
+-- "Computer Algebra and Symbolic Computation: Mathematical Methods"
+-- A.K. Peters, Ltd. ISBN 1-56881-159-4
+-- "Computer Algebra and Symbolic Computation: Elementary Algorithms"
+-- A.K. Peters, Ltd. ISBN 1-56881-158-6
+-- Description:
+-- This type represents symbolic algebra as detailed in Joel Cohen's
+-- books on Computer Algebra and Symbolic Computation
+--
+-- Elements of this domain are purely symbolic expressions represented
+-- in a tree structure. All operations on the elements are purely
+-- procedural so expressions are not simplified unless the simplification
+-- is requested. Thus, members of this domain can be used in program to
+-- explain the steps of a calculation, such as subtracting 1 from both
+-- sides of an equation and displaying the unsimplified result:\br
+-- \tab{5}x + 1 - 1 = 5 - 1
+
+CohenProgram(CNODE:CohenNode): SetCategory with
+
+
+ == add 
+
+  Rep:= Tree CNODE
+
+  coerce(term:%):OutputForm ==
+    term.coef = 0 => 0::CINT()::OutputForm    
+    typeof(term) = inttag => term.int::OutputForm
+    typeof(term) = symtag => term.sym::OutputForm
+    "failed"::OutputForm
+
+\end{verbatim}
+
+\begin{chunk}{cohen.spad}
+\getchunk{domain CINT CohenInteger}
+\getchunk{domain CID CohenIdentifier}
+\getchunk{domain COP CohenOperator}
+\getchunk{domain CFUNC CohenFunction}
+\getchunk{domain CNODE CohenNode}
+\getchunk{domain CREP CohenRepresentation}
+\getchunk{domain CPOLY CohenPolynomial}
+\getchunk{package CMANIP CohenManipulations}
+\end{chunk}
+
+\begin{chunk}{*}
+)set break resume
+)spool cohen.output
+)set message test on
+)set message auto off
+)clear all
+)sys cp $AXIOM/../../src/input/cohen.input.pamphlet .
+)lisp (tangle "cohen.input.pamphlet" "cohen.spad" "cohen.spad")
+
+-- CohenInteger
+--S 1 of 116
+)co cohen
+--E 1
+
+--S 2 of 116
+t1:CohenInteger:=1
+--R 
+--R
+--R   (1)  1
+--R                                                           Type: CohenInteger
+--E 2
+
+--S 3 of 116
+t2:Integer:=t1
+--R 
+--R
+--R   (2)  1
+--R                                                                Type: Integer
+--E 3
+
+)clear all
+ 
+-- CohenIdentifier
+
+--S 4 of 116
+t1:CohenIdentifier:=x
+--R 
+--R
+--R   (1)  x
+--R                                                        Type: CohenIdentifier
+--E 4
+
+--S 5 of 116
+t2:Symbol:=t1 
+--R 
+--R
+--R   (2)  x
+--R                                                                 Type: Symbol
+--E 5
+
+--S 6 of 116
+t3:=set(t1,5)
+--R 
+--R
+--R   (3)  x
+--R                                                        Type: CohenIdentifier
+--E 6
+
+--S 7 of 116
+t4:=get(t3)
+--R 
+--R
+--R   (4)  5
+--R                                                        Type: PositiveInteger
+--E 7
+
+--S 8 of 116
+t5:=name(t3)
+--R 
+--R
+--R   (5)  x
+--R                                                                 Type: Symbol
+--E 8
+
+)clear all
+ 
+--S 9 of 116
+t1:=copPrefix(x)
+--R 
+--R
+--R   (1)  x
+--R                                                          Type: CohenOperator
+--E 9
+
+--S 10 of 116
+t2:=copPrefix(x,1)
+--R 
+--R
+--R   (2)  x
+--R                                                          Type: CohenOperator
+--E 10
+
+--S 11 of 116
+t3:=copInfix(y)
+--R 
+--R
+--R   (3)  y
+--R                                                          Type: CohenOperator
+--E 11
+
+--S 12 of 116
+t4:=copInfix(y,1)
+--R 
+--R
+--R   (4)  y
+--R                                                          Type: CohenOperator
+--E 12
+
+--S 13 of 116
+t5:=copPostfix(z)
+--R 
+--R
+--R   (5)  z
+--R                                                          Type: CohenOperator
+--E 13
+
+--S 14 of 116
+t6:=copPostfix(z,3)
+--R 
+--R
+--R   (6)  z
+--R                                                          Type: CohenOperator
+--E 14
+
+--S 15 of 116
+t7:=name(t1)
+--R 
+--R
+--R   (7)  x
+--R                                                                 Type: Symbol
+--E 15
+
+--S 16 of 116
+t8:=arity(t6)
+--R 
+--R
+--R   (8)  3
+--R                                                        Type: PositiveInteger
+--E 16
+
+--S 17 of 116
+t9:=prefix?(t1)
+--R 
+--R
+--R   (9)  true
+--R                                                                Type: Boolean
+--E 17
+
+--S 18 of 116
+t10:=prefix?(t3)
+--R 
+--R
+--R   (10)  false
+--R                                                                Type: Boolean
+--E 18
+
+--S 19 of 116
+t11:=infix?(t1)
+--R 
+--R
+--R   (11)  false
+--R                                                                Type: Boolean
+--E 19
+
+--S 20 of 116
+t12:=infix?(t3)
+--R 
+--R
+--R   (12)  true
+--R                                                                Type: Boolean
+--E 20
+
+--S 21 of 116
+t13:=postfix?(t1)
+--R 
+--R
+--R   (13)  false
+--R                                                                Type: Boolean
+--E 21
+
+--S 22 of 116
+t14:=postfix?(t5)
+--R 
+--R
+--R   (14)  true
+--R                                                                Type: Boolean
+--E 22
+
+--S 23 of 116
+t15:Symbol:=t1
+--R 
+--R
+--R   (15)  x
+--R                                                                 Type: Symbol
+--E 23
+
+--S 24 of 116
+t16:CohenOperator:=x
+--R 
+--R
+--R   (16)  x
+--R                                                          Type: CohenOperator
+--E 24
+
+)clear all
+ 
+-- CohenFunction
+
+--S 25 of 116
+t1:=cFunction(x)
+--R 
+--R
+--R   (1)  x
+--R                                                          Type: CohenFunction
+--E 25
+
+--S 26 of 116
+t2:=cFunction(x,3)
+--R 
+--R
+--R   (2)  x
+--R                                                          Type: CohenFunction
+--E 26
+
+--S 27 of 116
+t3:=name(t2)
+--R 
+--R
+--R   (3)  x
+--R                                                                 Type: Symbol
+--E 27
+
+--S 28 of 116
+t4:=arity(t2)
+--R 
+--R
+--R   (4)  3
+--R                                                        Type: PositiveInteger
+--E 28
+
+--S 29 of 116
+t5:Symbol:=t1
+--R 
+--R
+--R   (5)  x
+--R                                                                 Type: Symbol
+--E 29
+
+--S 30 of 116
+t6:CohenFunction:=y
+--R 
+--R
+--R   (6)  y
+--R                                                          Type: CohenFunction
+--E 30
+
+)clear all
+ 
+-- CohenNode
+
+--S 31 of 116
+t1:=intnode(1)
+--R 
+--R
+--R   (1)  1
+--R                                                              Type: CohenNode
+--E 31
+
+--S 32 of 116
+t2:CohenNode:=coerce(1)
+--R 
+--R
+--R   (2)  1
+--R                                                              Type: CohenNode
+--E 32
+
+--S 33 of 116
+t3:CohenInteger:=retract(t1)
+--R 
+--R
+--R   (3)  1
+--R                                                           Type: CohenInteger
+--E 33
+
+--S 34 of 116
+t4:CohenInteger:=retractIfCan(t1)
+--R 
+--R
+--R   (4)  1
+--R                                                           Type: CohenInteger
+--E 34
+
+--S 35 of 116
+t5:=cohenInteger?(t1)
+--R 
+--R
+--R   (5)  true
+--R                                                                Type: Boolean
+--E 35
+
+--S 36 of 116
+t6:=idnode(x)
+--R 
+--R
+--R   (6)  x
+--R                                                              Type: CohenNode
+--E 36
+
+--S 37 of 116
+t7:CohenNode:=coerce(x)
+--R 
+--R
+--R   (7)  x
+--R                                                              Type: CohenNode
+--E 37
+
+--S 38 of 116
+t8:CohenIdentifier:=retract(t6)
+--R 
+--R
+--R   (8)  x
+--R                                                        Type: CohenIdentifier
+--E 38
+
+--S 39 of 116
+t9:CohenIdentifier:=retractIfCan(t6)
+--R 
+--R
+--R   (9)  x
+--R                                                        Type: CohenIdentifier
+--E 39
+
+--S 40 of 116
+t10:=cohenIdentifier?(t6)
+--R 
+--R
+--R   (10)  true
+--R                                                                Type: Boolean
+--E 40
+
+--S 41 of 116
+t11:=opnode(_+)
+--R 
+--R
+--R   (11)  +
+--R                                                              Type: CohenNode
+--E 41
+
+--S 42 of 116
+t12:CohenOperator:=coerce(_+)
+--R 
+--R
+--R   (12)  +
+--R                                                          Type: CohenOperator
+--E 42
+
+--S 43 of 116
+t13:CohenOperator:=retract(t11)
+--R 
+--R
+--R   (13)  +
+--R                                                          Type: CohenOperator
+--E 43
+
+--S 44 of 116
+t14:CohenOperator:=retractIfCan(t11)
+--R 
+--R
+--R   (14)  +
+--R                                                          Type: CohenOperator
+--E 44
+
+--S 45 of 116
+t15:=cohenOperator?(t11)
+--R 
+--R
+--R   (15)  true
+--R                                                                Type: Boolean
+--E 45
+
+--S 46 of 116
+t16:=fnnode(f)
+--R 
+--R
+--R   (16)  f
+--R                                                              Type: CohenNode
+--E 46
+
+--S 47 of 116
+t17:CohenFunction:=coerce(f)
+--R 
+--R
+--R   (17)  f
+--R                                                          Type: CohenFunction
+--E 47
+
+--S 48 of 116
+t18:CohenFunction:=retract(t16)
+--R 
+--R
+--R   (18)  f
+--R                                                          Type: CohenFunction
+--E 48
+
+--S 49 of 116
+t19:CohenFunction:=retractIfCan(t16)
+--R 
+--R
+--R   (19)  f
+--R                                                          Type: CohenFunction
+--E 49
+
+--S 50 of 116
+t20:=cohenFunction?(t16)
+--R 
+--R
+--R   (20)  true
+--R                                                                Type: Boolean
+--E 50
+
+)clear all
+ 
+-- CohenRepresentation
+
+--S 51 of 116
+t0:=empty()
+--R 
+--R
+--R   (1)  ()
+--R                                                    Type: CohenRepresentation
+--E 51
+
+--S 52 of 116
+explain? t0
+--R 
+--R
+--R   (2)  false
+--R                                                                Type: Boolean
+--E 52
+
+--S 53 of 116
+explain(t0,true)
+--R 
+--R
+--R   (3)  true
+--R                                                                Type: Boolean
+--E 53
+
+--S 54 of 116
+explain? t0
+--R 
+--R
+--R   (4)  true
+--R                                                                Type: Boolean
+--E 54
+
+--S 55 of 116
+t1:=makeNode(1::CINT) --> 1
+--R 
+--R
+--R   (5)  1
+--R                                                    Type: CohenRepresentation
+--E 55
+
+--S 56 of 116
+t2:=makeNode(2::CINT) --> 2
+--R 
+--R
+--R   (6)  2
+--R                                                    Type: CohenRepresentation
+--E 56
+
+--S 57 of 116
+t3:=makeNode(opnode(copInfix(-)),[t1,t2]) --> (1 - 2)
+--R 
+--R
+--R   (7)  1 - 2
+--R                                                    Type: CohenRepresentation
+--E 57
+
+--S 58 of 116
+t4:=makeNode(opnode(copPrefix(_+)),[t1,t2]) --> (+ 1 2)
+--R 
+--R
+--R   (8)  + 1 2
+--R                                                    Type: CohenRepresentation
+--E 58
+
+--S 59 of 116
+t5:=makeNode(opnode(copPrefix(_+)),[t3,t1]) --> (+ (1 - 2) 1)
+--R 
+--R
+--R   (9)  + 1 - 2 1
+--R                                                    Type: CohenRepresentation
+--E 59
+
+--S 60 of 116
+t5:=makeNode(opnode(copPostfix(_+)),[t1,t2]) --> (1 2 +)
+--R 
+--R
+--R   (10)  1 2 +
+--R                                                    Type: CohenRepresentation
+--E 60
+
+--S 61 of 116
+t6:=makeNode(opnode(copPostfix(_+)),[t3,t1]) --> ((1 - 2) 1 +)
+--R 
+--R
+--R   (11)  1 - 2 1 +
+--R                                                    Type: CohenRepresentation
+--E 61
+
+--S 62 of 116
+t7:=makeNode(opnode(copPostfix(_+)),[t4,t1]) --> ((+ 1 2) 1 +)
+--R 
+--R
+--R   (12)  + 1 2 1 +
+--R                                                    Type: CohenRepresentation
+--E 62
+
+--S 63 of 116
+t8:=makeNode(3::CINT) --> 3
+--R 
+--R
+--R   (13)  3
+--R                                                    Type: CohenRepresentation
+--E 63
+
+--S 64 of 116
+t9:=makeNode(opnode(copInfix(_+)),[t1,t2,t8]) --> (1 + 2 + 3)
+--R 
+--R
+--R   (14)  1 + 2 + 3
+--R                                                    Type: CohenRepresentation
+--E 64
+
+--S 65 of 116
+n1:=makeNode(x::CID) --> x
+--R 
+--R
+--R   (15)  x
+--R                                                    Type: CohenRepresentation
+--E 65
+
+--S 66 of 116
+n2:=makeNode(2::CINT) --> 2
+--R 
+--R
+--R   (16)  2
+--R                                                    Type: CohenRepresentation
+--E 66
+
+--S 67 of 116
+n3:=makeNode(opnode(copInfix(_^)),[n1,n2]) --> (x ^ 2)
+--R 
+--R
+--R   (17)  x ^ 2
+--R                                                    Type: CohenRepresentation
+--E 67
+
+--S 68 of 116
+n4:=makeNode(d::CID) --> d
+--R 
+--R
+--R   (18)  d
+--R                                                    Type: CohenRepresentation
+--E 68
+
+--S 69 of 116
+n5:=makeNode(opnode(copInfix(_*)),[n4,n3]) --> (d * (x ^ 2))
+--R 
+--R
+--R   (19)  d * x ^ 2
+--R                                                    Type: CohenRepresentation
+--E 69
+
+--S 70 of 116
+n6:=makeNode(c::CID) --> c
+--R 
+--R
+--R   (20)  c
+--R                                                    Type: CohenRepresentation
+--E 70
+
+--S 71 of 116
+n7:=makeNode(opnode(copInfix(_+)),[n6,n5]) --> (c + (d * (x ^ 2)))
+--R 
+--R
+--R   (21)  c + d * x ^ 2
+--R                                                    Type: CohenRepresentation
+--E 71
+
+--S 72 of 116
+m1:=makeNode(opnode(copPostfix(_!)),[t2]) --> (2 !)
+--R 
+--R
+--R   (22)  2 !
+--R                                                    Type: CohenRepresentation
+--E 72
+
+--S 73 of 116
+p1:=makeNode(x::CID) --> x
+--R 
+--R
+--R   (23)  x
+--R                                                    Type: CohenRepresentation
+--E 73
+
+--S 74 of 116
+p2:=makeNode(y::CID) --> y
+--R 
+--R
+--R   (24)  y
+--R                                                    Type: CohenRepresentation
+--E 74
+
+--S 75 of 116
+p3:=makeNode(z::CID) --> z
+--R 
+--R
+--R   (25)  z
+--R                                                    Type: CohenRepresentation
+--E 75
+
+--S 76 of 116
+p4:=makeNode(fnnode(cFunction(f)),[p1,p2,p3]) --> f(x,y,z)
+--R 
+--R
+--R   (26)  f(x,y,z)
+--R                                                    Type: CohenRepresentation
+--E 76
+
+--S 77 of 116
+q1:=makeNode(subscript(a,[1])::CID) --> a_1
+--R 
+--R
+--R   (27)  a
+--R          1
+--R                                                    Type: CohenRepresentation
+--E 77
+
+--S 78 of 116
+q2:=makeNode(subscript(a,[2])::CID) --> a_2
+--R 
+--R
+--R   (28)  a
+--R          2
+--R                                                    Type: CohenRepresentation
+--E 78
+
+--S 79 of 116
+q3:=makeNode(subscript(a,[3])::CID) --> a_3
+--R 
+--R
+--R   (29)  a
+--R          3
+--R                                                    Type: CohenRepresentation
+--E 79
+
+--S 80 of 116
+q4:=makeNode(fnnode(cFunction(f)),[q1,q2,q3]) --> f(a ,a ,a )
+--R 
+--R
+--R   (30)  f(a ,a ,a )
+--R            1  2  3
+--R                                                    Type: CohenRepresentation
+--E 80
+
+)clear all
+ 
+-- CohenPolynomial
+
+--S 81 of 116
+cohen2poly makeNode(intnode(1::CINT))
+--R 
+--R
+--R   (1)  1
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 81
+
+--S 82 of 116
+cohen2poly makeNode(idnode(x::CID))
+--R 
+--R
+--R   (2)  x
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 82
+
+--S 83 of 116
+poly:=1/3*x^4+1/2*x^2+15*x+18
+--R 
+--R
+--R        1  4   1  2
+--R   (3)  - x  + - x  + 15x + 18
+--R        3      2
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 83
+
+--S 84 of 116
+t1:=poly2cohen(poly)$CPOLY
+--R 
+--R
+--R   (4)  1 / 3 * x ^ 4 + 1 / 2 * x ^ 2 + 15 * x ^ 1 + 18
+--R                                                    Type: CohenRepresentation
+--E 84
+
+--S 85 of 116
+t2:=cohen2poly t1
+--R 
+--R
+--R        1  4   1  2
+--R   (5)  - x  + - x  + 15x + 18
+--R        3      2
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 85
+
+--S 86 of 116
+tt1:=makeNode(2::CINT)
+--R 
+--R
+--R   (6)  2
+--R                                                    Type: CohenRepresentation
+--E 86
+
+--S 87 of 116
+tt2:=makeNode(4::CINT)
+--R 
+--R
+--R   (7)  4
+--R                                                    Type: CohenRepresentation
+--E 87
+
+--S 88 of 116
+tt3:=makeNode(opnode(copInfix(_/)),[tt1,tt2])
+--R 
+--R
+--R   (8)  2 / 4
+--R                                                    Type: CohenRepresentation
+--E 88
+
+--S 89 of 116
+tt4:=cohen2poly tt3
+--R 
+--R
+--R        1
+--R   (9)  -
+--R        2
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 89
+
+--S 90 of 116
+poly:=1/3*x^4+1/2*x^2+15*x+18
+--R 
+--R
+--R         1  4   1  2
+--R   (10)  - x  + - x  + 15x + 18
+--R         3      2
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 90
+
+--S 91 of 116
+poly+poly
+--R 
+--R
+--R         2  4    2
+--R   (11)  - x  + x  + 30x + 36
+--R         3
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 91
+
+--S 92 of 116
+p1:=poly2cohen poly
+--R 
+--R
+--R   (12)  1 / 3 * x ^ 4 + 1 / 2 * x ^ 2 + 15 * x ^ 1 + 18
+--R                                                    Type: CohenRepresentation
+--E 92
+
+--S 93 of 116
+t1:=p1+p1
+--R 
+--R
+--R   (13)
+--R   1 / 3 * x ^ 4 + 1 / 2 * x ^ 2 + 15 * x ^ 1 + 18 + 1 / 3 * x ^ 4 + 1 / 2 * x
+--R     ^ 2 + 15 * x ^ 1 + 18
+--R                                                    Type: CohenRepresentation
+--E 93
+
+--S 94 of 116
+cohen2poly t1
+--R 
+--R
+--R         2  4    2
+--R   (14)  - x  + x  + 30x + 36
+--R         3
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 94
+
+--S 95 of 116
+parens(t1,true)
+--R 
+--R
+--R   (15)  true
+--R                                                                Type: Boolean
+--E 95
+
+--S 96 of 116
+t1
+--R 
+--R
+--R   (16)
+--R   ((((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18)))  +
+--R    (((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18))))
+--R                                                    Type: CohenRepresentation
+--E 96
+
+--S 97 of 116
+t2:=p1-p1
+--R 
+--R
+--R   (17)
+--R   ((((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18)))  -
+--R    (((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18))))
+--R                                                    Type: CohenRepresentation
+--E 97
+
+--S 98 of 116
+cohen2poly t2
+--R 
+--R
+--R   (18)  0
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 98
+
+--S 99 of 116
+t3:=p1*p1
+--R 
+--R
+--R   (19)
+--R   ((((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18)))  *
+--R    (((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18))))
+--R                                                    Type: CohenRepresentation
+--E 99
+
+--S 100 of 116
+cohen2poly t3
+--R 
+--R
+--R         1  8   1  6      5   49  4      3       2
+--R   (20)  - x  + - x  + 10x  + -- x  + 15x  + 243x  + 540x + 324
+--R         9      3              4
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 100
+
+)clear all
+ 
+-- CohenManipulations
+)clear all
+ 
+--S 101 of 116
+t1:=x+1
+--R 
+--R
+--R   (1)  x + 1
+--R                                                    Type: Polynomial(Integer)
+--E 101
+
+--S 102 of 116
+t2:=poly2cohen t1
+--R 
+--R
+--R   (2)  ((1 * (x ^ 1)) + 1)
+--R                                                    Type: CohenRepresentation
+--E 102
+
+--S 103 of 116
+parens(t2,false)
+--R 
+--R
+--R   (3)  false
+--R                                                                Type: Boolean
+--E 103
+
+--S 104 of 116
+t3:=t2+makeNode((-1)::CINT)
+--R 
+--R
+--R   (4)  1 * x ^ 1 + 1 + - 1
+--R                                                    Type: CohenRepresentation
+--E 104
+
+--S 105 of 116
+t4:=sum(monomials(t3))
+--R 
+--R
+--R   (5)  1 * x ^ 1
+--R                                                    Type: CohenRepresentation
+--E 105
+
+--S 106 of 116
+parens(t4,true)
+--R 
+--R
+--R   (6)  true
+--R                                                                Type: Boolean
+--E 106
+
+--S 107 of 116
+t3:=t2+makeNode((-1)::CINT)
+--R 
+--R
+--R   (7)  (((1 * (x ^ 1)) + 1) + - 1)
+--R                                                    Type: CohenRepresentation
+--E 107
+
+--S 108 of 116
+poly:=1/3*x^4+1/2*x^2+15*x+18
+--R 
+--R
+--R        1  4   1  2
+--R   (8)  - x  + - x  + 15x + 18
+--R        3      2
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 108
+
+--S 109 of 116
+p1:=poly2cohen poly
+--R 
+--R
+--R   (9)  (((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18)))
+--R                                                    Type: CohenRepresentation
+--E 109
+
+--S 110 of 116
+t1:=p1+makeNode(6::CINT)
+--R 
+--R
+--R   (10)
+--R   ((((1 / 3) * (x ^ 4)) + (((1 / 2) * (x ^ 2)) + ((15 * (x ^ 1)) + 18))) + 6)
+--R                                                    Type: CohenRepresentation
+--E 110
+
+--S 111 of 116
+cohen2poly t1
+--R 
+--R
+--R         1  4   1  2
+--R   (11)  - x  + - x  + 15x + 24
+--R         3      2
+--R                                          Type: Polynomial(Fraction(Integer))
+--E 111
+
+--S 112 of 116
+t2:=monomials t1
+--R 
+--R
+--R   (12)  [((1 / 3) * (x ^ 4)),((1 / 2) * (x ^ 2)),(15 * (x ^ 1)),18,6]
+--R                                              Type: List(CohenRepresentation)
+--E 112
+
+--S 113 of 116
+constantFold monomials t1
+--R 
+--R
+--R   (13)  [((1 / 3) * (x ^ 4)),((1 / 2) * (x ^ 2)),(15 * (x ^ 1)),24]
+--R                                              Type: List(CohenRepresentation)
+--E 113
+
+--S 114 of 116
+map(x+->constant? x,t2)
+--R 
+--R
+--R   (14)  [false,false,false,true,true]
+--R                                                          Type: List(Boolean)
+--E 114
+
+--S 115 of 116
+sum(t2)
+--R 
+--R
+--R   (15)  (((((1 / 3) * (x ^ 4)) + ((1 / 2) * (x ^ 2))) + (15 * (x ^ 1))) + 24)
+--R                                                    Type: CohenRepresentation
+--E 115
+
+--S 116 of 116
+sum(t2,"Now we sum all the constant terms")
+--R 
+--R   Now we sum all the constant terms
+--R
+--R   (16)  (((((1 / 3) * (x ^ 4)) + ((1 / 2) * (x ^ 2))) + (15 * (x ^ 1))) + 24)
+--R                                                    Type: CohenRepresentation
+--E 116
+
+)spool
+)lisp (bye)
+
+\end{chunk}
+
+\newpage
+\begin{thebibliography}{99}
+\bibitem{1} nothing
+\end{thebibliography}
+\end{document}
