用欧拉计划学Rust编程,第323题,随机整数按位或运算

用欧拉计划学Rust编程,第323题,随机整数按位或运算

题目描述:

理解题意,一个数初始值为0,不断与一个随机数进行‘或’运算,直到32个二进制位都为1,所需要步数记为N。进行无数数这样的试验,问N的平均值是多少。

第一步:

先用蒙特卡罗模拟方法试验一下,看看答案的大概范围。

fn trial() -> u32 {
    let mut x = 0_u32;
    for i in 1.. {
        x |= rand::random::<u32>();
        if x == 0xFFFFFFFF {
            return i;
        }
    }
    0
}

运行这个trail()函数10万次,取平均值,可以得到答案在6.35附近。

第二步:

原题的x是32位整数,可以从最简单的情况入手,假设x为1位的整数。
记此时的期望值为E(1)。
进行第一步操作之后,有两种可能:
第一种情况,1/2的概率,直接变成了1,即一步完成
第二种情况,1/2的概率,仍是0,回到了初始情况,即需要1+E(1)步完成。
这时有:E(1) = 1/2 + 1/2 * (1 + E(1))
可以解出:E(1) = 2

第三步:

假设x是2位整数,此时的期望值记为E(2)。
进行第一步操作之后,有三种可能:
第一种情况,1/4的概率,直接变成了(11)binary,即一步完成
第二种情况,2/4的概率,有一位变为1,即(01)binary或(10)binary,这时是E(1)时的情况,期望值为1+E(1)。
第三种情况,1/4的概率,仍为0,这时是1+E(2)。
列出式子,E(2) = 1/4 + 2/4 * (1+ E(1)) + 1/4 * (1 + E(2))
3/4 * E(2) = 8/4
E(2) = 8/3

第四步:

这时,你可能还没有总结出规律,再来一次,假设x是3位整数,此时的期望值记为E(3)。
进行第一步操作之后,有四种可能:
第一种情况,1/8的概率,直接变成了(111)binary,即一步完成
第二种情况,3/8的概率,有两位变为1,即(011)binary或(101)binary或(110)binary,这时是E(1)时的情况,期望值为1+E(1)。
第三种情况,3/8的概率,有一位变为1,即(001)binary或(010)binary或(100)binary,这时是E(2)时的情况,期望值为1+E(2)。
第四种情况,仍为0,这时是1+E(3)。
列出式子,E(3) = 1/8 + 3/8 * (1+ E(1)) + 3/8 * (1 + E(2)) + 1/8 * (1 + E(3))
7/8 * E(3) = 1/8 + 3/8 * (1+ E(1)) + 3/8 * (1 + E(2)) + 1/8
7/8 * E(3) = 1/8 + 3/8 * 3 + 3/8 * 11/3 + 1/8
E(3) = 22/7

第五步,现在可以总结出最一般的规律:

E(N) = 1/2^N + C(N,1)/2^N * (1+E(1)) + C(N,2)/2^N * (1+E(2)) + ... + C(N,N-1)/2^N * (1+E(N-1)) + 1/2^N (1+E(N))
这里的C(N, i)是排列组合数,从N个数里取i个数的组合数。
耐心地把右侧的E(N)挪到左侧:
(2N-1)/2N * E(N) = 1/2^N + C(N,1)/2^N * (1+E(1)) + C(N,2)/2^N * (1+E(2)) + ... + C(N,N-1)/2^N * (1+E(N-1)) + 1/2^N
最后可以得到:
E(N) = [1 + C(N,1) * (1+E(1)) + C(N,2) * (1+E(2)) + ... + C(N,N-1) * (1+E(N-1)) + 1] / (2^N-1)
我们知道E(1)=2,不断套用上面的公式,就可以得到最终的答案。
求组合数的函数:

fn comb(n: usize, r: usize) -> f64 {
    let mut c = 1_f64;
    for i in 0..r {
        c *= ((n - i) as f64) / (r - i) as f64;
    }
    c
}

主程序:

let mut e = vec![0.0, 2.0];
for n in 2..=32 {
    let mut f = 1.0 as f64;
    for r in 1..=n {
        f += comb(n, r) * (1.0 + e[n - r]);
    }
    f = f / (2_f64.powf(n as f64) - 1.0);
    e.push(f);
    println!("E({}) = {:.10}", n, f);
}


----==== Email: slofslb (GTD) qq.com 请将(GTD)换成@ ====----
版权声明:自由转载-非商用-非衍生-保持署名(创意共享3.0许可证
作者:申龙斌的程序人生

---- 魔方桥牌象棋、游戏人生...
---- BASIC、C++、JAVA、C#HaskellObjective-COpen Inventor、程序人生...
---- GTD伴我实现人生目标
---- 区块链生存训练
---- 用欧拉计划学Rust编程
---- 申龙斌的读书笔记(2011-2019)
----
posted @ 2020-09-07 22:13  申龙斌的程序人生  阅读(278)  评论(0编辑  收藏  举报