1 ! Copyright (C) 2018 Alexander Ilin.
2 ! See http://factorcode.org/license.txt for BSD license.
3 USING: combinators.smart kernel locals math math.bitwise
4 math.functions math.order math.parser ryu.data sequences
11 : mul-shift ( x mul shift -- y )
12 [ first2 rot [ * ] keep swapd * -64 shift + ] [ 64 - neg ] bi* shift ;
14 : mul-shift-all ( mmShift m mul shift -- vm vp vr )
16 [ [ 1 - swap - ] 2dip mul-shift ]
17 [ [ 2 + ] 2dip mul-shift ]
20 :: pow-5-factor ( x -- y )
21 x f 0 [ 2dup x > or ] [
22 [ 5 /mod ] 2dip rot zero? [ 1 + ] [ nip dupd ] if
25 : multiple-of-power-of-5 ( p value -- ? )
28 : double-pow-5-bits ( n -- m )
30 DOUBLE_LOG2_5_NUMERATOR * DOUBLE_LOG2_5_DENOMINATOR + 1 -
31 DOUBLE_LOG2_5_DENOMINATOR /i
34 : decimal-length ( m -- n )
54 } [ dupd >= ] find-last [ 2 + ] [ drop 1 ] if nip ; inline
56 CONSTANT: mantissaBits 52
57 CONSTANT: exponentBits 11
58 CONSTANT: offset 1023 ! (1 << (exponentBits - 1)) - 1
60 :: unpack-bits ( value -- e2 m2 acceptBounds ieeeExponent<=1? neg? string/f )
62 dup mantissaBits exponentBits + bit? :> sign
63 dup mantissaBits bits :> ieeeMantissa
64 mantissaBits neg shift exponentBits bits :> ieeeExponent
67 exponentBits on-bits ieeeExponent = [
68 ieeeMantissa zero? [ sign "-Inf" "Inf" ? ] [ "NaN" ] if
71 ieeeMantissa [ sign "-0e0" "0e0" ? ] [
73 -1 offset - mantissaBits - e2!
77 offset - mantissaBits - 2 - e2!
78 ieeeMantissa mantissaBits set-bit m2!
81 ] if [ e2 m2 dup even? ieeeExponent 1 <= sign ] dip ; inline
83 :: prepare-output ( vp! acceptBounds vmIsTrailingZeros! vrIsTrailingZeros! vr! vm! -- output )
84 ! vr is converted into the output
86 ! the if has this stack-effect: ( lastRemovedDigit -- lastRemovedDigit' output )
87 vmIsTrailingZeros vrIsTrailingZeros or [
89 [ vp 10 /i vm 10 /i 2dup > ] [
91 vmIsTrailingZeros [ vm 10 divisor? vmIsTrailingZeros! ] when
92 vrIsTrailingZeros [ dup zero? vrIsTrailingZeros! ] when
93 vr 10 /mod swap vr! nip ! lastRemovedDigit!
96 [ vm dup 10 /i dup 10 * swapd = ] [
98 vrIsTrailingZeros [ dup zero? vrIsTrailingZeros! ] when
99 vr 10 /mod swap vr! nip ! lastRemovedDigit!
101 ] while drop ! Drop (vm 10 /i) result from the while condition.
105 vr even? [ drop 4 ] when ! 4 lastRemovedDigit!
108 vr over 5 >= [ 1 + ] [
110 acceptBounds vmIsTrailingZeros and not [ 1 + ] when
115 [ vp 10 /i vm 10 /i 2dup > ] [
117 vr 10 /mod swap vr! nip ! lastRemovedDigit!
119 vr dup vm = [ 1 + ] [
120 over 5 >= [ 1 + ] when
124 :: produce-output ( exp sign output -- string )
127 output number>string 1 cut-slice dup empty? f "." ? swap
130 ] "" append-outputs-as ; inline
134 :: print-float ( value -- string )
135 value >float unpack-bits
136 :> ( e2 m2 acceptBounds ieeeExponent<=1 sign output )
139 mantissaBits 2^ m2 = not ieeeExponent<=1 or 1 0 ? :> mmShift
140 f f 0 0 0 :> ( vmIsTrailingZeros! vrIsTrailingZeros! e10! vr! vm! )
141 ! After the following loop vp is left on stack.
143 e2 DOUBLE_LOG10_2_NUMERATOR * DOUBLE_LOG10_2_DENOMINATOR /i 0 max :> q
145 q double-pow-5-bits DOUBLE_POW5_INV_BITCOUNT + 1 - :> k
147 mmShift m2 q DOUBLE_POW5_INV_SPLIT nth-unsafe i mul-shift-all vr! swap vm! ! vp on stack
150 q mv multiple-of-power-of-5 vrIsTrailingZeros!
153 q mv mmShift - 1 - multiple-of-power-of-5 vmIsTrailingZeros!
155 q mv 2 + multiple-of-power-of-5 1 0 ? - ! vp!
160 e2 neg DOUBLE_LOG10_5_NUMERATOR * DOUBLE_LOG10_5_DENOMINATOR /i 1 [-] :> q
163 i double-pow-5-bits DOUBLE_POW5_BITCOUNT - :> k
165 mmShift m2 i DOUBLE_POW5_SPLIT nth-unsafe j mul-shift-all vr! swap vm! ! vp on stack
167 mv 1 bitand bitnot q >= vrIsTrailingZeros!
169 mv 1 - mmShift - bitnot 1 bitand q >= vmIsTrailingZeros!
173 q 1 - on-bits mv bitand zero? vrIsTrailingZeros!
177 [ decimal-length e10 + 1 - sign ] keep ! exp sign vp
178 acceptBounds vmIsTrailingZeros vrIsTrailingZeros vr vm
179 prepare-output produce-output
182 ALIAS: d2s print-float