LISTING 1 ­ 64-BIT OPERATIONS

// 64-bit Arithmetic in C++
//
// This file contains C/C++ functions that allow the manipulation of
// 64-bit integers. It represents an alternative to the C++ class
// longlong. In this file, the type Longlong is treated as an ordinary
// C typedef, and the functions perform operations on the type without
// collection into a C++ class. For the sake of efficiency, the C++
// feature of operator overloading is not used. This implementation
// should be appreciated by those who are programming in C, since it is
// more easily translated into pure C.
//
// In this formulation, structs and unions are also used as sparingly
// as possible to maintain simplicity, while avoiding unnecessary
// data copying.
//
//
#include <iostream.h>

// macro definitions

#define min(x, y) (((x)<(y))?(x):(y))
#define max(x, y) (((x)>(y))?(x):(y))
#define FALSE 0
#define TRUE 1


// Definition of basic type
typedef struct{
	unsigned long lo, hi;
}Longlong;


union Long{
	unsigned long lng;
	struct{
		unsigned short lo, hi;
	}pr;
};


// function prototypes
void assign(Longlong x, Longlong & y);	// data assignment
void assign(long x, Longlong & y);

Longlong get_longlong(void);			// I/O functions (hex)
void put_longlong(const Longlong &a);

void increment(Longlong &x);			// increment/decrement
void decrement(Longlong &x);

void complement(Longlong &x);			// complement, negate, abs, sign
void negate(Longlong &x);
void make_abs(Longlong &x);

void shiftleft(Longlong &x, int n);		// shifts
void shiftright(Longlong &x, int n);

int iszero(const Longlong &x);			// tests
int isneg(const Longlong &x);

int equal(const Longlong &x, const Longlong &y);	// comparison operators
int not_equal(const Longlong &x, const Longlong &y);
int less(const Longlong &x, const Longlong &y);
int greater(const Longlong &x, const Longlong &y);
int less_or_equal(const Longlong &x, const Longlong &y);
int greater_or_equal(const Longlong &x, const Longlong &y);

Longlong and(const Longlong &x, const Longlong &y);	// logical operators
Longlong or(const Longlong &x, const Longlong &y);
Longlong xor(const Longlong &x, const Longlong &y);

Longlong add(const Longlong &x, const Longlong &y);	// arithmetic operators
Longlong add(const Longlong &x, const long y);
Longlong subtract(const Longlong &x, const Longlong &y);
Longlong subtract(const Longlong &x, const long y);
Longlong umult(const unsigned long x, const unsigned long y);
Longlong mult(const long x, const long y);
Longlong mult(const Longlong &x, const Longlong &y);
Longlong udiv (const Longlong &u, unsigned short a);
unsigned short div3(unsigned short u[3], Long v);
Longlong udiv(const Longlong &U, unsigned long v);
Longlong udiv(const Longlong &u, Longlong &a);
Longlong udiv(const Longlong &u, Longlong &a);
Longlong div (const Longlong &u, long a);


// Function definitions

// assign longlong (y = x)
void assign(Longlong x, Longlong & y){
	y.hi = x.hi;
	y.lo = x.lo;
}


// assign signed long to longlong (y = x)
void assign(long x, Longlong & y){
	y.lo = x;
	y.hi = 0;
	if(x < 0)
		--y.hi;
}


// input longlong in hex
Longlong get_longlong(void){
	Longlong retval;
	cin >> hex >> retval.hi >> retval.lo;
	cin >> dec;
	return retval;
}


// output longlong in hex
void put_longlong(const Longlong &a){
	cout.fill('0');
	cout.width(8);
	cout << hex << a.hi;
	cout.fill('0');
	cout.width(8);
	cout << hex << a.lo << ' ';
}


// increment a longlong
void increment(Longlong &x){
	++x.lo;
	if(x.lo == 0)
		++ x.hi;
}


