Discussion:
exercise in double number arithmetic
(too old to reply)
Krishna Myneni
2024-07-06 20:20:45 UTC
Permalink
I've been working on extending the kForth-64 User's Manual and, in
particular, illustrating double length arithmetic, which, being one of
the strengths of Forth, often does not get enough exposure. Here's an
exercise which you can do using only standard Forth words.

How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.

Please show your code.

--
Krishna Myneni
Ahmed
2024-07-06 21:59:35 UTC
Permalink
Why using double arithmetic when we can use simply Pascal triangle?


create table 50 50 * cells allot

table 50 50 * cells erase
: table.init 50 0 do 1 table i 50 * cells + ! loop ;
table.init
: table.calc 49 1 do 50 1 do table j 1 - 50 * i 1 - + cells + @
table j 1 - 50 * i + cells + @ +
table j 50 * i + cells + !
loop loop ;
table.calc

table 42 50 * 21 + cells + @ . 538257874440

Ahmed
Krishna Myneni
2024-07-06 23:50:28 UTC
Permalink
Post by Ahmed
Why using double arithmetic when we can use simply Pascal triangle?
...
Yes! I had forgotten about Pascal's triangle when I coded my solution. A
benefit of your double number solution is that it works on both 32-bit
and 64-bit Forth systems.

Here's my solution for the specific problem, which will only work on a
64-bit system:

: solution ( -- d)
1 s>d 43 22 do i 1 m*/ loop d>f
21 dup s>d rot 2 do I 1 m*/ loop d>f
f/ fround f>d ;

We don't have D/ in the standard, which necessitated the conversion to
float for dividing two double length numbers and then convert back to
double.

--
Krishna
Ahmed
2024-07-07 05:48:47 UTC
Permalink
...
Post by Krishna Myneni
: solution ( -- d)
1 s>d 43 22 do i 1 m*/ loop d>f
21 dup s>d rot 2 do I 1 m*/ loop d>f
f/ fround f>d ;
...
Post by Krishna Myneni
Krishna
You have not to use fround. (42!/21!) is an integer and is multiple of
21!.
in other words 42! is multiple of (21!)^2.

Ahmed
Krishna Myneni
2024-07-07 13:17:11 UTC
Permalink
Post by Ahmed
...
Post by Krishna Myneni
: solution ( -- d)
     1 s>d 43 22 do i 1 m*/ loop d>f
     21 dup s>d rot 2 do I 1 m*/ loop d>f
     f/ fround f>d ;
...
Post by Krishna Myneni
Krishna
You have not to use fround. (42!/21!) is an integer and is multiple of
21!.
in other words 42! is multiple of (21!)^2.
The FROUND is necessary because 42!/21! is too large of an integer to be
exactly represented by double-precision fp (IEEE format) and F>D is a
truncating operation, not rounding.

The whole D>F and floating point operations can be dispensed with if we
have D/MOD -- a double length version /MOD.

--
Krishna
Gerry Jackson
2024-07-07 14:33:37 UTC
Permalink
Post by Krishna Myneni
Post by Ahmed
...
Post by Krishna Myneni
: solution ( -- d)
     1 s>d 43 22 do i 1 m*/ loop d>f
     21 dup s>d rot 2 do I 1 m*/ loop d>f
     f/ fround f>d ;
...
Post by Krishna Myneni
Krishna
You have not to use fround. (42!/21!) is an integer and is multiple of
21!.
in other words 42! is multiple of (21!)^2.
The FROUND is necessary because 42!/21! is too large of an integer to be
exactly represented by double-precision fp (IEEE format) and F>D is a
truncating operation, not rounding.
The whole D>F and floating point operations can be dispensed with if we
have D/MOD -- a double length version /MOD.
--
Krishna
There is a +D/MOD in this link

http://www3.cs.stonybrook.edu/~algorith/implement/random-number/distrib/r250.seq
--
Gerry
Krishna Myneni
2024-07-08 12:30:11 UTC
Permalink
...
Post by Gerry Jackson
Post by Krishna Myneni
The whole D>F and floating point operations can be dispensed with if
we have D/MOD -- a double length version /MOD.
--
Krishna
There is a +D/MOD in this link
http://www3.cs.stonybrook.edu/~algorith/implement/random-number/distrib/r250.seq
Thanks, Gerry. It looks like D/MOD is not essential to computing the
binomial coefficients with double length numbers. We can use M*/ for
single computations, or simply D+ to build up a table of coefficients.
But, it's good to have a reference for D/MOD.

I remember that I wanted to update the FSL R250 code to generate 64-bit
random numbers, but I never got around to it. There is also a
permutations and combinations module in FSL which may benefit from what
has been demonstrated in this thread.

--
Krishna
Gerry Jackson
2024-07-10 07:51:02 UTC
Permalink
Post by Krishna Myneni
...
Post by Gerry Jackson
Post by Krishna Myneni
The whole D>F and floating point operations can be dispensed with if
we have D/MOD -- a double length version /MOD.
--
Krishna
There is a +D/MOD in this link
http://www3.cs.stonybrook.edu/~algorith/implement/random-number/distrib/r250.seq
Thanks, Gerry. It looks like D/MOD is not essential to computing the
binomial coefficients with double length numbers. We can use M*/ for
single computations, or simply D+ to build up a table of coefficients.
But, it's good to have a reference for D/MOD.
I remember that I wanted to update the FSL R250 code to generate 64-bit
random numbers, but I never got around to it. There is also a
permutations and combinations module in FSL which may benefit from what
has been demonstrated in this thread.
--
Krishna
I didn't realise the paper I linked to was in the FSL as I've never had
much to do with the FSL. I found I had another definition for D/MOD from
Forth Dimensions Vol 14 Issue 6, p 27 "Math - Who Needs It"
https://www.forth.org/fd/FD-V14N6.pdf
I've no idea how good it is
--
Gerry
dxf
2024-07-10 08:41:03 UTC
Permalink
[...] I found I had another definition for D/MOD from Forth Dimensions Vol 14 Issue 6, p 27 "Math - Who Needs It"
https://www.forth.org/fd/FD-V14N6.pdf
I've no idea how good it is
Here's the source file (32MATH.SEQ) for those interested:

https://pastebin.com/mfU8FZ1x
Krishna Myneni
2024-07-10 11:05:11 UTC
Permalink
Post by dxf
[...] I found I had another definition for D/MOD from Forth Dimensions Vol 14 Issue 6, p 27 "Math - Who Needs It"
https://www.forth.org/fd/FD-V14N6.pdf
I've no idea how good it is
https://pastebin.com/mfU8FZ1x
Nice. Thanks.

--
KM
dxf
2024-07-10 13:08:59 UTC
Permalink
Post by dxf
[...] I found I had another definition for D/MOD from Forth Dimensions Vol 14 Issue 6, p 27 "Math - Who Needs It"
https://www.forth.org/fd/FD-V14N6.pdf
I've no idea how good it is
https://pastebin.com/mfU8FZ1x
Here's one I had lying around. Not sure if I used it so may require checking
first! The 2variable is an eye-sore. Perhaps someone can eliminate it without
resorting to locals :)

2variable d

\ Divide quad by double. Unsigned.
: DUM/MOD ( uq ud -- udrem udquot )
d 2! [ 16 cells ] literal 0 do
dup >r 2swap dup >r d2* 2swap d2*
r> 0< dup d- 2dup d 2@ du< 0= r> 0< or
if d 2@ d- 2swap 1 0 d+ 2swap then
loop 2swap ;

\ Divide doubles. Unsigned.
: DU/MOD ( ud1 ud2 -- udrem udquot )
0 0 2swap dum/mod ;

\ Divide doubles. Signed. Symmetric.
: D/MOD ( d1 d2 -- drem dquot )
2 pick 2dup xor 2>r dabs 2swap dabs
2swap du/mod r> 0< if dnegate then
r> 0< if 2swap dnegate 2swap then ;
Krishna Myneni
2024-07-11 12:47:11 UTC
Permalink
Post by dxf
Post by dxf
[...] I found I had another definition for D/MOD from Forth Dimensions Vol 14 Issue 6, p 27 "Math - Who Needs It"
https://www.forth.org/fd/FD-V14N6.pdf
I've no idea how good it is
https://pastebin.com/mfU8FZ1x
Here's one I had lying around. Not sure if I used it so may require checking
first! The 2variable is an eye-sore. Perhaps someone can eliminate it without
resorting to locals :)
2variable d
\ Divide quad by double. Unsigned.
: DUM/MOD ( uq ud -- udrem udquot )
d 2! [ 16 cells ] literal 0 do
dup >r 2swap dup >r d2* 2swap d2*
loop 2swap ;
\ Divide doubles. Unsigned.
: DU/MOD ( ud1 ud2 -- udrem udquot )
0 0 2swap dum/mod ;
\ Divide doubles. Signed. Symmetric.
: D/MOD ( d1 d2 -- drem dquot )
2 pick 2dup xor 2>r dabs 2swap dabs
2swap du/mod r> 0< if dnegate then
r> 0< if 2swap dnegate 2swap then ;
This is interesting. I had recently implemented the non-standard word
UD/MOD based on the corresponding word in Gforth. Thus, I have the
double length division words,

UM/MOD ( ud u -- urem uquot )
standard word for dividing unsigned double length by single length;
return unsigned single length remainder and quotient.

UD/MOD ( ud u -- urem udquot )
non-standard word for dividing unsigned double length by single;
return unsigned single remainder and double quotient.

UM/MOD can have an overflow, while UD/MOD cannot overflow, but neither
of the two words above can divide two double length numbers. The DUM/MOD
and D/MOD words do what I'm looking for, in principle -- I don't know if
they actually work. We need some tests. Also, the word naming doesn't
seem to follow any real convention and the operation of each word can't
be discerned from the name.

--
Krishna
minforth
2024-07-11 15:18:47 UTC
Permalink
Had them somewhere on my drive, but have no access
to my system right now, so could not double-check
whether this is still correct:

\ Check double symmetric and floored division

1 8 CELLS 1- LSHIFT CONSTANT MSB
MSB INVERT CONSTANT MINT

