poj 2417: baby-step giant-step的模板题
hash用map实现会超时...所以我就手动模拟了一个模除hash...不会hash真受伤...
有一个纠结的地方是，我推出来应该是不用从0循环到m，只用到m-1就行了。但是run到m-1的代码挂了。自己测试了一下，发现“8 2 4”这组数据过不了。
算了一下，发现这组数据的特别在于它求出来的逆元是0. 所以我们需要特判逆元为0的情况。至于为什么run到m就行了，暂时还不清楚。
（这题是离散对数的简化版。等有时间了去真正的研究一下离散对数。）
1 //13440960 ooyyloo 2417 Accepted 1924K 63MS G++ 2107B 2014-09-14 16:27:24 2 #include <cstdio> 3 #include <iostream> 4 #include <cmath> 5 #include <algorithm> 6 #define pb push_back 7 #define mp make_pair 8 #define rep(i,n) for(int i=0;i<n;++i) 9 #define For(i,l,r) for(int i=l;i<=r;++i) 10 #define dbg(x) cout<<#x<<"="<<x<<endl 11 using namespace std; 12 typedef long long LL; 13 const int maxn = 100007; 14 LL ksm(LL a, LL n, LL p) 15 { 16 LL ret = 1; 17 for (a %= p; n; n >>= 1) 18 { 19 if (n & 1) ret = ret * a % p; 20 a = a * a % p; 21 } 22 return ret; 23 } 24 25 int ta[maxn], lin[maxn], pos; LL sd[maxn], ssd[maxn]; 26 void biu(LL s, LL t, LL i) { ++pos; lin[pos] = ta[s]; ta[s] = pos; sd[pos] = t; ssd[pos] = i; } 27 28 LL p, b, n; 29 LL ans; 30 void insert(LL tp, LL i) 31 { 32 for (LL i = ta[tp % maxn]; i; i=lin[i]) 33 if (sd[i] == tp) 34 { 35 ssd[i] = min(ssd[i], i); 36 return; 37 } 38 biu(tp % maxn, tp, i); 39 } 40 bool find(LL tp, LL &f) 41 { 42 for (LL i = ta[tp % maxn]; i; i=lin[i]) 43 if (sd[i] == tp) 44 { 45 f = ssd[i]; 46 return 1; 47 } 48 return 0; 49 } 50 void init() 51 { 52 pos = 0; 53 for (int i = 0; i < maxn; i++) ta[i] = 0; 54 } 55 int main() 56 { 57 for (LL bm, m, tp; scanf("%lld%lld%lld", &p, &b, &n) != EOF; ) 58 { 59 init(); 60 m = sqrt(p * 1.0); bm = ksm(b, m, p); 61 tp = 1; insert(tp, 0); 62 for (LL i = 1; i <= m; i++)//we must use "<=" here, just consider 8 2 4 63 tp = tp * b % p, insert(tp, i); 64 bm = ksm(ksm(b, m, p), p - 2, p); 65 tp = n % p; ans = 0x7fffffffll + 1ll; 66 for (LL i = 0, f; i <= m; i++) 67 { 68 if (find(tp, f)) ans = min(ans, i * m + f); 69 tp = tp * bm % p; 70 } 71 if (ans <= 0x7fffffffll) printf("%lld ", ans); 72 else puts("no solution"); 73 /*m = sqrt(p * 1.0); bm = ksm(b, m, p); 74 tp = 1; insert(tp, 0); 75 for (int i = 1; i <= m; i++) 76 tp = tp * bm % p, insert(tp, i); 77 bm = ksm(b, p-2, p); 78 tp = n * b % p; 79 ans = 0x7fffffffll + 1ll; 80 for (LL i = 0, f; i <= m; i++) 81 { 82 tp = tp * bm % p; 83 if (find(tp, f)) ans = min(ans, f * m + i); 84 } 85 if (ans <= 0x7fffffffll) printf("%lld ", ans); 86 else puts("no solution");*/ 87 } 88 return 0; 89 } 90 //when p==b, we can't get b^(-1)
poj 1284: 求原根个数
简单题。先求出2~65535里所有的原根，然后输入p，输出φ(p-1)就行了。
个数定理：设d|(p-1), 则模p的一个既约剩余系中对模p的指数等于d的元素个数等于φ(d).
另：模m有原根的充要条件是m=1, 2, 4, p^a, 2p^a. 其中，p为奇素数。
1 //13443437 ooyyloo 1284 Accepted 692K 16MS G++ 643B 2014-09-15 15:43:18 2 #include <cstdio> 3 using namespace std; 4 const int maxn = 65537; 5 6 int euler[maxn], plist[maxn], num; 7 void doeuler(int n = maxn - 1) 8 { 9 for (int i = 2; i <= n; i++) 10 { 11 if (!euler[i]) euler[i] = i - 1, plist[++num] = i; 12 for (int j = 1; j <= num && i * plist[j] <= n; j++) 13 { 14 if (i % plist[j] == 0) 15 { 16 euler[i * plist[j]] = euler[i] * plist[j]; 17 break; 18 } 19 euler[i * plist[j]] = euler[i] * (plist[j] - 1); 20 } 21 } 22 } 23 24 int main() 25 { 26 doeuler(); 27 for (int p; scanf("%d", &p) != EOF; ) printf("%d ", euler[p - 1]); 28 return 0; 29 }
poj 3696: 原根的应用
需要先公式推导出10^m % (9*L/gcd(8,9L)) == 1 这个结论。然后用欧拉函数和原根相关的定理求出最小的m。
要注意(9*L)^2会超long long，所以此题快速幂里的乘法应该用快速加来实现。
1 //13445399 ooyyloo 3696 Accepted 704K 360MS G++ 1897B 2014-09-16 00:58:50 2 #include <cstdio> 3 #include <algorithm> 4 using namespace std; 5 typedef long long LL; 6 const int maxn = 150001; 7 8 LL gcd(LL a, LL b) 9 { 10 LL temp; 11 while (b != 0) temp = b, b = a % b, a = temp; 12 return a; 13 } 14 LL mul(LL a, LL n, LL p) 15 { 16 LL ret = 0; 17 for (a %= p; n; n >>= 1) 18 { 19 if (n & 1) ret = (ret + a) % p; 20 a = (a + a) %p; 21 } 22 return ret; 23 } 24 LL ksm(LL a, LL n, LL p) 25 { 26 LL ret = 1; 27 for (a %= p; n; n >>= 1) 28 { 29 if (n & 1) ret = mul(ret, a, p); 30 a = mul(a, a, p); 31 } 32 return ret; 33 } 34 35 bool euler0[maxn]; 36 LL plist[maxn], num; 37 void doeuler(int n = maxn - 1) 38 { 39 for (int i = 2; i <= n; i++) 40 { 41 if (!euler0[i]) euler0[i] = 1, plist[++num] = i; 42 for (int j = 1; j <= num && i * plist[j] <= n; j++) 43 { 44 euler0[i * plist[j]] = 1; 45 if (i % plist[j] == 0) break; 46 } 47 } 48 } 49 50 LL zs[maxn], gs[maxn], cnt; 51 LL ys[maxn], po; 52 void div(LL d) 53 { 54 cnt = 0; 55 for (int i = 1; i <= num; i++) 56 if (d % plist[i] == 0) 57 { 58 zs[++cnt] = plist[i]; gs[cnt] = 0; 59 while (d % plist[i] == 0) gs[cnt]++, d /= plist[i]; 60 } 61 if (d > 1) zs[++cnt] = d, gs[cnt] = 1; 62 } 63 LL euler(LL d) 64 { 65 div(d); 66 LL ret = 1; 67 for (int i = 1; i <= cnt; i++) 68 { 69 ret *= (zs[i] - 1); 70 for (int j = 2; j <= gs[i]; j++) ret *= zs[i]; 71 } 72 return ret; 73 } 74 void dfs(int s, LL now) 75 { 76 if (s == cnt) return; 77 s++; 78 dfs(s, now); 79 for (int i = 1; i <= gs[s]; i++) 80 { 81 now *= zs[s]; 82 ys[++po] = now; 83 dfs(s, now); 84 } 85 } 86 87 LL L; 88 int main() 89 { 90 doeuler(); 91 for (int tt = 1; scanf("%lld", &L) != EOF; tt++) 92 { 93 if (!L) break; 94 LL tp = 9 * L / gcd(8, 9 * L); 95 div(euler(tp)); 96 ys[po = 1] = 1; dfs(0, 1); 97 sort(ys + 1, ys + 1 + po); 98 printf("Case %d: ", tt); 99 bool flag = 0; 100 for (int i = 1; i <= po; i++) 101 if (ksm(10, ys[i], tp) == 1) 102 { 103 printf("%lld ", ys[i]); 104 flag = 1; 105 break; 106 } 107 if (!flag) puts("0"); 108 } 109 return 0; 110 }
hdu 4992: 求模n的所有原根
首先，我们需要知道φ(n)的求法。我用的是筛法，把φ(n)用O(n)的时间预处理出来了。理论上直接用φ(n)=n*(1-1/p1)*……*(1-1/pm)也是可以的。
然后，需要知道原根存在定理。只有n=1, 2, 4, p^a, 2*p^a 的时候有原根。加上这个判断会快很多。
接着，就是正式的求解过程：
1. 求最小原根。因为最小原根一般都很小（不信自己输出来看），所以可以直接枚举求解。原根的判定应该用分解φ(n)求约数的方法做，否则速度会很慢。
2. 根据定理，若a为对模n的最小原根，那么a^x%n也是对模n的最小原根，其中(x,n)=1. 我们现在已经知道a，那么直接求解即可。
3. 排序输出答案。
第一步中原根判定也能用gφ(n)/pi != 1 (mod n) , pi | φ(n) 来做，而且据说这样做会快一些。经过试验，确实可行，从理论上来说，效率确实快一些，代码也少一些...
以上说明若有不懂之处，请找本数论书看看...
（某本书上说用随机法找原根的效率是多项式...有兴趣的尝试之后请联系窝...）
1 //11663376 2014-09-16 11:06:55 Accepted 4992 1015MS 6416K 2657 B G++ TofE 2 #include <cstdio> 3 #include <algorithm> 4 using namespace std; 5 typedef long long LL; 6 const int maxn = 1000001, maxx = 10000; 7 8 int gcd(int a, int b) 9 { 10 int temp; 11 while (b != 0) temp = b, b = a % b, a = temp; 12 return a; 13 } 14 int pow(int a, int n, int p) 15 { 16 int ret = 1; 17 for (a %= p; n; n >>= 1) 18 { 19 if (n & 1) ret = (LL)ret * a % p; 20 a = (LL)a * a % p; 21 } 22 return ret; 23 } 24 25 int euler[maxn], plist[maxn], num; 26 void doeuler(int n = maxn - 1) 27 { 28 euler[1] = 1; 29 for (int i = 2; i <= n; i++) 30 { 31 if (!euler[i]) euler[i] = i - 1, plist[++num] = i; 32 for (int j = 1; j <= num && i * plist[j] <= n; j++) 33 { 34 if (i % plist[j] == 0) 35 { 36 euler[i * plist[j]] = euler[i] * plist[j]; 37 break; 38 } 39 euler[i * plist[j]] = euler[i] * (plist[j] - 1); 40 } 41 } 42 } 43 44 int zs[maxx], gs[maxx], cnt; 45 int ys[maxx], yct; 46 void div(int d) 47 { 48 cnt = 0; 49 for (int i = 1; i <= num; i++) 50 if (d % plist[i] == 0) 51 { 52 zs[++cnt] = plist[i]; gs[cnt] = 0; 53 while (d % plist[i] == 0) gs[cnt]++, d /= plist[i]; 54 } 55 if (d > 1) zs[++cnt] = d, gs[cnt] = 1; 56 } 57 void dfs(int s, int now) 58 { 59 if (s == cnt) return; 60 s++; 61 dfs(s, now); 62 for (int i = 1; i <= gs[s]; i++) 63 { 64 now *= zs[s]; 65 ys[++yct] = now; 66 dfs(s, now); 67 } 68 } 69 70 int n; 71 int ans[maxn], po; 72 bool YES(int a) 73 { 74 for (int i = 1; i <= yct; i++) 75 if (pow(a, ys[i], n) == 1) return 0; 76 return pow(a, euler[n], n) == 1; 77 } 78 bool rich(int n) 79 { 80 if (n == 1 || n == 2 || n == 4) return 1; 81 div(n); 82 if (cnt > 2) return 0; 83 if (cnt == 2 && zs[1] != 2) return 0; 84 if (zs[1] == 2 && gs[1] > 1) return 0; 85 return 1; 86 } 87 int main() 88 { 89 doeuler(); 90 for (int a; scanf("%d", &n) != EOF; ) 91 { 92 if (n == 1) { puts("1"); continue; }//ke yi bu xie 93 if (!rich(n)) 94 { 95 puts("-1"); 96 continue; 97 } 98 div(euler[n]); 99 ys[yct = 1] = 1; dfs(0, 1); yct--; 100 for (int i = 1; i <= n; i++)//-- 101 if (YES(i)) { a = i; break; } 102 po = 0; 103 for (int i = 1, tp = a; i <= euler[n]; i++, tp = tp * a % n)// Is i == 0 possible ? 104 if (gcd(i, euler[n]) == 1) ans[++po] = tp; 105 sort(ans + 1, ans + 1 + po); 106 for (int i = 1, T = 0; i <= po; i++) printf("%s%d", T?" ":"", ans[i]), T = 1; 107 puts(""); 108 } 109 return 0; 110 }
sgu 261: 求 x^K = A (mod P) 的解
这道题是一道综合题。要用n多定理及模板。
解题思路：设r是对模P的原根，那么r^0, r^1, …… , r^(φ(P)-1)构成模P的既约剩余系。而P又是质数，所以它们就够成了模P的剩余系（除了0）。由上述定理知，可设r^j = A (mod P), (r^i)^K = x^K (mod P), 即 r^iK = r^j (mod P), 所以得到i*K = j (mod P-1)，然后解线性方程组求出i.
解题步骤：
1. 求出模P的最小原根r，再求出r^j = A (mod P).
2. 求出满足i*K = j (mod P-1) 的最小的i.
3. 枚举所有可能的i，求出对应的x，排序，输出。
ps. 注意到当A=0时，r^j = A (mod P)无解，所以这里要特判一下。
（我在写的时候，一开始枚举i的时候是for到P-1为止，WA了一次。后来发现，当最小的i为0时，i=P-1的答案和i=0的答案是重合的...改了就AC了）
1 //1593026 16.09.14 18:02 TofE 261 .CPP Accepted 93 ms 2022 kb 2 #include <cstdio> 3 #include <algorithm> 4 #include <cmath> 5 #include <map> 6 using namespace std; 7 typedef long long LL; 8 const int maxn = 40001; 9 10 int gcd(int a, int b) 11 { 12 int temp; 13 while (b != 0) temp = b, b = a % b, a = temp; 14 return a; 15 } 16 void _gcd(LL a, LL b, LL &x, LL &y) 17 { 18 if (b == 0) { x = 1, y = 0; return; } 19 _gcd(b, a % b, x, y); 20 LL temp = y; 21 y = x - a / b * y; 22 x = temp; 23 } 24 int pow(int a, int n, int p) 25 { 26 int ret = 1; 27 for (a %= p; n; n >>= 1) 28 { 29 if (n & 1) ret = (LL)ret * a % p; 30 a = (LL)a * a % p; 31 } 32 return ret; 33 } 34 35 int euler[maxn], plist[maxn], num; 36 void doeuler(int n = maxn - 1) 37 { 38 for (int i = 2; i <= n; i++) 39 { 40 if (!euler[i]) euler[i] = i - 1, plist[++num] = i; 41 for (int j = 1; j <= num && i * plist[j] <= n; j++) 42 if (i % plist[j] != 0) euler[i * plist[j]] = euler[i] * (plist[j] - 1); 43 else j = num + 1, euler[i * plist[j]] = euler[i] * plist[j]; 44 } 45 } 46 47 int zs[maxn], gs[maxn], cnt; 48 void div(int d) 49 { 50 cnt = 0; 51 for (int i = 1; i <= num; i++) 52 if (d % plist[i] == 0) 53 { 54 zs[++cnt] = plist[i]; gs[cnt] = 0; 55 while (d % plist[i] == 0) gs[cnt]++, d /= plist[i]; 56 } 57 if (d > 1) zs[++cnt] = d, gs[cnt] = 1; 58 } 59 60 int P, K, A; 61 int ans[maxn], act; //how many answers? 62 63 bool YES(int d) 64 { 65 for (int i = 1; i <= cnt; i++) 66 if (pow(d, (P - 1) / zs[i], P) == 1) return 0; 67 return 1; 68 } 69 int proot() 70 { 71 div(P - 1); 72 for (int i = 1; i <= P; i++) 73 if (YES(i)) return i; 74 } 75 int dlog(int r) 76 { 77 map<int, int> hashf; map<int, int>::iterator it; 78 int m = sqrt(P * 1.0), tp = 1; 79 for (int i = 0; i <= m; i++) 80 { 81 if (!hashf[tp]) hashf[tp] = i; 82 tp = (LL)tp * r % P; 83 } 84 tp = A; 85 int rm = pow(pow(r, m, P), P - 2, P); 86 for (int i = 0; i < m; i++) 87 { 88 if ((it = hashf.find(tp)) != hashf.end()) 89 return i * m + it->second; // only exist one answer. 90 tp = (LL)tp * rm % P; 91 } 92 } 93 void __gcd(int r, int j) 94 { 95 if (j % gcd(K, P - 1) != 0) return; 96 LL x, y, gbs = gcd(K, P - 1), a = K / gbs, b = (P - 1) / gbs; 97 _gcd(a, b, x, y); 98 x = (x * ((LL)j / gbs) % b + b) % b; 99 for (; x < P - 1; x += b) //-- 100 ans[++act] = pow(r, x, P); 101 } 102 int main() 103 { 104 doeuler(); 105 scanf("%d%d%d", &P, &K, &A); 106 if (A == 0) ans[++act] = 0; 107 else 108 { 109 int r = proot(); 110 int j = dlog(r); 111 __gcd(r, j); 112 } 113 sort(ans + 1, ans + 1 + act); 114 printf("%d ", act); 115 for (int i = 1; i <= act; i++) printf("%d ", ans[i]); 116 return 0; 117 }