随机算法(数值概率&舍伍德)

随机算法的分类
• 数值概率算法常用于数值问题的求解。将一个问题的计算与 某个概率分布已经确定的事件联系起来，求问题的近似解。 这类算法所得到的往往是近似解，且近似解的精度随计算时 间的增加而不断提高。在许多情况下，要计算出问题的精确 解是不可能或没有必要的，因此可以用数值随机化算法得到 相当满意的解。 • 蒙特卡罗算法用于求问题的准确解，但得到的解未必是正确 的。蒙特卡罗算法以正的概率给出正解，求得正确解的概率 依赖于算法所用的时间。算法所用的时间越多，得到正确解 的概率就越高。一般给定执行步骤的上界，给定一个输入，
选定也可用系统时间自动产生。函数Random的输入参数
n65536是一个无符号整形数，返回[0，n-1]范围内的一 个随机整数。函数fRandom()返回[0，1]内的一个随机实数。
//随机数类 Const unsigned long maxshort = 64436L; Const unsigned long multiplier = 1194211693L; Const unsigned long adder = 12345L; class RandomNumber { private: //当前种子 unsigned long randSeed; public: //构造函数，默认值0表示由系统自动产生种子 RandomNumber (unsigned long s = 0); //产生0到n-1之间的随机整数 unsigned short Random (unsigned long n); //产生[0,1)之间的随机实数 double fRandom (void); }
*
i 1
用 I 作为I的近似值。
b
假设设要计算的积分形式为： I
g ( x ) dx
，可以任选一个有简便
b
a 方法可以进行抽样的概率密度函数f(x)，使其满足：
（1） f ( x ) 0 , 当 g ( x ) 0时（ a x b ） ; （2） f ( x ) dx 1
7.2.2 计算定积分
（1）用随机投点法计算定积分 设f(x)是[0，1]上的连续函数，且0f(x)1。需要计算的积
1
分 I 为

0
f ( x ) dx
, 积分 I 等于图中的面积G。 在图所示单位正方形内均匀地作投点
试 验，则随机点落在曲线下面的概率为： 1 f (x) 1
Pr { y f ( x )}
有很坏的性能。由于随机性选择比最优选择省时间，因此
引入随机化算法可以在很大程度上降低算法的复杂度。
随机算法的基本特征
• 随机算法对所求解问题的同一个实例用同一随机算法求解 两次可能得到完全不同的效果。这两次求解所需要的时间， 甚至所得到的结果都可能会有相当大的差别。
• 包括
– 数值概率算法 – 蒙特卡罗（Monte Carlo）算法 – 拉斯维加斯（Las Vegas）算法 – 舍伍德（Sherwood）算法
f (x)
随机变量，而且
E ( g ( i ))
*

b
g ( x ) f ( x ) dx
*
a

b
g ( x ) dx I
a
由强大数定律：
Pr ( lim

n
1 n

n
g ( i ) I ) 1
*

i 1
若选 I

1 n

n
g ( i ) ，则 I 依概率1收敛于I，平均值法就是
如果所遇到的积分形式为： I
b

a
f ( x ) dx
其中，a和b为有限值，被积函数f(x)在区间[a，b]中 有界，并用M，L分别表示其最大值和最小值。此时可做 变量代换x=a+(b-a)z，将所求积分变为： cI * d I 式中： c ( M L )( b a )
1
d L (b a )

0 0
dydx