T{ 0. 1. d/rem -> 0. 0. }T
T{ 2. 3. d/rem -> 2. 0. }T
T{ 10. 7. d/rem -> 3. 1. }T
T{ -10. 7. d/rem -> -3. -1. }T
T{ 10. -7. d/rem -> 3. -1. }T
T{ -10. -7. d/rem -> -3. 1. }T
T{ 0 1 2. d/rem -> 0. msb 0 }T
T{ 0 mint 2. d/rem -> 0. msb mint 2/ }T
T{ -1 mint 2dup d/rem -> 0. 1. }T
T{ 1 mint 0 mint d/rem -> 1. 1. }T
T{ 0 mint 1 mint d/rem -> 0 mint 0. }T
T{ -1 mint 2. d/rem -> 1. -1 mint 1 rshift }T

T{ 10. 7. d/mod -> 3. 1. }T
T{ -10. 7. d/mod -> 4. -2. }T
T{ 10. -7. d/mod -> -4. -2. }T
T{ -10. -7. d/mod -> -3. 1. }T

T{ 0. 1. d/ -> 0. }T
T{ 20. 7. d/ -> 2. }T
T{ -20. 7. d/ -> -2. }T
T{ 20. -7. d/ -> -2. }T
T{ -20. -7. d/ -> 2. }T
Krishna Myneni
2024-07-12 12:04:17 UTC
Permalink
Post by minforth
Had them somewhere on my drive, but have no access
to my system right now, so could not double-check
\ Check double symmetric and floored division
1 8 CELLS 1- LSHIFT CONSTANT MSB
MSB INVERT CONSTANT MINT
T{ 0. 1. d/rem -> 0. 0. }T
T{ 2. 3. d/rem -> 2. 0. }T
T{ 10. 7. d/rem -> 3. 1. }T
T{ -10. 7. d/rem -> -3. -1. }T
T{ 10. -7. d/rem -> 3. -1. }T
T{ -10. -7. d/rem -> -3. 1. }T
T{ 0 1 2. d/rem -> 0. msb 0 }T
T{ 0 mint 2. d/rem -> 0. msb mint 2/ }T
T{ -1 mint 2dup d/rem -> 0. 1. }T
T{ 1 mint 0 mint d/rem -> 1. 1. }T
T{ 0 mint 1 mint d/rem -> 0 mint 0. }T
T{ -1 mint 2. d/rem -> 1. -1 mint 1 rshift }T
T{ 10. 7. d/mod -> 3. 1. }T
T{ -10. 7. d/mod -> 4. -2. }T
T{ 10. -7. d/mod -> -4. -2. }T
T{ -10. -7. d/mod -> -3. 1. }T
T{ 0. 1. d/ -> 0. }T
T{ 20. 7. d/ -> 2. }T
T{ -20. 7. d/ -> -2. }T
T{ 20. -7. d/ -> -2. }T
T{ -20. -7. d/ -> 2. }T
I ran your tests for D/REM and D/ against the source definitions
provided by dxf. dxf's D/MOD is your D/REM (I assume your D/MOD is a
floored division word). I renamed dxf's D/MOD to D/REM for your tests
and commented out your D/MOD tests.

On symmetric division systems, the definition of D/REM below should be
reverted to D/MOD as in his original posting

Below is the code and test results, under kForth-64.

--
Krishna

=== begin double-number-division-tests.4th ===
\ Check double symmetric and floored division
\ Posted by minforth on comp.lang.forth on 11 July 2024

include ans-words
include ttester

\ double number division words posted on comp.lang.forth
\ by dxf on 10 July 2024; renamed D/MOD to D/REM
\
\ DUM/MOD ( uq ud -- udrem udquot )
\ DU/MOD ( ud1 ud2 -- udrem udquot )
\ D/REM ( d1 d2 -- drem dquot ) ( originally D/MOD )
\ D/ ( d1 d2 -- dquot )

2variable d

\ Divide quad by double. Unsigned.
: DUM/MOD ( uq ud -- udrem udquot )
d 2! [ 16 cells ] literal 0 do
dup >r 2swap dup >r d2* 2swap d2*
r> 0< dup d- 2dup d 2@ du< 0= r> 0< or
if d 2@ d- 2swap 1 0 d+ 2swap then
loop 2swap ;

\ Divide doubles. Unsigned.
: DU/MOD ( ud1 ud2 -- udrem udquot )
0 0 2swap dum/mod ;

\ Divide doubles. Signed. Symmetric; originally D/MOD
: D/REM ( d1 d2 -- drem dquot )
2 pick 2dup xor 2>r dabs 2swap dabs
2swap du/mod r> 0< if dnegate then
r> 0< if 2swap dnegate 2swap then ;

\ Divide doubles. Signed symmetric
: D/ ( d1 d2 -- dquot )
D/REM 2swap 2drop ;

1 8 CELLS 1- LSHIFT CONSTANT MSB
MSB INVERT CONSTANT MINT

TESTING D/REM
T{ 0 s>d 1 s>d d/rem -> 0 s>d 0 s>d }T
T{ 2 s>d 3 s>d d/rem -> 2 s>d 0 s>d }T
T{ 10 s>d 7 s>d d/rem -> 3 s>d 1 s>d }T
T{ -10 s>d 7 s>d d/rem -> -3 s>d -1 s>d }T
T{ 10 s>d -7 s>d d/rem -> 3 s>d -1 s>d }T
T{ -10 s>d -7 s>d d/rem -> -3 s>d 1 s>d }T
T{ 0 1 2 s>d d/rem -> 0 s>d msb 0 }T
T{ 0 mint 2 s>d d/rem -> 0 s>d msb mint 2/ }T
T{ -1 mint 2dup d/rem -> 0 s>d 1 s>d }T
T{ 1 mint 0 mint d/rem -> 1 s>d 1 s>d }T
T{ 0 mint 1 mint d/rem -> 0 mint 0 s>d }T
T{ -1 mint 2 s>d d/rem -> 1 s>d -1 mint 1 rshift }T

0 [IF]
TESTING D/MOD
T{ 10 s>d 7 s>d d/mod -> 3 s>d 1 s>d }T
T{ -10 s>d 7 s>d d/mod -> 4 s>d -2 s>d }T
T{ 10 s>d -7 s>d d/mod -> -4 s>d -2 s>d }T
T{ -10 s>d -7 s>d d/mod -> -3 s>d 1 s>d }T
[THEN]

TESTING D/
T{ 0 s>d 1 s>d d/ -> 0 s>d }T
T{ 20 s>d 7 s>d d/ -> 2 s>d }T
T{ -20 s>d 7 s>d d/ -> -2 s>d }T
T{ 20 s>d -7 s>d d/ -> -2 s>d }T
T{ -20 s>d -7 s>d d/ -> 2 s>d }T

=== end double-number-division-tests.4th ===

=== begin execution output ===
include double-number-division-tests
TESTING D/REM
TESTING D/
ok
=== end execution output ===
minforth
2024-07-12 14:43:50 UTC
Permalink
DUM/MOD quoted by dxf is the classic looping shift-subtract
algorithm. OTOH the double division from FSL #23 quoted by
Gerry does it without looping and should be faster. But I
have no time to check whether both algorithms are really
equivalent.
Ahmed
2024-07-07 16:23:11 UTC
Permalink
Post by Krishna Myneni
Post by Ahmed
...
Post by Krishna Myneni
: solution ( -- d)
     1 s>d 43 22 do i 1 m*/ loop d>f
     21 dup s>d rot 2 do I 1 m*/ loop d>f
     f/ fround f>d ;
...
Post by Krishna Myneni
Krishna
You have not to use fround. (42!/21!) is an integer and is multiple of
21!.
in other words 42! is multiple of (21!)^2.
The FROUND is necessary because 42!/21! is too large of an integer to be
exactly represented by double-precision fp (IEEE format) and F>D is a
truncating operation, not rounding.
The whole D>F and floating point operations can be dispensed with if we
have D/MOD -- a double length version /MOD.
--
Krishna
Hi,
You are right.

Your solution with fround gives

: solution ( -- d)
1 s>d 43 22 do i 1 m*/ loop d>f
21 dup s>d rot 2 do I 1 m*/ loop d>f
f/ fround f>d ;
solution d. 538257874440 ok ( exact solution)

your solution without fround gives

: solution2 ( -- d)
1 s>d 43 22 do i 1 m*/ loop d>f
21 dup s>d rot 2 do I 1 m*/ loop d>f
f/ f>d ;
solution2 d. 538257874439 ok ( false)

And thanks for the explanation.

Ahmed
minforth
2024-07-07 06:32:24 UTC
Permalink
: BCOEF ( n k -- )
swap s>f 1e 1+ 1 ?DO
fover 1e f+ i s>f f- i s>f f/ f*
LOOP fswap fdrop fround f>d d. ;
Ahmed
2024-07-06 22:07:36 UTC
Permalink
Here is a version with double numbers
and Pascal triangle


create table 50 50 2* * cells allot


table 50 50 2* * cells erase
: table.init 50 0 do 1. table i 50 * cells + 2! loop ;
table.init
: table.calc 49 1 do 50 1 do table j 1 - 50 * i 1 - 2* + cells + 2@
table j 1 - 50 * i 2* + cells + 2@ d+
table j 50 * i 2* + cells + 2!
loop loop ;
table.calc

table 42 50 * 21 2* + cells + 2@ d.


It gives 538257874440. (the same result)

Ahmed
Krishna Myneni
2024-07-06 23:42:25 UTC
Permalink
Post by Ahmed
Here is a version with double numbers
and Pascal triangle
create table 50 50 2* * cells allot
table 50 50 2* * cells erase
: table.init 50 0 do 1. table i 50 * cells + 2! loop ;
table.init
                            table j     50 * i     2* + cells + 2!
loop loop ;
table.calc
It gives  538257874440. (the same result)
This version works on both 32-bit and 64-bit Forth systems.

--
Krishna
Ahmed
2024-07-07 08:00:37 UTC
Permalink
I noticed there was a problem in my double artihmetic version of the
program.
It is about table's elements addressing.
for example, it gives 2 for 26 0 istead of 1, and the error spreads for
the next calculations

