1/*
2 * Copyright (c) 1993 Martin Birgmeier
3 * All rights reserved.
4 *
5 * You may redistribute unmodified or modified versions of this source
6 * code provided that the above copyright notice and this and the
7 * following conditions are retained.
8 *
9 * This software is provided ``as is'', and comes with no warranties
10 * of any kind. I shall in no event be liable for anything that happens
11 * to anyone/anything when using this software.
12 *
13 * $FreeBSD: src/lib/libc/gen/rand48.h,v 1.2 2002/02/01 01:32:19 obrien Exp $
14 */
15
16#ifndef _RAND48_H_
17#define _RAND48_H_
18
19#include <math.h>
20#include <stdlib.h>
21
22#define	RAND48_SEED_0	(0x330e)
23#define	RAND48_SEED_1	(0xabcd)
24#define	RAND48_SEED_2	(0x1234)
25#define	RAND48_MULT_0	(0xe66d)
26#define	RAND48_MULT_1	(0xdeec)
27#define	RAND48_MULT_2	(0x0005)
28#define	RAND48_ADD	(0x000b)
29
30typedef unsigned long long uint48;
31
32extern uint48 _rand48_seed;
33extern uint48 _rand48_mult;
34extern uint48 _rand48_add;
35
36#define TOUINT48(x,y,z) \
37	((uint48)(x) + (((uint48)(y)) << 16) + (((uint48)(z)) << 32))
38
39#define RAND48_SEED	TOUINT48(RAND48_SEED_0, RAND48_SEED_1, RAND48_SEED_2)
40#define RAND48_MULT	TOUINT48(RAND48_MULT_0, RAND48_MULT_1, RAND48_MULT_2)
41
42#define LOADRAND48(l,x) \
43	(l) = TOUINT48((x)[0], (x)[1], (x)[2])
44
45#define STORERAND48(l,x) \
46	(x)[0] = (unsigned short)(l); \
47	(x)[1] = (unsigned short)((l) >> 16); \
48	(x)[2] = (unsigned short)((l) >> 32)
49
50#define _DORAND48(l) \
51	(l) = (l) * _rand48_mult + _rand48_add
52
53#define DORAND48(l,x) \
54	LOADRAND48(l, x); \
55	_DORAND48(l); \
56	STORERAND48(l, x)
57
58#include "fpmath.h"
59
60/*
61 * Optimization for speed: avoid int-to-double conversion.  Assume doubles
62 * are IEEE-754 and insert the bits directly.  To normalize, the (1 << 52)
63 * is the hidden bit, which the first set bit is shifted to.
64 */
65#define ERAND48_BEGIN \
66	union { \
67		union IEEEd2bits ieee; \
68		unsigned long long l; \
69	} u; \
70	int s
71
72#define ERAND48_END(x) \
73	u.l = ((x) & 0xffffffffffffULL); \
74	if (u.l == 0) \
75		return 0.0; \
76	u.l <<= 5; \
77	for(s = 0; !(u.l & (1LL << 52)); s++, u.l <<= 1) {} \
78	u.ieee.bits.exp = 1022 - s; \
79	return u.ieee.d
80
81#endif /* _RAND48_H_ */
82