Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

Rosetta code binomial coefficient

379 views
Skip to first unread message

krishna...@ccreweb.org

unread,
Nov 7, 2015, 1:05:04 AM11/7/15
to
The Forth entry on Rosetta Code for computing the binomial coefficient is beautiful. Compare its compactness of implementation with those in other languages. This is a very nice code example in Forth. Although it meets the requirements of the entry, single integers for 32 bit word size quickly lead to failure, e.g. "60 30 choose". The implementation is better suited for 64-bit word size.

: choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;

Good job, whoever submitted this to Rosetta Code!

Krishna

m...@iae.nl

unread,
Nov 7, 2015, 2:24:04 AM11/7/15
to
This algorithm must be quite common, in Forth or elsewhere --
I have it in Euler158.frt (2008) as

: binomial ( n k -- d )
LOCALS| k n |
1 k 1+ 1 DO n I - 1+ I */ LOOP ;

In Euler113.frt (2008) it is defined using arbitrary precision:

-- result in binom
: binomial ( n k -- binom )
LOCALS| k n |
1 binom V!
k 1+ 1 DO binom n I - 1+ VS*
binom I VS/MOD DROP
LOOP
binom ; PRIVATE

-marcel

WJ

unread,
Nov 7, 2015, 3:16:18 AM11/7/15
to
Since only stack operations are performed, it is very difficult
to discern the algorithm.

Ruby:

def choose n,k
(0...k).reduce(1){|t,i| t*(n-i)/(i+1)}
end

choose(8,3)
==>56
choose(22,9)
==>497420
choose(60,30)
==>118264581564861424
choose(99,50)
==>50445672272782096667406248628
choose(200,100)
==>90548514656103281165404177077484163874504589675413336841320

Paul Rubin

unread,
Nov 7, 2015, 3:23:46 AM11/7/15
to
krishna...@ccreweb.org writes:
> Compare its compactness of implementation with those in
> other languages.

choose n r = product [n-r+1 .. n] `div` product [1 .. r]

