zoukankan      html  css  js  c++  java
  • Min25筛求素数和

    复杂度n0.75/logn,实际上比n2/3的杜教筛要快一些。

      1 //#include <bits/stdc++.h>
      2 #include<cstdio>
      3 #include<cstring>
      4 #include<algorithm>
      5 #include<iostream>
      6 #include<string>
      7 #include<vector>
      8 #include<stack>
      9 #include<bitset>
     10 #include<cstdlib>
     11 #include<cmath>
     12 #include<set>
     13 #include<list>
     14 #include<deque>
     15 #include<map>
     16 #include<queue>
     17 
     18 using namespace std;
     19 
     20 const int N = 1000010;
     21 
     22 typedef long long LL;
     23 LL TT,nn,k;
     24 namespace Min25 {
     25 
     26     int prime[N], id1[N], id2[N], flag[N], ncnt, m;
     27 
     28     LL g[N], sum[N], a[N], T, n;
     29     inline void fff()
     30     {
     31         for(int i=0;i<=N;i++){
     32             prime[i]=0;
     33             id1[i]=0;
     34             id2[i]=0;
     35             flag[i]=0;
     36             g[i]=0;
     37             sum[i]=0;
     38             a[i]=0;
     39         }
     40         ncnt=0;
     41         m=0;
     42         T=0;
     43         n=0;
     44     }
     45     inline int ID(LL x) {
     46         return x <= T ? id1[x] : id2[n / x];
     47     }
     48 
     49     inline LL calc(LL x) {
     50         return x * (x + 1) / 2 - 1;
     51     }
     52 
     53     inline LL f(LL x) {
     54         return x;
     55     }
     56 
     57     inline void init() {
     58         T = sqrt(n + 0.5);
     59         for (int i = 2; i <= T; i++) {
     60             if (!flag[i]) prime[++ncnt] = i, sum[ncnt] = sum[ncnt - 1] + i;
     61             for (int j = 1; j <= ncnt && i * prime[j] <= T; j++) {
     62                 flag[i * prime[j]] = 1;
     63                 if (i % prime[j] == 0) break;
     64             }
     65         }
     66         for (LL l = 1; l <= n; l = n / (n / l) + 1) {
     67             a[++m] = n / l;
     68             if (a[m] <= T) id1[a[m]] = m; else id2[n / a[m]] = m;
     69             g[m] = calc(a[m]);
     70         }
     71         for (int i = 1; i <= ncnt; i++)
     72             for (int j = 1; j <= m && (LL)prime[i] * prime[i] <= a[j]; j++)
     73                 g[j] = g[j] - (LL)prime[i] * (g[ID(a[j] / prime[i])] - sum[i - 1]);
     74     }
     75 
     76     inline LL solve(LL x) {
     77         if (x <= 1) return x;
     78         return n = x, init(), g[ID(n)];
     79     }
     80 
     81 }
     82 
     83 void extend_gcd(LL a,LL b,LL &x,LL &y)
     84 {
     85     if(b==0) {
     86         x=1,y=0;
     87         return;
     88     }
     89     extend_gcd(b,a%b,x,y);
     90     LL tmp=x;
     91     x=y;
     92     y=tmp-(a/b)*y;
     93 }
     94 LL mod_inverse(LL a,LL m)
     95 {
     96     LL x,y;
     97     extend_gcd(a,m,x,y);
     98     return (m+x%m)%m;
     99 }
    100 /*int main()
    101 {
    102     LL sum;
    103     cin>>TT;
    104     while(TT--){
    105         sum=0;
    106 
    107         scanf("%lld%lld",&nn,&k);
    108         if(nn==1){
    109             printf("0
    ");
    110         }else if(nn==2){
    111             printf("6
    ");
    112         }else{
    113         LL a=mod_inverse(2,k);
    114         //cout<<a<<endl;
    115         sum=(((((nn%k)*(nn%k))%k+(nn*3)%k)%k)*a)%k;
    116         //cout<<sum<<endl;
    117         Min25::fff();
    118         sum=(sum-4+Min25::solve(nn+1))%k;
    119         printf("%lld
    ",sum%k);
    120         }
    121     }
    122     return 0;
    123 }*/
    124 int main()
    125 {
    126     int n;
    127     scanf("%d",&n);
    128     printf("%d
    ",Min25::solve(n));
    129 }
  • 相关阅读:
    ABAP Help Document(2):1.2 表达式
    ABAP Help Document(1):1.1关键字
    api——》将.doc文件转成.docx文件后缀,且仅需要输入单个文件绝对路径
    python 更改默认输出 解决编码常出错问题
    爬取法律法规代码(可直接使用)
    python datetime 模块详解
    python 获得日期列表中最大日期(能够剔出不是日期类型)
    博客园页面css
    日期大小比较令解决{strftime('%Y年%m月%d日')}出错问题
    CodeForces
  • 原文地址:https://www.cnblogs.com/St-Lovaer/p/13907017.html
Copyright © 2011-2022 走看看