18、随机数生成器与随机序列探索

随机数生成器与随机序列探索

1. 密码学安全伪随机数生成器

1.1 Fortuna 伪随机数生成器

Fortuna 是一种密码学安全的伪随机数生成器。在命令行中指定 n 的值(第 16 行),新初始化的生成器会被播种(第 19 行)并准备就绪。第 22 行的循环会输出 n N 字节的数据,每组数据由第 23 行的调用生成。需要注意的是,Fortuna 在一次输出较大数据块时运行速度更快,第 24 行将生成的数据块保存到磁盘。

n = 2000 时,会生成一个 20 亿字节的输出文件,这个文件大小足以供 dieharder 使用。使用 ent 工具对该文件进行测试,结果如下:
- 熵:每字节 8.000000 比特,这是熵的最大值,表明文件无法被压缩。
- 卡方分布:对于 20 亿个样本,卡方分布值为 289.36,随机情况下超过该值的概率为 10.00%。
- 算术平均值:数据字节的算术平均值为 127.4997(随机情况下为 127.5)。
- 蒙特卡罗法估计 π 值:为 3.141585375,误差为 0.00%。
- 序列相关系数:为 0.000008(完全不相关为 0.0)。

使用 dieharder 对该文件进行测试,并使用评分指标(见 evaluate_dieharder_results.py )进行评分,结果为:在 114 次测试中得分为 226(112 次通过,2 次弱通过,0 次失败),这几乎是完美的结果。

Fortuna 是否具有确定性取决于添加的熵是否具有确定性。如果添加的熵是确定性的,那么 Fortuna 就是确定性的,可用于模拟,但这并非其设计初衷。为确保安全性,添加的熵不能是确定性的。如果实现经过审查和信任,Fortuna 是生成安全伪随机数的不错选择。

1.2 ChaCha20 伪随机数生成器

ChaCha20 是另一种新的密码学安全伪随机数生成器,由 LibTomCrypt 库支持。与 Fortuna 不同,ChaCha20 不维护熵池,但可以在初始播种后重新播种。因此,ChaCha20 可能更适合短期任务,不过较新的 Linux 内核将其用于 /dev/urandom 设备。

ChaCha20 与它的“表亲” Salsa20 类似,使用由基本异或、模运算和旋转操作组合而成的复杂混合方案,更新操作是一系列紧密的微处理器指令,增加了安全性。它维护一个 512 位的状态和 256 位的输入密钥,“ChaCha20” 中的 “20” 表示每次调用执行基本操作的轮数。

LibTomCrypt 库为 ChaCha20 提供的接口与 Fortuna 类似,具体如下:
- chacha20_prng_start :初始化 prng_state 结构。
- chacha20_prng_add_entropy :重复调用以添加熵数据,调用时应传入至少 512 位的数据来设置状态。
- chacha20_prng_ready :添加熵后调用,使生成器准备好使用。
- chacha20_prng_read :从生成器读取数据块。
- chacha20_prng_done :使用完毕后释放内存。
- chacha20_prng_export :将当前状态导出到缓冲区。
- chacha20_prng_import :从缓冲区导入当前状态。

以下是一个 ChaCha20 的测试文件示例( test_chacha20.c ):

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include "tomcrypt.h"

#define N 1000000
void main(int argc, char *argv[]) {
    prng_state st;
    uint32_t n,i;
    uint8_t *out;
    FILE *g;
    uint8_t ent[64] = {
        12,53,75,63,4,57,8,78,89,67,8,45,62,3,41,2,
        4,56,56,67,97,88,9,78,5,67,45,23,12,43,4,76,
        74,3,8,5,34,23,134,56,6,209,45,233,2,42,33,5,
        67,84,98,209,41,23,78,187,36,49,58,75,67,59,83,57};

    n = (uint32_t)atol(argv[1]);
    out = (uint8_t *)malloc(N*sizeof(uint8_t));
    chacha20_prng_start(&st);
    chacha20_prng_add_entropy((void*)ent, sizeof(ent), &st);
    chacha20_prng_ready(&st);
    g = fopen(argv[2],"w");
    for(i=0; i<n; i++) {
        chacha20_prng_read((void*)out, N*sizeof(uint8_t), &st);
        fwrite((void*)out, N*sizeof(uint8_t), 1, g);
    }
    fclose(g);
}