0
f ( x ) dx
假设向单位正方形内随机地投入n个 点(xi，yi)，i = 1,2，…，n。随机点（xi，yi）落入G内，则yi f（xi）。 m 如果有m个点落入G内，则随机点落入G内的概率，即 I n
double Darts (int n) { // 用随机投点法计算定积分 static RandomNumber dart; int k=0; for (int i=1;i<=n;i++) { double x=dart.fRandom(); double y=dart.fRandom(); if (y<=f(x)) k++; } return k/double(n) }
第7章 随机算法
学习要点
了解随机算法的基本特征 理解产生伪随机数的算法 掌握数值随机化算法的设计思想 掌握舍伍德算法的设计思想 掌握拉斯维加斯算法的设计思想
掌握蒙特卡罗算法的设计思想
随机算法的基本特征
前面各章讨论的算法的每一个步骤都是确定的，而本 章讨论的随机算法允许算法在执行过程中随机地选择一下 计算步骤。 在许多情况下，一般算法比较复杂，性能较差，很 多具有很好平均运行时间的算法，在最坏的情况下，却具
得到反面，返回1表示得到正面。下面的算法
TossCoins模拟抛10次硬币这一事件50000次。用
head[i] (0 i 10)记录这50000此模拟恰好得到i次 正面的次数。最终输出模拟抛硬币事件得到正面事 件的频率图，如下图所示：
0 1 2 3 4 5 6 7 8 9 10
* *
*
*
* *
//产生种子 RandomNumber:: RandomNumber(unsigned long s) { if (s == 0) randSeed = time (0); //用系统时间产生种子 else randSeed = s; //由用户提供种子 } //产生0~n-1之间的随机数 Unsigned short RandomNumber:: Random (unsigned long n) { randSeed = multiplier * randSeed + adder; return (unsigned short) ( (randSeed >>16) % n ); }
int TossCoins (int numberCoins)
7.2 数值随机化算法
7.2 数值随机化算法
7.2.1 用随机投点法计算值
设有一半径为r的圆及其外切四边形。向该正
方形随机地投掷n个点。设落入圆内的点数为k。 由于所投入的点在正方形上均匀分布，因而所投 入的点落入圆内的概率为 所以当n足够大时，k 与n之比就逼近这一概率。 从而
I
*

0
f ( z ) dz
1 M L
*
f (z)
*
[ f ( a (b a ) z L ]
（2）用平均值法计算定积分
任取一组相互独立、同分布的随机变量{ξi},ξi在[a,b] 中服从分 布律f(x), 令
g (x)
*
g (x)
，则{g*(ξi)}也是一组互独立、同分布的
函数Random在每次计算时，用线性同余式计算 新的种子randSeed。它的高16位的随机性较好。将 randSeed右移16位得到一个0~65535间的随机整数， 然后再将此随机整数映射到0~（n-1）范围内。
对于函数fRandom，先用函数Random(maxshort) 产生一个0~（maxshort-1）之间的整形随机序列，将每 个整型随机数除以maxshort，就得到[0,1）区间中的随 机实数。
n 1, 2 ,
• 其中b0，c0，dm。d称为该随机序列的种子。如何选取该 方法中的常数b、c和m直接关系到所产生的随机序列的随机性 能。这是随机性理论研究的内容，已超出本书讨论的范围。从 直观上看，m应取得充分大，因此可取m为机器大数。
为了在设计概率算法时便于产生所需的随机数，建立 一个随机数类RandomNumber。该类包含一个需有用户初始 化的种子randseed。给定初始种子后，即可产生与之相应 的随机序列。种子seed 是一个无符号长整型数，可由用户
7.1随机数
7.1随机数
• 随机数在随机化算法设计中扮演着十分重要的角色。在现实计 算机上无法产生真正的随机数，因此在随机化算法中使用的随 机数都是一定程度上随机的，即伪随机数。 • 线性同余法是产生伪随机数的最常用的方法。由线性同余法产 生的随机序列a0,a1,…,an，满足：
a0 d a n ( ba n 1 c ) mod m
*
* * * *
模拟抛硬币得到的正面事件频率图
void main (void) { //随机抛硬币 { // 模拟随机抛硬币事件 static RandomNumber coinToss; const int NCOINS = 10; const long NTOSSES = 50000L; int i, tosses = 0; // heads[i]是得到i次正面的次数 int TossCoins (int numberCoins)0; i < numberCoins; i++) for (i = long i, heads[NCOINS+1]; // Random (2) = 1 表示正面 { //随机抛硬币 int j, position; // 初始化数组heads tosses += coinToss.Random (2); static RandomNumber coinToss; for (j=0;j< NCOIN+1; j++) return tosses; heads[j] = 0; int i, tosses = 0; } // 重复50000次模拟事件 for (i=0;i< NTOSSES; i++) for (i = 0; i < numberCoins; i++) heads[TossCoins(NCOINS)]++; // 输出频率图// Random (2) = 1 表示正面 for (i=0;i<=NCOINS; i++) tosses += coinToss.Random (2); { position = int (float (heads[i])/NTOSSES*72); return tosses; cout<<setw(6)<<i<<" "; } for (j=0;j<position-1;j++) cout<<" "; cout<<"*"<<endl; } }