「VMware校園挑戰賽」小V的和式
阿新 • • 發佈：2020-10-10
Description
給定 \(n,m\) ，求
\[\sum\limits_{x_1=1}^{n}\sum\limits_{x_2=1}^{n}\sum\limits_{y_1=1}^{m}\sum\limits_{y_2=1}^{m}\frac{(x_1y_2+x_2y_1)\ mod\ (max\{x_1,x_2\}max\{ y_1,y_2\})}{gcd(max\{x_1,x_2\},max\{ y_1,y_2\})} \]資料範圍 \(1\le n,m\le 10^9\)，答案對 \(998244353\) 取模。
Solution
約定
\(S_k(n)=\sum\limits_{i=1}^{n} i^k\)
若 \(n>m\) ，則交換 \(n,m\) ，保證 \(n\le m\)。
拆上下指標
先對上下指標做一下處理，得：
\[[1\le x_1\le n][1\le x_2\le n][1\le y_1\le m][1\le y_2\le m] \]\[=([x_1<x_2]+[x_1=x_2]+[x_1>x_2])([y_1<y_2]+[y_1=y_2]+[y_1>y_2]) \]\[=2[x_1<x_2][y_1<y_2]+2[x_1>x_2][y_1<y_2]+2([x_1<x_2][y_1=y_2]+[x_1=x_2][y_1<y_2])+[x_1=x_2][y_1=y_2] \]考慮拆開 \([x_1<x_2][y_1<y_2]\) ：
\[[1\le x_1<x_2][1\le y_1<y_2] \]\[=([0\le x_1<x_2]-[0=x_1<x_2])([0\le y_1<y_2]-[0=y_1<y_2]) \]迴帶到原來式子，得到：
\(2[0\le x_1<x_2][0\le y_1<y_2]\)
\(-2([0=x_1<x_2][0\le y_1<y_2]+[0\le x_1<x_2][0=y_1<y_2])\)
\(+2[0=x_1<x_2][0=y_1<y_2]\)
\(+2[x_1>x_2][y_1<y_2]\)
\(+2([x_1<x_2][y_1=y_2]+[x_1=x_2][y_1<y_2])\)
\(+[x_1=x_2][y_1=y_2]\)
所以分別求出這 \(6\) 種情形即可，下面開始推導。
1. \([0\le x_1<x_2][0\le y_1<y_2]\)
即求
\[f_1=\sum\limits_{x_1=0}^{n}\sum\limits_{x_2=x_1+1}^{n} \sum\limits_{y_1=0}^{m} \sum\limits_{y_2=y_1+1}^{m} \frac{(x_1y_2+x_2y_1)\ mod\ (x_2y_2)}{gcd(x_2,y_2)} \]\[=\sum\limits_{x_2=1}^{n} \sum\limits_{y_2=1}^{m} \frac{1}{gcd(x_2,y_2)} \sum\limits_{x_1=0}^{x_2-1}\sum\limits_{y_1=0}^{y_2-1} (x_1y_2+x_2y_1)\ mod\ (x_2y_2) \]- 引理1
對於所有的 \(a \in [0,y), b\in[0,x)\)，\((ax+by)\ mod\ (xy)\) 將集合 \(S=\{ (a,b) | 0\le a<y, 0\le b<x\}\) 對映到 \(T=\{ k\times gcd(a,b) | 0\le k< \frac{ab}{gcd(a,b)}\}\) ，並且對於 \(T\) 中的每一元素恰好有 \(gcd(a,b)\) 個原象。
感性理解一下即可，故：
\[f_1=\frac{1}{2}\sum\limits_{x_2=1}^{n} \sum\limits_{y_2=1}^{m} (\frac{x_2^2y_2^2}{gcd(x_2,y_2)}-x_2y_2) \]將前半部分拎出來：
\[\sum\limits_{x=1}^{n}\sum\limits_{y=1}^{m} \frac{x^2y^2}{gcd(x,y)} \]\[=\sum\limits_{d=1}^{n} d^3 \sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{m/d}i^2j^2 [gcd(i,j)=1] \]\[=\sum\limits_{d=1}^{n}d^3\sum\limits_{k=1}^{n/d} \mu(k)k^4 S_2(\frac{n}{dk})S_2(\frac{m}{dk}) \]\[=\sum\limits_{T=1}^{n} S_2(\frac{n}{T})S_2(\frac{m}{T})G_1(T) \]其中 \(G_1(n)=n^3 \sum\limits_{d|n} \mu(d)d\)
因此
\[f_1=\frac{1}{2}\sum\limits_{T=1}^{n} S_2(\frac{n}{T})S_2(\frac{m}{T})G_1(T)-\frac{n(n+1)m(m+1)}{8} \]2. \([0= x_1<x_2][0\le y_1<y_2]\)
即求
\[f_2=\sum\limits_{x_2=1}^{n} \sum\limits_{y_1=0}^{m}\sum\limits_{y_2=y_1+1}^{m} \frac{(x_2y_1)\ mod\ (x_2y_2)}{gcd(x_2,y_2)} \]\[=\sum\limits_{x_2=1}^{n} \sum\limits_{y_1=0}^{m}\sum\limits_{y_2=y_1+1}^{m} \frac{x_2y_1}{gcd(x_2,y_2)} \]其實還有對稱的 \([0\le x_1<x_2][0=y_1<y_2]\) ，推導類似。
3. \([0=x_1<x_2][0=y_1<y_2]\)
\[f_3=\sum\limits_{x_2=1}^{n}\sum\limits_{y_2=1}^{m}\frac{0}{gcd(n,m)}=0 \]4. \([x_1>x_2][y_1<y_2]\)
\[f_4=\sum\limits_{x_2=1}^{n}\sum\limits_{x_1=x_2+1}^{n}\sum\limits_{y_1=1}^{m}\sum\limits_{y_2=y_1+1}^{m} \frac{(x_1y_2+x_2y_1)\ mod\ (x_1y_2)}{gcd(x_1,y_2)} \]\[=\sum\limits_{x_1=1}^{n}\sum\limits_{y_2=1}^{m}\frac{1}{gcd(x_1,y_2)} \sum\limits_{x_2=1}^{x_1-1}\sum\limits_{y_1=1}^{y_2-1}x_2y_1 \]\[=\sum\limits_{x_1=1}^{n}\sum\limits_{y_2=1}^{m}\frac{S_1(x_1-1)S_1(y_2-1)}{gcd(x_1,y_2)} \]\[=\frac{1}{4}\sum\limits_{d=1}^{n}\frac{1}{d} \sum\limits_{x_1=1}^{n}\sum\limits_{y_2=1}^{m}x_1(x_1-1)y_2(y_2-1)[(x_1,y_2)=d] \]\[=\frac{1}{4}\sum\limits_{d=1}^{n} \frac{1}{d}\sum\limits_{x_1=1}^{n/d}\sum\limits_{y_2=1}^{m/d}(x_1^2y_2^2d^3-x_1y_2^2d^2-x_1^2y_2d^2+x_1y_2d)[gcd(x_1,y_2=1)] \]\[=\frac{1}{4}\sum\limits_{d=1}^{n}\sum\limits_{t=1}^{n/d}\mu(t)\sum\limits_{x_1=1}^{n/dt}\sum\limits_{y_2=1}^{m/dt}(x_1^2y_2^2d^3t^4-x_1y_2^2d^2t^3-x_1^2y_2d^2t^3+x_1y_1dt^2) \]\[=\frac{1}{4}\sum\limits_{T=1}^{n}(S_2(\frac{n}{T})S_2(\frac{m}{T})G_1(T)-(S_1(\frac{n}{T})S_2(\frac{m}{T})+S_2(\frac{n}{T})S_1(\frac{m}{T}))G_2(T)+S_1(\frac{n}{T})S_1(\frac{m}{T})G_3(T)) \]其中 \(G_2(n)=n^2\sum\limits_{d|n}\mu(d)d\)
\(G_3(n)=n\sum\limits_{d|n}\mu(d)d\)
5. \([x_1=x_2][y_1<y_2]\)
\[f_5=\sum\limits_{x_2=1}^{n}\sum\limits_{y_1=1}^{m}\sum\limits_{y_2=y_1+1}^{m} \frac{x_2(y_1+y_2)\ mod\ (x_2y_2)}{gcd(x2,y2)} \]\[=\sum\limits_{x_2=1}^{n}\sum\limits_{y_1=1}^{m}\sum\limits_{y_2=y_1+1}^{m} \frac{x_2y_1}{gcd(x_2,y_2)} \]可以發現 \(f_5\) 和 \(f_2\) 的式子完全一致，且符號相反，因此可以抵消。
顯然對於另一半對稱的式子也是如此。
6. \([x_1=x_2][y_1=y_2]\)
\[f_6=\sum\limits_{x_1=1}^{n}\sum\limits_{y_1=1}^{m} \frac{0}{gcd(x_2,y_2)}=0 \]合併
答案為
\[2f_1-2f_2+2f_3+2f_4+2f_5+f_6 \]\[=2f_1+2f_4 \]因此我們只需要快速求出 \(G_1(n),G_2(n),G_3(n)\) 的字首和即可。
杜教篩
考慮杜教篩，以 \(G_1(n)\) 為例。
\(f(n)=n^3\sum\limits_{d|n}\mu(d)d\)
記 \(S(n)=\sum\limits_{i=1}^{n}f(i)\)
則 \(S(n)=\sum\limits_{i=1}^{n}i^3\sum\limits_{d|i}\mu(d)d\)
\(=\sum\limits_{d=1}^{n}\mu(d)d^4S_3(\lfloor \frac{n}{d}\rfloor)\)
所以我們需要杜教篩求 \(f'(i)=\mu(i)i^4\) 的字首和。
考慮 \(g=ID^4\) ，則 \(h(n)=f'*g=[n==1]\)
故 \(S(n)=1-\sum\limits_{i=2}^{n}i^4 \times S(\lfloor \frac{n}{i}\rfloor)\)
實現
時間複雜度：\(O(n^{\frac{5}{6}})\)
空間複雜度：\(O(n^{\frac{2}{3}})\)
#pragma GCC optimize(2)
#pragma GCC optimize(3)
#pragma GCC optimize("Ofast")
#pragma GCC optimize("inline")
#pragma GCC optimize("-fgcse")
#pragma GCC optimize("-fgcse-lm")
#pragma GCC optimize("-fipa-sra")
#pragma GCC optimize("-ftree-pre")
#pragma GCC optimize("-ftree-vrp")
#pragma GCC optimize("-fpeephole2")
#pragma GCC optimize("-ffast-math")
#pragma GCC optimize("-fsched-spec")
#pragma GCC optimize("unroll-loops")
#pragma GCC optimize("-falign-jumps")
#pragma GCC optimize("-falign-loops")
#pragma GCC optimize("-falign-labels")
#pragma GCC optimize("-fdevirtualize")
#pragma GCC optimize("-fcaller-saves")
#pragma GCC optimize("-fcrossjumping")
#pragma GCC optimize("-fthread-jumps")
#pragma GCC optimize("-funroll-loops")
#pragma GCC optimize("-fwhole-program")
#pragma GCC optimize("-freorder-blocks")
#pragma GCC optimize("-fschedule-insns")
#pragma GCC optimize("inline-functions")
#pragma GCC optimize("-ftree-tail-merge")
#pragma GCC optimize("-fschedule-insns2")
#pragma GCC optimize("-fstrict-aliasing")
#pragma GCC optimize("-fstrict-overflow")
#pragma GCC optimize("-falign-functions")
#pragma GCC optimize("-fcse-skip-blocks")
#pragma GCC optimize("-fcse-follow-jumps")
#pragma GCC optimize("-fsched-interblock")
#pragma GCC optimize("-fpartial-inlining")
#pragma GCC optimize("no-stack-protector")
#pragma GCC optimize("-freorder-functions")
#pragma GCC optimize("-findirect-inlining")
#pragma GCC optimize("-fhoist-adjacent-loads")
#pragma GCC optimize("-frerun-cse-after-loop")
#pragma GCC optimize("inline-small-functions")
#pragma GCC optimize("-finline-small-functions")
#pragma GCC optimize("-ftree-switch-conversion")
#pragma GCC optimize("-foptimize-sibling-calls")
#pragma GCC optimize("-fexpensive-optimizations")
#pragma GCC optimize("-funsafe-loop-optimizations")
#pragma GCC optimize("inline-functions-called-once")
#pragma GCC optimize("-fdelete-null-pointer-checks")
#include<bits/stdc++.h>
using namespace std;
#define rint register int
#define rep(i,l,r) for(rint i=l;i<=r;i++)
#define per(i,l,r) for(rint i=l;i>=r;i--)
#define ll long long
#define ull unsigned long long
#define pii pair<int,int>
#define pll pair<ll,ll>
#define pb push_back
#define fir first
#define sec second
#define mset(s,t) memset(s,t,sizeof(s))
template<typename T1,typename T2>void ckmin(T1 &a,T2 b){if(a>b)a=b;}
template<typename T1,typename T2>void ckmax(T1 &a,T2 b){if(a<b)a=b;}
template<typename T>T gcd(T a,T b){return b?gcd(b,a%b):a;}
int read(){
int x=0,f=0;
char ch=getchar();
while(!isdigit(ch))f|=ch=='-',ch=getchar();
while(isdigit(ch))x=10*x+ch-'0',ch=getchar();
return f?-x:x;
}
const int N=10000005;
const int mod=998244353;
ll qpow(ll a,ll b=mod-2){
ll res=1;
a%=mod;
while(b>0){
if(b&1)res=res*a%mod;
a=a*a%mod;
b>>=1;
}
return res;
}
const ll inv2=qpow(2);
const ll inv4=qpow(4);
const ll inv6=qpow(6);
const ll inv8=qpow(8);
const ll inv30=qpow(30);
int mu[N],g1[N],g2[N],g3[N];
bool vis[N];
int pr[N>>3],len,L=1e7;
void init(int n){
mu[1]=1;
for(register int i=2;i<=n;i++){
if(!vis[i])pr[len++]=i,mu[i]=-1;
for(register int j=0;j<len&&pr[j]*i<=n;j++){
vis[pr[j]*i]=1;
if(i%pr[j]==0)break;
mu[pr[j]*i]=-mu[i];
}
}
rep(i,1,n){
g1[i]=(g1[i-1]+1ll*i*i%mod*i%mod*i%mod*mu[i])%mod;
g2[i]=(g2[i-1]+1ll*i*i%mod*i%mod*mu[i])%mod;
g3[i]=(g3[i-1]+1ll*i*i%mod*mu[i])%mod;
}
}
ll S1(ll x){
x%=mod;
return x*(x+1)%mod*inv2%mod;
}
ll S2(ll x){
x%=mod;
return x*(x+1)%mod*(2*x+1)%mod*inv6%mod;
}
ll S3(ll x){
x%=mod;
return S1(x)*S1(x)%mod;
}
ll S4(ll x){
x%=mod;
return x*(x+1)%mod*(2*x+1)%mod*(3*x*x%mod+3*x-1)%mod*inv30%mod;
}
int n,m;
unordered_map<int,ll>Map1,Map2,Map3;
ll GG1(int n){
if(n<=L)return g1[n];
if(Map1[n])return Map1[n];
ll ans=1;
for(int i=2,j;i<=n;i=j+1){
j=n/(n/i);
ans=(ans-GG1(n/i)*(S4(j)-S4(i-1)+mod))%mod;
}
return Map1[n]=(ans%mod+mod)%mod;
}
ll GG2(ll n){
if(n<=L)return g2[n];
if(Map2[n])return Map2[n];
ll ans=1;
for(int i=2,j;i<=n;i=j+1){
j=n/(n/i);
ans=(ans-GG2(n/i)*(S3(j)-S3(i-1)+mod))%mod;
}
return Map2[n]=(ans%mod+mod)%mod;
}
ll GG3(ll n){
if(n<=L)return g3[n];
if(Map3[n])return Map3[n];
ll ans=1;
for(int i=2,j;i<=n;i=j+1){
j=n/(n/i);
ans=(ans-GG3(n/i)*(S2(j)-S2(i-1)+mod))%mod;
}
return Map3[n]=(ans%mod+mod)%mod;
}
unordered_map<int,ll>map4,map5,map6;
ll G1(int n){
if(map4[n])return map4[n];
ll ans=0;
for(int i=1,j;i<=n;i=j+1){
j=n/(n/i);
ans=(ans+(GG1(j)-GG1(i-1)+mod)*S3(n/i))%mod;
}
return map4[n]=ans;
}
ll G2(int n){
if(map5[n])return map5[n];
ll ans=0;
for(int i=1,j;i<=n;i=j+1){
j=n/(n/i);
ans=(ans+(GG2(j)-GG2(i-1)+mod)*S2(n/i))%mod;
}
return map5[n]=ans;
}
ll G3(int n){
if(map6[n])return map6[n];
ll ans=0;
for(int i=1,j;i<=n;i=j+1){
j=n/(n/i);
ans=(ans+(GG3(j)-GG3(i-1)+mod)*S1(n/i))%mod;
}
return map6[n]=ans;
}
ll f1(){
ll ans=0;
for(int i=1,j;i<=n;i=j+1){
j=min(n/(n/i),m/(m/i));
ans=(ans+S2(n/i)*S2(m/i)%mod*(G1(j)-G1(i-1)+mod))%mod;
}
ans=ans*inv2%mod;
ans=(ans+mod-S1(n)*S1(m)%mod*inv2%mod)%mod;
return (ans+mod)%mod;
}
ll f4(){
ll ans=0;
for(int i=1,j;i<=n;i=j+1){
j=min(n/(n/i),m/(m/i));
ans=(ans+S2(n/i)*S2(m/i)%mod*(G1(j)-G1(i-1)+mod))%mod;
ans=(ans-(S1(n/i)*S2(m/i)+S2(n/i)*S1(m/i))%mod*(G2(j)-G2(i-1)+mod))%mod;
ans=(ans+S1(n/i)*S1(m/i)%mod*(G3(j)-G3(i-1)+mod))%mod;
}
ans=ans*inv4%mod;
return (ans+mod)%mod;
}
int main(){
init(L);
//n=1e9,m=1e9;
scanf("%d%d",&n,&m);
if(n>m)swap(n,m);
printf("%lld\n",2ll*(f1()+f4())%mod);
return 0;
}
關於對 \([0= x_1<x_2][0\le y_1<y_2]\) 的彩(推)蛋(導)
qwq 其實是因為我一開始沒意料到能相互抵消，所以推了不少，這裡就保留一下好了。
\[=\sum\limits_{d=1}^{n} \sum\limits_{x_2=1}^{n/d}\sum\limits_{y_2=1}^{m/d}x_2S_1(y_2d-1)[gcd(x_2,y_2)=1] \]\[=\sum\limits_{d=1}^{n}\sum\limits_{k=1}^{n/d}\mu(k)kS_1(\frac{n}{dk})\sum\limits_{y_2=1}^{m/dk}S_1(y_2dk-1) \]右邊 \(y_2\) 的式子，令 \(t=m/dk\) ，經化簡得
\[\frac{1}{2}(d^2k^2S_2(t)-dkS_1(t)) \]因此：
\[f_2=\frac{1}{2}\sum\limits_{d=1}^{n}\sum\limits_{k=1}^{n/d}\mu(k)kS_1(\frac{n}{dk})(d^2k^2S_2(\frac{m}{dk})-dkS_1(\frac{m}{dk})) \]然後令 \(T=dk\) 化簡一下就行了。。。