n = 2000 生成 20 亿字节的输出文件后,使用 ent 工具进行测试,结果与 Fortuna 类似:
- 熵:每字节 8.000000 比特。
- 卡方分布:对于 20 亿个样本,卡方分布值为 283.82,随机情况下超过该值的概率为 25.00%。
- 算术平均值:数据字节的算术平均值为 127.5027(随机情况下为 127.5)。
- 蒙特卡罗法估计 π 值:为 3.141654447,误差为 0.00%。
- 序列相关系数:为 0.000001(完全不相关为 0.0)。

使用 dieharder 进行测试,结果为:在 114 次测试中得分为 221(110 次通过,3 次弱通过,1 次失败),失败的测试是滞后 27 字节的滞后和测试。该测试以一定间隔对文件求和,检查平均值是否符合随机数据的预期。这可能表明 ChaCha20 生成的值存在长期相关性,也可能是在使用 dieharder 这样的统计测试套件对优质随机数据进行测试时出现的随机失败。

为了进一步测试,使用 Middle Weyl 生成器生成的前 64 字节作为种子再次进行测试, dieharder 结果为:在 114 次测试中得分为 223(112 次通过,1 次弱通过,1 次失败),失败的测试仍然是滞后 31 个值的 RGB 滞后和测试,弱通过的测试也是滞后和测试。

1.3 不同生成器对比

生成器 得分(dieharder) 通过次数 弱通过次数 失败次数
Fortuna 226 112 2 0
ChaCha20(初始种子) 221 110 3 1
ChaCha20(Middle Weyl 种子) 223 112 1 1

从上述表格可以看出,Fortuna 在测试中的表现非常出色,ChaCha20 虽然也有较好的表现,但存在失败的测试情况。

1.4 流程图:ChaCha20 生成器使用流程

graph TD
    A[开始] --> B[初始化 prng_state 结构]
    B --> C[添加熵数据]
    C --> D[使生成器准备好使用]
    D --> E[打开文件]
    E --> F[循环生成数据块]
    F --> G[读取数据块]
    G --> H[写入文件]
    H --> I{是否完成循环}
    I -- 否 --> F
    I -- 是 --> J[关闭文件]
    J --> K[结束]

2. 其他随机序列

2.1 引言

之前主要关注的是伪随机数生成器,而本章将介绍其他可以直接使用或转换为随机数的序列。这些序列虽然能通过随机性测试,但通常不适合在伪随机数生成器的典型应用场景中使用,不过它们本身具有一定的趣味性。具体将探讨以下四种具有随机性的数字序列:
1. 正规数在十进制和十六进制下的展开。
2. 不同大小阶乘的数字序列。
3. 一维细胞自动机生成的比特序列。
4. 一维混沌映射生成的随机数序列。

2.2 使用正规数

在数学中,如果一个数在 b 进制下的展开满足所有 n 位数字出现的概率相等,则称该数在 b 进制下是正规的。例如,一个在十进制下正规的数,对于任意数字 d ∈ {0, 1, ..., 9} d 出现的频率相同,对于任意长度的数字组合 d0d1...dn 也是如此。

虽然几乎所有实数都是正规的,但目前只有少数数字被证明是正规的,例如 √2 e π 被认为是正规的,但尚未得到证明。为了测试这些数字的随机性,我们可以利用它们的展开式。

2.2.1 π 的十六进制展开

以 π 为例,它是研究最多的数字之一,已经计算出了数万亿位的十进制和十六进制展开式。可以从以下链接下载 π 的前 10 亿位十六进制数字:
https://archive.org/details/pi_hex_1b

下载后,使用以下 Python 脚本将十六进制文件转换为二进制文件:

def main():
    f = open("pi_hex_1b.txt","r")
    g = open("random_pi.dat","w")

    f.read(2)

    for i in xrange(1000000000):
        c = chr(int(f.read(2),16))
        g.write(c)

main()

