题意
求小于等于M的数中有多少个可以使得A^x=B(mod C)成立。
思路
首先利用extended BSGS求出最小的x使得A^x=B(mod C)成立。
接下来我们要求出最小的len使得A^(x+len*k)=B(mod C)，我们知道模C的指数循环节是Phi(C)，对于某个特定的幂的最小循环节则一定是Phi(C)的约数。我们要枚举Phi(C)的各个素因数，不断试除试算，知道得到最小的循环节长度len。
得到了x和len之后，答案就是(M-x)/len+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 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 |
#include <iostream> #include <algorithm> #include <queue> #include <vector> #include <cstdlib> #include <cstdio> #include <string> #include <cstring> #include <ctime> #include <iomanip> #include <cmath> #include <set> #include <stack> #include <cmath> #include <map> #include <complex> using namespace std; #pragma mark - Extended Baby-Step-Gaint-Step Algorithm (exBSGS) struct Hashtable{//Mod Based Hash Table static const int HASHMOD = 1200611;//A Prime long long top, hash[HASHMOD+100], value[HASHMOD+100], stk[1 << 16]; int locate(const long long x){ int h = x % HASHMOD; while(hash[h] != -1 && hash[h] != x){ h++; } return h; } void insert(const long long key,const long long val){ int h = locate(key); if(hash[h] == -1){ hash[h] = key; value[h] = val; stk[++top] = h; } } long long find(const long long key){ int h = locate(key); return (hash[h] == key)?value[h]:-1; } void clear(){ while(top){ hash[stk[top--]] = -1; } } Hashtable(){ top = 0; memset(hash, -1, sizeof(hash)); } } hashtab; struct Triple{ long long x,y,z; Triple(const long long _x,const long long _y,const long long _z){ x = _x; y = _y; z = _z; } Triple(){ x = y = z = 0; } }; Triple exgcd(const long long a,const long long b){ //for inv use(A/B % C):let a = B,b = C,return.x = inv,return.z = gcd //for solve use:return.x = specx,return.y = specy,returnz = gcd if (b == 0) { return Triple(1,0,a); } Triple last = exgcd(b, a%b); return Triple(last.y, last.x - a / b * last.y, last.z); } long long exBSGS(long long A,long long B,long long C){ //A^x=B(mod C) //Ensure B is legal B %= C; if (B < 0) { B += C; } //Procceed A,B,C to normal BSGS if C is not prime long long tmp = 1 % C, cnt = 0, D = 1 % C; for (int i = 0; i < 64; i++) { if (tmp == B) { return i; } tmp = tmp * A % C; } for (Triple res; (res = exgcd(A, C)).z != 1; cnt++) { if (B % res.z) { return -1; } B /= res.z; C /= res.z; D = D * A / res.z % C; } //Normal BSGS long long sqrtn = (long long)(ceil(sqrt(C))); hashtab.clear(); long long base = 1 % C; for (int i = 0; i < sqrtn; i++) { hashtab.insert(base, i); base = base * A % C; } long long j = -1,i = 0; for (; i < sqrtn; i++) { Triple res = exgcd(D, C); long long c = C / res.z; res.x = (res.x * B/res.z % c + c) % c; j = hashtab.find(res.x); if (j != -1) { return i * sqrtn + j + cnt; } D = D * base % C; } return -1; } #pragma mark - long long fastpow(long long a,long long p,long long m){ long long ret = 1; long long fac = a; while(p){ if(p&1){ ret *= fac; ret %= m; } p >>= 1; fac *= fac; fac %= m; } return ret; } long long euler_phi(long long x){ long long ret = 1; for(long long i = 2; i * i <= x; i++){ if(x % i == 0){ x /= i; ret *= i-1; while(x % i == 0){ x /= i; ret *= i; } } } if(x > 1){ ret *= x-1; } return ret; } long long loop_len(long long start,long long A,long long D,long long P){ long long phi = euler_phi(P),ret = phi; for(long long i = 2; i * i <= phi; i++){ while(phi % i == 0){ phi /= i; } while(ret % i == 0 && fastpow(A,start + ret/i,P) == D){ ret /= i; } } if (phi > 1) { while(ret % phi == 0 && fastpow(A,start + ret/phi,P) == D){ ret /= phi; } } return ret; } int main(){ long long A,P,D; unsigned long long M; while(scanf(" %lld %lld %lld %llu",&A,&P,&D,&M) != EOF){ A %= P,D %= P; long long first = exBSGS(A,D,P); unsigned long long res; long long phi = euler_phi(P); if(first == -1 || first > M){ res = 0; }else if(fastpow(A, first + phi, P) != D%P){ res = 1; }else{ unsigned long long len = loop_len(first,A,D,P); res = (M - (unsigned long long)first)/len + 1; } printf("%llu\n",res); } return 0; } |