Here is the new version. It works fine for n<=49 (it is due to my choice
of 50 lines of Pascal's triangle)



create table 50 50 2* * cells allot


table 50 50 * 2* cells erase
: table.init 50 0 do 1. table i 50 * 2* cells + 2! loop ;
table.init
: table.calc 50 1 do 50 1 do table j 1 - 50 * i 1 - + 2* cells + 2@
table j 1 - 50 * i + 2* cells + 2@ d+
table j 50 * i + 2* cells + 2!
loop loop ;
table.calc


: bcoef_tab swap 50 * swap + 2* cells table + 2@ d. ;


Here a comaprison with results given table version and minforth verson
side by side, the results are identical.

: BCOEF ( n k -- )
swap s>f 1e 1+ 1 ?DO
fover 1e f+ i s>f f- i s>f f/ f*
LOOP fswap fdrop fround f>d d. ;

minforth version works for any n and k. It uses the fact that:

C(n,k) = prod((n+1-i)/(i)) for i = 1 to k
= ((n+1-k)/(k))...(n/1) = (n!/(n-k)!)/k!

: go 50 0 do i 1+ 0 ?do cr j . i . j i bcoef_tab j i bcoef loop loop ;


Ahmed
a***@spenarnc.xs4all.nl
2024-07-07 11:24:35 UTC
Permalink
Post by Krishna Myneni
I've been working on extending the kForth-64 User's Manual and, in
particular, illustrating double length arithmetic, which, being one of
the strengths of Forth, often does not get enough exposure. Here's an
exercise which you can do using only standard Forth words.
How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.
Please show your code.
It is in the special function screen:
0 ( CHS TRI PYR SQ CUB FIB ) \ AvdH B8jan21
2 \ For N M return "N OVER M " (N/M)
3 : CHS >R R@ - 1 R> 1+ 1 ?DO OVER I + I */ LOOP NIP ;
4 \ '(./.) ALIAS CHS
5
6 \ For x return its TRIANGLE, number
7 : TRI DUP 1+ * 2/ ;
8 \ For x return its PYRAMIDAL number
9 ( : PYR ; )
10 \ For x return its SQUARE number
11 : SQ DUP * ;
12 \ For x return its CUBE number
13 : CUB DUP DUP * * ;
14 \ For x return the xth Fibonacci .
15 : FIB >R 0 1 R> 0 ?DO SWAP OVER + LOOP DROP ;

S[ ] OK 42 21 chs
S[ 538257874440 ] OK
Post by Krishna Myneni
--
Krishna Myneni
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
mhx
2024-07-07 15:41:54 UTC
Permalink
Post by a***@spenarnc.xs4all.nl
Post by Krishna Myneni
How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.
Please show your code.
[..]
Post by a***@spenarnc.xs4all.nl
2 \ For N M return "N OVER M " (N/M)
[..]

Really great!

FORTH> 42 21 CHS . 538257874440 ok
FORTH> 66 33 CHS . 7219428434016265740 ok
FORTH> 68 34 CHS . ( ... )

-marcel
Ahmed
2024-07-07 16:13:33 UTC
Permalink
Post by mhx
Post by a***@spenarnc.xs4all.nl
Post by Krishna Myneni
How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.
Please show your code.
[..]
Post by a***@spenarnc.xs4all.nl
2 \ For N M return "N OVER M " (N/M)
[..]
Really great!
FORTH> 42 21 CHS . 538257874440 ok
FORTH> 66 33 CHS . 7219428434016265740 ok
FORTH> 68 34 CHS . ( ... )
-marcel
Hi,

100 value N

create table N N * 2* cells allot


table N N * 2* cells erase
: table.init N 0 do 1. table i N * 2* cells + 2! loop ;
table.init
: table.calc N 1 do N 1 do table j 1 - N * i 1 - + 2* cells + 2@
table j 1 - N * i + 2* cells + 2@ d+
table j N * i + 2* cells + 2!
loop loop ;
table.calc

table 68 N * 34 + 2* cells + 2@ d. 28453041475240576740 ok

Wolfram alpha gives C(68,43) = 28453041475240576740


minforth version gives: 68 34 bcoef 28453041475240574976 ok (it is
different)


Ahmed
mhx
2024-07-07 16:45:19 UTC
Permalink
[..]
Post by Ahmed
Post by a***@spenarnc.xs4all.nl
2 \ For N M return "N OVER M " (N/M)
[..]
Post by Ahmed
table.calc
Wolfram alpha gives C(68,43) = 28453041475240576740
minforth version gives: 68 34 bcoef 28453041475240574976 ok (it is
different)
Interesting! I don't immediately see where that goes wrong...
FORTH> : test 67 1 do i i 2/ CHS CR I . dup H. space . loop ;
FORTH> test
1 $00000001 1
2 $00000002 2
3 $00000003 3
4 $00000006 6
5 $0000000A 10
6 $00000014 20
7 $00000023 35
8 $00000046 70
9 $0000007E 126
10 $000000FC 252
11 $000001CE 462
12 $0000039C 924
13 $000006B4 1716
14 $00000D68 3432
15 $00001923 6435
16 $00003246 12870
17 $00005EF6 24310
18 $0000BDEC 48620
19 $000168DA 92378
20 $0002D1B4 184756
21 $000561CC 352716
22 $000AC398 705432
23 $0014A18E 1352078
24 $0029431C 2704156
25 $004F59AC 5200300
26 $009EB358 10400600
27 $013210BC 20058300
28 $02642178 40116600
29 $049F73E8 77558760
30 $093EE7D0 155117520
31 $11E9E123 300540195
32 $23D3C246 601080390
33 $458C00A6 1166803110
34 $8B18014C 2333606220
35 $000000010E75C9A2 4537567650
36 $000000021CEB9344 9075135300
37 $000000041D5EF65C 17672631900
38 $000000083ABDECB8 35345263800
39 $000000100C258D9A 68923264410
40 $00000020184B1B34 137846528820
41 $0000003EA955AF04 269128937220
42 $0000007D52AB5E08 538257874440
43 $000000F4F3092084 1052049481860
44 $000001E9E6124108 2104098963720
45 $000003BE7F5B5DD8 4116715363800
46 $0000077CFEB6BBB0 8233430727600
47 $00000EAA1D7B2F8E 16123801841550
48 $00001D543AF65F1C 32247603683100
49 $0000397C21A572BC 63205303218876
50 $000072F8434AE578 126410606437752
51 $0000E18483FF3844 247959266474052
52 $0001C30907FE7088 495918532948104
53 $0003755D946EB6F8 973469712824056
54 $0006EABB28DD6DF0 1946939425648112
55 $000D9638C720AA3C 3824345300380220
56 $001B2C718E415478 7648690600760440
57 $0035690281893C18 15033633249770520
58 $006AD20503127830 30067266499541040
59 $00D2148152D785F8 59132290782430712
60 $01A42902A5AF0BF0 118264581564861424
61 $033AC44F881661D0 232714176627630544
62 $0675889F102CC3A0 465428353255261088
63 $0CB764F927D82123 916312070471295267
64 $196EC9F24FB04246 1832624140942590534
65 $321847F48D727306 3609714217008132870
66 $64308FE91AE4E60C 7219428434016265740 ok
( higher gives integer divide by 0 )

-marcel
minforth
2024-07-07 18:15:14 UTC
Permalink
Post by Ahmed
minforth version gives: 68 34 bcoef 28453041475240574976 ok (it is
different)
No surprise here. Float64 can't represent more than about 16 decimal
significand digits (mantissa has 52 bits)
Ahmed
2024-07-07 18:27:51 UTC
Permalink
Ah,
I see.
Thanks for the explanation.

Ahmed
Ahmed
2024-07-07 18:27:02 UTC
Permalink
Hi again,

Here is another version, we haven't to save all the Pascal's triangle as
before.
Just two rows are sufficient.
Also, we can use locals (n, k) or values (n, k).


100 value M

create coeffs 2 M * 2* cells allot



: coeffs.init
coeffs 2 M * 2* cells erase
2 0 do 1. coeffs i M * 2* cells + 2! loop ;


0 value k
0 value n


: coef.calc { n k -- }
\ to k to n
coeffs.init
n 1+ 1 do
k 1+ 1 do
coeffs j 1- 2 mod M * i 1- + 2* cells + 2@
coeffs j 1- 2 mod M * i + 2* cells + 2@
d+
coeffs j 2 mod M * i + 2* cells + 2!
loop
loop
coeffs n 2 mod M * k + 2* cells + 2@ d.
;

42 21 coef.calc 538257874440 ok
68 34 coef.calc 28453041475240576740 ok


Ahmed
Ahmed
2024-07-07 18:59:30 UTC
Permalink
The last version I posted gives the same results as WolframAlpha to
C(n,k) for n up to 130.
130 65 coef.calc 95067625827960698145584333020095113100
ok
WolframAlpha gives C(130,65)=95067625827960698145584333020095113100

For n>130, my results are wrong.


Notice that I used M=100 the number of columns of coeffs array, so one
have to take M=130 as max value for M.
And we have
130 90 coef.calc 5334728207404302438396453128740400 ok
WolframAlpha gives C(130,90)=5334728207404302438396453128740400


When we take M=130,

100 100 coef.calc 1 ok (true)

130 128 coef.calc 8385 ok (true)
130 129 coef.calc 130 ok (true)
130 130 coef.calc 2 ok (false) this is the extreme case, must be 1.
And after this case begin the wrong results.

Ahmed
mhx
2024-07-07 20:52:43 UTC
Permalink
OK, let's extend it a little.

: CHS ( M N -- M!/{M-N}!/N! )
DUP local N
- local M-N
1. N 1+ 1 DO M-N I + I M*/ LOOP ;

FORTH> : test 131 2 do i i 2/ CHS CR I . 2dup H. H. space d. loop ;
ok
FORTH> test
2 $00000000$00000002 2
3 $00000000$00000003 3
4 $00000000$00000006 6
5 $00000000$0000000A 10
6 $00000000$00000014 20
7 $00000000$00000023 35
8 $00000000$00000046 70
9 $00000000$0000007E 126
10 $00000000$000000FC 252
11 $00000000$000001CE 462
12 $00000000$0000039C 924
13 $00000000$000006B4 1716
14 $00000000$00000D68 3432
15 $00000000$00001923 6435
16 $00000000$00003246 12870
17 $00000000$00005EF6 24310
18 $00000000$0000BDEC 48620
19 $00000000$000168DA 92378
20 $00000000$0002D1B4 184756
21 $00000000$000561CC 352716
22 $00000000$000AC398 705432
23 $00000000$0014A18E 1352078
24 $00000000$0029431C 2704156
25 $00000000$004F59AC 5200300
26 $00000000$009EB358 10400600
27 $00000000$013210BC 20058300
28 $00000000$02642178 40116600
29 $00000000$049F73E8 77558760
30 $00000000$093EE7D0 155117520
31 $00000000$11E9E123 300540195
32 $00000000$23D3C246 601080390
33 $00000000$458C00A6 1166803110
34 $00000000$8B18014C 2333606220
35 $00000000$000000010E75C9A2 4537567650
36 $00000000$000000021CEB9344 9075135300
37 $00000000$000000041D5EF65C 17672631900
38 $00000000$000000083ABDECB8 35345263800
39 $00000000$000000100C258D9A 68923264410
40 $00000000$00000020184B1B34 137846528820
41 $00000000$0000003EA955AF04 269128937220
42 $00000000$0000007D52AB5E08 538257874440
43 $00000000$000000F4F3092084 1052049481860
44 $00000000$000001E9E6124108 2104098963720
45 $00000000$000003BE7F5B5DD8 4116715363800
46 $00000000$0000077CFEB6BBB0 8233430727600
47 $00000000$00000EAA1D7B2F8E 16123801841550
48 $00000000$00001D543AF65F1C 32247603683100
49 $00000000$0000397C21A572BC 63205303218876
50 $00000000$000072F8434AE578 126410606437752
51 $00000000$0000E18483FF3844 247959266474052
52 $00000000$0001C30907FE7088 495918532948104
53 $00000000$0003755D946EB6F8 973469712824056
54 $00000000$0006EABB28DD6DF0 1946939425648112
55 $00000000$000D9638C720AA3C 3824345300380220
56 $00000000$001B2C718E415478 7648690600760440
57 $00000000$0035690281893C18 15033633249770520
58 $00000000$006AD20503127830 30067266499541040
59 $00000000$00D2148152D785F8 59132290782430712
60 $00000000$01A42902A5AF0BF0 118264581564861424
61 $00000000$033AC44F881661D0 232714176627630544
62 $00000000$0675889F102CC3A0 465428353255261088
63 $00000000$0CB764F927D82123 916312070471295267
64 $00000000$196EC9F24FB04246 1832624140942590534
65 $00000000$321847F48D727306 3609714217008132870
66 $00000000$64308FE91AE4E60C 7219428434016265740
67 $00000000$C56EC13C4B95E372 14226520737620288370
68 $00000001$8ADD8278972BC6E4 28453041475240576740
69 $00000003$0A72DCA497BCB3FC 56093138908331422716
70 $00000006$14E5B9492F7967F8 112186277816662845432
71 $0000000B$FE8C2D6CC84BE262 221256270138418389602
72 $00000017$FD185AD99097C4C4 442512540276836779204
73 $0000002F$5436F86EFAAEE514 873065282167813104916
74 $0000005E$A86DF0DDF55DCA28 1746130564335626209832
75 $000000BA$D329D4A89A2BA334 3446310324346630677300
76 $00000175$A653A95134574668 6892620648693261354600
77 $000002E1$B7FA82CE4684ED78 13608507434599516007800
78 $000005C3$6FF5059C8D09DAF0 27217014869199032015600
79 $00000B61$FD1D84AEC9C0439A 53753604366668088230810
80 $000016C3$FA3B095D93808734 107507208733336176461620
81 $00002CF9$CF237667B3042A54 212392290424395860814420
82 $000059F3$9E46ECCF660854A8 424784580848791721628840
83 $0000B1C2$F5BCEC5CE81CA74C 839455243105945545123660
84 $00016385$EB79D8B9D0394E98 1678910486211891090247320
85 $0002BEC7$3CA376D483CA9568 3318776542511877736535400
86 $00057D8E$7946EDA907952AD0 6637553085023755473070800
87 $000ADB2B$29FACA4866440904 13124252690842425594480900
88 $0015B656$53F59490CC881208 26248505381684851188961800
89 $002AF127$E4A16FC90BFC0CE8 51913710643776705684835560
90 $0055E24F$C942DF9217F819D0 103827421287553411369671120
91 $00A9E6A8$F7E2E6CD8875F048 205397724721029574666088520
92 $0153CD51$EFC5CD9B10EBE090 410795449442059149332177040
93 $02A05FCD$B450EDFC5D65CCB0 812850570172585125274307760
94 $0540BF9B$68A1DBF8BACB9960 1625701140345170250548615520
95 $0A657B38$E9C058B19C5D9F8E 3217533506933149454210801550
96 $14CAF671$D380B16338BB3F1C 6435067013866298908421603100
97 $29294B20$05F44F7B4682585C 12738806129490428451365214300
98 $52529640$0BE89EF68D04B0B8 25477612258980856902730428600
99 $A2FFAE9D$88381C06E4042AB4 50445672272782096667406248628
100 $0000000145FF5D3B$1070380DC8085568 100891344545564193334812497256
101 $00000002859A5942$C633922554ED5DD8 199804427433372226016001220056
102 $000000050B34B285$8C67244AA9DABBB0 399608854866744452032002440112
103 $00000009FD94B061$24DFFE0A0B84F3C4 791532924062974587678774064068
104 $00000013FB2960C2$49BFFC141709E788 1583065848125949175357548128136
105 $0000002795CF8F63$EDE1C7EDD6B304A8 3136262529306125724764953838760
106 $0000004F2B9F1EC7$DBC38FDBAD660950 6272525058612251449529907677520
107 $0000009CDFEAB382$88CA9D0D5C53AA28 12428892245768720464809261509160
108 $00000139BFD56705$11953A1AB8A75450 24857784491537440929618523018320
109 $0000026DCB4E7D0A$0B92CB96B3C4A270 49263609265046928387789436527216
110 $000004DB969CFA14$1725972D678944E0 98527218530093856775578873054432
111 $000009A0F8404B1E$ADE15DF0DAF0163C 195295022443578894680165266232892
112 $00001341F080963D$5BC2BBE1B5E02C78 390590044887157789360330532465784
113 $0000262D6385A799$143A3119491F38B8 774327632846470705223111406467256
114 $00004C5AC70B4F32$28746232923E7170
1548655265692941410446222812934512
115 $000097648AA81432$E647DD2F4E1AB4C8
3070609578529107968988200404956360
116 $00012EC915502865$CC8FBA5E9C356990
6141219157058215937976400809912720
117 $0002587062ABF954$B85E1ACCF9061F70
12178349853827309571919303301013360
118 $0004B0E0C557F2A9$70BC3599F20C3EE0
24356699707654619143838606602026720
119 $00094DBDCBAA29D0$0E86593E200FC0F8
48307454420181661301946569760686328
120 $00129B7B975453A0$1D0CB27C401F81F0
96614908840363322603893139521372656
121 $0024E8E02C2D90E5$7892E424A0C4CB30
191645966716130525165099506263706416
122 $0049D1C0585B21CA$F125C84941899660
383291933432261050330199012527412832
123 $009272B343EE99C0$07B22E5FC8361DF0
760401738905937245009910944207609328
124 $0124E56687DD3380$0F645CBF906C3BE0
1520803477811874490019821888415218656
125 $0245249EBC4D3D8C$4F4D3A0E5F91ABA0
3017467217880703353213932318284164000
126 $048A493D789A7B18$9E9A741CBF235740
6034934435761406706427864636568328000
127 $09026955FB528C44$DABA7E690B4A2123
11975573020964041433067793888190275875
128 $1204D2ABF6A51889$B574FCD216944246
23951146041928082866135587776380551750
129 $23C2ADEAF15F4854$40BCDA2EBA9873C6
47533812913980349072792166510047556550
130 $47855BD5E2BE90A8$8179B45D7530E78C
95067625827960698145584333020095113100 ok

FORTH> 131 dup 2/ chs cr ud.
188694833082770476622296176145946360850
( correct cf. WolframAlpha, but 132 66 CHS is incorrect )

-marcel
Ahmed
2024-07-07 21:41:27 UTC
Permalink
On Sun, 7 Jul 2024 20:52:43 (UTC), mhx wrote:

..
Post by mhx
FORTH> 131 dup 2/ chs cr ud.
188694833082770476622296176145946360850
( correct cf. WolframAlpha, but 132 66 CHS is incorrect )
-marcel
Same results when using ud. instead of d.

200 value M

create coeffs 2 M * 2* cells allot


: coeffs.init
coeffs 2 M * 2* cells erase
2 0 do 1. coeffs i M * 2* cells + 2! loop ;


0 value k
0 value n


: coef.calc { n k -- }
\ to k to n
coeffs.init
n 1+ 1 do
k 1+ 1 do
coeffs j 1- 2 mod M * i 1- + 2* cells + 2@
coeffs j 1- 2 mod M * i + 2* cells + 2@
d+
coeffs j 2 mod M * i + 2* cells + 2!
loop
loop
coeffs n 2 mod M * k + 2* cells + 2@ ud.
;



131 65 coef.calc 188694833082770476622296176145946360850 ok (correct)
132 66 coef.calc 37107299244602489781217744860124510244 ok (incorrect)


Ahmed
Krishna Myneni
2024-07-08 01:02:37 UTC
Permalink
Post by mhx
OK, let's extend it a little.
: CHS ( M N -- M!/{M-N}!/N! )
    DUP local  N
     -  local M-N
    1.  N 1+ 1 DO  M-N I +  I  M*/  LOOP ;
Great! My non-locals version is

\ Compute binomial coefficient as a double length number
: binom ( n k -- d )
dup >r - 1 s>d
r> 1+ 1 ?DO
2 pick I + I m*/
LOOP rot drop ;

130 65 binom d.
95067625827960698145584333020095113100 ok

This provides a nice example for use of double length number
computation, thanks to Albert.

As a side benefit to this exercise, I discovered that my division
overflow check in kForth's M*/ (and in the non-standard UTM/ used by
M*/) was not correct. This has been fixed in kForth-64.