looks more natural to me (that's Haskell fwiw).

Albert van der Horst

unread,
Nov 7, 2015, 5:58:06 AM11/7/15
to
m...@iae.nl writes:

>On Saturday, November 7, 2015 at 7:05:04 AM UTC+1, krishna...@ccreweb.org wrote:
>> The Forth entry on Rosetta Code for computing the binomial
>> coefficient is beautiful. Compare its compactness of
>> implementation with those in other languages. This is
>> a very nice code example in Forth. Although it meets
>> the requirements of the entry, single integers for 32
>> bit word size quickly lead to failure, e.g. "60 30
>> choose". The implementation is better suited for 64-bit
>> word size.
>>
>> : choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;
>>
>> Good job, whoever submitted this to Rosetta Code!

>This algorithm must be quite common, in Forth or elsewhere --
>I have it in Euler158.frt (2008) as

>: binomial ( n k -- d )
> LOCALS| k n |
> 1 k 1+ 1 DO n I - 1+ I */ LOOP ;

I do some more preparation, and save a 1+.

8 \ For N M return "N OVER M " (N/M)
9 : CHS >R R@ - 1 R> 1+ 1 ?DO OVER I + * I / LOOP NIP ;
10 \ '(./.) ALIAS CHS

Of course I do so many Euler problems that I have it in my library.
I overlooked the possibility to use */, which extends the
range of usable numbers considerably. Damn.

>In Euler113.frt (2008) it is defined using arbitrary precision:

>-- result in binom
>: binomial ( n k -- binom )
> LOCALS| k n |
> 1 binom V!
> k 1+ 1 DO binom n I - 1+ VS*
> binom I VS/MOD DROP
> LOOP
> binom ; PRIVATE

With k single precision.

>-marcel
--
Albert van der Horst, UTRECHT,THE NETHERLANDS
Economic growth -- being exponential -- ultimately falters.
albert@spe&ar&c.xs4all.nl &=n http://home.hccnet.nl/a.w.m.van.der.horst

WJ

unread,
Nov 7, 2015, 9:07:49 AM11/7/15
to
Shorter:

def choose n,k
(1..k).reduce(1){|t,i| t*(n-i+1)/i}
end

Explanation:

1..k creates a range of integers (k is the last).

"reduce" is sometimes called "fold" in Lispy languages;
in Ruby, it's also known as "inject".

It iterates over a sequence of items, building up a result
as it goes.

In the code above, the value we are generating is called "t"
(for temporary). Its initial value is set to 1.
(If we did not specify an initial value, the first item of
the sequence would be used.)

This part
{|t,i| t*(n-i+1)/i}
is an anonymouse code-block (lambda or "quotation").
It is called for each item in the sequence, receiving two
values. The first is the value that we have built so far,
and the second is the current item of the sequence.
t*(n-i+1)/i
represents the new value that we are generating.

When "reduce" finishes, it yields the value that we have
built up.

Bernd Paysan

unread,
Nov 7, 2015, 9:11:47 AM11/7/15
to
The Forth solution gets pretty far into the binomial coefficients without
overflowing integers. Your solution requires bignums for pretty small
coefficients.

With the concept of lazy evaluation, you could do much better in terms of
computation speed, as the building rule of Pascal's triangle only requires
the + operation.

choose n r = choose n-1 r + choose n-1 r+1
choose 0 r = r==0 ? 1 : 0

Lazy evaluation assumes that you keep what already has been computed.

I.e. just compute all of Pascal's triangle instead of trying to generate
every binomial coefficient on the fly; most of the numbers stay small, and
the operation is very simple, just +.

--
Bernd Paysan
"If you want it done right, you have to do it yourself"
net2o ID: kQusJzA;7*?t=uy@X}1GWr!+0qqp_Cn176t4(dQ*
http://bernd-paysan.de/

krishna...@ccreweb.org

unread,
Nov 7, 2015, 10:49:14 AM11/7/15
to
On Saturday, November 7, 2015 at 8:11:47 AM UTC-6, Bernd Paysan wrote:
> Paul Rubin wrote:
>
> > krishna...@ccreweb.org writes:
> >> Compare its compactness of implementation with those in
> >> other languages.
> >
> > choose n r = product [n-r+1 .. n] `div` product [1 .. r]
> >
> > looks more natural to me (that's Haskell fwiw).
>
> The Forth solution gets pretty far into the binomial coefficients without
> overflowing integers. Your solution requires bignums for pretty small
> coefficients.
>

The Forth solution permits up to n = 33, for all valid k, for 32-bit cell size. The growth of the factorials, thought, don't extend this range very far even for 64-bit cell size. I haven't determined where it fails, but it is for n < 100 on 64 bits.

Also, it would be interesting to compare speed of evaluation in some of the higher level languages such as Haskell, with that of the Forth version.

> With the concept of lazy evaluation, you could do much better in terms of
> computation speed, as the building rule of Pascal's triangle only requires
> the + operation.
>
> choose n r = choose n-1 r + choose n-1 r+1
> choose 0 r = r==0 ? 1 : 0
>
> Lazy evaluation assumes that you keep what already has been computed.
>
> I.e. just compute all of Pascal's triangle instead of trying to generate
> every binomial coefficient on the fly; most of the numbers stay small, and
> the operation is very simple, just +.
>

Yes. It's not clear to me that the entire tree needs to be computed either, up to the needed row and column.

> --
> Bernd Paysan
> "If you want it done right, you have to do it yourself"
> net2o ID: kQusJzA;7*?t=uy@X}1GWr!+0qqp_Cn176t4(dQ*
> http://bernd-paysan.de/


Krishna

krishna...@ccreweb.org

unread,
Nov 7, 2015, 11:10:55 AM11/7/15
to
On Saturday, November 7, 2015 at 8:11:47 AM UTC-6, Bernd Paysan wrote:
...
> With the concept of lazy evaluation, you could do much better in terms of
> computation speed, as the building rule of Pascal's triangle only requires
> the + operation.
>
> choose n r = choose n-1 r + choose n-1 r+1
> choose 0 r = r==0 ? 1 : 0
>

What is "lazy evaluation"? Isn't the above just recursion, as is implemented in the Julia version on Rosetta code?

function binom(n,k)
n >= k || return 0 #short circuit base cases
n == 1 && return 1
k == 0 && return 1

binom(n-1,k-1) + binom (n-1,k) #recursive call
end


Krishna

Paul Rubin

unread,
Nov 7, 2015, 11:50:58 AM11/7/15
to
Bernd Paysan <bernd....@gmx.de> writes:
> The Forth solution gets pretty far into the binomial coefficients without
> overflowing integers. Your solution requires bignums for pretty small
> coefficients.

True, maybe there's a simple way to improve that. Hmm.

> choose n r = choose n-1 r + choose n-1 r+1
> choose 0 r = r==0 ? 1 : 0
> Lazy evaluation assumes that you keep what already has been computed.

Unfortunately that memoization only happens if there's a data object
that stays accessible between calls. Otherwise the above takes
exponential time like a naive implementation. There are memoization
schemes in Haskell that do exploit lazy evaluation by keeping a data
object around (typically an infinite tree that's computed lazily) that
are clever and beautiful, but using one would have required loading a
library and complicating the code somewhat.

Bernd Paysan

unread,
Nov 7, 2015, 2:20:54 PM11/7/15
to
Here's a Forth program that creates Pascal's triangle up to the maximum
value for unsigned cells (4/8 bytes per cell for desktop Forths):

cell 4 = [IF] 34 [ELSE] 67 [THEN] constant maxchoose#

: tri ( n -- n' ) dup 1+ * 2/ ;
maxchoose# 1+ tri cells buffer: pascal-triangle
1 pascal-triangle !
: tri' ( k n -- addr ) tri + cells pascal-triangle + ;
: out-of ( k n -- result )
2dup u> IF 2drop 0 ELSE tri' @ THEN ;
: out-of+ ( k n -- result )
1- over 1- over out-of -rot out-of + ;
: fill-line ( n -- )
dup 1+ 0 DO I over 2dup out-of+ -rot tri' ! LOOP drop ;
: fill-triangle ( n -- )
dup 1+ 1 DO I fill-line LOOP drop ;
maxchoose# fill-triangle

Parameter order reversed (fits more naturally to the problem), and therefor
using OUT-OF as name. Test on a 64 bit machine:

6 49 out-of . 13983816

Looks good. If you need a lot of binomial coefficients within the range of
integers, this approach should be the fastest.

krishna...@ccreweb.org

unread,
Nov 7, 2015, 8:45:12 PM11/7/15
to
Nice. I'm adding this to my toolbox. In fact, I think it would make a good module for the FSL, perhaps named pascalt.x. Some additional words I would
add are shown below.

Krishna

============

\ Validate the triangle: each row must sum to 2^row
: tri-linesum ( n -- sum )
0 over 1+ 0 do over i swap tri' @ + loop nip ;

: validate-triangle ( -- )
cell 4 = if 32 else 64 then
0 do i tri-linesum 1 i lshift
<> ABORT" Invalid Pascal's Triangle!"
loop ;

validate-triangle

\ Display the triangle up to row n.
: tri. ( n -- )
0 ?do I 1+ 0 do I J tri' @ . 2 spaces loop cr loop ;

Paul Rubin

unread,
Nov 7, 2015, 10:50:52 PM11/7/15
to
Bernd Paysan <bernd....@gmx.de> writes:
> Here's a Forth program that creates Pascal's triangle up to the maximum
> value for unsigned cells (4/8 bytes per cell for desktop Forths): ...
> 6 49 out-of . 13983816
> Looks good. If you need a lot of binomial coefficients within the range of
> integers, this approach should be the fastest.

Nice. Pascal's triangle works ok with Haskell lazy evaluation too:

triangle = rows [1] where
rows lastRow = lastRow : rows (shiftAdd lastRow)
shiftAdd row = zipWith (+) (0 : row) (row ++ [0])

This gives an infinite list where each row is made by putting 0's around
the last row, then adding that result to itself shifted by 1. So for
example the 10th row is:

*Main> triangle !! 10
[1,10,45,120,210,252,210,120,45,10,1]

The 6th element of the 49th row:

*Main> triangle !! 49 !! 6
13983816

I guess it's all bignums even for the low entries, but the 100th element
of the 200th row is:

*Main> triangle !! 200 !! 100
90548514656103281165404177077484163874504589675413336841320

These all run pretty much instantly.
Message has been deleted

WJ

unread,
Nov 8, 2015, 10:37:41 AM11/8/15
to
Ruby:

def triangle n
(1..n).reduce([1]){|row| [0,*row,0].each_cons(2).map{|a,b| a+b}}
end

triangle 10
==>[1, 10, 45, 120, 210, 252, 210, 120, 45, 10, 1]

triangle(49)[6]
==>13983816

triangle(200)[100]
==>90548514656103281165404177077484163874504589675413336841320

These all run pretty much instantly.

--
The report card by the American Society of Civil Engineers showed the national
infrastructure a single grade above failure, a step from declining to the point
where everyday things simply stop working the way people expect them to. ---
washingtonpost.com/local/trafficandcommuting/us-infrastructure-gets-d-in-annual-report/2013/03/19/c48cb010-900b-11e2-9cfd-36d6c9b5d7ad_story.html

Bernd Paysan

unread,
Nov 8, 2015, 12:51:36 PM11/8/15
to
Paul Rubin wrote:

> Bernd Paysan <bernd....@gmx.de> writes:
>> Here's a Forth program that creates Pascal's triangle up to the maximum
>> value for unsigned cells (4/8 bytes per cell for desktop Forths): ...
>> 6 49 out-of . 13983816
>> Looks good. If you need a lot of binomial coefficients within the range
>> of integers, this approach should be the fastest.
>
> Nice. Pascal's triangle works ok with Haskell lazy evaluation too:
>
> triangle = rows [1] where
> rows lastRow = lastRow : rows (shiftAdd lastRow)
> shiftAdd row = zipWith (+) (0 : row) (row ++ [0])

Yes, that looks elegant.

> This gives an infinite list where each row is made by putting 0's around
> the last row, then adding that result to itself shifted by 1. So for
> example the 10th row is:
>
> *Main> triangle !! 10
> [1,10,45,120,210,252,210,120,45,10,1]
>
> The 6th element of the 49th row:
>
> *Main> triangle !! 49 !! 6
> 13983816
>
> I guess it's all bignums even for the low entries, but the 100th element
> of the 200th row is:
>
> *Main> triangle !! 200 !! 100
> 90548514656103281165404177077484163874504589675413336841320
>
> These all run pretty much instantly.

Can you benchmark that vs. the 200 multiplications which generate really big
intermediate results and one division for these really big numbers? These
are 10k additions on pretty small bignums (200*100/2), supposed you do only
that query (once done, it's precomputed, and much cheaper).

With bignum multiplications each going with O(n²), I suppose, the 10k
additions are actually faster.

Albert van der Horst

unread,
Nov 8, 2015, 1:51:52 PM11/8/15
to
In euler problems you often need the binomials modulo a prime.
Although it is possible to divide modulo a prime, the addition
method is much faster, sometimes to the point that a solution
is unfeasible without it. Often one needs all the intermediate
results anyway.

>--
>Bernd Paysan

Groetjes Albert

Paul Rubin

unread,
Nov 8, 2015, 10:55:09 PM11/8/15
to
Bernd Paysan <bernd....@gmx.de> writes:
> Can you benchmark that vs. the 200 multiplications

I'll try to do this soon.

krishna...@ccreweb.org

unread,
Nov 9, 2015, 9:57:43 AM11/9/15
to
On Saturday, November 7, 2015 at 12:05:04 AM UTC-6, krishna...@ccreweb.org wrote:
> The Forth entry on Rosetta Code for computing the binomial coefficient is beautiful. Compare its compactness of implementation with those in other languages. This is a very nice code example in Forth. Although it meets the requirements of the entry, single integers for 32 bit word size quickly lead to failure, e.g. "60 30 choose". The implementation is better suited for 64-bit word size.
>
> : choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;
>

The FSL module permcomb.x, FSL Algorithm #59 by Gordon Charlton, uses a double length number and M*/ for a triple length intermediate product, to compute the correct double length result for n up to 66 (all k) on a 32-bit machine. Pascal's triangle formed by addition should allow at least n = 67 when double length integers are used on a 32-bit system.

Pascal's triangle suggests a strategy for computing binomial coefficients for larger n:

1) Use single length integers to store a Pascal's triangle array of binomial coefficients, up to the maximum unsigned integer representation for 32/64 bit.

2) Store an extended portion of Pascal's triangle computed with double length integers, until the maximum double length representation is reached.

3) Use floating point calculations (with addition only) to obtain binomial coefficients up to some maximum n value, N. No coefficients are stored.

4) Use logarithmic forms for calculating binomial coefficients for n > N, up to the maximum IEEE double prec. fp value.

5) Provide words which can return single integer, double length integer, and floating point binomial coefficients: two variants, ones without a check on the upper limit of n for each type of output, and ones with a bounds check on n.

