|
| stwotox.sa 3.1 12/10/90
|
| stwotox --- 2**X
| stwotoxd --- 2**X for denormalized X
| stentox --- 10**X
| stentoxd --- 10**X for denormalized X
|
| Input: Double-extended number X in location pointed to
| by address register a0.
|
| Output: The function values are returned in Fp0.
|
| Accuracy and Monotonicity: The returned result is within 2 ulps in
| 64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
| result is subsequently rounded to double precision. The
| result is provably monotonic in double precision.
|
| Speed: The program stwotox takes approximately 190 cycles and the
| program stentox takes approximately 200 cycles.
|
| Algorithm:
|
| twotox
| 1. If |X| > 16480, go to ExpBig.
|
| 2. If |X| < 2**(-70), go to ExpSm.
|
| 3. Decompose X as X = N/64 + r where |r| <= 1/128. Furthermore
| decompose N as
| N = 64(M + M') + j, j = 0,1,2,...,63.
|
| 4. Overwrite r := r * log2. Then
| 2**X = 2**(M') * 2**(M) * 2**(j/64) * exp(r).
| Go to expr to compute that expression.
|
| tentox
| 1. If |X| > 16480*log_10(2) (base 10 log of 2), go to ExpBig.
|
| 2. If |X| < 2**(-70), go to ExpSm.
|
| 3. Set y := X*log_2(10)*64 (base 2 log of 10). Set
| N := round-to-int(y). Decompose N as
| N = 64(M + M') + j, j = 0,1,2,...,63.
|
| 4. Define r as
| r := ((X - N*L1)-N*L2) * L10
| where L1, L2 are the leading and trailing parts of log_10(2)/64
| and L10 is the natural log of 10. Then
| 10**X = 2**(M') * 2**(M) * 2**(j/64) * exp(r).
| Go to expr to compute that expression.
|
| expr
| 1. Fetch 2**(j/64) from table as Fact1 and Fact2.
|
| 2. Overwrite Fact1 and Fact2 by
| Fact1 := 2**(M) * Fact1
| Fact2 := 2**(M) * Fact2
| Thus Fact1 + Fact2 = 2**(M) * 2**(j/64).
|
| 3. Calculate P where 1 + P approximates exp(r):
| P = r + r*r*(A1+r*(A2+...+r*A5)).
|
| 4. Let AdjFact := 2**(M'). Return
| AdjFact * ( Fact1 + ((Fact1*P) + Fact2) ).
| Exit.
|
| ExpBig
| 1. Generate overflow by Huge * Huge if X > 0; otherwise, generate
| underflow by Tiny * Tiny.
|
| ExpSm
| 1. Return 1 + X.
|
| Copyright (C) Motorola, Inc. 1990
| All Rights Reserved
|
| For details on the license for this file, please see the
| file, README, in this same directory.
|STWOTOX idnt 2,1 | Motorola 040 Floating Point Software Package
fmovex %fp0,%fp1
fmuls #0x42800000,%fp1 | ...64 * X
fmovel %fp1,N(%a6) | ...N = ROUND-TO-INT(64 X)
movel %d2,-(%sp)
lea EXPTBL,%a1 | ...LOAD ADDRESS OF TABLE OF 2^(J/64)
fmovel N(%a6),%fp1 | ...N --> FLOATING FMT
movel N(%a6),%d0
movel %d0,%d2
andil #0x3F,%d0 | ...D0 IS J
asll #4,%d0 | ...DISPLACEMENT FOR 2^(J/64)
addal %d0,%a1 | ...ADDRESS FOR 2^(J/64)
asrl #6,%d2 | ...d2 IS L, N = 64L + J
movel %d2,%d0
asrl #1,%d0 | ...D0 IS M
subl %d0,%d2 | ...d2 IS M', N = 64(M+M') + J
addil #0x3FFF,%d2
movew %d2,ADJFACT(%a6) | ...ADJFACT IS 2^(M')
movel (%sp)+,%d2
|--SUMMARY: a1 IS ADDRESS FOR THE LEADING PORTION OF 2^(J/64),
|--D0 IS M WHERE N = 64(M+M') + J. NOTE THAT |M| <= 16140 BY DESIGN.
|--ADJFACT = 2^(M').
|--REGISTERS SAVED SO FAR ARE (IN ORDER) FPCR, D0, FP1, a1, AND FP2.
EXPBIG:
|--|X| IS LARGE, GENERATE OVERFLOW IF X > 0; ELSE GENERATE UNDERFLOW
|--REGISTERS SAVE SO FAR ARE FPCR AND D0
movel X(%a6),%d0
cmpil #0,%d0
blts EXPNEG
bclrb #7,(%a0) |t_ovfl expects positive value
bra t_ovfl
EXPNEG:
bclrb #7,(%a0) |t_unfl expects positive value
bra t_unfl
.global stentoxd
stentoxd:
|--ENTRY POINT FOR 10**(X) FOR DENORMALIZED ARGUMENT
.global stentox
stentox:
|--ENTRY POINT FOR 10**(X), HERE X IS FINITE, NON-ZERO, AND NOT NAN'S
fmovemx (%a0),%fp0-%fp0 | ...LOAD INPUT, do not set cc's
fmovel %fp1,N(%a6) | ...N=INT(X*64*LOG10/LOG2)
movel %d2,-(%sp)
lea EXPTBL,%a1 | ...LOAD ADDRESS OF TABLE OF 2^(J/64)
fmovel N(%a6),%fp1 | ...N --> FLOATING FMT
movel N(%a6),%d0
movel %d0,%d2
andil #0x3F,%d0 | ...D0 IS J
asll #4,%d0 | ...DISPLACEMENT FOR 2^(J/64)
addal %d0,%a1 | ...ADDRESS FOR 2^(J/64)
asrl #6,%d2 | ...d2 IS L, N = 64L + J
movel %d2,%d0
asrl #1,%d0 | ...D0 IS M
subl %d0,%d2 | ...d2 IS M', N = 64(M+M') + J
addil #0x3FFF,%d2
movew %d2,ADJFACT(%a6) | ...ADJFACT IS 2^(M')
movel (%sp)+,%d2
|--SUMMARY: a1 IS ADDRESS FOR THE LEADING PORTION OF 2^(J/64),
|--D0 IS M WHERE N = 64(M+M') + J. NOTE THAT |M| <= 16140 BY DESIGN.
|--ADJFACT = 2^(M').
|--REGISTERS SAVED SO FAR ARE (IN ORDER) FPCR, D0, FP1, a1, AND FP2.
expr:
|--FPCR, FP2, FP3 ARE SAVED IN ORDER AS SHOWN.
|--ADJFACT CONTAINS 2**(M'), FACT1 + FACT2 = 2**(M) * 2**(J/64).
|--FP0 IS R. THE FOLLOWING CODE COMPUTES
|-- 2**(M'+M) * 2**(J/64) * EXP(R)
fmovex %fp0,%fp1
fmulx %fp1,%fp1 | ...FP1 IS S = R*R
fmoved EXPA5,%fp2 | ...FP2 IS A5
fmoved EXPA4,%fp3 | ...FP3 IS A4
fmulx %fp1,%fp2 | ...FP2 IS S*A5
fmulx %fp1,%fp3 | ...FP3 IS S*A4
faddd EXPA3,%fp2 | ...FP2 IS A3+S*A5
faddd EXPA2,%fp3 | ...FP3 IS A2+S*A4
fmulx %fp1,%fp2 | ...FP2 IS S*(A3+S*A5)
fmulx %fp1,%fp3 | ...FP3 IS S*(A2+S*A4)
faddd EXPA1,%fp2 | ...FP2 IS A1+S*(A3+S*A5)
fmulx %fp0,%fp3 | ...FP3 IS R*S*(A2+S*A4)
fmulx %fp1,%fp2 | ...FP2 IS S*(A1+S*(A3+S*A5))
faddx %fp3,%fp0 | ...FP0 IS R+R*S*(A2+S*A4)
faddx %fp2,%fp0 | ...FP0 IS EXP(R) - 1
|--FINAL RECONSTRUCTION PROCESS
|--EXP(X) = 2^M*2^(J/64) + 2^M*2^(J/64)*(EXP(R)-1) - (1 OR 0)
Die Informationen auf dieser Webseite wurden
nach bestem Wissen sorgfältig zusammengestellt. Es wird jedoch weder Vollständigkeit, noch Richtigkeit,
noch Qualität der bereit gestellten Informationen zugesichert.
Bemerkung:
Die farbliche Syntaxdarstellung und die Messung sind noch experimentell.