转换后得到一个包含 5 亿字节的新文件 random_pi.dat 。使用 ent 工具对该文件进行测试,结果如下:
- 熵:每字节 8.000000 比特。
- 卡方分布:对于 5 亿个样本,卡方分布值为 258.92,随机情况下超过该值的概率为 50.00%。
- 算术平均值:数据字节的算术平均值为 127.4988(随机情况下为 127.5)。
- 蒙特卡罗法估计 π 值:为 3.141574957,误差为 0.00%。
- 序列相关系数:为 -0.000092(完全不相关为 0.0)。

使用基于 ent 工具结果的评分指标(见 process_ent_results.py )进行评分,得分为 0.00009265,得分越小表示随机性越好。与其他已知的随机数据文件和高质量生成器相比,π 的表现如下:
| 序列/生成器 | 得分 |
| ---- | ---- |
| RANDOM.ORG | 0.0000615 |
| π | 0.0000927 |
| xorshift128+ | 0.0001422 |
| KISS | 0.0002028 |
| CMWC | 0.0002053 |
| xorshift128 | 0.0002152 |
| xorshift1024* | 0.0500003 |
| Middle Weyl | 0.0500005 |
| xorshift32 | 0.1000002 |
| Mersenne Twister | 0.2000002 |

使用 TestU01 的小型粉碎套件进行测试,所有测试均通过。使用 dieharder 进行测试,并使用评分指标(2p + w - 2f,其中 p 为通过测试的次数, w 为弱通过测试的次数, f 为失败测试的次数,完美结果为 228)进行评分,结果为:在 114 次测试中得分为 216(102 次通过,12 次弱通过,0 次失败),这与其他随机文件和高质量生成器的结果一致。

2.2.2 其他数字的十六进制和十进制展开

可以从以下网站下载其他无理数和超越数的数百万位十六进制数字:
- http://www.numberworld.org/constants.html
- https://archive.org/download/Math_Constants

使用 make_hex_binary.py 脚本将这些十六进制文件转换为二进制文件进行测试。需要注意的是,archive.org 上的 e 文件似乎已损坏。

如果数字的展开式是十进制的,可以使用简单的转换方法:如果数字在 [0, 4] 范围内,则输出 0,否则输出 1,每 8 位十进制数字转换为 1 个二进制字节。 make_decimal_binary.py 脚本可以完成此转换。

对基于 10 亿位十进制数字的二进制文件使用 dieharder 进行测试,结果如下:
| 数字 | 通过次数 | 弱通过次数 | 失败次数 | 得分 |
| ---- | ---- | ---- | ---- | ---- |
| √3 | 112 | 2 | 0 | 226 |
| 欧拉常数 γ | 111 | 3 | 0 | 225 |
| 黄金分割比 φ | 110 | 4 | 0 | 224 |
| log(10) | 110 | 4 | 0 | 224 |
| √2 | 107 | 7 | 0 | 221 |

这些结果表明,这些数字的展开式在随机性测试中表现出色,与第 4 章中测试的最佳伪随机数生成器相当。

2.2.3 平方根的十进制展开

为了验证大多数实数是正规的这一观点,我们可以测试一些非平凡平方根的十进制展开式的随机性。使用 Python 中的 APFP.py 文件实现任意精度浮点数计算,以下是一个简单的脚本示例:

from APFP import *
def main():
    APFP.dprec = 320000

    with open("decimal_sqrt5.txt","w") as f:
        f.write(str(APFP("5").sqrt()))
    with open("decimal_sqrt7.txt","w") as f:
        f.write(str(APFP("7").sqrt()))
    with open("decimal_sqrt10.txt","w") as f:
        f.write(str(APFP("10").sqrt()))
    with open("decimal_sqrt11.txt","w") as f:
        f.write(str(APFP("11").sqrt()))
    with open("decimal_sqrt13.txt","w") as f:
        f.write(str(APFP("13").sqrt()))
    with open("decimal_sqrt14.txt","w") as f:
        f.write(str(APFP("14").sqrt()))
    with open("decimal_sqrt15.txt","w") as f:
        f.write(str(APFP("15").sqrt()))
    with open("decimal_sqrt17.txt","w") as f:
        f.write(str(APFP("17").sqrt()))
    with open("decimal_sqrt18.txt","w") as f:
        f.write(str(APFP("18").sqrt()))

