luogu 3768 简单的数学题
题目描述
由于出题人懒得写背景了，题目还是简单一点好。
输入一个整数n和一个整数p，你需要求出$(\sum_{i=1}^n\sum_{j=1}^n ijgcd(i,j))~mod~p$，其中gcd(a,b)表示a与b的最大公约数。
刚才题面打错了，已修改
输入输出格式
输入格式：一行两个整数p、n。
输出格式：一行一个整数$(\sum_{i=1}^n\sum_{j=1}^n ijgcd(i,j))~mod~p$。
输入输出样例
998244353 2000
883968974
说明
对于20%的数据，$n \leq 1000$。
对于30%的数据，$n \leq 5000$。
对于60%的数据，$n \leq 10^6$，时限1s。
对于另外20%的数据，$n \leq 10^9$，时限3s。
对于最后20%的数据，$n \leq 10^{10}$，时限6s。
对于100%的数据，$5 \times 10^8 \leq p \leq 1.1 \times 10^9$且p为质数。
$$\sum_{d=1}^nd\sum_{i=1}^n\sum_{j=1}^nij[gcd(i,j)==d]$$
$$\sum_{d=1}^nd^3\sum_{i=1}^{n/d}\sum_{j=1}^{n/d}ij[gcd(i,j)==1]$$
莫比乌斯反演得
$$\sum_{d=1}^nd^3\sum_{i=1}^{n/d}\mu(i)i^2(1+2+...[\frac{n}{id}])^2$$
设$sum(x)=1+2+3+...x$
令$T=id$
$$\sum_{T=1}^nsum(\frac{n}{T})^2\sum_{d|T}d^3\frac{T}{d}^2\mu(\frac{T}{d}$$
$$\sum_{T=1}^nsum(\frac{n}{T})^2T^2\sum_{d|T}d\mu(\frac{T}{d})$$
但是后面的$$T^2\sum_{d|T}d\mu(\frac{T}{d})$$没法O(n)求前缀
这时候就要用到O(n$\frac{2}{3}$)的杜教筛
首先$$(id*\mu)(i)=\varphi(i)$$
$\mu$和$id(x)=x$的狄利克雷卷积等于$\varphi(i)$
即$$\sum_{d|T}d\mu(\frac{T}{d})=\varphi(i)$$
$$T^2\sum_{d|T}d\mu(\frac{T}{d})=T^2\varphi(T)$$
令$$f(i)=i^2\varphi(i)$$
$$S(n)=\sum_{i=1}^nf(i)$$
根据狄利克雷卷积：
\begin{aligned} \sum_{i = 1}^N (f * g)(i) & = \sum_{i = 1}^N \sum_{d \mid i} g(d) f\left(\frac i d\right) \\ & = \sum_{d = 1}^N g(d) \sum_{1 \leq i \leq N, d | i} f\left(\frac i d\right) \\ & = \sum_{d = 1}^N g(d) \sum_{1 \leq i \leq \lfloor \frac n d \rfloor} f(i) \\ & = \sum_{d = 1}^N g(d) S\left(\left\lfloor \frac n d \right\rfloor\right) \end{aligned}
于是有$$g(1)S(n)=\sum_{i=1}^n(g*f)(i)-\sum_{i=2}^ng(i)S(\frac{n}{i})$$
构造经可能能约掉d2的g(x),发现g(x)=x2时最好
因为$$\sum_{d|i}\varphi(d)=i$$
$$(g*f)(i)=\sum_{d|i}f(d)g(\frac{i}{d})=\sum_{d|i}d^2\varphi(d)\frac{i}{d}^2$$ $$=i^2\sum_{d|i}\varphi(d)=i^3$$
所以 $$S(n)=\sum_{i=1}^ni^3-\sum_{i=2}^ni^2S(\frac{n}{i})$$
先预处理出S(x)的前N位(N值自定，一般为maxn$\frac{2}{3}$)
否则就内部分块递归
外部直接数论分块
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cmath> 5 #include<algorithm> 6 #include<map> 7 using namespace std; 8 typedef long long lol; 9 const lol N=8000000; 10 lol Mod,tot,n; 11 lol ans,phi[N+5]; 12 lol prime[N+5],inv2,inv6,MAX; 13 bool vis[N+5]; 14 map<lol,lol> M; 15 lol qpow(lol x,lol y) 16 { 17 lol res=1; 18 while (y) 19 { 20 if (y&1) res=res*x%Mod; 21 x=x*x%Mod; 22 y/=2; 23 } 24 return res; 25 } 26 void pre() 27 {lol i,j; 28 phi[1]=1; 29 for (i=2;i<=MAX;i++) 30 { 31 if (vis[i]==0) 32 { 33 tot++; 34 prime[tot]=i; 35 phi[i]=(i-1)%Mod; 36 } 37 for (j=1;j<=tot;j++) 38 { 39 if (i*prime[j]>MAX) break; 40 vis[i*prime[j]]=1; 41 if (i%prime[j]==0) 42 { 43 phi[i*prime[j]]=phi[i]*prime[j]%Mod; 44 break; 45 } 46 else phi[i*prime[j]]=phi[i]*(prime[j]-1)%Mod; 47 } 48 } 49 for (i=1;i<=MAX;i++) 50 phi[i]=(i*i%Mod*phi[i]%Mod+phi[i-1])%Mod; 51 } 52 lol query(lol x) 53 { 54 if (x<=N) return phi[x]; 55 if (M[x]) return M[x]; 56 lol s=0,i,pos,sum,p=x%Mod; 57 s=p*(p+1)%Mod*inv2%Mod; 58 s=s*s%Mod; 59 lol as=0; 60 for (i=2;i<=x;i=pos+1) 61 { 62 pos=x/(x/i); 63 p=pos%Mod; 64 sum=(p+1)*p%Mod*(2*p+1)%Mod*inv6%Mod; 65 p=(i-1)%Mod; 66 sum-=(p+1)*p%Mod*(2*p+1)%Mod*inv6%Mod; 67 sum=(sum+Mod)%Mod; 68 sum=sum*query(x/i)%Mod; 69 as+=sum;as%=Mod; 70 } 71 return M[x]=(s-as+Mod)%Mod; 72 } 73 int main() 74 {lol i,pos,x,sum; 75 cin>>Mod>>n; 76 MAX=min(N,n); 77 pre(); 78 inv2=qpow(2,Mod-2); 79 inv6=qpow(6,Mod-2); 80 ans=0; 81 for (i=1;i<=n;i=pos+1) 82 { 83 pos=(n/(n/i)); 84 x=n/i;x%=Mod; 85 sum=(x+1)*x%Mod*inv2%Mod; 86 sum=sum*sum%Mod; 87 ans=ans+(sum*((query(pos)-query(i-1)+Mod)%Mod))%Mod; 88 ans%=Mod; 89 } 90 printf("%lld",ans); 91 }