--
Krishna
Ahmed
2024-07-08 02:00:47 UTC
Permalink
On Mon, 8 Jul 2024 1:02:37 (UTC), Krishna Myneni wrote:
..
Post by Krishna Myneni
130 65 binom d.
95067625827960698145584333020095113100 ok
..
Post by Krishna Myneni
Krishna
And it works for n=131, k=65
131 65 binom ud. 188694833082770476622296176145946360850 ok (correct,
see WolframAlpha)
But for n>=132 it gives uncorrect results.

Good work. This thread was fruitful.

Ahmed
minforth
2024-07-08 11:49:50 UTC
Permalink
Post by Ahmed
And it works for n=131, k=65
131 65 binom ud. 188694833082770476622296176145946360850 ok
Good work. This thread was fruitful.
Good to know how many different ways one can choose 65 distinct
objects from a set of 131, at a time.

Given the age of the universe with 13.787 billion years, this
must be a tremendous show every second.

A bit repetitive, agreed, but always new. :o)
Krishna Myneni
2024-07-08 12:25:23 UTC
Permalink
Post by Ahmed
..
Post by Krishna Myneni
130 65 binom d.
95067625827960698145584333020095113100  ok
..
Post by Krishna Myneni
Krishna
And it works for n=131, k=65
131 65 binom ud. 188694833082770476622296176145946360850  ok (correct,
see WolframAlpha)
But for n>=132 it gives uncorrect results.
I reproduce the same results for those cases here.
Post by Ahmed
Good work. This thread was fruitful.
All of the different solutions have utility, particularly the Pascal
triangle method when a number of binomial coefficients must be
available. I recall that modern Fortran has ways to specify
non-rectangular arrays. We should be able to define a triangle array
structure in Forth as well, to reduce the storage required for n = 1 to
nmax array of coefficients.

