(FILECREATED "16-Sep-86 13:10:19" {ERIS}<LISPCORE>LIBRARY>CMLFLOAT.;39 81771  

      changes to:  (FNS %%EXPT-INTEGER)

      previous date: "20-Aug-86 20:27:43" {ERIS}<LISPCORE>LIBRARY>CMLFLOAT.;38)


(* Copyright (c) 1986 by Xerox Corporation. All rights reserved.)

(PRETTYCOMPRINT CMLFLOATCOMS)

(RPAQQ CMLFLOATCOMS 
       ((* * CMLFLOAT -- Covering sections 12.5-12.5.3 irrational, transcendental, exponential, 
           logarithmic, trigonometric, and hyperbolic functions. Section 12.10, implementation 
           parameters. -- By Kelly Roach. *)
        (DECLARE: EVAL@COMPILE DONTCOPY (FILES (LOADCOMP)
                                               LLFLOAT))
        (COMS (* Section 12.10, implementation parameters. The constants in this COMS are exported to 
                 the user. *)
              (* %%FLOAT allows us to recreate FLOATPs in a way that is independent of the ordinairy 
                 reading and printing FLOATPs to files which involves loss of the last couple bits of 
                 accuracy due to rounding effects. *)
              (* Using INITVARS instead of CONSTANTS in various places here because of problems with 
                 the way BYTECOMPILER stores FLOATPs in DCOM files. *)
              (FNS %%FLOAT)
              (CONSTANTS (MOST-POSITIVE-FIXNUM 65535)
                     (MOST-NEGATIVE-FIXNUM -65536))
              (INITVARS (MOST-POSITIVE-SINGLE-FLOAT (%%FLOAT 32639 65535))
                     (LEAST-POSITIVE-SINGLE-FLOAT (%%FLOAT 128 0))
                     (LEAST-NEGATIVE-SINGLE-FLOAT (%%FLOAT 32896 0))
                     (MOST-NEGATIVE-SINGLE-FLOAT (%%FLOAT 65407 65535))
                     (MOST-POSITIVE-SHORT-FLOAT MOST-POSITIVE-SINGLE-FLOAT)
                     (LEAST-POSITIVE-SHORT-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)
                     (LEAST-NEGATIVE-SHORT-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)
                     (MOST-NEGATIVE-SHORT-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)
                     (MOST-POSITIVE-DOUBLE-FLOAT MOST-POSITIVE-SINGLE-FLOAT)
                     (LEAST-POSITIVE-DOUBLE-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)
                     (LEAST-NEGATIVE-DOUBLE-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)
                     (MOST-NEGATIVE-DOUBLE-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)
                     (MOST-POSITIVE-LONG-FLOAT MOST-POSITIVE-SINGLE-FLOAT)
                     (LEAST-POSITIVE-LONG-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)
                     (LEAST-NEGATIVE-LONG-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)
                     (MOST-NEGATIVE-LONG-FLOAT MOST-NEGATIVE-SINGLE-FLOAT))
              (* EPSILON is the smallest positive floating point number such that
                 (NOT (= (FLOAT 1 EPSILON)
                         (+ (FLOAT 1 EPSILON)
                            EPSILON)))
                 *)
              (INITVARS (SINGLE-FLOAT-EPSILON (%%FLOAT 13312 0))
                     (SHORT-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)
                     (DOUBLE-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)
                     (LONG-FLOAT-EPSILON SINGLE-FLOAT-EPSILON))
              (* NEGATIVE-EPSILON is the smallest negative floating point number such that
                 (NOT (= (FLOAT 1 NEGATIVE-EPSILON)
                         (- (FLOAT 1 NEGATIVE-EPSILON)
                            NEGATIVE-EPSILON)))
                 *)
              (INITVARS (SINGLE-FLOAT-NEGATIVE-EPSILON (%%FLOAT 13184 0))
                     (SHORT-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)
                     (DOUBLE-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)
                     (LONG-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)))
        (COMS (* Miscellaneous. *)
              (DECLARE: EVAL@COMPILE DONTCOPY (FILES (LOADCOMP)
                                                     LLFLOAT))
              (INITVARS (PI (%%FLOAT 16457 4059)))
              (* Should be (DECLARE: EVAL@COMPILE DONTCOPY (CONSTANTS ...))
                 except that compiler does a poor job of compiling FLOATPs. Use an INITVARS to patch 
                 around this situation for now. *)
              (INITVARS (%%E (%%FLOAT 16429 63572))
                     (%%2PI (%%FLOAT 16585 4059))
                     (%%PI (%%FLOAT 16457 4059))
                     (%%2PI/3 (%%FLOAT 16390 2706))
                     (%%PI/2 (%%FLOAT 16329 4059))
                     (%%-PI/2 (%%FLOAT 49097 4059))
                     (%%PI/3 (%%FLOAT 16262 2706))
                     (%%PI/4 (%%FLOAT 16201 4059))
                     (%%-PI/4 (%%FLOAT 48969 4059))
                     (%%PI/6 (%%FLOAT 16134 2706))
                     (%%2/PI (%%FLOAT 16162 63875)))
              (FNS %%MAKE-ARRAY))
        (COMS (* EXP *)
              (COMS (INITVARS (%%LOG-BASE2-E (%%FLOAT 16312 43579)))
                    (* * %%EXP-POLY contains P and Q coefficients of Harris et al EXPB 1103 rational 
                       approximation to (EXPT 2 X)
                       in interval (0 .125)
                       %. %%EXP-TABLE contains values of powers (EXPT 2 (/ N 8)) . *)
                    (VARS (%%EXP-POLY (%%MAKE-ARRAY (LIST (%%FLOAT 15549 17659)
                                                          (%%FLOAT 16256 0)
                                                          (%%FLOAT 16801 38273)
                                                          (%%FLOAT 17257 7717)
                                                          (%%FLOAT 17597 11739)
                                                          (%%FLOAT 17800 30401))))
                          (%%EXP-TABLE (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                                           (%%FLOAT 16267 38338)
                                                           (%%FLOAT 16280 14320)
                                                           (%%FLOAT 16293 65239)
                                                           (%%FLOAT 16309 1267)
                                                           (%%FLOAT 16325 26410)
                                                           (%%FLOAT 16343 17661)
                                                           (%%FLOAT 16362 49351))))))
              (FNS EXP %%EXP-FLOAT %%EXP-COMPLEX))
        (COMS (* EXPT *)
              (FNS CL:EXPT %%EXPT-INTEGER %%EXPT-FLOAT %%EXPT-COMPLEX %%EXPT-COMPLEX-POWER))
        (COMS (* LOG *)
              (COMS (INITVARS (%%LOG2 (%%FLOAT 16177 29208))
                           (%%SQRT2 (%%FLOAT 16309 1267)))
                    (* * %%LOG-PPOLY and %%LOG-QPOLY contain P and Q coefficients of Harris et al 
                       LOGE 2707 rational approximation to (LOG X)
                       in interval ((SQRT .5)
                                    (SQRT 2)) . *)
                    (VARS (%%LOG-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16042 22803)
                                                           (%%FLOAT 49484 23590)
                                                           (%%FLOAT 17044 17982)
                                                           (%%FLOAT 49926 37153)
                                                           (%%FLOAT 17046 5367))))
                          (%%LOG-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                                           (%%FLOAT 49512 9103)
                                                           (%%FLOAT 16992 42274)
                                                           (%%FLOAT 49823 38048)
                                                           (%%FLOAT 16918 5367))))))
              (FNS CL:LOG %%LOG-FLOAT %%LOG-COMPLEX))
        (COMS (* SQRT *)
              (FNS CL:SQRT %%SQRT-FLOAT %%SQRT-COMPLEX))
        (COMS (* SIN *)
              (COMS (* * %%SIN-EPSILON is sufficiently small that (SIN X)
                       = X for X in interval (0 %%SIN-EPSILON)
                       %. It suffices to take %%SIN-EPSILON a little bit smaller than
                       (SQRT (CL:* 6 SINGLE-FLOAT-EPSILON))
                       which we get by the Taylor series expansion (SIN X)
                       =
                       (+ X (/ (EXPT X 3)
                               6)
                          ...)
                       (The relative error caused by ommitting (/ (EXPT X 3)
                                                                  6)
                            isn't observable.)
                       Comparison against %%SIN-EPSILON is used to avoid POLYEVAL microcode underflow 
                       when computing SIN. *)
                    (INITVARS (%%SIN-EPSILON (%%FLOAT 14720 0)))
                    (* * %%SIN-PPOLY and %%SIN-QPOLY contain adapted P and Q coefficients of Harris 
                       et al SIN 3374 rational approximation to (SIN X)
                       in interval (0 (/ PI 2))
                       %. The coefficients for %%SIN-PPOLY and %%SIN-QPOLY have been computed from 
                       Harris using extended precision routines and the relations %%SIN-PPOLY =
                       (REVERSE (for I from 0 as ENTRY in PS collect
                                     (/ (CL:* (EXPT (/ 2 PI)
                                                    (1+ (CL:* 2 I)))
                                              ENTRY)
                                        Q0)))
                       and %%SIN-QPOLY = (REVERSE (for I from 0 as ENTRY in QS collect
                                                       (/ (CL:* (EXPT (/ 2 PI)
                                                                      (CL:* 2 I))
                                                                ENTRY)
                                                          Q0)))
                       *)
                    (VARS (%%SIN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 45236 25611)
                                                           (%%FLOAT 13589 26148)
                                                           (%%FLOAT 47286 34797)
                                                           (%%FLOAT 15295 3306)
                                                           (%%FLOAT 48666 34805)
                                                           (%%FLOAT 16256 0))))
                          (%%SIN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 11384 52865)
                                                           (%%FLOAT 12553 9550)
                                                           (%%FLOAT 13604 38385)
                                                           (%%FLOAT 14593 18841)
                                                           (%%FLOAT 15489 5549)
                                                           (%%FLOAT 16256 0))))))
              (FNS CL:SIN %%SIN-FLOAT %%SIN-COMPLEX))
        (COMS (* COS *)
              (FNS CL:COS %%COS-COMPLEX))
        (COMS (* TAN *)
              (COMS (* * %%TAN-EPSILON is sufficiently small that (TAN X)
                       = X for X in interval (0 %%TAN-EPSILON)
                       %. It suffices to take %%TAN-EPSILON a little bit smaller than
                       (SQRT (CL:* 3 SINGLE-FLOAT-EPSILON))
                       which we get by the Taylor series expansion (TAN X)
                       =
                       (+ X (/ (EXPT X 3)
                               3)
                          ...)
                       (The relative error caused by ommitting (/ (EXPT X 3)
                                                                  3)
                            isn't observable.)
                       Comparison against %%TAN-EPSILON is used to avoid POLYEVAL microcode underflow 
                       when computing TAN. *)
                    (INITVARS (%%TAN-EPSILON (%%FLOAT 14720 0)))
                    (* * %%TAN-PPOLY and %%TAN-QPOLY contain adapted P and Q coefficients of Harris 
                       et al TAN 4288 rational approximation to (TAN X)
                       in interval (-PI/4 PI/4)
                       %. The coefficients for %%TAN-PPOLY and %%TAN-QPOLY have been computed from 
                       Harris using extended precision routines and the relations %%TAN-PPOLY =
                       (REVERSE (for I from 0 as ENTRY in PS collect
                                     (/ (CL:* (EXPT (/ 4 PI)
                                                    (1+ (CL:* 2 I)))
                                              ENTRY)
                                        Q0)))
                       and %%TAN-QPOLY = (REVERSE (for I from 0 as ENTRY in QS collect
                                                       (/ (CL:* (EXPT (/ 4 PI)
                                                                      (CL:* 2 I))
                                                                ENTRY)
                                                          Q0)))
                       *)
                    (VARS (%%TAN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 13237 21090)
                                                           (%%FLOAT 47141 15825)
                                                           (%%FLOAT 15246 8785)
                                                           (%%FLOAT 48655 48761)
                                                           (%%FLOAT 16256 0))))
                          (%%TAN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 45267 36947)
                                                           (%%FLOAT 13848 46875)
                                                           (%%FLOAT 47612 53738)
                                                           (%%FLOAT 15596 52854)
                                                           (%%FLOAT 48882 35303)
                                                           (%%FLOAT 16256 0))))))
              (FNS CL:TAN %%TAN-FLOAT %%TAN-COMPLEX))
        (COMS (* ASIN *)
              (COMS (* * %%ASIN-EPSILON is sufficiently small that (ASIN X)
                       = X for X in interval (0 %%ASIN-EPSILON)
                       %. It suffices to take %%ASIN-EPSILON a little bit smaller than
                       (CL:* 2 SINGLE-FLOAT-EPSILON)
                       which we get by the Taylor series expansion (ASIN X)
                       =
                       (+ X (/ (EXPT X 3)
                               6)
                          ...)
                       (The relative error caused by ommitting (/ (EXPT X 3)
                                                                  6)
                            isn't observable.)
                       Comparison against %%ASIN-EPSILON is used to avoid POLYEVAL microcode 
                       underflow when computing SIN. *)
                    (INITVARS (%%ASIN-EPSILON (%%FLOAT 14720 0)))
                    (* * %%ASIN-PPOLY and %%ASIN-QPOLY contain P and Q coefficients of Harris et al 
                       ARCSN 4671 rational approximation to (ASIN X)
                       in interval (0 (SQRT .5)) . *)
                    (VARS (%%ASIN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16007 50045)
                                                            (%%FLOAT 49549 8020)
                                                            (%%FLOAT 17236 15848)
                                                            (%%FLOAT 50285 63464)
                                                            (%%FLOAT 17650 31235)
                                                            (%%FLOAT 50403 62852)
                                                            (%%FLOAT 17440 39471))))
                          (%%ASIN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                                            (%%FLOAT 49672 25817)
                                                            (%%FLOAT 17308 55260)
                                                            (%%FLOAT 50326 38098)
                                                            (%%FLOAT 17674 22210)
                                                            (%%FLOAT 50417 22451)
                                                            (%%FLOAT 17440 39471))))))
              (FNS ASIN %%ASIN-FLOAT %%ASIN-COMPLEX))
        (COMS (* ACOS *)
              (FNS ACOS %%ACOS-COMPLEX))
        (COMS (* ATAN *)
              (INITVARS (%%SQRT3 (%%FLOAT 16349 46039))
                     (%%2-SQRT3 (%%FLOAT 16009 12451))
                     (%%INV-2-SQRT3 (%%FLOAT 16494 55788)))
              (FNS CL:ATAN %%ATAN-FLOAT1 %%ATAN-FLOAT2 %%ATAN-DOMAIN-CHECK %%ATAN-FLOAT 
                   %%ATAN-COMPLEX))
        (COMS (* CIS *)
              (FNS CIS))
        (COMS (* SINH COSH TANH *)
              (FNS SINH COSH TANH))
        (COMS (* ASINH ACOSH ATANH *)
              (FNS ASINH ACOSH ATANH %%ATANH-DOMAIN-CHECK))
        (PROP FILETYPE CMLFLOAT)
        (DECLARE: DONTEVAL@LOAD DOEVAL@COMPILE DONTCOPY COMPILERVARS
               (ADDVARS (NLAMA)
                      (NLAML)
                      (LAMA %%ATANH-DOMAIN-CHECK ATANH ACOSH ASINH TANH COSH SINH CIS %%ATAN-COMPLEX 
                            %%ATAN-FLOAT %%ATAN-DOMAIN-CHECK %%ATAN-FLOAT2 %%ATAN-FLOAT1 CL:ATAN 
                            %%ACOS-COMPLEX ACOS %%ASIN-COMPLEX ASIN %%TAN-COMPLEX CL:TAN 
                            %%COS-COMPLEX CL:COS %%SIN-COMPLEX CL:SIN CL:SQRT %%LOG-COMPLEX CL:LOG 
                            %%EXPT-COMPLEX-POWER %%EXPT-COMPLEX %%EXPT-INTEGER CL:EXPT %%EXP-COMPLEX 
                            %%EXP-FLOAT EXP)))))
