关于复杂度优于线性的质数筛
  • 板块学术版
  • 楼主donghanwen1225
  • 当前回复13
  • 已保存回复13
  • 发布时间2022/1/30 17:21
  • 上次更新2023/10/28 10:02:28
查看原帖
关于复杂度优于线性的质数筛
153687
donghanwen1225楼主2022/1/30 17:21

这篇帖子中,@wangchenyi 提出了一种质数筛法,但其复杂度并没有优于线性,实际上大致在欧拉筛和埃氏筛之间。而 @wkywkywky 给出的这篇论文中提出了一种时间复杂度为 O(nloglogn)O\left(\dfrac{n}{\log\log n}\right),最低空间复杂度为 O(N0.5+o(1))O(60N)O(N^{0.5+o(1)})\approx O(60\sqrt N) 的筛法。然后蒟蒻就实现了一下这个科技,并且其时间的确要优于一般的欧拉筛。
实际上我在实现该算法的时候,因为模板要求输出第 kk 个质数,所以需要开一个 int 数组去存储所有质数;并且为了不对这些质数排序而浪费时间,把空间开到了 O(16×N0.5+o(1))O(960N)O(16\times N^{0.5+o(1)})\approx O(960\sqrt N) 级别。
目前在未经任何卡常的情况下,线性筛模板单点时间在 350ms\text{350ms} 左右,空间在 33M33M 左右。相比一般的欧拉筛,无论是时间还是空间上都有较大的优化。

代码:

#include<iostream>
#include<cstdio>
#include<ctime>
#include<cstring>
using namespace std;
const int INF=2000000000;
const int N=100000000,B=560000;
int del1[9]={0,1,13,17,29,37,41,49,53};
int del2[5]={0,7,19,31,43};
int del3[5]={0,11,23,47,59};
int cnt,pri[6000001];
int pre[2001];bool ycl[11001];
bool bj1[9][B+5],bj2[5][B+5],bj3[5][B+5];
double be_time,en_time;

template<typename T> void Read(T &x)
{
	x=0;bool f=1;
	char ch=getchar();
	while(ch<'0'||ch>'9') ch=getchar();
	while(ch>='0'&&ch<='9')
	{
		x=x*10+(ch^48);
		ch=getchar();
	}
	x=f==1?x:-x;
}
template<typename T> void Put(T x)
{
	if(x<10) putchar(x+'0');
	else Put(x/10),putchar(x%10+'0');
}

void be_speed_test(){be_time=clock();}
void en_speed_test(){en_time=clock();}
void show_time(){printf("%.3lfms\n",en_time-be_time);}

void init()
{
	int n,tmp=0;Read(n);
	for(int i=2;i<=110;i++)
		if(!ycl[i])
			for(int j=i;i*j<=11000;j++)
				ycl[i*j]=1;
	for(int i=7;i<=11000;i++)
		if(!ycl[i]) pre[++tmp]=i;
	pri[1]=2;pri[2]=3;pri[3]=5;cnt=3;
}
void get_pri(int L)
{
	for(int i=1;i<=8;i++)
	{
		int d=del1[i];
		for(int f=1;f<=15;f++)
			for(int g=1;g<=30;g++)
				if((4*f*f+g*g)%60==d)
				{
					int x=f,y0=g,k0=(4*f*f+g*g-d)/60;
					while(k0<L+B) k0+=2*x+15,x+=15;
					while(x>15)
					{
						x-=15;k0-=2*x+15;
						while(k0<L) k0+=y0+15,y0+=30;
						int k=k0,y=y0;
						while(k<L+B) bj1[i][k-L]=1-bj1[i][k-L],k+=y+15,y+=30;
					}
				}
		for(int q=1;pre[q]*pre[q]<60*(L+B);q++)
		{
			int now=pre[q]*pre[q],be=60*L/now+(60*L%now>0);
			for(int j=be*now;j<=60*(L+B);j+=now)
				if((j-d)%60==0) bj1[i][(j-d)/60-L]=0;
		}
	}
	for(int i=1;i<=4;i++)
	{
		int d=del2[i];
		for(int f=1;f<=10;f++)
			for(int g=1;g<=30;g++)
				if((3*f*f+g*g)%60==d)
				{
					int x=f,y0=g,k0=(3*f*f+g*g-d)/60;
					while(k0<L+B) k0+=x+5,x+=10;
					while(x>10)
					{
						x-=10;k0-=x+5;
						while(k0<L) k0+=y0+15,y0+=30;
						int k=k0,y=y0;
						while(k<L+B) bj2[i][k-L]=1-bj2[i][k-L],k+=y+15,y+=30;
					}
				}
		for(int q=1;pre[q]*pre[q]<60*(L+B);q++)
		{
			int now=pre[q]*pre[q],be=60*L/now+(60*L%now>0);
			for(int j=be*now;j<=60*(L+B);j+=now)
				if((j-d)%60==0) bj2[i][(j-d)/60-L]=0;
		}
	}
	for(int i=1;i<=4;i++)
	{
		int d=del3[i];
		for(int f=1;f<=10;f++)
			for(int g=1;g<=30;g++)
				if((3*f*f-g*g+6000)%60==d)
				{
					int x=f,y0=g,k0=(3*f*f-g*g-d)/60;
					while(1)
					{
						int en=0;
						while(k0>=L+B)
						{
							if(x<=y0){en=1;break;}
							k0-=(y0+15),y0+=30;
						}
						if(en) break;
						int k=k0,y=y0;
						while(k>=L&&y<x) bj3[i][k-L]=1-bj3[i][k-L],k-=(y+15),y+=30;
						k0+=x+5,x+=10;
					}
				}
		for(int q=1;pre[q]*pre[q]<60*(L+B);q++)
		{
			int now=pre[q]*pre[q],be=60*L/now+(60*L%now>0);
			for(int j=be*now;j<=60*(L+B);j+=now)
				if((j-d)%60==0) bj3[i][(j-d)/60-L]=0;
		}
	}
	for(int i=0;i<B;i++)
	{
		int cur1=1,cur2=1,cur3=1;
		int cur=60*(L+i);
		while(cur1<=8||cur2<=4||cur3<=4)
		{
			int num1=INF,num2=INF,num3=INF;
			if(cur1<=8) num1=del1[cur1];
			if(cur2<=4) num2=del2[cur2];
			if(cur3<=4) num3=del3[cur3];
			if(num1<num2&&num1<num3)
			{
				if(bj1[cur1][i]) pri[++cnt]=cur+num1,bj1[cur1][i]=0;
				cur1++;
			}
			else if(num2<num1&&num2<num3)
			{
				if(bj2[cur2][i]) pri[++cnt]=cur+num2,bj2[cur2][i]=0;
				cur2++;
			}
			else
			{
				if(bj3[cur3][i]) pri[++cnt]=cur+num3,bj3[cur3][i]=0;
				cur3++;
			}
		}
	}
}
void cal()
{
	for(int i=0;i<=N/(60*B);i++) get_pri(i*B);
}
int main()
{
	init();
//	be_speed_test();
	cal();
//	en_speed_test();
//	show_time();
	int k,x;
	Read(k);
	while(k--)
	{
		Read(x);
		Put(pri[x]);putchar('\n');
	}
	return 0;
}
2022/1/30 17:21
加载中...