HDU6607Easy Math Problem（min_25筛+杜教筛+拉格朗日插值）

最新推荐文章于 2021-09-14 21:45:13 发布

原创最新推荐文章于 2021-09-14 21:45:13 发布 · 839 阅读

1 ·

CC 4.0 BY-SA版权

数论专栏收录该内容

49 篇文章

订阅专栏

本文详细解析了HDU6607题目，采用min_25筛、杜教筛及拉格朗日插值等算法，解决了一个涉及素数、最大公约数、最小公倍数的复杂数学问题。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

HDU6607Easy Math Problem（min_25筛+杜教筛+拉格朗日插值）

题目大意

求
$\sum_{i=1}^n\sum_{j=1}^ngcd(i,j)^klcm(i,j)[gcd(i,j)\in prime]$

其中 $1≤n≤1010,1≤k≤1001\le n\le 10^{10},1\le k\le 100$

解题思路

$\begin{aligned} \sum_{i=1}^n\sum_{j=1}^ngcd(i,j)^klcm(i,j)[gcd(i,j)\in prime]&=\sum_p\sum_{i=1}^{\left\lfloor\frac{n}{p}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{p}\right\rfloor}p^{k+1}\cdot(i\cdot j)[gcd(i,j)=1]\\ &=\sum_pp^{k+1}\sum_{i=1}^{\left\lfloor\frac{n}{p}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{p}\right\rfloor}\cdot(i\cdot j)[gcd(i,j)=1] \end{aligned}$

对于 $∑i=1⌊np⌋∑j=1⌊np⌋⋅(i⋅j)[gcd(i,j)=1]\sum_{i=1}^{\left\lfloor\frac{n}{p}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{p}\right\rfloor}\cdot(i\cdot j)[gcd(i,j)=1]$ 单独分析
$\sum_{i=1}^n\sum_{j=1}^nj\cdot i[gcd(i,j)=1]$

可以对原本的区间进行拆分
$[1\le i\le n][1\le j\le n]=[1\le i\le n][1\le j< i]+[1\le i<j][1\le j\le n]+[i=j][1\le i\le n]$
对第一部分求和
$\sum_{i=1}^ni\sum_{j=1}^{i-1}j [gcd(i,j)]$
根据公式 $KaTeX parse error: \tag works only in display equations$ 就有
$\sum_{i=1}^n\frac{i[i=1]+i^2\varphi(i)}{2}$
而区间的第二部分和第一部分实质上是一样的，那么前面两部分就是 $∑i=1ni[i=1]+i2φ(i)\sum_{i=1}^ni[i=1]+i^2\varphi(i)$ 而第三部分除了 $i = j = 1$ 的情况以外，其总为0，而当其为1时答案恰好为1，故最终的答案就是
$\sum_{i=1}^n\sum_{j=1}^nj\cdot i[gcd(i,j)=1]=\sum_{i=1}^ni^2\varphi(i)$
而对于 $id2⋅φid^2\cdot \varphi$ 这个函数可以通过与 $id^2$ 卷积
$(id^2\cdot\varphi)*(id^2)=\sum_{d|n}d^2\varphi(d)\frac{n^2}{d^2}=n^2\sum_{d|n}\varphi(d)=n^3$
其中 $n^3$ 和 $n^2$ 都有公式可以 $O (1)$ 求出，因此可以很方便地应用杜教筛求解

而对于 $∑ppk+1\sum_pp^{k+1}$ 可以通过min25筛中的子方案求出所有的整除分块中需要用到的前缀和

这个问题的主要复杂度就是通过拉格朗日插值处理出前 $n\sqrt n$ 的 $∑i=1nik+1\sum_{i=1}^ni^{k+1}$ 的前缀和，这部分的复杂度是 $O(kn)O(k\sqrt n)$

加上min25筛加上杜教筛的复杂度，总的复杂度就是 $O(kn+n34logn+n23)O(k\sqrt n+\frac{n^{\frac{3}{4}}}{logn}+n^{2\over 3})$

