bzoj 2820 YY的GCD

题目描述

神犇YY虐完数论后给傻×kAc出了一题

给定N, M,求1<=x<=N, 1<=y<=M且gcd(x, y)为质数的(x, y)有多少对

kAc这种傻×必然不会了，于是向你来请教……

多组输入

输入输出格式

输入格式：

第一行一个整数T 表述数据组数

接下来T行，每行两个正整数，表示N, M

输出格式：

T行，每行一个整数表示第i组数据的结果

输入输出样例

输入样例#1： 复制 
2
10 10
100 100
输出样例#1： 复制 
30
2791

说明

T = 10000

N, M <= 10000000

∑p为素数∑i∑j[gcd(i,j)==p]

=∑p为素数∑in/p∑jm/p[gcd(i/p,j/p)==1]

=∑p为素数∑in/p∑jm/p∑d|i,d|j μ(d)

=∑p∑μ(d)[n/dp][m/dp]

=∑dp[n/dp][m/dp]∑d|dp且p为素数μ(d)

F(dp)=∑d|dp且p为素数μ(d)

F可以nlnn求解

查询时分块，因为显然有很多n/i(或m/i)相等

然后常数优化：比如尽量写int，加register，减少循环次数

 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<algorithm>
 5 #include<cmath>
 6 using namespace std;
 7 typedef long long lol;
 8 int f[10000001],mu[10000001],n,m;
 9 int prime[10000001],tot;
10 bool vis[10000001];
11 void pre()
12 {register int i,j;
13   mu[1]=1;
14   tot=0;
15   for (i=2;i<=10000000;i++)
16     {
17       if (vis[i]==0)
18     {
19       prime[++tot]=i;
20       mu[i]=-1;
21     }
22       for (j=1;j<=tot;j++)
23     {
24       if (1ll*i*prime[j]>10000000)break;
25       vis[i*prime[j]]=1;
26       if (i%prime[j]==0)
27         {
28           mu[i*prime[j]]=0;
29           break;
30         }
31       else mu[i*prime[j]]=-mu[i];
32     }
33     }
34   for (i=1;i<=10000000;i++)
35     {
36       for (j=1;j<=tot;j++)
37     {
38       if (prime[j]*i>10000000) break;
39       f[prime[j]*i]+=mu[i];
40     }
41       f[i]+=f[i-1];
42     }
43 }
44 int main()
45 {lol T;
46   register int i,pos;
47   register lol ans;
48   pre();
49   cin>>T;
50   while (T--)
51     {
52       scanf("%d%d",&n,&m);
53       if (n>m) swap(n,m);
54       pos=0;
55       ans=0;
56       for (i=1;i<=n;i=pos+1)
57     {
58       if ((n/(n/i))<(m/(m/i))) pos=n/(n/i);
59       else pos=m/(m/i);
60       ans+=1ll*(n/i)*(m/i)*(f[pos]-f[i-1]);
61     }
62       printf("%lld\n",ans);
63     }
64 }

 

posted @ 2018-01-16 08:07  Z-Y-Y-S  阅读(233)  评论(0)  编辑  收藏  举报