题目:
分析:
首先用线性筛把\(A\)分解质因数,得到:
\[A=p_1^{a_1}*p_2^{a_2}...*p_n^{a_n} (p_i是质数且a_i>0) \] 则显然\(A^B\)分解质因数后为\[A=p_1^{a_1B}*p_2^{a_2B}...*p_n^{a_nB} (p_i是质数且a_i>0) \]接下来隆重推出约数和定理:(证明见)
\[Sum=\prod_{i=1}^n \sum_{j=0}^{a_i}p_i^j\] 那么很明显可以对于每一个\(p_i\)计算\(p_i^0+p_i^1...+p_1^{a_iB}\)然后乘起来就是答案。这就是一个等比数列求和了。 等比数列求和公式中含有除法,所以取模求和的时候不能直接用求和公式,否则如果除数刚好是模数的倍数就会出现逆元不存在的尴尬情况……例如POJ该题讨论区中的数据\(59407 \ 1\) (\(59407=9901*6+1\),求和公式中除数是\(59406\),此数在模\(9901\)意义下没有逆元)这里介绍一种二分等比数列求和的方法,思路和快速幂相似,即代码中的\(powersum\)函数
可以把这个等比数列平分成长度相等的两部分。 当\(n\)是偶数\[(p_i^0+p_i^1...+p_i^{n/2})+(p_i^{n/2+1}+p_i^{n/2+2}...+p_i^n)\] 然后从后半部分提出一个\(p_i^{n/2+1}\),它就和偶数部分一样了!得到\[(p_i^0+p_i^1...+p_i^{n/2})*(p_i^{n/2+1}+1)\] 显然左边可以递归地算下去,右边用快速幂求出。 当\(n\)是奇数,只要上述\(n/2\)均向下取整,算完以后加上\(p_i^n\)就可以了。这个也可以用快速幂解决。代码:
\(powersum\)函数求的是\(\sum_{j=1}^{a_i}p_i^j\),所以最后统计答案的时候要手动加上\(p_i^0\) (也就是\(1\))
筛质数时有一个小技巧。并不需要筛出\(5e7\)范围内的所有质数。x不可能含有两个或以上大于\(\sqrt x\)的质因数,所以\(x\)除以\(\sqrt{5e7}\)范围内的所有质数后如果仍不为\(1\),那么此时剩下的\(x\)一定是一个质数#includeusing namespace std;namespace zyt{ typedef long long ll; const int QQ_kotori = 9901; const int M = 7100; ll prime[M], index[M]; int cnt; void init_prime(ll x) { static bool mark[M]; for (int i = 2; i < M && x > 1; i++) { if (!mark[i]) { ll tmp = 0; prime[cnt] = i; while (x % i == 0) tmp++, x /= i; index[cnt++] = tmp; } for (int j = 0; j < cnt && i * prime[j] < M; j++) { ll k = i * prime[j]; mark[k] = true; if (i % prime[j] == 0) break; } } if (x > 1) { prime[cnt] = x; index[cnt++] = 1; } } ll power(ll a, ll b) { ll ans = 1; while (b) { if (b % 2) ans = ans * a % QQ_kotori; a = a * a % QQ_kotori; b /= 2; } return ans; } ll powersum(ll a, ll b) { if (b == 1) return a % QQ_kotori; ll ans = powersum(a, b / 2) * (1 + power(a, b / 2)) % QQ_kotori; if(b % 2) ans = (ans + power(a, b)) % QQ_kotori; return ans; } void work() { ll a, b, ans = 1; cin >> a >> b; if (a == 0) { cout << 0; return; } else if (b == 0) { cout << 1; return; } init_prime(a); for (int i = 0; i < cnt; i++) if (index[i]) ans = ans * (powersum(prime[i], index[i] * b) + 1) % QQ_kotori; cout << ans; }}int main(){ zyt::work(); return 0;}