--
Krishna
Ahmed
2024-07-08 14:00:53 UTC
Permalink
On Mon, 8 Jul 2024 12:25:23 +0000, Krishna Myneni wrote:
...
Post by Krishna Myneni
All of the different solutions have utility, particularly the Pascal
triangle method when a number of binomial coefficients must be
available. I recall that modern Fortran has ways to specify
non-rectangular arrays. We should be able to define a triangle array
structure in Forth as well, to reduce the storage required for n = 1 to
nmax array of coefficients.
--
Krishna
Here is a version that uses Pascal's triangle saved in a table (array)
with minimal size.
The table has just 8778 double cells (132*133/2 = 8778 double cells) and
not (132*132 double cells).
The table is structured as below:
the number of double cell: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
16 17 18 19 20 ...
The content of the double cell: 1 1 1 1 2 1 1 3 3 1 1 4 6 4 1 1
5 10 10 5 1 ...
n: 0 1 1 2 2 2 3 3 3 3 4 4 4 4 4 5
5 5 5 5 5
p: 0 0 1 0 1 2 0 1 2 3 0 1 2 3 4 0
1 2 3 4 5 ...

The code begins here.

132 value M
M dup 1+ * 2/ value table_size (132*133/2 = 8778 double cells)

create binomial_coeffs table_size 2* cells allot

: th ( n p -- r) swap dup 1+ * 2/ + ;

: bc_th th 2* cells binomial_coeffs + ;
: bc_th! bc_th 2! ;
: bc_th@ bc_th 2@ ;

: coeffs.init
binomial_coeffs table_size 2* cells erase
M 0 do 1. i 0 bc_th! loop ;

: coeffs.calc
coeffs.init
M 1 do
M 1 do
j 1- i 1- bc_th@
j 1- i < if 0. else j 1- i bc_th@ then d+
j i bc_th!
loop
loop
;


: get_binomial_coeff ( n p -- c) bc_th@ ;
: binomial bc_th@ ;

coeffs.calc

Here, the code ends.

Some calculations:

table_size . 8778 ok

42 21 binomial d. 538257874440 ok
68 34 binomial d. 28453041475240576740 ok
130 65 binomial d. 95067625827960698145584333020095113100 ok
131 65 binomial d. -151587533838167986841078431285821850606 ok
131 65 binomial ud. 188694833082770476622296176145946360850 ok (notice
ud. here)

Ahmed
Marc Olschok
2024-07-14 21:43:07 UTC
Permalink
Post by Krishna Myneni
I've been working on extending the kForth-64 User's Manual and, in
particular, illustrating double length arithmetic, which, being one of
the strengths of Forth, often does not get enough exposure. Here's an
exercise which you can do using only standard Forth words.
How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.
Yes, M*/ comes in handy for C(n,0) = 1 , C(n+1,k+1) = C(n,k)*n/k

42 21 binom d.
gives 538257874440

where

: binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;

Gforth SEE nicely replaced the original 'recurse' with 'binom' for
better readability.
--
M.O.
Marc Olschok
2024-07-14 21:49:10 UTC
Permalink
Post by Marc Olschok
[...]
Yes, M*/ comes in handy for C(n,0) = 1 , C(n+1,k+1) = C(n,k)*n/k
_________________________________________________________^^^^^^^^^^
of course this shoul[d read C(n,k)*(n+1)/(k+1)
--
M.O.
Krishna Myneni
2024-07-15 01:15:10 UTC
Permalink
...
Post by Marc Olschok
Post by Krishna Myneni
How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.
Yes, M*/ comes in handy for C(n,0) = 1 , C(n+1,k+1) = C(n,k)*n/k
42 21 binom d.
gives 538257874440
where
: binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
Gforth SEE nicely replaced the original 'recurse' with 'binom' for
better readability.
THank you for the recursive version. It's nice to have both looping and
recursive examples.

There's a reason why RECURSE (or equivalent) is preferable to having the
name of the word in the output of SEE in Forth. This is because it is
possible to have an earlier definition with the same name and to call it
from within the definition e.g.

: binom ... ;

: binom ... binom ... ;

In the later definition of BINOM Forth requires the call to BINOM be the
earlier definition if it exists in the search order. If it does not
exist in the search order, then I believe the standard would require an
error to occur.

Code such as shown in your output of SEE cannot be compiled when no
earlier definition of BINOM exists. For this reason, it is not helpful
for SEE to replace RECURSE with the word name.

The SEE in kForth shows the following:

see binom
5631665B9710 DUP
5631665B9711 0=
5631665B9712 JZ $1D
5631665B971B 2DROP
5631665B971C #1
5631665B9725 S>D
5631665B9726 JMP $1A
5631665B972F 2DUP
5631665B9730 1-
5631665B9731 SWAP
5631665B9732 1-
5631665B9733 SWAP
5631665B9734 $5631665B9710
5631665B973D EXECUTE-BC
5631665B973E 2SWAP
5631665B973F M*/
5631665B9740 RET
ok

A slightly smarter SEE would replace the address and the EXECUTE-BC
virtual machine instruction with RECURSE.

--
Krishna
minforth
2024-07-15 07:38:50 UTC
Permalink
Yet another fast variant without looping/recursion:

\ LOGBINOM

100 FLOATS BUFFER: LCS{ \ log sum array
: } ( a i -- f: n ) floats + f@ ;
: }! ( a i f: n -- ) floats + f! ;
: PREP-LCS ( -- ) \ prepare array
0e lcs{ 0 }! 100 1 DO
lcs{ i 1- } i s>f flog f+ lcs{ i }!
LOOP ;

: LBINOM ( n k -- binom ) \ calc binom coeff by array lookup
2dup = over 0= or
IF 0e
ELSE 2dup lcs{ swap } lcs{ swap } - lcs{ swap } f- fswap f-
THEN falog fround f>d ;

PREP-LCS
42 21 LBINOM D.
minforth
2024-07-15 07:58:51 UTC
Permalink
oversight: IF 2drop 0e
Ahmed
2024-07-15 08:12:48 UTC
Permalink
Hi, thanks for sharing.
I noticed that it works fine for C(n,p) where 0<=n<48 and 0<=p<=n.
But for n = 48, it works fine for p<=21. and when 21<p<=n, it gives
different results than WolframAlpha.
When n>=49 the results are false.

It's nice to know why and what are the limits.

Ahmed
minforth
2024-07-15 08:45:43 UTC
Permalink
The limitiation is within FALOG that works with float64
i.e. 52 significant bits ~ about 16 decimal digits only.

A Forth with DFLOAT=float80 (x87) works with 64 signif. bits
i.e. about 21 decimal digits.

If this is a problem, the Pascal triangle does not have
such limitations.
minforth
2024-07-15 08:55:00 UTC
Permalink
Another way to improve results for high n's and k's would be
to use Kahan summation when filling the CSL array.
Ahmed
2024-07-15 08:57:31 UTC
Permalink
Post by minforth
The limitiation is within FALOG that works with float64
i.e. 52 significant bits ~ about 16 decimal digits only.
A Forth with DFLOAT=float80 (x87) works with 64 signif. bits
i.e. about 21 decimal digits.
If this is a problem, the Pascal triangle does not have
such limitations.
Thanks for the explanation.

Ahmed
a***@spenarnc.xs4all.nl
2024-07-15 10:46:23 UTC
Permalink
Post by minforth
The limitiation is within FALOG that works with float64
i.e. 52 significant bits ~ about 16 decimal digits only.
A Forth with DFLOAT=float80 (x87) works with 64 signif. bits
i.e. about 21 decimal digits.
If this is a problem, the Pascal triangle does not have
such limitations.
For these large numbers you should use the floating point
approximation for factorials.
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
minforth
2024-07-15 10:56:42 UTC
Permalink
Is the Stirling approximation really precise enough in order
to compute binomial coefficients to the last decimnal digit?
a***@spenarnc.xs4all.nl
2024-07-15 17:46:28 UTC
Permalink
Post by minforth
Is the Stirling approximation really precise enough in order
to compute binomial coefficients to the last decimnal digit?
Better a mistake in the last digit than a totally bogus result.
And you expect unaccuracy.

Groetjes Albert
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
Anton Ertl
2024-07-15 14:00:40 UTC
Permalink
Post by Krishna Myneni
There's a reason why RECURSE (or equivalent) is preferable to having the
name of the word in the output of SEE in Forth. This is because it is
possible to have an earlier definition with the same name and to call it
from within the definition e.g.
: binom ... ;
: binom ... binom ... ;
Yes, in Gforth SEE produces the same output for

: foo ;
: foo foo ;
see foo

and for

: foo recursive foo ;
see foo

Gforth's SEE also does not tell you that BAR is shadowed in the
following example:

: bar 1 . ;
: foo bar ;
: bar 2 . ;
see foo

Nor does Gforth's SEE tell you that a word it calls is in a wordlist
that is not in a search order.

Sometimes I think about giving some indication of such issues, but up
to now these ideas are pretty low on my todo list, because these
things cause little pain.