Krishna


krishna...@ccreweb.org

unread,
Nov 10, 2015, 6:49:05 AM11/10/15
to
Below is the double length integer version of Bernd Paysan's code for computing Pascal's triangle. On 32-bit systems, it provides valid binomial coefficients up to n = 67, and on 64-bit systems, up to n = 130.

Now I need to flesh out the floating point code to deal with n > 130:

: fchoose ( k n -- r ) \ or ( n k -- ) ( F: -- r )
dup dmaxchoose# < IF
dout-of D>F
ELSE
\ Compute the combinations for n >= dmaxchoose#
0e
THEN ;

Krishna

============================================

kForth (32-bit system):
========================

33 66 dout-of d.
7219428434016265740 ok

33 67 dout-of d. \ kForth needs ud.
-4220223336089263246 ok


gforth (64-bit system):
=======================

65 130 dout-of ud. 95067625827960698145584333020095113100 ok

For comparison,

Maxima:
=======
binomial(66, 33);
(%o5) 7219428434016265740

binomial(67, 33);
(%o6) 14226520737620288370

binomial(130,65);
(%o3) 95067625827960698145584333020095113100


The double integer version of B. Paysan's code:

\ Double length Pascal's triangle
cell 4 = [IF] 67 [ELSE] 130 [THEN] constant dmaxchoose#
dmaxchoose# 1+ tri cells 2* buffer: dpascal-triangle
1 s>d dpascal-triangle 2!

