随机化数论算法总结

好吧这个名字很蠢（）

1 Miller_Rabin

作用试判断 \(10^{18}\) 级别或以上的数是否是质数，显然此时 \(O(\sqrt n)\) 的朴素算法是无法使用的。

1.1 费马小定理

当 \(p\) 为质数时，对于任意整数 \(a\)，有 \(a^{p-1}\equiv 1 \pmod p\)，证明如下：

引理 1

令 \(a,b,c\) 为任意整数， \(m\) 为正整数使 \(\gcd(m,c)=1\)，则当 \(a\times c\equiv b\times c\pmod m\) 时，有 \(a\equiv b\pmod m\)。

证明：\(a\times b\equiv b\times c\pmod m\) 可推出 \(a\times c-b\times c\equiv 0\pmod m\)，由于 \(\gcd(m,c)=1\)，可约去 \(c\)，即 \(a-b\equiv 0\pmod m\) 也就是 \(a\equiv b\pmod m\)

引理 2

令 \(m\) 为一个 \(>1\) 的正整数，\(b\) 为整数使 \(\gcd(m,b)=1\)，若 \(a_1,a_2,\ldots, a_m\) 为模式 \(m\) 的完全剩余系，则有 \(a_1\times b,a_2\times b,\ldots, a_m\times b\) 也为 \(m\) 的完全剩余系。

证明：若对于 \(\forall i,j\in [1,m]\) 有 \(a_i \times b\equiv a_j \times b \pmod m\)，则根据引理 1 有 \(a_i\equiv a_j \pmod m\)，与完全剩余系定义矛盾，故 \(a_1\times b,a_2\times b,\ldots, a_m\times b\) 也为 \(m\) 的完全剩余系。

于是构造对于 \(p\) 的完全剩余系：

\[P=\{1,2,3,\ldots, p-1\} \]

因为 \(\gcd(a,p)=1\)，由引理 2 得：

\[A=\{a,2a,3a,\ldots,(p-1)a\} \]

也为 \(p\) 的完全剩余系，由完全剩余系性质有：

\[1\times 2\times 3\times\ldots \times(p-1)\equiv a\times 2a\times 3a\times \ldots \times (p-1)a \pmod p \]

整理得：

\[(p-1)!\equiv (p-1)!\times a^{p-1}\pmod p \]

由于 \(\gcd ((p-1)!,p)=1\)，故有：

\[a^{p-1}\equiv 1\pmod p \]

证毕。

1.2 二次探测定理

当 \(p\) 为质数且 \(a^2\equiv 1 \pmod p\)，有 \(a\equiv \pm 1 \pmod p\)。

证明：由 \(a^2\equiv 1 \pmod p\) 得 \(a^2-1\equiv 0 \pmod p\)，即 \((a+1)(a-1)\equiv 0 \pmod p\)，则有 \(a\equiv \pm1\pmod p\)。

1.3 算法实现

若当前需判断的数为 \(p\)，先将 \(p-1\) 中的 \(2\) 都拆出拆成 \(2^k\times t\) 的形式，则当 \(p\) 为质数时有 \(a^{2^k\times t}\equiv 1\pmod p\) ，先算出 \(a^t\) 的值，再通过自乘弄 \(2^k\) 的系数，此时可以使用二次探测定理判断，自乘结束后再通过费马小定理判断一次，\(a\) 取一些小质数即可。

int pri[15]={0,2,3,5,7,11,13,17,19,23,29,31,37,41};
bool Miller_rabin(int x){
	if(x==2) return 1;
	int t=x-1,k=0;
	while(!(t&1)){
		t>>=1;
		k++;
	}
	for(int i=1;i<=13;i++){
		if(x==pri[i]) return 1;
		int a=qpow(pri[i]%x,t,x),nxt;
		for(int j=1;j<=k;j++){
			nxt=a*a%x;
			if(nxt==1&&a!=1&&a!=x-1) return 0;//二次探测定理判断
			a=nxt;
		}
		if(a!=1) return 0;//费马小定理判断
	}
	return 1;
}

2 Pollard-Rho 算法

