bzoj 4407
莫比烏斯反演
還是推式子：
設(shè)$f(n)=n^{k}$
那就是上一道題了
推的過程如下：
$\sum_{i=1}^{a}\sum_{j=1}^{b}f(gcd(i,j))$
$\sum_{i=1}^{a}\sum_{j=1}^{b}\sum_{d=1}^{min(a,b)}[gcd(i,j)\equiv d]f(d)$
$\sum_{d=1}^{min(a,b)}f(d)\sum_{i=1}^{a}\sum_{j=1}^{b}[gcd(i,j)\equiv d]$
$\sum_{d=1}^{min(a,b)}f(d)\sum_{i=1}^{\frac{a}w0obha2h00}\sum_{j=1}^{\frac{b}w0obha2h00}[gcd(i,j)\equiv 1]$
$\sum_{d=1}^{min(a,b)}f(d)\sum_{i=1}^{\frac{a}w0obha2h00}\sum_{j=1}^{\frac{b}w0obha2h00}\sum_{t=1}^{min(\frac{a}w0obha2h00,\frac{b}w0obha2h00)}\mu(t)$
$\sum_{d=1}^{min(a,b)}f(d)\sum_{t=1}^{min(a,b)}\mu(t)\frac{a}{dt}\frac{b}{dt}$
令$T=dt$，得到：
$\sum_{T=1}^{min(a,b)}\frac{a}{T}\frac{b}{T}\sum_{d|T}f(d)\mu(\frac{T}w0obha2h00)$
也就是：
$\sum_{T=1}^{min(a,b)}\frac{a}{T}\frac{b}{T}\sum_{d|T}d^{k}\mu(\frac{T}w0obha2h00)$
考慮線性篩后面那堆東西，仍然分類討論：
①.篩到的$p$與$i$互質(zhì)：
此時(shí)我們考慮增加一個(gè)$p$的貢獻(xiàn)，如果增加到$\mu$里，則原先那些直接取反
如果增加到$d^{k}$里，則相當(dāng)于原先那些乘$p^{k}$
因此$g(ip)=(p^{k}-1)g(i)$
②.篩到的$p$與$i$不互質(zhì)：
此時(shí)我們考慮增加一個(gè)$p$的貢獻(xiàn)，如果增加到$\mu$里，則原先那些仍然取反
如果增加到$f$里，則原先那些多一個(gè)$p^{k}$的貢獻(xiàn)
可...等等！
還有一種可能！
再考慮如果原先$\mu$里有一個(gè)$p$，然后增加到$f$里，此時(shí)會(huì)抵消掉取反的效果！
因此只需乘一個(gè)$p^{k}$即可
貼代碼：
#include <cstdio> #include <cmath> #include <cstring> #include <cstdlib> #include <iostream> #include <algorithm> #include <queue> #include <stack> #define ll long long using namespace std; const ll mode=1000000007; int mu[5000005]; int pri[5000005]; ll f[5000005]; ll pow_mul(ll x,ll y) { ll ret=1; while(y) { if(y&1)ret=ret*x%mode; x=x*x%mode,y>>=1; } return ret; } bool used[10000005]; int cnt=0; ll T,x,y,k; void init() { mu[1]=1; f[1]=1; for(int i=2;i<=5000000;i++) { if(!used[i])mu[i]=-1,pri[++cnt]=i,f[i]=(pow_mul(i,k)+mode-1)%mode; for(int j=1;j<=cnt&&i*pri[j]<=5000000;j++) { used[i*pri[j]]=1; if(i%pri[j]==0) { mu[i*pri[j]]=0; f[i*pri[j]]=f[i]*(f[pri[j]]+1)%mode; break; } mu[i*pri[j]]=-mu[i],f[i*pri[j]]=f[i]*f[pri[j]]%mode; } } for(int i=2;i<=5000000;i++)f[i]+=f[i-1],f[i]%=mode; } ll solve(ll a,ll b) { ll las=1,ans=0; for(int i=1;i<=a&&i<=b;i=las+1) { las=min(a/(a/i),b/(b/i)); ans+=(f[las]-f[i-1]+mode)*(a/i)%mode*(b/i)%mode; ans%=mode; } return ans; } template <typename T>inline void read(T &x) { T f=1,c=0;char ch=getchar(); while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();} while(ch>='0'&&ch<='9'){c=c*10+ch-'0';ch=getchar();} x=c*f; } int main() { read(T),read(k); init(); while(T--) { read(x),read(y); printf("%lld\n",solve(x,y)); } return 0; }
浙公網(wǎng)安備 33010602011771號(hào)