--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/plain;
  charset="iso-8859-1"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline

[Resending with somewhat trimmed source...]

I enclose several revisions (without the arc4random() stirring removed, but it 
can be easily readded) demonstrating various optimizations.  The last 
version, j.c, is about 7.8x as fast as the original.

Things I note:

* Reducing the number of modmults in pmod2() by "unrolling" definitely helps.  
The level to which you do this affects both the final speed and the amount of 
memory required.

* Using 32-bit values is always faster on x86 than 16-bit values.  It's also 
faster on some other platforms.

* The result of this optimization suggests that there is no way to make 
generation of 20-bit and larger numbers efficient enough to do in the packet 
forwarding path, if a modular exponentiation is involved.  However, we might 
be able to reach a compromise where some bits are pure LCG output (or e.g. 
the output of one LCG stuff into a power-of-2 LCG) and the rest are generated 
by exponentiation.



--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="d.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="d.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[16];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t s;
	u_int32_t u;

	s = 1;
	u = expo;

	while (u) {
		if (u & 1) {
			if (s == 1)
				s = t[0];
			else
				s = (s * t[0]) % RU_N;
		}
		if (u & 2) {
			if (s == 1)
				s = t[1];
			else
				s = (s * t[1]) % RU_N;
		}
		u >>= 2, t += 2;
	}
	return (s);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n;
	int noprime = 1;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	while (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % RU_N;
	tmp = tmp >> 16;

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N
	 */

	while (noprime) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			noprime = 0;
		else
			j = (j + 1) % RU_N;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	for (n = 0; n < 16; n++)
		ru_t[n] = pmod(ru_g, 1 << n, RU_N);
	ru_counter = 0;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter + 0 >= RU_MAX)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter += 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; n; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="f.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="f.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[22];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t s;
	u_int32_t u;

	s = 1;
	u = expo;

	t--;
	while (u) {
		if (u & 3) {
			if (s == 1)
				s = t[u & 3];
			else
				s = (s * t[u & 3]) % RU_N;
		}
		u >>= 2, t += 3;
	}
	return (s);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n;
	int noprime = 1;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	while (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % RU_N;
	tmp = tmp >> 16;

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N
	 */

	while (noprime) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			noprime = 0;
		else
			j = (j + 1) % RU_N;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	ru_t[0] = ru_g;
	for (n = 0; n < 21; n += 3) {
		ru_t[n+1] = (ru_t[n+0] * ru_t[n+0]) % RU_N;
		ru_t[n+2] = (ru_t[n+0] * ru_t[n+1]) % RU_N;
		ru_t[n+3] = (ru_t[n+1] * ru_t[n+1]) % RU_N;
	}
	ru_counter = 0;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter + 0 >= RU_MAX)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter += 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; n; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="g.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="g.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[36];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t s;
	u_int32_t u;

	s = 1;
	u = expo;

	t--;
	while (u) {
		if (u & 7) {
			if (s == 1)
				s = t[u & 7];
			else
				s = (s * t[u & 7]) % RU_N;
		}
		u >>= 3, t += 7;
	}
	return (s);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n;
	int noprime = 1;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	while (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % RU_N;
	tmp = tmp >> 16;

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N
	 */

	while (noprime) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			noprime = 0;
		else
			j = (j + 1) % RU_N;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	ru_t[0] = ru_g;
	for (n = 0; n < 35; n += 7) {
		ru_t[n+1] = (ru_t[n+0] * ru_t[n+0]) % RU_N;
		ru_t[n+2] = (ru_t[n+0] * ru_t[n+1]) % RU_N;
		ru_t[n+3] = (ru_t[n+1] * ru_t[n+1]) % RU_N;
		ru_t[n+4] = (ru_t[n+1] * ru_t[n+2]) % RU_N;
		ru_t[n+5] = (ru_t[n+2] * ru_t[n+2]) % RU_N;
		ru_t[n+6] = (ru_t[n+2] * ru_t[n+3]) % RU_N;
		ru_t[n+7] = (ru_t[n+3] * ru_t[n+3]) % RU_N;
	}
	ru_counter = 0;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter + 0 >= RU_MAX)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter += 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; n; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="e.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="e.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[24];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t s;
	u_int32_t u;

	s = 1;
	u = expo;

	t--;
	while (u) {
		if (u & 3) {
			if (s == 1)
				s = t[u & 3];
			else
				s = (s * t[u & 3]) % RU_N;
		}
		u >>= 2, t += 3;
	}
	return (s);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n;
	int noprime = 1;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	while (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % RU_N;
	tmp = tmp >> 16;

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N
	 */

	while (noprime) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			noprime = 0;
		else
			j = (j + 1) % RU_N;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	for (n = 0; n < 24; n += 3) {
		ru_t[n+0] = pmod(ru_g, 1 << ((n/3) * 2), RU_N);
		ru_t[n+1] = pmod(ru_t[n+0], 2, RU_N);
		ru_t[n+2] = pmod(ru_t[n+0], 3, RU_N);
	}
	ru_counter = 0;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter + 0 >= RU_MAX)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter += 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; n; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="h.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="h.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[40];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t a, b, c;

	a = (t[ 0 + ((expo >> 0) & 7)] * t[ 8 + ((expo >> 3) & 7)]) % RU_N;
	b = (t[16 + ((expo >> 6) & 7)] * t[24 + ((expo >> 9) & 7)]) % RU_N;
	c = (a * b) % RU_N;
	c = (c * t[32 + ((expo >> 12) & 7)]) % RU_N;

	return (c);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n;
	int noprime = 1;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	while (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % RU_N;
	tmp = tmp >> 16;

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N
	 */

	while (noprime) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			noprime = 0;
		else
			j = (j + 1) % RU_N;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	u_int32_t x = ru_g, y;
	for (n = 0; n < 40; n += 8) {
		ru_t[n+0] = 1;				/*0*/
		ru_t[n+1] = x;				/*1*/
		y = ru_t[n+2] = (x * x) % RU_N;		/*2*/
		y = ru_t[n+3] = (x * y) % RU_N;		/*3*/
		y = ru_t[n+4] = (x * y) % RU_N;		/*4*/
		y = ru_t[n+5] = (x * y) % RU_N;		/*5*/
		y = ru_t[n+6] = (x * y) % RU_N;		/*6*/
		y = ru_t[n+7] = (x * y) % RU_N;		/*7*/
		x = (x * y) % RU_N;			/*8*/
	}
	ru_counter = 0;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter + 0 >= RU_MAX)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter += 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; n; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="i.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="i.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[64];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t a, b, c;

	a = (t[ 0 + ((expo >>  0) & 15)] * t[16 + ((expo >>  4) & 15)]) % RU_N;
	b = (t[32 + ((expo >>  8) & 15)] * t[48 + ((expo >> 12) & 15)]) % RU_N;
	c = (a * b) % RU_N;

	return (c);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n, m;
	int noprime = 1;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	while (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % RU_N;
	tmp = tmp >> 16;

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N
	 */

	while (noprime) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			noprime = 0;
		else
			j = (j + 1) % RU_N;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	u_int32_t x = ru_g, y;
	for (n = 0; n < 64; n += 16) {
		ru_t[n+0] = 1;				/*0*/
		ru_t[n+1] = x;				/*1*/
		y = ru_t[n+2] = (x * x) % RU_N;		/*2*/
		for (m = 3; m < 16; m++)
			y = ru_t[n+m] = (x * y) % RU_N;	/*m*/
		x = (x * y) % RU_N;			/*16*/
	}
	ru_counter = 0;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter + 0 >= RU_MAX)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter += 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; n; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr
