51nod 1965 奇怪的式子 —— min_25筛

题目:http://www.51nod.com/Challenge/Problem.html#!#problemId=1965

推式子就同这里:https://www.cnblogs.com/yoyoball/p/9196092.html

一开始想设 \( g(n,j) = \sum\limits_{i=1}^{n} [ min(i) >= p_{j} ] f(i) \),其中 \( f(i) = d(i)  \mu(i) \) 或 \( f(i) = mu(i) \),\( d(i) \) 是质因子个数;

结果发现枚举最小质因子还好(但是太慢了会TLE啊!),改成循环来筛的话很不好赋初值啊!而且筛的过程也感觉怪怪的,弄了半天还是不对...

于是就也把 \( g(n,j) \) 设成 \( g(n,j) = \sum\limits_{i=1}^{n} [ i \in prime || min(i) >= p_{j} ] f(i) \) 了,果然就对了...

这里因为 \( d(i) \mu(i) \) 不是积性函数,所以要考虑 \( d(i) \) 少加的部分,再补回来即可(就是加一个 \( \sum \mu(i) \) );

long double 快速乘;

注意多组数据!虽然也不用清空什么数组,但因为我在过程中用到了 \( s[m+1] \) 和 \( sum[m+1] \),所以要把这两个清零。

代码如下:

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;
typedef long long ll;
typedef long double ld;
int const xn=1e6+5;
ll const mod=1e12+39,md=mod-1;
int m,pri[xn],cnt,sqr;
ll n,w[xn],h[xn],s[xn],ps[xn],sum[xn],f[xn],g[xn];
bool vis[xn];
ll upt(ll x){while(x>=mod)x-=mod; while(x<0)x+=mod; return x;}
ll upt2(ll x){while(x>=md)x-=md; while(x<0)x+=md; return x;}
/*
ll mul(ll a,ll b){ll ret=0; for(;b;b>>=1,a=upt(a+a))if(b&1)ret=upt(ret+a); return ret;}
ll mul2(ll a,ll b){ll ret=0; for(;b;b>>=1,a=upt2(a+a))if(b&1)ret=upt2(ret+a); return ret;}
*/
ll mul(ll a,ll b){ll r=(ld)a*b/mod; return a*b-r*mod;}
ll mul2(ll a,ll b){ll r=(ld)a*b/md; return a*b-r*md;}
ll pw(ll a,ll b){ll ret=1; for(;b;b>>=1,a=mul(a,a))if(b&1)ret=mul(ret,a); return ret;}
ll pw2(ll a,ll b){ll ret=1; for(;b;b>>=1,a=mul2(a,a))if(b&1)ret=mul2(ret,a); return ret;}
void init(int mx)
{
  for(int i=2;i<=mx;i++)
    {
      if(!vis[i])pri[++cnt]=i,ps[cnt]=upt2(ps[cnt-1]+i);
      for(int j=1;j<=cnt&&(ll)i*pri[j]<=mx;j++)
    {
      vis[i*pri[j]]=1;
      if(i%pri[j]==0)break;
    }
    }
}
int Id(ll x)
{
  if(x>sqr)return n/x;
  return m-x+1;
}
ll F(ll x,int y)
{
  if(pri[y]>x)return 0;
  ll ret=upt2(-h[Id(x)]+y-1);
  for(int i=y;i<=cnt&&(ll)pri[i]*pri[i]<=x;i++)
    ret=upt2(ret-F(n/pri[i],i+1));
  return ret;
}
ll G(ll x,int y)
{
  if(pri[y]>x)return 0;
  ll ret=upt2(-h[Id(x)]+y-1);
  for(int i=y;i<=cnt&&(ll)pri[i]*pri[i]<=x;i++)
    ret=upt2(ret-G(n/pri[i],i+1)-F(n/pri[i],i+1));
  return ret;
}
int main()
{
  int T; scanf("%d",&T); init(1e6);
  while(T--)
    {
      scanf("%lld",&n); sqr=sqrt(n); m=0;//
      int up; for(up=cnt;(ll)pri[up]*pri[up]>n;up--);
      for(ll i=1,j;i<=n;i=j+1)
    {
      w[++m]=n/i; j=n/w[m];
      h[m]=w[m]-1;
      if(w[m]&1)sum[m]=mul2((w[m]+1)/2,w[m]);
      else sum[m]=mul2(w[m]/2,(w[m]+1));
      s[m]=upt2(sum[m]-1);
    }
      s[m+1]=0; sum[m+1]=0;//!!!
      for(int j=1;j<=up;j++)
    for(int i=1;i<=m&&(ll)pri[j]*pri[j]<=w[i];i++)
      {
        int k=Id(w[i]/pri[j]);
        h[i]=upt2(h[i]-h[k]+(j-1));
        s[i]=upt2(s[i]-mul2(s[k],pri[j]));
        s[i]=upt2(s[i]+mul2(ps[j-1],pri[j]));//pri[j]
      }
      //for(int i=1;i<=m;i++)printf("sum[%d]=%lld ",i,sum[i]); puts("");
      //for(int i=1;i<=m;i++)printf("s[%d]=%lld ",i,s[i]); puts("");
      ll res=0;
      for(ll i=1,j;i<=n;i=j+1)
    {
      j=n/(n/i);
      res=upt2(res+mul2(sum[Id(n/i)],upt2(s[Id(j)]-s[Id(i-1)])));//s[m+1]!
    }
      for(int i=1;i<=up;i++)
    res=upt2(res-mul2(pri[i],sum[Id(n/pri[i])]));
      ll ans=pw(2,res);
      for(int i=1;i<=up;i++)
    for(ll k=1,p0=pri[i];p0<=n;p0*=pri[i],k++)
      {
        ll r=upt2(mul2(p0,sum[Id(n/p0)])-mul2(p0*pri[i],sum[Id(n/p0/pri[i])]));//sum[m+1]!
        ans=mul(ans,pw(k+1,r));
      }
      //ans=mul(ans,pw(2,G(n,1)));
      /*
      for(int i=1;i<=m;i++)
    {
      //if(w[i]>sqr)f[i]=g[i]=upt2(-(h[Id(n)]-h[i]+1));
      //else f[i]=g[i]=upt2(-(h[Id(n)]-up));
      if(h[i]<=up)f[i]=g[i]=upt2(-(h[Id(n)]-up));
      else f[i]=g[i]=upt2(-(h[Id(n)]-h[i]+1));
    }
      for(int j=up;j;j--)
    for(int i=1;i<=m&&(ll)pri[j]*pri[j]<=w[i];i++)
      {
        int k=Id(w[i]/pri[j]);
        f[i]=upt2(f[i]-f[k]-1);//mu[pj]
        g[i]=upt2(g[i]-g[k]-1-f[k]-1);
      }
      */
      for(int i=1;i<=m;i++)f[i]=g[i]=upt2(-h[i]);
      for(int j=up;j;j--)
    for(int i=1;i<=m&&(ll)pri[j]*pri[j]<=w[i];i++)
      {
        int k=Id(w[i]/pri[j]);
        f[i]=upt2(f[i]-f[k]-j);
        g[i]=upt2(g[i]-g[k]-j-f[k]-j);
      }
      ans=mul(ans,pw(2,g[Id(n)]));
      printf("%lld\n",ans);
    }
  return 0;
}

 

posted @ 2019-01-18 11:24  Zinn  阅读(187)  评论(0编辑  收藏  举报