// decrement a longlong
void decrement(Longlong &x){
	--x.lo;
	if((long)x.lo == -1)
		-- x.hi;
}


// complement a longlong
void complement(Longlong &x){
	x.hi = ~x.hi;
	x.lo = ~x.lo;
}


// negate a longlong
void negate(Longlong &x){
	complement(x);
	increment(x);
}


// convert a longlong to its absolute value
void make_abs(Longlong &x){
	if(isneg(x))
		negate(x);
}


// shift a longlong left n bits
void shiftleft(Longlong &x, int n){
	x.hi <<= n;
	x.hi |= (x.lo >> (32-n));
	x.lo <<= n;
}


// shift a longlong left n bits
void shiftright(Longlong &x, int n){
	x.lo >>= n;
	x.lo |= (x.hi << (32-n));
	x.hi >>= n;
}


// test a longlong for zero
int iszero(const Longlong &x){
	return ((x.hi | x.lo) == 0);
}


// test a longlong for sign
int isneg(const Longlong &x){
	return ((long)x.hi < 0);
}


/* compare two longlongs
 * return true if x == y
 */
int equal(const Longlong &x, const Longlong &y){
	return ((x.hi == y.hi) && (x.lo == y.lo));
}


/* compare two longlongs
 * return true if x != y
 */
int not_equal(const Longlong &x, const Longlong &y){
	return ((x.hi != y.hi) || (x.lo != y.lo));
}


/* compare two longlongs
 * return true if x < y
 */
int less(const Longlong &x, const Longlong &y){
	if(x.hi == y.hi)
		return (x.lo < y.lo);
	else
		return (x.hi < y.hi);
}


/* compare two longlongs
 * return true if x > y
 */
int greater(const Longlong &x, const Longlong &y){
	if(x.hi == y.hi)
		return (x.lo > y.lo);
	else
		return (x.hi > y.hi);
}


/* compare two longlongs
 * return true if x <= y
 */
int less_or_equal(const Longlong &x, const Longlong &y){
	if(x.hi == y.hi)
		return (x.lo <= y.lo);
	else
		return (x.hi < y.hi);
}


/* compare two longlongs
 * return true if x >= y
 */
int greater_or_equal(const Longlong &x, const Longlong &y){
	if(x.hi == y.hi)
		return (x.lo >= y.lo);
	else
		return (x.hi > y.hi);
}


// logical and two longlongs
Longlong and(const Longlong &x, const Longlong &y){
	Longlong retval = x;
	retval.hi &= y.hi;
	retval.lo &= y.lo;
	return retval;
}


// logical and two longlongs
Longlong or(const Longlong &x, const Longlong &y){
	Longlong retval = x;
	retval.hi |= y.hi;
	retval.lo |= y.lo;
	return retval;
}


// logical exclusive or two longlongs
Longlong xor(const Longlong &x, const Longlong &y){
	Longlong retval = x;
	retval.hi ^= y.hi;
	retval.lo ^= y.lo;
	return retval;
}


// add two longlongs
Longlong add(const Longlong &x, const Longlong &y){
	Longlong retval = x;
	retval.lo += y.lo;
	retval.hi += y.hi;
	if(retval.lo < x.lo)
		++ retval.hi;
	return retval;
}


// subtract two longlongs
Longlong subtract(const Longlong &x, const Longlong &y){
	Longlong retval = x;
	retval.lo -= y.lo;
	retval.hi -= y.hi;
	if(retval.lo > x.lo)
		-- retval.hi;
	return retval;
}


// add a (signed) long to a longlong
Longlong add(const Longlong &x, const long y){
	Longlong retval = x;
	retval.lo += y;
	if(y > 0){
		if(retval.lo < x.lo)
			++ retval.hi;
	}
	else{
		if(retval.lo > x.lo)
			-- retval.hi;
	}
	return retval;
}


