Luogu P2522 [HAOI2011]Problem b

如果你做过[Luogu P3455 POI2007]ZAP-Queries就很好办了，我们发现那一题求的是(sum_{i=1}^asum_{j=1}^b[gcd(i,j)=d])，就是这道题的特殊情况。

因此我们直接令(operatorname{calc}(x,y,d))表示(sum_{i=1}^xsum_{j=1}^y[gcd(i,j)=d])，然后直接容斥即可：

[ans=operatorname{calc}(b,d,k)-operatorname{calc}(a-1,d,k)-operatorname{calc}(b,c-1,k)+operatorname{calc}(a-1,c-1,k) ]

关于(operatorname{calc}(x,y,d))的求法可以看上面那题的Sol，这里不再赘述。

CODE

#include<cstdio>
#include<cctype>
#define RI register int
using namespace std;
const int P=50005;
int t,a,b,c,d,k,prime[P+5],cnt,mu[P+5],sum[P+5]; bool vis[P+5];
class FileInputOutput
{
    private:
        #define tc() (A==B&&(B=(A=Fin)+fread(Fin,1,S,stdin),A==B)?EOF:*A++)
        #define pc(ch) (Ftop<S?Fout[Ftop++]=ch:(fwrite(Fout,1,S,stdout),Fout[(Ftop=0)++]=ch))
        #define S 1<<21
        char Fin[S],Fout[S],*A,*B; int Ftop,pt[25];
    public:
        FileInputOutput() { A=B=Fin; Ftop=0; } 
        inline void read(int &x)
        {
            x=0; char ch; while (!isdigit(ch=tc()));
            while (x=(x<<3)+(x<<1)+(ch&15),isdigit(ch=tc()));
        }
        inline void write(long long x)
        {
            if (!x) return (void)(pc(48),pc('
')); RI ptop=0;
            while (x) pt[++ptop]=x%10,x/=10; while (ptop) pc(pt[ptop--]+48); pc('
');
        }
        inline void Fend(void)
        {
            fwrite(Fout,1,Ftop,stdout);
        }
        #undef tc
        #undef pc
        #undef S
}F;
#define Pi prime[j]
inline void Euler(void)
{
    vis[1]=mu[1]=1; RI i,j; for (i=2;i<=P;++i)
    {
        if (!vis[i]) prime[++cnt]=i,mu[i]=-1;
        for (j=1;j<=cnt&&i*Pi<=P;++j)
        {
            vis[i*Pi]=1; if (i%Pi) mu[i*Pi]=-mu[i]; else break;
        }
    }
    for (i=1;i<=P;++i) sum[i]=sum[i-1]+mu[i];
}
#undef Pi
inline int min(int a,int b)
{
    return a<b?a:b;
}
inline long long calc(int n,int m,int d)
{
    long long ans=0; int lim=min(n/d,m/d);
    for (RI l=1,r;l<=lim;l=r+1)
    {
        r=min(n/(n/l),m/(m/l)); ans+=1LL*(n/(l*d))*(m/(l*d))*(sum[r]-sum[l-1]);
    }
    return ans;
}
int main()
{
    //freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
    for (Euler(),F.read(t);t;--t)
    {
        F.read(a); F.read(b); F.read(c); F.read(d); F.read(k);
        F.write(calc(b,d,k)-calc(a-1,d,k)-calc(b,c-1,k)+calc(a-1,c-1,k));
    }
    return F.Fend(),0;
}