AC代码

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int size=2e5+5;
const int mod=1e9+7;
typedef long long LL;
bool prime[size]; int p[size];
int tot=0;
int inv6,inv2;
int phi[size];
int sumiiphi[size];
int quick_pow(int a,int b)
{
    register int ans=1;
    while(b)
    {
        if(b&1) ans=1LL*ans*a%mod;
        a=1LL*a*a%mod;
        b>>=1;
    }
    return ans;
}
int inv[205];
void init()
{
    tot=0;
    phi[1]=1;
    for(int i=1;i<size;i++) prime[i]=true;
    for(int i=2;i<size;i++)
    {
        if(prime[i])
            p[++tot]=i,phi[i]=i-1;
        for(int j=1;j<=tot&&p[j]*i<size;j++) 
        {
            prime[i*p[j]]=false;
            if(i%p[j]==0) 
            {
                phi[i*p[j]]=phi[i]*p[j];
                break;
            }
            else
            phi[i*p[j]]=phi[i]*phi[p[j]];                
        }
    }
    sumiiphi[0]=0;
    for(int i=1;i<size;i++) sumiiphi[i]=(sumiiphi[i-1]+1LL*i*i%mod*phi[i]%mod)%mod;
    inv6=quick_pow(6,mod-2);
    inv2=quick_pow(2,mod-2);
    inv[0]=inv[1]=1;
    for(int i=2;i<=200;i++)
	inv[i] = (mod - (mod / i)) * inv[mod % i] % mod;
	for (int i = 1; i <= 200; ++i)inv[i] = inv[i] * inv[i - 1] % mod;
}
int tol;
int hk[size<<1];LL num[size<<1];
int pre[size];
int coeff[205];
void pre_lage(int k)
{
    for(int i=1;i<=k+2;i++) coeff[i]=quick_pow(i,k);
    coeff[0]=0;
    for(int i=1;i<=k+2;i++) coeff[i]=(coeff[i-1]+coeff[i])%mod;
}
LL suf[205],bef[205];
int lage(int n,int k)
{
    if(n<=k+2) return coeff[n];
    bef[0] = suf[0] = 1;
	for (int i = 1; i <= k + 2; ++i) {
		bef[i] = bef[i - 1] * ((n - i) % mod) % mod;
		suf[i] = suf[i - 1] * ((n + i - k - 3) % mod) % mod;
	}
	LL res = 0;
	for (int i = 1; i <= k + 2; ++i) {
		LL s = coeff[i] * bef[i - 1] % mod * suf[k - i + 2] % mod * inv[i - 1] % mod * inv[k + 2 - i] % mod;
		if ((k + 2 - i) & 1)s = -s;
		res += s;
		res = (res % mod + mod) % mod;
	}
	return res;
}
int s;LL n;
int id1[size],id2[size];
inline int ID(LL x){return x<=s?x:tol-n/x+1;}
void get_h(int k)
{
    s=sqrt(n);
    while(1LL*s*s<=n) s++;
    while(1LL*s*s>n) s--;
    pre[0]=0;
    for(register int i=1;p[i]<=s;i++)
    {
        pre[i]=(1LL*pre[i-1]+quick_pow(p[i],k+1))%mod;
    }
    tol=0;
    for(int i=1;i<=s;i++) num[++tol]=i;
    for(int i=s;i>=1;i--) if(n/i>s) num[++tol]=n/i;
    pre_lage(k+1);
    for(register int i=1;i<=tol;i++) hk[i]=(lage(num[i]%mod,k+1)-1+mod)%mod;
    hk[0]=0;
    int x=1;
    for(int j=1;j<=tot&&p[j]<=s;j++)
    {
        while(num[x]<p[j]*p[j]) x++;
        for(int i=tol;i>=x;i--)
        {
            hk[i]=((hk[i]-1LL*(pre[j]-pre[j-1]+mod)%mod*(hk[ID(num[i]/p[j])]-pre[j-1]))%mod+mod)%mod;
        }
    }
}
unordered_map<LL,int> mp;
inline int sum2(int n)
{
    return 1LL*n*(n+1)%mod*(2*n+1)%mod*inv6%mod;
}
inline int sum3(int n)
{
    int ans=1LL*n*(n+1)%mod*inv2%mod;
    return 1LL*ans*ans%mod;
}
inline int S(LL n)
{
    if(n<size) return sumiiphi[n];
    if(mp.count(n)) return mp[n];
    LL ans=sum3(n%mod);//n取模 
    for(LL l=2,r;l<=n;l=r+1)
    {
        r=n/(n/l);
        ans=((ans-1LL*S(n/r)*(sum2(r%mod)-sum2((l-1)%mod)))%mod+mod)%mod;
    }
    return mp[n]=ans;
}
int32_t main()
{
    int t;
    scanf("%lld",&t);
    int k;
    init();
    while(t--)
    {
        scanf("%lld%lld",&n,&k);
        mp.clear();
        get_h(k);
        int ans=0;
        for(int i=1;i<=tol;i++)
        {
            ans=(ans+1LL*(hk[i]-hk[i-1])%mod*S(n/num[i])%mod+mod)%mod;
        }
        printf("%lld\n",ans);
    }
}