// subtract a (signed) long from a longlong
Longlong subtract(const Longlong &x, const long y){
	Longlong retval = x;
	retval.lo -= y;
	if(y < 0){
		if(retval.lo < x.lo)
			++ retval.hi;
	}
	else{
		if(retval.lo > x.lo)
			-- retval.hi;
	}
	return retval;
}


// Multiply two (unsigned) longs to yield a longlong
Longlong umult(const unsigned long x, const unsigned long y){
	Longlong retval;
	Long X, Y;
	Long mid;	// middle two words
	unsigned long ad, bc;

	// pointer to address individual words
	unsigned short *p = (unsigned short *)(&retval);

	 // copy is as quick as pointer assignment
	X.lng = x;
	Y.lng = y;

	// generate partial products
	retval.hi = (unsigned long)X.pr.hi * Y.pr.hi;	// ac
	ad    = (unsigned long)X.pr.hi * Y.pr.lo;
	bc    = (unsigned long)X.pr.lo * Y.pr.hi;
	retval.lo = (unsigned long)X.pr.lo * Y.pr.lo;	// bd

	mid.lng = ad + p[1];		// ad, high half of bd
	if(mid.lng < ad)      	// if carry,
		++p[3];       	// bump high word
	mid.lng += bc;       	// all of bc
	if(mid.lng < bc)      	// if carry,
		++p[3];       	// bump high word

	// results into place
	retval.hi += mid.pr.hi;		// mid.hi to hi.lo
	p[1] = mid.pr.lo;			// mid.lo to lo.hi
	return retval;
}


/* Multiply two (signed) longs to yield a longlong
 * Uses Booth's algorithm
 */
Longlong mult(const long x, const long y){
	Longlong retval = umult((unsigned long)x, (unsigned long)y);
	if(x < 0)
		retval.hi -= y;
	if(y < 0)
		retval.hi -= x;
	return retval;
}


/* multiply two longlongs, retain lower half
 * (works for signed or unsigned numbers)
 */
Longlong mult(const Longlong &x, const Longlong &y){
	Longlong retval;

	// pointers to address individual words
	unsigned short *r = (unsigned short *)(&retval);
	unsigned short *X = (unsigned short *)(&x);
	unsigned short *Y = (unsigned short *)(&y);

	Long mid;	// middle two words
	Long high;	// top two of five words
	 unsigned long prod;

	/* generate partial products,
	 * add like terms
	 */
	// words 3-4
	high.lng  = (unsigned long)X[3] * Y[0];
	high.lng += (unsigned long)X[2] * Y[1];
	high.lng += (unsigned long)X[1] * Y[2];
	high.lng += (unsigned long)X[0] * Y[3];

	// words 2-3
	retval.hi = (unsigned long)X[2] * Y[0];
	retval.hi += (unsigned long)X[1] * Y[1];
	retval.hi += (unsigned long)X[0] * Y[2];

	// words 1-2
	mid.lng  = (unsigned long)X[1] * Y[0];
	prod   = (unsigned long)X[0] * Y[1];
	mid.lng += prod;
	if(mid.lng < prod)
		++high.lng;

	// words 0-1
	retval.lo = (unsigned long)X[0] * Y[0];

	// ripple carry term to term

	prod = (unsigned long)r[1];
	mid.lng += prod;
	if(mid.lng < prod)
		++high.lng;

	retval.hi += (unsigned long)mid.pr.hi;
	r[3] += high.pr.lo;
	r[1] = mid.pr.lo;
	return retval;
}


