bntest16.c

包含哈希,对称以及非对称的经典算法 包含经典事例

			gettime(&stop);

			subtime(stop, start);
			twoexpsec += cursec = sec(stop);
			twoexpms += curms = msec(stop);

			printf("2^<%d>:%3lu.%03u   ",
			       (int)expbits, cursec, curms);
#else
			printf("2^<%d>    ", (int)expbits);
#endif
			fflush(stdout);

#if CLOCK_AVAIL
			gettime(&start);
#endif
			lbnExpMod_16(A, A, i, B, k, C, i);
#if CLOCK_AVAIL
			gettime(&stop);

			subtime(stop, start);
			expsec += cursec = sec(stop);
			expms += curms = msec(stop);
			printf("<%d>^<%d>:%3lu.%03u   ",
			       (int)modbits, (int)expbits, cursec, curms);
#else
			printf("<%d>^<%d>    ", (int)modbits, (int)expbits);
#endif
			fflush(stdout);

#if CLOCK_AVAIL
			gettime(&start);
#endif
			lbnDoubleExpMod_16(D, A, i, B, k, D, i, E, l, C, i);
#if CLOCK_AVAIL
			gettime(&stop);

			subtime(stop, start);
			dblexpsec += cursec = sec(stop);
			dblexpms += curms = msec(stop);
			printf("<%d>^<%d>*<%d>^<%d>:%3lu.%03u\n",
			       (int)modbits, (int)expbits,
			       (int)modbits, (int)expbits2,
			       cursec, curms);
#else
			printf("<%d>^<%d>*<%d>^<%d>\n",
			       (int)modbits, (int)expbits,
			       (int)modbits, (int)expbits2);
#endif
		}
#if CLOCK_AVAIL
		putchar('\n');
		twoexpms += 1000 * (unsigned)(twoexpsec % j);
		while (twoexpms > 1000*j) {
			twoexpsec += j;
			twoexpms -= 1000*j;
		}
		printf("2^<%d> mod <%d> bits AVERAGE: %lu.%03u s\n",
		       (int)expbits, (int)modbits,
		       twoexpsec/j, twoexpms/j);

		expms += 1000 * (unsigned)(expsec % j);
		while (expms > 1000*j) {
			expsec += j;
			expms -= 1000*j;
		}
		printf("<%d>^<%d> mod <%d> bits AVERAGE: %lu.%03u s\n",
		       (int)modbits, (int)expbits, (int)modbits,
		       expsec/j, expms/j);

		dblexpms += 1000 * (unsigned)(dblexpsec % j);
		while (dblexpms > 1000*j) {
			dblexpsec += j;
			dblexpms -= 1000*j;
		}
		printf("<%d>^<%d> * <%d>^<%d> mod <%d> bits AVERAGE:"
		       " %lu.%03u s\n",
		       (int)modbits, (int)expbits, (int)modbits, (int)expbits2,
		       (int)modbits, dblexpsec/j, dblexpms/j);
