导读
--------------------------------------------------------------------------------------------------------------------------------------------------------------
Redis源码中有一个rand.c的源文件，很明显这是一个和（伪）随机数有关的文件。细看该文件代码只有寥寥50行，不过涉及到的算法原理却不简单，读起来虽然有些晦涩，但对于深入理解48位空间中的伪随机数算法是不可多得的范本。作者在该文件的注释中写道：这个伪随机数生成函数是从pysam源码中的drand48()派生过来的。关于pysam是什么项目，并不是重点，其实很多Unix系统中都存在drand48这个函数（SVr4，POSIX.1-2001），我们可在终端中man一下drand48。
可以看到在头文件stdlib.h中定义了一系列和随机数有关的函数，除drand48()外，还有erand48()、lrand48()、nrand48()等等。所谓的48指的是随机数种子的二进制位数。它们的主要操作使用的都是同一个随机数生成算法，但是后续的处理不同，比如drand48()、erand48()返回的是 [0.0, 1.0)之间的浮点型，而lrand48()和nrand48()返回的是[0,
2^31)之间的整数。rand.c文件中有一个重要的函数——redisLrand48()。虽然作者是受了drand48()源码的启发编写的该函数，但实际上redisLrand48()和lrand48()更像，拥有几乎一样的编程接口，并且设置相同的随机数种子，生成的伪随机数序列相同。
继续读rand.c开头的注释，可以发现作者编写这个文件的目的是用于替换Lua中实现的math.random()函数。作者认为Lua随机数生成算法的默认实现中使用了libc的rand()函数，而libc的rand()并不能保证在设置了相同种子的情况下，不同平台能生成相同的随机数序列。因此作者重写这一文件，保证跨平台性。
--------------------------------------------------------------------------------------------------------------------------------------------------------------
算法原理
线性同余方程
该伪随机数生成算法使用了线性同余方程，顾名思义：所谓“线性”，指的是自变量（方程中的x）和因变量（方程中的y）是线性关系；所谓“同余”，表示它们与同一个数（比如m）相除，余数相同。可以表示为：y = f(x)%m 或 y = f(x) mod m 【f(x)是线性函数，即x的次数为1】
在Linux中，通过lrand48的手册页，可以看到lrand48()所使用的线性同余方程的具体内容（Redis使用的是同一个方程）：
Xn+1 = (aXn + c) mod m, where n >= 0 //默认情况下 a = 0x5DEECE66D c = 0xB m = 2^48
X0、X1、X2……Xn就是伪随机数序列。它们每一个都要由上一个数字通过运算来生成。
种子
随机数生成算法中都有种子的概念，无论是普通的rand()还是lrand48()。我们前面多次提到伪随机数这一术语，之所以说“伪”，是因为当种子确定之后，此后生成到随机数序列也就确定了。换句话说只要是相同的种子，生成的随机数序列总是固定的（相同的平台下），并非完全随机。
通过上一小节中的线性同余方程，我们观察到随机数序列的生成受到X0、a、c的影响（m总是2^48，表示48位空间中的运算，所以这些函数后缀才会有48）。因而X0、a、c其本质就是种子。POSIX的lrand48()的实现中和Redis的实现中给这三个参数提供了相同默认值。但不同之处是POSIX标准提供了lcong48()函数来修改这三个参数的默认值，而Redis中仅仅能修改X0（表示随机数序列的第一个数）的值，所以在Redis中真正的只有一个——X0.
源码实现
变量
rand.c中定义了三个静态的全局变量：x、a、c 与线性同余方程中的参数相对应
static uint32_t x[3] = { X0, X1, X2 }, a[3] = { A0, A1, A2 }, c = C;
这里涉及到几个宏定义：
#define X0 0x330E #define X1 0xABCD #define X2 0x1234 #define A0 0xE66D #define A1 0xDEEC #define A2 0x5 #define C 0xB
可以看出数组a中的三个元素保存的就是线性同余方程中的常数 0x5DEECE66D的三个部分。注意的是a中每个元素是uint32_t类型，就是说4个字节，实际上我们把0x5DEECE66D拆成三个部分（0x5、0xDEEC、0xE66D），每个部分最多只占用了2个字节。另外要注意的是宏定义中的X0和前文（上一节）中的X0语义是不同的，这里的X0表示的是数组X默认初始化时的第一个元素的内容，而上一节中提到的X0表示的是第一个随机数，相当于本节中的X2*2^32
+ X1*2^16 + X0。
对外函数
rand.c对外只提供了两个函数：
- redisLrand48()：返回一个随机数
- redisSrand48()：设置随机数种子
另外还有一个静态函数next()，它是为redisLrand48()服务的，对外不可见。
redisSrand48()函数
void redisSrand48(int32_t seedval) {
SEED(X0, LOW(seedval), HIGH(seedval));
}
设置随机数的种子，用到了宏函数SEED。看一下定义：
#define SET3(x, x0, x1, x2) ((x)[0] = (x0), (x)[1] = (x1), (x)[2] = (x2)) #define SEED(x0, x1, x2) (SET3(x, x0, x1, x2), SET3(a, A0, A1, A2), c = C)
可以看成宏函数SEED进行了三个赋值操作，分别是给数组x、数组a和变量c。不过a和c赋的都是默认值，所以SEED(X0,LOW(seedval),HIGH(seedval))实际进行的操作就是用X0（0x330E）、seedval的低16位，seedval的高16位分别赋值给x[0]、x[1]、x[2]。LOW和HIGH也是两个宏函数
#define N 16 #define MASK ((1 << (N - 1)) + (1 << (N - 1)) - 1) //MASK=65535(16个二进制位1) #define LOW(x) ((unsigned)(x) & MASK) //LOW(x)获得3x的低16位（2个字节） #define HIGH(x) LOW((x) >> N) //HIGH(x)获得x的高16位（2个字节）
注意上面所指的高16位和低16位指的是32位整型数据的高低。最终每次要保留的结果是48位的（存储在数组x[ ]中），它分为高16位(x[2])、中16位(x[1])和低16位(x[0])。不要混淆。
redisLrand48()函数
int32_t redisLrand48() {
next();
return (((int32_t)x[2] << (N - 1)) + (x[1] >> 1));
}
next函数容后再禀，先看一下返回值，把N=16带进去，就是：
((int32_t)x[2] << 15) + (x[1] >> 1)
就是x[2]*2^15 + x[1]/2。这里没有什么道理可讲，随机数嘛，这只是一种方案而已（lrand48()采用的方案）。
next()函数
公式推导与化简
next()是为redisLrand48()服务的static函数，也就是说对外不可见的，但是next()函数才是伪随机数生成算法的精髓所在。再回顾一遍线性同余公式：
Xn+1 = (aXn + c) mod m, where n >= 0,m = 2^48
注意：n和n+1都是X的下标。从数学角度来看很简单啊，先乘再加，最后取模。但是计算机做起来却不尽然，因为会溢出，所以我们是用数组存储的x和a。所以要做的是用数组模拟两个数的乘法运算。先把方程中的X和a表示出来：
对线性同余方程进行一下换算：
Xn+1 = (aXn + c) mod m Xn+1 = ((aXn) mod m + c mod m) mod m Xn+1 = ((aXn) mod m + c) mod m
先计算一下a*Xn：
这个多项式共有9个部分组成。再计算一下(a*Xn) mod m。因为m = 2^48，所以上述多项式中，2的幂次大于等于48的项一定是2^48的倍数，取模之后就是0，故可去掉这些项。所以(a*Xn)mod m的结果为：
只剩下6项了。接着合并同类项，并加上c，很容易写出 a*Xn + c的结果：
上述算式的每个括号中的内容就是在一次随机数生成之后x[2]、x[1]、x[0] 里面保存的新的内容（即用来表示Xn+1）。当然了编程去计算的时候还要考虑进位的问题。先看一个宏定义：
用到的宏
#define MUL(x, y, z) { int32_t l = (long)(x) * (long)(y); (z)[0] = LOW(l); (z)[1] = HIGH(l); }
#define CARRY(x, y) ((int32_t)(x) + (long)(y) > MASK)
#define ADDEQU(x, y, z) (z = CARRY(x, (y)), x = LOW(x + (y)))
很明显宏函数MUL(x,y,z)表示的是乘法运算，x和y相乘，结果用z来存储。z[0]存储低16位，z[1]存储高16位。
CARRY(x,y)就是判断两个数相加是否会进位”，这里的“进位”指的超出2个字节（16位）的大小。因为我们无论是数组a还是数组x，其中每个元素保存的都是2个字节的有效数据。也就是说，如果x[0]里面的数字超过2个字节，就要把多出的部分“进位”到x[1]中。每个元素其实都是4个字节的大小，之所以没有用到全部的4字节来存储，是因为当4字节全部用来存储数据时，相加以后可能就会溢出，编译器会弹出警告，其结果也会变成负数。
ADDEQU(x,y,z)执行的操作是：x和y相加，其结果存储到x中。z中保存是否“进位”。
next()源码
static void next(void) {
uint32_t p[2], q[2], r[2], carry0, carry1;
MUL(a[0], x[0], p);
ADDEQU(p[0], c, carry0);
ADDEQU(p[1], carry0, carry1);
MUL(a[0], x[1], q);
ADDEQU(p[1], q[0], carry0);
MUL(a[1], x[0], r);
x[2] = LOW(carry0 + carry1 + CARRY(p[1], r[0]) + q[1] + r[1] +
a[0] * x[2] + a[1] * x[1] + a[2] * x[0]);
x[1] = LOW(p[1] + r[0]);
x[0] = LOW(p[0]);
}
我们来一行一行的分析一下。声明语句不看，从第4行开始读。
MUL(a[0], x[0], p);
ADDEQU(p[0], c, carry0);
ADDEQU(p[1], carry0, carry1);
- MUL运算：a[0]和x[0]相乘，结果保存到数组p中。
- 第一个ADDEQU：然后p[0]再和c相加，carry0标识是否有“进位”。这里完成的就是运算。
- 第二个ADDEQU之后就是把低16位加上c之后是否进位，加到p[1]，并且判断这次相加之后有没有产生新的进位（指2^16次多项式向2^32次多项式的进位）并保存到carry1。
至此一次多项式
相关的运算计数完毕，接下来是2^16次多项式的运算。
MUL(a[0], x[1], q);
ADDEQU(p[1], q[0], carry0);
MUL(a[1], x[0], r);
- 第一个MUL运算：a[0]和x[1]相乘，结果保存到数组q中。
- ADDEQU运算：将p[1]（存储了“低16位”运算之后向“中16位”进位的值）和q[0]相加，结果保存到p[1]，carry0标识“中16位”（2^16次多项式）对高16位（2^32次多项式）是否有新的进位。
- 第二个MUL运算：a[1]和x[0]相乘，结果保存到数组r中。
至此完成了
系数相关运算。
x[2] = LOW(carry0 + carry1 + CARRY(p[1], r[0]) + q[1] + r[1] +
a[0] * x[2] + a[1] * x[1] + a[2] * x[0]);
x[1] = LOW(p[1] + r[0]);
x[0] = LOW(p[0]);
倒着读吧，先看x[0]，它保存的应该是一次多项式的结果。
x[1]存储2^16次多项式的系数。p[1]里面存储着a[0]x[1]的结果以及低16位的进位。r[0]存储着a[1]x[0]的结果（不包括进位）。
x[2]最复杂。回顾一下2^32次多项式的系数：
再看一眼，那行代码，这几个两两直接相乘的部分读懂了吧，那么可以去掉他们再看其他部分：
carry0 + carry1 +CARRY(p[1],
r[0])
+ q[1]
+ r[1]
剩余的多项式就是2^16次多项式在运算过程中向上的进位了。
画了一个图给大家表示一下：
这个图的上面两行是表头，下面彩色部分是所存储的内容。比如第三列表示x[0]里面存储的是p[0]和c的和。
从第一幅图到第二幅图的转变是因为执行了 ADDEQU(p[1],
q[0],
carry0);这条语句，将q[0]加到了p[1]上。然后再去看看x[2]的那段代码，就清楚多了吧，重点是考虑进位。