lip.c

应用密码学这本书的源代码

	if (ALLOCATE && !d)
	{
		zzero(bb);
		return;
	}
	if ((sa = a[0]) < 0)
	{
		sa = (-sa);
		if (d < 0)
			d = (-d);
		else
			bnegative = 1;
	}
	else if ((bnegative = (d < 0)))
		d = (-d);
	zsetlength(&b, sa + 1, "in zsmul, third argument");
        if (a == *bb) a = b;
	*bb = b;
	pb = &b[0];
	for (i = sa; i >= 0; i--)
		*pb++ = *a++;
	b[0] = sa;
	sa++;
	*pb = 0;
	zaddmul(d - 1, &b[1], &b[0]);
	while ((sa > 1) && (!(b[sa])))
		sa--;
	b[0] = sa;
	if (bnegative && (b[1] || b[0] != 1))
		b[0] = (-b[0]);
}

static void
kar_mul(
	verylong a,
	verylong b,
	verylong *c,
	long shi
	)
{
	register long al;
	register long hal;
	register long i;
	register long restoreb0 = b[0];
	register verylong pc;
	register long bbig = 1;
	verylong *a0;
	verylong *a1;
	verylong *a2;
	verylong *a3;
	verylong *a4;

	zsetlength(c, (hal = (al = a[0]) + (i = b[0])), "in kar_mul, third argument");
	if ((shi >= (5 * KAR_DEPTH)) || (al < KAR_MUL_CROV) || (i < KAR_MUL_CROV))
	{
		pc = &(*c)[1];
		for (i = hal; i > 0; i--)
			*pc++ = 0;
		pc = &(*c)[1];
		if (al <= *b)
			for (i = al; i; i--)
			{
				zaddmul(*(++a), pc++, b);
			}
		else
			for (i = *b; i; i--)
			{
				zaddmul(*(++b), pc++, a);
			}
		while ((hal > 1) && (!((*c)[hal])))
			hal--;
		(*c)[0] = hal;
		return;
	}
	a0 = &(kar_mem[shi]);
	a1 = &(kar_mem[shi + 1]);
	a2 = &(kar_mem[shi + 2]);
	a3 = &(kar_mem[shi + 3]);
	a4 = &(kar_mem[shi + 4]);
	hal = ((al + 1) >> 1);
	zsetlength(a0, al, "in kar_mul, locals\n");
	zsetlength(a1, al, "");
	zsetlength(a2, al, "");
	zsetlength(a3, al + hal, "");
	zsetlength(a4, al + 2, "");
	i = hal;
	while ((i > 1) && (!(a[i])))
		i--;
	a[0] = i;
	if (hal >= b[0])
		bbig = 0;
	else
	{
		i = hal;
		while ((i > 1) && (!(b[i])))
			i--;
		b[0] = i;
	}
	for (i = hal + 1; i <= al; i++)
		(*a1)[i - hal] = a[i];
	(*a1)[0] = al - hal;
	if (bbig)
	{
		for (i = hal + 1; i <= restoreb0; i++)
			(*a3)[i - hal] = b[i];
		(*a3)[0] = restoreb0 - hal;
	}
	kar_mul(a, b, a4, shi + 5);
	zadd(a, (*a1), a0);
	a[0] = al;
	if (bbig)
	{
		kar_mul((*a1), (*a3), c, shi + 5);
		zadd(b, (*a3), a2);
		b[0] = restoreb0;
		kar_mul((*a0), (*a2), a3, shi + 5);
	}
	else
		kar_mul((*a0), b, a3, shi + 5);
	zsubpos((*a3), (*a4), a3);
	if (bbig)
		zsubpos((*a3), *c, a3);
	zlshift((*a3), hal * NBITS, a3);
	hal <<= 1;
	if (bbig)
	{
		for (i = (*c)[0]; i; i--)
			(*c)[i + hal] = (*c)[i];
		for (i = hal; i > (*a4)[0]; i--)
			(*c)[i] = 0;
		for (; i; i--)
			(*c)[i] = (*a4)[i];
		(*c)[0] += hal;
	}
	else
	{
		for (i = (*a4)[0]; i >= 0; i--)
			(*c)[i] = (*a4)[i];
	}
	zadd(*c, (*a3), c);
}