: dtri' ( k n -- addr ) tri + cells 2* dpascal-triangle + ;

: dout-of ( k n -- dresult )
2dup u> IF 2drop 0 s>d ELSE dtri' 2@ THEN ;

: dout-of+ ( k n -- dresult )
1- over 1- over dout-of 2>r dout-of 2r> d+ ;

: dfill-line ( n -- )
dup 1+ 0 DO I over 2dup dout-of+ 2over dtri' 2! 2drop LOOP drop ;

: dfill-triangle ( n -- )
1+ 1 DO I dfill-line LOOP ;

dmaxchoose# dfill-triangle

Ron Aaron

unread,
Nov 10, 2015, 7:01:20 AM11/10/15
to
On 11/10/2015 13:49, krishna...@ccreweb.org wrote:
> =======
> binomial(66, 33);
> (%o5) 7219428434016265740
>
> binomial(67, 33);
> (%o6) 14226520737620288370
>
> binomial(130,65);
> (%o3) 95067625827960698145584333020095113100

In 8th no need for special number words, just implement the algo:

: perm \ p k -- P(p,k)
over swap \ p p k
n:- n:1+ \ p p-k+1
over n:1-
' n:* -rot loop ;

: comb \ p k -- C(p,k) == P(p,k) / fact(k)
tuck perm \ k P(p,k)
swap 2 swap
' n:/ -rot loop ;