可以解决 \(n\) 级别在 \(10^{18}\) 次方或更大的寻找任意质因子的问题

2.1 算法思想

考虑小学就听过的经典悖论——生日悖论，悖论的内容大体是对于一个生成值域在 \([1,n]\) 的随机数生成器，期望上会在生成 \(\sqrt{\frac{\pi n}{2}}\) 个数后得到两个相同的数。

类似应用到这个问题上就是设 \(n\) 的最小质因数为 \(p\)（显然有 \(p\le \sqrt n\)）最小随机生成 \(\sqrt [4]{n}\) 个数后就会有两个数 \(a,b\) 使得 \(a\equiv b \pmod p\)，然后我们求一下 \(\gcd(|a-b|,n)\) 就可以求出 \(p\) 了。

理想很美好，但当你用 c++ 原配的随机数生成后重新暴力找两个数匹配的复杂度又乘回去了,甚至因为要求 \(\gcd\) 的原因复杂度更是来到了 \(O(\sqrt n \log_2n)\)，沦落到不如暴力🤦‍♂️。

Pollard-Rho 算法给出了一个精细构造的伪随机数生成器：\(a_{i}=({a_{i-1}}^2+c)\bmod n\)，其中 \(a_1\) 和 \(c\) 是你用原配随机数生成的。

首先考虑这个序列的随机性，这里引用 OI-wiki 上的证明：

经过 OI-wiki 上的严谨证明，看来这个构造满足生日悖论的🎉。

作者许多资料也没有严谨的证明，所以或许看作一个正确性奇高的乱搞？

而且它还有一个性质：考虑若有 \(|a_i-a_j|\equiv 0\pmod p\)，那么同时就有：

\[\begin{aligned} &|a_{i+1}-a_{j+1}|\\ =\ &|(a_i^2+c)-(a_j^2+c)|\\ =\ &|a_i^2-a_j^2|\\ =\ &|(a_i+a_j)\textcolor{red}{(a_i-a_j)}|\equiv 0\pmod p \end{aligned} \]

看来我们只需要检查不同距离的 \(a_i\) 和 \(a_j\) 即可🧐。

这个序列的另外一个性质是由于 \(a_i\) 的生成只依赖 \(a_{i-1}\) 和 \(c\)，而后者又是固定的，所以当生成出重复的数字时序列会进入循环，也就是说我们需要判环。

2.2 算法实现

我们给出两个指针 \(A\) 和 \(B\)，让 \(A\) 每次动 1 步，\(B\) 每次动 2 步，这里的动几步相当于调用几次伪随机数生成器，那么 \(A\) 和 \(B\) 在“相遇”（相等）的时候就是有环了，此时如果还没有找到 \(|a_i-a_j|\equiv 0\pmod p\) 就退出换 \(c\)。

但是注意到这个东西的的复杂度期望有高额的 \(\sqrt [4]n \log_2n\)，注意到这个 \(\log_2\) 是 \(\gcd\) 给的，所以考虑每隔约 \(\log_2n\) 次再进行一次 \(\gcd\)，在此前把所有的差值乘起来，这样就消掉 \(\log\) 了。

还没完，这篇论文描述了一种加速方式不太好描述，下面代码中有。

论文中似乎是有精细证明该方法能比传统方法快 \(24\%\)，但是原论文

导致我看不了😭，有兴趣的自己看一下。

int Pollard_Rho(int n,int c){
	int a,b,i=1,k=2,val=1;
	a=ran()%(n-1)+1;
	b=a;
	while(1){
		a=f(a,c,n);
		val=val*abs(a-b)%n;
		if(i%127==0){//用来去 log 的东西
			int d=__gcd(val,n);
			if(d>1){
				return d;
			}
		}
		if(a==b){
			return n;
		}
		i++;
		if(i==k){//刚刚说的神秘优化
			k<<=1;
			b=a;
			int d=__gcd(val,n);
			if(d>1){
				return d;
			}
		}
	}
}

哦补充一点，求任意一个质因数没啥用，但这个东西配合 Miller Rabin 可以在极优的复杂度分解质因数，这应该是不难实现的就不贴代码了。