void
zmul(
	verylong a,
	verylong b,
	verylong *c
	)
{	/* output not input */
	register long aneg;
	register long bneg;
	verylong olda;
	verylong oldb;

	if (ALLOCATE && (!a || !b))
	{
		zzero(c);
		return;
	}
	if (a == b)
	{
		zsq(a, c);
		return;
	}
	if (!kar_mem_initialized)
	{
		kar_mem_initialized = 1;
		for (aneg = (5 * KAR_DEPTH) - 1; aneg >= 0; aneg--)
			kar_mem[aneg] = (verylong) 0;
	}
	olda = a;
	oldb = b;
	if ((aneg = (*a < 0)))
		a[0] = -a[0];
	if ((bneg = (*b < 0)))
		b[0] = -b[0];
	if (*a > *b)
		kar_mul(a, b, c, (long) 0);
	else
		kar_mul(b, a, c, (long) 0);
	if (aneg != bneg && ((*c)[1] || (*c)[0] != 1))
		(*c)[0] = -(*c)[0];
	if (aneg)
		olda[0] = -olda[0];
	if (bneg)
		oldb[0] = -oldb[0];
}

static void
kar_sq(
	verylong a,
	verylong *c,
	long shi
	)
{
	register long al;
	register long hal;
	register long i;
	register verylong pc;
	verylong *a0;
	verylong *a1;
	verylong *a2;

	zsetlength(c, (i = ((al = a[0]) << 1)), "in kar_sq, second argument");
	if ((shi >= (3 * KAR_DEPTH)) || (al < KAR_SQU_CROV))
	{
		register unsigned long uncar;
		long carry = 0;
		pc = &(*c)[1];
		for (; i > 0; i--)
			*pc++ = 0;
		for (hal = 1; hal <= al; hal++)
		{
			i += 2;
			{
				zaddmulsq(al - hal, &((*c)[i]), &(a[hal]));
			}
			uncar = ((*c)[i - 1] << 1) + carry;
			(*c)[i - 1] = uncar & RADIXM;
			uncar = ((*c)[i] << 1) + (uncar >> NBITS);
			{
				zaddmulpsq(&(*c)[i - 1], a[hal], &carry);
			}
			uncar += carry;
			carry = uncar >> NBITS;
			(*c)[i] = uncar & RADIXM;
		}
		while ((i > 1) && (!((*c)[i])))
			i--;
		(*c)[0] = i;
		return;
	}
	a0 = &(kar_mem[shi]);
	a1 = &(kar_mem[shi + 1]);
	a2 = &(kar_mem[shi + 2]);
	hal = ((al + 1) >> 1);
	zsetlength(a0, al + hal + 2, "in kar_sq, locals\n");
	zsetlength(a1, al + 2, "");
	zsetlength(a2, al, "");
	i = hal;
	while ((i > 1) && (!(a[i])))
		i--;
	a[0] = i;
	for (i = hal + 1; i <= al; i++)
		(*a0)[i - hal] = a[i];
	(*a0)[0] = al - hal;
	kar_sq(a, a1, shi + 3);
	zadd(a, (*a0), a2);
	kar_sq((*a0), c, shi + 3);
	a[0] = al;
	kar_sq((*a2), a0, shi + 3);
	zsubpos((*a0), (*a1), a0);
	zsubpos((*a0), *c, a0);
	zlshift((*a0), hal * NBITS, a0);
	hal <<= 1;
	for (i = (*c)[0]; i; i--)
		(*c)[i + hal] = (*c)[i];
	for (i = hal; i > (*a1)[0]; i--)
		(*c)[i] = 0;
	for (; i; i--)
		(*c)[i] = (*a1)[i];
	(*c)[0] += hal;
	zadd(*c, (*a0), c);
}

