2019牛客多校第7场 K.Function（min_25筛）

最新推荐文章于 2023-03-10 23:03:32 发布

Flyppy_White

最新推荐文章于 2023-03-10 23:03:32 发布

阅读量551

点赞数 1

分类专栏：数论

本文链接：https://blog.csdn.net/baiyifeifei/article/details/98894560

版权

数论专栏收录该内容

49 篇文章 0 订阅

订阅专栏

2019牛客多校第7场 K.Function（min_25筛）

题目大意

定义一个函数
$csl(p,x)=\begin{cases} 3e+1&x=p^e\ and\ \exist a,b\ p=a^2+b^2\\ 1&x=p^e\ and\ \not\exist a,b\ p=a^2+b^2\\ 0&other \end{cases}$
再定义 $tl(p,x)=max_{d|x}(csl(p,d))$

现在要求求
$\sum_{d=1}^n\prod_ptl(p,d)$

解题思路

很显然 $f(d)=\prod_ptl(p,d)$ 是个积性函数。

这个问题看数据范围显然要用亚线性筛，且函数函数显然也很难构造卷积，因此我们选择使用min25筛。为此这个问题的难点就剩下了求函数在为素数时的前缀和需要解决。根据二平方定理，我们可以知道一个素数如果想要可以分解成两个平方数之和，需要满足其为2或mod4余1.据此做DP求出函数在素数部分的前缀和即可

AC代码

#include<bits/stdc++.h>
using namespace std;
#define int long long
typedef long long LL;
int s,n;
const int size=1e5+5;
int pre[size][4],p[size];
bool prime[size];
int tot;
int num[size];
int h[size][4];
int tol;
inline int ID(int x){return x<=s?x:tol-n/x+1;}
inline int f(int p,int e)
{
    if(p%4==1) return 3*e+1;
    return 1;
}
LL get_g(int x,int k)
{
    if(x<=1||p[k]>x) return 0;
    LL ans=h[ID(x)][3]-pre[k-1][3]+(h[ID(x)][1]-pre[k-1][1])*4;
    if(k==1) ans++;//'2' is specil
    while(1LL*p[k]*p[k]<=x&&k<=tot)
    {
        int pk=p[k],t=p[k];
        for(int e=1;t*pk<=x;e++,t=t*pk)
            ans=ans+f(pk,e)*(get_g(x/t,k+1))+f(pk,e+1);
        k++;
    }
    return ans;
}  
void get_h(int n)
{
    s=sqrt(n);
    tot=0;
    while(s*s<=n) s++;
    while(s*s>n) s--;
    for(int j=0;j<4;j++) pre[0][j]=0;
    for(int i=1;i<=s;i++) prime[i]=true;
    for(int i=2;i<=s;i++)
    {
        if(prime[i])
        {
            p[++tot]=i;
            for(int j=0;j<4;j++) pre[tot][j]=pre[tot-1][j]+(i%4==j);
        }
        for(int j=1;j<=tot&&p[j]*i<=s;j++)
        {
            prime[i*p[j]]=false;
            if(i%p[j]==0) break;
        }
    }
    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;
    for(int i=1;i<=tol;i++)
    {
        h[i][1]=num[i]/4+(num[i]%4>=1)-1;
        h[i][3]=num[i]/4+(num[i]%4>=3);
        h[i][2]=num[i]/4+(num[i]%4>=2);
        h[i][0]=num[i]/4;
    }
    int x=1;
    for(int j=1;j<=tot;j++)
    {
        while(num[x]<p[j]*p[j]) x++;
        for(int i=tol;i>=x;i--)
        {
            for(int r=0;r<4;r++)
            {
                h[i][r*p[j]%4]=h[i][r*p[j]%4]-(h[ID(num[i]/p[j])][r]-pre[j-1][r]);
            }
        }
    }
}
int32_t main()
{
    int t;
    //freopen("k.in","r",stdin);
    scanf("%lld",&t);
    while(t--)
    {
        scanf("%lld",&n);
        get_h(n);
        printf("%lld\n",get_g(n,1)+1);
    }
}