【国家集训队2012】【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 <=5*10^3 14 <=3*10^8

5 <=2*10^4 15 <=5*10^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=mx,b=my则a+b=m(x+y),a×b=xym2
m为a和b的最大公约数
有a+b|a×b即m(x+y)|xym2
故x+y|xym
假设x<y
有gcd(x,y)=1,所以gcd(x+y,x)=gcd(x+y,y)=1
那么只能是x+y|m
所以我们可以设m=k(x+y)
既然有1≤a<b≤n那么得到yk(x+y)≤n
所以对每组确定的x,y,符合条件的k有⌊ny(x+y)⌋个
显然可以有y≤n√−1,所以枚举y的取值.这一步是O(n√)的.
这时x的取值必然是1→y−1
直接枚举这样的x,y，当gcd(x,y)=1时计算答案，ans+=ny(x+y)是正确的.可以得到50%的分数.
对于剩下的50%,我们可以考虑使用Mobius反演.
化到上面我们可以发现自己的问题转化成了求

∑i=lr[gcd(i,y)=1]

这样感觉就可以用反演了.
令F(k)=⌊lk⌋−⌊rk⌋
可以知道答案是

∑k|yμ(k)∗F(k)

处理出y的因子然后用反演常用的枚举除法就行了.这一步是O(n√−−−√)
所以总复杂度O(n34).
这么久没写过数论感觉自己都变成傻逼了…
UPD.2015-9-2
我的时间复杂度被Hack啦真是抱歉…确实是分析错了QAQ
好像时间复杂度看起来就像是某种…接近暴力的?…
但是实践出真知!他真的跑得很快!

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>#include<algorithm>#define MAXN 100010#define LL long longusing namespace std;int n;bool not_prime[MAXN];int prime[MAXN],top;int phi[MAXN]={0,1},mu[MAXN]={0,1};int divs[MAXN],Top;LL ans,temp;inline void check_prime(){    for (int i=2;i<MAXN;i++)    {        if (!not_prime[i])  prime[++top]=i,phi[i]=i-1,mu[i]=-1;        for (int j=1;j<=top&&i*prime[j]<MAXN;j++)        {            not_prime[i*prime[j]]=1;            if (i%prime[j]==0)            {                mu[i*prime[j]]=0;                phi[i*prime[j]]=phi[i]*prime[j];                break;            }            else            {                mu[i*prime[j]]=-mu[i];                phi[i*prime[j]]=phi[i]*(prime[j]-1);            }        }    }}inline void solve(int x)//处理因子 {    Top=0;    for (int i=1;i*i<=x;i++)        if (x%i==0) divs[++Top]=i,divs[++Top]=x/i;    if (divs[Top]==divs[Top-1]) Top--;    sort(divs+1,divs+Top+1);}int main(){    check_prime();    scanf("%d",&n);    int last;    for (LL i=1;i*(i+1)<=n;i++)    {        solve(i);        for (LL j=1;j<i&&i*(i+j)<=n;j=last+1)        {            last=min(n/(n/i/(i+j))/i-i,i-1);//kx(x+y)<=n x<y 已知y=i,求x            temp=0;            for (LL k=1;divs[k]<=last;k++)  temp+=mu[divs[k]]*(last/divs[k]-(j-1)/divs[k]);            ans+=n/i/(i+j)*temp;//ans+=k*temp        }    }    printf("%lld\n",ans);}