If you want to see the source code, use LOCATE (in the development
version).

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2024: https://euro.theforth.net
Gerry Jackson
2024-07-15 22:39:37 UTC
Permalink
Post by Krishna Myneni
...
Post by Krishna Myneni
How many different ways can you choose 42 distinct objects, 21 at a
time? This is "n choose k" or the binomial coefficent.
Yes, M*/ comes in handy for  C(n,0) = 1 , C(n+1,k+1) = C(n,k)*n/k
42 21 binom d.
gives  538257874440
where
: binom  ( n1 n2 -- nd ) \ n k --> C(n,k)
   dup 0=
   IF     2drop 1 s>d
   ELSE   2dup 1- swap 1- swap binom 2swap m*/
   THEN ;
Gforth SEE nicely replaced the original 'recurse' with 'binom' for
better readability.
THank you for the recursive version. It's nice to have both looping and
recursive examples.
There's a reason why RECURSE (or equivalent) is preferable to having the
name of the word in the output of SEE in Forth. This is because it is
possible to have an earlier definition with the same name and to call it
from within the definition e.g.
: binom ... ;
: binom ... binom ... ;
In the later definition of BINOM Forth requires the call to BINOM be the
earlier definition if it exists in the search order. If it does not
exist in the search order, then I believe the standard would require an
error to occur.
That's why the simplest way to achieve recursion by name is:

synonym binom recurse
: binom ... binom ... ;
but Gforth's RECURSIVE is better
--
Gerry
Marc Olschok
2024-07-31 21:59:09 UTC
Permalink
Post by Krishna Myneni
[...]
THank you for the recursive version. It's nice to have both looping and
recursive examples.
On second thought it might be nicer to make it tail-recursive.
This way one can see the comparison with the equivalent loop.
Of course one needs to accumulate the intermediate products in
ascending order to preserve integrality.

\ d , k , n-k , n --> d*(n+1)/(k+1) , k+1 , n-k+1 , n
: mstep ( d n1 n2 n3 -- d n1 n2 n3 )
Post by Krishna Myneni
r 1+ swap 1+ 2dup >r >r m*/ r> r> swap r> ;
\ d , 0 , n-k , n --> d*C(n,k) , k , n , n
: binom2_rec ( d n1 n2 -- d n1 n2 ) recursive
2dup < if mstep binom2_rec then ;

\ d , 0 , n-k , n --> d*C(n,k) , k , n , n
: binom2_loop ( d n1 n2 -- d n1 n2 )
2dup swap - 0 ?do mstep loop ;

: binom2 binom2_rec ;

\ n , k --> C(n,k)
: binom ( n1 n2 -- n )
Post by Krishna Myneni
r >r 1 s>d 0 r> dup r> - swap binom2 2drop drop ;
There's a reason why RECURSE (or equivalent) is preferable to having the
name of the word in the output of SEE in Forth. This is because it is
possible to have an earlier definition with the same name and to call it
from within the definition e.g.
I remember this feature from long ago (I think it was mentioned in
'Starting Forth'), but it seems that since the edit/compile/run cycle
became so fast I had never used it.
--
M.O.
a***@spenarnc.xs4all.nl
2024-07-15 10:41:54 UTC
Permalink
In article <v71gpb$jsug$***@solani.org>,
Marc Olschok <***@nowhere.invalid> wrote:
<SNIP>
Post by Marc Olschok
: binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
Gforth SEE nicely replaced the original 'recurse' with 'binom' for
better readability.
That is not smartness. It is quite an effort to discover that it
is actually RECURSE.
The ciforth decompiler does the same.
You can no longer paste the source recovered with impunity.
So I prefer:

:F binom ;

:R binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;

