【51nod 1847】奇怪的数学题

最新推荐文章于 2021-07-22 20:55:56 发布

原创最新推荐文章于 2021-07-22 20:55:56 发布 · 1.1k 阅读

1 ·

CC 4.0 BY-SA版权

题解同时被 2 个专栏收录

177 篇文章

订阅专栏

数论

31 篇文章

订阅专栏

该博客讨论了一道数学题目，涉及计算特定条件下两数的最大公约数的第二大公约数的K次方之和，并对2^32取模的结果。解题思路提及了最小质因子的概念在解决此类问题中的应用。

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

题目描述

给出 N,K ,请计算下面这个式子:
$∑_{i=1}^N∑_{j=1}^Nsgcd(i,j)^k$
其中,sgcd(i, j)表示(i, j)的所有公约数中第二大的,特殊地,如果gcd(i, j) = 1, 那么sgcd(i, j) = 0。
考虑到答案太大,请输出答案对2^32取模的结果.
1≤N≤109,1≤K≤50
样例解释:
因为gcd(i, j)=1时sgcd(i,j)=0对答案没有贡献,所以我们只考虑gcd(i,j)>1的情况.
当i是2时,j是2时,sgcd(i,j)=1,它的K次方是1
当i是2时,j是4时,sgcd(i,j)=1,它的K次方是1
当i是3时,j是3时,sgcd(i,j)=1,它的K次方是1
当i是4时,j是2时,sgcd(i,j)=1,它的K次方是1
当i是4时,j是4时,sgcd(i,j)=2,它的K次方是8
当i是5时,j是5时,sgcd(i,j)=1,它的K次方是1

解题思路

设minp(x)表示x最小的质因子（当x等于1时，minp(x)为0，当x质数时，minp(x)为1）。
于是

\sum n i = 1 \sum n j = 1 s g c d (i, j) k

$∑_{i=1}^n∑_{j=1}^nsgcd(i,j)^k$

= \sum d = 2 n d m i n p ( d ) \sum ⌊ n d ⌋ i = 1 \sum ⌊ n d ⌋ j = 1 [s g c d (i, j) = = 1]

$=\sum_{d=2}^{n}\frac{d}{minp(d)}∑^{\lfloor {\frac{n}{d}} \rfloor}_{i=1}∑^{\lfloor {\frac{n}{d}} \rfloor}_{j=1}[sgcd(i,j)==1]$

= \sum d = 2 n d m i n p ( d ) (2 \sum ⌊ n d ⌋ i = 1 φ (i) - 1)

$=\sum_{d=2}^{n}\dfrac{d}{minp(d)}(2∑^{\lfloor {\frac{n}{d}} \rfloor}_{i=1}φ(i)−1)$
对于

phi(i)phi(i) $phi(i)$ 的前缀和就可以直接杜教筛。
至于如何求出

dminp(d)dminp(d) $\dfrac{d}{minp(d)}$ ，
当

d>n−−√且d为质数d>n且d为质数 $d>\sqrt n 且 d为质数$ 时，

dminp(d)dminp(d) $\dfrac{d}{minp(d)}$ 为1。
这个就可以通过求

>n−−√>n $>\sqrt n$ 的质数来得出。
我们设

F(i,j)表示在[1,j]中，不能被前i个质数整除的数的K次方和F(i,j)表示在[1,j]中，不能被前i个质数整除的数的K次方和 $F(i,j)表示在[1,j]中，不能被前i个质数整除的数的K次方和$

H(i,j)表示在[1,j]中，不能被前i个质数整除的数的个数H(i,j)表示在[1,j]中，不能被前i个质数整除的数的个数 $H(i,j)表示在[1,j]中，不能被前i个质数整除的数的个数$

G(i,j)表示所有x∈[1,j],minp(x)≤pi,dminp(d)K的和G(i,j)表示所有x∈[1,j],minp(x)≤pi,dminp(d)K的和 $G(i,j)表示所有x∈[1,j],minp(x)≤p_i,\frac{d}{minp(d)}^K的和$
于是，我们可以的出一个递推式

F(i,j)=F(i−1,j)−pkiF(i−1,⌊jpi⌋)F(i,j)=F(i−1,j)−pikF(i−1,⌊jpi⌋) $F(i,j)=F(i-1,j)-p_i^kF(i-1,\lfloor\dfrac{j}{p_i}\rfloor)$

H(i,j)=H(i−1,j)−H(i−1,⌊jpi⌋)H(i,j)=H(i−1,j)−H(i−1,⌊jpi⌋) $H(i,j)=H(i-1,j)-H(i-1,\lfloor\dfrac{j}{p_i}\rfloor)$

G(i,j)=G(i−1,j)+F(i−1,⌊jpi⌋)G(i,j)=G(i−1,j)+F(i−1,⌊jpi⌋) $G(i,j)=G(i-1,j)+F(i-1,\lfloor\dfrac{j}{p_i}\rfloor)$
答案就是

H(n−−√以内质数个数,n)+G(n−−√以内质数个数,n)H(n以内质数个数,n)+G(n以内质数个数,n) $H(\sqrt n以内质数个数,n)+G(\sqrt n以内质数个数,n)$
然后我们发现，当

pi>jpi>j $p_i>j$ 时，F(i,j)=1,H(i,j)=1;
当

