题目概述

求 $\sum_{i=1}^{n}\sum_{j=1}^{n}f^K[gcd(i,j)]$ ，其中 $f(n)$ 表示 $n$ 的次大质因子（相同质因子算多次），特殊的，$f(1)=0,f(p)=1$ 。

解题报告

我菜死了……刚开始写了个除法分块套除法分块+杜教筛+Min_25筛然后就果断TLE了。

先用莫比乌斯函数常规操作一波，可以得到：

$$ \sum_{T=1}^{n}\lfloor{n\over T}\rfloor^2\sum_{d|T}f^K(d)\mu({T\over d}) $$

令 $g=f^K*\mu$ ，只要能快速求 $g$ 的前缀和，就可以进行除法分块了。

因为带 $\mu$ ，所以可以考虑卷上 $1$ 然后杜教筛：

$$ \sum_{i=1}^{n}f^K(i)=\sum_{i=1}^{n}(g*1)(i)=\sum_{i=1}^{n}S(\lfloor{n\over i}\rfloor)\\ S(n)=\sum_{i=1}^{n}f^K(i)-\sum_{i=2}^{n}S(\lfloor{n\over i}\rfloor) $$

但是 $\sum_{i=1}^{n}f^K(i)$ 怎么求啊？和这题一样，可以用Min_25筛。

由于进行的是 $n$ 的除法分块，所以Min_25筛只需要预处理一次，不需要每次杜教筛都重新预处理。

把杜教筛和Min_25筛都记忆化一下复杂度大概就对了，但是我还是算不来……

示例程序

#include<cmath>
#include<cstdio>
using namespace std;
typedef unsigned int uint;typedef long long LL;
const int maxs=89443,maxt=2e6,HA=19260817;

int N,K,p[maxs+5];uint pp[maxs+5];bool pri[maxs+5];
int n,S,m,ID[2][maxs+5],a[maxs+5];uint g[maxs+5],ans;
struct Hashmap{
    int E,lnk[HA],nxt[maxt+5];uint w[maxt+5];LL son[maxt+5];
    inline void Add(int x,LL y) {son[++E]=y;w[E]=0;nxt[E]=lnk[x];lnk[x]=E;}
    inline int Count(LL x) {for (int j=lnk[x%HA];j;j=nxt[j]) if (son[j]==x) return 1;return 0;}
    inline uint& operator [] (const LL &x){
        for (int j=lnk[x%HA];j;j=nxt[j]) if (son[j]==x) return w[j];
        Add(x%HA,x);return w[E];
    }
}F,G;

inline uint Pow(uint w,int b) {uint s;for (s=1;b;b>>=1,w=w*w) if (b&1) s=s*w;return s;}
inline void Make(int n){
    for (int i=2;i<=n;i++){
        if (!pri[i]) p[++p[0]]=i,pp[p[0]]=Pow(i,K);
        for (int j=1,t;j<=p[0]&&(t=i*p[j])<=n;j++)
            {pri[t]=true;if (!(i%p[j])) break;}
    }
}
#define ID(x) ((x)<=S?ID[0][x]:ID[1][n/(x)])
uint Sum(int a,int b){
    if (a<=1||p[b]>a) return 0;if (G.Count((LL)a*p[0]+b-1)) return G[(LL)a*p[0]+b-1];
    int k=ID(a);uint ans=(b>1?pp[b-1]*(g[k]-(b-1)):0);
    for (int i=b;i<=p[0]&&(LL)p[i]*p[i]<=a;i++){
        LL A=p[i],B=A*A;
        for (int j=1;B<=a;j++,A=B,B*=p[i])
            ans+=Sum(a/A,i+1)+pp[i];
    }return G[(LL)a*p[0]+b-1]=ans;
}
inline void Min25(int x){
    n=x;S=sqrt(n);m=0;
    for (int l=1,r;l<=n;l=r+1){
        r=n/(n/l);a[++m]=n/l;g[m]=a[m]-1;
        a[m]<=S?ID[0][a[m]]=m:ID[1][n/a[m]]=m;
    }
    for (int j=1;j<=p[0]&&p[j]<=S;j++)
        for (int i=1;i<=m&&(LL)p[j]*p[j]<=a[i];i++)
            g[i]+=j-1-g[ID(a[i]/p[j])];
}
uint Sum(int x){
    if (F.Count(x)) return F[x];uint ans=Sum(x,1)+g[ID(x)];
    for (int l=2,r;l<=x;l=r+1) r=x/(x/l),ans-=Sum(x/l)*(r-l+1);return F[x]=ans;
}
int main(){
    freopen("program.in","r",stdin);freopen("program.out","w",stdout);
    scanf("%d%d",&N,&K);Make(maxs);Min25(N);
    for (int l=1,r;l<=N;l=r+1) r=N/(N/l),ans+=(Sum(r)-Sum(l-1))*(N/l)*(N/l);
    printf("%u\n",ans);return 0;
}