void
zsq(
	verylong a,
	verylong *c
	)
{	/* output is not input */
	register long aneg;

	if (ALLOCATE && !a)
	{
		zzero(c);
		return;
	}
        if (!kar_mem_initialized)
        {
                kar_mem_initialized = 1;
                for (aneg = (5 * KAR_DEPTH) - 1; aneg >= 0; aneg--)
                        kar_mem[aneg] = (verylong) 0;
        }
	if ((aneg = (*a < 0)))
		a[0] = -a[0];
	kar_sq(a, c, (long) 0);
	if (aneg)
		a[0] = -a[0];
}

void
zmul_plain(
	verylong a,
	verylong b,
	verylong *cc
	)
{
	register long i;
	register verylong pc;
	register long sc;
	verylong c = *cc;
	long anegative;
	long bnegative;
	verylong olda;
	verylong oldb;

	if (ALLOCATE && (!a || !b))
	{
		zzero(cc);
		return;
	}
	olda = a;
	oldb = b;
	if ((anegative = (*a < 0)))
		a[0] = -a[0];
	if (olda == oldb)
		bnegative = anegative;
	else if ((bnegative = (*b < 0)))
		b[0] = -b[0];
	zsetlength(&c, (sc = *a + *b), "in zmul_plain, third argument");
	*cc = c;
	pc = &c[1];
	for (i = sc; i > 0; i--)
		*pc++ = 0;
	pc = &c[1];
	if (*a <= *b)
		for (i = *a; i; i--)
		{
			zaddmul(*(++a), pc++, b);
		}
	else
		for (i = *b; i; i--)
		{
			zaddmul(*(++b), pc++, a);
		}
	while ((sc > 1) && (!(c[sc])))
		sc--;
	c[0] = sc;
	if (anegative != bnegative && (c[1] || c[0] != 1))
		c[0] = -c[0];
	if (anegative)
		olda[0] = -olda[0];
	if (bnegative && oldb != olda)
		oldb[0] = -oldb[0];
}

void
zmulin(
	verylong b,
	verylong *a
	)
{
	STATIC verylong mem = 0;

	if (ALLOCATE && (!*a || !b))
	{
		zzero(a);
		return;
	}
	zcopy(*a, &mem);
	zmul(mem, b, a);
	FREESPACE(mem);
}

void
zsq_plain(
	verylong a,
	verylong *cc
	)
{
	register long i;
	register long sc;
	register verylong pc;
	register unsigned long uncar;
	long carry = 0;
	verylong c = *cc;
	long anegative;

	if (ALLOCATE && !a)
	{
		zzero(cc);
		return;
	}
	if ((anegative = (*a)) < 0)
		anegative = -anegative;
	zsetlength(&c, (sc = 2 * anegative), "in zsq_plain, second argument");
	*cc = c;
	pc = &c[1];
	for (i = sc; i > 0; i--)
		*pc++ = 0;
	for (sc = 1; sc <= anegative; sc++)
	{
		i += 2;
		{
			zaddmulsq(anegative - sc, &(c[i]), &(a[sc]));
		}
		uncar = (c[i - 1] << 1) + carry;
		c[i - 1] = uncar & RADIXM;
		uncar = (c[i] << 1) + (uncar >> NBITS);
		{
			zaddmulpsq(&c[i - 1], a[sc], &carry);
		}
		uncar += carry;
		carry = uncar >> NBITS;
		c[i] = uncar & RADIXM;
	}
	while ((i > 1) && (!(c[i])))
		i--;
	c[0] = i;
}

void
zsqin(
	verylong *a
	)
{
	STATIC verylong mem = 0;

	if (ALLOCATE && !*a)
	{
		zzero(a);
		return;
	}
	zcopy(*a, &mem);
	zsq(mem, a);
	FREESPACE(mem);
}