ok> 130 65 comb . cr
95067625827960698145584333020095113100
ok> 200 65 comb . cr
355410124430498052252434369613867238863579093582797880

krishna...@ccreweb.org

unread,
Nov 10, 2015, 10:13:53 PM11/10/15
to
On Tuesday, November 10, 2015 at 6:01:20 AM UTC-6, Ron Aaron wrote:
...
> In 8th no need for special number words, just implement the algo:
>
> : perm \ p k -- P(p,k)
> over swap \ p p k
> n:- n:1+ \ p p-k+1
> over n:1-
> ' n:* -rot loop ;
>
> : comb \ p k -- C(p,k) == P(p,k) / fact(k)
> tuck perm \ k P(p,k)
> swap 2 swap
> ' n:/ -rot loop ;
...

Congrats on your 8th language. I haven't looked at it, but from your examples above, it appears to be a Forth like language that supports a higher level of abstraction, and includes arbitrary precision arithmetic. I use a variety of low and high level languages/environments for my computing needs, from assembly language to Forth to scientific computing environments such as R and Maxima. At present, I'm not looking to learn another environment. For reasons of convenience, familiarity, and ability to write code over a range of levels, from bare metal to modular programming, Forth suits quite a few of my computing needs, and I occasionally want to port functions/library routines from higher level systems to Forth. Sometimes this is just for a check on the output of higher level systems where I can't see what's going on under the hood very easily.

