游戏用幸运算法及其优化

KaaPexei

一般游戏中带有随机概率的元素，很容易就会出现一个幸运的属性，这个属性越高，则获取高级物品的概率越高，低级物品的概率越低。
在开始算法之前，我们先精确描述一下这个问题：
我们需要在区间 \mathbb V = [0,1) 上工作，有一个幸运度入参 n\in\Bbb N^+ ，设计一个分布 D(n) ，让随机变量 X\sim D(n) 随着 n 的增大，值逐渐向 1 偏移，品级是在 \Bbb V 上的区间集 \mathbf R = \{[0,a_1),[a_1,a_2)\cdots[a_n,1)\} ，其中有 0<a_1<a_2<\cdots<a_n<1 .
基础方法

实现这个的一种简单方法是掷骰子，掷多次骰子取最好的结果，幸运度越高，掷骰子的次数越多，抽到高级物品的概率就越高。
也就是 X=\max\{r_n\},\quad r_n= \underset{x\in \Bbb V}{\rm rand} . 我们可以很直接地把它转化成代码：
use rand::Rng;
use num_enum::TryFromPrimitive;

use std::convert::TryFrom;

#[derive(Debug, Eq, PartialEq, TryFromPrimitive)]
#[repr(u8)]
enum Rarity {
    Common, // 普通
    Great, // 优秀
    Execellent, // 精良
    Epic, // 史诗
    Legend, // 传说
}

/// 稀有度区间
const RARITY_SECTIONS: [(f64,f64);5] = [
    (0.00, 0.50),
    (0.50, 0.75),
    (0.75, 0.90),
    (0.90, 0.98),
    (0.98, 1.00),
];

pub fn generate_random_rarity(luck: usize) -> Rarity {
    let mut rd: f64 = 0.0;
    for _ in 0..luck {
        let r = rand::random(); // 生成 [0,1) 之间的浮点数
        rd = max(rd,r);
    }
    for (i, (from, to)) in (0..RARITY_SECTIONS.len()).zip(RARITY_SECTIONS) {
        // 随机到的数在这个区间内
        if from <= rd && rd < to {
            return Rarity::try_from(i as u8).unwrap_unchecked();
        }
    }
    unreachable!();
}
但是这种方法在幸运度较高时需要连续投掷多次骰子，可以考虑进行优化。
优化方法

首先，我们知道随机变量 X 根据稀有度被划分成多个区间，而对于我们来说，只需要知道随机变量在某个区间内即可，不需要知道这个变量确切的值。如果我们能提前计算出在某一幸运度下的区间，然后提前储存起来，这样调用时就只需要使用一次随机+查表就能解决。
考虑将 \Bbb V 分为三个区间 [0, a_1), [a_1, a_2), [a_2,1) 的情况，很容易知道：

注：由于区间划分是任意的，对于超过三个区间的情况，我们可以从某个区间开始，把左边的看作一整个区间，右边的看作一整个区间，这样仍然会回归到三个区间的情况。

\begin{align} P\Big(X\in[0,a_1)\Big) &=  P\Big(\max \{r_n\}\in[0,a_1)\Big) \\ &= P\Big(\forall r \in \{r_n\}, r\in[0,a_1)\Big) \\ &=P^n\Big(r\in[0,a_1)\Big) \\ &=P^n\left(\underset{x\in\Bbb V}{\rm rand \,} \in [0,a_1)\right) \end{align}  
其中，由于 \underset {x\in \Bbb V}{\rm rand} 是均匀的，所以 P\left (\underset{x\in\Bbb V}{\rm rand \,} \in [0, a_1)\right)=\cfrac{a_1-0}{1-0}=a_1 ，所以有：
P(X\in[0,a_1))=P^n\left (\underset{x\in\Bbb V}{\rm rand \,} \in [0, a_1)\right)=a_1^n \tag{1}  我们先把第二个区间放在一边，先来处理第三个区间：
\begin{align} P\big(X\in[a_2,1)\big)  &= 1-P\big(X\in[0,a_2)\big) \\ &= 1-a_2^n \tag{2} \end{align}  
于是再来处理第二个区间：
\begin{align} P\big(X\in[a_1,a_2)\big) &= 1-P\big(X\in[0,a_1) \vee X\in [a_2,1)\big) \\ &= 1-\Big(P\big(X\in[0,a_1)\big) + P\big(X\in[0,a_1)\big)\Big) \\ \end{align}  
由 (1),(2) 得：
P\big(X\in[a_1,a_2)\big) = 1-(a_1^n+1-a_2^n) = a_2^n-a_1^n \tag{3}  
接下来进行验证：

• 将 a_1=0,a_2=1 代入，得 P (X\in[0,1))=1-0=1 ;
  
• 将 a_1=0 代入，得 P (X\in[0, a_2))=a_2^n-0=a_2^n ，满足 (1) ；
  
• 将 a_2=1 代入，得 P (X\in[a_1,1))=1-a_1^n ，满足 (2)
  

因此， (3) 实际上推广到了任意的稀有度区间，而不局限于中间区间。

