From Newsgroup: comp.lang.forth

Implementers of floating point output have likely encountered the exponent function which returns the exponent of a given fp number. It's typically defined:

: exp# ( f -- exp )
fdup f0= if fdrop 0 exit then
fabs flog
floor f>s ;

It works well enough until applied to the x87 FPU where the FLOG FLOOR calculation is done at 80-bit precision:

1e-1 exp# . -1 ok
1e-2 exp# . -2 ok
1e-3 exp# . -4 ok
1e-4 exp# . -4 ok
1e-5 exp# . -5 ok

In the third case -3 was expected but -4 was produced. Needless to say this can/has resulted in problems down the track. Testing over the full range of exponents (0..4932) 235 such duplicated exponents were found. Below is the test program used, along with a potential fix.

Enjoy

-----------------

create aa 8 allot

: exp# ( f -- exp )
fdup f0= if fdrop 0 exit then
fabs flog
\ aa df! aa df@ \ fix
floor f>s ;

4932 value size \ assume 80 bit extended precision

variable last
variable start
variable cnt

: .duplicates ( -- )
start @ size 0 do dup @ dup last @ = if dup . 1 cnt +! then last !
1 cells + loop drop ;

: test+ ( -- )
-10000 last ! here start !
1e size 1+ 0 do fdup exp# , 10e f* loop fdrop .duplicates ;

: test- ( -- )
-10000 last ! here start !
0.1e size 1+ 0 do fdup exp# , 10e f/ loop fdrop .duplicates ;

: test ( -- )
cr cr ." Testing for duplicated exponents ... "
0 cnt ! test+ test- cr cr cnt @ . ." duplicates " ;

test

--------------------
--- Synchronet 3.22a-Linux NewsLink 1.2