Cheers,
Krishna

Ron Aaron

unread,
Nov 11, 2015, 12:38:02 AM11/11/15
to


On 11/11/2015 05:13, krishna...@ccreweb.org wrote:

>
> Congrats on your 8th language. I haven't looked at it, but from your examples above, it appears to be a Forth like language that supports a higher level of abstraction, and includes arbitrary precision arithmetic. I use a variety of low and high level languages/environments for my computing needs, from assembly language to Forth to scientific computing environments such as R and Maxima. At present, I'm not looking to learn another environment. For reasons of convenience, familiarity, and ability to write code over a range of levels, from bare metal to modular programming, Forth suits quite a few of my computing needs, and I occasionally want to port functions/library routines from higher level systems to Forth. Sometimes this is just for a check on the output of higher level systems where I can't see what's going on under the hood very easily.

Hi, Krishna -

Indeed, it's exactly as you've described it. The tack we took with 8th
regarding numbers is that the end-user usually just wants to get the
work done, and doesn't really care about 'integers' vs 'floats' vs 'big
numbers' etc. So 8th takes care of the low-level details and as far as
the user is concerned there is (usually) no reason to care.

8th is well suited to applications development (and the reason it was
written to begin with); however, it's not suitable for bare-metal stuff
though there are accessor words for things like Bluetooth etc.

Best regards,
Ron

glid...@gmail.com

unread,
Nov 14, 2015, 8:26:53 AM11/14/15
to
On 11/7/15 1:05 AM, krishna...@ccreweb.org wrote:
> The Forth entry on Rosetta Code for computing the binomial
> coefficient is beautiful. Compare its compactness of implementation
> with those in other languages. This is a very nice code example in
> Forth.

> : choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;

An Oforth solution can be pretty good I think:

: choose // n k -- nCk
| i | 1 swap loop: i [ over i - 1 + * i / ] nip ;
200 100 choose .
90548514656103281165404177077484163874504589675413336841320

For those not familiar with Oforth:

1) | i | is an uninitialized local variable declaration, not Forth's i

2) " ( limit ) loop: y [ ... ]" is a looping construct
- loop: removes one stack item which becomes the loop limit
- y must be a declared local variable
- y steps from 1 through loop limit inclusive
- the code in the [ ... ] block is executed loop limit times

