1. 程式人生 > 實用技巧 >杜教篩&Min_25篩學習筆記

杜教篩&Min_25篩學習筆記

0 - 前置

莫比烏斯反演

詳見 Link ,這裡摘錄一些式子。

Dirichlet卷積:

\[(f*g)(x)=\sum_{d|x}f(d)g(\frac{x}{d}) \]

莫比烏斯反演:

\[f=g*1=>g=f*\mu\\ f(x)=\sum_{d|x}g(d)=>g(x)=\sum_{d|x}f(d)\times \mu(\frac{x}{d}). \]

尤拉函式

性質:\(\sum_{d|n} \varphi(d)=n\)

表示成 Dirichlet卷積 的形式:\(\varphi*I=id\) .

捲上 \(\mu\)

\[\varphi*I=id\\ =>\varphi*I*\mu=id*\mu\\ =>\varphi * \epsilon=id*\mu \]

所以有

\[\varphi(n)=\sum_{d|n}\mu(d)\cdot\frac{n}{d} \]

\[\frac{\varphi(n)}{n}=\sum_{d|n}\frac{\mu(d)}{d} \]

1 - 杜教篩

Method

用途:以低於線性的複雜度計算積性函式字首和 ,即計算:\(\sum_{i=1}^n f(i)\) .

推導:

構造兩個積性函式 \(h,g\) ,使得 \(h=f*g\) .

\(S(n)=\sum_{i=1}^n f(i)\) ,那麼有

\[\begin{aligned} \sum_{i=1}^nh(i) &=\sum_{i=1}^n\sum_{d|i}g(d)\cdot f(\frac{i}{d})\\\\ &=\sum_{d=1}^ng(d)\cdot \sum_{i=1}^{\lfloor n/d\rfloor} f(i)\\\\ &=\sum_{d=1}^ng(d)\cdot S(\lfloor \frac{n}{d}\rfloor )\\\\ &=g(1)\cdot S(n)+\sum_{d=2}^{n}g(d)\cdot S(\lfloor\frac{n}{d}\rfloor)\\\\ g(1)S(n)&=\sum_{i=1}^nh(i)-\sum_{d=2}^ng(d)\cdot S(\lfloor \frac{n}{d}\rfloor ) \end{aligned} \]

所以關鍵就在於找一個容易得到的 \(h\) .求得之後對後面整除分塊就好了,複雜度 \(\mathcal{O}(n^{\frac 23})\) .

關於複雜度 (@George1123)

直接遞迴求是 \(\Theta(n^{\frac{3}{4}})\)

\[\begin{split} T(n)=&\sqrt n+\sum\limits_{i=2}^{\sqrt n}\left(T(i)\cdot T(\lfloor\frac{n}{i}\rfloor)\right)\\ \ge&\sqrt n+\sum\limits_{i=2}^{\sqrt n}\left(\sqrt i\cdot \sqrt{\lfloor\frac{n}{i}\rfloor)}\right)\\ \ge&\sqrt n+\sum\limits_{i=2}^{\sqrt n}2\sqrt{\sqrt n}\\ =&n^{\frac{3}{4}}\\ \end{split} \]

線性篩求出前 \(n^{\frac23}\) 個,然後再杜教篩,是 \(\Theta(n^{\frac{2}{3}})\) .

\[\begin{split} T(n)=&m+\sum\limits_{i=2}^{\lfloor\frac nm\rfloor}\sqrt{\lfloor\frac ni\rfloor}\\ =&m+\frac{n}{\sqrt m}\\ \ge&2n^{\frac{2}{3}}(=: \operatorname{while} m=n^{\frac{2}{3}})\\ \end{split} \]

Example

莫比烏斯函式

\[S(n)=\sum_{i=1}^n\mu(i). \]

由於 \(\mu*I=\epsilon\) ,所以 令 \(I=>g,\epsilon=>h\) ,根據上面的式子

\[g(1)S(n)=\sum_{i=1}^nh(i)-\sum_{d=2}^ng(d)\cdot S(\lfloor \frac{n}{d}\rfloor ) \]

\[\begin{aligned} 1\times S(n)&=\sum_{i=1}^n\epsilon(i)-\sum_{d=2}^nS(\lfloor \frac{n}{d}\rfloor )\\\\ S(n) &= 1-\sum_{d=2}^nS(\lfloor \frac nd\rfloor) \end{aligned} \]

尤拉函式

\[S(n)=\sum_{i=1}^n\varphi(i) \]

由於 \(\varphi *I=ID\) ,所以

\[\begin{aligned} S(n) &=\sum_{i=1}^n i-\sum_{d=2}^n S(\lfloor \frac{n}{d}\rfloor)\\\\ &=\frac{n(n+1)}{2}-\sum_{d=2}^nS(\lfloor \frac{n}{d}\rfloor) \end{aligned} \]

?函式

\[S(n)=\sum_{i=1}^ni\cdot \varphi(i) \]

