games@tekred.TEK.COM (09/01/88)
Submitted by: uunet.uu.net!unido!nixpbe!kebsch
Comp.sources.games: Volume 5, Issue 58
Archive-name: pi-perl
[Is this a game or not? I'm not sure. In any event, it was sent
to me, so here it is. -br]
[[The following two programs are for your entertainment. pi.pl is a perl script
an pi.c is the equivalent written in C. Now you may compute pi (3,14....)
about thousands of digits. Yes, the programs are just for fun, especially
the pi.pl (Perl version)! :-)]]
#! /bin/sh
# This is a shell archive, meaning:
# 1. Remove everything above the #! /bin/sh line.
# 2. Save the resulting text in a file.
# 3. Execute the file with /bin/sh (not csh) to create:
# pi.pl
# pi.c
# This archive created: Sat Aug 6 12:33:15 1988
export PATH; PATH=/bin:/usr/bin:$PATH
if test -f 'pi.pl'
then
echo shar: "will not over-write existing file 'pi.pl'"
else
sed 's/^X//' << \SHAR_EOF > 'pi.pl'
Xeval "exec nice -19 perl $0 $*"
X if $running_under_some_shell;
X# ---------------------------------------------------------------------------
X# pi.perl computes pi (3.14...) about 5120 Digits
X#
X# W. Kebsch, July-1988 {uunet!mcvax}!unido!nixpbe!kebsch
X
X$my_name = `basename $0`; chop($my_name);
X$version = $my_name . "-1.2";
X
X# some working parameter
X
X$smax = 5120; # max digits
X$lmax = 4; # digits per one array element
X$hmax = 10000; # one array element contains: 0..9999
X$smin = $lmax; # min digits
X$mag = 7; # magic number
X
X# subroutines
X
Xsub mul_tm # multiply the tm array with a long value
X{
X $cb = pop(@_); # elements(array)
X $x = pop(@_); # value
X
X $c = 0;
X for($i = 1; $i <= $cb; $i++)
X {
X $z = $tm[$i] * $x + $c;
X $c = int($z / $hmax);
X $tm[$i] = $z - $c * $hmax;
X }
X}
X
Xsub mul_pm # multiply the pm array with a long value
X{
X $cb = pop(@_); # elements(array)
X $x = pop(@_); # value
X
X $c = 0;
X for($i = 1; $i <= $cb; $i++)
X {
X $z = $pm[$i] * $x + $c;
X $c = int($z / $hmax);
X $pm[$i] = $z - $c * $hmax;
X }
X}
X
Xsub divide # divide the tm array by a long value
X{
X $cb = pop(@_); # elements(array)
X $x = pop(@_); # value
X
X $c = 0;
X for($i = $cb; $i >= 1; $i--)
X {
X $z = $tm[$i] + $c;
X $q = int($z / $x);
X $tm[$i] = $q;
X $c = ($z - $q * $x) * $hmax;
X }
X}
X
Xsub add # add tm array to pm array
X{
X $cb = pop(@_); # elements(array)
X
X $c = 0;
X for($i = 1; $i <= $cb; $i++)
X {
X $z = $pm[$i] + $tm[$i] + $c;
X if($z >= $hmax)
X {
X $pm[$i] = $z - $hmax;
X $c = 1;
X }
X else
X {
X $pm[$i] = $z;
X $c = 0;
X }
X }
X}
X
X$m0 = 0; $m1 = 0; $m2 = 0;
X
Xsub check_xb # reduce current no. of elements (speed up!)
X{
X $cb = pop(@_); # current no. of elements
X
X if(($pm[$cb] == $m0) && ($pm[$cb - 1] == $m1) && ($pm[$cb - 2] == $m2))
X {
X $cb--;
X }
X $m0 = $pm[$cb];
X $m1 = $pm[$cb - 1];
X $m2 = $pm[$cb - 2];
X $cb;
X}
X
Xsub display # show the result
X{
X $cb = pop(@_); # elements(array);
X
X printf("\n%3d.", $pm[$cb]);
X $j = $mag - $lmax;
X for($i = $cb - 1; $i >= $j; $i--)
X {
X printf(" %04d", $pm[$i]);
X }
X print "\n";
X}
X
Xsub the_job # let's do the job
X{
X $s = pop(@_); # no. of digits
X
X $s = int(($s + $lmax - 1) / $lmax) * $lmax;
X $b = int($s / $lmax) + $mag - $lmax;
X $xb = $b;
X $t = int($s * 5 / 3);
X
X for($i = 1; $i <= $b; $i++) # init arrays
X {
X $pm[$i] = 0;
X $tm[$i] = 0;
X }
X $pm[$b - 1] = $hmax / 2;
X $tm[$b - 1] = $hmax / 2;
X
X printf("digits:%5d, terms:%5d, elements:%5d\n", $s, $t, $b);
X for($n = 1; $n <= $t; $n++)
X {
X printf("\r\t\t\t term:%5d", $n);
X if($n < 200)
X {
X do mul_tm((4 * ($n * $n - $n) + 1), $xb);
X }
X else
X {
X do mul_tm((2 * $n - 1), $xb);
X do mul_tm((2 * $n - 1), $xb);
X }
X if($n < 100)
X {
X do divide(($n * (16 * $n + 8)), $xb);
X }
X else
X {
X do divide((8 * $n), $xb);
X do divide((2 * $n + 1), $xb);
X }
X do add($xb);
X if($xb > $mag)
X {
X $xb = do check_xb($xb);
X }
X }
X do mul_pm(6, $b);
X do display($b);
X ($user,$sys,$cuser,$csys) = times;
X printf("\n[u=%g s=%g cu=%g cs=%g]\n",$user, $sys, $cuser, $csys);
X}
X
X# main block ----------------------------------------------------------------
X
X$no_of_args = $#ARGV + 1;
Xprint("$version, ");
Xdie("usage: $my_name <no. of digits>") unless($no_of_args == 1);
X$digits = int($ARGV[0]);
Xdie("no. of digits out of range [$smin\..$smax]")
X unless(($digits >= $smin) && ($digits <= $smax));
Xdo the_job($digits);
Xexit 0;
X
X# That's all ----------------------------------------------------------------
SHAR_EOF
if test 3802 -ne "`wc -c < 'pi.pl'`"
then
echo shar: "error transmitting 'pi.pl'" '(should have been 3802 characters)'
fi
fi
if test -f 'pi.c'
then
echo shar: "will not over-write existing file 'pi.c'"
else
sed 's/^X//' << \SHAR_EOF > 'pi.c'
X/* ------------------------------------------------------------------
X * Function: Compute a long PI-constant (4..15,000 see "smax")
X *
X * Made by Waldemar Kebsch, Salzkotten, Germany - December 1985
X * UUCP: {uunet,mcvax}!unido!nixpbe!kebsch
X *
X * Usage: pi [number] (number of digits on the right hand)
X */
X
X#include <stdio.h>;
X
X#define VERSION (" V-1.10 ")
X
X#define hmax 10000l /* one element contains: 0..9999 */
X#define lmax 4 /* max digits of an element */
X#define smin lmax /* min digits on the right hand */
X#define smax 15000 /* max digits .. */
X#define ml (1+smax/lmax+7-lmax) /* max. necessary matrix length */
X
X
Xmain(argc, argv) /* Parameter: Number of digits */
Xint argc;
Xchar **argv;
X{ /* 1. Check Parameter */
X int s;
X
X if(argc < 2)
X {
X printf("Missing argument.\n");
X exit(1);
X }
X if(argc > 2)
X {
X printf("Too many arguments.\n");
X exit(1);
X }
X if(sscanf(*++argv, "%d", &s) != 1)
X {
X printf("Syntax error.\n");
X exit(1);
X }
X if(s < smin || s > smax)
X {
X printf("Sorry, out of range.\n");
X exit(1);
X }
X
X doit(s); /* 2. Let's do it */
X}
X
X
Xdoit(s) /* Work */
Xint s;
X{
X char *malloc();
X long *pm; /* primarily matrix */
X long *tm; /* secondary matrix */
X int b, xb; /* number of necessary blocks (elem.) in the matrix */
X long n, t; /* n=current term, t=max. terms */
X
X s = ((s + lmax - 1) / lmax) * lmax; /* digits */
X b = xb = s / lmax + 7 - lmax; /* blocks */
X t = s + 2l * s / 3l; /* terms */
X pm = (long*)malloc(b * sizeof(long)); /* array pm */
X tm = (long*)malloc(b * sizeof(long)); /* array tm */
X
X initial(pm, tm, b);
X
X printf("\n");
X printf VERSION;
X printf("digits: %5d, terms: %5d, ", s, t);
X printf("blocks: %5d\n", b);
X nice(19);
X
X for(n = 1l; n <= t; n++) /* time eating loop */
X {
X printf("\rterm: %5d ", n); fflush(stdout);
X if(n < 232l)
X multiply(tm, (4l * (n * n - n) + 1l), xb);
X else
X {
X multiply(tm, (2l * n - 1l), xb );
X multiply(tm, (2l * n - 1l), xb );
X }
X if(n < 115l)
X divide(tm, (n * (16l * n + 8l)), xb);
X else
X {
X divide(tm, (8l * n), xb );
X divide(tm, (2l * n + 1l), xb );
X }
X
X add(pm, tm, xb);
X xb = checkxb(pm, xb);
X }
X
X multiply(pm, 6l, b);
X display(pm, b);
X}
X
Xinitial(cm1, cm2, cb) /* cm1 = cm2 = 0.5 */
Xregister long *cm1, *cm2; /* current matrix */
Xint cb;
X{
X long *cmm;
X
X for(cmm = cm1 + cb; cm1 < cmm; )
X *++cm1 = *++cm2 = 0;
X
X *--cm1 = *--cm2 = hmax / 2l;
X}
X
Xcheckxb(cm, cb) /* Reduce current block number (speed up) */
Xregister long *cm;
Xint cb;
X{
X static long mm = 0, mn = 0, nn = 0;
X register long *cmm;
X
X cmm = cm;
X cm += cb;
X if(*cm == mm && *--cm == mn && *--cm == nn)
X cb--;
X
X cmm += cb;
X mm = *cmm;
X mn = *--cmm;
X nn = *--cmm;
X return(cb);
X}
X
Xmultiply(cm, x, cb) /* cm *= x */
Xregister long *cm; /* current matrix */
Xlong x; /* multiply with .. */
Xint cb; /* current blocks */
X{
X long *cmm;
X long c, z; /* carry */
X
X for(c = 0, cmm = cm + cb; cm < cmm; )
X {
X z = *++cm * x + c;
X c = z / hmax;
X *cm = z - c * hmax; /* the same as "*cm = z % hmax;" */
X } /* but faster! */
X}
X
Xdivide(cm, x, cb) /* cm /= x */
Xregister long *cm; /* current matrix */
Xlong x; /* divide by .. */
Xint cb; /* current blocks */
X{
X long *cmm;
X long c, z, q; /* c=carry, z, q=temporary variable */
X
X for(c = 0, cmm = cm, cm += cb; cm > cmm; )
X {
X z = *cm + c;
X *cm-- = q = z / x;
X c = (z - q * x) * hmax; /* the same as "c = (z % x) * hmax;" */
X } /* but faster! */
X}
X
Xadd(cm1, cm2, cb) /* cm1 += cm2 */
Xregister long *cm1, *cm2; /* current matrix */
Xint cb; /* current blocks */
X{
X long *cmm;
X long c; /* carry */
X
X for(c = 0, cmm = cm1 + cb; cm1 < cmm; )
X {
X if((*++cm1 += (*++cm2 + c)) >= hmax)
X {
X *cm1 -= hmax;
X c = 1;
X }
X else
X c = 0;
X }
X}
X
Xdisplay(cm, cb)
Xregister long *cm; /* current matrix */
Xint cb;
X{
X int i;
X
X printf("\n%1d. \n", *(cm += cb)); /* one digit on the left hand */
X for(i = cb-1; i >= (7-lmax); i--)
X printf(" %04d", *--cm); /* print all digits on the rigth hand */
X
X printf("\n");
X}
X
X/* That's all */
SHAR_EOF
if test 4336 -ne "`wc -c < 'pi.c'`"
then
echo shar: "error transmitting 'pi.c'" '(should have been 4336 characters)'
fi
fi
# End of shell archive
exit 0