The same algorithm in standard Forth would be:
: choose \ n k -- nCk
1 swap 1+ 1 do over i - 1 + * i / loop nip ;

Note use of the same stack manipulation words as in Oforth
and exactly the same code inside the loop.

-Doug

krishna...@ccreweb.org

unread,
Nov 14, 2015, 9:19:42 AM11/14/15
to
On Saturday, November 14, 2015 at 7:26:53 AM UTC-6, glid...@gmail.com wrote:
> On 11/7/15 1:05 AM, krishna...@ccreweb.org wrote:
> > The Forth entry on Rosetta Code for computing the binomial
> > coefficient is beautiful. Compare its compactness of implementation
> > with those in other languages. This is a very nice code example in
> > Forth.
>
> > : choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;
>
> An Oforth solution can be pretty good I think:
>
> : choose // n k -- nCk
> | i | 1 swap loop: i [ over i - 1 + * i / ] nip ;
> 200 100 choose .
> 90548514656103281165404177077484163874504589675413336841320
>
...
> The same algorithm in standard Forth would be:
> : choose \ n k -- nCk
> 1 swap 1+ 1 do over i - 1 + * i / loop nip ;
>
...

Interesting example of Oforth. However, please note that your standard Forth version, above, does not use */ for an intermediate double length, as does the Rosetta code example, so it overflows earlier, e.g.

30 15 choose u.

gives the correct answer in the Rosetta code example, but not in your standard Forth version.

Regards,
Krishna


krishna...@ccreweb.org

unread,
Nov 14, 2015, 9:22:14 AM11/14/15
to
I mean on a 32-bit system, for the particular example I gave. The same problem would exist on a 64-bit system, but for a larger value of n.

glid...@gmail.com

unread,
Nov 14, 2015, 11:16:28 AM11/14/15
to
On Saturday, November 14, 2015 at 9:19:42 AM UTC-5, krishna...@ccreweb.org wrote:
> On Saturday, November 14, 2015 at 7:26:53 AM UTC-6, glid...@gmail.com wrote:
> > On 11/7/15 1:05 AM, krishna...@ccreweb.org wrote:
> > > The Forth entry on Rosetta Code for computing the binomial
> > > coefficient is beautiful. Compare its compactness of implementation
> > > with those in other languages. This is a very nice code example in
> > > Forth.
> >
> > > : choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;
> >
> > An Oforth solution can be pretty good I think:
> >
> > : choose // n k -- nCk
> > | i | 1 swap loop: i [ over i - 1 + * i / ] nip ;
> > 200 100 choose .
> > 90548514656103281165404177077484163874504589675413336841320
> >
> ...
> > The same algorithm in standard Forth would be:
> > : choose \ n k -- nCk
> > 1 swap 1+ 1 do over i - 1 + * i / loop nip ;
> >
> ...
>
> Interesting example of Oforth. However, please note that your standard Forth version, above, does not use */ for an intermediate double length, as does the Rosetta code example, so it overflows earlier, e.g.

Yes. I knew that. The idea was to illustrate how Forth-like Oforth is.
The Oforth solution has little, if any, issues with overflow as it will automatically use bignums as required. Nothing like */ required. Not saying Oforth is better than standard Forth necessarily.
But it can do this kind of calculation in a straightforward manner.

-Doug

krishna...@ccreweb.org

unread,
Nov 14, 2015, 2:16:46 PM11/14/15
to
On Saturday, November 14, 2015 at 10:16:28 AM UTC-6, glid...@gmail.com wrote:
> ... The idea was to illustrate how Forth-like Oforth is.
> The Oforth solution has little, if any, issues with overflow as it will automatically use bignums as required. Nothing like */ required. Not saying Oforth is better than standard Forth necessarily.
> But it can do this kind of calculation in a straightforward manner.
>

