[bzoj2671][莫比乌斯反演][数论]Calc

Description

给出N，统计满足下面条件的数对(a,b)的个数：
　　1.1<=a<b<=N
　　2.a+b整除a*b

Input

一行一个数N

Output

一行一个数表示答案

Sample Input

15

Sample Output

4

HINT

数据规模和约定

Test N Test N

1 <=10 11 <=5*10^7

2 <=50 12 <=10^8

3 <=10^3 13 <=2*10^8

4 <=510^3 14 <=310^8

5 <=210^4 15 <=510^8

6 <=2*10^5 16 <=10^9

7 <=2*10^6 17 <=10^9

8 <=10^7 18 <=2^31-1

9 <=2*10^7 19 <=2^31-1

10 <=3*10^7 20 <=2^31-1

题解

这题复杂度不会证的情况下还是不敢写的啊…
首先要求的就是$(a+b)|(a?b)$的数量
提出来一个$gcd=d$，可以知道就是$d?(u+v)|d^2?(u?v)$
化简一下就是$(u+v)|d?(u?v)$，这个时候$u,v$互质
随便反证法一下就可以知道$(u+v)$和$u?v$没有共同因子，所以只可能有$(u+v)|d$
所以$d=k?(u+v)$，其中$k=1...y$
扔进去就可以知道你需要找到满足$(u+v)?u?k<=n$的$k$的数量
显然这个$k$只有$?\frac{n}{(u+v)u}?$个
然后就是计数
$$\sum _{i=1}^n\sum _{j=1}^{i?1}[gcd(i,j)=1]?\frac{n}{(i+j)i}?$$
然后可以发现这个$i$有 $\sqrt{n}$的上界，所以不妨设$m=\sqrt{n}$然后继续做
随便反演一下就可以知道
$$\sum _{d=1}^m\mu (d)?\sum _{i=1}^{?\frac{m}{d}?}\sum _{j=1}^{i?1}?\frac{n}{d^2(i+j)i}?$$
然后前面是根号的复杂度，后面可以数论分块然后也是根号的复杂度
所以这大概是个分块套分块
所以复杂度大概是$n^{\frac{3}{4}}$的也就是$2^{24}$左右…
所以就可以跑的飞快了

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
#include<cmath>
#include<queue>
#include<vector>
#include<ctime>
#include<map>
#include<bitset>
#include<set>
#define LL long long
#define mp(x,y) make_pair(x,y)
#define pll pair<long long,long long>
#define pii pair<int,int>
using namespace std;
inline int read()
{
    int f=1,x=0;char ch=getchar();while(ch<'0'||ch>'9'){
    if(ch=='-')f=-1;ch=getchar();}while(ch>='0'&&ch<='9'){
    x=x*10+ch-'0';ch=getchar();}return x*f;
}
int stack[20];
inline void write(LL x)
{
    if(x<0){
    putchar('-');x=-x;}if(!x){
    putchar('0');return;}int top=0;while(x)stack[++top]=x%10,x/=10;while(top)putchar(stack[top--]+'0');
}
inline void pr1(int x){
    write(x);putchar(' ');}
inline void pr2(LL x){
    write(x);putchar('\n');}
const int MAXN=70005;
int mu[MAXN],pr[MAXN],plen;bool vis[MAXN];
void getmu(int N)
{
    mu[1]=1;for(int i=2;i<=N;i++){
    if(!vis[i])mu[i]=-1,pr[++plen]=i;for(int j=1;j<=plen&&i*pr[j]<=N;j++){
    vis[i*pr[j]]=true;if(!(i%pr[j])){
    mu[i*pr[j]]=0;break;}else mu[i*pr[j]]=-mu[i];}}
}
int n,m;
LL qry(int x)
{
    LL ret=0;LL temp=n/x/x;for(int i=1;i<=m/x;i++){
    LL T=temp/i;for(int j=i+1,last;j<i*2&&j<=T;j=last+1){
    last=min(T/(T/j),(LL)i*2-1);ret+=(last-j+1)*(T/j);}}return mu[x]*ret;
}
int gcd(int a,int b){
    return a==0?b:gcd(b%a,a);}
int main()
{
    getmu(MAXN-5);n=read();m=sqrt(n);LL ans=0;for(int i=1;i<=m;i++)if(mu[i]){
    ans+=qry(i);
// for(int j=1;j<=i-1;j++)if(gcd(i,j)==1)ans+=n/i/(i+j);}pr2(ans);return 0;
}