洛谷P4466 [国家集训队] 和与积

所求即为：

[large sum_{a=1}^nsum_{b=a+1}^nleft[ a+b mid ab ight] ]

设 (gcd(a,b)=d,a=id,b=jd)，得：

[largeegin{aligned} d(i+j)&mid ijd^2\ (i+j)&mid ijd\ end{aligned} ]

发现因为 (gcd(i,j)=1)，所以 ((i+j) ot mid ij)，这是因为有：

[large gcd(i+j,i)=gcd(i+j,j)=1 ]

因此得 ((i+j)mid d)。原式变为：

[largeegin{aligned} &sum_{i=1}^{sqrt n}sum_{j=i+1}^{sqrt n}left[ gcd(i,j)=1 ight]leftlfloor frac{n}{j(i+j)} ight floor \ =&sum_{i=1}^{sqrt n}sum_{j=i+1}^{sqrt n}sum_{dmid i and d mid j}mu(d)leftlfloor frac{n}{j(i+j)} ight floor \ =&sum_{d=1}^{sqrt n}mu(d)sum_{i=1}^{sqrt n}left[ dmid i ight]sum_{j=i+1}^{sqrt n}left[ dmid j ight]leftlfloor frac{n}{j(i+j)} ight floor \ =&sum_{d=1}^{sqrt n}mu(d)sum_{i=1}^{leftlfloorfrac{sqrt n}{d} ight floor}sum_{j=i+1}^{leftlfloorfrac{sqrt n}{d} ight floor}leftlfloor frac{n}{d^2j(i+j)} ight floor \ =&sum_{d=1}^{sqrt n}mu(d)sum_{i=2}^{leftlfloorfrac{sqrt n}{d} ight floor}sum_{j=i+1}^{2i-1}leftlfloor frac{leftlfloor frac{n}{d^2i} ight floor}{j} ight floor \ end{aligned} ]

最后一步是分别枚举 (j) 和 (i+j)。枚举 (d,j) 后，(leftlfloor frac{n}{d^2i} ight floor) 就为定值了，然后就可以数论分块了。

大致分析一下复杂度，得总枚举次数为：

[large sum_{i=1}^{sqrt n}frac{sqrt n}{i}sqrt{frac{sqrt n}{i}}=n^{frac{3}{4}}sum_{i=1}^{sqrt n}frac{1}{i^{frac{3}{2}}} ]

后一项为黎曼函数 (zeta(x))，当 (x=frac{3}{2}) 时，其取值约为 (2.6)，因此复杂度为 (O(n^{frac{3}{4}}))。

#include<bits/stdc++.h>
#define maxn 47350
using namespace std;
typedef long long ll;
template<typename T> inline void read(T &x)
{
    x=0;char c=getchar();bool flag=false;
    while(!isdigit(c)){if(c=='-')flag=true;c=getchar();}
    while(isdigit(c)){x=(x<<1)+(x<<3)+(c^48);c=getchar();}
    if(flag)x=-x;
}
int n,m,tot;
ll ans;
int p[maxn],mu[maxn];
bool tag[maxn];
void init()
{
    mu[1]=1;
    for(int i=2;i<=m;++i)
    {
        if(!tag[i]) p[++tot]=i,mu[i]=-1;
        for(int j=1;j<=tot;++j)
        {
            int k=i*p[j];
            if(k>m) break;
            tag[k]=true;
            if(i%p[j]) mu[k]=mu[i]*mu[p[j]];
            else
            {
                mu[k]=0;
                break;
            }
        }
    }
}
ll calc(ll d)
{
    ll v=0;
    for(int i=2;i<=m/d;++i)
    {
        ll val=n/(d*d*i);
        if(!val) continue;
        for(int l=i+1,r;l<=2*i-1;l=r+1)
        {
            if(val/l==0) break;
            r=min(val/(val/l),(ll)2*i-1),v+=val/l*(r-l+1);
        }
    }
    return v;
}
int main()
{
    read(n),m=sqrt(n),init();
    for(int i=1;i<=m;++i)
        if(mu[i])
            ans+=mu[i]*calc(i);
    printf("%lld",ans);
    return 0;
}