(* * CMLFLOAT -- Covering sections %12.5-12.5.3 irrational, transcendental, exponential, 
logarithmic, trigonometric, and hyperbolic functions. Section 12.10, implementation parameters.
 -- By Kelly Roach. *)

(DECLARE: EVAL@COMPILE DONTCOPY 
(FILESLOAD (LOADCOMP)
       LLFLOAT)
)



(* Section 12.10, implementation parameters. The constants in this COMS are exported to the 
user. *)




(* %%FLOAT allows us to recreate FLOATPs in a way that is independent of the ordinairy reading 
and printing FLOATPs to files which involves loss of the last couple bits of accuracy due to 
rounding effects. *)




(* Using INITVARS instead of CONSTANTS in various places here because of problems with the way 
BYTECOMPILER stores FLOATPs in DCOM files. *)

(DEFINEQ

(%%FLOAT
  (LAMBDA (HIWORD LOWORD)                                    (* kbr: "14-May-86 16:43")
    (\FLOATBOX (\VAG2 HIWORD LOWORD))))
)
(DECLARE: EVAL@COMPILE 

(RPAQQ MOST-POSITIVE-FIXNUM 65535)

(RPAQQ MOST-NEGATIVE-FIXNUM -65536)

(CONSTANTS (MOST-POSITIVE-FIXNUM 65535)
       (MOST-NEGATIVE-FIXNUM -65536))
)

(RPAQ? MOST-POSITIVE-SINGLE-FLOAT (%%FLOAT 32639 65535))

(RPAQ? LEAST-POSITIVE-SINGLE-FLOAT (%%FLOAT 128 0))

(RPAQ? LEAST-NEGATIVE-SINGLE-FLOAT (%%FLOAT 32896 0))

(RPAQ? MOST-NEGATIVE-SINGLE-FLOAT (%%FLOAT 65407 65535))

(RPAQ? MOST-POSITIVE-SHORT-FLOAT MOST-POSITIVE-SINGLE-FLOAT)

(RPAQ? LEAST-POSITIVE-SHORT-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)

(RPAQ? LEAST-NEGATIVE-SHORT-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)

(RPAQ? MOST-NEGATIVE-SHORT-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)

(RPAQ? MOST-POSITIVE-DOUBLE-FLOAT MOST-POSITIVE-SINGLE-FLOAT)

(RPAQ? LEAST-POSITIVE-DOUBLE-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)

(RPAQ? LEAST-NEGATIVE-DOUBLE-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)

(RPAQ? MOST-NEGATIVE-DOUBLE-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)

(RPAQ? MOST-POSITIVE-LONG-FLOAT MOST-POSITIVE-SINGLE-FLOAT)

(RPAQ? LEAST-POSITIVE-LONG-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)

(RPAQ? LEAST-NEGATIVE-LONG-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)

(RPAQ? MOST-NEGATIVE-LONG-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)



(* EPSILON is the smallest positive floating point number such that (NOT (= (FLOAT 1 EPSILON) (
+ (FLOAT 1 EPSILON) EPSILON))) *)


(RPAQ? SINGLE-FLOAT-EPSILON (%%FLOAT 13312 0))

(RPAQ? SHORT-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)

(RPAQ? DOUBLE-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)

(RPAQ? LONG-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)



(* NEGATIVE-EPSILON is the smallest negative floating point number such that (NOT (= (FLOAT 1 
NEGATIVE-EPSILON) (- (FLOAT 1 NEGATIVE-EPSILON) NEGATIVE-EPSILON))) *)


(RPAQ? SINGLE-FLOAT-NEGATIVE-EPSILON (%%FLOAT 13184 0))

(RPAQ? SHORT-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)

(RPAQ? DOUBLE-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)

(RPAQ? LONG-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)



(* Miscellaneous. *)

(DECLARE: EVAL@COMPILE DONTCOPY 
(FILESLOAD (LOADCOMP)
       LLFLOAT)
)

(RPAQ? PI (%%FLOAT 16457 4059))



(* Should be (DECLARE: EVAL@COMPILE DONTCOPY (CONSTANTS ...)) except that compiler does a poor 
job of compiling FLOATPs. Use an INITVARS to patch around this situation for now. *)