考慮 Dirichlet卷積的形式:

\[\sum_{d|n}(d\cdot \varphi(d))\cdot g(\frac{n}{d}) \]

\(ID=>g\) ,有

\[\sum_{d|n}(d\cdot \varphi(d))\cdot \frac{n}{d}=\sum_{d|n}n\cdot \varphi(d)=n\sum_{d|n}\varphi(d)=n^2 \]

所以有

\[S(n)=\sum_{i=1}^ni^2-\sum_{d=2}^nd\cdot S(\lfloor \frac{n}{d}\rfloor) \]

??函式

\[\sum_{i=1}^n\mu(i)i^2 \]

\(g(i)=i^2\) ,有

\[g*f(n)=\sum_{d|n}\mu(d)d^2\cdot\Big(\frac{n}{d}\Big)^2=\sum_{d|n}\mu(d)=[n=1] \]

所以

\[S(n)=\sum_{i=1}^nh(i)-\sum_{d=2}^nd^2\cdot S(\lfloor \frac{n}{d}\rfloor )=1-\sum_{d=2^n}d^2\cdot S(\lfloor \frac{n}{d}\rfloor) \]

Code

大常數選手打的 模板

//Author: RingweEH
const int N=5e5+10;
int pri[N],tot=0,n;
ll mu[N],phi[N];
bool is[N];

void Sieve()
{
    mu[1]=phi[1]=1; is[1]=1;
    for ( int i=2; i<=N-10; i++ )
    {
        if ( !is[i] ) { pri[++tot]=i,mu[i]=-1,phi[i]=i-1; }
        for ( int j=1; j<=tot && (i*pri[j]<=(N-10)); j++ )
        {
            is[i*pri[j]]=1;
            if ( i%pri[j]==0 ) { mu[i*pri[j]]=0,phi[i*pri[j]]=phi[i]*pri[j]; break; }
            mu[i*pri[j]]=-mu[i]; phi[i*pri[j]]=phi[i]*phi[pri[j]];
        }
    }
    for ( int i=2; i<=(N-10); i++ )
        mu[i]+=mu[i-1],phi[i]+=phi[i-1];
}

map<int,ll> mp_mu,mp_phi;

ll Sieve_Phi( int n )
{
    if ( n<=N-10 ) return phi[n];
    if ( mp_phi[n] ) return mp_phi[n];
    ll _res=0;
    for ( int l=2,r; l<=n; l=r+1 )
    {
        r=n/(n/l); _res+=(r-l+1)*Sieve_Phi(n/l);
        if ( r>=2147483647 ) break;
    }
    ll res=(ull)n*(n+1ll)/2ll-_res;
    mp_phi[n]=res; return res;
}

ll Sieve_Mu( int n )
{
    if ( n<=N-10 ) return mu[n];
    if ( mp_mu[n] ) return mp_mu[n];
    ll res=1;
    for ( int l=2,r; l<=n; l=r+1 )
    {
        r=n/(n/l); res-=(r-l+1)*Sieve_Mu(n/l);
        if ( r>=2147483647 ) break;
    }
    return mp_mu[n]=res;
}

int main()
{
    Sieve();
    int T=read();
    while ( T-- )
    {
        int n=read();
        printf("%lld %lld\n",Sieve_Phi(n),Sieve_Mu(n) );
    }
    return 0;
}

例題

毒瘤的數學題

\[\left(\sum\limits_{i=1}^n\sum\limits_{j=1}^nij\gcd(i,j)\right)\bmod p \]

\(1\le n\le 10^{10},5\times10^{8}\le p\le 1.1\times 10^{9}\) .

Solution

難得用上五級標題

\[\begin{aligned} \sum_{i=1}^n\sum_{j=1}^nij\gcd(i,j) &=\sum_{d=1}^n d^3\sum_{i=1}^{\lfloor n/d\rfloor}\sum_{j=1}^{\lfloor n/d\rfloor} ij[\gcd(i,j)=1]\\\\ &=\sum_{d=1}^nd^3\sum_{i=1}^{\lfloor n/d\rfloor}\sum_{j=1}^{\lfloor n/d\rfloor}ij\sum_{k|d}\mu(k)\\\\ 令S(n)=\frac{n*(n+1)}2,&=\sum_{d=1}^nd^3\sum_{k|d}\mu(k)k^2S(\lfloor \frac{n}{dk}\rfloor)^2\\\\ &=\sum_{d=1}^nd^3\sum_{k=1}^{\lfloor n/d\rfloor}\mu(k)k^2S(\lfloor \frac{n}{dk}\rfloor)^2 \end{aligned} \]

\(T=dk\) 進行代換,得到:

\[\begin{aligned} \sum_{d=1}^nd^3\sum_{k=1}^{\lfloor n/d\rfloor}\mu(k)k^2S(\lfloor \frac{n}{dk}\rfloor)^2 &=\sum_{d=1}^nd^3\sum_{k=1}^{\lfloor n/d\rfloor}\mu(k)k^2S(\lfloor \frac{n}{T}\rfloor)^2\\\\ &=\sum_{T=1}^nS(\lfloor \frac nT\rfloor)^2\sum_{d|T}d^3\mu\Big(\frac{T}{d}\Big)\Big(\frac{T}{d}\Big)^2\\\\ &=\sum_{T=1}^nS(\lfloor \frac nT\rfloor)^2T^2\sum_{d|T}\mu\Big(\frac{T}{d}\Big)\cdot d\\\\ \end{aligned} \]

根據莫反中的前置知識,有

\[\varphi(x)=(\mu*ID)(x)=\sum\limits_{d|x}d\cdot\mu(\frac xd)\\\\ \]

所以:

\[=\sum_{T=1}^nS(\lfloor \frac{n}{T}\rfloor )^2T^2\varphi(T)\\\\ \]

然後就可以愉快地杜教篩了。把模板掏下來以便觀察:

\[g(1)S(n)=\sum_{i=1}^nh(i)-\sum_{d=2}^ng(d)\cdot S(\lfloor \frac{n}{d}\rfloor )(杜教篩)\\\\ (f*g)(x)=\sum_{d|x}f(d)g(\frac{x}{d})(Dirichlet卷積)\\\\ \varphi*1=ID(和\varphi 有關的前置式子) \]

\(sum(n)=\sum_{i=1}^ni^2\varphi(i)\) ,顯然我們需要把 \(sum(n)\) 中這個 \(i^2\) 消掉。看著卷積式子,考慮令 \(g(n)=n^2\) .

於是有:

\[(f*g)(n)=\sum_{d|n}d^2\varphi(d)\Big(\frac{n}{d}\Big)^2=n^2\sum_{d|n}\varphi(d) \]

然後再把最後一個式子給摁上去:

\[(f*g)(n)=n^3. \]

那麼就有:

\[sum(n)=\sum_{i=1}^ni^3-\sum_{d=2}^nd^2sum(\lfloor \frac{n}{d}\rfloor)\\\\ \]

對於原式:

\[\sum_{T=1}^nS(\lfloor \frac{n}{T}\rfloor )^2T^2\varphi(T)\\\\ \]

就可以杜教篩解決了。

Code
//Author: RingweEH
const int N=5e6+10;
int pri[N],tot=0,phi[N],sum[N],Mod;
ll n,inv6,inv2;
bool is[N];

int power( int a,int b )
{
	int res=1;
	for ( ; b; b>>=1,a=1ll*a*a%Mod )
		if ( b&1 ) res=1ll*res*a%Mod;
	return res;
}

void Sieve()
{
	phi[1]=1; is[1]=1;
	for ( int i=2; i<=(N-10); i++ )
	{
		if ( !is[i] ) pri[++tot]=i,phi[i]=i-1;
		for ( int j=1; j<=tot && (i*pri[j]<=(N-10)); j++ )
		{
			int pos=i*pri[j]; is[pos]=1;
			if ( i%pri[j]==0 ) { phi[pos]=1ll*phi[i]*pri[j]%Mod; break; }
			else phi[pos]=1ll*phi[i]*phi[pri[j]]%Mod;
		}
	}
	for ( int i=1; i<=(N-10); i++ )
		sum[i]=(1ll*i*i%Mod*phi[i]%Mod+sum[i-1])%Mod;
}

int pow1( ll x )
{
	x%=Mod; int res=x*(x+1)%Mod*inv2%Mod;
	return res;
}

int pow2( ll x )
{
	x%=Mod; int res=1ll*x*(x+1)%Mod*(2*x+1)%Mod*inv6%Mod;
	return res;
}

int pow3( ll x )
{
	int res=pow1(x); res=1ll*res*res%Mod;
	return res;
}

unordered_map<ll,int> mp;
int Sieve_Mu( ll x )
{
	if ( x<=(N-10) ) return sum[x];
	if ( mp[x] ) return mp[x];
	int res=pow3(x);
	for ( ll l=2,r; l<=x; l=r+1 )
	{
		r=x/(x/l); res=(res-1ll*(pow2(r)-pow2(l-1))*Sieve_Mu(x/l)%Mod)%Mod;
	}
	mp[x]=(res+Mod)%Mod; return mp[x];
}

void init()
{
	Sieve(); inv2=power( 2ll,Mod-2 ); inv6=power( 6ll,Mod-2 );
}

int main()
{
	Mod=read(),n=read(); init();

	ll ans=0;
	for ( ll l=1,r; l<=n; l=r+1 )
	{
		r=n/(n/l); ans=(ans+1ll*(Sieve_Mu(r)-Sieve_Mu(l-1))*pow3(n/l)%Mod)%Mod;
	}

	printf( "%lld\n",(ans+Mod)%Mod );

	return 0;
}

迴圈之美和Min_25 篩還在咕咕咕(