Advertisement

2019年ACM-ICPC - 南京网络赛E:K Sum【反演+杜教筛】

阅读量:

题目:

2019年ICPC南京网络赛E:K Sum

题意:

给定

mall N

mall K

;求出下面又长又难看的式子的值:

mall um_{i=2}{k}f_n(i)=\sum_{i=2}{k}um_{l_1}{n}\sum_{l_2}{n}...um_{l_i}{n}(gcd(l_1,l_2,...,l_i))2

分析:

首先从简单的推一下(纯属自己娱乐):

mall um_{i=1}{N}\sum_{j=1}{N}^2
mall =um_{t=1}{N}\sum_{i=1}{eft floor rac{N}{t} ight floor}um_{j=1}^{eft floor rac{N}{t} ight floor}um_{d|gcd}t^2u
mall =um_{t=1}{N}t2um_{d=1}^{eft floor rac{N}{t} ight floor}ueft floor rac{N}{td} ight floor^2
mall =um_{l=1}^{N}eft floor rac{N}{l} ight floor2\sum_{d|l}\mu(d)\frac{l2}{d^2}

那么原式就可以化简成如下形式:

mall =um_{i=2}{k}\sum_{l=1}{n}eft floor rac{n}{l} ight floori\sum_{d|l}\mu(d)\frac{l2}{d^2}
mall =um_{l=1}^{n}rac{l2}{d2}
mall =um_{l=1}^{n}g

令最后面的一系列数值定义为函数g(L);观察该式子中的关键部分,则可以看出其内部结构构成一个等比数列;这些取整后的值被划分为若干相等的部分;可据此特性进行优化处理;进而需要迅速计算得到g(L_1)-g(L_2);令S(n)表示gL处前n项之和,则上式即等于S(L_1)-S(L_2);需要注意的是,在此过程中需特别关注这一操作所带来的影响;值得注意的是,在后续步骤中我们需要利用杜教筛这一技巧来快速求解该函数的前缀和;通过深入分析可以看出:该特定类型的积性函数与其相关联的是著名的莫比乌斯函数

u

通过与单位函数id²进行狄利克雷卷积的形式,则莫比乌斯函数与其对应的狄利克雷卷积恒等运算结果为元函数。定义h(n)为g与I的狄利克雷卷积,则有h(n)=id²(n)=n²。(这部分较为基础,在深入理解后会更加清楚)即得到:

S = um_{i=1}{N}i2-um_{d=2}^{N}S

通过上式即可实现杜教筛的部分;而等比数列的求和可通过欧拉降幂公式来简化;特别地,在公比等于1的情况下(此时求和结果为k-1),需进行模运算处理。

代码:

复制代码
 #include <bits/stdc++.h>

    
  
    
 #define sz(x) (int)(x).size()
    
 using namespace std;
    
 typedef long long LL;
    
 const int mod = 1e9+7;
    
 const int maxn = 1e6+16;
    
 const int inv6 = 166666668;
    
 map<int,int> mp;
    
 string s;
    
 int ans[maxn],mob[maxn],prime[maxn],vis[maxn];
    
 void init(){                     //O(nlogn)预处理g函数的部分前缀和
    
     mob[1] = 1;int cnt = 0;
    
     for(int i = 2;i < maxn; ++i){
    
     if(!vis[i]) prime[cnt++]=i,mob[i]=-1;
    
     for(int j = 0;j < cnt&&1ll*i*prime[j]<maxn; ++j){
    
         vis[i*prime[j]] = 1;
    
         if(i%prime[j] == 0){
    
             mob[i*prime[j]] = 0;
    
             break;
    
         }
    
         else mob[i*prime[j]] = -1*mob[i];
    
     }
    
     }
    
     for(int i = 1;i < maxn; ++i){
    
     for(int j = i;j < maxn; j+=i){
    
         ans[j] += ((1ll*j/i)*(1ll*j/i)%mod*mob[i]%mod+mod)%mod;
    
         if(ans[j] >= mod) ans[j] -= mod;
    
     }
    
     }
    
     for(int i = 1;i < maxn; ++i){
    
     ans[i] += ans[i-1];
    
     if(ans[i] >= mod) ans[i] -= mod;
    
     }
    
 }
    
 LL qpow(LL a,LL x,LL mod){
    
     LL res = 1; a%=mod;
    
     while(x){
    
     if(x & 1) res = res*a % mod;
    
     a = a*a % mod;
    
     x >>= 1;
    
     }
    
     return res;
    
 }
    
 LL inv(LL a,LL mod){
    
     return qpow(a,mod-2,mod);
    
 }
    
 LL cal1(int n){               //杜教筛
    
     if(n < maxn) return ans[n];
    
     if(mp.count(n)) return mp[n];
    
     LL res = 0;
    
     for(int i = 2,last;i <= n; i=last+1){
    
     int tep = n/i; last = n/tep;
    
     res += 1ll*(last-i+1)*cal1(tep)%mod;
    
     if(res>=mod) res -= mod; 
    
     }
    
     return mp[n] = (1ll*n*(n+1)%mod*(2*n+1)%mod*inv6%mod-res+mod)%mod;
    
 }
    
 LL cal2(LL q,LL n,LL nn){     //等比数列求和
    
     if(q == 1) return (nn-1+mod)%mod;
    
     return inv(q-1,mod)*(qpow(q,n+1,mod)-q*q%mod+mod)%mod; 
    
 }
    
 LL solve(LL n){
    
     LL res = 0,k = 0,kk = 0;
    
     for(int i = 0;i < sz(s); ++i){
    
     k = (k*10+(s[i]-'0'))%(mod-1);
    
     kk = (kk*10+(s[i]-'0'))%mod;
    
     }
    
     for(int i = 1,last;i <= n; i=last+1){
    
     int tep = n/i; last = n/tep;
    
     res += cal2(tep,k,kk)*(cal1(last)-cal1(i-1)+mod)%mod;
    
     if(res >= mod) res -= mod;
    
     }
    
     return res%mod;
    
 }
    
 int main(){
    
     int T,n; 
    
     init(); cin >> T;
    
     while(T--){
    
     cin >> n >> s;
    
     cout << solve(n) << '\n';
    
     }
    
     return 0;
    
 }

全部评论 (0)

还没有任何评论哟~