Contributors

• 二、概率算法
  • 2.1 数字算法
  • 2.2 Sherwood 算法
  • 2.3 Las Vegas 算法
    • 2.3.1 八皇后问题
    • 2.3.2 模 p 平方根
    • 2.3.3 整数的因数分解
  • 2.4 Monte Carlo 算法
    • 2.4.1 主元素问题
    • 2.4.2 素数测定
    • 2.4.3 矩阵乘法验证

转载自 cnblogs/学科笔记/中科大-算法分析与设计-概率算法复习知识点.md at main · BamLubi/cnblogs 中科大-算法分析与设计-概率算法复习知识点 - 小陆斑比 - 博客园

二、概率算法

概率算法分类：

• 数字算法：求近似解
• Monte Carlo 算法：解判定问题或只有一个正确解的问题。答案未必正确。
• Las Vegas 算法：决不返回错误答案。
• Sherwood 算法：总是给出正确答案。

2.1 数字算法

找到数字问题的近似解

$\Pi$值计算：投飞镖模拟，单位圆和外切正方形。

$$\begin{gathered} >\frac{落在圆内}{总数}=\frac{k}{n}=\frac{\pi r^2}{4r^2}=\frac{\pi }{4} \\ >\Rightarrow \pi =4\frac{k}{n}> \end{gathered}$$

Darts(n){
    k=0;
    for i=1 to n do {
        x = uniform(0,1); y = uniform(0,1);
        if(x*x+y*y<=1) then k++;
    }
    return 4k/n;
}

数字积分：Monte Carlo 积分，但不是本课中的 MC 算法。

计算积分$S={\int }_0^{!}f(x)dx$

HitorMiss 算法：单位面积正方形投飞镖，落入积分函数下方的为 k，$\frac{k}{n}=\frac{S}{1}$

HitorMiss(f, n){
    k=0;
    for i=1 to n do {
        x = uniform(0,1); y = uniform(0,1);
        if(y<=f(x)) then k++;
    }
    return k/n;
}

Crude 算法：积分区间上随机取点，取算术平均值乘区间宽度。${\int }_a^{b}f(x)dx=(b-1)\frac{1}{n}{\sum }_{i=1}^{n}f(x_i)$

Crude(f,n,a,b){
    sum = 0;
    for i=1 to n do{
        x = uniform(a,b);
        sum = sum + f(x);
    }
    return (b-a)*sum/n;
}

相对而言，Crude 算法方差更小，但是 HitorMiss 算法在相同时间迭代次数更多。

解题思路：$\pi =4{\int }_0^1\sqrt{1-x^2}dx$

Trapezoid 确定性算法：梯形算法，将区间划分为 n-1 个子区间

$S=\delta \ast (f(a+\delta )+f(a+2\delta )+\cdots +\frac{f(a)+f(b)}{2})$

Trapezoid(f,n,a,b){
    // 假设n>=2
    delta = (b-a)/(n-1);
    sum = (f(a)+f(b))/2;
    for(x=a+delta; x<=b-delta; x+=delta){
        sum = sum + f(x);
    }
    return sum * delta;
}

确定性算法有时求不出解，概率算法适用于高维定积分。

概率计数：估计整数值，如集合的大小或者集合中不同对象数目的估计。

求集合的大小：有放回的随机抽取元素，k 是第 1 次出现重复元素之前所选出的数目。

$n足够大时,k=\beta \sqrt{n}=\sqrt{\frac{\pi }{2}n},即n=\frac{2k^2}{\pi }$

SetCount(X){
    k = 0; S = {}; a = uniform(X);
    do{
        k++;
        S = S U {a};
        a = uniform(X);
    }while(a不属于S)
    return 2k**2/pi
}

