|Re: Allow me some questions?|
Message #15 Posted by Eamonn on 29 July 2003, 3:12 p.m.,
in response to message #14 by Eamonn
Sorry about the last post - used the wrong formatting code for the listing. This should be easier to read.
Here is a program for the HP-25 that calculates either C(n,r) or P(n,r) for all cases where C(n,r) and P(n,r) are less than 1e100. The program fits comfortably in the 49 steps of program memory available on the HP-25. I haven't tested these routines on an actual HP-25, but I have tested them on the simulator found on this museums HP-25 page.
The C(n,r) routine is based on the algorithm described by Victor. It makes use of the identity C(n,r) = C(n,n-r), selecting the minimum value of r and n-r for the calculations. In the loop body, the division is done first, to avoid overflows in cases where the intermediate calculations exceed the dynamic range of the HP-25. Since doing division first can result in non-integer results, the final result is rounded to the nearest integer. This routine does work for the case C(355,167), returning a value of 3.0443588e99. Incidentially, the HP-11C also works for this case, returning a value of 3.044358e99. There is a possibility that the final rounding does not produce the correct result, but all the cases I've tested it with seem to work fine so far.
The P(n,r) routine uses the identity P(n,r) = n x (n-1) x (n-2) x ... x (n-r+1). No divisions required here, so the result is always an integer. There is no limitation on the values of n and r, as long as P(n,r) is less than 1e100. For example, P(5000,27) is calculated to be 6.9446225e99.
To use the C(n,r) routine, position the program counter at step 0 if it's not there already, key in n, hit enter, key in r and hit R/S.
To use the P(n,r) routine, position the program counter at step 30 (GTO 30), key in n, hit enter, key in r and hit R/S.
01 sto 0 ; X=r, Y=n. Store r in R0
03 sto 1 ; X=n, Y=r. Store n in R1
04 x<>y ; X=r, Y=n
05 - ; X=n-r
06 rcl 0 ; X=r, Y=n-r
07 x<y ; Want to use the minimum of r and (n-r). If r < n-r ...
08 gto 11 ; then skip ahead 2 steps
09 x<>y ; X=n-r, Y=r
10 sto 0 ; store n-r in R0. This will be used as the loop counter
11 1 ; X=1
12 sto 2 ; R2 = result
13 Rv ; Let's preserve the Z register while we're at it
14 rcl 0 ; X = loop counter
15 x=0 ; have we reached the end of the loop yet?
16 gto 24 ; if so, then exit loop (pc + 8)
17 sto/ 2 ; update the result
18 rcl 1 ; recall next numerator factor
19 sto* 2 ; update the result
21 sto- 0 ; decrement loop count
22 sto- 1 ; decrement numerator factor
23 gto 14 ; loop back (pc-9)
24 rcl 2 ; place the result in X
27 + ; round result to nearest integer
28 int ; by adding 0.5 and taking the integer part
29 gto 00 ; done
30 x<>y ; P(n,r) routine. X=n, Y=r
31 sto 0 ; Store n in R0
32 x<>y ; X=r, Y=n
33 1 ; X=1, Y=r, Z=n
34 - ; X=r-1, Y=n
35 - ; X=n-r+1
36 rcl 0 ; X=n, y=n-r+1
38 - ; subtract 1 from X
39 x<y ; are we done yet?
40 gto 43 ; if so, then display the result
41 sto* 0 ; update result
42 gto 37 ; loop
43 rcl 0 ; display the result