(RPAQ? %%E (%%FLOAT 16429 63572))

(RPAQ? %%2PI (%%FLOAT 16585 4059))

(RPAQ? %%PI (%%FLOAT 16457 4059))

(RPAQ? %%2PI/3 (%%FLOAT 16390 2706))

(RPAQ? %%PI/2 (%%FLOAT 16329 4059))

(RPAQ? %%-PI/2 (%%FLOAT 49097 4059))

(RPAQ? %%PI/3 (%%FLOAT 16262 2706))

(RPAQ? %%PI/4 (%%FLOAT 16201 4059))

(RPAQ? %%-PI/4 (%%FLOAT 48969 4059))

(RPAQ? %%PI/6 (%%FLOAT 16134 2706))

(RPAQ? %%2/PI (%%FLOAT 16162 63875))
(DEFINEQ

(%%MAKE-ARRAY
  (LAMBDA (LIST)                                             (* kbr: "15-May-86 18:09")
                                                             (* Function to build Interlisp arrays.
                                                             I would prefer to build Common Lisp 
                                                             arrays, but I don't think POLYEVAL 
                                                             opcode will understand them.
                                                             *)
    (PROG (ARRAY)
          (SETQ ARRAY (ARRAY (LENGTH LIST)
                             (QUOTE FLOATP)
                             0.0 0))
          (for ELEMENT in LIST as I from 0 do (SETA ARRAY I ELEMENT))
          (RETURN ARRAY))))
)



(* EXP *)


(RPAQ? %%LOG-BASE2-E (%%FLOAT 16312 43579))
(* * %%EXP-POLY contains P and Q coefficients of Harris et al EXPB 1103 rational approximation 
to (EXPT 2 X) in interval (0 .125) %. %%EXP-TABLE contains values of powers (EXPT 2 (/ N 8)) . 
*)


(RPAQ %%EXP-POLY (%%MAKE-ARRAY (LIST (%%FLOAT 15549 17659)
                                     (%%FLOAT 16256 0)
                                     (%%FLOAT 16801 38273)
                                     (%%FLOAT 17257 7717)
                                     (%%FLOAT 17597 11739)
                                     (%%FLOAT 17800 30401))))

(RPAQ %%EXP-TABLE (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                      (%%FLOAT 16267 38338)
                                      (%%FLOAT 16280 14320)
                                      (%%FLOAT 16293 65239)
                                      (%%FLOAT 16309 1267)
                                      (%%FLOAT 16325 26410)
                                      (%%FLOAT 16343 17661)
                                      (%%FLOAT 16362 49351))))
(DEFINEQ

(EXP
  (CL:LAMBDA (NUMBER)                                        (* lmm "27-Jun-86 22:40")
                                                             (* Calculates e to the given power, 
                                                             where e is the base of natural 
                                                             logarithms. *)
         (CTYPECASE NUMBER (COMPLEX (%%EXP-COMPLEX NUMBER))
                (NUMBER (%%EXP-FLOAT (FLOAT NUMBER))))))

(%%EXP-FLOAT
  (CL:LAMBDA (X)                                             (* kbr: "23-Jul-86 13:45")
          
          (* * (EXP X) for float X calculated via EXPB 1103 rational approximation of 
          Harris et al. *)

         (PROG (RECIPFLG M N R ANSWER)
               (DECLARE (GLOBALVARS %%EXP-TABLE %%EXP-POLY))
          
          (* * First, arrange X to be in interval (0 infinity) via identity
          (EXP (minus X)) = (/ (EXP X)) *)

               (SETQ R X)
               (COND
                  ((< R 0)
                   (SETQ R (- R))
                   (SETQ RECIPFLG T)))
          
          (* * Next, the problem of (EXP X) is converted into a problem
          (EXPT 2 Y) where Y = (CL:* %%LOG-BASE2-E X)%.
          Then range reduction is accomplished via
          (EXPT 2 Y) = (CL:* (EXPT 2 M) (EXPT 2 (/ N 8))
          (EXPT 2 R)) where M and N are integers and R is a float in the interval
          (0.0 .125) After M N R are determined, (EXPT 2 M) is effected by scaling,
          (EXPT 2 (/ N 8)) is found by table lookup, and
          (EXPT 2 R) is calculated by rational approximation EXPB 1103 of Harris et al.
          *)

               (MULTIPLE-VALUE-SETQ (M R)
                      (TRUNCATE (CL:* %%LOG-BASE2-E R)))
               (MULTIPLE-VALUE-SETQ (N R)
                      (TRUNCATE R .125))
               (SETQ ANSWER (SCALE-FLOAT (CL:* (ELT %%EXP-TABLE N)
                                               (/ (POLYEVAL R %%EXP-POLY 5)
                                                  (POLYEVAL (- R)
                                                         %%EXP-POLY 5)))
                                   M))
               (COND
                  (RECIPFLG (SETQ ANSWER (/ ANSWER))))
               (RETURN ANSWER))))

(%%EXP-COMPLEX
  (CL:LAMBDA (Z)                                             (* compute exp (z) where z is complex 
                                                             is called by exp *)
         (LET* ((X (REALPART Z))
                (Y (IMAGPART Z)))
               (CL:* (EXP X)
                     (CIS Y)))))
)



(* EXPT *)

(DEFINEQ

(CL:EXPT
  (CL:LAMBDA (BASE-NUMBER POWER-NUMBER)                      (* lmm "27-Jun-86 22:52")
          
          (* This function calculates BASE-NUMBER raised to the nth power.
          It separates the cases by the type of POWER-NUMBER for efficiency reasons, as 
          powers can be calculated more efficiently if POWER-NUMBER is a positive 
          integer, Therefore, All integer values of POWER-NUMBER are calculated as 
          positive integers, and inverted if negative.
          *)

         (COND
            ((AND (RATIONALP BASE-NUMBER)
                  (INTEGERP POWER-NUMBER))
             (%%EXPT-INTEGER BASE-NUMBER POWER-NUMBER))
            ((= BASE-NUMBER 0)
             (CL:IF (PLUSP POWER-NUMBER)
                    BASE-NUMBER
                    (CL:ERROR "~A to a non-positive power ~A." BASE-NUMBER POWER-NUMBER)))
            ((AND (COMPLEXP BASE-NUMBER)
                  (INTEGERP POWER-NUMBER))
             (%%EXPT-INTEGER BASE-NUMBER POWER-NUMBER))
            ((COMPLEXP POWER-NUMBER)
             (CL:IF (OR (COMPLEXP BASE-NUMBER)
                        (PLUSP BASE-NUMBER))
                    (%%EXPT-COMPLEX-POWER BASE-NUMBER POWER-NUMBER)
                    (CL:ERROR "~A negative number to a complex power ~A." BASE-NUMBER POWER-NUMBER)))
            ((COMPLEXP BASE-NUMBER)
             (%%EXPT-COMPLEX BASE-NUMBER POWER-NUMBER))
            ((AND (NOT (INTEGERP POWER-NUMBER))
                  (MINUSP BASE-NUMBER))
             (CL:ERROR "Negative number ~A to non-integral power ~A." BASE-NUMBER POWER-NUMBER))
            (T (%%EXPT-FLOAT (FLOAT BASE-NUMBER)
                      (FLOAT POWER-NUMBER))))))

(%%EXPT-INTEGER
  (CL:LAMBDA (BASE POWER)                                    (* kbr: "16-Sep-86 13:03")
          
          (* * (EXPT BASE POWER) where BASE is rational and POWER is an integer.
          *)

         (COND
            ((EQ BASE MIN.INTEGER)
             (COND
                ((INTEGERP POWER)
                 (COND
                    ((< POWER 0)
                     0)
                    ((EQ POWER 0)
                     1)
                    ((EQ POWER MAX.INTEGER)
                     (ERROR "Can't raise negative infinity to infinite power."))
                    ((EVENP POWER)
                     MAX.INTEGER)
                    (T                                       (* Odd integer POWER *)
                       MIN.INTEGER)))
                (T                                           (* POWER is a noninteger rational *)
                   (ERROR "Can't raise negative infinity to noninteger power."))))
            ((EQ BASE MAX.INTEGER)
             (COND
                ((EQ POWER 0)
                 1)
                (T MAX.INTEGER)))
            ((EQ POWER MIN.INTEGER)
             (COND
                ((EQ BASE 0)
                 (ERROR "Can't expt 0 to a negative power."))
                (T 0)))
            ((EQ POWER MAX.INTEGER)
             (COND
                ((EQ BASE 0)
                 0)
                ((> BASE 0)
                 MAX.INTEGER)
                (T                                           (* Undefined because sign oscillates 
                                                             as you take the limit.
                                                             *)
                   (ERROR "Can't expt negative number to infinite power."))))
            ((MINUSP POWER)
             (/ (%%EXPT-INTEGER BASE (- POWER))))
            ((AND (RATIONALP BASE)
                  (NOT (INTEGERP BASE)))
             (%%MAKE-RATIO (%%EXPT-INTEGER (NUMERATOR BASE)
                                  POWER)
                    (%%EXPT-INTEGER (DENOMINATOR BASE)
                           POWER)))
            ((AND (INTEGERP BASE)
                  (= BASE 2))
             (ASH 1 POWER))
            (T (CL:DO ((NEXTN (ASH POWER -1)
                              (ASH POWER -1))
                       (TOTAL (CL:IF (ODDP POWER)
                                     BASE 1)
                              (CL:IF (ODDP POWER)
                                     (CL:* BASE TOTAL)
                                     TOTAL)))
                      ((ZEROP NEXTN)
                       TOTAL)
                      (SETQ BASE (CL:* BASE BASE))
                      (SETQ POWER NEXTN))))))

(%%EXPT-FLOAT
  (LAMBDA (X Y)                                              (* kbr: "29-May-86 00:03")
          
          (* * (EXPT X Y) where X is a nonnegative float and Y is a float.
          *)

    (COND
       ((= X 0.0)
        0.0)
       (T (EXP (CL:* Y (CL:LOG X)))))))

(%%EXPT-COMPLEX
  (CL:LAMBDA (Z N)                                           (* compute (complex) ↑n where n is an 
                                                             integer some round-off error if n is 
                                                             not a fixnum *)
         (CL:* (CL:EXPT (%%COMPLEX-ABS Z)
                      N)
               (CIS (CL:* N (PHASE Z))))))

(%%EXPT-COMPLEX-POWER
  (CL:LAMBDA (Z W)                                           (* this function computes z ↑ w where 
                                                             w is a complex number it can also be 
                                                             used for any positive real number.
                                                             *)
         (EXP (CL:* W (CL:LOG Z)))))
)