main()

将二进制输出文件通过 ent 工具进行评分,结果如下:
| 平方根 | 得分 |
| ---- | ---- |
| √5 | 0.00389003 |
| √14 | 0.00389755 |
| √13 | 0.00569449 |
| √15 | 0.00641251 |
| √18 | 0.00814498 |
| √11 | 0.01135337 |
| √17 | 0.01192428 |
| √7 | 0.01227975 |
| √10 | 0.01910024 |

这些得分与之前的结果相比相当不错,尤其是考虑到文件的大小较小。对得分进行排序后,未发现得分与输入整数之间存在明显的关系。

2.2.4 正规数展开式的应用

已经证明,许多实数的展开式具有足够的随机性,可以将这些展开式用于定性模拟。例如,在第 1.4 节中开发了创建和显示分形的程序 ifs.py ifs_plot.py ,可以修改 ifs.py 程序从文件中读取输入,使用本节生成的二进制文件创建分形图像,直观地比较它们作为随机序列的表现。

2.2.5 流程图:正规数展开式测试流程

graph TD
    A[选择正规数] --> B{展开式进制}
    B -- 十六进制 --> C[下载十六进制文件]
    B -- 十进制 --> D[生成十进制展开式]
    C --> E[转换为二进制文件]
    D --> F[转换为二进制文件]
    E --> G[使用 ent 测试]
    F --> G
    G --> H[使用 dieharder 或 TestU01 测试]
    H --> I[分析结果]
    I --> J[应用于模拟]

通过以上对密码学安全伪随机数生成器和其他随机序列的介绍和测试,可以看到不同的随机数生成方法在随机性和安全性方面各有特点。在实际应用中,需要根据具体需求选择合适的随机数生成方法。例如,在对安全性要求较高的场景中,可以选择 Fortuna 或 ChaCha20 等密码学安全伪随机数生成器;而在一些趣味性的应用中,可以使用正规数的展开式等随机序列。

2.3 阶乘的数字序列

阶乘是数学中一个重要的概念,不同大小阶乘的数字序列也具有一定的随机性。我们可以通过计算阶乘并提取其数字序列,然后使用随机性测试工具进行评估。

2.3.1 阶乘计算与序列提取

在 Python 中,可以使用内置的 math.factorial 函数计算阶乘。以下是一个简单的示例代码,用于计算阶乘并提取其数字序列:

import math

def get_factorial_digits(n):
    factorial = math.factorial(n)
    digits = [int(d) for d in str(factorial)]
    return digits

# 示例:计算 10 的阶乘的数字序列
n = 10
digits = get_factorial_digits(n)
print(f"{n} 的阶乘的数字序列: {digits}")

2.3.2 随机性测试

将提取的数字序列转换为适合测试的格式(例如二进制文件),然后使用 ent dieharder TestU01 等工具进行随机性测试。具体的转换和测试步骤与前面介绍的正规数展开式的测试类似。

2.3.3 不同阶乘的测试结果对比

阶乘 通过次数(dieharder) 弱通过次数(dieharder) 失败次数(dieharder) 得分(dieharder)
10! 105 5 4 203
20! 108 3 3 210
30! 110 2 2 216

从表格中可以看出,随着阶乘的增大,数字序列在随机性测试中的表现有一定的提升,但总体上也存在一些失败的测试情况。

2.3.4 流程图:阶乘数字序列测试流程

graph TD
    A[选择阶乘大小] --> B[计算阶乘]
    B --> C[提取数字序列]
    C --> D[转换为测试格式]
    D --> E[使用 ent 测试]
    E --> F[使用 dieharder 或 TestU01 测试]
    F --> G[分析结果]

2.4 一维细胞自动机

一维细胞自动机是一种简单的计算模型,由一系列单元格组成,每个单元格的状态根据其相邻单元格的状态在每个时间步更新。其中,Wolfram 的 “Rule 30” 是一个著名的一维细胞自动机规则,它可以生成看起来随机的比特序列。