#endif
		putchar('\n');
	}

        puts("Beginning 1000 interations of sanity checking.\n"
             "Any output indicates a bug.  No output is very strong\n"
             "evidence that all the important low-level bignum routines\n"
             "are working properly.\n");

	/*
	 * If you change this loop to have an iteration 0, all results
	 * are primted on that iteration.  Useful to see what's going
	 * on in case of major wierdness, but it produces a *lot* of
	 * output.
	 */
        for (j = 1; j <= 1000; j++) {
		/* Dop the tests for lots of different number sizes. */
		for (i = 1; i <= SIZE/2; i++) {
			/* Make a random number i words long */
			do {
				randnum(A,i);
			} while (lbnNorm_16(A,i) < i);

			/* Checl lbnCmp - does a == a? */
			if (lbnCmp_16(A,A,i) || !j) {
				bnput16("a = ", A, i);
				printf("(a <=> a) = %d\n", lbnCmp_16(A,A,i));
			}

			memcpy(c, a, sizeof(a));

			/* Check that the difference, after copy, is good. */
			if (lbnCmp_16(A,C,i) || !j) {
				bnput16("a = ", A, i);
				bnput16("c = ", C, i);
				printf("(a <=> c) = %d\n", lbnCmp_16(A,C,i));
			}

			/* Generate a non-zero random t */
			do {
				t = rand16();
			} while (!t);

			/*
			 * Add t to A.  Check that:
			 * - lbnCmp works in both directions, and
			 * - A + t is greater than A.  If there was a carry,
			 *   the result, less the carry, should be *less*
			 *   than A.
			 */
			carry = lbnAdd1_16(A,i,t);
			if (lbnCmp_16(A,C,i) + lbnCmp_16(C,A,i) != 0 ||
			    lbnCmp_16(A,C,i) != (carry ? -1 : 1) || !j)
			{
				bnput16("c       = ", C, i);
				printf("t = %lX\n", (unsigned long)t);
				bnput16("a = c+t = ", A, i);
				printf("carry = %lX\n", (unsigned long)carry);
				printf("(a <=> c) = %d\n", lbnCmp_16(A,C,i));
				printf("(c <=> a) = %d\n", lbnCmp_16(C,A,i));
			}

			/* Subtract t again */
			memcpy(d, a, sizeof(a));
			borrow = lbnSub1_16(A,i,t);

			if (carry != borrow || lbnCmp_16(A,C,i) || !j) {
				bnput16("a = ", C, i);
				printf("t = %lX\n", (unsigned long)t);
				lbnAdd1_16(A,i,t);
				bnput16("a += t = ", A, i);
				printf("Carry = %lX\n", (unsigned long)carry);
				lbnSub1_16(A,i,t);
				bnput16("a -= t = ", A, i);
				printf("Borrow = %lX\n", (unsigned long)borrow);
				printf("(a <=> c) = %d\n", lbnCmp_16(A,C,i));
			}

			/* Generate a random B */
			do {
				randnum(B,i);
			} while (lbnNorm_16(B,i) < i);

			carry = lbnAddN_16(A,B,i);
			memcpy(d, a, sizeof(a));
			borrow = lbnSubN_16(A,B,i);

			if (carry != borrow || lbnCmp_16(A,C,i) || !j) {
				bnput16("a = ", C, i);
				bnput16("b = ", B, i);
				bnput16("a += b = ", D, i);
				printf("Carry = %lX\n", (unsigned long)carry);
				bnput16("a -= b = ", A, i);
				printf("Borrow = %lX\n", (unsigned long)borrow);
				printf("(a <=> c) = %d\n", lbnCmp_16(A,C,i));
			}

			/* D = B * t */
			lbnMulN1_16(D, B, i, t);
			memcpy(e, d, sizeof(a));
			carry = *(BIGLITTLE(D-i-1,D+i));
			/* D = A + B * t */
			carry = lbnAddN_16(D,A,i);

			/* For comparison with the result of lbnMulAdd1_16 */
			borrow = carry + *(BIGLITTLE(D-i-1,D+i));

			*(BIGLITTLE(A-i-1,A+i)) = 0;
			carry = lbnMulAdd1_16(A, B, i, t);

			/* Did MulAdd get the same answer as mul then add? */
			if (carry != borrow || lbnCmp_16(A, D, i) || !j) {
				bnput16("a = ", C, i);
				bnput16("b = ", B, i);
				printf("t = %lX\n", (unsigned long)t);
				bnput16("e = b * t = ", E, i+1);
				bnput16("a + e = ", D, i);
				printf("carry = %lX\n", (unsigned long)borrow);
				bnput16("a + b * t = ", A, i);
				printf("carry = %lX\n", (unsigned long)carry);
			}

			memcpy(d, a, sizeof(a));
			borrow = lbnMulSub1_16(A, B, i, t);

			/* Did MulSub perform the inverse of MulAdd */
			if (carry != borrow || lbnCmp_16(A,C,i) || !j) {
				bnput16("a = ", C, i);
				bnput16("b = ", B, i);
				bnput16("a += b*t = ", D, i+1);
				printf("Carry = %lX\n", (unsigned long)carry);
				bnput16("a -= b*t = ", A, i+1);
				printf("Borrow = %lX\n", (unsigned long)borrow);
				printf("(a <=> c) = %d\n", lbnCmp_16(A,C,i));
				bnput16("b*t = ", E, i+1);
			}

			/* At this point we're done with t, so it's scratch */

			/* Does Mul work both ways symmetrically */
			lbnMul_16(C,A,i,B,i);
			lbnMul_16(D,B,i,A,i);
			if (lbnCmp_16(C,D,i+i) || !j) {
				bnput16("a = ", A, i);
				bnput16("b = ", B, i);
				bnput16("a * b = ", C, i+i);
				bnput16("b * a = ", D, i+i);
				printf("(a*b <=> b*a) = %d\n",
				       lbnCmp_16(C,D,i+i));
			}

			/* Generate an F less than A and B */
			do {
				randnum(F,i);
			} while (lbnCmp_16(F,A,i) >= 0 ||
			         lbnCmp_16(F,B,i) >= 0);

			/* Add F to D */
			lbnAdd1_16(BIGLITTLE(D-i,D+i), i, lbnAddN_16(D, F, i));
			memcpy(c, d, sizeof(d));

			/*
			 * Divide by A and check that quotient and remainder
			 * match
			 */
			t = lbnDiv_16(E,C,i+i,A,i);
			if (t || lbnCmp_16(E,B,i) || lbnCmp_16(C, F, i) || !j) {
				bnput16("a = ", A, i);
				bnput16("b = ", B, i);
				bnput16("f = ", F, i);
				bnput16("a * b + f = ", D, i+i);
				printf("qhigh = %lX\n", (unsigned long)t);
				bnput16("(a*b+f) / a = ", E, i);
				bnput16("(a*b+f) % a = ", C, i);
			}

			memcpy(c, d, sizeof(d));

			/* Divide by B and check similarly */
			t = lbnDiv_16(E,C,i+i,B,i);
			if (lbnCmp_16(E,A,i) || lbnCmp_16(C, F, i) || !j) {
				bnput16("a = ", A, i);
				bnput16("b = ", B, i);
				bnput16("f = ", F, i);
				bnput16("a * b + f = ", D, i+i);
				printf("qhigh = %lX\n", (unsigned long)t);
				bnput16("(a*b+f) / b = ", E, i);
				bnput16("(a*b+f) % b = ", C, i);
			}

			/* Checl that multiply and square do the same thing */
			lbnMul_16(C,A,i,A,i);
			lbnSquare_16(D,A,i);
			if (lbnCmp_16(C,D,i+i) || !j) {
				bnput16("a*a = ", C, i+i);
				bnput16("a^2 = ", D, i+i);
				printf("(a * a == a^2) = %d\n",
				       lbnCmp_16(C,D,i+i));
			}

			/* Compute a GCD */
			lbnCopy_16(C,A,i);
			lbnCopy_16(D,B,i);
			z = lbnGcd_16(C,i,D,i);

			/* Approximate check that the GCD came out right */
			for (l = 0;
			     l < sizeof(smallprimes)/sizeof(*smallprimes);
			     l++)
			{
				t = lbnModQ_16(z > 0 ? C : D,
				               (unsigned)(z > 0 ? z : -z),
					       smallprimes[l]);
				carry = lbnModQ_16(A, i, smallprimes[l]);
				borrow = lbnModQ_16(B, i, smallprimes[l]);
				if (!t != (!carry && !borrow)) {
					bnput16("a = ", A, i);
					printf("a mod %u = %lu\n",
					       smallprimes[l],
					       (unsigned long)carry);
					bnput16("b = ", B, i);
					printf("b mod %u = %lu\n",
					       smallprimes[l],
					       (unsigned long)borrow);
					bnput16("gcd(a,b) = ", z > 0 ? C : D,
						(unsigned)(z > 0 ? z : -z));
					printf("gcd(a,b) mod %u = %lu\n",
					       smallprimes[l],
					       (unsigned long)t);
				}
			}

			/*
			 * Do some Montgomery operations
			 * Start with A > B, and also place a copy of B
			 * into C.
			 */
			if (lbnCmp_16(A, B, i) < 0) {
				memcpy(c, a, sizeof(c));
				memcpy(a, b, sizeof(a));
				memcpy(b, c, sizeof(b));
			} else {
				memcpy(c, b, sizeof(c));
			}
			
			/* Convert to and from */
			BIGLITTLE(A[-1],A[0]) |= 1;
			lbnToMont_16(B, i, A, i);
			lbnFromMont_16(B, A, i);
			if (lbnCmp_16(B, C, i)) {
				memcpy(b, c, sizeof(c));
				bnput16("mod = ", A, i);
				bnput16("input = ", B, i);
				lbnToMont_16(B, i, A, i);
				bnput16("mont = ", B, i);
				lbnFromMont_16(B, A, i);
				bnput16("output = ", B, i);
			}
			/* E = B^5 (mod A), the simple way */
			lbnSquare_16(E, B, i);
			(void)lbnDiv_16(BIGLITTLE(E-i,E+i),E,i+i,A,i);
			lbnSquare_16(D, E, i);
			(void)lbnDiv_16(BIGLITTLE(D-i,D+i),D,i+i,A,i);
			lbnMul_16(E, D, i, B, i);
			(void)lbnDiv_16(BIGLITTLE(E-i,E+i),E,i+i,A,i);

			/* D = B^5, using ExpMod */
			BIGLITTLE(F[-1],F[0]) = 5;
			lbnExpMod_16(D, B, i, F, 1, A, i);
			if (lbnCmp_16(D, E, i)  || !j) {
				bnput16("mod = ", A, i);
				bnput16("input = ", B, i);
				bnput16("input^5 = ", E, i);
				bnput16("input^5 = ", D, i);
				printf("a>b (x <=> y) = %d\n",
				       lbnCmp_16(D,E,i));
			}
			/* TODO: Test lbnTwoExpMod, lbnDoubleExpMod */
		} /* for (i) */
		printf("\r%d ", j);
		fflush(stdout);
        } /* for (j) */
        printf("%d iterations of up to %d 16-bit words completed.\n",
	       j-1, i-1);
        return 0;
}