(* LOG *)


(RPAQ? %%LOG2 (%%FLOAT 16177 29208))

(RPAQ? %%SQRT2 (%%FLOAT 16309 1267))
(* * %%LOG-PPOLY and %%LOG-QPOLY contain P and Q coefficients of Harris et al LOGE 2707 
rational approximation to (LOG X) in interval ((SQRT .5) (SQRT 2)) . *)


(RPAQ %%LOG-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16042 22803)
                                      (%%FLOAT 49484 23590)
                                      (%%FLOAT 17044 17982)
                                      (%%FLOAT 49926 37153)
                                      (%%FLOAT 17046 5367))))

(RPAQ %%LOG-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                      (%%FLOAT 49512 9103)
                                      (%%FLOAT 16992 42274)
                                      (%%FLOAT 49823 38048)
                                      (%%FLOAT 16918 5367))))
(DEFINEQ

(CL:LOG
  (CL:LAMBDA (NUMBER &OPTIONAL (BASE NIL BASE-SUPPLIED))     (* lmm "27-Jun-86 22:40")
         (COND
            (BASE-SUPPLIED (/ (CL:LOG NUMBER)
                              (CL:LOG BASE)))
            (T (CTYPECASE NUMBER (COMPLEX (%%LOG-COMPLEX NUMBER))
                      (NUMBER (CL:IF (MINUSP NUMBER)
                                     (COMPLEX (%%LOG-FLOAT (FLOAT (- NUMBER)))
                                            %%PI)
                                     (%%LOG-FLOAT (FLOAT NUMBER)))))))))

(%%LOG-FLOAT
  (LAMBDA (X)                                                (* kbr: "30-May-86 00:04")
          
          (* * (EXP X) for float X calculated via EXPB 1103 rational approximation of 
          Harris et al. *)
          
          (* * (LOG X) for nonnegative float X. *)

    (PROG (R EXP Z Z2 ANSWER)
          (DECLARE (GLOBALVARS %%LOG-PPOLY %%LOG-QPOLY))
          
          (* * NOTE: Probably best not to declare type of variables inside this routine 
          for now. I've found that FLOATP smashing combined with FLOATP declarations 
          compiles wrong. *)

          (COND
             ((NOT (FGREATERP X 0.0))
              (ERROR "LOG OF ZERO:" X)))
          
          (* * Range reduce to an R in interval ((SQRT .5)
          (SQRT 2)) via identity (LOG X) = (+ (LOG R)
          (CL:* %%LOG-2 EXP)) for a suitable integer EXP.
          This reduction depends on and smashes the FLOATP representation of R.
          *)

          (SETQ R (\FLOAT.BOX X))
          (SETQ EXP (- (fetch (FLOATP EXPONENT) of R)
                       \EXPONENT.BIAS))
          (replace (FLOATP EXPONENT) of R with \EXPONENT.BIAS)
          (COND
             ((> R %%SQRT2)
              (SETQ EXP (1+ EXP))
              (SETQ R (/ R 2.0))))
          
          (* * (LOG R) is calculated by rational approximation LOGE 2707 of Harris et al.
          *)

          (SETQ Z (/ (1- R)
                     (1+ R)))
          (SETQ Z2 (CL:* Z Z))
          (SETQ ANSWER (+ (CL:* Z (/ (POLYEVAL Z2 %%LOG-PPOLY 4)
                                     (POLYEVAL Z2 %%LOG-QPOLY 4)))
                          (CL:* %%LOG2 EXP)))
          (RETURN ANSWER))))

(%%LOG-COMPLEX
  (CL:LAMBDA (Z)                                             (* natural log of a complex number *)
         (COMPLEX (CL:LOG (%%COMPLEX-ABS Z))
                (PHASE Z))))
)



(* SQRT *)

(DEFINEQ

(CL:SQRT
  (CL:LAMBDA (X)                                             (* kbr: "28-Jul-86 14:06")
         (COND
            ((COMPLEXP X)
             (%%SQRT-COMPLEX X))
            ((MINUSP X)                                      (* Negative real axis maps into 
                                                             positive imaginary axis.
                                                             *)
             (COMPLEX 0 (CL:SQRT (- X))))
            (T (%%SQRT-FLOAT (FLOAT X))))))

(%%SQRT-FLOAT
  (LAMBDA (X)                                                (* kbr: "21-Jul-86 16:38")
          
          (* * (SQRT X) for nonnegative float X. *)

    (PROG (V)
          (DECLARE (TYPE FLOATP X V))
          (COND
             ((NOT (FGREATERP X 0.0))                        (* Trichotomy ==> X = 0.0 *)
              (RETURN 0.0)))
          (SETQ V (create FLOATP
                         EXPONENT ← (LOGAND (IPLUS \EXPONENT.BIAS
                                                   (LRSH (LOGAND (IDIFFERENCE (fetch (FLOATP EXP)
                                                                                 of X)
                                                                        \EXPONENT.BIAS)
                                                                (MASK.1'S 0 BITSPERWORD))
                                                         1))
                                           \MAX.EXPONENT)
                         HIFRACTION ← (fetch (FLOATP HIFRAC) of X)))
          
          (* Exponent is stored as excess \EXPONENT.BIAS and although the LRSH doesn't 
          really do division by 2 (e.g., when the arg is negative) at least the low-order 
          8 bits will be right. It doesn't even matter that it may be off-by-one, due to 
          the infamous "Arithmetic Shifting Considered Harmful" since it is only an 
          estimate. *)

          (FRPTQ 4 (SETQ V (CL:* .5 (+ V (/ X V)))))
          (RETURN V))))

(%%SQRT-COMPLEX
  (LAMBDA (X)                                                (* kbr: "28-Jul-86 14:04")
          
          (* * (SQRT X) for complex X. *)

    (PROG (ABS PHASE A B C D E ANSWER)
          (DECLARE (TYPE FLOATP ABS PHASE A B C D E))
          (SETQ A (FLOAT (REALPART X)))
          (SETQ B (FLOAT (IMAGPART X)))
          
          (* * Make initial guess. *)

          (SETQ ABS (CL:SQRT (ABS X)))
          (SETQ PHASE (/ (PHASE X)
                         2.0))
          (SETQ C (CL:* ABS (CL:COS PHASE)))
          (SETQ D (CL:* ABS (CL:SIN PHASE)))
          
          (* * Newton's method. *)

          (for I from 1 to 4 do (SETQ E (+ (CL:* C C)
                                           (CL:* D D)))
                                (SETQ C (CL:* .5 (+ C (/ (+ (CL:* A C)
                                                            (CL:* B D))
                                                         E))))
                                (SETQ D (CL:* .5 (+ D (/ (- (CL:* B C)
                                                            (CL:* A D))
                                                         E)))))
          (SETQ ANSWER (COMPLEX C D))
          (RETURN ANSWER))))
)



(* SIN *)

(* * %%SIN-EPSILON is sufficiently small that (SIN X) = X for X in interval (0 %%SIN-EPSILON) 
%. It suffices to take %%SIN-EPSILON a little bit smaller than (SQRT (CL:* 6 
SINGLE-FLOAT-EPSILON)) which we get by the Taylor series expansion (SIN X) = (+ X (/ (EXPT X 3)
 6) ...) (The relative error caused by ommitting (/ (EXPT X 3) 6) isn't observable.) Comparison
 against %%SIN-EPSILON is used to avoid POLYEVAL microcode underflow when computing SIN. *)


(RPAQ? %%SIN-EPSILON (%%FLOAT 14720 0))
(* * %%SIN-PPOLY and %%SIN-QPOLY contain adapted P and Q coefficients of Harris et al SIN 3374 
rational approximation to (SIN X) in interval (0 (/ PI 2)) %. The coefficients for %%SIN-PPOLY 
and %%SIN-QPOLY have been computed from Harris using extended precision routines and the 
relations %%SIN-PPOLY = (REVERSE (for I from 0 as ENTRY in PS collect (/ (CL:* (EXPT (/ 2 PI) (
1+ (CL:* 2 I))) ENTRY) Q0))) and %%SIN-QPOLY = (REVERSE (for I from 0 as ENTRY in QS collect (/
 (CL:* (EXPT (/ 2 PI) (CL:* 2 I)) ENTRY) Q0))) *)


(RPAQ %%SIN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 45236 25611)
                                      (%%FLOAT 13589 26148)
                                      (%%FLOAT 47286 34797)
                                      (%%FLOAT 15295 3306)
                                      (%%FLOAT 48666 34805)
                                      (%%FLOAT 16256 0))))

(RPAQ %%SIN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 11384 52865)
                                      (%%FLOAT 12553 9550)
                                      (%%FLOAT 13604 38385)
                                      (%%FLOAT 14593 18841)
                                      (%%FLOAT 15489 5549)
                                      (%%FLOAT 16256 0))))
