[bzoj5332] [SDOI2018]旧试题

题目描述

时光匆匆，转眼间又是一年省选季……

这是小 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) 表示 i × j × k 的约数个数。

输入输出格式

输入格式：

第一行包含一个正整数 T，表示有 T 组测试数据。

接下来 T 行，每行描述一组测试数据，包含三个整数 A,B 和 C，含义见题目描述。

输出格式：

对于每组测试数据，输出一行，包含一个整数，表示所求表达式的值。

输入输出样例

输入样例：

5
10 10 10
100 100 100
1000 1000 1000
10000 10000 10000
100000 100000 100000

输出样例：

11536
51103588
165949340
19234764
176764584

Solution

关于约数个数(d(x))有这样一个性质：

[d(ijk)=sum_{r|i}sum_{s|j}sum_{t|k}[(r,s)=1][(s,t)=1][(r,t)=1] ]

证明：对于每一个质因子进行考虑，设(i,j,k)当前质因子的指数分别是(a,b,c)，显然这个质因子的贡献是((a+b+c+1))，那么由于限制了(r,s,t)两两互质，最多这三者只有一个含有当前质因子：

1. (r)含有(x)个，则表示当前质因子用了(x)个
2. (s)含有(x)个，则表示用了(a+x)个
3. (t)含有(x)个，则表示用了(a+b+x)个
4. 若都不含当前质因子，则表示没用这个质因子。

所以总情况是(a+b+c+1)个，这样正好就可以算对。

根据上面的证明，推广到(n)个数相乘也是成立的，只要枚举的数两两互质就好。

所以带进原式可得：

[ans=sum_{i=1}^asum_{j=1}^bsum_{k=1}^csum_{r|i}sum_{s|j}sum_{t|k}epsilon(r,s)epsilon(s,t)epsilon(r,t) ]

这里将([gcd(i,j)=1])简记成(epsilon(i,j))。

然后把(r,s,t)提前：

[ans=sum_{r=1}^{a}sum_{s=1}^bsum_{t=1}^clfloorfrac{a}{r} floorlfloorfrac{b}{s} floorlfloorfrac{c}{t} floorepsilon(r,s)epsilon(s,t)epsilon(r,t) ]

然后对于每一个(epsilon(i,j))都进行莫比乌斯反演，整理化简可得：

[ans=sum_{u=1}^{min(a,b)}sum_{v=1}^{min(b,c)}sum_{w=1}^{min(a,c)}mu(u)mu(v)mu(w)s(frac{a}{{ m{lcm}}(u,v)})s(frac{b}{{ m{lcm}}(v,w)})s(frac{c}{{ m{lcm}}(u,w)}) ]

其中：

[s(n)=sum_{i=1}^{n}lfloorfrac{n}{i} floor=sum_{i=1}^{n}d(i) ]

这一步由约数个数的定义可得。

然后，，上面那个长的像(O(n^3))的东西就是正解，，

考虑下如何让式子后面那一块有值，当且仅当：

[mu(u) e 0,mu(v) e 0,mu(w) e 0,{{ m{lcm}}(u,v)}leqslant a,{{ m{lcm}}(v,w)}leqslant b,{{ m{lcm}}(u,w)}leqslant c ]

所以对于每一个质数，在(u,v,w)的质因数分解中只可能至多出现一次，所以我们爆搜有没有出现的八种情况，然后剪枝。。

• 可以预处理出前面一部分的答案，当当前质数大于限制的时候，且范围较小，可以直接用预处理出来的答案。
• 考虑到超过两个含有当前质数时，除开当前的(mu)的影响，所有情况答案都是一样的，所以可以少算三种情况。
• 可以证明从大到小枚举质数，上面那个剪枝会快些。。。

然后尽量不要取模，(ans)不会爆(long\,\,long)，最后取一次模就行了。

代码很好写(是真的比三元环计数好些啊，3k的码量差距啊)

#include<bits/stdc++.h>
using namespace std;
 
void read(int &x) {
    x=0;int f=1;char ch=getchar();
    for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
    for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
}
  
void print(int x) {
    if(x<0) putchar('-'),x=-x;
    if(!x) return ;print(x/10),putchar(x%10+48);
}
void write(int x) {if(!x) putchar('0');else print(x);putchar('
');}
 
const int maxn = 1e5+10;
const int lim = 40;
const int mod = 1e9+7;
 
#define ll long long 
 
int f[maxn],pri[maxn],tot,d[maxn],ff[lim+1][lim+1][lim+1],vis[maxn];
 
void prepare() {
    int n=1e5;for(int i=1;i<=n;i++) d[i]=1;
    for(int i=2;i<=n;i++) {
        if(!vis[i]) pri[++tot]=i;
        for(int j=i;j<=n;j+=i) d[j]++,vis[j]=1;
    }
    for(int i=1;i<=n;i++) f[i]=(f[i-1]+d[i])%mod;
    for(int i=1;i<=lim;i++)
        for(int j=1;j<=lim;j++)
            for(int k=1;k<=lim;k++)
                ff[i][j][k]=(0ll+ff[i-1][j][k]+ff[i][j-1][k]+ff[i][j][k-1]-ff[i-1][j-1][k]-ff[i-1][j][k-1]-ff[i][j-1][k-1]+ff[i-1][j-1][k-1]+d[i*j*k])%mod;
}
 
ll dfs(int now,int a,int b,int c) {
    if(!now) return 1ll*f[a]*f[b]%mod*f[c]%mod;
    if((!a)&&(!b)&&(!c)) return 0;
    int p=pri[now],mx=max(max(a,b),c);
    if(p>mx&&mx<=lim) return ff[a][b][c];
    ll ans=dfs(now-1,a,b,c);
    if(a>=p&&b>=p) ans-=dfs(now-1,a/p,b/p,c);
    if(a>=p&&c>=p) ans-=dfs(now-1,a/p,b,c/p);
    if(b>=p&&c>=p) ans-=dfs(now-1,a,b/p,c/p);
    if(a>=p&&b>=p&&c>=p) ans+=dfs(now-1,a/p,b/p,c/p)<<1ll;
    return ans;
}
 
void solve() {
    int a,b,c;read(a),read(b),read(c);
    int mx=max(max(a,b),c);
    int now=tot;
    while(pri[now]>mx) now--;
    write((int)((dfs(now,a,b,c)%mod+mod)%mod));
}
 
signed main() {
    prepare();int t;read(t);
    while(t--) solve();
    return 0;
}