// divide a longlong by an unsigned short
Longlong udiv (const Longlong &u, unsigned short a){
	int i;
	Longlong retval;
	Long rem;
	unsigned long temp;
	unsigned short *up = (unsigned short *)(&u.hi);
	unsigned short *qp = (unsigned short *)(&retval.hi);

	// test for divide by zero
	if(a == 0){
		retval.lo = 0xffffffffL;
		retval.hi = 0x7fffffffL;
		if((long)u.hi < 0)
			retval.hi |= 0x80000000L;
		return retval;
	}

	// do first division (result may be long)

	retval.hi = u.hi / a;
	 rem.lng = u.hi - retval.hi * a;

	 // loop for next two digits
	 for(i=0; i<2; i++){
		--up;
		--qp;
		rem.pr.hi = rem.pr.lo;
		rem.pr.lo = *up;
		temp = rem.lng / a;
		*qp = (unsigned short)temp;
		rem.lng -= *qp * a;
	}
	return retval;
}


/* divide three digits by two
 * (used by multiword division)
 * inputs: dividend u (leftmost word = 0)
 *     divisor v,
 * return value: quotient word
 * 		 u <- remainder
 */
unsigned short div3(unsigned short u[3], Long v){
	unsigned long q;
	unsigned short temp;
	Long low;
	unsigned long *p1 = (unsigned long *)(&u[1]);
	unsigned long *p0 = (unsigned long *)(&u[0]);

	/* get initial guess of q
	 * using high digit of v
	 */
	q = *p1 / v.pr.hi;
	q = min(q, 0xffff);

	// compute first remainder
	*p1 -= q * v.pr.hi;

	// then second remainder
	low.lng = q * v.pr.lo;
	temp = u[0];
	u[0] -= low.pr.lo;
	if(u[0] > temp)
		--(*p1);
	*p1 -= low.pr.hi;

	// refine as needed
	while((long)(*p1) < 0){
		--q;
		*p0 += v.lng;
		if(*p0 < v.lng)
			++u[2];
	}
	return (unsigned short)q;
}


// divide unsigned longlong by unsigned long
Longlong udiv(const Longlong &U, unsigned long v){
	 int i, j;

	// structure to hold dividend
	union{
		unsigned short u[5];
		struct{
			Longlong u0u3;
			unsigned short u4;
			//(number matches index)
		}s;
	}div;

	// structure to hold quotient
	union{
		Longlong Q;
		unsigned short q[4];
	}Q;


	// pointers to the structures
	unsigned short *up = &div.u[2];
	unsigned short *qp = &Q.q[3];
	Long *V = (Long *)(&v);

	// first, make sure v is not a short
	if((long)v <= 0xffff)
		return udiv(U, (unsigned short)v);


	// normalize the divisor
	for(i=0; (long)v > 0; i++)
		v <<= 1;

	// and the dividend
	div.s.u4 = (unsigned short)(U.hi >> (32-i));
	div.s.u0u3.hi = U.hi << i;
	div.s.u0u3.hi += (U.lo >> (32-i));
	div.s.u0u3.lo = U.lo << i;

	// start the division loop
	*qp-- = 0;
	for(j=0; j<3; j++)
		*qp-- = div3(up--, *V);
	return Q.Q;
}


/* divide a longlong by another longlong
 * (using shift-and-subtract. Slower than the
 * other functions, so use it sparingly.
 */
Longlong udiv(const Longlong &u, Longlong &a){
	Longlong V = {0, 0};
	Longlong retval = u;
	for(int i=0; i< 64; i++){
		shiftleft(V, 1);
		if(isneg(retval))
			increment(V);
		shiftleft(retval, 1);
		if(greater_or_equal(V,a)){
			increment(retval);
			V = subtract(V, a);
		}
	}
	return retval;
}


// signed division of longlong by long
Longlong div (const Longlong &u, long a){
	Longlong retval;
	Longlong U = u;	// local copy to preserve u
	int sign = FALSE;

	// check signs of arguments, force positive
	if(isneg(U)){
		negate(U);
		sign = !sign;
	}
	 if(a < 0){
		a = -a;
		sign = !sign;
		}

		// do unsigned divide
		retval = udiv(U, (unsigned long)a);

		// restore sign
		if(sign)
			negate(retval);
		return retval;
}