(DEFINEQ

(CL:SIN
  (CL:LAMBDA (RADIANS)                                       (* kbr: " 7-Apr-86 22:57")
         (CTYPECASE RADIANS (COMPLEX (%%SIN-COMPLEX RADIANS))
                (FLOAT (%%SIN-FLOAT RADIANS NIL))
                (NUMBER (%%SIN-FLOAT (FLOAT RADIANS)
                               NIL)))))

(%%SIN-FLOAT
  (LAMBDA (X COS-FLAG)                                       (* kbr: "23-Jul-86 15:53")
          
          (* * SIN of a FLOAT X calculated via SIN 3374 rational approximation of Harris 
          et al. *)

    (PROG (R SIGN R2 ANSWER)
          (DECLARE (GLOBALVARS %%SIN-PPOLY %%SIN-QPOLY)
                 (TYPE FLOATP X R SIGN R2 ANSWER))
          
          (* * If this function called by COS then use
          (COS X) = (SIN (minus %%PI/2 X)) = (SIN (+ %%PI/2 X)) Case out on sign of X for 
          improved numerical stability. Avoids unnecessary rounding and promotes 
          symmetric properties. (COS X) = (COS (minus X)) is guaranteed by this strategy.
          *)

          (SETQ R (COND
                     ((NOT COS-FLAG)
                      X)
                     ((> X 0)
                      (- %%PI/2 X))
                     (T (+ %%PI/2 X))))
          
          (* * First range reduce to (0 infinity) by
          (SIN (minus X)) = (minus (SIN X)) This strategy guarantees
          (SIN (minus X)) = (minus (SIN X)) *)

          (COND
             ((< R 0)
              (SETQ SIGN -1.0)
              (SETQ R (- R)))
             (T (SETQ SIGN 1.0)))
          
          (* * Next range reduce to interval (0 %%2PI) by
          (SIN X) = (SIN (MOD X %%2PI)) . *)

          (SETQ R (REM R %%2PI))
          
          (* * Next range reduce to interval (0 PI) by
          (SIN (+ X PI)) = (minus (SIN X)) *)

          (COND
             ((> R %%PI)
              (SETQ SIGN (- SIGN))
              (SETQ R (- R %%PI))))
          
          (* * Next range reduce to interval (0 %%PI/2) by
          (SIN (+ X %%PI/2)) = (SIN (minus %%PI/2 X)) *)

          (COND
             ((> R %%PI/2)
              (SETQ R (- %%PI R))))
          (COND
             ((< R %%SIN-EPSILON)
          
          (* * If R is in the interval (0 %%SIN-EPSILON) then
          (SIN R) = R to the precision that we can offer.
          Return R because (1) it is desirable that
          (SIN R) = R exactly for small R and (2) microcode POLYEVAL will underflow on 
          sufficiently small positive R. *)

              (RETURN (CL:* SIGN R))))
          
          (* * Now use SIN 3374 rational approximation of Harris et al.
          which works on interval (0 %%PI/2) *)

          (SETQ R2 (CL:* R R))
          (SETQ ANSWER (CL:* SIGN R (/ (POLYEVAL R2 %%SIN-PPOLY 5)
                                       (POLYEVAL R2 %%SIN-QPOLY 5))))
          (RETURN ANSWER))))

(%%SIN-COMPLEX
  (CL:LAMBDA (Z)                                             (* sin of a complex number *)
         (LET* ((X (REALPART Z))
                (Y (IMAGPART Z)))
               (COMPLEX (CL:* (CL:SIN X)
                              (COSH Y))
                      (CL:* (CL:COS X)
                            (SINH Y))))))
)



(* COS *)

(DEFINEQ

(CL:COS
  (CL:LAMBDA (RADIANS)                                       (* lmm " 1-Jul-86 15:55")
         (CTYPECASE RADIANS (FLOAT (%%SIN-FLOAT RADIANS T))
                (COMPLEX (%%COS-COMPLEX RADIANS))
                (NUMBER (%%SIN-FLOAT (FLOAT RADIANS)
                               T)))))

(%%COS-COMPLEX
  (CL:LAMBDA (Z)                                             (* cosine of a complex number *)
         (LET* ((X (REALPART Z))
                (Y (IMAGPART Z)))
               (COMPLEX (CL:* (CL:COS X)
                              (COSH Y))
                      (- (CL:* (CL:SIN X)
                               (SINH Y)))))))
)



(* TAN *)

(* * %%TAN-EPSILON is sufficiently small that (TAN X) = X for X in interval (0 %%TAN-EPSILON) 
%. It suffices to take %%TAN-EPSILON a little bit smaller than (SQRT (CL:* 3 
SINGLE-FLOAT-EPSILON)) which we get by the Taylor series expansion (TAN X) = (+ X (/ (EXPT X 3)
 3) ...) (The relative error caused by ommitting (/ (EXPT X 3) 3) isn't observable.) Comparison
 against %%TAN-EPSILON is used to avoid POLYEVAL microcode underflow when computing TAN. *)


(RPAQ? %%TAN-EPSILON (%%FLOAT 14720 0))
(* * %%TAN-PPOLY and %%TAN-QPOLY contain adapted P and Q coefficients of Harris et al TAN 4288 
rational approximation to (TAN X) in interval (-PI/4 PI/4) %. The coefficients for %%TAN-PPOLY 
and %%TAN-QPOLY have been computed from Harris using extended precision routines and the 
relations %%TAN-PPOLY = (REVERSE (for I from 0 as ENTRY in PS collect (/ (CL:* (EXPT (/ 4 PI) (
1+ (CL:* 2 I))) ENTRY) Q0))) and %%TAN-QPOLY = (REVERSE (for I from 0 as ENTRY in QS collect (/
 (CL:* (EXPT (/ 4 PI) (CL:* 2 I)) ENTRY) Q0))) *)


(RPAQ %%TAN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 13237 21090)
                                      (%%FLOAT 47141 15825)
                                      (%%FLOAT 15246 8785)
                                      (%%FLOAT 48655 48761)
                                      (%%FLOAT 16256 0))))

(RPAQ %%TAN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 45267 36947)
                                      (%%FLOAT 13848 46875)
                                      (%%FLOAT 47612 53738)
                                      (%%FLOAT 15596 52854)
                                      (%%FLOAT 48882 35303)
                                      (%%FLOAT 16256 0))))
(DEFINEQ

(CL:TAN
  (CL:LAMBDA (RADIANS)                                       (* lmm "27-Jun-86 22:41")
         (CTYPECASE RADIANS (COMPLEX (%%TAN-COMPLEX RADIANS))
                (NUMBER (%%TAN-FLOAT (FLOAT RADIANS))))))

(%%TAN-FLOAT
  (LAMBDA (X)                                                (* kbr: "20-Aug-86 20:13")
          
          (* * TAN of a FLOAT X calculated via TAN 4288 rational approximation of Harris 
          et al. *)

    (PROG (R SIGN RECIPFLG R2 ANSWER)
          (DECLARE (GLOBALVARS %%TAN-PPOLY %%TAN-QPOLY)
                 (TYPE FLOATP X R R2 ANSWER))
          (SETQ R X)
          
          (* * First range reduce to (0 infinity) by
          (TAN (minus X)) = (minus (TAN X)) *)

          (COND
             ((< R 0)
              (SETQ SIGN -1.0)
              (SETQ R (- R)))
             (T (SETQ SIGN 1.0)))
          
          (* * Next range reduce to (0 PI) *)

          (SETQ R (REM R %%PI))
          
          (* * Next, range reduce to (-PI/4 PI/4) using
          (TAN X) = (TAN (minus X PI)) to get into interval
          (-PI/2 PI/2) and then (TAN X) = (/ (TAN (minus PI/2 X))) to get into interval
          (-PI/4 PI/4) *)

          (COND
             ((> R %%PI/2)
              (SETQ R (- R %%PI))
              (COND
                 ((< R %%-PI/4)
                  (SETQ RECIPFLG T)
                  (SETQ R (- %%-PI/2 R)))))
             (T (COND
                   ((> R %%PI/4)
                    (SETQ RECIPFLG T)
                    (SETQ R (- %%PI/2 R))))))
          (COND
             ((< (ABS R)
               %%TAN-EPSILON)
          
          (* * If R is in the interval (0 %%TAN-EPSILON) then
          (TAN R) = R to the precision that we can offer.
          Return R because (1) it is desirable that
          (TAN R) = R exactly for small R and (2) microcode POLYEVAL will underflow on 
          sufficiently small positive R. *)

              (SETQ ANSWER (CL:* SIGN R))
              (COND
                 (RECIPFLG (SETQ ANSWER (/ ANSWER))))
              (RETURN ANSWER)))
          
          (* * Now use TAN 4288 rational approximation of Harris et al.
          which works on interval (0 %%PI/4) *)

          (SETQ R2 (CL:* R R))
          (SETQ ANSWER (CL:* SIGN R (/ (POLYEVAL R2 %%TAN-PPOLY 4)
                                       (POLYEVAL R2 %%TAN-QPOLY 5))))
          (COND
             (RECIPFLG (SETQ ANSWER (/ ANSWER))))
          (RETURN ANSWER))))

(%%TAN-COMPLEX
  (CL:LAMBDA (Z)                                             (* tan of a complex number there was a 
                                                             nicer algorithm but it turned out not 
                                                             to work so well. *)
         (LET* ((NUM (CL:SIN Z))
                (DENOM (CL:COS Z)))
               (CL:IF (ZEROP DENOM)
                      (CL:ERROR "~S undefined tangent." Z)
                      (/ NUM DENOM)))))
)



(* ASIN *)

(* * %%ASIN-EPSILON is sufficiently small that (ASIN X) = X for X in interval (0 %%ASIN-EPSILON
) %. It suffices to take %%ASIN-EPSILON a little bit smaller than (CL:* 2 SINGLE-FLOAT-EPSILON)
 which we get by the Taylor series expansion (ASIN X) = (+ X (/ (EXPT X 3) 6) ...) (The 
relative error caused by ommitting (/ (EXPT X 3) 6) isn't observable.) Comparison against 
%%ASIN-EPSILON is used to avoid POLYEVAL microcode underflow when computing SIN. *)


(RPAQ? %%ASIN-EPSILON (%%FLOAT 14720 0))
(* * %%ASIN-PPOLY and %%ASIN-QPOLY contain P and Q coefficients of Harris et al ARCSN 4671 
rational approximation to (ASIN X) in interval (0 (SQRT .5)) . *)


(RPAQ %%ASIN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16007 50045)
                                       (%%FLOAT 49549 8020)
                                       (%%FLOAT 17236 15848)
                                       (%%FLOAT 50285 63464)
                                       (%%FLOAT 17650 31235)
                                       (%%FLOAT 50403 62852)
                                       (%%FLOAT 17440 39471))))

(RPAQ %%ASIN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                       (%%FLOAT 49672 25817)
                                       (%%FLOAT 17308 55260)
                                       (%%FLOAT 50326 38098)
                                       (%%FLOAT 17674 22210)
                                       (%%FLOAT 50417 22451)
                                       (%%FLOAT 17440 39471))))
