projects / pforth / blame
? search:
summary | tags | clone url | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge pull request #75 from SeekingMeaning/0BSD
[pforth] / fth / floats.fth
CommitLineData
8e9db35f
PB		 1\ @(#) floats.fth 98/02/26 1.4 17:51:40
2\ High Level Forth support for Floating Point
3\
4\ Author: Phil Burk and Darren Gibbs
1a088514 5\ Copyright 1994 3DO, Phil Burk, Larry Polansky, David Rosenboom
8e9db35f 6\
1f99f95d
S		 7\ Permission to use, copy, modify, and/or distribute this
8\ software for any purpose with or without fee is hereby granted.
9\
10\ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL
11\ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED
12\ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL
13\ THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR
14\ CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING
15\ FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF
16\ CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
17\ OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
8e9db35f
PB		 18\
19\ 19970702 PLB Drop 0.0 in REPRESENT to fix 0.0 F.
20\ 19980220 PLB Added FG. , fixed up large and small formatting
21\ 19980812 PLB Now don't drop 0.0 in REPRESENT to fix 0.0 F. (!!!)
22\ Fixed F~ by using (F.EXACTLY)
23
24ANEW TASK-FLOATS.FTH
25
26: FALIGNED ( addr -- a-addr )
27 1 floats 1- +
28 1 floats /
29 1 floats *
30;
31
32: FALIGN ( -- , align DP )
33 dp @ faligned dp !
34;
35
36\ account for size of create when aligning floats
37here
38create fp-create-size
39fp-create-size swap - constant CREATE_SIZE
40
41: FALIGN.CREATE ( -- , align DP for float after CREATE )
42 dp @
43 CREATE_SIZE +
44 faligned
45 CREATE_SIZE -
46 dp !
47;
48
49: FCREATE ( <name> -- , create with float aligned data )
50 falign.create
51 CREATE
52;
53
54: FVARIABLE ( <name> -- ) ( F: -- )
55 FCREATE 1 floats allot
56;
57
58: FCONSTANT
59 FCREATE here 1 floats allot f!
60 DOES> f@
61;
62
63: F0SP ( -- ) ( F: ? -- )
64 fdepth 0 max 0 ?DO fdrop LOOP
65;
66
67\ Convert between single precision and floating point
68: S>F ( s -- ) ( F: -- r )
69 s>d d>f
70;
71: F>S ( -- s ) ( F: r -- )
72 f>d d>s
73;
74
75: (F.EXACTLY) ( r1 r2 -f- flag , return true if encoded equally ) { | caddr1 caddr2 fsize fcells }
76 1 floats -> fsize
77 fsize cell 1- + cell 1- invert and \ round up to nearest multiple of stack size
78 cell / -> fcells ( number of cells per float )
79\ make room on data stack for floats data
80 fcells 0 ?DO 0 LOOP
81 sp@ -> caddr1
82 fcells 0 ?DO 0 LOOP
83 sp@ -> caddr2
84\ compare bit representation
85 caddr1 f!
86 caddr2 f!
87 caddr1 fsize caddr2 fsize compare 0=
88 >r fcells 2* 0 ?DO drop LOOP r> \ drop float bits
89;
90
91: F~ ( -0- flag ) ( r1 r2 r3 -f- )
92 fdup F0<
93 IF
94 frot frot ( -- r3 r1 r2 )
95 fover fover ( -- r3 r1 r2 r1 r2 )
96 f- fabs ( -- r3 r1 r2 |r1-r2| )
97 frot frot ( -- r3 |r1-r2| r1 r2 )
98 fabs fswap fabs f+ ( -- r3 |r1-r2| |r1|+|r2| )
99 frot fabs f* ( -- |r1-r2| |r1|+|r2|*|r3| )
100 f<
101 ELSE
102 fdup f0=
103 IF
104 fdrop
105 (f.exactly) \ f- f0= \ 19980812 Used to cheat. Now actually compares bit patterns.
106 ELSE
107 frot frot ( -- r3 r1 r2 )
108 f- fabs ( -- r3 |r1-r2| )
109 fswap f<
110 THEN
111 THEN
112;
113
114\ FP Output --------------------------------------------------------
115fvariable FVAR-REP \ scratch var for represent
116: REPRESENT { c-addr u | n flag1 flag2 -- n flag1 flag2 , FLOATING } ( F: r -- )
117 TRUE -> flag2 \ FIXME - need to check range
118 fvar-rep f!
119\
120 fvar-rep f@ f0<
121 IF
122 -1 -> flag1
123 fvar-rep f@ fabs fvar-rep f! \ absolute value
124 ELSE
125 0 -> flag1
126 THEN
127\
128 fvar-rep f@ f0=
129 IF
130\ fdrop \ 19970702 \ 19980812 Remove FDROP to fix "0.0 F."
131 c-addr u [char] 0 fill
132 0 -> n
133 ELSE
134 fvar-rep f@
135 flog
136 fdup f0< not
137 IF
138 1 s>f f+ \ round up exponent
139 THEN
140 f>s -> n
141\ ." REP - n = " n . cr
142\ normalize r to u digits
143 fvar-rep f@
144 10 s>f u n - s>f f** f*
145 1 s>f 2 s>f f/ f+ \ round result
146\
147\ convert float to double_int then convert to text
148 f>d
149\ ." REP - d = " over . dup . cr
150 <# u 1- 0 ?DO # loop #s #> \ ( -- addr cnt )
151\ Adjust exponent if rounding caused number of digits to increase.
152\ For example from 9999 to 10000.
153 u - +-> n
154 c-addr u move
155 THEN
156\
157 n flag1 flag2
158;
159
160variable FP-PRECISION
161
162\ Set maximum digits that are meaningful for the precision that we use.
1631 FLOATS 4 / 7 * constant FP_PRECISION_MAX
164
165: PRECISION ( -- u )
166 fp-precision @
167;
168: SET-PRECISION ( u -- )
169 fp_precision_max min
170 fp-precision !
171;
1727 set-precision
173
17432 constant FP_REPRESENT_SIZE
17564 constant FP_OUTPUT_SIZE
176
177create FP-REPRESENT-PAD FP_REPRESENT_SIZE allot \ used with REPRESENT
178create FP-OUTPUT-PAD FP_OUTPUT_SIZE allot \ used to assemble final output
179variable FP-OUTPUT-PTR \ points into FP-OUTPUT-PAD
180
181: FP.HOLD ( char -- , add char to output )
182 fp-output-ptr @ fp-output-pad 64 + <
183 IF
184 fp-output-ptr @ tuck c!
185 1+ fp-output-ptr !
186 ELSE
187 drop
188 THEN
189;
190: FP.APPEND { addr cnt -- , add string to output }
191 cnt 0 max 0
192 ?DO
193 addr i + c@ fp.hold
194 LOOP
195;
196
197: FP.STRIP.TRAILING.ZEROS ( -- , remove trailing zeros from fp output )
198 BEGIN
199 fp-output-ptr @ fp-output-pad u>
200 fp-output-ptr @ 1- c@ [char] 0 =
201 and
202 WHILE
203 -1 fp-output-ptr +!
204 REPEAT
205;
206
207: FP.APPEND.ZEROS ( numZeros -- )
208 0 max 0
209 ?DO [char] 0 fp.hold
210 LOOP
211;
212
213: FP.MOVE.DECIMAL { n prec -- , append with decimal point shifted }
214 fp-represent-pad n prec min fp.append
215 n prec - fp.append.zeros
216 [char] . fp.hold
217 fp-represent-pad n +
218 prec n - 0 max fp.append
219;
220
221: (EXP.) ( n -- addr cnt , convert exponent to two digit value )
222 dup abs 0
223 <# # #s
224 rot 0<
225 IF [char] - HOLD
226 ELSE [char] + hold
227 THEN
228 #>
229;
230
231: FP.REPRESENT ( -- n flag1 flag2 ) ( r -f- )
232;
233
234: (FS.) ( -- addr cnt ) ( F: r -- , scientific notation )
235 fp-output-pad fp-output-ptr ! \ setup pointer
236 fp-represent-pad precision represent
237\ ." (FS.) - represent " fp-represent-pad precision type cr
238 ( -- n flag1 flag2 )
239 IF
240 IF [char] - fp.hold
241 THEN
242 1 precision fp.move.decimal
243 [char] e fp.hold
244 1- (exp.) fp.append \ n
245 ELSE
246 2drop
247 s" <FP-OUT-OF-RANGE>" fp.append
248 THEN
249 fp-output-pad fp-output-ptr @ over -
250;
251
252: FS. ( F: r -- , scientific notation )
253 (fs.) type space
254;
255
256: (FE.) ( -- addr cnt ) ( F: r -- , engineering notation ) { | n n3 -- }
257 fp-output-pad fp-output-ptr ! \ setup pointer
258 fp-represent-pad precision represent
259 ( -- n flag1 flag2 )
260 IF
261 IF [char] - fp.hold
262 THEN
263\ convert exponent to multiple of three
264 -> n
265 n 1- s>d 3 fm/mod \ use floored divide
266 3 * -> n3
267 1+ precision fp.move.decimal \ amount to move decimal point
268 [char] e fp.hold
269 n3 (exp.) fp.append \ n
270 ELSE
271 2drop
272 s" <FP-OUT-OF-RANGE>" fp.append
273 THEN
274 fp-output-pad fp-output-ptr @ over -
275;
276
277: FE. ( F: r -- , engineering notation )
278 (FE.) type space
279;
280
281: (FG.) ( F: r -- , normal or scientific ) { | n n3 ndiff -- }
282 fp-output-pad fp-output-ptr ! \ setup pointer
283 fp-represent-pad precision represent
284 ( -- n flag1 flag2 )
285 IF
286 IF [char] - fp.hold
287 THEN
288\ compare n with precision to see whether we do scientific display
289 dup precision >
290 over -3 < OR
291 IF \ use exponential notation
292 1 precision fp.move.decimal
293 fp.strip.trailing.zeros
294 [char] e fp.hold
295 1- (exp.) fp.append \ n
296 ELSE
297 dup 0>
298 IF
299\ POSITIVE EXPONENT - place decimal point in middle
300 precision fp.move.decimal
301 ELSE
302\ NEGATIVE EXPONENT - use 0.000????
303 s" 0." fp.append
304\ output leading zeros
305 negate fp.append.zeros
306 fp-represent-pad precision fp.append
307 THEN
308 fp.strip.trailing.zeros
309 THEN
310 ELSE
311 2drop
312 s" <FP-OUT-OF-RANGE>" fp.append
313 THEN
314 fp-output-pad fp-output-ptr @ over -
315;
316
317: FG. ( F: r -- )
318 (fg.) type space
319;
320
321: (F.) ( F: r -- , normal or scientific ) { | n n3 ndiff prec' -- }
322 fp-output-pad fp-output-ptr ! \ setup pointer
323 fp-represent-pad \ place to put number
324 fdup flog 1 s>f f+ f>s precision max
325 fp_precision_max min dup -> prec'
326 represent
327 ( -- n flag1 flag2 )
328 IF
329\ add '-' sign if negative
330 IF [char] - fp.hold
331 THEN
332\ compare n with precision to see whether we must do scientific display
333 dup fp_precision_max >
334 IF \ use exponential notation
335 1 precision fp.move.decimal
336 fp.strip.trailing.zeros
337 [char] e fp.hold
338 1- (exp.) fp.append \ n
339 ELSE
340 dup 0>
341 IF
342 \ POSITIVE EXPONENT - place decimal point in middle
343 prec' fp.move.decimal
344 ELSE
345 \ NEGATIVE EXPONENT - use 0.000????
346 s" 0." fp.append
347 \ output leading zeros
348 dup negate precision min
349 fp.append.zeros
350 fp-represent-pad precision rot + fp.append
351 THEN
352 THEN
353 ELSE
354 2drop
355 s" <FP-OUT-OF-RANGE>" fp.append
356 THEN
357 fp-output-pad fp-output-ptr @ over -
358;
359
360: F. ( F: r -- )
361 (f.) type space
362;
363
364: F.S ( -- , print FP stack )
365 ." FP> "
366 fdepth 0>
367 IF
368 fdepth 0
369 DO
370 cr?
371 fdepth i - 1- \ index of next float
372 fpick f. cr?
373 LOOP
374 ELSE
375 ." empty"
376 THEN
377 cr
378;
379
380\ FP Input ----------------------------------------------------------
381variable FP-REQUIRE-E \ must we put an E in FP numbers?
382false fp-require-e ! \ violate ANSI !!
383
384: >FLOAT { c-addr u | dlo dhi u' fsign flag nshift -- flag }
385 u 0= IF false exit THEN
386 false -> flag
387 0 -> nshift
388\
389\ check for minus sign
390 c-addr c@ [char] - = dup -> fsign
391 c-addr c@ [char] + = OR
392 IF 1 +-> c-addr -1 +-> u \ skip char
393 THEN
394\
395\ convert first set of digits
396 0 0 c-addr u >number -> u' -> c-addr -> dhi -> dlo
397 u' 0>
398 IF
399\ convert optional second set of digits
400 c-addr c@ [char] . =
401 IF
402 dlo dhi c-addr 1+ u' 1- dup -> nshift >number
403 dup nshift - -> nshift
404 -> u' -> c-addr -> dhi -> dlo
405 THEN
406\ convert exponent
407 u' 0>
408 IF
409 c-addr c@ [char] E =
410 c-addr c@ [char] e = OR
411 IF
412 1 +-> c-addr -1 +-> u' \ skip E char
413 u' 0>
414 IF
415 c-addr c@ [char] + = \ ignore + on exponent
416 IF
417 1 +-> c-addr -1 +-> u' \ skip char
418 THEN
419 c-addr u' ((number?))
420 num_type_single =
421 IF
422 nshift + -> nshift
423 true -> flag
424 THEN
425 ELSE
426 true -> flag \ allow "1E"
427 THEN
428 THEN
429 ELSE
430\ only require E field if this variable is true
431 fp-require-e @ not -> flag
432 THEN
433 THEN
434\ convert double precision int to float
435 flag
436 IF
437 dlo dhi d>f
438 10 s>f nshift s>f f** f* \ apply exponent
439 fsign
440 IF
441 fnegate
442 THEN
443 THEN
444 flag
445;
446
4473 constant NUM_TYPE_FLOAT \ possible return type for NUMBER?
448
449: (FP.NUMBER?) ( $addr -- 0 | n 1 | d 2 | r 3 , convert string to number )
450\ check to see if it is a valid float, if not use old (NUMBER?)
451 dup count >float
452 IF
453 drop NUM_TYPE_FLOAT
454 ELSE
455 (number?)
456 THEN
457;
458
459defer fp.old.number?
460variable FP-IF-INIT
461
462: FP.TERM ( -- , deinstall fp conversion )
463 fp-if-init @
464 IF
465 what's fp.old.number? is number?
466 fp-if-init off
467 THEN
468;
469
470: FP.INIT ( -- , install FP converion )
471 fp.term
472 what's number? is fp.old.number?
473 ['] (fp.number?) is number?
474 fp-if-init on
475 ." Floating point numeric conversion installed." cr
476;
477
478FP.INIT
479if.forgotten fp.term
480
481
4820 [IF]
483
48423.8e-9 fconstant fsmall
4851.0 fsmall f- fconstant falmost1
486." Should be 1.0 = " falmost1 f. cr
487
488: TSEGF ( r -f- , print in all formats )
489." --------------------------------" cr
490 34 0
491 DO
492 fdup fs. 4 spaces fdup fe. 4 spaces
493 fdup fg. 4 spaces fdup f. cr
494 10.0 f/
495 LOOP
496 fdrop
497;
498
499: TFP
500 1.234e+22 tsegf
501 1.23456789e+22 tsegf
502 0.927 fsin 1.234e+22 f* tsegf
503;
504
505[THEN]
Fork of 'github/philburk/pforth', a portable ANS-like Forth written in ANSI C.