#ifndef ALPHA50
long 
zsdiv(
	verylong a,
	long d,
	verylong *bb
	)
{
	register long sa;
	verylong b = *bb;

#ifndef START
	if (fudge < 0)
		zstart();
#endif
	if (ALLOCATE && !a)
	{
		zzero(bb);
		return (0);
	}
	if (!d)
        {
		zhalt("division by zero in zsdiv");
		return (0);
	}
	if ((sa = a[0]) < 0)
		sa = (-sa);
	zsetlength(&b, sa, "in zsdiv, third argument");
	if (a == *bb) a = b;
	*bb = b;
	if ((d >= RADIX) || (d <= -RADIX))
	{
		STATIC verylong zd = 0;
		STATIC verylong zb = 0;
		zintoz(d, &zb);
		zdiv(a, zb, &b, &zd);
		*bb = b;
		sa = ztoint(zd);
		FREE2SPACE(zd,zb);
		return (sa);
	}
	else
	{
		register long den = d;
		register double deninv;
		register long carry = 0;
		register long i;
		long flag = (*a < 0 ? 2 : 0) | (den < 0 ? 1 : 0);

		if (den < 0)
			den = -den;
		deninv = (double)1/den;
		if (a[sa] < den && sa > 1)
			carry = a[sa--];

		for (i = sa; i; i--)
			zdiv21(carry, a[i], den, deninv, b[i]);

		while ((sa > 1) && (!(b[sa])))
			sa--;
		b[0] = sa;
		if (flag)
		{
			if (flag <= 2)
			{
				if (!carry)
					znegate(&b);
				else
				{
					zadd(b, one, &b);
					b[0] = -b[0];
					if (flag == 1)
						carry = carry - den;
					else
						carry = den - carry;
					*bb = b;
				}
			}
			else
				carry = -carry;
		}
		return (carry);
	}
}

#endif


long 
zsmod(
	verylong a,
	long b
	)
{
	STATIC verylong q = 0;
	long y;

	y= zsdiv(a, b, &q);
	FREESPACE(q);
	return y;
}


#define correct( q, x1, x2, x3, y1, y2, btopinv) { \
	register long ix1=x1,ix2=x2,ix3=x3,iy1=y1,iy2=y2; \
	long lprodlow=0,lprodhigh=0; \
\
	zaddmulp(&lprodlow, q, iy2, &lprodhigh); \
	if ((ix3 -= lprodlow) < 0) { \
		ix3+= RADIX; \
		ix2--; \
	} \
	lprodlow=0; \
	zaddmulp(&lprodlow, q, iy1, &lprodhigh); \
	if ((ix2 -= lprodlow) < 0) { \
		ix2 += RADIX; \
		--ix1; \
	} \
	if (ix1 < lprodhigh) q--; \
	else if ((ix1 -= lprodhigh)) { \
		q += (long) ((fradix * (fradix * ix1 + ix2))*btopinv ); \
	} \
	else { \
		q += (long) ((fradix * ix2 + ix3)*btopinv ); \
	} \
}


#ifndef ALPHA50
void
zdiv(
	verylong in_a,
	verylong in_b,
	verylong *qqq,
	verylong *rrr
	)

