IeeeFloatA.mesa - Mesa implementation of floating point ops +, -, *, /, Rem
Copyright © 1985 by Xerox Corporation. All rights reserved.
Stewart January 14, 1984 4:49 pm
Rovner May 4, 1983 9:53 am
Russ Atkinson (RRA) May 28, 1985 6:44:40 pm PDT
DIRECTORY
Ieee USING [ADC2, ADC3, BitOn, DeNormalize, Ext, funny, HiBit, HiddenBit, LongNumber, LowSignificandMask, Mul32, NormalType, Pack, PostNormalize, RShift1in1, SetDivisionByZero, SetInvalidOperation, SignBit, StepTwo, TNnn, TNnz, TNzn, TNzz, Unpack],
PrincOpsUtils USING [BITAND, BITOR, BITXOR],
Real USING [AddInfinityNaN, DivideInfinityNaN, Fix, MinusZero, MultiplyInfinityNaN, NumberType, PlusZero],
RealExceptions USING [ClearExceptions],
RealOps USING [GetMode, Mode, SetMode];
IeeeFloatA:
PROGRAM
IMPORTS Ieee, PrincOpsUtils, Real, RealExceptions, RealOps
EXPORTS RealOps
= BEGIN OPEN Ieee, PrincOpsUtils, Real, RealOps;
NormalizedCardinal:
PROC [g:
CARDINAL]
RETURNS [
BOOL] =
INLINE {
RETURN [PrincOpsUtils.BITAND[g, Ieee.HiddenBit] # 0];
};
AddNormal:
PROC [x, y: Ext]
RETURNS [Ext, Ext] = {
ediff: INTEGER ← x.exp - y.exp;
normalized: BOOL;
IF ediff > 0
THEN DeNormalize[@y, -ediff]
ELSE {DeNormalize[@x, ediff]; x.exp ← y.exp};
normalized ← NormalizedCardinal[BITOR[x.frac.highbits, y.frac.highbits]];
IF x.det.sign = y.det.sign
THEN {
cy: CARDINAL;
[cy, x.frac.lowbits] ← ADC2[x.frac.lowbits, y.frac.lowbits];
[cy, x.frac.highbits] ← ADC3[x.frac.highbits, y.frac.highbits, cy];
IF cy # 0
THEN {
IF BitOn[x.frac.lowbits, 1] THEN x.frac.lowbits ← BITOR[x.frac.lowbits, 2];
x.frac ← RShift1in1[x.frac];
x.exp ← x.exp + 1;
};
}
ELSE {
IF x.frac.lc >= y.frac.lc
THEN x.frac.lc ← x.frac.lc - y.frac.lc
ELSE {x.frac.lc ← y.frac.lc - x.frac.lc; x.det.sign ← y.det.sign; };
IF x.frac.highbits = 0
AND
NOT BitOn[x.frac.lowbits, LowSignificandMask]
THEN {
x.det.sign ← RealOps.GetMode[].round = rm;
x.det.type ← zero;
}
ELSE {IF normalized THEN PostNormalize[@x]; };
};
IF y.det.sticky THEN x.det.sticky ← TRUE;
RETURN [x, y];
};
FAdd:
PUBLIC
SAFE
PROC [a, b:
REAL, m: Mode]
RETURNS [r:
REAL] =
TRUSTED {
inner:
PROC = {
x: Ext ← Unpack[a];
y: Ext ← Unpack[b];
{
SELECT
TRUE
FROM
Basics.
BITAND[funny,
BITOR[
LOOPHOLE[x.det],
LOOPHOLE[y.det]]] = 0 => {
Normal case: Neither operand is infinity or NaN
SELECT NormalType[x.det, y.det]
FROM
TNnn => NULL;
TNnz => GOTO RetA;
TNzn => GOTO RetB;
TNzz => {
a ← PlusZero;
IF m.round = rm AND x.det.sign AND y.det.sign THEN a ← MinusZero;
GO TO RetA;
};
ENDCASE => ERROR;
};
x.det.type = nan => GOTO RetX;
y.det.type = nan => GOTO RetY;
x.det.type = infinity
AND y.det.type = infinity =>
IF m.im = affine
AND x.det.sign = y.det.sign
THEN GOTO RetA
ELSE {
SetInvalidOperation[];
x.det.type ← nan;
x.frac.lc ← AddInfinityNaN;
GOTO RetX;
};
x.det.type = infinity => GOTO RetA;
y.det.type = infinity => GOTO RetB;
ENDCASE => ERROR;
[x, y] ← AddNormal[x, y];
StepTwo[@x];
GOTO RetX;
EXITS
RetA => r ← a;
RetB => r ← b;
RetX => r ← Pack[@x];
RetY => r ← Pack[@y];
};
};
DoWithMode[inner, m];
};
FSub:
PUBLIC
SAFE
PROC [a, b:
REAL, m: Mode]
RETURNS [
REAL] =
TRUSTED {
Negate b
LOOPHOLE[b, LongNumber].highbits ← BITXOR[LOOPHOLE[b, LongNumber].highbits, SignBit];
RETURN[FAdd[a, b, m]];
};
FMul:
PUBLIC
SAFE
PROC [a, b:
REAL, m: Mode]
RETURNS [r:
REAL] =
TRUSTED {
inner:
PROC = {
x: Ext ← Ieee.Unpack[a];
y: Ext ← Ieee.Unpack[b];
{
lo: LongNumber;
x.det.sign ← x.det.sign # y.det.sign; -- due to different bias
x.exp ← x.exp + y.exp + 1;
SELECT
TRUE
FROM
Basics.
BITAND[funny,
BITOR[
LOOPHOLE[x.det],
LOOPHOLE[y.det]]] = 0 => {
Normal case: Neither operand is infinity or NaN
SELECT NormalType[x.det, y.det]
FROM
TNnn => NULL;
TNnz => GOTO RetZero;
TNzn => GOTO RetZero;
TNzz => GOTO RetZero;
ENDCASE => ERROR;
};
x.det.type = nan => GOTO RetX;
y.det.type = nan => GOTO RetY;
x.det.type = infinity
OR y.det.type = infinity => {
IF x.det.type = zero
OR y.det.type = zero
THEN {
SetInvalidOperation[];
x.det.type ← nan;
x.frac.lc ← MultiplyInfinityNaN;
}; -- Possible test for unnormal zero later.
x.det.type ← infinity;
GOTO RetX;
};
ENDCASE => NULL;
[x.frac, lo] ← Mul32[x.frac, y.frac]; -- normalize (double normalize??)
IF x.frac.lc=0
THEN {
x.frac.lc ← lo.lc;
x.exp ← x.exp-32;
}
ELSE {
IF lo.lc # 0 THEN x.det.sticky ← TRUE;
};
PostNormalize[@x];
StepTwo[@x];
GOTO RetX;
EXITS
RetX => r ← Pack[@x];
RetY => r ← Pack[@y];
RetZero => r ← IF x.det.sign THEN MinusZero ELSE PlusZero;
};
};
DoWithMode[inner, m];
};
FDiv:
PUBLIC
SAFE
PROC [a, b:
REAL, m: Mode]
RETURNS [c:
REAL] =
TRUSTED {
inner:
PROC = {
x: Ext ← Unpack[a];
y: Ext ← Unpack[b];
{
x.det.sign ← x.det.sign # y.det.sign;
{
SELECT
TRUE
FROM
Basics.
BITAND[funny,
BITOR[
LOOPHOLE[x.det],
LOOPHOLE[y.det]]] = 0 => {
Neither of the number types is infinity or NaN
SELECT NormalType[x.det, y.det]
FROM
TNnn => {};
this is usual flow
RRA: this used to be the following bogus code:
IF NormalizedCardinal[y.frac.highbits] THEN NULL ELSE GOTO ZeroDivide
TNnz => GOTO ZeroDivide;
TNzn => GOTO RetZero;
TNzz => GOTO Undef;
ENDCASE => ERROR;
};
x.det.type = nan => GOTO RetX;
y.det.type = nan => GOTO RetY;
x.det.type = infinity
AND y.det.type = infinity =>
GO
TO Undef;
both are some kind of infinity, so results are undefined
y.det.type = infinity => GOTO RetZero;
x.det.type = infinity => GOTO RetX;
ENDCASE =>
ERROR;
It should not be possible to get here!
EXITS
ZeroDivide => {
Divide by 0, so result is infnity
SetDivisionByZero[];
x.det.type ← infinity;
GOTO RetX;
};
Undef => {
The results are undefined (inf/inf or 0/0)
SetInvalidOperation[];
x.det.type ← nan;
x.frac.lc ← DivideInfinityNaN;
GOTO RetX;
};
};
Both numbers need to be normalized for the Divide routine to work.
IF NOT NormalizedCardinal[x.frac.highbits] THEN PostNormalize[@x];
IF NOT NormalizedCardinal[y.frac.highbits] THEN PostNormalize[@y];
x ← Divide[x, y];
StepTwo[@x];
GOTO RetX;
EXITS
RetZero => c ← IF x.det.sign THEN MinusZero ELSE PlusZero;
RetX => c ← Pack[@x];
RetY => c ← Pack[@y];
};
};
DoWithMode[inner, m];
};
Divide:
PROC [a, b: Ext]
RETURNS [y: Ext] = {
cy: BOOL;
Step:
PROC =
INLINE {
y.frac.lc ← y.frac.lc + y.frac.lc;
IF cy
OR a.frac.lc >= b.frac.lc
THEN {
a.frac.lc ← a.frac.lc - b.frac.lc;
y.frac.lowbits ← y.frac.lowbits + 1; -- can't carry!
};
cy ← BitOn[a.frac.highbits, HiBit];
a.frac.lc ← a.frac.lc + a.frac.lc;
};
IF b.frac.lc = 0 THEN ERROR;
y ← a;
y.det ← [sign: a.det.sign, sticky: FALSE, blank: 0, type: Real.NumberType[normal]];
y.exp ← a.exp - b.exp;
y.frac.lc ← 0;
cy ← FALSE;
THROUGH [0..32) DO Step[]; ENDLOOP;
WHILE PrincOpsUtils.
BITAND[y.frac.highbits, HiBit] = 0
DO
normalize
Step[];
y.exp ← y.exp - 1;
ENDLOOP;
IF a.frac.lc # 0 THEN y.frac.lowbits ← BITOR[y.frac.lowbits, 1];
};
FRem:
PUBLIC
SAFE
PROC [a, b:
REAL, m: Mode]
RETURNS [c:
REAL] =
TRUSTED {
n: INT ← Real.Fix[a / b];
c ← a - n * b;
};
DoWithMode:
PROC [inner:
PROC, mode: Mode] = {
old: Mode ← RealOps.SetMode[mode];
RealExceptions.ClearExceptions[];
inner[ ! UNWIND => [] ← RealOps.SetMode[old]];
[] ← RealOps.SetMode[old];
};
END.
July 8, 1980 6:36 PM, L. Stewart; FRem dummy added
August 25, 1980 10:07 AM, L. Stewart; new Divide
January 28, 1981 11:50 AM, L. Stewart; Fix to Multiply (denormalized numbers)
August 27, 1982 11:24 am, L. Stewart, added TRUSTED everywhere
January 14, 1984 4:49 pm, Stewart, change to Ieee