HDU 5525 Product 数论

题意:

给出一个长度为\(n(1 \leq n \leq 10^5)\)的序列\(A_i\)\(N=\prod\limits_{i=1}^{n}i^{A_i}\)。求\(N\)的所有约数的乘积。

分析:

首先还是将\(N\)质因数分解,\(N=\prod\limits_{i=1}^{k}p_i^{e_i}\)
考虑素因子\(p_a\)对结果的贡献,含有且仅含有\(p_a^x\)的约数有\(\frac{\prod\limits_{i=1}^{k}e_i+1}{e_a+1}\)个。
所以\(p_a\)对答案贡献了\(\frac{e_a(e_a+1)}{2} \cdot \frac{\prod\limits_{i=1}^{k}e_i+1}{e_a+1}\)次,也可以继续化简得到\(\frac{e_a\prod\limits_{i=1}^{k}e_i+1}{2}\)

指数太大,根据费马小定理,我们可以对\(MOD-1\)取模。
指数第一部分有个分母\(2\),由于不能计算\(MOD-1\)的逆元,所以可以先对\(2(MOD-1)\)取模。
第二部分可以通过计算前缀后缀积来避免除法。

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long LL;
const LL MOD = 1000000007;
const int maxn = 100000 + 10;

LL pow_mod(LL a, LL n) {
	LL ans = 1;
	while(n) {
		if(n & 1) ans = ans * a % MOD;
		a = a * a % MOD;
		n >>= 1;
	}
	return ans;
}

int pcnt, prime[maxn], pos[maxn];
bool vis[maxn];

void preprocess() {
	for(int i = 2; i < maxn; i++) {
		if(!vis[i]) prime[++pcnt] = i;
		for(int j = 1; j <= pcnt && i * prime[j] < maxn; j++) {
			vis[i * prime[j]] = true;
			if(i % prime[j] == 0) break;
		}
	}
	for(int i = 1; i <= pcnt; i++) pos[prime[i]] = i;
}

int n;
LL p[maxn];
LL pre[maxn], suf[maxn];

int main()
{
	preprocess();

	while(scanf("%d", &n) == 1) {
		memset(p, 0, sizeof(p));

		for(int i = 1; i <= n; i++) {
			int a; scanf("%d", &a);
			if(i == 1 || a == 0) continue;
			int t = i;
			for(int j = 1; j <= pcnt; j++) {
				if((LL) prime[j] * prime[j] > t) break;
				if(t % prime[j] == 0) {
					int cnt = 0;
					while(t % prime[j] == 0) {
						t /= prime[j];
						cnt++;
					}
					p[j] += (LL)cnt * a;
					p[j] %= 2 * (MOD - 1);
				}
			}
			if(t > 1) { p[pos[t]] += a; p[pos[t]] %= 2 * (MOD - 1); }
		}

		int tot = pcnt;
		while(tot > 1 && p[tot] == 0) tot--;
		pre[0] = suf[tot+1] = 1;
		for(int i = 1; i <= tot; i++) pre[i] = pre[i-1] * (p[i] + 1) % (MOD - 1);
		for(int i = tot; i >= 1; i--) suf[i] = suf[i+1] * (p[i] + 1) % (MOD - 1);

		LL ans = 1LL;
		for(int i = 1; i <= tot; i++) {
			LL ta, tb;
			if(p[i] % 2 == 0) ta = p[i]/2%(MOD-1)*((p[i]+1)%(MOD-1))%(MOD-1);
			else ta = (p[i]+1)/2%(MOD-1)*(p[i]%(MOD-1))%(MOD-1);
			tb = pre[i-1] * suf[i+1] % (MOD - 1);
			ta = ta * tb % (MOD - 1);
			ans = ans * pow_mod(prime[i], ta) % MOD;
		}

		printf("%lld\n", ans);
	}

	return 0;
}
posted @ 2015-11-03 17:14  AOQNRMGYXLMV  阅读(265)  评论(0编辑  收藏  举报