Lucas的數論：杜教篩，莫比烏斯反演

Description：

求$\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} d(i \times j)$

$d(i)$表示$i$的約數個數和。$n \leq 10^9$

 

廢話：

很久沒有作反演了感受本身都不會了。。。

作了一遍發現本身真的不會了

手推了不知道多久終於推出了式子中間還錯了一遍打了一半發現過不去樣例

 

題解：

標籤都知道了那也就沒什麼好說的了，直接上式子（類比《約數個數和》那道題）（如下分數皆表示整除向下取整）

$\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} d(i \times j)$

$=\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} \sum\limits_{a|i} \sum\limits_{b|j} [gcd(a,b)==1] $

$=\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{n} \sum\limits_{a|i} \sum\limits_{b|j} \sum\limits_{d|a \  and \  d|b} \mu (d) $

$=\sum\limits_{d=1}^{n} \mu (d) \sum\limits_{a=1}^{\frac{n}{d}} \sum\limits_{b=1}^{\frac{n}{d}} \frac{\frac{n}{d}}{a} \frac{\frac{n}{d}}{b}$

$=\sum\limits_{d=1}^{n} \mu (d) (\sum\limits_{i=1}^{\frac{n}{d}} \frac{\frac{n}{d}}{i})^2$

而後第二個$\sum$後面的東西視做函數$f(\frac{n}{d})$，能夠整除分塊，前面就是杜教篩獲得$\mu$的前綴和

而後函數$f(\frac{n}{d})$的求值其實也就是一個整除分塊。

其實$f(x)=\sum\limits_{i=1}^{x} d(i)$即$\sum\limits_{i=1}^{x} d(i) =\sum\limits_{i=1}^{x} \frac{x}{i}$

這裏的$d(x)$函數也是約數個數的意思。%%%LrefraiNC

從含義上也能夠理解。因此$f(x)$函數也能夠杜教篩。可是我尚未作過$d(i)$的杜教篩因此我沒有打。

杜教篩線篩預處理的範圍應該在$n^{\frac{2}{3}}=10^6$否則慢的要死（雖然說沒有T）

 1 #include<cstdio>
 2 #define mod 1000000007
 3 struct hash_map{
 4     int w[10000005],fir[3000005],l[10000005],to[10000005],cnt;
 5     int &operator[](int x){
 6         int r=x%3000001;
 7         for(int i=fir[r];i;i=l[i])if(to[i]==x)return w[i];
 8         l[++cnt]=fir[r];fir[r]=cnt;to[cnt]=x;return w[cnt]=-123456789;
 9     }
10 }M;
11 int p[1000005],mu[1000005],pc;char np[1000005];
12 int sum(int n){
13     if(n<=1000000)return mu[n];
14     if(M[n]!=-123456789)return M[n];
15     int a=1;
16     for(int l=2,r,A;l<=n;l=r+1)A=n/l,r=n/A,a-=sum(A)*(r-l+1);
17     return M[n]=a;
18 }
19 long long cal(int x,long long ans=0){
20     for(int l=1,r,a;l<=x;l=r+1)a=x/l,r=x/a,ans=(ans+a*(r-l+1ll))%mod;
21     return ans*ans%mod;
22 }
23 int main(){
24     mu[1]=1;
25     for(int i=2;i<=1000000;++i){
26         if(!np[i])p[++pc]=i,mu[i]=-1;
27         for(int j=1;j<=pc&&i*p[j]<=1000000;++j)
28             if(i%p[j])mu[i*p[j]]=-mu[i],np[i*p[j]]=-1;
29             else np[i*p[j]]=-1;
30         mu[i]+=mu[i-1];
31     }
32     int n;long long ans=0;scanf("%d",&n);
33     for(int l=1,r,a;l<=n;l=r+1)a=n/l,r=n/a,ans+=1ll*cal(a)*(sum(r)-sum(l-1))%mod;
34     printf("%lld\n",(ans%mod+mod)%mod);
35 }

View Code