Here is the lhmb.c source file. Clip and save.


/*
 * lhmb.c - long hand multiplication (binary version).
 *	Calculates ridiculously large powers of 13.
 *
 *	Scott Coburn  June 1992
 */

#include 
#include "lhmb.h"

#if SHOWTIME
#include 
#include 
#include 
#endif

/*
 * See notes in lhmb.h BEFORE COMPILING...
 */

/*
 * approximate run times:
 *
 *	1:  Zeos 386SX (80386SX) (16 MHz, 2 Meg, Minix bcc compiler)
 *	2:  Intel 302-25 (80386) (25 MHz, 4 Meg, 386/ix C compiler)
 *	3:  Gateway 486/33 (80486) (33 MHz, 8Meg, 386/ix C compiler)
 *	4a: SUN SPARCstation IPC (? MHz, 24 Meg, GNU C compiler ver 2.4.5)
 *	4b: SUN SPARCstation IPC (? MHz, 24 Meg, SUN C compiler)
 *	5:  Gateway 486/66DX2 (80486) (66 MHz, 16Meg, 386/ix C compiler)
 *	6:  IBM RS/6000 model 320H (Power1, 25? MHz, 64 Meg, IBM XLC compiler)
 *	7:  DEC Alpha AXP (VMS) (Alpha, 150 MHz, VMS Alpha C compiler)
 *	8:  Gateway P5-120 (120 MHz, 32Meg, Linux/GCC compiler)
 *	9:  Generic P5-133 (133 MHz, 32Meg, Linux/GCC compiler)
 *	10:  IBM RS/6000 model 380  (Power2, 59 MHz, 64 Meg, IBM XLC compiler)
 *
 *	11a: IBM RS/6000 model 39H  (Power?, ?? MHz, 64 Meg, IBM C, 32 bit)
 *	11b: IBM RS/6000 model 39H  (Power?, ?? MHz, 64 Meg, IBM C, 64 bit)
 *
 *	12a: IBM RS/6000 model 43P  (PPC604, 133 MHz, 128 Meg, AIX/GCC, 32 bit)
 *	12b: IBM RS/6000 model 43P  (PPC604, 133 MHz, 128 Meg, AIX/GCC, 64 bit)
 *
 *	13a: Dell Precision 410 (450Mhz PIII, 128M, Linux/GCC, 32 bit)
 *	13b: Dell Precision 410 (450Mhz PIII, 128M, Linux/GCC, 64 bit)
 *
 *	14a: ? (XEON 500 MHz, 512 Meg, Linux/GCC compiler, 32 bit)
 *	14b: ? (XEON 500 MHz, 512 Meg, Linux/GCC compiler, 64 bit)
 *
 *	machine	 13^1	 13^2	 13^3	 13^4 	  13^5 	   13^6
 *	-------	-------	-------	-------	-------	-------- --------
 *	1	~0	~0	00:03.7	08:37
 *	2	~0	~0	00:00.9	02:11
 *	3	~0	~0	00:00.6	01:31
 *	4a	~0	~0	00:00.5	01:06
 *	4b	~0	~0	00:00.5	01:06	03:20:52
 *	5	~0	~0	00:00.3	00:26	01:13:55
 *	6	~0	~0	00:00.2	00:16.8	00:46:32
 *	7	~0	~0	~0	00:08.5	00:21:45
 *	8	~0	~0	~0	00:07.2	00:20:51
 *	9	~0	~0	~0	00:06.4	00:17:35
 *	10	~0	~0	~0	00:04.5	00:12:16
 *
 *	11a	~0	~0	~0	00:03.2	00:09:24
 *	11b	~0	~0	~0	00:02.2	00:05:40 16:04:30
 *
 *	12a	~0	~0	~0	00:02.1	00:07:37
 *	12b	~0	~0	~0	00:01.1	00:03:37 10:53:09
 *
 *	13a	~0	~0	~0	00:00.9	00:02:42 09:37:10
 *	13b	~0	~0	~0	~0	00:01:42 09:37:10
 *
 *	14a	~0	~0	~0	00:00.8	00:02:21
 *	14b	~0	~0	~0	~0	00:01:30 05:29:22
 *
 *	Times shown are min:sec or hr:min:sec.  All times gotten
 *	from the unix timex or time command.  Times shown as ~0 were too
 *	small to be relied upon with any accuracy.  Other times > 10 seconds
 *	or so should be accurate to a second or so.
 *
 *	On unix systems:
 *		cc -O -o lhmb lhmb.c
 *		timex nice -15 lhmb > stdout 2> stderr
 *
 *	To get a decimal version of the output:
 *
 *		lhdb < stdout > decimal
 */

