37 lines
3.6 KiB
Rexx
37 lines
3.6 KiB
Rexx
/*REXX program computes the average loop length mapping a random field 1···N ───► 1···N */
|
|
parse arg runs tests seed . /*obtain optional arguments from the CL*/
|
|
if runs =='' | runs =="," then runs = 40 /*Not specified? Then use the default.*/
|
|
if tests =='' | tests =="," then tests= 1000000 /* " " " " " " */
|
|
if datatype(seed, 'W') then call random ,, seed /*Is integer? For RAND repeatability.*/
|
|
!.=0; !.0=1 /*used for factorial (!) memoization.*/
|
|
numeric digits 100000 /*be able to calculate 25k! if need be.*/
|
|
numeric digits max(9, length( !(runs) ) ) /*set the NUMERIC DIGITS for !(runs). */
|
|
say right( runs, 24) 'runs' /*display number of runs we're using.*/
|
|
say right( tests, 24) 'tests' /* " " " tests " " */
|
|
say right( digits(), 24) 'digits' /* " " " digits " " */
|
|
say
|
|
say " N average exact % error " /* ◄─── title, header ►────────┐ */
|
|
hdr=" ═══ ═════════ ═════════ ═════════"; pad=left('',3) /* ◄────────┘ */
|
|
say hdr
|
|
do #=1 for runs; av=fmtD( exact(#) ) /*use four digits past decimal point. */
|
|
xa=fmtD( exper(#) ) /* " " " " " " */
|
|
say right(#,9) pad xa pad av pad fmtD( abs(xa-av) * 100 / av) /*show values.*/
|
|
end /*#*/
|
|
say hdr /*display the final header (some bars).*/
|
|
exit /*stick a fork in it, we're all done. */
|
|
/*──────────────────────────────────────────────────────────────────────────────────────*/
|
|
!: procedure expose !.; parse arg z; if !.z\==0 then return !.z
|
|
!=1; do j=2 for z -1; !=!*j; !.j=!; end; /*compute factorial*/ return !
|
|
/*──────────────────────────────────────────────────────────────────────────────────────*/
|
|
exact: parse arg x; s=0; do j=1 for x; s=s + !(x) / !(x-j) / x**j; end; return s
|
|
/*──────────────────────────────────────────────────────────────────────────────────────*/
|
|
exper: parse arg n; k=0; do tests; $.=0 /*do it TESTS times.*/
|
|
do n; r=random(1, n); if $.r then leave
|
|
$.r=1; k=k + 1 /*bump the counter. */
|
|
end /*n*/
|
|
end /*tests*/
|
|
return k/tests
|
|
/*──────────────────────────────────────────────────────────────────────────────────────*/
|
|
fmtD: parse arg y,d; d=word(d 4, 1); y=format(y, , d); parse var y w '.' f
|
|
if f=0 then return w || left('', d +1); return y
|