注：上面的推导实际上用了两个假设，即 rand() 生成的随即数是均匀的，以及多次 rand() 之间是独立的，互不干扰。

根据上述事实，我们可以从稀有度区间集 \bf R 和幸运度 n 计算出随机变量属于任何一个区间的概率，用这些概率来划分出新的区间集 \bf R' ，在 \bf R' 上做均匀随机，只用随机一次，就可以得到与原先在 \bf R 上随机 n 次等价的结果。
use std::convert::TryFrom;

#[derive(Debug, Eq, PartialEq, TryFromPrimitive)]
#[repr(u8)]
enum Rarity {
    Common, // 普通
    Great, // 优秀
    Execellent, // 精良
    Epic, // 史诗
    Legend, // 传说
}

/// 稀有度区间 R
const RARITY_SECTIONS: [(f64,f64);5] = [
    (0.00, 0.50),
    (0.50, 0.75),
    (0.75, 0.90),
    (0.90, 0.98),
    (0.98, 1.00),
];

/// 经过计算后的区间集 R'
static mut G_CALCULATED_SECTIONS: [(f64,f64),5] = [
    (0.00, 0.50),
    (0.50, 0.75),
    (0.75, 0.90),
    (0.90, 0.98),
    (0.98, 1.00),
];

pub fn calculate_sections(luck: usize) {
    let mut previous_pow: f64 = 0.0;
    let it = G_CALCULATED_SECTIONS.iter_mut();
    for (from, to) in RARITY_SECTIONS {
        let section = it.next();
        let next_pow = to.powi(luck as u32);
        section = (previous_pow, next_pow);
        // 由于两个连续区间的端点值是共享的，所以只用计算一次 powi
        previous_pow = next_pow;
    }
}

pub fn generate_random_rarity() -> Rarity {
    let rd = rand::random();
    let it = (0..G_CALCULATED_SECTIONS.len()).zip(G_CALCULATED_SECTIONS);
    for (i, (from, to)) in it {
        // 随机到的数在这个区间内
        if from <= rd && rd < to {
            return Rarity::try_from(i as u8).unwrap_unchecked();
        }
    }
    unreachable!();
}
在实际使用中，由于品级数量非常有限，幸运值往往也在一定的范围内（例如 0~999），所以打表是很好的选择。本例中 5 个品级打表只有 5 * 64 * 2 * 1000 bit = 640Kib = 80KiB 的内存占用，而如果去除冗余（这里为了编码方便，存入的数据包括两个端点，造成了冗余，实际上只存末端点即可，并且首尾的 0 和 1 不需要存储），则只需要存储4 * 64 * 1000 bit = 32KiB，因此可以在进入战斗之前提前计算好所有的幸运值对应的 \bf R' ，或者直接在编译期计算打表，这样在战斗中就可以查表调用，减小爆一地金装时的 CPU 负载。

后续

X 服从的分布 D(n) 究竟是什么呢？我们考虑一个无穷小区间 s = [x, x+\mathrm dx) ， X 在这个小区间内概率的导数 \mathrm dP 满足：

\frac{\mathrm d}{\mathrm dx}P(X\in s) = \frac{\mathrm d}{\mathrm dx}\big((x+\mathrm dx)^n-x^n\big)

二项式定理展开，并舍去高阶无穷小量，得：

\begin{align} \frac{\mathrm d}{\mathrm dx}P(X\in s)  &= \frac{\mathrm d}{\mathrm dx}\big((x+\mathrm dx)^n-x^n\big) \\ &= \frac{\mathrm d}{\mathrm dx}\Big(\sum_i^n C_n^ix^{n-i}\mathrm d^{i}x-x^n\Big) \\ &= \frac{\mathrm d}{\mathrm dx}\Big(x^n+C_n^1x^{n-1}\mathrm dx  +\sum_{i=2}^nC_n^ix^{n-i}\mathrm dx - x^n\Big) \\ &= \frac{\mathrm d}{\mathrm dx}\Big(C_n^1x^{n-1}\mathrm dx + o\Big) \\&= C_n^1x^{n-1} \\ &=nx^{n-1}\end{align}

也就是说， X 的概率密度分布函数 f_X(x)=nx^{n-1},\;(x\in[0,1)) 。我们做个图看看这个分布随着 n 的增大，变化合不合理（下面是 n 取 1~100 之间整数时概率密度函数的图像）：

n 取 1~100 之间整数时图像
https://www.zhihu.com/video/1595074434385203200
对这个概率密度函数积分，得到累计密度函数 F_X(x)=x^n ，当 x=1 时，有 F_X(1)=1 恒成立，此时我们发现，就算我们对 n 取小数，也不会影响到这个分布的成立，要知道我们是从离散的掷骰子问题开始的！因此，我们把这个幸运计算公式中的 n 推广到了 [0,+\infty) 。
同时，我们发现这个函数的变化在 n 较小的时候比较大，而 n 变大之后，变化就不明显了。因此，我们可以在实现幸运算法时在幸运值前预乘一个较小的常系数 \lambda ，这样让 n 在合理的取值范围内都可以让概率有较为明显的变化，这个就看工程里实际应用如何自行判断。