时光匆匆,转眼间又是一年寒暑……
这是小 Q 同学第二次参加省队选拔赛。
今年,小 Q 痛定思痛,不再冒险偷取试题,而是通过练习旧 试题提升个人实力。可是旧试题太多了,小 Q 没日没夜地做题,却看不到前方的光明在哪里。
一天,因做题过度而疲惫入睡的小 Q 梦到自己在考场上遇到了一道好像做过的题目,却怎么也想不 起曾经自己是怎么解决它的,直到醒来还心有余悸。
小 Q 眉头一皱,感觉事情不妙,于是他找到了你,希望你能教他解决这道题目。小 Q 依稀记得题目 要计算如下表达式的值
$({sum_{i=1}^{A}}{sum_{j=1}^{B}}{sum_{k=1}^{C}}d(ijk))mod(10^9+7)$
其中 d(ijk) 表示 ijk 的约数个数。
Input
第一行包含一个正整数 T,表示有 T 组测试数据。 接下来 T 行,每行描述一组测试数据,包含三个整数 A, B 和 C,含义见题目描述。
Output
对于每组测试数据,输出一行,包含一个整数,表示所求表达式的值。
解题思路:
极限卡常题,具体卡常的方法主要是使用vector代替链式前向星,如果还过不去,还是看看dalao们的博客吧。
这里主要讲算法。
首先是喜闻乐见的莫比乌斯反演,三个$sum$有点中暑,慢慢推吧其实都差不多。
不妨设$Aleq Bleq C$,则:
$Ans={sum_{i=1}^{A}}{sum_{j=1}^{B}}{sum_{k=1}^{C}}d(ijk)$
$={sum_{i=1}^{A}}{sum_{x|i}}{sum_{j=1}^{B}}{sum_{y|j}}{sum_{k=1}^{C}}{sum_{z|k}}{varepsilon(gcd(x,y))}{varepsilon(gcd(y,z))}{varepsilon(gcd(x,z))}$
$={sum_{x=1}^{A}}{sum_{x|i}^{A}}{sum_{y=1}^{B}}{sum_{y|j}^{B}}{sum_{z=1}^{C}}{sum_{z|k}^{C}}{varepsilon(gcd(x,y))}{varepsilon(gcd(y,z))}{varepsilon(gcd(x,z))}$
$={sum_{x=1}^{A}}{sum_{y=1}^{B}}{sum_{z=1}^{C}}{leftlfloor{frac{A}{x}} ight floor}{leftlfloor{frac{B}{y}} ight floor}{leftlfloor{frac{C}{z}} ight floor}{varepsilon(gcd(x,y))}{varepsilon(gcd(y,z))}{varepsilon(gcd(x,z))}$
$={sum_{x=1}^{A}}{sum_{y=1}^{B}}{sum_{z=1}^{C}}{leftlfloor{frac{A}{x}} ight floor}{leftlfloor{frac{B}{y}} ight floor}{leftlfloor{frac{C}{z}} ight floor}{sum_{i|gcd(x,y)}}{sum_{j|gcd(y,z)}}{sum_{k|gcd(x,z)}}{mu(i))}{mu(j))}{mu(k))}$
$={sum_{i=1}^{A}}{sum_{j=1}^{B}}{sum_{k=1}^{A}}{mu(x)}{mu(y)}{mu(z)}{sum_{x|lcm(i,k)}}{sum_{y|lcm(i,j)}}{sum_{z|lcm(j,k)}}{left lfloor {frac{A}{i}} ight floor}{left lfloor {frac{B}{j}} ight floor}{left lfloor {frac{C}{k}} ight floor}$
将$mu$不等于0的连边,双向改单向枚举三元环就好了。
代码:
1 #include<cstdio> 2 #include<vector> 3 #include<cstring> 4 #include<algorithm> 5 typedef long long lnt; 6 const int N=100010; 7 const lnt mod=(lnt)(1e9+7); 8 struct edge{ 9 int from; 10 int to; 11 int lcm; 12 }ed[N*100]; 13 struct ent{ 14 int twd; 15 int lcm; 16 }; 17 int prime[N]; 18 int miu[N]; 19 int ind[N]; 20 bool vis[N]; 21 lnt Lcm[N]; 22 lnt F[3][N]; 23 int cnt; 24 std::vector<ent>hd[N]; 25 int A,B,C; 26 void gtp(void) 27 { 28 miu[1]=1; 29 for(int i=2;i<N;i++) 30 { 31 if(!vis[i]) 32 { 33 prime[++cnt]=i; 34 miu[i]=-1; 35 } 36 for(int j=1;j<=cnt&&prime[j]*i<N;j++) 37 { 38 vis[i*prime[j]]=true; 39 if(i%prime[j]==0) 40 break; 41 miu[i*prime[j]]=-miu[i]; 42 } 43 } 44 return ; 45 } 46 lnt gcd(lnt x,lnt y) 47 { 48 if(!y) 49 return x; 50 return gcd(y,x%y); 51 } 52 int main() 53 { 54 int T; 55 gtp(); 56 scanf("%d",&T); 57 while(T--) 58 { 59 int a,b,c; 60 scanf("%d%d%d",&a,&b,&c); 61 A=std::min(a,std::min(b,c)); 62 C=std::max(a,std::max(b,c)); 63 B=a+b+c-A-C; 64 cnt=0; 65 for(int i=1;i<=C;i++) 66 hd[i].clear(),ind[i]=0; 67 for(int i=1;i<=C;i++) 68 F[0][i]=F[1][i]=F[2][i]=0; 69 for(int i=1;i<=A;i++) 70 for(int j=i;j<=A;j+=i) 71 F[0][i]+=A/j; 72 for(int i=1;i<=B;i++) 73 for(int j=i;j<=B;j+=i) 74 F[1][i]+=B/j; 75 for(int i=1;i<=C;i++) 76 for(int j=i;j<=C;j+=i) 77 F[2][i]+=C/j; 78 lnt ans=0; 79 for(int i=1;i<=A;i++) 80 if(miu[i]) 81 ans+=miu[i]*F[0][i]*F[1][i]*F[2][i]; 82 ans%=mod; 83 for(int d=1;d<=C;d++) 84 { 85 for(int i=1;i*d<=C;i++) 86 { 87 if(!miu[i*d]) 88 continue; 89 for(int j=i+1;1ll*j*i*d<=C;j++) 90 { 91 if(miu[j*d]==0||gcd(j,i)!=1) 92 continue; 93 int x=i*d,y=j*d,lcm=i*j*d; 94 ans+=miu[x]*(F[0][y]*F[1][lcm]*F[2][lcm]+F[1][y]*F[0][lcm]*F[2][lcm]+F[2][y]*F[1][lcm]*F[0][lcm]); 95 ans=ans%mod; 96 ans+=miu[y]*(F[0][x]*F[1][lcm]*F[2][lcm]+F[1][x]*F[0][lcm]*F[2][lcm]+F[2][x]*F[1][lcm]*F[0][lcm]); 97 ans=ans%=mod; 98 ind[x]++; 99 ind[y]++; 100 ed[++cnt]=(edge){x,y,lcm}; 101 } 102 } 103 } 104 for(int i=1;i<=cnt;i++) 105 { 106 int f=ed[i].from,t=ed[i].to,v=ed[i].lcm; 107 if(ind[f]<ind[t]||(ind[f]==ind[t]&&f<t)) 108 hd[f].push_back((ent){t,v}); 109 else 110 hd[t].push_back((ent){f,v}); 111 } 112 for(int i=1;i<=C;i++) 113 { 114 for(int j=0;j<hd[i].size();j++) 115 Lcm[hd[i][j].twd]=hd[i][j].lcm; 116 for(int j=0;j<hd[i].size();j++) 117 { 118 int a=i,b=hd[i][j].twd; 119 for(int k=0;k<hd[b].size();k++) 120 { 121 int c=hd[b][k].twd;; 122 int lab=hd[a][j].lcm; 123 int lbc=hd[b][k].lcm; 124 int lac=Lcm[c]; 125 if(!lac) 126 continue; 127 lnt tmp=miu[a]*miu[b]*miu[c]; 128 lnt ttt=F[0][lab]*F[1][lbc]*F[2][lac]+ 129 F[0][lab]*F[1][lac]*F[2][lbc]+ 130 F[0][lbc]*F[1][lac]*F[2][lab]+ 131 F[0][lbc]*F[1][lab]*F[2][lac]+ 132 F[0][lac]*F[1][lab]*F[2][lbc]+ 133 F[0][lac]*F[1][lbc]*F[2][lab]; 134 ans=(ans+ttt*tmp%mod)%mod; 135 } 136 } 137 for(int j=0;j<hd[i].size();j++) 138 Lcm[hd[i][j].twd]=0; 139 } 140 printf("%lld ",(ans%mod+mod)%mod); 141 } 142 return 0; 143 }