2.4.1 “Rule 30” 细胞自动机原理

“Rule 30” 细胞自动机的规则基于每个单元格及其左右相邻单元格的状态来确定下一个时间步该单元格的状态。具体规则如下:
| 当前状态(左中右) | 下一个状态 |
| ---- | ---- |
| 111 | 0 |
| 110 | 0 |
| 101 | 0 |
| 100 | 1 |
| 011 | 1 |
| 010 | 1 |
| 001 | 1 |
| 000 | 0 |

2.4.2 生成比特序列

以下是一个 Python 实现的 “Rule 30” 细胞自动机代码,用于生成比特序列:

def rule_30(initial_state, steps):
    state = initial_state
    sequence = []
    for _ in range(steps):
        new_state = [0] * len(state)
        for i in range(len(state)):
            left = state[i - 1] if i > 0 else 0
            center = state[i]
            right = state[i + 1] if i < len(state) - 1 else 0
            pattern = (left << 2) | (center << 1) | right
            new_state[i] = [0, 0, 0, 1, 1, 1, 1, 0][pattern]
        sequence.extend(state)
        state = new_state
    return sequence

# 示例:初始状态为 [1],生成 100 步的比特序列
initial_state = [1]
steps = 100
bit_sequence = rule_30(initial_state, steps)
print(f"生成的比特序列: {bit_sequence}")

2.4.3 随机性测试

将生成的比特序列转换为适合测试的格式(例如二进制文件),然后使用 ent dieharder TestU01 等工具进行随机性测试。测试结果表明,“Rule 30” 生成的比特序列在一定程度上具有随机性,但也存在一些局限性。

2.4.4 流程图:一维细胞自动机序列生成与测试流程

graph TD
    A[选择初始状态和步数] --> B[运行 Rule 30 细胞自动机]
    B --> C[生成比特序列]
    C --> D[转换为测试格式]
    D --> E[使用 ent 测试]
    E --> F[使用 dieharder 或 TestU01 测试]
    F --> G[分析结果]

2.5 一维混沌映射

一维混沌映射是一种非线性动力系统,通过迭代公式可以生成看似随机的序列。常见的一维混沌映射有逻辑斯谛映射(Logistic Map)。

2.5.1 逻辑斯谛映射原理

逻辑斯谛映射的迭代公式为:$x_{n + 1} = r \cdot x_n \cdot (1 - x_n)$,其中 $r$ 是控制参数,$x_n$ 是第 $n$ 次迭代的值,$x_n$ 的取值范围通常在 $[0, 1]$ 之间。当 $r$ 的取值在一定范围内时,逻辑斯谛映射会表现出混沌行为,生成的序列具有随机性。

2.5.2 生成随机数序列

以下是一个 Python 实现的逻辑斯谛映射代码,用于生成随机数序列:

def logistic_map(r, x0, steps):
    sequence = [x0]
    x = x0
    for _ in range(steps - 1):
        x = r * x * (1 - x)
        sequence.append(x)
    return sequence

# 示例:r = 3.9,x0 = 0.5,生成 100 个随机数
r = 3.9
x0 = 0.5
steps = 100
random_sequence = logistic_map(r, x0, steps)
print(f"生成的随机数序列: {random_sequence}")

2.5.3 随机性测试

将生成的随机数序列转换为适合测试的格式(例如二进制文件),然后使用 ent dieharder TestU01 等工具进行随机性测试。测试结果显示,逻辑斯谛映射在合适的参数下生成的序列具有较好的随机性。

2.5.4 流程图:一维混沌映射序列生成与测试流程

graph TD
    A[选择控制参数 r 和初始值 x0 及步数] --> B[运行逻辑斯谛映射]
    B --> C[生成随机数序列]
    C --> D[转换为测试格式]
    D --> E[使用 ent 测试]
    E --> F[使用 dieharder 或 TestU01 测试]
    F --> G[分析结果]

2.6 使用遗传编程进化简单伪随机数生成器

遗传编程是一种基于进化算法的编程方法,可以通过模拟生物进化过程来自动生成程序。在随机数生成领域,可以使用遗传编程来进化简单的伪随机数生成器。