Content-Type: text/x-csrc;
  charset="iso-8859-1";
  name="j.c"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
	filename="j.c"

#include <sys/cdefs.h>
__KERNEL_RCSID(0, "$NetBSD: ip_id.c,v 1.4 2003/11/25 18:13:55 itojun Exp $");

#include <sys/types.h>
#include <sys/param.h>

#include <net/if.h>
#include <netinet/in.h>
#include <netinet/ip_var.h>

#define RU_OUT  180		/* Time after wich will be reseeded */
#define RU_MAX	30000		/* Uniq cycle, avoid blackjack prediction */
#define RU_GEN	2		/* Starting generator */
#define RU_N	32749		/* RU_N-1 = 2*2*3*2729 */
#define RU_AGEN	7		/* determine ru_a as RU_AGEN^(2*rand) */
#define RU_M	31104		/* RU_M = 2^7*3^5 - don't change */

#define PFAC_N 3
const static u_int32_t pfacts[PFAC_N] = {
	2,
	3,
	2729
};

static u_int32_t ru_x;
static u_int32_t ru_seed = 0;
static u_int32_t ru_a, ru_b;
static u_int32_t ru_g, ru_t[96];
static u_int32_t ru_counter = 0;
#if 0
static long ru_reseed;
#endif
static u_int32_t tmp;		/* Storage for unused random */