(DEFINEQ

(ASIN
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
         (CTYPECASE NUMBER (FLOAT (%%ASIN-FLOAT NUMBER))
                (COMPLEX (%%ASIN-COMPLEX NUMBER))
                (NUMBER (%%ASIN-FLOAT (FLOAT NUMBER))))))

(%%ASIN-FLOAT
  (LAMBDA (X ACOS-FLAG)                                      (* kbr: "30-May-86 16:46")
          
          (* * (ASIN X) for float X calculated via ARCSN 4671 rational approximation of 
          Harris et al. *)

    (PROG (NEGATIVE REDUCED R R2 ANSWER)
          (DECLARE (GLOBALVARS %%ASIN-PPOLY %%ASIN-QPOLY))
          (SETQ R X)
          (COND
             ((OR (< R -1.0)
                  (> R 1.0))
              (ERROR "ARCSIN: arg not in range" R)))
          (COND
             ((< R 0.0)
          
          (* * Range reduce to (0 1) via identity (ASIN
          (minus X)) = (minus (ASIN X)) *)

              (SETQ NEGATIVE T)
              (SETQ R (- R))))
          (COND
             ((> R .5)
          
          (* * Range reduce to (0 .5) via identity
          (ASIN X) = (minus %%PI/2 (CL:* 2.0 (ASIN
          (CL:SQRT (CL:* .5 (minus 1.0 R)))))) Avoids numerical instability calculating
          (ASIN X) for X near one. SIN is horizontally flat near %%PI/2 so calculating
          (ASIN X) by rational approximation wouldn't work well for X near
          (SIN %%PI/2) = 1 *)

              (SETQ REDUCED T)
              (SETQ R (CL:SQRT (CL:* .5 (- 1.0 R))))))
          
          (* * R is now in range (0 .5) Use ARCSN 4671 rational approximation to 
          calculate (ASIN R) *)

          (SETQ ANSWER (COND
                          ((< R %%ASIN-EPSILON)
          
          (* * If R is in the interval (0 %%SIN-EPSILON) then
          (ASIN R) = R to the precision that we can offer.
          Return R because (1) it is desirable that
          (ASIN R) = R exactly for small R and (2) microcode POLYEVAL will underflow on 
          sufficiently small positive R. *)

                           R)
                          (T (SETQ R2 (CL:* R R))
                             (CL:* R (/ (POLYEVAL R2 %%ASIN-PPOLY 6)
                                        (POLYEVAL R2 %%ASIN-QPOLY 6))))))
          (COND
             (REDUCED (SETQ ANSWER (- %%PI/2 (CL:* 2.0 ANSWER)))))
          (COND
             (NEGATIVE (SETQ ANSWER (- ANSWER))))
          (COND
             (ACOS-FLAG 
          
          (* * In case we want (ACOS X) then use identity
          (ACOS X) = (minus %%PI/2 (ASIN X)) *)

                    (SETQ ANSWER (- %%PI/2 ANSWER))))
          (RETURN ANSWER))))

(%%ASIN-COMPLEX
  (CL:LAMBDA (Z)                                             (* kbr: "27-May-86 22:38")
         (%%COMPLEX-MINUS (%%COMPLEX-TIMESI (CL:LOG (+ (%%COMPLEX-TIMESI Z)
                                                       (CL:SQRT (- 1 (CL:* Z Z)))))))))
)



(* ACOS *)

(DEFINEQ

(ACOS
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
         (CTYPECASE NUMBER (FLOAT (%%ASIN-FLOAT NUMBER T))
                (COMPLEX (%%ACOS-COMPLEX NUMBER))
                (NUMBER (%%ASIN-FLOAT (FLOAT NUMBER)
                               T)))))

(%%ACOS-COMPLEX
  (CL:LAMBDA (Z)                                             (* kbr: "27-May-86 22:31")
         (%%COMPLEX-MINUS (%%COMPLEX-TIMESI (CL:LOG (+ Z (%%COMPLEX-TIMESI (CL:SQRT
                                                                            (- 1 (CL:* Z Z))))))))))
)



(* ATAN *)


(RPAQ? %%SQRT3 (%%FLOAT 16349 46039))

(RPAQ? %%2-SQRT3 (%%FLOAT 16009 12451))

(RPAQ? %%INV-2-SQRT3 (%%FLOAT 16494 55788))
(DEFINEQ

(CL:ATAN
  (CL:LAMBDA (X &OPTIONAL Y)                                 (* lmm "27-Jun-86 22:25")
         (COND
            (Y (%%ATAN-FLOAT (FLOAT X)
                      (FLOAT Y)))
            ((COMPLEXP X)
             (%%ATAN-COMPLEX X))
            (T (%%ATAN-FLOAT (FLOAT X))))))

(%%ATAN-FLOAT1
  (CL:LAMBDA (X)                                             (* kbr: "23-Jul-86 18:26")
         (COND
            ((MINUSP X)                                      (* (ATAN (minus X)) =
                                                             (minus (ATAN X)) *)
             (- (%%ATAN-FLOAT1 (- X))))
            ((< X %%2-SQRT3)                                 (* (ATAN X) = (minus %%PI/2
                                                             (ATAN (/ X))) *)
             (%%ATAN-FLOAT2 X))
            ((< X 1)
             (+ %%PI/6 (%%ATAN-FLOAT2 (/ (1- (CL:* X %%SQRT3))
                                         (+ %%SQRT3 X)))))
            ((< X %%INV-2-SQRT3)
             (- %%PI/3 (%%ATAN-FLOAT2 (/ (- %%SQRT3 X)
                                         (1+ (CL:* X %%SQRT3))))))
            (T (- %%PI/2 (%%ATAN-FLOAT2 (/ X)))))))

(%%ATAN-FLOAT2
  (CL:LAMBDA (X)                                             (* kbr: "23-Jul-86 18:26")
         (CL:DO* ((SQR (- (CL:* X X)))
                  (INT 1 (+ 2 INT))
                  (OLD 0 VAL)
                  (POW X (CL:* POW SQR))
                  (VAL POW (+ VAL (/ POW INT))))
                ((= OLD VAL)
                 VAL))))

(%%ATAN-DOMAIN-CHECK
  (CL:LAMBDA (Z)                                             (* kbr: "27-May-86 22:42")
                                                             (* Return T if Z is in the domain of 
                                                             ATAN. ATAN is singular at i and -i.
                                                             *)
         (NOT (AND (ZEROP (REALPART Z))
                   (= (ABS (IMAGPART Z))
                      1)))))

(%%ATAN-FLOAT
  (CL:LAMBDA (Y &OPTIONAL X)
         (COND
            ((NOT X)
             (CL:IF (ZEROP Y)
                    0.0
                    (%%ATAN-FLOAT1 Y)))
            ((= Y X 0)
             (CL:ERROR "Error in double entry atan both 0."))
            ((= X 0)
             (CL:* (%%SIGNUM Y)
                   %%PI/2))
            ((= Y 0)
             (CL:IF (PLUSP X)
                    0 %%PI))
            ((AND (PLUSP X)
                  (PLUSP Y))
             (%%ATAN-FLOAT1 (/ Y X)))
            ((PLUSP Y)
             (- %%PI (%%ATAN-FLOAT1 (/ (- Y)
                                       X))))
            ((PLUSP X)
             (- (%%ATAN-FLOAT1 (/ Y (- X)))))
            (T (- (%%ATAN-FLOAT1 (/ Y X))
                  %%PI)))))

(%%ATAN-COMPLEX
  (CL:LAMBDA (Z)
         (CL:IF (%%ATAN-DOMAIN-CHECK Z)
                (LET ((I (%%COMPLEX-TIMESI Z)))
                     (%%COMPLEX-MINUS (%%COMPLEX-TIMESI (CL:* .5 (CL:LOG (/ (+ 1 I)
                                                                            (- 1 I)))))))
                (CL:ERROR "Argument not in domain for atan. ~S" Z))))
)



(* CIS *)

(DEFINEQ

(CIS
  (CL:LAMBDA (RADIANS)                                       (* kbr: "12-Jul-86 18:05")
                                                             (* Return cos (Theta) + i sin
                                                             (Theta)%. *)
         (CL:IF (COMPLEXP RADIANS)
                (CL:ERROR "Argument to CIS is complex: ~S" RADIANS)
                (COMPLEX (CL:COS RADIANS)
                       (CL:SIN RADIANS)))))
)



(* SINH COSH TANH *)

(DEFINEQ

(SINH
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
          
          (* Hyperbolic trig functions. Each of the hyperbolic trig functions is 
          calculated directly from their definition.
          Exp (x) is calculated only once for efficiency.
          They all work with complex arguments without modification.
          *)

         (LET ((Z (EXP NUMBER)))
              (/ (- Z (/ Z))
                 2))))

(COSH
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
         (LET ((Z (EXP NUMBER)))
              (/ (+ Z (/ Z))
                 2))))

(TANH
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
                                                             (* Different form than in the manual.
                                                             Does much better. *)
         (LET* ((Z (EXP (CL:* 2 NUMBER)))
                (Y (/ Z)))
               (- (/ (1+ Y))
                  (/ (1+ Z))))))
)



(* ASINH ACOSH ATANH *)

(DEFINEQ

(ASINH
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
         (CL:LOG (+ NUMBER (CL:SQRT (+ (CL:* NUMBER NUMBER)
                                       1))))))

(ACOSH
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
         (CL:LOG (+ NUMBER (CL:SQRT (- (CL:* NUMBER NUMBER)
                                       1))))))

(ATANH
  (CL:LAMBDA (NUMBER)                                        (* kbr: "12-Jul-86 18:05")
         (CL:IF (%%ATANH-DOMAIN-CHECK NUMBER)
                (CL:* .5 (CL:LOG (/ (1+ NUMBER)
                                    (- 1 NUMBER))))
                (CL:ERROR "~S argument out of range." NUMBER))))