In my efforts to Make A Lisp (a github https://github.com/kanaka/mal )
I discovered that using recurse is an ugly cludge that present
a lot of problems in refactoring code, if not prevent it.
Forward and resolve definitions is the more sane method, cf. c.
It is hardly more complicated.
Post by Marc Olschok
M.O.
Groetjes Albert
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
minforth
2024-07-15 13:19:59 UTC
Permalink
Post by a***@spenarnc.xs4all.nl
:F binom ;
:R binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
In my efforts to Make A Lisp (a github https://github.com/kanaka/mal )
I discovered that using recurse is an ugly cludge that present
a lot of problems in refactoring code, if not prevent it.
Forward and resolve definitions is the more sane method, cf. c.
It is hardly more complicated.
IIRC gforth has RECURSIVE to avoid duplicating definitions.
Anton Ertl
2024-07-15 13:29:17 UTC
Permalink
Post by minforth
Post by a***@spenarnc.xs4all.nl
:F binom ;
:R binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
In my efforts to Make A Lisp (a github https://github.com/kanaka/mal )
I discovered that using recurse is an ugly cludge that present
a lot of problems in refactoring code, if not prevent it.
Forward and resolve definitions is the more sane method, cf. c.
It is hardly more complicated.
IIRC gforth has RECURSIVE to avoid duplicating definitions.
Yes. So for a direct recursion like this one you can write

: binom ( n1 n2 -- nd ) recursive \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;

RECURSIVE also allows you to tick the word in its own definition (not
possible with RECURSE), a feature that I actually have used; something
along the lines:

: walk ( node -- ) recursive
dup is-leaf? if
... \ do something for the leaf
else \ the node has subnodes
node ['] walk for-each-subnode
then ;

Gforth (development) also has FORWARD, which also allows you to do
mutual recursion, e.g.

forward foo

: bar ( n -- ) dup . 1- foo ;

: foo ( n -- ) dup 0> if bar else drop then ;

5 foo \ prints "5 4 3 2 1 "

let's see what the decompiler says:

simple-see foo
$7F54E256E870 dup 0->0
$7F54E256E878 0> 0->0
$7F54E256E880 ?branch 0->0
$7F54E256E888 <foo+$40>
$7F54E256E890 call 0->0
$7F54E256E898 bar
$7F54E256E8A0 branch 0->0
$7F54E256E8A8 <foo+$48>
$7F54E256E8B0 drop 0->0
$7F54E256E8B8 ;s 0->0 ok
simple-see bar
$7F54E256E810 dup 0->0
$7F54E256E818 call 0->0
$7F54E256E820 .
$7F54E256E828 1- 0->0
$7F54E256E830 call 0->0
$7F54E256E838 foo
$7F54E256E840 ;s 0->0 ok

$7F54E256E838 @ hex. \ prints $7F54E256E870

The last like shows that the call to FOO inside BAR directly calls the
FOO defined above, no DEFER or somesuch involved in this case.

The definition of FORWARD is quite intricate and uses several recent
features of Gforth. Read all about it in
<***@mips.complang.tuwien.ac.at>.

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2024: https://euro.theforth.net
Marc Olschok
2024-07-15 14:45:28 UTC
Permalink
Post by Anton Ertl
Post by minforth
Post by a***@spenarnc.xs4all.nl
:F binom ;
:R binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
In my efforts to Make A Lisp (a github https://github.com/kanaka/mal )
I discovered that using recurse is an ugly cludge that present
a lot of problems in refactoring code, if not prevent it.
Forward and resolve definitions is the more sane method, cf. c.
It is hardly more complicated.
IIRC gforth has RECURSIVE to avoid duplicating definitions.
Yes. So for a direct recursion like this one you can write
: binom ( n1 n2 -- nd ) recursive \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
Oh that is nice, I did not know about that. And it also avoids the
source paste problem that Albert noted.
Post by Anton Ertl
RECURSIVE also allows you to tick the word in its own definition (not
possible with RECURSE), a feature that I actually have used; something
: walk ( node -- ) recursive
dup is-leaf? if
... \ do something for the leaf
else \ the node has subnodes
node ['] walk for-each-subnode
then ;
Gforth (development) also has FORWARD, which also allows you to do
mutual recursion, e.g.
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
simple-see foo
$7F54E256E870 dup 0->0
$7F54E256E878 0> 0->0
$7F54E256E880 ?branch 0->0
$7F54E256E888 <foo+$40>
$7F54E256E890 call 0->0
$7F54E256E898 bar
$7F54E256E8A0 branch 0->0
$7F54E256E8A8 <foo+$48>
$7F54E256E8B0 drop 0->0
$7F54E256E8B8 ;s 0->0 ok
simple-see bar
$7F54E256E810 dup 0->0
$7F54E256E818 call 0->0
$7F54E256E820 .
$7F54E256E828 1- 0->0
$7F54E256E830 call 0->0
$7F54E256E838 foo
$7F54E256E840 ;s 0->0 ok
The last like shows that the call to FOO inside BAR directly calls the
FOO defined above, no DEFER or somesuch involved in this case.
The definition of FORWARD is quite intricate and uses several recent
features of Gforth. Read all about it in
Will these parts eventually go into future standards?
--
M.O.
Ruvim
2024-07-15 19:37:34 UTC
Permalink
Post by Marc Olschok
Post by Anton Ertl
Post by minforth
Post by a***@spenarnc.xs4all.nl
:F binom ;
:R binom ( n1 n2 -- nd ) \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
In my efforts to Make A Lisp (a github https://github.com/kanaka/mal )
I discovered that using recurse is an ugly cludge that present
a lot of problems in refactoring code, if not prevent it.
Forward and resolve definitions is the more sane method, cf. c.
It is hardly more complicated.
IIRC gforth has RECURSIVE to avoid duplicating definitions.
Yes. So for a direct recursion like this one you can write
: binom ( n1 n2 -- nd ) recursive \ n k --> C(n,k)
dup 0=
IF 2drop 1 s>d
ELSE 2dup 1- swap 1- swap binom 2swap m*/
THEN ;
Oh that is nice, I did not know about that. And it also avoids the
source paste problem that Albert noted.
I think, the word "recursive" is confusing. It looks like an ordinary
word, but behaves like an immediate word that changes the current word
list, and it does not add any behavior to the current definition.


Some better variants:

: binom [recursive] ... binom ... ;

recursive
: binom ... binom ... ;


rec: binom ... binom ... ;



An even more better and more general approach:

: binom ... forward binom ... ;


: binom ... forward:binom ... ;

The later one can also be used as:

: foo ... ['] forward:foo ... ;


An advantage of this approach is that it is more general than the
"recursive" word, and you also don't have to repeat the definition name.
Post by Marc Olschok
Post by Anton Ertl
RECURSIVE also allows you to tick the word in its own definition (not
possible with RECURSE), a feature that I actually have used;
I think, there should be a standard method to get the xt of the current
definition (regardless whether it is a named definition, or nameless
definition).
Post by Marc Olschok
Post by Anton Ertl
: walk ( node -- ) recursive
dup is-leaf? if
... \ do something for the leaf
else \ the node has subnodes
node ['] walk for-each-subnode
then ;
Gforth (development) also has FORWARD, which also allows you to do
mutual recursion, e.g.
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
simple-see foo
$7F54E256E870 dup 0->0
$7F54E256E878 0> 0->0
$7F54E256E880 ?branch 0->0
$7F54E256E888 <foo+$40>
$7F54E256E890 call 0->0
$7F54E256E898 bar
$7F54E256E8A0 branch 0->0
$7F54E256E8A8 <foo+$48>
$7F54E256E8B0 drop 0->0
$7F54E256E8B8 ;s 0->0 ok
simple-see bar
$7F54E256E810 dup 0->0
$7F54E256E818 call 0->0
$7F54E256E820 .
$7F54E256E828 1- 0->0
$7F54E256E830 call 0->0
$7F54E256E838 foo
$7F54E256E840 ;s 0->0 ok
The last like shows that the call to FOO inside BAR directly calls the
FOO defined above, no DEFER or somesuch involved in this case.
The definition of FORWARD is quite intricate and uses several recent
features of Gforth. Read all about it in
Will these parts eventually go into future standards?
--
Ruvim
Gerry Jackson
2024-07-15 22:41:44 UTC
Permalink
Post by Ruvim
Post by Anton Ertl
RECURSIVE also allows you to tick the word in its own definition (not
possible with RECURSE), a feature that I actually have used;
I think, there should be a standard method to get the xt of the current
definition (regardless whether it is a named definition, or nameless
definition).
It can be done by using DEFER as a forward definition
e.g.
defer foo
:noname ... ['] foo defer@ ... ; is foo

using DEFER@ gives the xt of the code, omittimg it gives the xt of the name.

as FOO can be called by name by executing
synoname foo recurse
I would guess that your suggestion of FORWARD FOO could be defined using
that and something like EXECUTE-PARSING e.g. copying "FOO RECURSE" to a
buffer and doing:
BUF COUNT ' SYNONYM EXECUTE-PARSING
followed by IMMEDIATE of course

As DEFER can be used as a forward definition it can also be used for
mutual recursion
--
Gerry
Gerry Jackson
2024-07-15 22:44:05 UTC
Permalink
Post by Gerry Jackson
Post by Ruvim
Post by Anton Ertl
RECURSIVE also allows you to tick the word in its own definition (not
possible with RECURSE), a feature that I actually have used;
I think, there should be a standard method to get the xt of the
current definition (regardless whether it is a named definition, or
nameless definition).
It can be done by using DEFER as a forward definition
e.g.
defer foo
as FOO can be called by name by executing
   synoname foo recurse
^^^^^^^^
Sorry a typo
synonym foo recurse
Post by Gerry Jackson
I would guess that your suggestion of FORWARD FOO could be defined using
that and something like EXECUTE-PARSING e.g. copying "FOO RECURSE" to a
  BUF COUNT ' SYNONYM EXECUTE-PARSING
followed by IMMEDIATE of course
As DEFER can be used as a forward definition it can also be used for
mutual recursion
--
Gerry
Gerry Jackson
2024-07-16 06:36:44 UTC
Permalink
Post by Gerry Jackson
I would guess that your suggestion of FORWARD FOO could be defined using
that and something like EXECUTE-PARSING e.g. copying "FOO RECURSE" to a
   BUF COUNT ' SYNONYM EXECUTE-PARSING
followed by IMMEDIATE of course
Sorry, ignore this as a standard system won't permit a definition within
a colon definition. My excuse - it was late last night and I was tired
when I posted it, but it would work in my Forth system.
--
Gerry
Ruvim
2024-07-16 07:40:07 UTC
Permalink
Post by Gerry Jackson
Post by Ruvim
Post by Anton Ertl
RECURSIVE also allows you to tick the word in its own definition (not
possible with RECURSE), a feature that I actually have used;
I think, there should be a standard method to get the xt of the
current definition (regardless whether it is a named definition, or
nameless definition).
It can be done by using DEFER as a forward definition
e.g.
defer foo
Or this way:

variable _foo
:noname ... _foo @ ... ; _foo !


Of course, in some particular cases this can be done in a standard way.

But what if we are creating the definition programmatically? Things
become more cumbersome without intrinsic necessity.


And this is not an API — this is not a method of getting the xt of the
current definition in the general case.

In the general case, you don't know how compilation of the current
definition was started, and you cannot require a special way of starting
compilation of the definition.
In this use case the named definition (created by "defer") is only
auxiliary, its xt is not needed.


[...]
Post by Gerry Jackson
As DEFER can be used as a forward definition it can also be used for
mutual recursion
This is possible but, awkward. Why else do some systems provide the
words like "recursive" and "forward"?


--
Ruvim
dxf
2024-07-16 08:08:18 UTC
Permalink
[...]
As DEFER can be used as a forward definition it can also be used for mutual recursion
This is possible but, awkward. Why else do some systems provide the words like "recursive" and "forward"?
"In Forth if you don't like it you can change it." ?
mhx
2024-07-16 07:40:41 UTC
Permalink
Post by Ruvim
I think, there should be a standard method to get the xt of the current
definition (regardless whether it is a named definition, or nameless
definition).
It would be nice if there was a repository somewhere with
non-standard tricks available in existing systems. To make it
more than a curiosity it should show examples of use (for the
specific system).

In general I am not surprised that anything imaginable can be
done in Forth. *Why* somebody would want do implement it is
much more interesting. That is still a long way from needing
standardization.

-marcel
Ruvim
2024-07-16 07:55:29 UTC
Permalink
Post by Ruvim
Post by Marc Olschok
Post by minforth
Post by a***@spenarnc.xs4all.nl
:F binom ;
:R binom  ( n1 n2 -- nd ) \ n k --> C(n,k)
   dup 0=
   IF     2drop 1 s>d
   ELSE   2dup 1- swap 1- swap binom 2swap m*/
   THEN ;
In my efforts to Make A Lisp (a github https://github.com/kanaka/mal )
I discovered that using recurse is an ugly cludge that present
a lot of problems in refactoring code, if not prevent it.
Forward and resolve definitions is the more sane method, cf. c.
It is hardly more complicated.
IIRC gforth has RECURSIVE to avoid duplicating definitions.
Yes.  So for a direct recursion like this one you can write
: binom  ( n1 n2 -- nd ) recursive \ n k --> C(n,k)
   dup 0=
   IF     2drop 1 s>d
   ELSE   2dup 1- swap 1- swap binom 2swap m*/
   THEN ;
Oh that is nice, I did not know about that. And it also avoids the
source paste problem that Albert noted.
I think, the word "recursive" is confusing. It looks like an ordinary
word, but behaves like an immediate word that changes the current word
list, and it does not add any behavior to the current definition.
  : binom [recursive] ... binom ... ;
  recursive
  : binom ... binom ... ;
  rec: binom ... binom ... ;
  : binom ... forward binom ... ;
  : binom ... forward:binom ... ;
  : foo ... ['] forward:foo ... ;
An advantage of this approach is that it is more general than the
"recursive" word, and you also don't have to repeat the definition name.
One more advantage of this approach is that it works as expected when
the compilation word list is absent in the search order. Because the
reference in a definition is resolved once the forwarded word name is
placed in the same compilation word list (regardless of the search order).


In the former approach (via variants of "recursive") you can get bugs
that seem very strange at the first glance.


--
Ruvim
Ruvim
2024-07-16 11:52:20 UTC
Permalink
Post by Ruvim
   : binom ... forward binom ... ;
   : binom ... forward:binom ... ;
   : foo ... ['] forward:foo ... ;
An advantage of this approach is that it is more general than the
"recursive" word, and you also don't have to repeat the definition name.
One more advantage of this approach is that it works as expected when
the compilation word list is absent in the search order. Because the
reference in a definition is resolved once the forwarded word name is
placed in the same compilation word list (regardless of the search order).
A low-level API (porcelain) to implement this functionality can be
represented by the following words:

forward-lit, ( -- x.lit-sys )
forward-call, ( -- x.call-sys )

compile-lit ( x x.lit-sys -- )
compile-call ( xt x.call-sys -- )


Or which names are better?


Conceptually, the standard words can be defined as:

: compile, ( xt -- ) forward-call, compile-call ;
: lit, ( x -- ) forward-lit, compile-lit ;
: literal state @ if lit, then ; immediate



--
Ruvim
a***@spenarnc.xs4all.nl
2024-07-15 18:16:11 UTC
Permalink
In article <***@mips.complang.tuwien.ac.at>,
Anton Ertl <***@mips.complang.tuwien.ac.at> wrote:
<SNIP>
Post by Anton Ertl
Gforth (development) also has FORWARD, which also allows you to do
mutual recursion, e.g.
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
<SNIP>
Post by Anton Ertl
The definition of FORWARD is quite intricate and uses several recent
features of Gforth. Read all about it in
In a simple Forth like ciforth this is much easier:

\ Could use plain : , :F marks the definition as forward.
': ALIAS :F
\ Resolve an earlier dummy definition for recursion.
: :R >IN @ NAME FOUND >R >IN ! : LATEST >DFA @ R> >DFA ! ;

We look up the previous definition and store it in the return stack.
After compilation, the data-field-address (pointing to high
level interpreted code) of the latest definition is copied to
the data-field-address of the forward definition header (dea).

[ NAME FOUND is approximately WORD FIND . ]
Post by Anton Ertl
- anton
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
sjack
2024-07-16 12:02:39 UTC
Permalink
Post by Anton Ertl
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
Toad code using Wil Baden type local macros ( /MM: MM ):

"dup . 1- myself" /mm: bar OK

: foo dup 0> if mm bar else drop then ; OK

5 foo 5 4 3 2 1 OK

Note: The macro BAR does not remain in the system image.
--
me
sjack
2024-07-16 12:53:00 UTC
Permalink
sjack <***@dontemail.me> wrote:

Shortly after I posted I recalled that some have an aversion to
the use of macros and the main point of the example code, I
assume, is the use of a hard coded fragment that includes some
variable that can somehow bind to the entry of the word(s) in
which it will be used.

--
me
a***@spenarnc.xs4all.nl
2024-07-17 10:11:33 UTC
Permalink
Post by sjack
Shortly after I posted I recalled that some have an aversion to
the use of macros and the main point of the example code, I
assume, is the use of a hard coded fragment that includes some
variable that can somehow bind to the entry of the word(s) in
which it will be used.
I like macro's, especially if the return stack is involved,
that otherwise prevents factoring.
In fact Forth words are macro's, `:I is a word that inlines its
code (different from `: ) but it makes no sense to invent weird
syntax around it.

WANT :I >IN

:I save >IN @ >R ;

:I restore R> >IN ! ;

: VARIABLE save NAME ." defining " TYPE CR restore
VARIABLE ;
: CONSTANT save NAME ." defining " TYPE CR restore
CONSTANT ;
VARIABLE aap
12 CONSTANT mies

output:
defining aap
defining mies

[ :2 :I :F :R fits comfortably in one screen. ]
Post by sjack
--
me
Groetjes Albert
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
minforth
2024-07-16 13:02:41 UTC
Permalink
Pretzel coding. I wondered what that could be useful for.
FWIW:
https://stackoverflow.com/questions/2725038/are-there-any-example-of-mutual-recursion
sjack
2024-07-16 19:39:11 UTC
Permalink
Post by minforth
Pretzel coding. I wondered what that could be useful for.
Kind of wonder why myself; more focused on ways to do the given
example rather than why it's needed to be a forward case. By the
way, the macro example I gave wasn't right; see example 5 for
correct way. Example 4 seems to me the most proper way as it
uses all basic Forth elements. Example 3 seems the best but it's
not a forwarding case so it may not be applicable to what's
being sought.

Example 1
: bar dup . 1- [ here tmp! ] noop ;
: foo dup 0> if bar [ latest pfa cfa tmp@ ! ] else drop then ;
5 foo5 4 3 2 1

Example 2
: forward here tmp! 0 , ; immediate
: resolve latest pfa cfa tmp@ ! ; immediate
: bar dup . 1- forward ;
: foo dup 0> if bar resolve else drop then ;
5 foo5 4 3 2 1

Example 3
: bar dup . 1- ;
: foo dup 0> if bar myself else drop then ;
5 foo5 4 3 2 1

Example 4
defer foo
: bar dup . 1- foo ;
anon dup 0> if bar else drop then ; is foo
5 foo5 4 3 2 1

Example 5
"here tmp! 0 ," /mm: forward
"latest pfa cfa tmp@ !" /mm: resolve
: bar dup . 1- [ mm forward ] ;
: foo dup 0> if bar [ mm resolve ] else drop then ;
5 foo5 4 3 2 1
-- remove all macros
mm.clear
-- do again
5 foo5 4 3 2 1

--
me
minforth
2024-07-16 20:15:55 UTC
Permalink
With a 2-pass Forth compiler, it becomes a trivial matter.
mhx
2024-07-16 20:33:04 UTC
Permalink
Post by minforth
With a 2-pass Forth compiler, it becomes a trivial matter.
Wow, that is more like a quantum jump!

(Put it on my to-do list.)

-marcel
dxf
2024-07-17 07:47:45 UTC
Permalink
Post by sjack
Post by Anton Ertl
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
"dup . 1- myself" /mm: bar OK
: foo dup 0> if mm bar else drop then ; OK
5 foo 5 4 3 2 1 OK
Note: The macro BAR does not remain in the system image.
I would use:

defer foo

: bar ( n -- ) dup . 1- foo ;

:noname ( n -- ) dup 0> if bar else drop then ; is foo

5 foo \ prints "5 4 3 2 1 "

DEFER may not be as fast as a directly patched definition but neither
has that prevented a generation from using it.
a***@spenarnc.xs4all.nl
2024-07-17 10:17:54 UTC
Permalink
Post by Gerry Jackson
Post by sjack
Post by Anton Ertl
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
"dup . 1- myself" /mm: bar OK
: foo dup 0> if mm bar else drop then ; OK
5 foo 5 4 3 2 1 OK
Note: The macro BAR does not remain in the system image.
defer foo
: bar ( n -- ) dup . 1- foo ;
:noname ( n -- ) dup 0> if bar else drop then ; is foo
5 foo \ prints "5 4 3 2 1 "
DEFER may not be as fast as a directly patched definition but neither
has that prevented a generation from using it.
I applaud a sensible solution with standard words.

Groetjes Albert
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
Stephen Pelc
2024-07-22 13:00:22 UTC
Permalink
Post by dxf
DEFER may not be as fast as a directly patched definition but neither
has that prevented a generation from using it.
At least on x64 and CISC CPUs, calling a deferred word is just
CALL [] foo
rather than
CALL foo

The difference on x64 is one byte and a a few (hardware and cache
dependent) cycles.

IMHO Forward referencing and resolving words are likely just to be
wrappers for syntactic sugar around DEFER and IS.


Stephen
--
Stephen Pelc, ***@vfxforth.com
MicroProcessor Engineering, Ltd. - More Real, Less Time
133 Hill Lane, Southampton SO15 5AF, England
tel: +44 (0)78 0390 3612, +34 649 662 974
http://www.mpeforth.com
MPE website
http://www.vfxforth.com/downloads/VfxCommunity/
downloads
Anton Ertl
2024-07-22 14:44:41 UTC
Permalink
Post by Stephen Pelc
Post by dxf
DEFER may not be as fast as a directly patched definition but neither
has that prevented a generation from using it.
At least on x64 and CISC CPUs, calling a deferred word is just
CALL [] foo
rather than
CALL foo
The difference on x64 is one byte and a a few (hardware and cache
dependent) cycles.
It's a bit more in Gforth (see below), but indeed, the difference is
small.
Post by Stephen Pelc
IMHO Forward referencing and resolving words are likely just to be
wrappers for syntactic sugar around DEFER and IS.
As I demonstrated in <***@mips.complang.tuwien.ac.at>,
Gforth's FORWARD uses direct calls in compiled words. Maybe the
example was unclear. Here is a simpler one:

forward x1
: y1 x1 ;
: x1 ;
y1

defer x2
: y2 x2 ;
: z2 ;
' z2 is x2
z2

cr see-code y1
cr see-code y2

Here's what the SEE-CODE calls produce:

cr see-code y1
$7F411D19D810 call 1->1
$7F411D19D818 x1
0x00007f411ce3f2c0: mov 0x8(%rbx),%rax
0x00007f411ce3f2c4: sub $0x8,%r14
0x00007f411ce3f2c8: add $0x10,%rbx
0x00007f411ce3f2cc: mov %rbx,(%r14)
0x00007f411ce3f2cf: mov %rax,%rbx
0x00007f411ce3f2d2: mov (%rbx),%rax
0x00007f411ce3f2d5: jmp *%rax
$7F411D19D820 ;s 1->1
0x00007f411ce3f2d7: mov (%r14),%rbx
0x00007f411ce3f2da: add $0x8,%r14
0x00007f411ce3f2de: mov (%rbx),%rax
0x00007f411ce3f2e1: jmp *%rax
ok
cr see-code y2
$7F411D19D8B0 lit-perform 1->1
$7F411D19D8B8 x2
0x00007f411ce3f2f0: mov 0x8(%rbx),%rax
0x00007f411ce3f2f4: add $0x8,%rbx
0x00007f411ce3f2f8: mov (%rax),%rdx
0x00007f411ce3f2fb: mov -0x10(%rdx),%rax
0x00007f411ce3f2ff: jmp *%rax
$7F411D19D8C0 ;s 1->1
0x00007f411ce3f301: mov (%r14),%rbx
0x00007f411ce3f304: add $0x8,%r14
0x00007f411ce3f308: mov (%rbx),%rax
0x00007f411ce3f30b: jmp *%rax

So the LIT-PERFORM in Y2 looks shorter, but it calls
the docol of Z2:

docol: 19 discode
0x00005647dfff8a08 <gforth_engine+104>: add $0x8,%rbx
0x00005647dfff8a0c <gforth_engine+108>: sub $0x8,%r14
0x00005647dfff8a10 <gforth_engine+112>: mov %rbx,(%r14)
0x00005647dfff8a13 <gforth_engine+115>: mov %rdx,%rbx
0x00005647dfff8a16 <gforth_engine+118>: mov (%rbx),%rax
0x00005647dfff8a19 <gforth_engine+121>: jmp *%rax

before it executes the code of Z2. By contrast, the CALL in Y1
directly continues with the code of X1. So overall, the
LIT-PERFORM+DOCOL cost 5+6=11 instructions, while the CALL costs 7
instructions.

And before we forget, here we see clearly that the FORWARD is resolved
to a direct call (CALL), not the code generated for a deferred word.

- anton
--
M. Anton Ertl http://www.complang.tuwien.ac.at/anton/home.html
comp.lang.forth FAQs: http://www.complang.tuwien.ac.at/forth/faq/toc.html
New standard: https://forth-standard.org/
EuroForth 2024: https://euro.theforth.net
a***@spenarnc.xs4all.nl
2024-07-17 09:48:12 UTC
Permalink
Post by sjack
Post by Anton Ertl
forward foo
: bar ( n -- ) dup . 1- foo ;
: foo ( n -- ) dup 0> if bar else drop then ;
5 foo \ prints "5 4 3 2 1 "
"dup . 1- myself" /mm: bar OK
: foo dup 0> if mm bar else drop then ; OK
5 foo 5 4 3 2 1 OK
Note: The macro BAR does not remain in the system image.
I really doubt it is useful to avoid adding one definition
to the 1000+ that is already in gforth.
It is dangerous to use `myself in this context, obviously
`recurse doesn't cut it. So you are obliged to do a ton of
reading to avoid pitfalls.
By the way, I thought we are talking about the usefulness
slash unavoidability of forward definitions.
Post by sjack
--
me
Groetjes Albert
--
Don't praise the day before the evening. One swallow doesn't make spring.
You must not say "hey" before you have crossed the bridge. Don't sell the
hide of the bear until you shot it. Better one bird in the hand than ten in
the air. First gain is a cat purring. - the Wise from Antrim -
Loading...