bzoj P3309 DZY Loves Math——solution

对于正整数n，定义f(n)为n所含质因子的最大幂指数。例如f(1960)=f(2^3 * 5^1 * 7^2)=3, f(10007)=1, f(1)=0。
给定正整数a,b，求:

$$\sum_{i=1}^{i<=a}\sum_{j=1}^{j<=b}f(gcd(i,j))$$

bzojP3309

http://www.lydsy.com/JudgeOnline/problem.php?id=3309

---

---

 

化式子：

$$\sum_{i=1}^{i<=a}\sum_{j=1}^{j<=b}f(gcd(i,j))$$

$$\sum_{x=1}^{min(a,b)}f(x)·\sum_{x|d}\mu({d \over x})({a \over d})({b \over d})$$

$$\sum_{x=1}^{min(a,b)}f(x)·\sum_{d=1}^{min(a,b) \over x}\mu(d)({a \over dx})({b \over dx})$$

$$\sum_{dx=1}^{min(a,b)}({a \over dx})({b \over dx})\sum_{d=1}^{d<=dx}f(x)*\mu(d)$$

$$\sum_{x=1}^{min(a,b)}({a \over dx})({b \over dx})(f*\mu)(x)$$

这个化式子的过程旨在尽可能地把无关a,b的，可预处理的部分提出来；

现在只要预处理出f与$\mu$的卷积即可通过枚举除法做到单次$O(\sqrt{n})的效率$

f和$\mu$可以通过线性筛求出，

如何在较好的时间内求出他们的卷积呢；

不会，

不过打表可以发现：

$$当\mu(i)=1时，(f*\mu)(i)=-1$$

$$当\mu(i)=-1时，(f*\mu)(i)=1$$

$$当\mu(i)=0时，若i的所有质因子齐次，次数为k，则(f*\mu)(i)=(f*\mu)(^{k}\sqrt{i})，否则(f*\mu)(i)=0$$

这样处理好$f*\mu$，然后枚举除法即可；

由于笔者太弱，于是直接上了单次处理$O(log_2log_2n)$的线性筛

2018.01.25 upd:我原来以为暴力计算卷积函数是$O(n\sqrt{n})$的，后来才发现，这其实是$O(n*ln(n))$的，所以暴力计算f*mu应该也能过

2018.01.25 upd:一开始不觉得数列求和分析复杂度可以积分（好像在项数少的时候不太拟合），现在看来好像可以（项数多的时候比较好）

2018.07.12 upd:这里顺便记录一下枚举除法套枚举除法是$O(N^{3\over 4})$的，积分证得

代码：

 1 #include<cmath>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<algorithm>
 5 #define LL long long 
 6 using namespace std;
 7 const int MAXN=1e7;
 8 LL f_mu[MAXN+5];
 9 int pri[MAXN],f[MAXN+5],mu[MAXN+5],p[MAXN+5],cnt;
10 bool vis[MAXN];
11 LL sum[MAXN+5];
12 void prime();
13 LL Sqr(LL ,int );
14 LL work(int ,int );
15 int main()
16 {
17     int i,j,k,n,a,b;
18     prime();
19     sum[0]=f_mu[1]=sum[1]=0;
20     for(i=2;i<=MAXN;i++){
21         if(mu[i]==1)
22             f_mu[i]=-1;
23         if(mu[i]==-1)
24             f_mu[i]=1;
25         if(mu[i]==0){
26             if(i==Sqr(p[i],f[i]))
27                 f_mu[i]=f_mu[p[i]];
28             else
29                 f_mu[i]=0;
30         }
31         sum[i]=sum[i-1]+f_mu[i];
32     }
33     scanf("%d",&n);
34     for(i=1;i<=n;i++){
35         scanf("%d%d",&a,&b);
36         printf("%lld\n",work(a,b));
37     }
38     return 0;
39 }
40 void prime(){
41     int i,j;
42     vis[1]=true,mu[1]=1;
43     for(i=2;i<=MAXN;i++){
44         if(!vis[i]){
45             f[i]=1,p[i]=i,mu[i]=-1;
46             pri[++cnt]=i;
47         }
48         for(j=1;j<=cnt&&pri[j]*i<=MAXN;j++){
49             vis[i*pri[j]]=true;
50             if(i%pri[j]){
51                 mu[i*pri[j]]=-mu[i];
52                 f[i*pri[j]]=f[i],p[i*pri[j]]=p[i]*pri[j];
53             }
54             else{
55                 mu[i*pri[j]]=0;
56                 f[i*pri[j]]=f[i]+(i%Sqr(pri[j],f[i])==0);
57                 p[i*pri[j]]=p[i];
58                 break;
59             }
60         }
61     }
62 }
63 LL Sqr(LL x,int n){
64     LL ret=1;
65     while(n){
66         if(n&1)
67             ret*=x;
68         n>>=1,x*=x;
69     }
70     return ret;
71 }
72 LL work(int a,int b){
73     int MIN=min(a,b),i,las;
74     LL ret=0;
75     for(las=0,i=1;i<=MIN;las=i,i=min(a/(a/(las+1)),b/(b/(las+1)))){
76         ret+=(sum[i]-sum[las])*(1ll*a/i)*(1ll*b/i);
77         if(i==MIN)break;
78     }
79     return ret;
80 }

 

转载于:https://www.cnblogs.com/nietzsche-oier/p/8337950.html