(%%ATANH-DOMAIN-CHECK
  (CL:LAMBDA (Z)                                             (* kbr: "27-May-86 22:49")
                                                             (* Return T if Z is in the domain of 
                                                             ATANH. ATANH is singular at 1 and -1.0 
                                                             *)
         (NOT (AND (ZEROP (IMAGPART Z))
                   (= (ABS (REALPART Z))
                      1)))))
)

(PUTPROPS CMLFLOAT FILETYPE COMPILE-FILE)
(DECLARE: DONTEVAL@LOAD DOEVAL@COMPILE DONTCOPY COMPILERVARS 

(ADDTOVAR NLAMA )

(ADDTOVAR NLAML )

(ADDTOVAR LAMA 
          %%ATANH-DOMAIN-CHECK ATANH ACOSH ASINH TANH COSH SINH CIS %%ATAN-COMPLEX %%ATAN-FLOAT 
                %%ATAN-DOMAIN-CHECK %%ATAN-FLOAT2 %%ATAN-FLOAT1 CL:ATAN %%ACOS-COMPLEX ACOS 
                %%ASIN-COMPLEX ASIN %%TAN-COMPLEX CL:TAN %%COS-COMPLEX CL:COS %%SIN-COMPLEX CL:SIN 
                CL:SQRT %%LOG-COMPLEX CL:LOG %%EXPT-COMPLEX-POWER %%EXPT-COMPLEX %%EXPT-INTEGER 
                CL:EXPT %%EXP-COMPLEX %%EXP-FLOAT EXP)
)
(PRETTYCOMPRINT CMLFLOATCOMS)

(RPAQQ CMLFLOATCOMS ((* * CMLFLOAT -- Covering sections 12.5-12.5.3 irrational, transcendental, 
                        exponential, logarithmic, trigonometric, and hyperbolic functions. Section 
                        12.10, implementation parameters. -- By Kelly Roach. *)
                     (DECLARE: EVAL@COMPILE DONTCOPY (FILES (LOADCOMP)
                                                            LLFLOAT))
                     (COMS (* Section 12.10, implementation parameters. The constants in this COMS 
                              are exported to the user. *)
                           (* %%FLOAT allows us to recreate FLOATPs in a way that is independent of 
                              the ordinairy reading and printing FLOATPs to files which involves loss 
                              of the last couple bits of accuracy due to rounding effects. *)
                           (* Using INITVARS instead of CONSTANTS in various places here because of 
                              problems with the way BYTECOMPILER stores FLOATPs in DCOM files. *)
                           (FNS %%FLOAT)
                           (CONSTANTS (MOST-POSITIVE-FIXNUM 65535)
                                  (MOST-NEGATIVE-FIXNUM -65536))
                           (INITVARS (MOST-POSITIVE-SINGLE-FLOAT (%%FLOAT 32639 65535))
                                  (LEAST-POSITIVE-SINGLE-FLOAT (%%FLOAT 128 0))
                                  (LEAST-NEGATIVE-SINGLE-FLOAT (%%FLOAT 32896 0))
                                  (MOST-NEGATIVE-SINGLE-FLOAT (%%FLOAT 65407 65535))
                                  (MOST-POSITIVE-SHORT-FLOAT MOST-POSITIVE-SINGLE-FLOAT)
                                  (LEAST-POSITIVE-SHORT-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)
                                  (LEAST-NEGATIVE-SHORT-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)
                                  (MOST-NEGATIVE-SHORT-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)
                                  (MOST-POSITIVE-DOUBLE-FLOAT MOST-POSITIVE-SINGLE-FLOAT)
                                  (LEAST-POSITIVE-DOUBLE-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)
                                  (LEAST-NEGATIVE-DOUBLE-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)
                                  (MOST-NEGATIVE-DOUBLE-FLOAT MOST-NEGATIVE-SINGLE-FLOAT)
                                  (MOST-POSITIVE-LONG-FLOAT MOST-POSITIVE-SINGLE-FLOAT)
                                  (LEAST-POSITIVE-LONG-FLOAT LEAST-POSITIVE-SINGLE-FLOAT)
                                  (LEAST-NEGATIVE-LONG-FLOAT LEAST-NEGATIVE-SINGLE-FLOAT)
                                  (MOST-NEGATIVE-LONG-FLOAT MOST-NEGATIVE-SINGLE-FLOAT))
                           (* EPSILON is the smallest positive floating point number such that
                              (NOT (= (FLOAT 1 EPSILON)
                                      (+ (FLOAT 1 EPSILON)
                                         EPSILON)))
                              *)
                           (INITVARS (SINGLE-FLOAT-EPSILON (%%FLOAT 13312 0))
                                  (SHORT-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)
                                  (DOUBLE-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)
                                  (LONG-FLOAT-EPSILON SINGLE-FLOAT-EPSILON))
                           (* NEGATIVE-EPSILON is the smallest negative floating point number such 
                              that (NOT (= (FLOAT 1 NEGATIVE-EPSILON)
                                           (- (FLOAT 1 NEGATIVE-EPSILON)
                                              NEGATIVE-EPSILON)))
                              *)
                           (INITVARS (SINGLE-FLOAT-NEGATIVE-EPSILON (%%FLOAT 13184 0))
                                  (SHORT-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)
                                  (DOUBLE-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)
                                  (LONG-FLOAT-NEGATIVE-EPSILON SINGLE-FLOAT-NEGATIVE-EPSILON)))
                     (COMS (* Miscellaneous. *)
                           (DECLARE: EVAL@COMPILE DONTCOPY (FILES (LOADCOMP)
                                                                  LLFLOAT))
                           (INITVARS (PI (%%FLOAT 16457 4059)))
                           (* Should be (DECLARE: EVAL@COMPILE DONTCOPY (CONSTANTS ...))
                              except that compiler does a poor job of compiling FLOATPs. Use an 
                              INITVARS to patch around this situation for now. *)
                           (INITVARS (%%E (%%FLOAT 16429 63572))
                                  (%%2PI (%%FLOAT 16585 4059))
                                  (%%PI (%%FLOAT 16457 4059))
                                  (%%2PI/3 (%%FLOAT 16390 2706))
                                  (%%PI/2 (%%FLOAT 16329 4059))
                                  (%%-PI/2 (%%FLOAT 49097 4059))
                                  (%%PI/3 (%%FLOAT 16262 2706))
                                  (%%PI/4 (%%FLOAT 16201 4059))
                                  (%%-PI/4 (%%FLOAT 48969 4059))
                                  (%%PI/6 (%%FLOAT 16134 2706))
                                  (%%2/PI (%%FLOAT 16162 63875)))
                           (FNS %%MAKE-ARRAY))
                     (COMS (* EXP *)
                           (COMS (INITVARS (%%LOG-BASE2-E (%%FLOAT 16312 43579)))
                                 (* * %%EXP-POLY contains P and Q coefficients of Harris et al EXPB 
                                    1103 rational approximation to (EXPT 2 X)
                                    in interval (0 .125)
                                    %. %%EXP-TABLE contains values of powers (EXPT 2 (/ N 8)) . *)
                                 (VARS (%%EXP-POLY (%%MAKE-ARRAY (LIST (%%FLOAT 15549 17659)
                                                                       (%%FLOAT 16256 0)
                                                                       (%%FLOAT 16801 38273)
                                                                       (%%FLOAT 17257 7717)
                                                                       (%%FLOAT 17597 11739)
                                                                       (%%FLOAT 17800 30401))))
                                       (%%EXP-TABLE (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                                                        (%%FLOAT 16267 38338)
                                                                        (%%FLOAT 16280 14320)
                                                                        (%%FLOAT 16293 65239)
                                                                        (%%FLOAT 16309 1267)
                                                                        (%%FLOAT 16325 26410)
                                                                        (%%FLOAT 16343 17661)
                                                                        (%%FLOAT 16362 49351))))))
                           (FNS EXP %%EXP-FLOAT %%EXP-COMPLEX))
                     (COMS (* EXPT *)
                           (FNS CL:EXPT %%EXPT-INTEGER %%EXPT-FLOAT %%EXPT-COMPLEX 
                                %%EXPT-COMPLEX-POWER))
                     (COMS (* LOG *)
                           (COMS (INITVARS (%%LOG2 (%%FLOAT 16177 29208))
                                        (%%SQRT2 (%%FLOAT 16309 1267)))
                                 (* * %%LOG-PPOLY and %%LOG-QPOLY contain P and Q coefficients of 
                                    Harris et al LOGE 2707 rational approximation to (LOG X)
                                    in interval ((SQRT .5)
                                                 (SQRT 2)) . *)
                                 (VARS (%%LOG-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16042 22803)
                                                                        (%%FLOAT 49484 23590)
                                                                        (%%FLOAT 17044 17982)
                                                                        (%%FLOAT 49926 37153)
                                                                        (%%FLOAT 17046 5367))))
                                       (%%LOG-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                                                        (%%FLOAT 49512 9103)
                                                                        (%%FLOAT 16992 42274)
                                                                        (%%FLOAT 49823 38048)
                                                                        (%%FLOAT 16918 5367))))))
                           (FNS CL:LOG %%LOG-FLOAT %%LOG-COMPLEX))
                     (COMS (* SQRT *)
                           (FNS CL:SQRT %%SQRT-FLOAT %%SQRT-COMPLEX))
                     (COMS (* SIN *)
                           (COMS (* * %%SIN-EPSILON is sufficiently small that (SIN X)
                                    = X for X in interval (0 %%SIN-EPSILON)
                                    %. It suffices to take %%SIN-EPSILON a little bit smaller than
                                    (SQRT (CL:* 6 SINGLE-FLOAT-EPSILON))
                                    which we get by the Taylor series expansion (SIN X)
                                    =
                                    (+ X (/ (EXPT X 3)
                                            6)
                                       ...)
                                    (The relative error caused by ommitting (/ (EXPT X 3)
                                                                               6)
                                         isn't observable.)
                                    Comparison against %%SIN-EPSILON is used to avoid POLYEVAL 
                                    microcode underflow when computing SIN. *)
                                 (INITVARS (%%SIN-EPSILON (%%FLOAT 14720 0)))
                                 (* * %%SIN-PPOLY and %%SIN-QPOLY contain adapted P and Q 
                                    coefficients of Harris et al SIN 3374 rational approximation to
                                    (SIN X)
                                    in interval (0 (/ PI 2))
                                    %. The coefficients for %%SIN-PPOLY and %%SIN-QPOLY have been 
                                    computed from Harris using extended precision routines and the 
                                    relations %%SIN-PPOLY =
                                    (REVERSE (for I from 0 as ENTRY in PS collect
                                                  (/ (CL:* (EXPT (/ 2 PI)
                                                                 (1+ (CL:* 2 I)))
                                                           ENTRY)
                                                     Q0)))
                                    and %%SIN-QPOLY =
                                    (REVERSE (for I from 0 as ENTRY in QS collect
                                                  (/ (CL:* (EXPT (/ 2 PI)
                                                                 (CL:* 2 I))
                                                           ENTRY)
                                                     Q0)))
                                    *)
                                 (VARS (%%SIN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 45236 25611)
                                                                        (%%FLOAT 13589 26148)
                                                                        (%%FLOAT 47286 34797)
                                                                        (%%FLOAT 15295 3306)
                                                                        (%%FLOAT 48666 34805)
                                                                        (%%FLOAT 16256 0))))
                                       (%%SIN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 11384 52865)
                                                                        (%%FLOAT 12553 9550)
                                                                        (%%FLOAT 13604 38385)
                                                                        (%%FLOAT 14593 18841)
                                                                        (%%FLOAT 15489 5549)
                                                                        (%%FLOAT 16256 0))))))
                           (FNS CL:SIN %%SIN-FLOAT %%SIN-COMPLEX))
                     (COMS (* COS *)
                           (FNS CL:COS %%COS-COMPLEX))
                     (COMS (* TAN *)
                           (COMS (* * %%TAN-EPSILON is sufficiently small that (TAN X)
                                    = X for X in interval (0 %%TAN-EPSILON)
                                    %. It suffices to take %%TAN-EPSILON a little bit smaller than
                                    (SQRT (CL:* 3 SINGLE-FLOAT-EPSILON))
                                    which we get by the Taylor series expansion (TAN X)
                                    =
                                    (+ X (/ (EXPT X 3)
                                            3)
                                       ...)
                                    (The relative error caused by ommitting (/ (EXPT X 3)
                                                                               3)
                                         isn't observable.)
                                    Comparison against %%TAN-EPSILON is used to avoid POLYEVAL 
                                    microcode underflow when computing TAN. *)
                                 (INITVARS (%%TAN-EPSILON (%%FLOAT 14720 0)))
                                 (* * %%TAN-PPOLY and %%TAN-QPOLY contain adapted P and Q 
                                    coefficients of Harris et al TAN 4288 rational approximation to
                                    (TAN X)
                                    in interval (-PI/4 PI/4)
                                    %. The coefficients for %%TAN-PPOLY and %%TAN-QPOLY have been 
                                    computed from Harris using extended precision routines and the 
                                    relations %%TAN-PPOLY =
                                    (REVERSE (for I from 0 as ENTRY in PS collect
                                                  (/ (CL:* (EXPT (/ 4 PI)
                                                                 (1+ (CL:* 2 I)))
                                                           ENTRY)
                                                     Q0)))
                                    and %%TAN-QPOLY =
                                    (REVERSE (for I from 0 as ENTRY in QS collect
                                                  (/ (CL:* (EXPT (/ 4 PI)
                                                                 (CL:* 2 I))
                                                           ENTRY)
                                                     Q0)))
                                    *)
                                 (VARS (%%TAN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 13237 21090)
                                                                        (%%FLOAT 47141 15825)
                                                                        (%%FLOAT 15246 8785)
                                                                        (%%FLOAT 48655 48761)
                                                                        (%%FLOAT 16256 0))))
                                       (%%TAN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 45267 36947)
                                                                        (%%FLOAT 13848 46875)
                                                                        (%%FLOAT 47612 53738)
                                                                        (%%FLOAT 15596 52854)
                                                                        (%%FLOAT 48882 35303)
                                                                        (%%FLOAT 16256 0))))))
                           (FNS CL:TAN %%TAN-FLOAT %%TAN-COMPLEX))
                     (COMS (* ASIN *)
                           (COMS (* * %%ASIN-EPSILON is sufficiently small that (ASIN X)
                                    = X for X in interval (0 %%ASIN-EPSILON)
                                    %. It suffices to take %%ASIN-EPSILON a little bit smaller than
                                    (CL:* 2 SINGLE-FLOAT-EPSILON)
                                    which we get by the Taylor series expansion (ASIN X)
                                    =
                                    (+ X (/ (EXPT X 3)
                                            6)
                                       ...)
                                    (The relative error caused by ommitting (/ (EXPT X 3)
                                                                               6)
                                         isn't observable.)
                                    Comparison against %%ASIN-EPSILON is used to avoid POLYEVAL 
                                    microcode underflow when computing SIN. *)
                                 (INITVARS (%%ASIN-EPSILON (%%FLOAT 14720 0)))
                                 (* * %%ASIN-PPOLY and %%ASIN-QPOLY contain P and Q coefficients of 
                                    Harris et al ARCSN 4671 rational approximation to (ASIN X)
                                    in interval (0 (SQRT .5)) . *)
                                 (VARS (%%ASIN-PPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16007 50045)
                                                                         (%%FLOAT 49549 8020)
                                                                         (%%FLOAT 17236 15848)
                                                                         (%%FLOAT 50285 63464)
                                                                         (%%FLOAT 17650 31235)
                                                                         (%%FLOAT 50403 62852)
                                                                         (%%FLOAT 17440 39471))))
                                       (%%ASIN-QPOLY (%%MAKE-ARRAY (LIST (%%FLOAT 16256 0)
                                                                         (%%FLOAT 49672 25817)
                                                                         (%%FLOAT 17308 55260)
                                                                         (%%FLOAT 50326 38098)
                                                                         (%%FLOAT 17674 22210)
                                                                         (%%FLOAT 50417 22451)
                                                                         (%%FLOAT 17440 39471))))))
                           (FNS ASIN %%ASIN-FLOAT %%ASIN-COMPLEX))
                     (COMS (* ACOS *)
                           (FNS ACOS %%ACOS-COMPLEX))
                     (COMS (* ATAN *)
                           (INITVARS (%%SQRT3 (%%FLOAT 16349 46039))
                                  (%%2-SQRT3 (%%FLOAT 16009 12451))
                                  (%%INV-2-SQRT3 (%%FLOAT 16494 55788)))
                           (FNS CL:ATAN %%ATAN-FLOAT1 %%ATAN-FLOAT2 %%ATAN-DOMAIN-CHECK %%ATAN-FLOAT 
                                %%ATAN-COMPLEX))
                     (COMS (* CIS *)
                           (FNS CIS))
                     (COMS (* SINH COSH TANH *)
                           (FNS SINH COSH TANH))
                     (COMS (* ASINH ACOSH ATANH *)
                           (FNS ASINH ACOSH ATANH %%ATANH-DOMAIN-CHECK))
                     (PROP FILETYPE CMLFLOAT)
                     (DECLARE: DONTEVAL@LOAD DOEVAL@COMPILE DONTCOPY COMPILERVARS
                            (ADDVARS (NLAMA)
                                   (NLAML)
                                   (LAMA CL:LOG %%EXPT-INTEGER CL:EXPT)))))