多重集合不同对象数目估计：将单词映射到长度为 m 的位串，$\pi (y,b)=\{\begin{array}{l}i,y[i]=b的最小i\\ k+1,y[i]!=b\end{array}$

WordCount(){
    y[1...m+1] = 0;
    for 磁带上的每个单词 x do {
        i = pi(h(x), 1); // x的散列值中等于1的最小位置，表示x是以 00...01....开头的
        y[i] = 1;
    }
    return pi(y, 0); // 返回y中等于0的最小位置
}
2^(k-2)<=n<=2^k，更准确的n=2^k/1.54703

2.2 Sherwood 算法

总是给出正确的答案，平滑不同输入实例的执行时间。

用途：平均性能较优，但是最坏性能较差时才使用。比如选择和排序，对元素的相对顺序敏感。

基本过程：

1. 将被解实例变换到一个随机实例。（预处理）
2. 用确定算法解随机实例，得到一个解。
3. 将此解变换为原实例的解。（后处理）

RH(x){
    // 用Sherwood算法计算f(x)
    n = length(x); // x的size为n
    r = uniform(A[n]); // 随机选取一个元素
    y = u(x,r); // 将原实例x转化为随机实例y
    s = f(y); // 用确定算法求y的解s
    return v(r,s); // 将s的解还原为x的解
}

选择和排序的 Sherwood 算法：将输入实例打乱，无需后处理。在调用确定性算法前，先洗牌

Shuffle(T){
    n = length(T);
    for i=1 to n-1 do {
        // 在T[1...n]中随机选取一个元素放在T[i]上
        j = uniform(1...n);
        sqap(T[i],T[j]);
    }
}

离散对数计算：$a=g^{x}modp,计算x=lo{g}_{g,p}a$

简单算法：

1. 计算$g^{x}modp,0\le x\le p-1$，x 有限个取值，因为离散对数$<g>$是个循环群
2. 验证$a=g^{x}modp$是否成立

dlog(g,a,p){ // 当这样的算法不存在时，算法返回p
    x = 0; y = 1;
    do{
        x++;
        y = y*g; // 计算y=g**x
    }while(a!=y%p && x!=p);
    return x;
}

定理：

• $lo{g}_{g,p}(st\%p)=(lo{g}_{g,p}s+lo{g}_{g,p}t)\%(p-1)$
• $lo{g}_{g,p}(g^r\%p)=r,0\le r\le p-2$

Sherwood 改造算法：

dlogRH(g, a, p){
    r = uniform(0..p-2);
    b = g^r mod p; // 求幂模b=g**r%p
    c = ba mod p; // 定理1
    y = dlog(g,c,p); // 使用确定性算法求log_{g,p}c, y=r+x
    return (y-r) mod (p-1); // 求x
}

$$>\begin{array}{rl}>y& =lo{g}_{g,p}(ba\%p)=(lo{g}_{g,p}b+lo{g}_{g,p}a)\%(p-1)\\ >& =(lo{g}_{g,p}(g^r\%p)+lo{g}_{g,p}a)\%(p-1)\\ >& =(r+x)\%(p-1)>\end{array}>$$

搜索有序表：数值数组和指针数组构成有序的静态链表。

时间为$O(n)$的确定算法：直接顺序搜索，平均时间$O(\frac{n}{2})$

Search(x, i){
    while x > val[i] do
     i = ptr[i];
     return i;
}
A(x){
    return Search(x, head);
}

时间为$O(n)$的概率算法：从随机的一个点左右开始顺序搜索，平均时间$O(\frac{n}{3})$

D(x){
    i = uniform(1..n);
    y = val[i];
    case{
        x<y: return Search(x, head);   // 从头开始搜索
        x>y: return Search(x, ptr[i]); // 从下一个位置开始搜索
        otherwise: return i; // case3 x==y
    }
}

平均时间为$O(\sqrt{n})$的确定算法：在前$\lfloor \sqrt{n}\rfloor$个数中找不大于 x 的最大整数的下标，从该数开始向后查找。

B(x){ // 设x在val[1..n]中
    i = head;
    max = val[i]; // max初值是表val中的最小值
    for j=1 to [sqrt(n)] do { // 在前根号n个数中找不大于x的最大整数y的下标
        y = val[j];
        if max < y <= x then {
            i = y; max = y;
        }
    }
    return Search(x, i);
}

2.3 Las Vegas 算法

特点：

• 更有效率的算法
• 时间上界可能不存在
• 可能找不到解，但如果找到一定是正确的
• 成功的概率随执行时间的增加而增加

算法的一般形式：$LV(x,y,success),x为输入,y为输出$

约定记号：

• $p(x),算法成功的概率$
• $s(x),算法成功时的期望时间$
• $e(x),算法失败时的期望时间$

一般会多次运行 LV 算法，知道成功

Obstinate(x){
    repeat
     LV(x,y,success);
    until success;
    return y;
}

找到正确解的期望时间：

$$>\begin{array}{rl}>t(x)& =p(x)s(x)+(1-p(x))(e(x)+t(x))\\ >& =s(x)+\frac{1-p(x)}{p(x)}e(x)>\end{array}>$$

2.3.1 八皇后问题

要求：

• 行不冲突：行号相等
• 列不冲突：列号相等
• 主对角线不冲突：行列号之差相等
• 辅对角线不冲突：行列号之和相等

传统回溯法：分别决定每一行放在哪个位置

i = j = 1;
while i ≤ 8 do { // 当前行号i≤ 8
    检查当前行i：从当前列j起向后逐列试探，寻找安全列号；
    if 找到安全列号 then {
        放置皇后，将列号记入栈，并将下一行置为当前行(i++);
        j=1; //当前列置为1
    } else {
        退栈回溯到上一行，即i--;
        移去该行已放置的皇后，以该皇后所在列的下一列作为当前列;
    }
}

Las Vegas 算法：

约定：

• $try[1\dots 8]表示第i个皇后放在(i,try[i])上$
• $try[1\dots k]称为k-promising$，当且仅当满足 k 皇后要求
• $\forall i!=j\in [1,k],try[i]-try[j]\notin \{i-j,0,j-i\}$

基本步骤：

1. 遍历改行所有的列位置是否可行，$\frac{1}{nb}$概率下选用该列
2. 如果有解则放，没解则退出

QueensLV(success){
    // 贪心的LV算法，所有皇后都是随机放置
    col, diag45, diag135 = {}; // 列、两个对角线初始化为空集
    k = 0; // 行号
    repeat
        nb = 0; // 计数器，nb值为(k+1)th皇后的open位置总数
     // 遍历列号，试探(k+1,i)是否安全
        for i=1 to 8 do {
            if (i!∈col) and (i-k-1!∈diag45) and (i+k+1!∈diag135) then {
                // 列i对(i+1)th皇后可用，但不马上放置
                nb = nb + 1;
                if uniform(1..nb)=1 then // 或许放在第i列
                 j = i;
            }
        }
        if(nb>0) then {
            // nb=0时无安全位置
            // 在所有nb个安全位置上，(k+1)th皇后选择位置j列的概率为1/nb
            k = k+1;
            try[k] = j;
            col = col ∪ {j};
            diag45 = diag45 ∪ {j-k};
            diag135 = diag135 ∪ {j+k};
        }
    until nb=0 or k=8; // 当找不到合适的位置，或者try是8-promising时，退出
    success = (nb>0);
}

算法 12 行，在所有可能的列位置上，选中某一列的概率都是$\frac{1}{nb}$，从最后一项往前倒推。

改进算法：联合回溯法和LV 算法，指定使用 LV 算法只确定前几行的位置。

...上述算法...
until nb=0 or k=stepVegas; // 已随机放好stepVegas个皇后
if nb>0 then
    backtrace(k, col, diag45, diag135, success);
else
    success = false;

2.3.2 模 p 平方根

$x\equiv y^2(modp),x\in [1,p-1],p为奇素数$，则

• x 为模 p 的二次剩余
• y 为 x 模 p 的平方根

定理：

• 任何一个模 p 的二次剩余有两个不同的平方根 $b,(p-b)$
• 1 到 p-1 之间的整数恰有一半时模 p 的二次剩余
• $\forall x\in [1,p-1],p为奇素数,有x^{(p-1)/2}\equiv \pm 1(modp)$
• $x是模p的二次剩余，当且仅当x^{(p-1)/2}\equiv 1(modp)$

如何判定x 是否为模 p 的二次剩余？计算$x^{\frac{p-1}{2}}(modp)$是否为 1。

如何计算x 模 p 的两个平方根？

• $ifp\equiv 3(mod4)$，则平方根为 $\pm x^{(p+1)/4}$

• $ifp\equiv 1(mod4)$，只能使用 LV 算法

例子：$p=53\equiv 1(mod4),x=7$，求 7 模 53 的平方根。

$$\begin{array}{rl}& \frac{p-1}{2}=26\\ & (2+\sqrt{7}{)}^{26}\equiv 0+41\sqrt{7}(mod53)\\ & 找到唯一的y,使得41y\equiv 1mod53,1\le y\le 52\\ & 使用Euclid算法,41\times 22\equiv 1(mod53)\\ & \therefore 22为一个平方根,53-22=31为另一个平方根\end{array}$$

基本步骤：

1. 随机取$a\in [1,p-1]$
2. 计算 c 和 d，$(a+\sqrt{x}{)}^{(p-1)/2}\equiv c+d\sqrt{x}(modp)$
3. 若$d=0$，则返回假
4. 若$c=0$，则找到 y，$dy\equiv 1(modp)$

rootLV(x, p, y, success){ // 计算y
    a = uniform(1..p-1); // 并不知道a取多少
    if a**2 % p == x % p then { // 可能性很小
        success = true; y =a;
    }else{
        计算c和d，使得0<=c,d<=p-1，(a+sqrt(x))**((p-1)/2) % p == c+d*sqrt(x) % p;
        if d==0 then success = false;
        else{ // c=0
            success = true;
            计算y使得，d*y % p == 1, 1<=y<=p-1 // 修改Euclid算法可求
        }
    }
}

算法成功的概率大于 0.5，因此平均调用两次即可求得 x 的平方根。

2.3.3 整数的因数分解

整数的素因数分解，$n=p_1^{m_1}{p}_2^{m_2}\cdots p_k^{m_k}$

朴素 split 算法：找到 n 的一个非平凡因数

split(n){
    // n是素数返回1，否则返回找到的n的最小非平凡因数
    for i=2 to [sqrt(n)] do
        if n % i == 0 then return i;
    return 1;
}

合数的二次剩余：上节的二次剩余要求 p 为奇素数，则恰有两个不同的平方根。若$n=pq,p和q为不同素数$，则二次剩余恰有4 个平方根。

Dixon 因数分解算法：

基本步骤：

1. 找到与 n 互素的整数 a 和 b，使得 $a^2\equiv b^2(modn)且a!\equiv \pm b(modn)$
2. 则 n 和 a+b 的最大公因子是一个非平凡因子，$x=gcd(n,a+b)$。

Dixon(n,x,success){ // 找合数n的某一非平凡因子
    if n是偶数 then {
        x = 2; success = true;
    }else{
        for i=2 to [log_{2}n] do {
            if n**(1/i) 是整数 then {
                x = n**(1/i); success = true; return;
            }
        }
        找到a, b，使得 a**2 % n == b**2 % n;
        if a % n == +-b % n then success = false;
        else {
            x = gcd(a+b, n); // 最大公因子
            success = true;
        }
    }
}

k-平滑：一个整数 x 的所有素因子均在前 k 个素数中。

确定 a 和 b 的取值：提前选定 n 和 k。

1. 随机选择$x\in [1,n-1]$。

   • 若 x 不与 n 互素，则找到了非平凡因子
   • 否则，$y=x^{2}modn$，若 y 是 k 平滑的，则将 x 和 y 的因数分解保存在表里。
   • 一直重复，直到选择了k+1 个互不相同的整数 x。

2. 在 k+1 个因式分解式中找到一个非空子集$\{y_1{y}_2\cdots y_k\}$，使得相应的因式分解的积中前 k 个素数的指数均为偶数（包含 0），如$y_1{y}_2{y}_4{y}_8=2^6{3}^4{5}^4{7}^{0}11^{2}13^{2}17^2$

3. 令$S=\{y_{s_1},y_{s_2},\cdots ,y_{s_k}\},a={\prod }_{y_i\in S}{x}_i,b=\sqrt{{\prod }_{y_i\in S}{y}_i}$，$ifa!\equiv \pm b(modn)$，则只需要求 a+b 和 n 的最大公因子即可。

时间分析：k 的取值将会影响 $x^{2}modn$是 k 平滑的可能性，以及因数分解 y 时的时间。通常取如下时间

$$\begin{gathered} L=e^{\sqrt{ln(n)\ast ln(ln(n))}} \\ k=\sqrt{L},期望时间O(L^2),成功概率至少0.5 \end{gathered}$$

2.4 Monte Carlo 算法

存在某些问题无法找到有效的算法获得正确解。MC 以高概率找到正确解。

约定：

• p-正确：一个 MC 算法以不小于 p 的概率返回一个正确的解
• 相容的、一致的：一个 MC 算法对同一实例返回相同的正确解。

多次调用同一个算法，选择出现次数最多的解，以增加一致的，p-正确算法成功的概率。

有偏算法：

• 偏真算法：MC 算法返回 true 时的解一定正确，只有返回 false 时才有可能产生错误的解。
• 偏$y_0$算法：算法返回$y_0$时总是正确的，返回非$y_0$时 p 概率正确。

重复调用一个一致的，p-正确的，偏$y_0$的 MC 算法 k 次，得到一个 $1-(1-p{)}^k$ -正确的算法。

2.4.1 主元素问题

一个具有 n 个元素的数组，若某一元素的个数大于$\frac{n}{2}$，则其为主元素。

判断 T 是否有主元素：随机设置一个 x，然后看看有多少个

maj(T){
    i = uniform(1..n);
    x = T[i];
    k = 0;
    for j=1 to n do
        if T[j] == x then
            k = k+1;
    return k>n/2;
}
偏真，0.5-正确的算法

重复调用降低错误率：

maj2(T){
    if maj(T) then return true; // 1次成功
    else return maj(T); // 调用2次
}
偏真，0.75-正确的算法

算法改进：控制算法出错概率小于$\epsilon >0$，应调用次数为$k=\lceil lg\frac{1}{\epsilon }\rceil$

// e 为错误概率
majMC(T, e){
    k = [lg(1/e)]; // 上取整
    for i=1 to k do
        if maj(T) then return true;
    return false;
}
// O(nlg(1/e))

2.4.2 素数测定

判定一个整数是否为素数，尚未找到有效的确定性算法或 LV 算法

简单的概率算法：

prime(n){
    d = uniform(2..[sqrt(n)]); // 下取整
    return (n mod d)!=0;
}

Fermat 小定理：$ifn是素数,则\forall a\in [1,n-1],有a^{n-1}modn=1$

变换命题：$if\exists a\in [1,n-1],使a^{n-1}modn!=1,则n不是素数$

素性测定算法：使用 Fermat 小定理，偏假的。即算法返回假，一定为假。

Fermat(n){
    a = uniform(1..n-1);
    if a**(n-1) mod n == 1 then return true; // 未必正确
    else return false; // 一定正确
}

伪素数和素性伪证据：符合 Fermat 小定理的合数，称为以 a 为底的伪素数，a 称为 n 的素性伪证据。

将 Fermat 测试重复多次，也不能降低误差。

强伪素数：n 是大于 4 的奇数，s 和 t 是整数，使得$n-1=2^{s}t,t为奇数$，集合$B(n)$满足如下条件：

$$\begin{array}{rl}& a\in B(n),2\le a\le n-2,且满足下述两个条件之一:\\ & 1.a^{t}modn=1\\ & 2.\exists 0\le i<s,a^{2^{i}t}modn=n-1\end{array}$$

• n 为素数，则 $\forall a\in [2,n-2],a\in B(n)$
• n 为合数，$ifa\in B(n)$，则称 n 为一个以 a 为底的强伪素数，a 为 n 素性的强伪证据。

Btest(a, n){
    // n为奇数，2<=a<=n-2, 返回true则代表a属于B(n)。则n为强伪素数或素数
    s=0; t=n-1; // t开始为偶数
    repeat
        s++;t/=2;
    until t mod 2 = 1; //n-1=(2**s)t, t为奇数
    x = a**t mod n;
    if x=1 or x=n-1 then return true; // 满足条件1或2
    for i=1 to s-1 do{
        x = x**2 mod n;
        if x=n-1 then return true; // 满足条件2
    }
    return false;
}
0.75-正确

强伪证据的数目比伪证据数目少很多。

$强伪证据\Rightarrow 伪证据$

Miller-Rabin 测试：0.75-正确，偏假的

MillRab(n){ // 奇n>4，返回真时表示素数，假时表示合数
    a = uniform(2..n-2);
    return Btest(a,n); // 测试n是否为强伪素数
}

RepeatMillRob：重复实验 k 次，$1-4^{-k}$-正确的 MC 算法

RepeatMillRob(n,k){
    for i=1 to k do
        if MillRob(n) == false then return false; // 一定是合数
    return true;
}

若给出出错概率不超过$\epsilon$，则重复次数为 $k=\lceil \frac{1}{2}lg\frac{1}{\epsilon }\rceil$

确定 n 素性时间为：$O(l{g}^{3}nlg\frac{1}{\epsilon })$

2.4.3 矩阵乘法验证

验证$矩阵AB=C$

设 $X_{1\times n}=(a_1,a_2,\cdots ,a_n),a_i\in \{0,1\}$，则问题转化为 $XAB=XC$

基本步骤：

1. 计算$XA$
2. 计算$XA$与$B$的乘积
3. 计算$XC$

GoodProduct(A,B,C,n){
    for i=1 to n do
        x[i] = uniform(0..1);
    if XAB=XC then return true;
    else return false;
}

偏假的，0.5-正确

返回假一定为假

改进算法：重复 k 次

RepeatGoodProduct(A,B,C,n,k){
    for i=1 to k do // 重复k次
        if GoodProduct(A,B,C,n) == false then return false;
    return true;
}

偏假的，$1-2^{-k}$-正确

给出错误概率：即$2^{-k}=\epsilon$

GP(A,B,C,n,e){
    k = [lg(1/e)]; // 上取整
    return RepeatGoodProduct(A,B,C,n,k);
}
// O(n**2log(1/e))