(FILECREATED "29-Aug-86 14:32:38" {ERIS}<LISPCORE>SOURCES>LLFLOAT.;47 72460
changes to: (FNS \UNBOXFLOAT3)
(VARS LLFLOATCOMS)
previous date: "23-Jun-86 09:11:45" {ERIS}<LISPCORE>SOURCES>LLFLOAT.;46)
(* Copyright (c) 1982, 1984, 1985, 1986 by Xerox Corporation. All rights reserved.)
(PRETTYCOMPRINT LLFLOATCOMS)
(RPAQQ LLFLOATCOMS
[(DECLARE: DONTCOPY (MACROS \HAND.FLOATUNBOX)
(EXPORT (MACROS POLYEVAL)))
(COMS (FNS \PUTBASEFLOATP \GETBASEFLOATP)
(MACROS \PUTBASEFLOATP \GETBASEFLOATP (* " the following deal with raw 32 bit numbers")
\.PUTBASE32 \.GETBASE32))
[COMS (FNS FTIMES FPLUS FQUOTIENT FDIFFERENCE FGREATERP FABS)
(* UFNs)
(FNS \SLOWFDIFFERENCE \SLOWFPLUS2 \SLOWFTIMES2 \SLOWFQUOTIENT \SLOWFGREATERP)
(FNS \FZEROP FEQP FLOAT \FLOAT \FIXP.FROM.FLOATP FIXR \BOXFPLUSDIF \BOXFQUOTIENT
\BOXFTIMES2 \INFINITY \MAKEFLOAT MAKEFLOATNUMBER PutFloat \FZEROP MAKERATIONAL)
(PROP DMACRO ZEROP)
(FNS SQRT)
(DECLARE: EVAL@COMPILE DONTCOPY (EXPORT (RECORDS FLOATP)
(CONSTANTS (MAX.DIGITS.ACCURACY 9)))
(CONSTANTS (\8BITS 255)
(\MAX.HI.FRAC 127)
(\SIGNBIT 32768)
(\EXPONENT.BIAS 127)
(\HIDDENBIT 128)
(\MAX.EXPONENT 255))
(MACROS .FLOATUNBOX. .LLSH1. .LLSH8. .LRSH1. .LRSH8. .LRSHSTICKY. .ADDSMALL2.
.ADDSMALL3. .SUBSMALL. .POWEROF2.)
(LOCALVARS . T))
(DECLARE: DONTEVAL@LOAD DOCOPY (VARS (\UNDERFLOW)
(MAX.FLOAT (\INFINITY 0))
(MIN.FLOAT (\INFINITY 1)))
(P (MOVD? (QUOTE FGREATERP)
(QUOTE FGTP]
[COMS (* * "unboxed ufns")
(FNS \UNBOXFLOAT1 \UNBOXFLOAT2 \UNBOXFLOAT3)
(FNS \MATMULT133 \MATMULT144 \MATMULT331 \MATMULT333 \MATMULT441 \MATMULT444)
(* unboxed arg handling)
(DECLARE: DONTCOPY (EXPORT (MACROS \CALLER.ARGS]
[COMS (FNS \FLOATINGSCALE \INIT.POWERS.OF.TEN \CONVERT.FLOATING.NUMBER \PRODUCE.FDIGIT)
(DECLARE: DONTCOPY (GLOBALVARS \POWERS.OF.TEN)
(MACROS \POWER.OF.TEN))
(DECLARE: DONTEVAL@LOAD DOCOPY (P (\INIT.POWERS.OF.TEN]
(COMS (FNS FLOATP.TO.BCPL BCPL.TO.FLOATP)
(DECLARE: EVAL@COMPILE DONTCOPY (RECORDS BCPLNUM)))
(PROP ARGNAMES \UNBOXFLOAT1 \UNBOXFLOAT2 \UNBOXFLOAT3)
(DECLARE: DONTEVAL@LOAD DOEVAL@COMPILE DONTCOPY COMPILERVARS (ADDVARS (NLAMA)
(NLAML)
(LAMA FPLUS FTIMES])
(DECLARE: DONTCOPY
(DECLARE: EVAL@COMPILE
[PUTPROPS \HAND.FLOATUNBOX MACRO (LAMBDA (X)
(* this doesn't call \FLOATUNBOX because it's used by the UFN
case of \FLOATUNBOX. - takes a FLOATP and returns the raw
unboxed bits of the value. Must be used with great caution
as raw unboxed bits are not allowed in many places.)
(\VAG2 (fetch (FLOATP HIWORD)
of
(SETQ X (FLOAT X)))
(fetch (FLOATP LOWORD)
of X]
)
(* FOLLOWING DEFINITIONS EXPORTED)
(DECLARE: EVAL@COMPILE
[PUTPROPS POLYEVAL DMACRO ((X COEFFS DEGREE)
(* execute the POLYEVAL opcode on the value X, the array COEFFS with
degree DEGREE)
(\FLOATBOX ((OPCODES UBFLOAT3 0)
(\FLOATUNBOX X)
(fetch (ARRAYP BASE)
of COEFFS)
DEGREE]
)
(* END EXPORTED DEFINITIONS)
)
(DEFINEQ
(\PUTBASEFLOATP
[LAMBDA (BASE OFFST VAL) (* JonL " 7-FEB-83 19:42")
(\MACRO.MX (\PUTBASEFLOATP BASE OFFST VAL))])
(\GETBASEFLOATP
[LAMBDA (BASE OFFST) (* JonL " 7-FEB-83 19:41")
(\MACRO.MX (\GETBASEFLOATP BASE OFFST))])
)
(DECLARE: EVAL@COMPILE
[PUTPROPS \PUTBASEFLOATP DMACRO ((BASE OFFST VAL)
(* put the floatp VAL at offset OFFST from BASE. Used by
REPLACEFIELD of floatp fields)
(\FLOATBOX (\.PUTBASE32 BASE OFFST (\FLOATUNBOX VAL]
[PUTPROPS \GETBASEFLOATP DMACRO ((BASE OFFST)
(* get the floatp at OFFST from BASE)
(\FLOATBOX (\.GETBASE32 BASE OFFST]
(PUTPROPS \.PUTBASE32 DMACRO (= . \PUTBASEPTR))
(PUTPROPS \.GETBASE32 DMACRO (APPLY* COMP.GETBASE NIL GETBASE.32))
)
(DEFINEQ
(FTIMES
(LAMBDA N (* JonL "17-May-84 18:35")
(PROG (R (J 1))
(COND
((EQ 0 N)
(RETURN 1.0))
((EQ N 1)
(RETURN (FLOAT (ARG N 1)))))
(SETQ R (ARG N 1))
LP (COND
((NEQ J N)
(SETQ J (ADD1 J)) (* assumes that FTIMES compiles into opcode that punts
into \FTIMES.UFN)
(SETQ R (FTIMES R (ARG N J)))
(GO LP)))
(RETURN R))))
(FPLUS
(LAMBDA N (* JonL "17-May-84 18:35")
(PROG (R (J 1))
(COND
((EQ 0 N)
(RETURN 0.0))
((EQ N 1)
(RETURN (FLOAT (ARG N 1)))))
(SETQ R (ARG N 1))
LP (COND
((NEQ J N)
(SETQ J (ADD1 J))
(SETQ R (FPLUS R (ARG N J)))
(GO LP)))
(RETURN R))))
(FQUOTIENT
[LAMBDA (X Y) (* lmm "11-FEB-82 14:02")
((OPCODES FQUOTIENT)
X Y])
(FDIFFERENCE
[LAMBDA (X Y) (* lmm "14-MAR-84 22:20")
((OPCODES FDIFFERENCE)
X Y])
(FGREATERP
[LAMBDA (X Y) (* lmm "17-Oct-84 15:45")
(* to compare two floats, compare signbits, and if they
are equal compare the remaining 31 bits of each number
as unsigned integers)
((OPCODES FGREATERP)
X Y])
(FABS
[LAMBDA (X)
(\MACRO.MX (FABS X))])
)
(* UFNs)
(DEFINEQ
(\SLOWFDIFFERENCE
[LAMBDA (X Y) (* lmm "17-Oct-84 15:42")
(\CALLME (QUOTE FDIFFERENCE))
(\BOXFPLUSDIF X Y T])
(\SLOWFPLUS2
[LAMBDA (X Y) (* lmm "17-Oct-84 15:42")
(* UFN for FPLUS)
(\CALLME (QUOTE FPLUS))
(\BOXFPLUSDIF X Y])
(\SLOWFTIMES2
[LAMBDA (X Y) (* lmm "17-Oct-84 15:43")
(\CALLME (QUOTE FTIMES))
(\BOXFTIMES2 X Y])
(\SLOWFQUOTIENT
[LAMBDA (X Y) (* lmm "17-Oct-84 15:43")
(\CALLME (QUOTE FQUOTIENT)) (* UFN for FQUOTIENT)
(\BOXFQUOTIENT X Y NIL])
(\SLOWFGREATERP
[LAMBDA (X Y) (* JonL "17-May-84 18:34")
(* to compare two floats, compare signbits, and if they
are equal compare the remaining 31 bits of each number
as unsigned integers)
(COND
[(AND (FLOATP X)
(FLOATP Y)) (* Can speed this up by not unpacking--check signs,
then compare remaining 31d bits as unsigned numbers)
(PROG ((HX (fetch (FLOATP HIWORD) of X))
(HY (fetch (FLOATP HIWORD) of Y))
SIGNX)
(RETURN (COND
((NEQ (SETQ SIGNX (LOGAND HX \SIGNBIT))
(LOGAND HY \SIGNBIT))
(EQ 0 SIGNX))
[(EQ 0 SIGNX) (* numbers are positive)
(OR (IGREATERP HX HY)
(AND (EQ HX HY)
(IGREATERP (fetch LOWORD of X)
(fetch LOWORD of Y]
(T (* Numbers are negative, so compare in other direction)
(OR (IGREATERP HY HX)
(AND (EQ HX HY)
(IGREATERP (fetch LOWORD of Y)
(fetch LOWORD of X]
(T (PROG (HX LX SIGNX EXPX HY LY SIGNY EXPY)
(.FLOATUNBOX. X SIGNX EXPX HX LX)
(.FLOATUNBOX. Y SIGNY EXPY HY LY)
(RETURN (COND
((NEQ SIGNX SIGNY)
(EQ 0 SIGNX))
[(EQ 0 SIGNX) (* numbers are positive)
(OR (IGREATERP EXPX EXPY)
(AND (EQ EXPX EXPY)
(OR (IGREATERP HX HY)
(AND (EQ HX HY)
(IGREATERP LX LY]
(T (* Numbers are negative, so compare in other direction)
(OR (IGREATERP EXPY EXPX)
(AND (EQ EXPY EXPX)
(OR (IGREATERP HY HX)
(AND (EQ HY HX)
(IGREATERP LY LX])
)
(DEFINEQ
(\FZEROP
(LAMBDA (X) (* JonL "27-Sep-84 22:20")
(* Support for generic ZEROP macro)
(AND (EQ 0 (fetch LOWORD of (\DTEST X (QUOTE FLOATP))))
(EQ 0 (fetch HIWORDNOSIGNBIT of X)))))
(FEQP
(LAMBDA (X Y) (* JonL "17-May-84 20:26")
(COND
((AND (FLOATP X)
(FLOATP Y)) (* If they're both floatp already, can essentially
compare contents, since floatps are generally
normalized.)
(* Last OR clause is to check for comparing a negative
zero to a positive zero.)
(AND (EQ (fetch LOWORD of X)
(fetch LOWORD of Y))
(OR (EQ (fetch HIWORD of X)
(fetch HIWORD of Y))
(AND (EQ 0 (fetch HIWORDNOSIGNBIT of X))
(EQ 0 (fetch HIWORDNOSIGNBIT of Y))))))
(T (PROG (SIGNX EXPX HX LX SIGNY EXPY HY LY)
(.FLOATUNBOX. X SIGNX EXPX HX LX T)
(.FLOATUNBOX. Y SIGNY EXPY HY LY T)
(RETURN (AND (EQ HX HY)
(EQ LX LY)
(EQ EXPX EXPY)
(OR (EQ SIGNX SIGNY)
(AND (EQ 0 EXPX)
(EQ 0 HX)
(EQ 0 LX))))))))))
(FLOAT
[LAMBDA (X) (* lmm "19-Jun-86 15:18")
(* common lisp says (FLOAT N &OPTIONAL
...) but since Interlisp only has one
float precision, it can just ignore
the rests)
(* compiles this way, too)
(\DTEST X (QUOTE FLOATP])
(\FLOAT
[LAMBDA (X) (* lmm "23-Jun-86 09:04")
(CTYPECASE X (FLOAT X)
[INTEGER (ETYPECASE X ((DATATYPE SMALLP)
(LET* [(HI 0)
(SIGN 0)
(LO (COND
((IGEQ X 0)
X)
(T (SETQ SIGN 1)
(* X is negative--negate it)
(COND
((EQ 0 (\LOLOC X))
(* Min small integer)
(SETQ HI 1)
0)
(T (ADD1 (IDIFFERENCE MAX.SMALL.INTEGER
(\LOLOC X]
(\MAKEFLOAT SIGN (IPLUS \EXPONENT.BIAS 31)
HI LO T)))
((DATATYPE FIXP)
(LET ((HI (fetch (FIXP HINUM) of X))
(LO (fetch (FIXP LONUM) of X))
(SIGN 0))
(COND
((IGREATERP HI MAX.POS.HINUM)
(.NEGATE. HI LO)
(SETQ SIGN 1)))
(\MAKEFLOAT SIGN (IPLUS \EXPONENT.BIAS 31)
HI LO T)))
((DATATYPE BIGNUM)
(\BIGNUM.TO.FLOAT X]
(RATIO (FQUOTIENT (RATIO-NUMERATOR X)
(RATIO-DENOMINATOR X])
(\FIXP.FROM.FLOATP
[LAMBDA (X) (* lmm "19-Sep-85 15:08")
(PROG (SIGN EXP HI LO)
(.FLOATUNBOX. X SIGN EXP HI LO (GO RETZERO))
(SETQ EXP (IDIFFERENCE EXP (SUB1 \EXPONENT.BIAS)))
(* number of bits to left of binary point)
[COND
((ILESSP EXP 0)
(RETURN 0))
([OR (IGREATERP EXP 32)
(AND (EQ EXP 32)
(OR (EQ 0 SIGN)
(NEQ HI \SIGNBIT)
(NEQ LO 0]
(RETURN (LSH (LET ((VAL (PLUS (LSH HI 16)
LO)))
(if (EQ SIGN 1)
then (MINUS VAL)
else VAL))
(DIFFERENCE EXP 32]
[COND
((IGEQ (SETQ EXP (IDIFFERENCE 32 EXP))
16)
(SETQ LO (LRSH HI (IDIFFERENCE EXP 16)))
(SETQ HI 0))
(T (* large integer, have to manipulate both halves)
(FRPTQ EXP (.LRSH1. HI LO]
(COND
((EQ SIGN 1)
(.NEGATE. HI LO)))
(RETURN (\MAKENUMBER HI LO))
RETZERO
(RETURN 0])
(FIXR
[LAMBDA (X) (* lmm "22-JUL-84 20:48")
(OR (FIXP X)
(PROG (SIGN EXP HI LO ROUNDINGBITS)
(.FLOATUNBOX. X SIGN EXP HI LO (GO RETZERO))
(SETQ EXP (IDIFFERENCE EXP (SUB1 \EXPONENT.BIAS)))
(* number of bits to left of binary point)
[COND
((ILESSP EXP 0)
(RETURN 0))
[(IGEQ EXP 40Q) (* FIX handles this)
(RETURN (FIX (FPLUS X .5]
([OR (IGREATERP EXP 40Q)
(AND (EQ EXP 40Q)
(OR (EQ 0 SIGN)
(NEQ HI \SIGNBIT)
(NEQ LO 0] (* Overflow: number is larger than MAX.FIXP)
(RETURN (SELECTQ \OVERFLOW
(T (LISPERROR "OVERFLOW" X T))
(COND
((EQ 0 SIGN)
MAX.FIXP)
(T MIN.FIXP]
[COND
((IGEQ EXP 30Q) (* No decimal places to worry about, so no rounding,
just shift into place)
(FRPTQ (IDIFFERENCE 40Q EXP)
(.LRSH1. HI LO)))
(T
(* Shift right until binary point is in the middle of LO, as per \MAKEFLOAT; then decide how to round, and shift
right once more)
[COND
((IGEQ (SETQ EXP (IDIFFERENCE 30Q EXP))
20Q) (* shifting all the way out of the high word)
(SETQ LO (LRSH [LOGOR HI (COND
((EQ 0 LO)
0)
(T (* Sticky bits)
(LRSH \8BITS 1]
(IDIFFERENCE EXP 20Q)))
(SETQ HI 0))
(T (* Shift both halves, keeping sticky bits in LO)
(FRPTQ EXP (.LRSHSTICKY. HI LO]
(SETQ ROUNDINGBITS (LOGAND LO \8BITS))
(.LRSH8. HI LO) (* Shift the rest of the way)
(COND
((OR (IGREATERP ROUNDINGBITS 200Q)
(AND (EQ ROUNDINGBITS 200Q)
(ODDP LO))) (* Round up if greater than .5, or exactly .5 and
rounding up will make number even)
(COND
((EQ LO MAX.SMALL.INTEGER)
(SETQ LO 0)
(add HI 1))
(T (add LO 1]
(COND
((EQ SIGN 1)
(.NEGATE. HI LO)))
(RETURN (\MAKENUMBER HI LO))
RETZERO
(RETURN 0])
(\BOXFPLUSDIF
(LAMBDA (X Y SUBTRACT BOX) (* JonL "17-May-84 18:56")
(* Does X-Y if SUBTRACT is true)
(PROG (SIGNX EXPX HX LX SIGNY EXPY HY LY EXPDIFF PLEASENORMALIZE CARRY)
(.FLOATUNBOX. Y SIGNY EXPY HY LY)
(COND
(SUBTRACT (SETQ SIGNY (IDIFFERENCE 1 SIGNY))))
(.FLOATUNBOX. X SIGNX EXPX HX LX (GO RESULTISY))
(COND
((AND (EQ 0 HY)
(EQ 0 LY))
(GO DONE))
((EQ EXPX \MAX.EXPONENT) (* X = infinity, so result is infinity.
This is not quite right if Y is infinity of opposite
sign, though)
(RETURN (\INFINITY SIGNX BOX)))
((EQ EXPY \MAX.EXPONENT)
(RETURN (\INFINITY SIGNY BOX))))
(SETQ EXPDIFF (IDIFFERENCE EXPX EXPY)) (* first align the binary points by right-shifting the
smaller guy)
(COND
((IGREATERP EXPDIFF 0)
(COND
((IGREATERP EXPDIFF 31) (* Y would get shifted into oblivion)
(GO DONE))
(T (FRPTQ EXPDIFF (.LRSHSTICKY. HY LY)))))
((NEQ EXPDIFF 0)
(COND
((ILESSP EXPDIFF -31)
(GO RESULTISY))
(T (FRPTQ (IMINUS EXPDIFF)
(.LRSHSTICKY. HX LX))
(SETQ EXPX EXPY)))))
(COND
((EQ SIGNX SIGNY) (* same sign, add magnitudes)
(SETQ CARRY (.ADDSMALL2. LX LY))
(COND
((EQ (.ADDSMALL3. HX HY CARRY)
1) (* there was a carry out of HX, so shift everyone right
and stick it back in)
(.LRSHSTICKY. HX LX)
(add HX \SIGNBIT)
(add EXPX 1))))
(T (* subtract magnitudes, smaller from larger)
(COND
((OR (ILESSP HX HY)
(AND (EQ HX HY)
(ILESSP LX LY))) (* Y is bigger, so swap)
(swap HX HY)
(swap LX LY)
(SETQ SIGNX SIGNY)))
(SETQ PLEASENORMALIZE (NEQ (LOGAND HX \SIGNBIT)
0)) (* thus if neither operand is normalized, we won't waste
time normalizing and denormalizing the result)
(SETQ HX (IDIFFERENCE (IDIFFERENCE HX HY)
(.SUBSMALL. LX LY)))))
DONE(RETURN (\MAKEFLOAT SIGNX EXPX HX LX PLEASENORMALIZE BOX))
RESULTISY
(RETURN (\MAKEFLOAT SIGNY EXPY HY LY NIL BOX)))))
(\BOXFQUOTIENT
[LAMBDA (X Y BOX) (* lmm "18-DEC-80 13:40")
(PROG (SIGNX EXPX HX LX (SIGNY 0)
(EXPY 0)
HY LY BORROW (HZ 0)
(LZ 0))
(.FLOATUNBOX. X SIGNX EXPX HX LX (GO DONE))
(.FLOATUNBOX. Y SIGNY EXPY HY LY (GO DIVZERO))
(COND
((EQ EXPX \MAX.EXPONENT)
(* X is infinity)
(RETURN (\INFINITY SIGNX BOX)))
((EQ EXPY \MAX.EXPONENT)
(* Y = infinity, result is zero)
(GO DONE)))
(* * Divide X -- double length, implicitly extended with zeros -- by Y. At each step, Y is subtracted from X if
possible, putting a one bit in the quotient, and then X and the quotient are shifted left. Result is a 32-bit
quotient.)
(.LRSH1. HX LX)
(.LRSH1. HY LY) (* shift these right one so that we never have to worry about carrying out of the
high bit)
(FRPTQ 31 (PROGN (.LLSH1. HZ LZ)
(* shift quotient left one as we accumulate it)
(COND
((OR (AND (EQ HX HY)
(IGEQ LX LY))
(IGREATERP HX HY))
(* X GE Y, so subtract Y)
(SETQ HX (IDIFFERENCE (IDIFFERENCE HX HY)
(.SUBSMALL. LX LY)))
(SETQ LZ (ADD1 LZ))
(* note that this never overflows, because of the left shift we did above)
))
(* now shift dividend left one. After the subtraction the high-order bit must be
off, so this works okay)
(.LLSH1. HX LX)))
(.LLSH1. HZ LZ) (* left shift result 1 to compensate for the earlier right shifts)
[COND
((OR (NEQ HX 0)
(NEQ LX 0)) (* set sticky bit)
(SETQ LZ (LOGOR LZ 1]
DONE(RETURN (\MAKEFLOAT (LOGXOR SIGNX SIGNY)
(IPLUS (IDIFFERENCE EXPX EXPY)
\EXPONENT.BIAS)
HZ LZ T BOX))
DIVZERO
(RETURN (COND
((EQ \OVERFLOW T)
(ERROR "FLOATING DIVIDE BY ZERO" Y))
(T (\INFINITY SIGNX BOX])
(\BOXFTIMES2
(LAMBDA (X Y BOX) (* JonL "17-May-84 18:56")
(PROG (SIGNX EXPX HX LX (SIGNY 0)
(EXPY 0)
HY LY (HHY 0)
(HHZ 0)
(HZ 0)
(LZ 0)
SAVEHY SAVELY CARRY)
(.FLOATUNBOX. X SIGNX EXPX HX LX (GO DONE)
T)
(.FLOATUNBOX. Y SIGNY EXPY HY LY (GO DONE)
T)
(COND
((EQ EXPX \MAX.EXPONENT) (* X = infinity)
(RETURN (\INFINITY SIGNX BOX)))
((EQ EXPY \MAX.EXPONENT)
(RETURN (\INFINITY SIGNY BOX))))
(* * Multiply the significands. We have two 24-bit integers, so have a 48-bit, 3-word product, stored as
{HHZ,HZ,LZ}. Multiplication will be in two steps: multiply LX by {HY,LY}, storing in result, and then multiply HX
by {HY,LY}, storing in the top two words. The first multiplication can be omitted in the not uncommon case of a
zero low fraction, and the second multiplication is a little bit simpler, since result fits in two words.)
(COND
((EQ 0 LX)
(GO LP2))
((EQ 0 LY) (* swap operands to make life easier)
(swap HX HY)
(swap LX LY)
(GO LP2)))
(SETQ SAVEHY HY) (* we'll need these for second step)
(SETQ SAVELY LY)
LP1 (* multiply LX times HY,LY)
(COND
((NEQ (LOGAND LX 1)
0)
(SETQ CARRY (.ADDSMALL2. LZ LY))
(SETQ CARRY (.ADDSMALL3. HZ HY CARRY))
(SETQ HHZ (IPLUS HHZ HHY CARRY))))
(COND
((EQ 0 (SETQ LX (LRSH LX 1))) (* done with this step)
(SETQ HY SAVEHY)
(SETQ LY SAVELY)
(GO LP2)))
(SETQ HHY (LLSH HHY 1)) (* left shift Y by one)
(SETQ HY (LLSH (COND
((IGREATERP HY MAX.POS.HINUM)
(add HHY 1)
(LOGAND HY MAX.POS.HINUM))
(T HY))
1))
(SETQ LY (LLSH (COND
((IGREATERP LY MAX.POS.HINUM)
(add HY 1)
(LOGAND LY MAX.POS.HINUM))
(T LY))
1))
(GO LP1)
LP2
(* multiply HX times HY,LY, adding into high two words of Z. No overflow here, since HX has at most
(and usually exactly) 8 bits)
(COND
((NEQ (LOGAND HX 1)
0)
(SETQ CARRY (.ADDSMALL2. HZ LY))
(SETQ HHZ (IPLUS HHZ HY CARRY))))
(COND
((NEQ (SETQ HX (LRSH HX 1))
0)
(.LLSH1. HY LY)
(GO LP2)))
DONE
(* * We now have a 48-bit result in HHZ,HZ,LZ. \MAKEFLOAT can handle it from here. Note that the exponent we give
is bumped by 1, because the "binary point", which was between the first and second bits, was moved one to the
right by multiplying)
(RETURN (\MAKEFLOAT (LOGXOR SIGNX SIGNY)
(IPLUS EXPX EXPY (IDIFFERENCE 1 \EXPONENT.BIAS))
HHZ HZ T BOX)))))
(\INFINITY
[LAMBDA (SIGN BOX) (* lmm "17-DEC-80 20:32")
(* Returns "infinity" of the appropriate SIGN
(0 or 1), reusing floating BOX if given)
(* * For now, don't return true infinity, but rather the largest representable finite number, so that miscellaneous
floating-point routines don't die)
(OR (FLOATP BOX)
(SETQ BOX (create FLOATP)))
(replace (FLOATP SIGNBIT) of BOX with SIGN)
(replace (FLOATP EXPONENT) of BOX with (SUB1 \MAX.EXPONENT))
(replace (FLOATP HIFRACTION) of BOX with \MAX.HI.FRAC)
(replace (FLOATP LOFRACTION) of BOX with 65535)
BOX])
(\MAKEFLOAT
(LAMBDA (SIGN EXP HI LO NORMALIZE BOX) (* JonL "17-May-84 18:56")
(* * packs up the pieces of a floating point result into a single number box, n the process checking for
underflow, rounding, checking overflow. BOX is optional box to reuse. NORMALIZE is true if we should normalize the
result first (make sign bit of HI 1); otherwise we assume result is already normalized)
(PROG (ROUNDINGBITS)
(OR (FLOATP BOX)
(SETQ BOX (create FLOATP)))
TOP (COND
((AND (EQ 0 HI)
(EQ 0 LO))
(replace HIWORD of BOX with (replace LOWORD of BOX with 0))
(RETURN BOX)))
(COND
(NORMALIZE (COND
((EQ 0 HI)
(SETQ HI LO)
(SETQ LO 0)
(SETQ EXP (IDIFFERENCE EXP 16))))
(while (EQ 0 (LOGAND HI \SIGNBIT))
do (.LLSH1. HI LO)
(SETQ EXP (SUB1 EXP)))))
(COND
((ILEQ EXP 0) (* underflow. Scale by 2↑\Exponentbias in order to
deliver a useful value to the error handler)
(SELECTQ \UNDERFLOW
(T (RETURN (LISPERROR "FLOATING UNDERFLOW" (\MAKEFLOAT SIGN (IPLUS EXP
\EXPONENT.BIAS)
HI LO NIL BOX)
T)))
NIL) (* If we have to return a result, we must "denormalize"
this number. This gives us a little more time before
vanishing to zero)
(COND
((ILESSP EXP -24) (* too small even as denormalized number)
(SETQ HI (SETQ LO 0))
(GO TOP))
(T (* denormalize by shifting right until the exponent is
logically 1; final result will have exponent zero,
hidden bit zero)
(FRPTQ (IDIFFERENCE 1 EXP)
(.LRSHSTICKY. HI LO))
(SETQ EXP 0)))))
(SETQ ROUNDINGBITS (LOGAND LO \8BITS)) (* round result. low order 8 bits are used for rounding)
(.LRSH8. HI LO)
(COND
((OR (IGREATERP ROUNDINGBITS 128)
(AND (EQ ROUNDINGBITS 128)
(NOT (EQ 0 (LOGAND LO 1))))) (* round up if the left over fraction was greater than
1/2; if it was equal to a half, round to the even
result)
(COND
((EQ LO MAX.SMALL.INTEGER) (* can't add 1 directly)
(SETQ LO 0)
(SETQ HI (ADD1 HI))
(COND
((IGREATERP HI (LOGOR \HIDDENBIT \MAX.HI.FRAC))
(* "1.11111--" became "10.000--")
(SETQ HI (LRSH HI 1))
(add EXP 1))))
(T (SETQ LO (ADD1 LO))))))
(COND
((AND (EQ HI 0)
(EQ LO 0)) (* result is zero. This could have snuck in if we
denormalized a number that didn't have enough digits to
survive)
(GO TOP))
((IGEQ EXP \MAX.EXPONENT)
(* overflow. If trap enabled, wrap the exponent around to middle of range (divide by 2↑\exponentbias) to provide a
number of possible use to error handler)
(SELECTQ \OVERFLOW
(T (RETURN (LISPERROR "FLOATING OVERFLOW" (\MAKEFLOAT SIGN (IDIFFERENCE EXP
\EXPONENT.BIAS)
HI LO NIL BOX)
T)))
NIL)
(RETURN (\INFINITY SIGN BOX))))
(replace SIGNBIT of BOX with SIGN)
(replace EXPONENT of BOX with EXP)
(replace HIFRACTION of BOX with HI)
(replace LOFRACTION of BOX with LO)
(RETURN BOX))))
(MAKEFLOATNUMBER
[LAMBDA (N0 N1) (* lmm "12-Apr-85 18:19")
(* CALLED FROM FETCHFIELD)
(LET [(VAL (NCREATE (QUOTE FLOATP]
(replace (FLOATP HIWORD) of VAL with N0)
(replace (FLOATP LOWORD) of VAL with N1)
VAL])
(PutFloat
[LAMBDA (PTR N) (* lmm "29-Dec-84 11:32")
(* used by REPLACEFIELD)
(\PUTBASEFLOATP PTR 0 N)
N])
(\FZEROP
(LAMBDA (X) (* JonL "27-Sep-84 22:20")
(* Support for generic ZEROP macro)
(AND (EQ 0 (fetch LOWORD of (\DTEST X (QUOTE FLOATP))))
(EQ 0 (fetch HIWORDNOSIGNBIT of X)))))
(MAKERATIONAL
[LAMBDA (X EPS) (* lmm "24-Sep-85 16:10")
(if (OR (FIXP X)
(type? RATIONAL X))
then (* Integers are already rational)
X
elseif (FLOATP X)
then [OR EPS (SETQ EPS (CONSTANT (EXPT 2.0 -24]
[if (ZEROP X)
then 0
else (LET (A Y)
(LET ([SIGN (if (MINUSP X)
then (SETQ X (MINUS X]
(XX X)
(NUM (SETQ A (FIX X)))
(DEN 1)
(ONUM 1)
(ODEN 0))
(eachtime (PSETQ XX (SETQ Y (QUOTIENT 1.0
(LET ((DIF (DIFFERENCE XX A)))
(if (ZEROP DIF)
then (GO $$OUT))
DIF)))
NUM
(PLUS (TIMES (SETQ A (FIX Y))
NUM)
ONUM)
DEN
(PLUS (TIMES A DEN)
ODEN)
ONUM NUM ODEN DEN)
until (AND (NOT (ZEROP DEN))
(NOT (GREATERP (ABS (QUOTIENT (DIFFERENCE X
(FQUOTIENT
NUM DEN))
X))
EPS)))
finally (RETURN (\RATIONALIZE (if SIGN
then (MINUS NUM)
else NUM)
DEN]
else (LISPERROR "NON-NUMERIC ARG" X])
)
(PUTPROPS ZEROP DMACRO [OPENLAMBDA (X)
(COND ((EQ X 0))
((FLOATP X)
(\FZEROP X])
(DEFINEQ
(SQRT
[LAMBDA (N) (* lmm "24-Jan-85 19:13")
(PROG ((X (FLOAT N))
V)
(DECLARE (TYPE FLOATP X V))
(if (FLESSP X 0.0)
then (ERROR "SQRT OF NEGATIVE VALUE" N)
elseif (NOT (FGREATERP X 0.0))
then (* 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 (FTIMES .5 (FPLUS V (FQUOTIENT X V]
(RETURN V])
)
(DECLARE: EVAL@COMPILE DONTCOPY
(* FOLLOWING DEFINITIONS EXPORTED)
[DECLARE: EVAL@COMPILE
(BLOCKRECORD FLOATP ((SIGNBIT BITS 1)
(EXPONENT BITS 8)
(HIFRACTION BITS 7)
(LOFRACTION BITS 16))
(BLOCKRECORD FLOATP ((HIWORD WORD)
(LOWORD WORD)))
(BLOCKRECORD FLOATP ((NIL BITS 9)
(LONGFRACTION BITS 23)))
(BLOCKRECORD FLOATP ((FLOATCONTENTS BITS 32)))
(BLOCKRECORD FLOATP ((NIL BITS 1)
(HIWORDNOSIGNBIT BITS 15)))
(CREATE (\FLOATBOX (\VAG2 (LOGOR (LLSH SIGNBIT (PLUS 7 8))
(LLSH EXPONENT 7)
HIFRACTION)
LOFRACTION)))
LOFRACTION ← 0 HIFRACTION ← 0 EXPONENT ← 0 SIGNBIT ← 0
[ACCESSFNS FLOATP ((EXP (LOGAND (LRSH (\HILOC (\FLOATUNBOX DATUM))
7)
255))
(HIFRAC (LOGAND (\HILOC (\FLOATUNBOX DATUM))
127])
]
(DECLARE: EVAL@COMPILE
(RPAQQ MAX.DIGITS.ACCURACY 9)
(CONSTANTS (MAX.DIGITS.ACCURACY 9))
)
(* END EXPORTED DEFINITIONS)
(DECLARE: EVAL@COMPILE
(RPAQQ \8BITS 255)
(RPAQQ \MAX.HI.FRAC 127)
(RPAQQ \SIGNBIT 32768)
(RPAQQ \EXPONENT.BIAS 127)
(RPAQQ \HIDDENBIT 128)
(RPAQQ \MAX.EXPONENT 255)
(CONSTANTS (\8BITS 255)
(\MAX.HI.FRAC 127)
(\SIGNBIT 32768)
(\EXPONENT.BIAS 127)
(\HIDDENBIT 128)
(\MAX.EXPONENT 255))
)
(DECLARE: EVAL@COMPILE
[PUTPROPS .FLOATUNBOX. MACRO
((FLONUM SIGN EXP HI LO ZEROFORM DONTSHIFT RESTARTIFINTEGER)
(* Unpacks a floating point number FLONUM into its components. ZEROFORM is evaluated if the
number is true zero. The fraction is unpacked into HI and LO, with the binary point
implicitly between bits 0 and 1 of HI. If DONTSHIFT is true, the fraction is left in its
original state, with 8 bits in HI and 16 in LO. If FLONUM is not floating, it is coerced.)
(PROG NIL RETRY [COND ((NOT (FLOATP FLONUM))
(* Float and normalize the non-floatp)
(COND
((QUOTE RESTARTIFINTEGER)
(SETQ FLONUM (LISPERROR "NON-NUMERIC ARG" FLONUM T))
(GO RESTARTIFINTEGER))
(T (SELECTC (NTYPX FLONUM)
[\FIXP (SETQ HI (fetch (FIXP HINUM)
of FLONUM))
(SETQ LO (fetch (FIXP LONUM)
of FLONUM))
(SETQ SIGN (COND ((IGREATERP HI MAX.POS.HINUM)
(.NEGATE. HI LO)
1)
(T 0]
[\SMALLP (SETQ HI 0)
(SETQ LO
(COND
((SMALLPOSP FLONUM)
(SETQ SIGN 0)
FLONUM)
(T (SETQ SIGN 1)
(* FLONUM is negative--negate it)
(COND ((EQ 0 (LOLOC FLONUM))
(* Min small integer)
(SETQ HI 1)
0)
(T (ADD1 (IDIFFERENCE
MAX.SMALL.INTEGER
(LOLOC FLONUM]
(PROGN (SETQ FLONUM (FLOAT FLONUM))
(GO RETRY)))
[COND [(EQ 0 HI)
(COND ((EQ 0 LO)
(SETQ EXP 0)
(PROGN ZEROFORM (RETURN)))
(T (SETQ HI LO)
(SETQ LO 0)
(SETQ EXP (IPLUS \EXPONENT.BIAS 15]
((IGREATERP HI 255)
(* Not exact, punt)
(SETQ FLONUM (FLOAT FLONUM))
(GO UNPACK))
(T (SETQ EXP (IPLUS \EXPONENT.BIAS 31]
[COND ((ILEQ HI 255)
(* Do a big shift first.)
(.LLSH8. HI LO)
(SETQ EXP (IDIFFERENCE EXP 8]
(while (EQ 0 (LOGAND HI \SIGNBIT))
do
(.LLSH1. HI LO)
(SETQ EXP (SUB1 EXP)))
(COND (DONTSHIFT (.LRSH8. HI LO)))
(RETURN]
UNPACK
(SETQ SIGN (fetch (FLOATP SIGNBIT)
of FLONUM))
(SETQ LO (fetch (FLOATP LOFRACTION)
of FLONUM))
(SETQ HI (fetch (FLOATP HIFRACTION)
of FLONUM))
[COND [(EQ 0 (SETQ EXP (fetch (FLOATP EXPONENT)
of FLONUM)))
(* zero or a de-normalized number from underflow)
(COND ((AND (EQ 0 HI)
(EQ 0 LO))
(* A zero, regardless of the sign bit zero)
ZEROFORM)
(T (* need bias adjust to account for lack of hidden bit)
(SETQ EXP 1]
((NEQ EXP \MAX.EXPONENT)
(* might want to check for NaN's here if EXP = \MAX.EXPONENT)
(* OR in the implicit high bit of fraction)
(SETQ HI (IPLUS HI \HIDDENBIT]
(COND ((NOT DONTSHIFT)
(.LLSH8. HI LO]
[PUTPROPS .LLSH1. MACRO ((HI LO)
(* shift the pair left one, assuming no overflow)
(SETQ HI (LLSH HI 1))
(SETQ LO (LLSH (COND ((IGREATERP LO MAX.POS.HINUM)
(add HI 1)
(LOGAND LO MAX.POS.HINUM))
(T LO))
1]
[PUTPROPS .LLSH8. MACRO ((HI LO)
(* shift pair left 8, assuming no overflow)
(SETQ HI (IPLUS (LLSH HI 8)
(LRSH LO 8)))
(SETQ LO (LLSH (LOGAND LO \8BITS)
8]
[PUTPROPS .LRSH1. MACRO ((HI LO)
(SETQ LO (LRSH LO 1))
[COND ((NEQ (LOGAND HI 1)
0)
(SETQ LO (IPLUS LO \SIGNBIT]
(SETQ HI (LRSH HI 1]
[PUTPROPS .LRSH8. MACRO ((HI LO)
(SETQ LO (IPLUS (LRSH LO 8)
(LLSH (LOGAND HI \8BITS)
8)))
(SETQ HI (LRSH HI 8]
[PUTPROPS .LRSHSTICKY. MACRO ((HI LO)
(* shifts pair right one, but low-order bit is sticky -- if it ever
becomes 1, it stays 1)
(SETQ LO (LOGOR (LRSH LO 1)
(LOGAND LO 1)))
[COND ((NEQ (LOGAND HI 1)
0)
(SETQ LO (IPLUS LO \SIGNBIT]
(SETQ HI (LRSH HI 1]
[PUTPROPS .ADDSMALL2. MACRO ((X Y)
(PROGN (* does X ← X+Y, returning the carry bit)
(COND ((IGREATERP X (IDIFFERENCE MAX.SMALL.INTEGER Y))
[SETQ X (IDIFFERENCE X (IDIFFERENCE MAX.SMALL.INTEGER
(SUB1 Y]
1)
(T (SETQ X (IPLUS X Y))
0]
[PUTPROPS .ADDSMALL3. MACRO ((X Y CARRY)
(PROGN (* X ← X+Y+CARRY, returning the new carry bit)
(COND ((IGREATERP X (IDIFFERENCE (IDIFFERENCE MAX.SMALL.INTEGER Y
)
CARRY))
(SETQ X (IDIFFERENCE
X
(IDIFFERENCE [IDIFFERENCE
MAX.SMALL.INTEGER
(SUB1 (COND ((EQ Y 0)
(PROG1 CARRY
(SETQ CARRY 0))
)
(T Y]
CARRY)))
1)
(T (SETQ X (IPLUS X Y CARRY))
0]
[PUTPROPS .SUBSMALL. MACRO ((X Y)
(* Subtract Y from X, returning the borrow out of the next word)
(COND ((ILEQ Y X)
(SETQ X (IDIFFERENCE X Y))
0)
(T [SETQ X (ADD1 (IDIFFERENCE MAX.SMALL.INTEGER (IDIFFERENCE Y X]
1]
[PUTPROPS .POWEROF2. MACRO (OPENLAMBDA (X)
(COND ((ILESSP X 16)
(LSH 1 X))
(T (LSH (LSH 1 (IDIFFERENCE X 16))
16]
)
(DECLARE: DOEVAL@COMPILE DONTCOPY
(LOCALVARS . T)
)
)
(DECLARE: DONTEVAL@LOAD DOCOPY
(RPAQQ \UNDERFLOW NIL)
(RPAQ MAX.FLOAT (\INFINITY 0))
(RPAQ MIN.FLOAT (\INFINITY 1))
(MOVD? (QUOTE FGREATERP)
(QUOTE FGTP))
)
(* * "unboxed ufns")
(DEFINEQ
(\UNBOXFLOAT1
[LAMBDA (OP) (* lmm
" 5-Mar-86 11:35")
(* UFN for the
unboxed floating 1-arg
cases)
(\SLOWRETURN)
(SELECTQ OP
(0 (* BOX)
(\CALLER.ARGS (X)
(LET [(VAL (NCREATE (QUOTE FLOATP]
(replace (FLOATP HIWORD) of VAL with (\HILOC X))
(replace (FLOATP LOWORD) of VAL with (\LOLOC X))
VAL)))
(1 (* UNBOX)
(\CALLER.ARGS (X)
(\HAND.FLOATUNBOX X)))
(2 (* UFABS)
(\CALLER.ARGS ((X FLOATP))
(\FLOATUNBOX (ABS X))))
(3 (* UFNEGATE)
(\CALLER.ARGS ((X FLOATP))
(\FLOATUNBOX (FMINUS X))))
(4 (* UFIX)
(\CALLER.ARGS ((X FLOATP))
(FIX X)))
(HELP "\UNBOXFLOAT1 called with illegal op " OP])
(\UNBOXFLOAT2
[LAMBDA (OP) (* lmm " 7-Mar-85 16:48")
(* UFN for the 2-arg floating cases)
(\CALLER.ARGS ((X FLOATP)
(Y FLOATP))
(SELECTQ OP
(0 (* UFADD)
(\HAND.FLOATUNBOX (FPLUS X Y)))
(1 (* UFSUB)
(\HAND.FLOATUNBOX (FDIFFERENCE X Y)))
(2 (* UFISUB)
(\HAND.FLOATUNBOX (FDIFFERENCE Y X)))
(3 (* UFMULT)
(\HAND.FLOATUNBOX (FTIMES X Y)))
(4 (* UFDIV)
(\HAND.FLOATUNBOX (FQUOTIENT X Y)))
(5 (* UFGREAT)
(FGREATERP X Y))
(6 (* UFMAX)
(\HAND.FLOATUNBOX (FMAX X Y)))
(7 (* UFMIN)
(\HAND.FLOATUNBOX (FMIN X Y)))
(8 (* UFREM)
(\HAND.FLOATUNBOX (FREMAINDER X Y)))
(HELP "\UNBOXFLOAT2 called with illegal op " OP])
(\UNBOXFLOAT3
[LAMBDA (OP) (* jop: "29-Aug-86 14:30")
(if (EQ 0 OP)
then [\CALLER.ARGS ((X FLOATP)
COEFFICIENTS DEGREE) (* "Polynomial evaluation")
(bind (RESULT ← (\GETBASEFLOATP COEFFICIENTS 0)) declare (TYPE FLOATP RESULT)
for I from 2 to (LLSH DEGREE 1) by 2 do (SETQ RESULT (FPLUS (FTIMES RESULT X)
(\GETBASEFLOATP
COEFFICIENTS I)))
finally (RETURN (\FLOATUNBOX RESULT]
else (\CALLER.ARGS (MATRIX1 MATRIX2 RESULT)
(SELECTQ OP
(1 (* "3 x 3 matrix multiply")
(\MATMULT333 MATRIX1 MATRIX2 RESULT))
(2 (* "4 x 4 matrix multiply")
(\MATMULT444 MATRIX1 MATRIX2 RESULT))
(3 (* "(1,3) * (3,3) => (1,3)")
(\MATMULT133 MATRIX1 MATRIX2 RESULT))
(4 (* "(3,3) * (3,1) => (3,1)")
(\MATMULT331 MATRIX1 MATRIX2 RESULT))
(5 (* "(1,4) * (4,4) => (1,4)")
(\MATMULT144 MATRIX1 MATRIX2 RESULT))
(6 (* "(4,4) * (4,1) => (4,1)")
(\MATMULT441 MATRIX1 MATRIX2 RESULT))
(HELP "\UNBOXFLOAT3 called with illegal op " OP])
)
(DEFINEQ
(\MATMULT133
[LAMBDA (ABASE BBASE CBASE) (* jop: "29-Aug-86 12:20")
(* * "Multiply a 3 vector times a 3 by 3 array")
(for K from 0 to 4 by 2
do (bind (PRODUCT ← 0.0) declare (TYPE FLOATP PRODUCT) for I from 0 to 4 by 2 as J
from K by 6 do [SETQ PRODUCT (FPLUS PRODUCT (FTIMES (\GETBASEFLOATP ABASE I)
(\GETBASEFLOATP BBASE J]
finally (\PUTBASEFLOATP CBASE K PRODUCT)))
CBASE])
(\MATMULT144
[LAMBDA (ABASE BBASE CBASE) (* jop: "29-Aug-86 12:20")
(* * "Multiply a 4 vector times a 4 by 4 array")
(for K from 0 to 6 by 2
do (bind (PRODUCT ← 0.0) declare (TYPE FLOATP PRODUCT) for I from 0 to 6 by 2 as J
from K by 8 do [SETQ PRODUCT (FPLUS PRODUCT (FTIMES (\GETBASEFLOATP ABASE I)
(\GETBASEFLOATP BBASE J]
finally (\PUTBASEFLOATP CBASE K PRODUCT)))
CBASE])
(\MATMULT331
[LAMBDA (ABASE BBASE CBASE) (* jop: "29-Aug-86 12:20")
(* * "Multiply a 3 by 3 array by a 3 vector")
(for K from 0 to 4 by 2
do (bind (PRODUCT ← 0.0) declare (TYPE FLOATP PRODUCT) for I from (ITIMES K 3) by 2
as J from 0 to 4 by 2 do [SETQ PRODUCT (FPLUS PRODUCT (FTIMES (\GETBASEFLOATP ABASE I)
(\GETBASEFLOATP BBASE J]
finally (\PUTBASEFLOATP CBASE K PRODUCT)))
CBASE])
(\MATMULT333
[LAMBDA (ABASE BBASE CBASE) (* jop: "29-Aug-86 12:31")
(* * "Multiply two 3 by 3 arrays")
[bind (K ← 0) for ASTART from 0 to 12 by 6
do (for BSTART from 0 to 4 by 2
do (bind (PRODUCT ← 0.0) declare (TYPE FLOATP PRODUCT) for I from ASTART
to (IPLUS ASTART 4) by 2 as J from BSTART by 6
do [SETQ PRODUCT (FPLUS PRODUCT (FTIMES (\GETBASEFLOATP ABASE I)
(\GETBASEFLOATP BBASE J]
finally (\PUTBASEFLOATP CBASE K PRODUCT))
(SETQ K (IPLUS K 2]
CBASE])
(\MATMULT441
[LAMBDA (ABASE BBASE CBASE) (* jop: "29-Aug-86 12:20")
(* * "Multiply a 4 by 4 array by a 4 vector")
(for K from 0 to 6 by 2
do (bind (PRODUCT ← 0.0) declare (TYPE FLOATP PRODUCT) for I from (ITIMES K 4) by 2
as J from 0 to 6 by 2 do [SETQ PRODUCT (FPLUS PRODUCT (FTIMES (\GETBASEFLOATP ABASE I)
(\GETBASEFLOATP BBASE J]
finally (\PUTBASEFLOATP CBASE K PRODUCT)))
CBASE])
(\MATMULT444
[LAMBDA (ABASE BBASE CBASE) (* jop: "29-Aug-86 12:31")
(* * "Multiply two 4 by 4 arrays")
[bind (K ← 0) for ASTART from 0 to 24 by 8
do (for BSTART from 0 to 6 by 2
do (bind (PRODUCT ← 0.0) declare (TYPE FLOATP PRODUCT) for I from ASTART
to (IPLUS ASTART 6) by 2 as J from BSTART by 8
do [SETQ PRODUCT (FPLUS PRODUCT (FTIMES (\GETBASEFLOATP ABASE I)
(\GETBASEFLOATP BBASE J]
finally (\PUTBASEFLOATP CBASE K PRODUCT))
(SETQ K (IPLUS K 2]
CBASE])
)
(* unboxed arg handling)
(DECLARE: DONTCOPY
(* FOLLOWING DEFINITIONS EXPORTED)
(DECLARE: EVAL@COMPILE
[PUTPROPS
\CALLER.ARGS MACRO
(X (LET ((ARGS (CAR X))
(FORMS (CDR X)))
(BQUOTE (PROGN (\SLOWRETURN)
(LET ((AL (\MYALINK))
NEXT ., (for VAR in ARGS collect (if (LISTP VAR)
then
(LIST (CAR VAR)
0)
else VAR)))
[DECLARE ., (for VAR in ARGS when (LISTP VAR)
collect
(BQUOTE (TYPE %, (SELECTQ (CADR VAR)
((FLOATING FLOATP)
(CADR VAR))
(HELP))
%,
(CAR VAR]
(SETQ NEXT (fetch (FX NEXTBLOCK)
of AL))
.,
[for X in (REVERSE ARGS)
collect
(LET [(FORMS (BQUOTE (\.GETBASE32 \STACKSPACE (SETQ NEXT
(IDIFFERENCE
NEXT
WORDSPERCELL]
(if (LISTP X)
then
(BQUOTE (SETQ %, (CAR X)
(\FLOATBOX %, FORMS)))
else
(BQUOTE (SETQ %, X %, FORMS]
(\MAKEFREEBLOCK NEXT (TIMES %, (LENGTH ARGS)
WORDSPERCELL))
(replace (FX NEXTBLOCK)
of AL with NEXT)
(PROGN ., FORMS]
)
(* END EXPORTED DEFINITIONS)
)
(DEFINEQ
(\FLOATINGSCALE
(LAMBDA (X EXP10 BOX) (* JonL "17-May-84 19:04")
(* Scales X by indicated power of 10, i.e. does
X*10↑EXP10, reusing BOX if supplied.)
(COND
((NOT (FLOATP X))
(SETQ X (FLOAT X))
(OR (FLOATP BOX)
(SETQ BOX X))))
(COND
((EQ 0 EXP10)
(COND
((AND (FLOATP BOX)
(NEQ BOX X))
(replace HIWORD of BOX with (fetch HIWORD of X))
(replace LOWORD of BOX with (fetch LOWORD of X))
BOX)
(T X)))
((IGREATERP EXP10 0)
(COND
((IGREATERP EXP10 29)
(\FLOATINGSCALE (SETQ BOX (\BOXFTIMES2 X (\POWER.OF.TEN 29)
BOX))
(IDIFFERENCE EXP10 29)
BOX))
(T (\BOXFTIMES2 X (\POWER.OF.TEN EXP10)
BOX))))
((ILESSP EXP10 -29)
(\FLOATINGSCALE (SETQ BOX (\BOXFQUOTIENT X (\POWER.OF.TEN 29)
BOX))
(IPLUS EXP10 29)
BOX))
(T (\BOXFQUOTIENT X (\POWER.OF.TEN (IMINUS EXP10))
BOX)))))
(\INIT.POWERS.OF.TEN
[LAMBDA NIL (* bvm: "29-JAN-81 16:33")
(* Initialize array \POWERS.OF.TEN to values 10↑-29 thru 10↑+29. I suppose I could have the array cover the entire
range of floats, but the range is asymmetric and the numbers start losing significance at the ends, so it's not
really worth it)
(SETQ \POWERS.OF.TEN (ARRAY 59 (QUOTE POINTER)))
(SETA \POWERS.OF.TEN 30 1.0)
(for I from 1 to 29 bind POWTEN←1.0 do (SETA \POWERS.OF.TEN (IPLUS I 30)
(SETQ POWTEN (FTIMES POWTEN 10.0)))
(SETA \POWERS.OF.TEN (IDIFFERENCE 30 I)
(FQUOTIENT 1.0 POWTEN)))
\POWERS.OF.TEN])
(\CONVERT.FLOATING.NUMBER
(LAMBDA (X STR STRPTR FORMAT) (* DECLARATIONS: (RECORD FLOATFMT
(WIDTH DECPART EXPPART PAD SIGDIGITS)))
(DECLARE (SPECVARS X CNT POWTEN STR)) (* JonL "17-May-84 19:04")
(PROG ((SIGN (fetch (FLOATP SIGNBIT) of X))
(CNT 0)
(MAXWIDTH (fetch (STRINGP LENGTH) of STR))
EXP POWTEN BOX GUESS DIF DOEXP SIGDIGITS DECDIGITS EXPDIGITS WIDTH)
(* MAXWIDTH is currently 38, I think)
(OR FORMAT (SETQ FORMAT \FLOATFORMAT))
(COND
((AND FORMAT (NEQ FORMAT T))
(SETQ FORMAT (CDR FORMAT))
(COND
((SETQ WIDTH (fetch WIDTH of FORMAT)) (* Currently we use WIDTH only to decide when we must go
to free format)
(COND
((IGREATERP WIDTH MAXWIDTH) (* Can't produce anything larger than our scratch
string)
(SETQ WIDTH MAXWIDTH)))))
(SETQ SIGDIGITS (fetch SIGDIGITS of FORMAT))
(COND
((SETQ EXPDIGITS (fetch EXPPART of FORMAT))
(COND
((EQ 0 EXPDIGITS)
(SETQ EXPDIGITS NIL))
((IGREATERP EXPDIGITS 6) (* avoid scratch string overflow problems)
(SETQ EXPDIGITS 6)))))
(SETQ DECDIGITS (fetch DECPART of FORMAT))
(COND
((AND DECDIGITS (IGREATERP DECDIGITS 20)) (* avoid scratch string overflow)
(SETQ DECDIGITS 20)))
(COND
((AND DECDIGITS EXPDIGITS (NOT SIGDIGITS)) (* Exponent form always has exactly one digit before
point, so...)
(SETQ SIGDIGITS (ADD1 DECDIGITS))))))
(COND
((FEQP X 0.0)
(SETQ X 0)
(SETQ EXP 1)
(SETQ DOEXP EXPDIGITS)
(GO DOIT)))
(COND
((EQ SIGN 1)
(SETQ BOX (SETQ X (FDIFFERENCE 0.0 X)))))
(COND
((NOT SIGDIGITS)
(SETQ SIGDIGITS 7))
((IGREATERP SIGDIGITS 9) (* 9 is the most we can possibly squeeze out of our
floating-point numbers, and that's even pushing it)
(SETQ SIGDIGITS 9))
((ILESSP SIGDIGITS 1)
(SETQ SIGDIGITS 1)))
(SETQ EXP SIGDIGITS) (* EXP is to be the number of digits to the left of the
decimal point in the original number)
(* Scale X to have desired number of digits to left of
decimal point)
(SETQ GUESS (IPLUS (IQUOTIENT (IPLUS (ITIMES (IDIFFERENCE (fetch EXPONENT of X)
\EXPONENT.BIAS)
3)
105)
10)
-9))
(* Guess at the exponent of 10, which is about .3 times the exponent of 2 (LOG10 of 2 = .3010293). The contorted
expression is because IQUOTIENT truncates toward zero, and I want floor. GUESS is guess at the next larger
exponent of 10, hence 1 greater than our estimate of the exponent. GUESS is usually right on, but in case it
isn't, we test below)
(COND
((IGREATERP GUESS 25)
(SETQ X (SETQ BOX (\BOXFQUOTIENT X (\POWER.OF.TEN 25)
BOX)))
(SETQ EXP (IPLUS EXP 25))
(SETQ GUESS (IDIFFERENCE GUESS 25)))
((ILESSP GUESS -25)
(SETQ X (SETQ BOX (\BOXFTIMES2 X (\POWER.OF.TEN 25)
BOX)))
(SETQ EXP (IDIFFERENCE EXP 25))
(SETQ GUESS (IPLUS GUESS 25)))) (* Bring X safely into range of \POWERS.OF.TEN array)
(COND
((FGREATERP (\POWER.OF.TEN GUESS)
X)
(GO SMALLX)))
BIGX(add GUESS 1)
(COND
((NOT (FGREATERP (\POWER.OF.TEN GUESS)
X))
(GO BIGX)))
(GO GOTRANGE)
SMALLX (* exponent of X is less than guess)
(COND
((FGREATERP (\POWER.OF.TEN (SUB1 GUESS))
X)
(SETQ GUESS (SUB1 GUESS))
(GO SMALLX)))
GOTRANGE (* X is between powers GUESS-1 and GUESS)
(SETQ EXP (IPLUS EXP (SETQ DIF (IDIFFERENCE GUESS SIGDIGITS))))
(SETQ DOEXP (OR EXPDIGITS (COND
((NOT DECDIGITS)
(OR (IGREATERP EXP 8)
(ILESSP EXP -4)))
((IGREATERP EXP (IDIFFERENCE (OR WIDTH MAXWIDTH)
(IPLUS DECDIGITS (COND
((EQ SIGN 1)
2)
(T 1)))))
(* Respect DECDIGITS by not switching to exponent format
unless forced by WIDTH)
(SETQ DECDIGITS NIL) (* switch to free format)
T)))) (* Use exponent notation if requested, or if number
would have too many extra digits
(zeros) in normal notation.)
(COND
((NOT (EQ 0 DIF)) (* -DIF is the number of digits to right of decimal
point)
(COND
((AND DECDIGITS (NOT DOEXP)
(ILESSP DECDIGITS (IMINUS DIF))) (* Caller requested fewer decimal places than we would
by default provide, so round sooner)
(SETQ SIGDIGITS (IPLUS SIGDIGITS DIF DECDIGITS))
(SETQ DIF (IMINUS DECDIGITS))))
(COND
((NOT (EQ 0 DIF)) (* Check a second time because DECDIGITS could have been
zero in previous cond)
(SETQ X (SETQ BOX (\FLOATINGSCALE X (IMINUS DIF)
BOX)))))))
(SETQ X (\BOXFPLUSDIF X .5 NIL BOX)) (* round to an integer)
(COND
((NOT (FGREATERP (\POWER.OF.TEN SIGDIGITS)
X)) (* in rounding, we overflowed into next power of 10;
be careful not to round twice now...)
(SETQ X (\BOXFQUOTIENT (\BOXFPLUSDIF X 4.5 NIL X)
10.0 X))
(add EXP 1)))
(SETQ X (FIX X)) (* X now is SIGDIGITS long)
(SETQ POWTEN (NTH (QUOTE (100000000 10000000 1000000 100000 10000 1000 100 10 1))
(IDIFFERENCE 10 SIGDIGITS)))
(* CAR of POWTEN will be integer power of 10 with the 1 in the same place as the high digit of X.
Divide POWTEN by 10 by CDRing)
(COND
((EQ SIGN 1)
(\RPLCHARCODE STR (SETQ CNT 1)
(CHARCODE -))))
DOIT(COND
(DOEXP (* Have to go to exponent notation: n.nnnnEmm)
(\PRODUCE.FDIGIT)
(SETQ EXP (SUB1 EXP)) (* Account for the digit to left of dec pt)
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE %.))
(COND
(DECDIGITS (while (IGEQ (SETQ DECDIGITS (SUB1 DECDIGITS))
0)
do (\PRODUCE.FDIGIT)))
(T (do (\PRODUCE.FDIGIT) repeatuntil (EQ 0 X))))
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE E))
(COND
((ILESSP EXP 0)
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE -))
(SETQ EXP (IMINUS EXP)) (* Count the sign against the exponent count)
(AND EXPDIGITS (SETQ EXPDIGITS (SUB1 EXPDIGITS))))
((AND EXPDIGITS (IGREATERP EXPDIGITS 3))
(* We have room for + in any case)
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE +))
(SETQ EXPDIGITS (SUB1 EXPDIGITS))))
(COND
(EXPDIGITS (FRPTQ (IDIFFERENCE EXPDIGITS (COND
((ILESSP EXP 10)
(* count the chars we know we're printing: E and low
digit)
2)
(T 3)))
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE 0)))))
(COND
((IGEQ EXP 10)
(\RPLCHARCODE STR (add CNT 1)
(IPLUS (CHARCODE 0)
(IQUOTIENT EXP 10)))
(SETQ EXP (IREMAINDER EXP 10))))
(\RPLCHARCODE STR (add CNT 1)
(IPLUS EXP (CHARCODE 0))))
(T (FRPTQ EXP (\PRODUCE.FDIGIT))
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE %.))
(COND
((ILESSP EXP 0)
(COND
(DECDIGITS (COND
((ILESSP EXP (IMINUS DECDIGITS))
(* First sig digit is beyond the allowed dec places)
(SETQ EXP (SETQ X 0)))
(T (SETQ DECDIGITS (IPLUS DECDIGITS EXP))))))
(FRPTQ (IMINUS EXP)
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE 0)))))
(COND
(DECDIGITS (while (IGEQ (SETQ DECDIGITS (SUB1 DECDIGITS))
0)
do (\PRODUCE.FDIGIT)))
(T (do (\PRODUCE.FDIGIT) repeatuntil (EQ 0 X))))))
(RETURN (SUBSTRING STR 1 CNT STRPTR)))))
(\PRODUCE.FDIGIT
(LAMBDA NIL (* JonL "27-Sep-84 22:27")
(DECLARE (USEDFREE CNT STR X POWTEN))
(COND
((EQ 0 X)
(\RPLCHARCODE STR (add CNT 1)
(CHARCODE 0)))
(T (\RPLCHARCODE STR (add CNT 1)
(IPLUS (CHARCODE 0)
(IQUOTIENT X (CAR POWTEN))))
(SETQ X (IREMAINDER X (CAR POWTEN)))
(SETQ POWTEN (CDR POWTEN))))))
)
(DECLARE: DONTCOPY
(DECLARE: DOEVAL@COMPILE DONTCOPY
(GLOBALVARS \POWERS.OF.TEN)
)
(DECLARE: EVAL@COMPILE
[PUTPROPS \POWER.OF.TEN MACRO ((N)
(ELT \POWERS.OF.TEN (IPLUS N 30]
)
)
(DECLARE: DONTEVAL@LOAD DOCOPY
(\INIT.POWERS.OF.TEN)
)
(DEFINEQ
(FLOATP.TO.BCPL
[LAMBDA (FLONUM) (* bvm: "22-OCT-81 22:31")
(* * Converts a floating point number in IEEE format to an integer in BCPL floating-point format)
(OR (FLOATP FLONUM)
(SETQ FLONUM (FLOAT FLONUM)))
(PROG (RESULT FRAC (EXP (IPLUS (fetch EXPONENT of FLONUM)
2)))
(COND
((FEQP FLONUM 0.0)
(RETURN 0)))
[COND
((IGREATERP EXP 255) (* Overflow, so just return BCPL infinity)
(SETQ EXP 255)
(SETQ FRAC 4194303))
(T (SETQ FRAC (LRSH (fetch LONGFRACTION of FLONUM)
1]
(SETQ RESULT (create BCPLNUM
BCPLEXPONENT ← EXP
SIGNIFICANTBIT ← 1
BCPLHIFRACTION ←(LRSH FRAC 16)
BCPLLOFRACTION ←(LOGAND FRAC MAX.SMALL.INTEGER)))
(RETURN (COND
((EQ (fetch SIGNBIT of FLONUM)
1)
(IMINUS RESULT))
(T RESULT])
(BCPL.TO.FLOATP
[LAMBDA (BCPLNUM) (* bvm: "22-OCT-81 22:34")
(* Converts BCPLNUM, an integer in BCPL floating-point
format, to a FLOATP, which is IEEE standard)
(PROG (SIGN EXP FRAC)
(COND
((ILESSP BCPLNUM 0)
(SETQ BCPLNUM (IMINUS BCPLNUM)) (* In a negative BCPL format, whole number is
complemented)
(SETQ SIGN 1))
((IEQP BCPLNUM 0) (* Canonical form for 0.0)
(RETURN (FPLUS 0.0)))
(T (SETQ SIGN 0)))
(COND
((OR (SMALLP BCPLNUM)
(NEQ (fetch SIGNIFICANTBIT of BCPLNUM)
1))
(ERROR "Not a valid BCPL flonum" BCPLNUM)))
[COND
((ILESSP (SETQ EXP (IDIFFERENCE (fetch BCPLEXPONENT of BCPLNUM)
2))
0) (* Underflow. IEEE exponent is off by 2 because the bias
is one smaller in IEEE format and we shift the mantissa
left one)
(RETURN (FPLUS 0.0]
(SETQ FRAC (LLSH (fetch RESTOFFRACTION of BCPLNUM)
1))
(RETURN (create FLOATP
SIGNBIT ← SIGN
EXPONENT ← EXP
HIFRACTION ←(LRSH FRAC 16)
LOFRACTION ←(LOGAND FRAC MAX.SMALL.INTEGER])
)
(DECLARE: EVAL@COMPILE DONTCOPY
[DECLARE: EVAL@COMPILE
(BLOCKRECORD BCPLNUM ((BCPLSIGNBIT BITS 1)
(BCPLEXPONENT BITS 8) (* exponent, biased by 128)
(SIGNIFICANTBIT BITS 1) (* Always 1 in a bcpl num;
binary point is to left)
(RESTOFFRACTION BITS 22))
(BLOCKRECORD BCPLNUM ((NIL BITS 10)
(BCPLHIFRACTION BITS 6)
(BCPLLOFRACTION BITS 16)))
(CREATE (CREATECELL \FIXP)))
]
)
(PUTPROPS \UNBOXFLOAT1 ARGNAMES (X OP))
(PUTPROPS \UNBOXFLOAT2 ARGNAMES (X Y OP))
(PUTPROPS \UNBOXFLOAT3 ARGNAMES (X Y Z OP))
(DECLARE: DONTEVAL@LOAD DOEVAL@COMPILE DONTCOPY COMPILERVARS
(ADDTOVAR NLAMA )
(ADDTOVAR NLAML )
(ADDTOVAR LAMA FPLUS FTIMES)
)
(PUTPROPS LLFLOAT COPYRIGHT ("Xerox Corporation" 1982 1984 1985 1986))
(DECLARE: DONTCOPY
(FILEMAP (NIL (4468 4820 (\PUTBASEFLOATP 4478 . 4649) (\GETBASEFLOATP 4651 . 4818)) (5443 7074 (FTIMES
5453 . 5974) (FPLUS 5976 . 6372) (FQUOTIENT 6374 . 6515) (FDIFFERENCE 6517 . 6636) (FGREATERP 6638 .
7010) (FABS 7012 . 7072)) (7092 10021 (\SLOWFDIFFERENCE 7102 . 7282) (\SLOWFPLUS2 7284 . 7534) (
\SLOWFTIMES2 7536 . 7704) (\SLOWFQUOTIENT 7706 . 7940) (\SLOWFGREATERP 7942 . 10019)) (10022 32559 (
\FZEROP 10032 . 10374) (FEQP 10376 . 11438) (FLOAT 11440 . 12054) (\FLOAT 12056 . 14100) (
\FIXP.FROM.FLOATP 14102 . 15368) (FIXR 15370 . 17912) (\BOXFPLUSDIF 17914 . 20425) (\BOXFQUOTIENT
20427 . 22614) (\BOXFTIMES2 22616 . 25649) (\INFINITY 25651 . 26378) (\MAKEFLOAT 26380 . 30119) (
MAKEFLOATNUMBER 30121 . 30509) (PutFloat 30511 . 30746) (\FZEROP 30748 . 31090) (MAKERATIONAL 31092 .
32557)) (32749 34015 (SQRT 32759 . 34013)) (45964 50702 (\UNBOXFLOAT1 45974 . 47560) (\UNBOXFLOAT2
47562 . 48799) (\UNBOXFLOAT3 48801 . 50700)) (50703 54895 (\MATMULT133 50713 . 51349) (\MATMULT144
51351 . 51987) (\MATMULT331 51989 . 52643) (\MATMULT333 52645 . 53440) (\MATMULT441 53442 . 54096) (
\MATMULT444 54098 . 54893)) (57493 68738 (\FLOATINGSCALE 57503 . 58584) (\INIT.POWERS.OF.TEN 58586 .
59307) (\CONVERT.FLOATING.NUMBER 59309 . 68250) (\PRODUCE.FDIGIT 68252 . 68736)) (69009 71402 (
FLOATP.TO.BCPL 69019 . 70009) (BCPL.TO.FLOATP 70011 . 71400)))))
STOP