p2i>jpi2>j $p_i^2>j$ 时

F(i,j)=F(i−1,j)−pkiF(i,j)=F(i−1,j)−pik $F(i,j)=F(i-1,j)-p_i^k$

H(i,j)=H(i−1,j)−1H(i,j)=H(i−1,j)−1 $H(i,j)=H(i-1,j)-1$

G(i,j)=G(i−1,j)+1G(i,j)=G(i−1,j)+1 $G(i,j)=G(i-1,j)+1$
因为j只有

n−−√n $\sqrt n$ 种取值，我们直接预处理出每种取值，然后直接递推。

#include <cmath>
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
#include <queue>
#include <map>
const long long inf=2147483647;
const int mo=1e9+7;
const int N=1000005;
const int M=55;
using namespace std;
typedef unsigned int uint;
#define sqr(x) ((x)*(x))
long long _n[N],num,p[N],qn,n,m,SP[N],v1[N],phi[N],nphi[N],nn;
uint F[N],G[N],H[N],lu[N],mi[N],smi[N],s[M][M],S[M],ans,Smi[N];
bool bz[N];
uint poww(uint x,uint y)
{
    uint s=1;
    for(;y;y>>=1,x=x*x)
        if(y&1) s=s*x;
    return s;
}
void pre_P()
{
    mi[1]=phi[1]=1;
    for(int i=2;i<=qn;i++)
    {
        if(!bz[i])
        {
            p[++p[0]]=i,phi[i]=i-1,mi[i]=poww(i,m),smi[p[0]]=smi[p[0]-1]+mi[i];
        }
        for(int j=1;j<=p[0];j++)
        {
            int k=i*p[j];
            if(k>qn) break;
            bz[k]=true;
            mi[k]=mi[i]*mi[p[j]];
            if(i%p[j]==0)
            {
                phi[k]=phi[i]*p[j];
                break;
            }
            else phi[k]=phi[i]*(p[j]-1);
        }
    }
    for(int i=1;i<=qn;i++) Smi[i]=Smi[i-1]+mi[i],phi[i]+=phi[i-1];
}
void pre_STR()
{
    for(int i=0;i<=m;i++) s[i][i]=1;
    for(int i=1;i<=m;i++)
        for(int j=1;j<i;j++) s[i][j]=s[i-1][j-1]+s[i-1][j]*(i-1);
}
uint get_P(uint n,uint k)
{
    uint val=1;
    for(uint i=n;i>=n-k+1;i--)
        if(i%k==0) val*=i/k; 
        else val*=i;
    return val;
}
void pre()
{
    for(uint i=1,last=1;i<=n;i=last+1) last=n/(n/i),_n[++num]=n/i;
    reverse(_n+1,_n+1+num);

    pre_P(),pre_STR();
    for(int i=1;i<=num;i++)
    {
        H[i]=_n[i];
        if(_n[i]<=qn)
        {
            F[i]=Smi[_n[i]];
            continue;
        }
        S[0]=_n[i];
        for(int k=1;k<=m;k++)
        {
            S[k]=get_P(_n[i]+1,k+1);
            for(int j=0;j<k;j++) S[k]-=((k-j)&1?-1:1)*s[k][j]*S[j];
        }
        F[i]=S[m];
    }
}
uint get_F(long long i,long long j)
{
    long long k;
    k=j<=qn?j:(num-n/j+1);
    return (!i || sqr(p[i])<=j)?F[k]:(p[i]<=j?(F[k]-smi[i]+smi[lu[k]]):1);
}
uint get_H(long long i,long long j)
{
    long long k;
    k=j<=qn?j:(num-n/j+1);
    return (!i || sqr(p[i])<=j)?H[k]:(p[i]<=j?(H[k]-i+lu[k])%mo:1);
}
uint Sphi(int n)
{
    if(n<=qn) return phi[n];
    if(nphi[nn/n]) return nphi[nn/n];
    uint ans=(n&1)?((n+1)>>1)*n:(n>>1)*(n+1);
    for(int i=2,last=1;i<=n;i=last+1)
    {
        last=n/(n/i);
        ans-=Sphi(n/i)*(last-i+1);
    }
    return nphi[nn/n]=ans;
}
int main()
{
    scanf("%lld%lld",&n,&m),qn=sqrt(n);
    pre();
    for(int i=1;i<=p[0];i++)
    {
        for(int j=num;j>=1;j--)
        {
            F[j]-=mi[p[i]]*get_F(i-1,_n[j]/p[i]);
            G[j]+=get_F(i-1,_n[j]/p[i]);
            H[j]-=get_H(i-1,_n[j]/p[i]);
            lu[j]=i;
            if(i==p[0] || sqr(p[i+1])>_n[i])
                v1[j]=G[j]+H[j]-1;
            if(sqr(p[i])>_n[j-1])
            {
                v1[j]=G[j]+H[j]-1;
                break;
            }
        }
    }
    v1[2]=1;
    if(_n[3]==3) v1[3]=2;
    ans=0;
    nn=n;
    for(int i=2;i<=num;i++)
        ans=ans+(v1[i]-v1[i-1])*(2*Sphi(n/_n[i])-1);
    cout<<ans<<endl;
    return 0;
}