bzoj 2693: jzptab （反演）

题目描述

传送门

题目大意：
∑i=1n∑j=1mlcm(i,j) mod 100000009

题解

哎，虽然式子很恶心，但是还是要化的 ,设n<=m
∑i=1n∑j=1mlcm(i,j)
∑i=1n∑j=1mi∗jgcd(i,j)
∑k=1n∑i=1n∑j=1mi∗jk[gcd(i,j)=k]
∑k=1n1k∑i=1⌊nk⌋∑j=1⌊mk⌋i∗k∗j∗k[gcd(i,j)=1]
∑k=1nk∑i=1⌊nk⌋∑j=1⌊mk⌋i∗j∑d|(i,j)μ(d)
∑k=1nk∑nd=1μ(d)∑i=1⌊nk⌋[d|i]∗i∑j=1⌊mk⌋[d|j]∗j
∑k=1nk∑nd=1μ(d)∑i=1⌊nk∗d⌋d∗i∑j=1⌊mk∗d⌋d∗j
设sum(x,y)=∑i=1xi∑j=1yj 这个式子直接用等差数列求和公式计算即可。
那么上面的式子就是

∑k=1nk∑d=1nμ(d)∗d2sum(⌊nk∗d⌋,⌊mk∗d⌋)

然后设T=k∗d ,转变枚举顺序

∑T=1nsum(⌊nT⌋,⌊mT⌋)∗∑d|nμ(Td)∗Td2∗d

那么问题的关键就转换成了如何快速的求解h(n)=∑d|nμ(Td)∗Td2∗d
这个函数是一个积性函数，可以线性筛啊。a,b互质的时候，直接相乘。
a,b不互质的时候，直接乘上质数。h(pixi)=μ(pi)∗pi2∗pixi−1+μ(1)∗pixi

代码

#include<iostream>#include<cstdio>#include<algorithm>#include<cstring>#include<cmath>#define LL long long #define p 100000009#define N 10000000using namespace std;int n,m,T,pd[N+3],prime[N+3],g[N+3];LL f[N+3],inv;void init(){    f[1]=1;     for (int i=2;i<=N;i++) {        if (!pd[i]) {            prime[++prime[0]]=i;            f[i]=(LL)(i-(LL)i*i%p)%p;         }        for (int j=1;prime[j]*i<=N;j++) {            pd[i*prime[j]]=1;            if (i%prime[j]==0) {                f[i*prime[j]]=(LL)f[i]*prime[j]%p;                break;            }            f[i*prime[j]]=(LL)f[i]*f[prime[j]]%p;        }    }    for (int i=1;i<=N;i++) f[i]=(f[i-1]+f[i])%p;}LL calc(LL x,LL y){    LL t1=((x+1)*x>>1)%p; LL t2=((y+1)*y>>1)%p;    return t1*t2%p;}int main(){    freopen("a.in","r",stdin);    freopen("my.out","w",stdout);    init(); scanf("%d",&T);    inv=g[4];    while (T--) {        scanf("%d%d",&n,&m);        if (n>m) swap(n,m);        int j=0; LL ans=0;        for (int i=1;i<=n;i=j+1) {            j=min(n/(n/i),m/(m/i));            ans=(ans+calc(n/i,m/i)*(f[j]-f[i-1])%p)%p;        }        printf("%lld\n",(ans%p+p)%p);    }}