Initial population.

1 files changed, 914 insertions, 0 deletions

diff --git a/src/algebra/elemntry.spad.pamphlet b/src/algebra/elemntry.spad.pamphlet
new file mode 100644
index 00000000..8ddc8e52
--- /dev/null
+++ b/src/algebra/elemntry.spad.pamphlet
@@ -0,0 +1,914 @@
+\documentclass{article}
+\usepackage{axiom}
+\begin{document}
+\title{\$SPAD/src/algebra elemntry.spad}
+\author{Manuel Bronstein}
+\maketitle
+\begin{abstract}
+\end{abstract}
+\eject
+\tableofcontents
+\eject
+\section{package EF ElementaryFunction}
+<<package EF ElementaryFunction>>=
+)abbrev package EF ElementaryFunction
+++ Author: Manuel Bronstein
+++ Date Created: 1987
+++ Date Last Updated: 10 April 1995
+++ Keywords: elementary, function, logarithm, exponential.
+++ Examples:  )r EF INPUT
+++ Description: Provides elementary functions over an integral domain.
+ElementaryFunction(R, F): Exports == Implementation where
+  R: Join(OrderedSet, IntegralDomain)
+  F: Join(FunctionSpace R, RadicalCategory)
+
+  B   ==> Boolean
+  L   ==> List
+  Z   ==> Integer
+  OP  ==> BasicOperator
+  K   ==> Kernel F
+  INV ==> error "Invalid argument"
+
+  Exports ==> with
+    exp     : F -> F
+	++ exp(x) applies the exponential operator to x
+    log     : F -> F
+	++ log(x) applies the logarithm operator to x
+    sin     : F -> F
+	++ sin(x) applies the sine operator to x
+    cos     : F -> F
+	++ cos(x) applies the cosine operator to x 
+    tan     : F -> F
+	++ tan(x) applies the tangent operator to x
+    cot     : F -> F
+	++ cot(x) applies the cotangent operator to x
+    sec     : F -> F
+	++ sec(x) applies the secant operator to x
+    csc     : F -> F
+	++ csc(x) applies the cosecant operator to x
+    asin    : F -> F
+	++ asin(x) applies the inverse sine operator to x 
+    acos    : F -> F
+	++ acos(x) applies the inverse cosine operator to x
+    atan    : F -> F
+	++ atan(x) applies the inverse tangent operator to x
+    acot    : F -> F
+	++ acot(x) applies the inverse cotangent operator to x
+    asec    : F -> F
+	++ asec(x) applies the inverse secant operator to x
+    acsc    : F -> F
+	++ acsc(x) applies the inverse cosecant operator to x
+    sinh    : F -> F
+	++ sinh(x) applies the hyperbolic sine operator to x 
+    cosh    : F -> F
+	++ cosh(x) applies the hyperbolic cosine operator to x
+    tanh    : F -> F
+	++ tanh(x) applies the hyperbolic tangent operator to x
+    coth    : F -> F
+	++ coth(x) applies the hyperbolic cotangent operator to x
+    sech    : F -> F
+	++ sech(x) applies the hyperbolic secant operator to x
+    csch    : F -> F
+	++ csch(x) applies the hyperbolic cosecant operator to x
+    asinh   : F -> F
+	++ asinh(x) applies the inverse hyperbolic sine operator to x
+    acosh   : F -> F
+	++ acosh(x) applies the inverse hyperbolic cosine operator to x
+    atanh   : F -> F
+	++ atanh(x) applies the inverse hyperbolic tangent operator to x
+    acoth   : F -> F
+	++ acoth(x) applies the inverse hyperbolic cotangent operator to x
+    asech   : F -> F
+	++ asech(x) applies the	inverse hyperbolic secant operator to x
+    acsch   : F -> F
+	++ acsch(x) applies the inverse hyperbolic cosecant operator to x
+    pi      : () -> F
+	++ pi() returns the pi operator
+    belong? : OP -> Boolean
+	++ belong?(p) returns true if operator p is elementary
+    operator: OP -> OP
+	++ operator(p) returns an elementary operator with the same symbol as p
+    -- the following should be local, but are conditional
+    iisqrt2   : () -> F
+	++ iisqrt2() should be local but conditional
+    iisqrt3   : () -> F
+	++ iisqrt3() should be local but conditional
+    iiexp     : F -> F
+	++ iiexp(x) should be local but conditional
+    iilog     : F -> F
+	++ iilog(x) should be local but conditional
+    iisin     : F -> F
+	++ iisin(x) should be local but conditional
+    iicos     : F -> F
+	++ iicos(x) should be local but conditional
+    iitan     : F -> F
+	++ iitan(x) should be local but conditional
+    iicot     : F -> F
+	++ iicot(x) should be local but conditional
+    iisec     : F -> F
+	++ iisec(x) should be local but conditional
+    iicsc     : F -> F
+	++ iicsc(x) should be local but conditional
+    iiasin    : F -> F
+	++ iiasin(x) should be local but conditional
+    iiacos    : F -> F
+	++ iiacos(x) should be local but conditional
+    iiatan    : F -> F
+	++ iiatan(x) should be local but conditional
+    iiacot    : F -> F
+	++ iiacot(x) should be local but conditional
+    iiasec    : F -> F
+	++ iiasec(x) should be local but conditional
+    iiacsc    : F -> F
+	++ iiacsc(x) should be local but conditional
+    iisinh    : F -> F
+	++ iisinh(x) should be local but conditional
+    iicosh    : F -> F
+	++ iicosh(x) should be local but conditional
+    iitanh    : F -> F
+	++ iitanh(x) should be local but conditional
+    iicoth    : F -> F
+	++ iicoth(x) should be local but conditional
+    iisech    : F -> F
+	++ iisech(x) should be local but conditional
+    iicsch    : F -> F
+	++ iicsch(x) should be local but conditional
+    iiasinh   : F -> F
+	++ iiasinh(x) should be local but conditional
+    iiacosh   : F -> F
+	++ iiacosh(x) should be local but conditional
+    iiatanh   : F -> F
+	++ iiatanh(x) should be local but conditional
+    iiacoth   : F -> F
+	++ iiacoth(x) should be local but conditional
+    iiasech   : F -> F
+	++ iiasech(x) should be local but conditional
+    iiacsch   : F -> F
+	++ iiacsch(x) should be local but conditional
+    specialTrigs:(F, L Record(func:F,pole:B)) -> Union(F, "failed")
+	++ specialTrigs(x,l) should be local but conditional
+    localReal?: F -> Boolean
+	++ localReal?(x) should be local but conditional
+
+  Implementation ==> add
+    ipi      : List F -> F
+    iexp     : F -> F
+    ilog     : F -> F
+    iiilog   : F -> F
+    isin     : F -> F
+    icos     : F -> F
+    itan     : F -> F
+    icot     : F -> F
+    isec     : F -> F
+    icsc     : F -> F
+    iasin    : F -> F
+    iacos    : F -> F
+    iatan    : F -> F
+    iacot    : F -> F
+    iasec    : F -> F
+    iacsc    : F -> F
+    isinh    : F -> F
+    icosh    : F -> F
+    itanh    : F -> F
+    icoth    : F -> F
+    isech    : F -> F
+    icsch    : F -> F
+    iasinh   : F -> F
+    iacosh   : F -> F
+    iatanh   : F -> F
+    iacoth   : F -> F
+    iasech   : F -> F
+    iacsch   : F -> F
+    dropfun  : F -> F
+    kernel   : F -> K
+    posrem   :(Z, Z) -> Z
+    iisqrt1  : () -> F
+    valueOrPole : Record(func:F, pole:B) -> F
+
+    oppi  := operator("pi"::Symbol)$CommonOperators
+    oplog := operator("log"::Symbol)$CommonOperators
+    opexp := operator("exp"::Symbol)$CommonOperators
+    opsin := operator("sin"::Symbol)$CommonOperators
+    opcos := operator("cos"::Symbol)$CommonOperators
+    optan := operator("tan"::Symbol)$CommonOperators
+    opcot := operator("cot"::Symbol)$CommonOperators
+    opsec := operator("sec"::Symbol)$CommonOperators
+    opcsc := operator("csc"::Symbol)$CommonOperators
+    opasin := operator("asin"::Symbol)$CommonOperators
+    opacos := operator("acos"::Symbol)$CommonOperators
+    opatan := operator("atan"::Symbol)$CommonOperators
+    opacot := operator("acot"::Symbol)$CommonOperators
+    opasec := operator("asec"::Symbol)$CommonOperators
+    opacsc := operator("acsc"::Symbol)$CommonOperators
+    opsinh := operator("sinh"::Symbol)$CommonOperators
+    opcosh := operator("cosh"::Symbol)$CommonOperators
+    optanh := operator("tanh"::Symbol)$CommonOperators
+    opcoth := operator("coth"::Symbol)$CommonOperators
+    opsech := operator("sech"::Symbol)$CommonOperators
+    opcsch := operator("csch"::Symbol)$CommonOperators
+    opasinh := operator("asinh"::Symbol)$CommonOperators
+    opacosh := operator("acosh"::Symbol)$CommonOperators
+    opatanh := operator("atanh"::Symbol)$CommonOperators
+    opacoth := operator("acoth"::Symbol)$CommonOperators
+    opasech := operator("asech"::Symbol)$CommonOperators
+    opacsch := operator("acsch"::Symbol)$CommonOperators
+
+    -- Pi is a domain...
+    Pie, isqrt1, isqrt2, isqrt3: F
+
+    -- following code is conditionalized on arbitraryPrecesion to recompute in
+    -- case user changes the precision
+
+    if R has TranscendentalFunctionCategory then
+      Pie := pi()$R :: F
+    else
+      Pie := kernel(oppi, nil()$List(F))
+
+    if R has TranscendentalFunctionCategory and R has arbitraryPrecision then
+      pi() == pi()$R :: F
+    else
+      pi() == Pie
+
+    if R has imaginary: () -> R then
+      isqrt1 := imaginary()$R :: F
+    else isqrt1 := sqrt(-1::F)
+
+    if R has RadicalCategory then
+      isqrt2 := sqrt(2::R)::F
+      isqrt3 := sqrt(3::R)::F
+    else
+      isqrt2 := sqrt(2::F)
+      isqrt3 := sqrt(3::F)
+
+    iisqrt1() == isqrt1
+    if R has RadicalCategory and R has arbitraryPrecision then
+      iisqrt2() == sqrt(2::R)::F
+      iisqrt3() == sqrt(3::R)::F
+    else
+      iisqrt2() == isqrt2
+      iisqrt3() == isqrt3
+
+    ipi l == pi()
+    log x == oplog x
+    exp x == opexp x
+    sin x == opsin x
+    cos x == opcos x
+    tan x == optan x
+    cot x == opcot x
+    sec x == opsec x
+    csc x == opcsc x
+    asin x == opasin x
+    acos x == opacos x
+    atan x == opatan x
+    acot x == opacot x
+    asec x == opasec x
+    acsc x == opacsc x
+    sinh x == opsinh x
+    cosh x == opcosh x
+    tanh x == optanh x
+    coth x == opcoth x
+    sech x == opsech x
+    csch x == opcsch x
+    asinh x == opasinh x
+    acosh x == opacosh x
+    atanh x == opatanh x
+    acoth x == opacoth x
+    asech x == opasech x
+    acsch x == opacsch x
+    kernel x == retract(x)@K
+
+    posrem(n, m)    == ((r := n rem m) < 0 => r + m; r)
+    valueOrPole rec == (rec.pole => INV; rec.func)
+    belong? op      == has?(op, "elem")
+
+    operator op ==
+      is?(op, "pi"::Symbol)    => oppi
+      is?(op, "log"::Symbol)   => oplog
+      is?(op, "exp"::Symbol)   => opexp
+      is?(op, "sin"::Symbol)   => opsin
+      is?(op, "cos"::Symbol)   => opcos
+      is?(op, "tan"::Symbol)   => optan
+      is?(op, "cot"::Symbol)   => opcot
+      is?(op, "sec"::Symbol)   => opsec
+      is?(op, "csc"::Symbol)   => opcsc
+      is?(op, "asin"::Symbol)  => opasin
+      is?(op, "acos"::Symbol)  => opacos
+      is?(op, "atan"::Symbol)  => opatan
+      is?(op, "acot"::Symbol)  => opacot
+      is?(op, "asec"::Symbol)  => opasec
+      is?(op, "acsc"::Symbol)  => opacsc
+      is?(op, "sinh"::Symbol)  => opsinh
+      is?(op, "cosh"::Symbol)  => opcosh
+      is?(op, "tanh"::Symbol)  => optanh
+      is?(op, "coth"::Symbol)  => opcoth
+      is?(op, "sech"::Symbol)  => opsech
+      is?(op, "csch"::Symbol)  => opcsch
+      is?(op, "asinh"::Symbol) => opasinh
+      is?(op, "acosh"::Symbol) => opacosh
+      is?(op, "atanh"::Symbol) => opatanh
+      is?(op, "acoth"::Symbol) => opacoth
+      is?(op, "asech"::Symbol) => opasech
+      is?(op, "acsch"::Symbol) => opacsch
+      error "Not an elementary operator"
+
+    dropfun x ==
+      ((k := retractIfCan(x)@Union(K, "failed")) case "failed") or
+        empty?(argument(k::K)) => 0
+      first argument(k::K)
+
+    if R has RetractableTo Z then
+      specialTrigs(x, values) ==
+        (r := retractIfCan(y := x/pi())@Union(Fraction Z, "failed"))
+          case "failed" => "failed"
+        q := r::Fraction(Integer)
+        m := minIndex values
+        (n := retractIfCan(q)@Union(Z, "failed")) case Z =>
+          even?(n::Z) => valueOrPole(values.m)
+          valueOrPole(values.(m+1))
+        (n := retractIfCan(2*q)@Union(Z, "failed")) case Z =>
+--          one?(s := posrem(n::Z, 4)) => valueOrPole(values.(m+2))
+          (s := posrem(n::Z, 4)) = 1 => valueOrPole(values.(m+2))
+          valueOrPole(values.(m+3))
+        (n := retractIfCan(3*q)@Union(Z, "failed")) case Z =>
+--          one?(s := posrem(n::Z, 6)) => valueOrPole(values.(m+4))
+          (s := posrem(n::Z, 6)) = 1 => valueOrPole(values.(m+4))
+          s = 2 => valueOrPole(values.(m+5))
+          s = 4 => valueOrPole(values.(m+6))
+          valueOrPole(values.(m+7))
+        (n := retractIfCan(4*q)@Union(Z, "failed")) case Z =>
+--          one?(s := posrem(n::Z, 8)) => valueOrPole(values.(m+8))
+          (s := posrem(n::Z, 8)) = 1 => valueOrPole(values.(m+8))
+          s = 3 => valueOrPole(values.(m+9))
+          s = 5 => valueOrPole(values.(m+10))
+          valueOrPole(values.(m+11))
+        (n := retractIfCan(6*q)@Union(Z, "failed")) case Z =>
+--          one?(s := posrem(n::Z, 12)) => valueOrPole(values.(m+12))
+          (s := posrem(n::Z, 12)) = 1 => valueOrPole(values.(m+12))
+          s = 5 => valueOrPole(values.(m+13))
+          s = 7 => valueOrPole(values.(m+14))
+          valueOrPole(values.(m+15))
+        "failed"
+
+    else specialTrigs(x, values) == "failed"
+
+    isin x ==
+      zero? x => 0
+      y := dropfun x
+      is?(x, opasin) => y
+      is?(x, opacos) => sqrt(1 - y**2)
+      is?(x, opatan) => y / sqrt(1 + y**2)
+      is?(x, opacot) => inv sqrt(1 + y**2)
+      is?(x, opasec) => sqrt(y**2 - 1) / y
+      is?(x, opacsc) => inv y
+      h  := inv(2::F)
+      s2 := h * iisqrt2()
+      s3 := h * iisqrt3()
+      u  := specialTrigs(x, [[0,false], [0,false], [1,false], [-1,false],
+                         [s3,false], [s3,false], [-s3,false], [-s3,false],
+                          [s2,false], [s2,false], [-s2,false], [-s2,false],
+                           [h,false], [h,false], [-h,false], [-h,false]])
+      u case F => u :: F
+      kernel(opsin, x)
+
+    icos x ==
+      zero? x => 1
+      y := dropfun x
+      is?(x, opasin) => sqrt(1 - y**2)
+      is?(x, opacos) => y
+      is?(x, opatan) => inv sqrt(1 + y**2)
+      is?(x, opacot) => y / sqrt(1 + y**2)
+      is?(x, opasec) => inv y
+      is?(x, opacsc) => sqrt(y**2 - 1) / y
+      h  := inv(2::F)
+      s2 := h * iisqrt2()
+      s3 := h * iisqrt3()
+      u  := specialTrigs(x, [[1,false],[-1,false], [0,false], [0,false],
+                             [h,false],[-h,false],[-h,false],[h,false],
+                              [s2,false],[-s2,false],[-s2,false],[s2,false],
+                               [s3,false], [-s3,false],[-s3,false],[s3,false]])
+      u case F => u :: F
+      kernel(opcos, x)
+
+    itan x ==
+      zero? x => 0
+      y := dropfun x
+      is?(x, opasin) => y / sqrt(1 - y**2)
+      is?(x, opacos) => sqrt(1 - y**2) / y
+      is?(x, opatan) => y
+      is?(x, opacot) => inv y
+      is?(x, opasec) => sqrt(y**2 - 1)
+      is?(x, opacsc) => inv sqrt(y**2 - 1)
+      s33 := (s3 := iisqrt3()) / (3::F)
+      u := specialTrigs(x, [[0,false], [0,false], [0,true], [0,true],
+                      [s3,false], [-s3,false], [s3,false], [-s3,false],
+                       [1,false], [-1,false], [1,false], [-1,false],
+                        [s33,false], [-s33, false], [s33,false], [-s33, false]])
+      u case F => u :: F
+      kernel(optan, x)
+
+    icot x ==
+      zero? x => INV
+      y := dropfun x
+      is?(x, opasin) => sqrt(1 - y**2) / y
+      is?(x, opacos) => y / sqrt(1 - y**2)
+      is?(x, opatan) => inv y
+      is?(x, opacot) => y
+      is?(x, opasec) => inv sqrt(y**2 - 1)
+      is?(x, opacsc) => sqrt(y**2 - 1)
+      s33 := (s3 := iisqrt3()) / (3::F)
+      u := specialTrigs(x, [[0,true], [0,true], [0,false], [0,false],
+                         [s33,false], [-s33,false], [s33,false], [-s33,false],
+                          [1,false], [-1,false], [1,false], [-1,false],
+                           [s3,false], [-s3, false], [s3,false], [-s3, false]])
+      u case F => u :: F
+      kernel(opcot, x)
+
+    isec x ==
+      zero? x => 1
+      y := dropfun x
+      is?(x, opasin) => inv sqrt(1 - y**2)
+      is?(x, opacos) => inv y
+      is?(x, opatan) => sqrt(1 + y**2)
+      is?(x, opacot) => sqrt(1 + y**2) / y
+      is?(x, opasec) => y
+      is?(x, opacsc) => y / sqrt(y**2 - 1)
+      s2 := iisqrt2()
+      s3 := 2 * iisqrt3() / (3::F)
+      h  := 2::F
+      u  := specialTrigs(x, [[1,false],[-1,false],[0,true],[0,true],
+                           [h,false], [-h,false], [-h,false], [h,false],
+                            [s2,false], [-s2,false], [-s2,false], [s2,false],
+                             [s3,false], [-s3,false], [-s3,false], [s3,false]])
+      u case F => u :: F
+      kernel(opsec, x)
+
+    icsc x ==
+      zero? x => INV
+      y := dropfun x
+      is?(x, opasin) => inv y
+      is?(x, opacos) => inv sqrt(1 - y**2)
+      is?(x, opatan) => sqrt(1 + y**2) / y
+      is?(x, opacot) => sqrt(1 + y**2)
+      is?(x, opasec) => y / sqrt(y**2 - 1)
+      is?(x, opacsc) => y
+      s2 := iisqrt2()
+      s3 := 2 * iisqrt3() / (3::F)
+      h  := 2::F
+      u  := specialTrigs(x, [[0,true], [0,true], [1,false], [-1,false],
+                            [s3,false], [s3,false], [-s3,false], [-s3,false],
+                              [s2,false], [s2,false], [-s2,false], [-s2,false],
+                                 [h,false], [h,false], [-h,false], [-h,false]])
+      u case F => u :: F
+      kernel(opcsc, x)
+
+    iasin x ==
+      zero? x => 0
+--      one? x =>   pi() / (2::F)
+      (x = 1) =>   pi() / (2::F)
+      x = -1 => - pi() / (2::F)
+      y := dropfun x
+      is?(x, opsin) => y
+      is?(x, opcos) => pi() / (2::F) - y
+      kernel(opasin, x)
+
+    iacos x ==
+      zero? x => pi() / (2::F)
+--      one? x => 0
+      (x = 1) => 0
+      x = -1 => pi()
+      y := dropfun x
+      is?(x, opsin) => pi() / (2::F) - y
+      is?(x, opcos) => y
+      kernel(opacos, x)
+
+    iatan x ==
+      zero? x => 0
+--      one? x =>   pi() / (4::F)
+      (x = 1) =>   pi() / (4::F)
+      x = -1 => - pi() / (4::F)
+      x = (r3:=iisqrt3()) => pi() / (3::F)
+--      one?(x*r3)          => pi() / (6::F)
+      (x*r3) = 1          => pi() / (6::F)
+      y := dropfun x
+      is?(x, optan) => y
+      is?(x, opcot) => pi() / (2::F) - y
+      kernel(opatan, x)
+
+    iacot x ==
+      zero? x =>   pi() / (2::F)
+--      one? x  =>   pi() / (4::F)
+      (x = 1)  =>   pi() / (4::F)
+      x = -1  =>   3 * pi() / (4::F)
+      x = (r3:=iisqrt3())  =>  pi() / (6::F)
+      x = -r3              =>  5 * pi() / (6::F)
+--      one?(xx:=x*r3)       =>  pi() / (3::F)
+      (xx:=x*r3) = 1      =>  pi() / (3::F)
+      xx = -1           =>     2* pi() / (3::F)
+      y := dropfun x
+      is?(x, optan) => pi() / (2::F) - y
+      is?(x, opcot) => y
+      kernel(opacot, x)
+
+    iasec x ==
+      zero? x => INV
+--      one? x => 0
+      (x = 1) => 0
+      x = -1 => pi()
+      y := dropfun x
+      is?(x, opsec) => y
+      is?(x, opcsc) => pi() / (2::F) - y
+      kernel(opasec, x)
+
+    iacsc x ==
+      zero? x => INV
+--      one? x =>   pi() / (2::F)
+      (x = 1) =>   pi() / (2::F)
+      x = -1 => - pi() / (2::F)
+      y := dropfun x
+      is?(x, opsec) => pi() / (2::F) - y
+      is?(x, opcsc) => y
+      kernel(opacsc, x)
+
+    isinh x ==
+      zero? x => 0
+      y := dropfun x
+      is?(x, opasinh) => y
+      is?(x, opacosh) => sqrt(y**2 - 1)
+      is?(x, opatanh) => y / sqrt(1 - y**2)
+      is?(x, opacoth) => - inv sqrt(y**2 - 1)
+      is?(x, opasech) => sqrt(1 - y**2) / y
+      is?(x, opacsch) => inv y
+      kernel(opsinh, x)
+
+    icosh x ==
+      zero? x => 1
+      y := dropfun x
+      is?(x, opasinh) => sqrt(y**2 + 1)
+      is?(x, opacosh) => y
+      is?(x, opatanh) => inv sqrt(1 - y**2)
+      is?(x, opacoth) => y / sqrt(y**2 - 1)
+      is?(x, opasech) => inv y
+      is?(x, opacsch) => sqrt(y**2 + 1) / y
+      kernel(opcosh, x)
+
+    itanh x ==
+      zero? x => 0
+      y := dropfun x
+      is?(x, opasinh) => y / sqrt(y**2 + 1)
+      is?(x, opacosh) => sqrt(y**2 - 1) / y
+      is?(x, opatanh) => y
+      is?(x, opacoth) => inv y
+      is?(x, opasech) => sqrt(1 - y**2)
+      is?(x, opacsch) => inv sqrt(y**2 + 1)
+      kernel(optanh, x)
+
+    icoth x ==
+      zero? x => INV
+      y := dropfun x
+      is?(x, opasinh) => sqrt(y**2 + 1) / y
+      is?(x, opacosh) => y / sqrt(y**2 - 1)
+      is?(x, opatanh) => inv y
+      is?(x, opacoth) => y
+      is?(x, opasech) => inv sqrt(1 - y**2)
+      is?(x, opacsch) => sqrt(y**2 + 1)
+      kernel(opcoth, x)
+
+    isech x ==
+      zero? x => 1
+      y := dropfun x
+      is?(x, opasinh) => inv sqrt(y**2 + 1)
+      is?(x, opacosh) => inv y
+      is?(x, opatanh) => sqrt(1 - y**2)
+      is?(x, opacoth) => sqrt(y**2 - 1) / y
+      is?(x, opasech) => y
+      is?(x, opacsch) => y / sqrt(y**2 + 1)
+      kernel(opsech, x)
+
+    icsch x ==
+      zero? x => INV
+      y := dropfun x
+      is?(x, opasinh) => inv y
+      is?(x, opacosh) => inv sqrt(y**2 - 1)
+      is?(x, opatanh) => sqrt(1 - y**2) / y
+      is?(x, opacoth) => - sqrt(y**2 - 1)
+      is?(x, opasech) => y / sqrt(1 - y**2)
+      is?(x, opacsch) => y
+      kernel(opcsch, x)
+
+    iasinh x ==
+      is?(x, opsinh) => first argument kernel x
+      kernel(opasinh, x)
+
+    iacosh x ==
+      is?(x, opcosh) => first argument kernel x
+      kernel(opacosh, x)
+
+    iatanh x ==
+      is?(x, optanh) => first argument kernel x
+      kernel(opatanh, x)
+
+    iacoth x ==
+      is?(x, opcoth) => first argument kernel x
+      kernel(opacoth, x)
+
+    iasech x ==
+      is?(x, opsech) => first argument kernel x
+      kernel(opasech, x)
+
+    iacsch x ==
+      is?(x, opcsch) => first argument kernel x
+      kernel(opacsch, x)
+
+    iexp x ==
+      zero? x => 1
+      is?(x, oplog) => first argument kernel x
+      x < 0 and empty? variables x => inv iexp(-x)
+      h  := inv(2::F)
+      i  := iisqrt1()
+      s2 := h * iisqrt2()
+      s3 := h * iisqrt3()
+      u  := specialTrigs(x / i, [[1,false],[-1,false], [i,false], [-i,false],
+            [h + i * s3,false], [-h + i * s3, false], [-h - i * s3, false],
+             [h - i * s3, false], [s2 + i * s2, false], [-s2 + i * s2, false],
+              [-s2 - i * s2, false], [s2 - i * s2, false], [s3 + i * h, false],
+               [-s3 + i * h, false], [-s3 - i * h, false], [s3 - i * h, false]])
+      u case F => u :: F
+      kernel(opexp, x)
+
+-- THIS DETERMINES WHEN TO PERFORM THE log exp f -> f SIMPLIFICATION
+-- CURRENT BEHAVIOR:
+--     IF R IS COMPLEX(S) THEN ONLY ELEMENTS WHICH ARE RETRACTABLE TO R
+--     AND EQUAL TO THEIR CONJUGATES ARE DEEMED REAL (OVERRESTRICTIVE FOR NOW)
+--     OTHERWISE (e.g. R = INT OR FRAC INT), ALL THE ELEMENTS ARE DEEMED REAL
+
+    if (R has imaginary:() -> R) and (R has conjugate: R -> R) then
+         localReal? x ==
+            (u := retractIfCan(x)@Union(R, "failed")) case R
+               and (u::R) = conjugate(u::R)
+
+    else localReal? x == true
+
+    iiilog x ==
+      zero? x => INV
+--      one? x => 0
+      (x = 1) => 0
+      (u := isExpt(x, opexp)) case Record(var:K, exponent:Integer) =>
+           rec := u::Record(var:K, exponent:Integer)
+           arg := first argument(rec.var);
+           localReal? arg => rec.exponent * first argument(rec.var);
+           ilog x
+      ilog x
+
+    ilog x ==
+--      ((num1 := one?(num := numer x)) or num = -1) and (den := denom x) ^= 1
+      ((num1 := ((num := numer x) = 1)) or num = -1) and (den := denom x) ^= 1
+        and empty? variables x => - kernel(oplog, (num1 => den; -den)::F)
+      kernel(oplog, x)
+
+    if R has ElementaryFunctionCategory then
+      iilog x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iiilog x
+        log(r::R)::F
+
+      iiexp x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iexp x
+        exp(r::R)::F
+
+    else
+      iilog x == iiilog x
+      iiexp x == iexp x
+
+    if R has TrigonometricFunctionCategory then
+      iisin x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => isin x
+        sin(r::R)::F
+
+      iicos x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => icos x
+        cos(r::R)::F
+
+      iitan x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => itan x
+        tan(r::R)::F
+
+      iicot x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => icot x
+        cot(r::R)::F
+
+      iisec x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => isec x
+        sec(r::R)::F
+
+      iicsc x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => icsc x
+        csc(r::R)::F
+
+    else
+      iisin x == isin x
+      iicos x == icos x
+      iitan x == itan x
+      iicot x == icot x
+      iisec x == isec x
+      iicsc x == icsc x
+
+    if R has ArcTrigonometricFunctionCategory then
+      iiasin x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iasin x
+        asin(r::R)::F
+
+      iiacos x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iacos x
+        acos(r::R)::F
+
+      iiatan x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iatan x
+        atan(r::R)::F
+
+      iiacot x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iacot x
+        acot(r::R)::F
+
+      iiasec x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iasec x
+        asec(r::R)::F
+
+      iiacsc x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iacsc x
+        acsc(r::R)::F
+
+    else
+      iiasin x == iasin x
+      iiacos x == iacos x
+      iiatan x == iatan x
+      iiacot x == iacot x
+      iiasec x == iasec x
+      iiacsc x == iacsc x
+
+    if R has HyperbolicFunctionCategory then
+      iisinh x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => isinh x
+        sinh(r::R)::F
+
+      iicosh x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => icosh x
+        cosh(r::R)::F
+
+      iitanh x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => itanh x
+        tanh(r::R)::F
+
+      iicoth x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => icoth x
+        coth(r::R)::F
+
+      iisech x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => isech x
+        sech(r::R)::F
+
+      iicsch x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => icsch x
+        csch(r::R)::F
+
+    else
+      iisinh x == isinh x
+      iicosh x == icosh x
+      iitanh x == itanh x
+      iicoth x == icoth x
+      iisech x == isech x
+      iicsch x == icsch x
+
+    if R has ArcHyperbolicFunctionCategory then
+      iiasinh x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iasinh x
+        asinh(r::R)::F
+
+      iiacosh x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iacosh x
+        acosh(r::R)::F
+
+      iiatanh x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iatanh x
+        atanh(r::R)::F
+
+      iiacoth x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iacoth x
+        acoth(r::R)::F
+
+      iiasech x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iasech x
+        asech(r::R)::F
+
+      iiacsch x ==
+        (r:=retractIfCan(x)@Union(R,"failed")) case "failed" => iacsch x
+        acsch(r::R)::F
+
+    else
+      iiasinh x == iasinh x
+      iiacosh x == iacosh x
+      iiatanh x == iatanh x
+      iiacoth x == iacoth x
+      iiasech x == iasech x
+      iiacsch x == iacsch x
+
+    evaluate(oppi, ipi)$BasicOperatorFunctions1(F)
+    evaluate(oplog, iilog)
+    evaluate(opexp, iiexp)
+    evaluate(opsin, iisin)
+    evaluate(opcos, iicos)
+    evaluate(optan, iitan)
+    evaluate(opcot, iicot)
+    evaluate(opsec, iisec)
+    evaluate(opcsc, iicsc)
+    evaluate(opasin, iiasin)
+    evaluate(opacos, iiacos)
+    evaluate(opatan, iiatan)
+    evaluate(opacot, iiacot)
+    evaluate(opasec, iiasec)
+    evaluate(opacsc, iiacsc)
+    evaluate(opsinh, iisinh)
+    evaluate(opcosh, iicosh)
+    evaluate(optanh, iitanh)
+    evaluate(opcoth, iicoth)
+    evaluate(opsech, iisech)
+    evaluate(opcsch, iicsch)
+    evaluate(opasinh, iiasinh)
+    evaluate(opacosh, iiacosh)
+    evaluate(opatanh, iiatanh)
+    evaluate(opacoth, iiacoth)
+    evaluate(opasech, iiasech)
+    evaluate(opacsch, iiacsch)
+    derivative(opexp, exp)
+    derivative(oplog, inv)
+    derivative(opsin, cos)
+    derivative(opcos, - sin #1)
+    derivative(optan, 1 + tan(#1)**2)
+    derivative(opcot, - 1 - cot(#1)**2)
+    derivative(opsec, tan(#1) * sec(#1))
+    derivative(opcsc, - cot(#1) * csc(#1))
+    derivative(opasin, inv sqrt(1 - #1**2))
+    derivative(opacos, - inv sqrt(1 - #1**2))
+    derivative(opatan, inv(1 + #1**2))
+    derivative(opacot, - inv(1 + #1**2))
+    derivative(opasec, inv(#1 * sqrt(#1**2 - 1)))
+    derivative(opacsc, - inv(#1 * sqrt(#1**2 - 1)))
+    derivative(opsinh, cosh)
+    derivative(opcosh, sinh)
+    derivative(optanh, 1 - tanh(#1)**2)
+    derivative(opcoth, 1 - coth(#1)**2)
+    derivative(opsech, - tanh(#1) * sech(#1))
+    derivative(opcsch, - coth(#1) * csch(#1))
+    derivative(opasinh, inv sqrt(1 + #1**2))
+    derivative(opacosh, inv sqrt(#1**2 - 1))
+    derivative(opatanh, inv(1 - #1**2))
+    derivative(opacoth, inv(1 - #1**2))
+    derivative(opasech, - inv(#1 * sqrt(1 - #1**2)))
+    derivative(opacsch, - inv(#1 * sqrt(1 + #1**2)))
+
+@
+\section{License}
+<<license>>=
+--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
+--All rights reserved.
+--
+--Redistribution and use in source and binary forms, with or without
+--modification, are permitted provided that the following conditions are
+--met:
+--
+--    - Redistributions of source code must retain the above copyright
+--      notice, this list of conditions and the following disclaimer.
+--
+--    - Redistributions in binary form must reproduce the above copyright
+--      notice, this list of conditions and the following disclaimer in
+--      the documentation and/or other materials provided with the
+--      distribution.
+--
+--    - Neither the name of The Numerical ALgorithms Group Ltd. nor the
+--      names of its contributors may be used to endorse or promote products
+--      derived from this software without specific prior written permission.
+--
+--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
+--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+@
+<<*>>=
+<<license>>
+
+-- SPAD files for the functional world should be compiled in the
+-- following order:
+--
+--   op  kl  fspace  algfunc  ELEMNTRY  expr
+<<package EF ElementaryFunction>>
+@
+\eject
+\begin{thebibliography}{99}
+\bibitem{1} nothing
+\end{thebibliography}
+\end{document}