#define SQLEN	100		/* size of array square[] (imperical)	*/

#define	ZEROCNT	3		/* ZEROCNT sequential zero values in	*/
				/*  a2 will stop loop...		*/
/*
 * static data
 */

#if SHOWTIME
static struct timeval tv;
static struct timezone tz;
#endif

static int square[SQLEN];		/* square or don't square	*/
static USINT a1[2*SIZE];		/* scratch array		*/
static USINT a2[2*SIZE];		/* scratch array		*/

/*
 * main
 */

main()
{
    register ULINT t1, t2;		/* scratch			*/
    register ULINT mcarry;		/* carry from multiply		*/
    register ULINT acarry;		/* carry from addition		*/
    register int i, j;			/* loop counters (what else...)	*/
    register int stop;			/* sequential zero value count	*/
    register int k;			/* index for square[] array	*/

    fprintf(stderr,"start\n");
    fflush(stderr);

    for (i = 0 ; i < 2*SIZE ; i++)
	 a1[i] = 0;			/* init scratch arrays		*/
	 a2[i] = 0;			/* init scratch arrays		*/
    a1[0] = BASE;			/* init scratch array with base	*/

    fprintf(stderr,"init done\n");
    fflush(stderr);

    i = POWER;				/* determine when to square...	*/
    k = 0;				/* init loop counter		*/
    do  {
	if ((i % 2) == 0)		/* if is even number...		*/
	    {
	    i = i / 2;			/* ...half...			*/
	    square[k++] = TRUE;		/* ...and set flag		*/
	    }
	else				/* if is odd number...		*/
	    {
	    i = i - 1;			/* ...subtract one...		*/
	    square[k++] = FALSE;	/* ...and reset flag		*/
	    }
	if (k > SQLEN)
	    {
	    fprintf(stderr,"Array square[] too small.\n");
	    exit(0);
	    }
	}
    while (i != 1);
    k--;				/* set k back to last one...	*/

    fprintf(stderr,"# of loops: %d\n",k);
    fflush(stderr);

    for ( ; k >= 0 ; k--)
	{
#if SHOWTIME
	gettimeofday(&tv, &tz);		/* get current time...		*/
	fprintf(stderr,"k %d %d   %s\n", k, square[k], ctime(&tv.tv_sec));
#else
	fprintf(stderr,"k %d %d\n", k, square[k]);
#endif
	fflush(stderr);

	if (square[k])
	    {
	    for (i = 0 ; i < SIZE ; i++)
		{
		acarry = 0;			/* init carries to 0	*/
		mcarry = 0;			/* init carries to 0	*/
		stop = 0;
		t2 = a1[i];
		for (j = 0 ; j < SIZE ; j++)
		    {
		    t1 = (a1[j] * t2) + mcarry;	/* square...		*/
		    mcarry = t1 >> SHIFT;	/* carry		*/
		    t1 = t1 & MASK;		/* remainder		*/
		    t1 = a2[i+j] + t1 + acarry;	/* add to running total	*/
		    acarry = t1 >> SHIFT;	/* carry		*/
		    t1 = t1 & MASK;		/* remainder		*/
		    a2[i+j] = t1;
		    if (t1 == 0)
			stop++;
		    else
			stop = 0;
		    if (stop == ZEROCNT) break;
		    }
		}
	    }
	else
	    {
	    acarry = 0;				/* init carries to 0	*/
	    mcarry = 0;				/* init carries to 0	*/
	    stop = 0;
	    t2 = BASE;
	    for (j = 0 ; j < SIZE ; j++)
		{
		t1 = (a1[j] * t2) + mcarry;	/* * BASE ...		*/
		mcarry = t1 >> SHIFT;		/* carry		*/
		t1 = t1 & MASK;			/* remainder		*/
		t1 = a2[j] + t1 + acarry;	/* add to running total	*/
		acarry = t1 >> SHIFT;		/* carry		*/
		t1 = t1 & MASK;			/* remainder		*/
		a2[j] = t1;
		if (t1 == 0)
		    stop++;
		else
		    stop = 0;
		if (stop == ZEROCNT) break;
		}
	    }

	stop = 0;
	for (i = 0 ; i < SIZE ; i++)
	    {
	    if (a2[i] == 0)
		stop++;
	    else
		stop = 0;
	    if (stop == ZEROCNT) break;
	    a1[i] = a2[i];
	    a2[i] = 0;
	    }
	}

    k = 0;
    for (i = SIZE-1 ; i >= 0 ; i--)
	{
	if (k == 0 && a1[i] != 0) k = 1;	/* do not print leading zeros */
	if (k == 1) printf(HFORMAT1,a1[i]);
	}
}