Yes. Oforth appears to be a higher level language than standard Forth, and frees the user from having to worry about the mundane details of arithmetic precision. The application space for Oforth is likely a bit different than the one for standard Forth, or for other lower level languages. For applications where efficiency is a big concern, the coder must bound the problem and use the least computationally intensive solution, or at least have full control of the computational complexity, rather than the system deciding the complexity. On the other hand, where the necessary bounds are not known in advance, it is very useful to have a system that handles such details.

Cheers,
Krishna

WJ

unread,
Nov 15, 2015, 1:21:15 AM11/15/15
to
WJ wrote:

> Shorter:
>
> def choose n,k
> (1..k).reduce(1){|t,i| t*(n-i+1)/i}
> end

Using recursion:

def choose n,k
k<1?1:(n-k+1)*choose(n,k-1)/k
end

WJ

unread,
Nov 25, 2015, 3:31:57 PM11/25/15
to
In Ocaml, ints don't automatically transition into bigints.
So we have to use the Num module.
The addition operator is "+/".

#load "nums.cma" ;;
open Num ;;

let triangle n =
let rec newrow = function [n] -> [n] | m::n::ns -> m+/n :: newrow (n::ns)
and go row = function
0 -> row
| n -> go (newrow (Int 0::row)) (n-1)
in go [Int 1] n ;;

# triangle 10;;
- : Num.num list =
[Int 1; Int 10; Int 45; Int 120; Int 210; Int 252; Int 210; Int 120;
Int 45; Int 10; Int 1]

# string_of_num (List.nth (triangle 200) 100);;
- : string = "90548514656103281165404177077484163874504589675413336841320"


For timing, compile this file:


(* compile with
ocamlopt nums.cmxa Pascal_triangle.ml
*)

open Num ;;

let triangle n =
let rec newrow = function [n] -> [n] | m::n::ns -> m+/n :: newrow (n::ns)
and go row = function
0 -> row
| n -> go (newrow (Int 0::row)) (n-1)
in go [Int 1] n ;;

print_endline ( string_of_num (List.nth (triangle 200) 100));;


let start = Sys.time() in
for i = 0 to 99 do
ignore (List.nth (triangle 200) 100)
done;
Printf.printf "%.3f seconds\n" (Sys.time() -. start);;


----------

0.343 seconds on an old winXP laptop.

Paul Rubin

unread,
Nov 25, 2015, 3:43:09 PM11/25/15
to
"WJ" <w_a_...@yahoo.com> writes:
> let start = Sys.time() in
> for i = 0 to 99 do
> ignore (List.nth (triangle 200) 100)
> done;
> Printf.printf "%.3f seconds\n" (Sys.time() -. start);;

This stuff should go in comp.lang.functional. It would be much less
annoying there.

WJ

unread,
Jun 27, 2016, 8:02:20 PM6/27/16
to
krishna...@ccreweb.org wrote:

> The Forth entry on Rosetta Code for computing the binomial coefficient is
> beautiful. Compare its
> compactness of implementation with those in other languages
> . This is a very nice code example in Forth. Although
> it meets the requirements of the entry
> , single integers for 32 bit word size quickly lead to failure, e.
> g. "60 30 choose". The implementation is better suited for 64-bit word size.
>
> : choose ( n k -- nCk ) 1 swap 0 ?do over i - i 1+ */ loop nip ;

OCaml:

#load "nums.cma";;
open Num;;

let rec choose n = function
0 -> Int 1
| k -> Int (n-k+1) */ (choose n (k-1)) // Int k ;;

string_of_num (choose 200 100);;
===>
"90548514656103281165404177077484163874504589675413336841320"

--
If confirmed, Garland would be the fourth Jewish justice on the nation's
highest court, which is comprised entirely of Jews and Catholics. The three
current Jewish members of the Supreme Court are Ruth Bader Ginsburg, Elana
Kagan, and Stephen Breyer.
www.jta.org/2016/03/16/news-opinion/united-states/obama-to-name-jewish-federal-judge-to-supreme-court
0 new messages