在 0..1 之间将 u64 转换为 f64

Question

我需要一个非常快速的伪随机数生成器来用于我一直在从事的项目。到目前为止，我已经实现了 xoshiro 算法，并且可以生成伪随机 u64。但是，我需要将这些 u64 转换为 0 到 1 范围内的浮点值。

我主要使用this和this作为参考。

出于某种原因，我无法接近我想要的行为；这让我很困惑，因为我使用了与here完全相同的方法。尽管我看到实现没有任何差异，但我得到了不同的结果。

    let seeds: [u64; 64] = core::array::from_fn(|i| i as u64);

    let bitshift12 = u64x64::splat(12);
    let bitshift25 = u64x64::splat(25);
    let bitshift27 = u64x64::splat(27);
    
    let bitshift52 = u64x64::splat(52);
    
    let mut random_states = Simd::from(seeds);
    
    random_states ^= random_states >> bitshift12;
    random_states ^= random_states << bitshift25;
    random_states ^= random_states >> bitshift27;
    
    random_states = random_states | ((u64x64::splat(1023) + u64x64::splat(0)) << bitshift52);
    
    let mut generated = Simd::<f64, 64>::from_bits(random_states);
    
    println!("{:?}", generated);

输出：

[1.0, 1.0000000074505808, 1.0000000149011616, 1.0000000223517425, 1.0000000298023235, 1.0000000372529039, ...]

显然我没有做正确的事情，因为最后几位小数是根据需要“随机”的。为什么我不能正确地将它们向上移动？

如果有人指出我的错误，我将不胜感激。

Answer 1

该序列看起来就像将小整数填充到指数从 1.0 开始的

f64

位模式的尾数中所得到的，因此你得到 1.0 加上少量的值。不像

0, 1, 2, 3, ...

那么小； https://www.binaryconvert.com/result_double.html?decimal=049046048048048048048048048048055052053048053056048056显示数字由

f64

位模式

0x3FF0000002000001

表示，其中仅设置了2位尾数。

不过，这看起来像是在 xoshiro 迭代之后从 seed = 1 开始得到的位模式。请注意，第一个移位是向右移动，移出唯一留下 0 的位。下一步是向左移动，导致两个设置位。然后最后右移 27 将它们都移出，再次与 0 进行异或，使它们保持不变。

因此，在 xoshiro 的一步之后，

seeds[i] = i

的极其非随机种子就会导致这些非随机尾数。（并且

seeds[0]

永远不会变成非零；xoshiro 需要非零种子，因为移位和异或永远不能从零创建非零位。）

如果您确实有统一的随机

u64

值（例如，使用真实种子，或让生成器对非零种子运行多次迭代），则将它们与

1.0

中的指数进行“或”运算也会使指数随机化，导致巨大的价值。但其数量级始终大于 1.0，除了具有全 1 指数的 NaN（如果尾数为零，则为无穷大）。也是随机的标志。 OR 无法清除位，并且由于指数偏差，IEEE 浮点数大小随着整数位模式的增加而单调增加。 https://en.wikipedia.org/wiki/Double- precision_floating-point_format

如果您屏蔽随机

u64

以仅保留低 52 位，这样您只需随机化尾数，您就可以轻松地在

[1.0, 2.0)

中获得统一的随机数。正如 Chux 所说，在您链接的问答中（如何将多个伪随机位转换为 0 和 1 之间的统计随机浮点值？），从中减去

1.0

是获得

[0.0, 1.0)

数字。

越接近 0.0（指数越小），在减去两个附近的数字后，尾数的尾随零就越多：指数越小，可表示的值越接近，但我们想要均匀分布。此方法只有 52 位熵。这可能没问题，但理论上您可以检查指数字段并使用可变计数移位 + OR 来随机化低尾数位。

Chux 的其他方法（保值转换，如 C 类型转换）然后除法（实际上是乘以逆）如果没有 AVX-512 来进行从

u64

到

f64

的打包转换，则无法在 x86 上高效地完成。 如何使用 SSE/AVX 高效执行 double/int64 转换？ - 它需要多条指令，不仅仅是替换指数和减法。 （使用 AVX-512，替换指数字段也变得更加高效，只需一个

vpternlogd

以及覆盖指数+符号字段的位掩码。）

---

顺便说一句，除非编译器将

u64x64 bitshift12

优化回标量立即数，否则使用移位计数的 SIMD 向量看起来效率不高。至少在 x86 和 AArch64 上，向量移位可以使用标量计数，因此我希望

random_states >> 12

能够编译为

vpsrlq ymm, ymm, 12

（使用 AVX2），而不需要 AVX2 变量计数移位以及像

 这样的计数的向量常量vpsrlvq ymm, ymm, ymm

。 （Zen 2 上的每 2 个周期吞吐量与立即计数移位的每周期 1 个吞吐量：https://uops.info/。但在 Zen 3 及更高版本以及 Intel Skylake 及更高版本上，吞吐量是相同的。但是如果编译器实际上必须从 64x

u64

数组中加载计数向量，那就太糟糕了）。

我假设

u64x64::splat(1023) + u64x64::splat(0)

是用来处理不同的指数域的，但为什么要进行向量相加呢？只需

u64x64::splat((1023 + offset) << 52)

即可为您提供 1.0 的指数字段，使用标量常量进行所有数学运算，甚至不会诱使编译器在运行时执行此操作。