hdu 4344 Mark the Rope pollard rho

缘起

最近（【1】）学习了 pollard rho 这种基于概率的快速分解超大整数素因数的算法. 但是【1】中还是没有显式的搞出所有的素因数，遂感觉不爽，于是找到 hdu 4344 Mark the Rope 来玩玩.

分析

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
Eric有一根长为N的绳子，他想在绳子上用不同的颜色做标记. 他标记颜色的方式是
1. 每次标记用的颜色不同
2. 每次标记会选择一个数L，L是N的非平凡因子（即 1<L<N, 且N能被L整除）
3. 每次标记，从0开头开始每隔L用1中选择的颜色标记一下
Eric希望选择K种标记，每种标记对应的L两两互质.Eric希望取得最大的K之后，求出L们的和S.

【输入】
多样例，第一行是样例数
后面每一行一个样例，每个样例就是一个数N（1<N<2^63）代表绳子的长度

【输出】
对每个样例，输出K和S，用空格分隔

【样例输入】
2
180
198

【样例输出】
3 18
3 22

【限制】
Time limit 10000 ms
Memory limit 32768 kB

其实和【1】是类似的. 唯一的区别就是【1】是求最小素因子，而这里是分解素因数之后，求所有素因数的和.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
//#include "stdafx.h"

#include <cstdio>
#include <cstdlib>
#include <ctime>
#include <map>
//#define LOCAL
using namespace std;

typedef long long LL;

map<LL, int> coutmap; // <k,v>键值对表示k这个素因子在分解素因数中出现了几次, 例如108=2^2*3^3, 所以2出现了2次,3出现了3次
LL n, A= 240,test=12; // A是pollard rho随机发生函数中的a, Miller-Rabin测试12次（毕竟范围已经到了整个long long）
typedef map<LL,int>::iterator mit;

LL modu_mul(LL a, LL b, LL mod) // a*b 除以mod的余数
{
	LL ans = 0;
	a%=mod, b%=mod;
	while(b)
	{
		if (b&1)
		{
			ans = (ans+a)%mod;
		}
		if (b>1)
		{
			a = (a<<1)%mod;
		}
		b>>=1;
	}
	return ans;
}

LL modu_pow(LL a, LL b, LL mod) // a^b 除以 mod的余数
{
	LL ans = 1;
	while (b)
	{
		if (b&1)
		{
			ans = modu_mul(ans, a, mod);
		}
		if (b>1)
		{
			a = modu_mul(a,a,mod);
		}
		b>>=1;
	}
	return ans;
}

bool mr(LL p) // Miller-Rabin测试p是不是素数
{
	if(p == 2 || p == 3 || p == 5 || p == 7 || p == 11)    return true;
	if(p == 1 || !(p%2) || !(p%3) || !(p%5) || !(p%7) || !(p%11))    return false;//先用无赖手段跑一波~
	if (p<2)
	{
		return false;
	}
	if (!(p&1))
	{
		return p==2;
	}
	LL s = 0, t = p-1;
	while(!(t&1)) s++,t>>=1;
	for (LL i = 0; i<test; i++)
	{
		LL a = rand()%(p-1)+1; // 随机选择mr测试的底数([1,p-1]之间的)
		LL x = modu_pow(a,t,p);
		for (LL j = 0; j<s; j++)
		{
			LL y = modu_mul(x,x,p);
			if (y==1 && x!=1 && x!=p-1) // Miller-Rabin测试
			{
				return false;
			}
			x=y;
		}
		if (x!=1) // 费马测试
		{
			return false;
		}
	}
	return true; // 通过12次Miller-Rabin测试，认定其为素数
}

LL gcd(LL a, LL b)
{
	return b?gcd(b,a%b):a;
}

LL pollard_rho(LL n, LL a) // pollard-rho算法, a是它的随机发生函数中的a
{
	LL i = 1,k=2;
	LL x = rand()%n; // x随机选取[0,n-1]中一个数作为起点即可
	LL y = x;
	while(1)
	{
		i++;
		x = (modu_mul(x,x,n)+a)%n; // f(x)=x^2+a mod n
		LL d = gcd(y-x+n, n);
		if (1<d&&d<n)
		{
			return d; // 找到了n的一个因子——未必是素因子
		}
		if (x==y)
		{
			return n;
		}
		if (i==k)
		{
			y = x;
			k<<=1;
		}
	}
}

void factorize(LL n, LL a) // 对n进行递归分解素因数
{
	if (mr(n)) // 找到n的一个素因子
	{
		coutmap[n]++;
		return;
	}
	LL d = n;
	while(d==n)
	{
		d = pollard_rho(n, a--);
	}
	factorize(d,a);
	factorize(n/d,a);
}


int main() 
{
#ifdef LOCAL
	freopen("d:\\data.in", "r", stdin);
	//freopen("d:\\my.out", "w", stdout);
#endif
	srand((unsigned int)(time(0)));
	LL t;
	scanf("%lld", &t);
	while(t--)
	{
		coutmap.clear();
		scanf("%lld", &n);
		factorize(n, A);
		if (coutmap.size()==1 && coutmap.begin()->first==n) // 如果n本身是素数
		{
			printf("1 %lld\n", n);
			continue;
		}
		if (coutmap.size()==1) // 如果n是单个素数因子的乘幂的话
		{
			printf("1 %lld\n",n/coutmap.begin()->first);
			continue;
		}
		LL ans = 0;
		for (mit i = coutmap.begin(); i!=coutmap.end(); i++)
		{
			ans += modu_pow(i->first, i->second, n);
		}
		printf("%d %lld\n", coutmap.size(), ans);
	}
	return 0;
}

ac情况（很奇怪，在hdu上ac，但是在vjudge上提交总是T）

30325811 2019-08-18 07:57:17 Accepted 4344 5382MS 1248K 3050 B G++ yfsyfsyfs

参考

【1】https://yfsyfs.github.io/2019/08/17/poj-1811-Prime-Test-Pollard-rho%E7%AE%97%E6%B3%95-%E5%88%86%E8%A7%A3%E7%B4%A0%E5%9B%A0%E6%95%B0-%E5%BF%AB%E9%80%9F%E5%B9%82/