(DECLARE: DONTEVAL@LOAD DOEVAL@COMPILE DONTCOPY COMPILERVARS 

(ADDTOVAR NLAMA )

(ADDTOVAR NLAML )

(ADDTOVAR LAMA CL:LOG %%EXPT-INTEGER CL:EXPT)
)
(PUTPROPS CMLFLOAT COPYRIGHT ("Xerox Corporation" 1986))
(DECLARE: DONTCOPY
  (FILEMAP (NIL (18284 18444 (%%FLOAT 18294 . 18442)) (21256 22110 (%%MAKE-ARRAY 21266 . 22108)) (23219 
25943 (EXP 23229 . 23736) (%%EXP-FLOAT 23738 . 25595) (%%EXP-COMPLEX 25597 . 25941)) (25963 31623 (
CL:EXPT 25973 . 27699) (%%EXPT-INTEGER 27701 . 30435) (%%EXPT-FLOAT 30437 . 30747) (%%EXPT-COMPLEX 
30749 . 31174) (%%EXPT-COMPLEX-POWER 31176 . 31621)) (32502 35036 (CL:LOG 32512 . 33066) (%%LOG-FLOAT 
33068 . 34833) (%%LOG-COMPLEX 34835 . 35034)) (35056 38408 (CL:SQRT 35066 . 35602) (%%SQRT-FLOAT 35604
 . 37136) (%%SQRT-COMPLEX 37138 . 38406)) (40190 43517 (CL:SIN 40200 . 40524) (%%SIN-FLOAT 40526 . 
43155) (%%SIN-COMPLEX 43157 . 43515)) (43536 44237 (CL:COS 43546 . 43865) (%%COS-COMPLEX 43867 . 44235
)) (45961 49053 (CL:TAN 45971 . 46205) (%%TAN-FLOAT 46207 . 48525) (%%TAN-COMPLEX 48527 . 49051)) (
50610 53636 (ASIN 50620 . 50902) (%%ASIN-FLOAT 50904 . 53345) (%%ASIN-COMPLEX 53347 . 53634)) (53656 
54295 (ACOS 53666 . 53983) (%%ACOS-COMPLEX 53985 . 54293)) (54452 57736 (CL:ATAN 54462 . 54772) (
%%ATAN-FLOAT1 54774 . 55692) (%%ATAN-FLOAT2 55694 . 56061) (%%ATAN-DOMAIN-CHECK 56063 . 56565) (
%%ATAN-FLOAT 56567 . 57356) (%%ATAN-COMPLEX 57358 . 57734)) (57755 58241 (CIS 57765 . 58239)) (58271 
59411 (SINH 58281 . 58782) (COSH 58784 . 58978) (TANH 58980 . 59409)) (59444 60736 (ASINH 59454 . 
59675) (ACOSH 59677 . 59898) (ATANH 59900 . 60225) (%%ATANH-DOMAIN-CHECK 60227 . 60734)))))
STOP