static u_int32_t pmod(u_int32_t, u_int32_t, u_int32_t);
static void ip_initid(void);

/*
 * Do a fast modular exponation, returned value will be in the range
 * of 0 - (mod-1)
 */

static u_int32_t
pmod(u_int32_t gen, u_int32_t expo, u_int32_t mod)
{
	u_int32_t s, t;
	u_int32_t u;

	s = 1;
	t = gen;
	u = expo;

	while (u) {
		if (u & 1)
			s = (s * t) % mod;
		u >>= 1;
		t = (t * t) % mod;
	}
	return (s);
}

static inline u_int32_t
pmod2(u_int32_t t[], u_int32_t expo)
{
	u_int32_t a;

	a = (t[ 0 + ((expo >>  0) & 31)] * t[32 + ((expo >>  5) & 31)]) % RU_N;
	a = (t[64 + ((expo >> 10) & 31)] * a) % RU_N;

	return (a);
}

/*
 * Initalizes the seed and chooses a suitable generator. Also toggles
 * the msb flag. The msb flag is used to generate two distinct
 * cycles of random numbers and thus avoiding reuse of ids.
 *
 * This function is called from id_randomid() when needed, an
 * application does not have to worry about it.
 */
static void
ip_initid(void)
{
	u_int32_t j, i, t;
	int n, m;

	ru_x = ((tmp = arc4random()) & 0xFFFF) % RU_M;

	/* 15 bits of random seed, plus alternating MSB */
	ru_seed ^= (tmp >> 16) | 0x8000;

	/* Determine the LCG we use */
	ru_b = ((tmp = arc4random()) & 0xfffe) | 1;
	ru_a = pmod(RU_AGEN, (tmp >> 16) & 0xfffe, RU_M);
	if (ru_b % 3 == 0)
		ru_b += 2;

	j = (tmp = arc4random()) % (RU_N - 1);

	/*
	 * Do a fast gcd(j,RU_N-1), so we can find a j with
	 * gcd(j, RU_N-1) == 1, giving a new generator for
	 * RU_GEN^j mod RU_N.  gcd(RU_N-2, RU_N-1) == 1 always,
	 * so we don't need to check for wraparound.
	 */

	for (;;) {
		for (i = 0; i < PFAC_N; i++)
			if (j % pfacts[i] == 0)
				break;

		if (i >= PFAC_N)
			break;
		j++;
	}

	ru_g = pmod(RU_GEN, j, RU_N);
	u_int32_t x = ru_g, y;
	for (n = 0; n < 96; n += 32) {
		ru_t[n+0] = 1;				/*0*/
		ru_t[n+1] = x;				/*1*/
		y = ru_t[n+2] = (x * x) % RU_N;		/*2*/
		for (m = 3; m < 32; m++)
			y = ru_t[n+m] = (x * y) % RU_N;	/*m*/
		x = (x * y) % RU_N;			/*32*/
	}
	ru_counter = RU_MAX;

#if 0
	ru_reseed = time.tv_sec + RU_OUT;
#endif
}

u_int16_t
ip_randomid(void)
{

	if (ru_counter <= 0)
		ip_initid();
#if 0
	else if (time.tv_sec > ru_reseed)
		ip_initid();
#endif

	/* Linear Congruential Generator */
	ru_x = (ru_a * ru_x + ru_b) % RU_M;

	ru_counter -= 1;

	return (ru_seed ^ pmod2(ru_t, ru_x));
}

int
main()
{
	static int foo[65536];
	int n, min, diff;
	u_int16_t r;

	for (n = 0; n < 65536; n++)
		foo[n] = -1;

	ip_initid();
	min = 65537;
	for (n = 100000000; ; n--) {
		r = ip_randomid();
		if (foo[r] == -1)
			foo[r] = n;
		else {
			diff = foo[r] - n;
			foo[r] = n;
			if (diff < min) {
				printf("new min: %d\n", diff);
				min = diff;
			}
		}
	}
}

--Boundary-00=_Civx/qQxt/bIJlr--