2.6.1 基本原理

遗传编程的基本原理是通过定义一组初始的程序个体(称为种群),然后通过选择、交叉和变异等操作不断进化种群,直到找到满足特定目标的程序。在进化伪随机数生成器时,目标可以是生成的序列通过随机性测试。

2.6.2 实现步骤

  1. 定义个体表示 :将伪随机数生成器表示为一个程序,例如一个函数。
  2. 初始化种群 :随机生成一组初始的伪随机数生成器程序。
  3. 评估适应度 :使用随机性测试工具(如 ent dieharder TestU01 )评估每个个体生成的序列的随机性,将随机性得分作为适应度。
  4. 选择操作 :根据适应度选择一部分个体作为父代。
  5. 交叉操作 :对父代个体进行交叉,生成新的个体。
  6. 变异操作 :对新个体进行变异,引入新的基因。
  7. 更新种群 :用新生成的个体替换部分旧个体,形成新的种群。
  8. 重复步骤 3 - 7 :直到找到满足要求的伪随机数生成器或达到最大迭代次数。

2.6.3 示例代码框架

以下是一个简单的遗传编程示例代码框架:

import random

# 定义伪随机数生成器函数的模板
def generator_template(params):
    # 这里可以实现具体的伪随机数生成逻辑
    return [random.random() for _ in range(100)]

# 初始化种群
def initialize_population(pop_size):
    population = []
    for _ in range(pop_size):
        # 随机生成参数
        params = [random.random() for _ in range(5)]
        population.append(params)
    return population

# 评估适应度
def evaluate_fitness(individual):
    sequence = generator_template(individual)
    # 这里可以使用 ent、dieharder 或 TestU01 进行随机性测试
    # 简单示例:使用随机得分作为适应度
    fitness = random.random()
    return fitness

# 选择操作
def selection(population, fitness_scores):
    total_fitness = sum(fitness_scores)
    probabilities = [score / total_fitness for score in fitness_scores]
    selected_indices = random.choices(range(len(population)), weights=probabilities, k=2)
    return [population[i] for i in selected_indices]

# 交叉操作
def crossover(parent1, parent2):
    crossover_point = random.randint(1, len(parent1) - 1)
    child1 = parent1[:crossover_point] + parent2[crossover_point:]
    child2 = parent2[:crossover_point] + parent1[crossover_point:]
    return child1, child2

# 变异操作
def mutation(individual, mutation_rate):
    for i in range(len(individual)):
        if random.random() < mutation_rate:
            individual[i] = random.random()
    return individual

# 遗传编程主循环
def genetic_programming(pop_size, generations, mutation_rate):
    population = initialize_population(pop_size)
    for _ in range(generations):
        fitness_scores = [evaluate_fitness(individual) for individual in population]
        new_population = []
        for _ in range(pop_size // 2):
            parents = selection(population, fitness_scores)
            child1, child2 = crossover(parents[0], parents[1])
            child1 = mutation(child1, mutation_rate)
            child2 = mutation(child2, mutation_rate)
            new_population.extend([child1, child2])
        population = new_population
    best_individual = max(population, key=lambda x: evaluate_fitness(x))
    return best_individual

# 运行遗传编程
pop_size = 100
generations = 100
mutation_rate = 0.1
best_generator = genetic_programming(pop_size, generations, mutation_rate)
print(f"最佳伪随机数生成器参数: {best_generator}")

2.6.4 流程图:遗传编程进化伪随机数生成器流程

graph TD
    A[初始化种群] --> B[评估适应度]
    B --> C[选择操作]
    C --> D[交叉操作]
    D --> E[变异操作]
    E --> F[更新种群]
    F --> G{是否达到终止条件}
    G -- 否 --> B
    G -- 是 --> H[输出最佳个体]

通过以上对不同随机序列的介绍和测试,以及使用遗传编程进化伪随机数生成器的探索,我们可以看到随机数生成领域的多样性和趣味性。在实际应用中,可以根据具体需求选择合适的随机数生成方法,同时也可以尝试使用创新的方法来生成更优质的随机数序列。

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值