LOJ#572. 「LibreOJ Round #11」Misaka Network 与求和

IV.LOJ#572. 「LibreOJ Round #11」Misaka Network 与求和

首先先考虑莫反推一波式子。

\begin{aligned} \sum_{i = 1}^{n} \sum_{j = 1}^{n} f^{k} (gcd (i, j)) \\ = & \sum_{d = 1}^{n} f^{k} (d) \sum_{i = 1}^{n} \sum_{j = 1}^{n} [gcd (i, j) = d] \\ = & \sum_{d = 1}^{n} f^{k} (d) \sum_{p = 1}^{n / d} μ (p) {⌊ \frac{n}{d p} ⌋}^{2} \\ = & \sum_{T = 1}^{n} {⌊ \frac{n}{T} ⌋}^{2} (f^{k} * μ) (T) \end{aligned}

$\begin{aligned}&\sum\limits_{i=1}^n\sum\limits_{j=1}^nf^k\Big(\gcd(i,j)\Big)\\=&\sum\limits_{d=1}^nf^k(d)\sum\limits_{i=1}^n\sum\limits_{j=1}^n[\gcd(i,j)=d]\\=&\sum\limits_{d=1}^nf^k(d)\sum\limits_{p=1}^{n/d}\mu(p)\left\lfloor\dfrac{n}{dp}\right\rfloor^2\\=&\sum\limits_{T=1}^n\left\lfloor\dfrac{n}{T}\right\rfloor^2(f^k*\mu)(T)\end{aligned}$

考虑整除分块掉 $T$ ，问题转换为求 $f^k*\mu$ 的前缀和，其中 $*$ 是狄利克雷卷积符号。

明显 $f^k$ 不会是积性函数，故 $f^k*\mu$ 亦不可能是。但我们的确可以类杜教筛地计算其前缀和。具体而言，即使 $f^k$ 不是积性函数，但它仍然是数论函数，则狄利克雷卷积中什么单位元与逆元的理论仍然可以应用，因而在杜教筛中我们选取的函数是 $1$ 。于是就有

S (n) = \sum_{d = 1}^{n} f^{k} (d) - \sum_{d = 2}^{n} S (\frac{n}{d})

$S(n)=\sum\limits_{d=1}^nf^k(d)-\sum\limits_{d=2}^nS(\dfrac{n}{d})$

$f^k$ 的前缀和仍可以如上min25筛掉。时间复杂度就是 $O(\text{min25筛})$ ，因为整除分块调用的 $S(x)$ 也都会是 $\dfrac{n}{k}$ 的形式，故直接非递归预处理一下就行了。

需要注意的是，本题 $f(p)=1$ 。这意味着需要在DP完 $S$ 后统一加上 $g$ 。

代码：

#include<bits/stdc++.h>
using namespace std;
typedef unsigned int ui;
int n,k,m,pri[100100],num,kth[100100],tot,sml[100100],lar[100100];
ui g[100100],f[100100],h[100100],prik[100100],res;
ui ksm(ui x,int y){ui z=1;for(;y;y>>=1,x*=x)if(y&1)z*=x;return z;}
void sieve(){for(int i=2;i<=m;i++){if(!pri[i])pri[++num]=i,prik[num]=ksm(i,k);for(int j=1;j<=num&&i*pri[j]<=m;j++){pri[i*pri[j]]=true;if(!(i%pri[j]))break;}}}
bool vis[100100];
#define ID(x) (x<=m?sml[x]:lar[n/(x)])
ui S(int x){if(!x)return 0;int X=ID(x);if(vis[X])return h[X];vis[X]=true,h[X]=f[X]+g[X];for(int l=2,r;l<=x;l=r+1)r=x/(x/l),h[X]-=S(x/l)*(r-l+1);return h[X];}
int main(){
	scanf("%d%d",&n,&k),m=sqrt(n),sieve();
	for(int i=1,I;i<=n;i=n/(n/i)+1)I=kth[++tot]=n/i,g[tot]=I-1,ID(I)=tot;
	for(int i=1;i<=num;i++)for(int j=1;j<=tot&&pri[i]*pri[i]<=kth[j];j++)g[j]-=g[ID(kth[j]/pri[i])]-(i-1);
	for(int i=num;i;i--)for(int j=1;j<=tot&&pri[i]*pri[i]<=kth[j];j++)for(int X=kth[j]/pri[i];X>=pri[i];X/=pri[i])f[j]+=f[ID(X)]+prik[i]*(g[ID(X)]-i+1);
	for(int l=1,r;l<=n;l=r+1)r=n/(n/l),res+=(S(r)-S(l-1))*(n/l)*(n/l);
	printf("%u\n",res);
	return 0;
}