src/float.4th

   1 \ Floating point arithmetic
   2
   3 ( Cheating for now by using forth.jl CODE/END-CODE to
   4   access Julia's floating point support.  This isn't
   5   at all portable.  That said, the year is 2016 and most
   6   CPUs implement these operations - even the trig functions,
   7   so I don't feel too bad! )
   8
   9 CODE f+
  10     b = reinterpret(Float64, popPS())
  11     a = reinterpret(Float64, popPS())
  12     pushPS(reinterpret(Int64, a+b))
  13 END-CODE
  14
  15 CODE f-
  16     b = reinterpret(Float64, popPS())
  17     a = reinterpret(Float64, popPS())
  18     pushPS(reinterpret(Int64, a-b))
  19 END-CODE
  20
  21 CODE f*
  22     b = reinterpret(Float64, popPS())
  23     a = reinterpret(Float64, popPS())
  24     pushPS(reinterpret(Int64, a*b))
  25 END-CODE
  26
  27 CODE f/
  28     b = reinterpret(Float64, popPS())
  29     a = reinterpret(Float64, popPS())
  30     pushPS(reinterpret(Int64, a/b))
  31 END-CODE
  32
  33 CODE f^
  34     b = reinterpret(Float64, popPS())
  35     a = reinterpret(Float64, popPS())
  36     pushPS(reinterpret(Int64, a^b))
  37 END-CODE
  38
  39 CODE f>
  40     b = reinterpret(Float64, popPS())
  41     a = reinterpret(Float64, popPS())
  42     if a > b
  43         pushPS(-1)
  44     else
  45         pushPS(0)
  46     end
  47 END-CODE
  48
  49 CODE f<
  50     b = reinterpret(Float64, popPS())
  51     a = reinterpret(Float64, popPS())
  52     if a < b
  53         pushPS(-1)
  54     else
  55         pushPS(0)
  56     end
  57 END-CODE
  58
  59 CODE f=
  60     b = reinterpret(Float64, popPS())
  61     a = reinterpret(Float64, popPS())
  62     if a == b
  63         pushPS(-1)
  64     else
  65         pushPS(0)
  66     end
  67 END-CODE
  68
  69 : f<=
  70     f> invert ;
  71
  72 : f>=
  73     f< invert ;
  74
  75 CODE fmod
  76     b = reinterpret(Float64, popPS())
  77     a = reinterpret(Float64, popPS())
  78     pushPS(reinterpret(Int64, a%b))
  79 END-CODE
  80
  81 CODE flog
  82     a = reinterpret(Float64, popPS())
  83     pushPS(reinterpret(Int64, log(a)))
  84 END-CODE
  85
  86 CODE fexp
  87     a = reinterpret(Float64, popPS())
  88     pushPS(reinterpret(Int64, exp(a)))
  89 END-CODE
  90
  91 CODE fnan?
  92     a = reinterpret(Float64, popPS())
  93     if isnan(a)
  94         pushPS(-1)
  95     else
  96         pushPS(0)
  97     end
  98 END-CODE
  99
 100 CODE finf?
 101     a = reinterpret(Float64, popPS())
 102     if isinf(a)
 103         pushPS(-1)
 104     else
 105         pushPS(0)
 106     end
 107 END-CODE
 108
 109 CODE i->f
 110     pushPS(reinterpret(Int64, Float64(popPS())))
 111 END-CODE
 112
 113 CODE f->i
 114     a = reinterpret(Float64, popPS())
 115     pushPS(Int64(round(a)))
 116 END-CODE
 117
 118 : f/mod
 119     2dup fmod -rot f/ ;
 120
 121 : 0.0
 122     [ 0 i->f ] literal ;
 123
 124 : 1.0
 125     [ 1 i->f ] literal ;
 126
 127 : -1.0
 128     [ -1 i->f ] literal ;
 129
 130 : 10.0
 131     [ 10 i->f ] literal ;
 132
 133 : flog10
 134     flog [ 10 i->f flog ] literal f/ ;
 135
 136 : fabs
 137     dup 0.0 f< if
 138         -1.0 f*
 139     then
 140 ;
 141
 142 : floor
 143     dup 0.0 f>= if
 144         dup 1.0 fmod f-
 145     else
 146         dup 1.0 fmod dup 0.0 <> if
 147             f- 1.0 f-
 148         else
 149             drop
 150         then
 151     then
 152 ;
 153
 154 : fhead ( float -- )
 155     floor f->i
 156     0 .R  ;
 157
 158 : ftail ( float prec -- )
 159     dup 0<= if 2drop exit then
 160
 161     swap
 162
 163     1.0 fmod 10.0 f*
 164
 165     dup floor f->i 0 .R
 166
 167     1.0 fmod dup 0.0 f> if
 168         swap 1- recurse
 169     else
 170         2drop
 171     then
 172 ;
 173
 174 variable precision
 175 16 precision !
 176
 177 : f.plain ( float -- )
 178
 179     dup 0.0 = if
 180         ." 0.0"
 181         drop exit
 182     then
 183
 184     dup 0.0 f< if
 185         [char] - emit
 186         -1.0 f*
 187     then
 188
 189     dup fhead
 190
 191     [char] . emit
 192
 193     precision @ over flog10 floor f->i -
 194     ftail
 195 ;
 196
 197 : f.scientific ( float -- )
 198     dup 0.0 = if
 199         ." 0.0"
 200         drop exit
 201     then
 202
 203     dup 0.0 f< if
 204         [char] - emit
 205         -1.0 f*
 206     then
 207
 208     dup flog10 floor dup -rot
 209     10.0 swap f^ f/ f.plain
 210     ." e" f->i 0 .R
 211 ;
 212
 213 : f. ( float -- )
 214     dup fabs dup 1000000 i->f f>= swap 1 i->f 10000 i->f f/ f< or if
 215         f.scientific
 216     else
 217         f.plain
 218     then
 219
 220     space
 221 ;