hdu5228 Count a*b

链接

http://acm.hdu.edu.cn/showproblem.php?pid=5528

题解

今天做完这道题才发现自己数论是如此薄弱
通过观察可以发现，每一行的$0$的个数恰好等于$gcd(i,m)$
打个表发现也确实是如此
于是题目就是要求$$ans=\sum _{m|n}{m}^2-\sum _{m|n}\sum _{i=0}^{m-1}gcd(i,m)$$
下面开始推式子
$$\begin{gathered} \sum _{m|n}\sum _{i=1}^{m}gcd(i,m) \\ =\sum _{m|n}\sum _{d|m}^{m}d\sum _{d|i}^m[gcd(i,m)=d] \\ =\sum _{m|n}\sum _{d|m}^{m}d\phi (\frac{m}{d}) \\ =\sum _{d|n}d\sum _{d|m,m|d}\phi (\frac{m}{d}) \\ =\sum _{d|n}d\sum _{\frac{m}{d}=k|\frac{n}{d}}\phi (k) \\ =\sum _{d|n}d\frac{n}{d} \\ =\sum _{d|n}n \\ =n{\sigma }_0(n) \end{gathered}$$
所以答案就是
$${\sigma }_2(n)-n{\sigma }_0(n)$$
其中${\sigma }_k(n)$表示$n$的约数的$k$次方和
最坑爹的在于，这道题$O(\sqrt{n})$的暴力求会被卡常数
所以用分解质因数来做
设$n=\prod p_i^{q_i}$
$${\sigma }_k(n)=\prod (p_i^0+p_i^1+...+p_i^{q_i})$$

代码

#include <bits/stdc++.h>
#define maxn 100010
#define ull unsigned long long
#define sqr(x) ((x)*(x))
using namespace std;
ull f[maxn], d[maxn], prime[maxn], p[maxn], q[maxn];
bool mark[maxn];
void init()
{
	ull i, j, t, cnt;
	f[1]=1, d[1]=1;
	for(i=2;i<maxn;i++)
	{
		if(!mark[i])prime[++*prime]=i, f[i]=i*i+1, d[i]=2;
		for(j=1;j<=*prime and i*prime[j]<maxn;j++)
		{
			mark[i*prime[j]]=1;
			if(i%prime[j]==0)
			{
				f[i*prime[j]]=f[i]*prime[j]*prime[j]+1;
				for(t=i*prime[j],cnt=0;t%prime[j]==0;t/=prime[j],cnt++);
				d[i*prime[j]]=(cnt+1)*d[t];
				break;
			}
			else d[i*prime[j]]=d[i]*2, f[i*prime[j]]=f[i]*f[prime[j]];
		}
	}
}
ull calc(ull n)
{
	ull i, j, ans1=1, ans2=1, t, tot=0, tmp, N=n;
	for(i=1;sqr(prime[i])<=n;i++)
		if(n%prime[i]==0)for(tot++,p[tot]=prime[i],q[tot]=0;n%prime[i]==0;n/=prime[i],q[tot]++);
	if(n>1)p[++tot]=n, q[tot]=1;
	for(i=1;i<=tot;i++)
	{
		for(j=0,t=1,tmp=0;j<=q[i];j++,t*=sqr(p[i]))tmp+=t;
		ans1*=tmp;
		ans2*=j;
	}
	return ans1-N*ans2;
}
int main()
{
	ull T, n, i;
	scanf("%llu",&T);
	init();
	while(T--)
	{
		scanf("%llu",&n);
		printf("%llu\n",calc(n));
	}
	return 0;
}