BZOJ3994: [SDOI2015]约数个数和

BZOJ3994

d(x)x
ni=1mj=1d(i,j)
有个很神奇的结论。。

=i=1nj=1m[gcd(i,j)==1]nimj

QnQ 像我这样脑残看不懂证明的就只能记结论了哎。
有了上面这个式子以后很容易往莫比乌斯函数上想。
n<=m
i=1nj=1mnimjdμ(d)

=d=1nμ(d)i=1n/dj=1m/dn/dim/dj

f(x)=xi=1xi
O(nn) 预处理出 f .(还可以线性筛,但是这个题这个复杂度也可以接受)
然后每次查询就下底函数分块辣。总复杂度(O(N+T)N)

【代码】

#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cmath>
#include <queue>
#define N 50005
#define M 50000005
#define INF 0x7fffffff
using namespace std;
typedef long long ll;
typedef pair<ll,int> pa;

ll read()
{
    ll x=0,f=1;char ch=getchar();
    while(!isdigit(ch)){if(ch=='-') f=-1;ch=getchar();}
    while(isdigit(ch)){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
    return x*f;
}

int n,m,T;
int f[N],Miu[N],p[5135],sum[N];
bool Not_Prime[N]; 

void Pre()
{
    Miu[1]=sum[1]=1;
    for(int i=2;i<N;i++) {
        if(!Not_Prime[i]) p[++p[0]]=i,Miu[i]=-1;
        for(int j=1;j<=p[0]&&i*p[j]<N;j++)
        {
            Not_Prime[i*p[j]]=1;
            if(i%p[j]==0) break;
            Miu[i*p[j]]=-Miu[i];
        }
        sum[i]=sum[i-1]+Miu[i];
    }
    for(int i=1;i<N;i++)
    {
        int Next;
        for(int j=1;j<=i;j=Next+1)
        {
            Next=i/(i/j);
            f[i]+=(Next-j+1)*(i/j);
        }
    }
}

void Solve()
{
    if(n>m) swap(n,m);
    int Next;ll ans=0;
    for(int i=1;i<=n;i=Next+1)
    {
        Next=min(n/(n/i),m/(m/i));
        ans+=1LL*(sum[Next]-sum[i-1])*f[n/i]*f[m/i];
    }
    printf("%lld\n",ans);
}

int main()
{
    Pre();
    T=read();
    while(T--)
    {
        n=read(),m=read();
        Solve();
    }
    return 0;
}
  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值