{
	register long i;
	register long qq;
	long sa;
	long sb;
	long sq;
	verylong p;
	verylong pc;
	long sign;
	STATIC verylong a = 0;
	STATIC verylong b = 0;
	STATIC verylong c = 0;
	STATIC verylong d = 0;
	double btopinv;
#ifdef ALPHA
	double btopinv2;
#endif
	double aux;
	verylong q = *qqq;
	verylong r = *rrr;

#ifndef START
	if (fudge < 0)
		zstart();
#endif

	if (ALLOCATE && !in_a)
	{
		zzero(qqq);
		zzero(rrr);
		return;
	}
	if ((!in_b) || (((sb = in_b[0]) == 1) && (!in_b[1])))
	{
		zhalt("division by zero in zdiv");
		return;
	}
	zcopy(in_a,&a);

	zcopy(in_b,&b);

	sign = (*a < 0 ? 2 : 0) | (*b < 0 ? 1 : 0);
	if (*a < 0)
		(*a) = (-(*a));
	sa = (*a);
	if (*b < 0)
		(*b) = (-(*b));
	sb = (*b);
	zsetlength(&c, (i = (sa > sb ? sa : sb) + 1), "in zdiv, locals\n");
	zsetlength(&d, i, "");
	zsetlength(&q, i, "in zdiv, third argument");
	*qqq = q;
	zsetlength(&r, sb + 1, "in zdiv, fourth argument");
	*rrr = r;
	p = &b[sb];
	if ((sb == 1) && (*p < RADIX))
		zintoz(zsdiv(a, *p, &q), &r);
	else
	{
		sq = sa - sb;	/* size of quotient */
		btopinv = (double) (*p) * fradix;
		if (sb > 1)
			btopinv += (*(--p));
		btopinv *= fradix;
#ifdef ALPHA
		btopinv2=btopinv;
		btopinv2 = fudge2/btopinv2;
#else
		if (sb > 2)
			btopinv += (*(p - 1));
#endif
		btopinv = fudge / btopinv;
		p = &a[1];
		pc = &c[0];
		for (i = sa; i > 0; i--)
			*pc++ = *p++;
		p = pc;
		sa = 1 - sb;
		for (i = (-sq); i > 0; i--)
			*p++ = 0;
		*pc = 0;
		d[1] = 0;
		p = &d[sq + 1];
		for (i = sq; i >= 0; i--)
		{

			aux = fradix * (fradix * (*pc) + (*(pc - 1))) + 1.0;
#ifndef ALPHA
			if (i > sa)
				aux += (*(pc - 2));
#endif
#ifdef ALPHA
			qq = (long) (btopinv2 * aux + 1.0);
#else
			qq = (long) (btopinv * aux + 0.5);
		/* dirty, but safe. On most machines (long)(btopinv*aux) */
		/* is correct, or one too big; on some however it becomes */
		/* too small. Could change zstart, but +0.5 and a while */
		/* instead of one single if is safer */
#endif
			if (qq > 0)
			{
#ifdef ALPHA
				if (qq > (1L<<48)) {
					correct(qq,*pc,*(pc-1),(i>sa) ? *(pc-2):0,b[sb],b[sb-1],btopinv);
					if (qq >= RADIX)
						qq = RADIX-1;
				}
#else
				if (qq >= RADIX)
					qq = RADIX-1;
#endif
				zsubmul(qq, &c[i], &b[0]);
				while (*pc < 0)
				{
					qq--;
					{
						zaddmulone(&c[i], &b[0]);
					}
				}
			}
			pc--;
			*p-- = qq;
		}
		sb--;
		while ((sb > 0) && (!(c[sb])))
			sb--;
		sb++;
		r[0] = sb;
		p = &r[1];
		pc = &c[0];
		for (i = sb; i > 0; i--)
			*p++ = *pc++;
		if (sq < 0)
		{
			q[0] = 1;
			q[1] = d[1];
		}
		else
		{
			sq++;
			while ((sq > 1) && (!(d[sq])))
				sq--;
			q[0] = sq;
			p = &q[1];
			pc = &d[1];
			for (i = sq; i > 0; i--)
				*p++ = *pc++;
		}
	}	/* non-trivial case */

	if (sign)
	{
		if (sign <= 2)
		{
			if (!(r[1]) && (r[0] == 1))
				q[0] = -(q[0]);
			else
			{
				zadd(q, one, &q);
				q[0] = -q[0];
				if (sign == 1)
					zsub(r, b, &r);
				else
					zsub(b, r, &r);
			}
		}
		else
			znegate(&r);
	}
	*qqq = q;
	*rrr = r;
	FREE2SPACE(a,b); FREE2SPACE(c,d);
}


void
zmod(
	verylong in_a,
	verylong in_b,
	verylong *rr
	)
{
	register long i;
	register long qq;
	verylong r = *rr;
	STATIC verylong a = 0;
	STATIC verylong b = 0;
	STATIC verylong c = 0;
	long sa;
	long sb;
	long sq;
	verylong p;
	verylong pc;
	long sign;