[国家集训队] Crash的数字表格 / JZPTAB

题目描述

今天的数学课上，Crash 小朋友学习了最小公倍数（Least Common Multiple）。对于两个正整数 \(a\) 和 \(b\)，\(\text{lcm}(a,b)\) 表示能同时被 \(a\) 和 \(b\) 整除的最小正整数。例如，\(\text{lcm}(6, 8) = 24\)。

回到家后，Crash 还在想着课上学的东西，为了研究最小公倍数，他画了一张 $ n \times m$ 的表格。每个格子里写了一个数字，其中第 \(i\) 行第 \(j\) 列的那个格子里写着数为 \(\text{lcm}(i, j)\)。

看着这个表格，Crash 想到了很多可以思考的问题。不过他最想解决的问题却是一个十分简单的问题：这个表格中所有数的和是多少。当 \(n\) 和 \(m\) 很大时，Crash 就束手无策了，因此他找到了聪明的你用程序帮他解决这个问题。由于最终结果可能会很大，Crash 只想知道表格里所有数的和对 \(20101009\) 取模后的值。

输入格式

输入包含一行两个整数，分别表示 \(n\) 和 \(m\)。

输出格式

输出一个正整数，表示表格中所有数的和对 \(20101009\) 取模后的值。

样例 #1

样例输入 #1

4 5

样例输出 #1

122

提示

样例输入输出 1 解释

该表格为：

数据规模与约定

• 对于 \(30\%\) 的数据，保证 \(n, m \le 10^3\)。
• 对于 \(70\%\) 的数据，保证 \(n, m \le 10^5\)。
• 对于 \(100\%\) 的数据，保证 \(1\le n,m \le 10^7\)。

链接

\(\displaystyle\sum_{i=1}^n\displaystyle\sum_{j=1}^m\frac{ij}{gcd(i,j)}\)

\(\displaystyle\sum_{d=1}^{min(n,m)}\displaystyle\sum_{i=1}^{\frac{n}{d}}\displaystyle\sum_{j=1}^{\frac{m}{d}}i*j*d[gcd(i,j)=1]\)，\(\displaystyle\sum_{d=1}^{min(n,m)}d\displaystyle\sum_{i=1}^{\frac{n}{d}}\displaystyle\sum_{j=1}^{\frac{m}{d}}i*j*\epsilon(gcd(i,j)=1)\)

\(\displaystyle\sum_{d=1}^{min(n,m)}d\displaystyle\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\displaystyle\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i*j*\displaystyle\sum_{e|gcd(i,j)}\mu(e)\)
\(\displaystyle\sum_{d=1}^{min(n,m)}d\displaystyle\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\displaystyle\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i*j*\displaystyle\sum_{e=1}^{min(i,j)}\mu(e)[e|i][e|j]\)

\(\displaystyle\sum_{d=1}^{min(n,m)}d\displaystyle\sum_{e=1}^{min(n,m)}\mu(e)\displaystyle\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}i*[e|i]\displaystyle\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}j*[e|j]\)

\(\displaystyle\sum_{d=1}^{min(n,m)}d\displaystyle\sum_{e=1}^{min(n,m)}\mu(e)\times e^2\displaystyle\sum_{i=1}^{\lfloor\frac{n}{de}\rfloor}i\displaystyle\sum_{j=1}^{\lfloor\frac{m}{de}\rfloor}j\)

把后面的部分拿出来。设\(sum(n,m)=\displaystyle\sum_{e=1}^{min(n,m)}\mu(e)\times e^2\displaystyle\sum_{i=1}^{\lfloor\frac{n}{e}\rfloor}i\displaystyle\sum_{j=1}^{\lfloor\frac{m}{e}\rfloor}j\)
前面的部分的\(\mu(e)*e^2\)的前缀和可以用埃式筛\(O(n\log\log n)\)预处理，然后可以发现，后面的部分在\(e\)属于一个相邻的连续范围内的时候是不会变化的，这个范围就是整除分块的范围。这样整个\(sum\)函数就可以用整除分块在\(O(\sqrt n)\)的时间处理。
所以式子变为\(\displaystyle\sum_{d=1}^{min(n,m)}d\times sum(\lfloor\frac {n}{d}\rfloor,\lfloor\frac {m}{d}\rfloor)\)，又可以发现，\(\lfloor\frac {n}{d}\rfloor\)在\(d\)变化时在一定范围内不变，可以整除分块。所以两个整除分块就好了。

预处理的时间复杂度应该没问题。。但是两个整除分块的嵌套的复杂度真的没事情吗...

#include<bits/stdc++.h>
#define ll long long
using namespace std;
inline ll read()
{
    char c=getchar();ll a=0,b=1;
    for(;c<'0'||c>'9';c=getchar())if(b=='-')b=-1;
    for(;c>='0'&&c<='9';c=getchar())a=a*10+c-'0';
    return a*b;
}
const ll N=1e7,P=20101009;
ll prime[N+1],vis[N+1],cnt,mu[N+1],sum[N+1],n,m;
void init()
{
    mu[1]=1;
    for(ll i=2;i<=min(n,m);i++)
    {
        if(!vis[i])prime[++cnt]=i,mu[i]=-1;
        for(ll j=1;j<=cnt&&i*prime[j]<=min(n,m);j++)
        {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=mu[i]*mu[prime[j]];
        }
    }
    for(ll i=1;i<=min(n,m);i++)
    {
        sum[i]=(sum[i-1]%P+i*i%P*(mu[i]+P)%P)%P;
    }
}
ll count(ll n,ll m)
{
    return (((1LL*(n+1)*(n)/2) %P) *(1LL*(m+1)*(m)/2%P)%P) %P;
}
ll Sum(ll n,ll m)
{
    ll ans=0;
    for(ll l=1,r=0;l<=min(n,m);l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        ans=(ans+((sum[r]-sum[l-1])+P)%P*(count(n/l,m/l)%P))%P;
    }
    return ans;
}
int main()
{
    n=read(),m=read();
    init();
    ll ans=0;
    for(ll l=1,r=0;l<=min(n,m);l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        ans=(ans+(((r+l)*(r-l+1)/2)%P*Sum(n/l,m/l)%P)%P)%P;
    }
    cout<<ans%P<<endl;
    return 0;
}

跑的挺快...
玛德跳了2两个小时，就因为一个部分的式子写错了，然后刚刚好的没有算错，而是导致了大数据时取模会爆炸。我草我草我草我草受不了了。