O(1) 查询gcd

时间:2023-03-09 18:22:51
O(1) 查询gcd

我们来安利一个黑科技。(其实是Claris安利来的

比如我现在有一坨询问,每次询问两个不超过n的数的gcd。

n大概1kw,询问大概300w(怎么输入就不是我的事了,大不了交互库

http://mimuw.edu.pl/~kociumaka/files/stacs2013_slides.pdf

http://drops.dagstuhl.de/opus/volltexte/2013/3938/pdf/26.pdf

我们定义一个数k的一种因式分解k=k1*k2*k3为“迷之分解”当且仅当k1、k2、k3为质数或小于等于$\sqrt{k}$ 。

我们发现线筛的时候对于一个数x,设x最小的质因子为p,x/p=g,那么x的“迷之分解”可以通过g的“迷之分解”中三个数最小的一个乘上p得到。

证明似乎可以用数学归纳法证(然而我证不出来啊

然后对于每两个小于等于$\sqrt{n}$ 的数我们可以打一张gcd表出来。

最后如果我们要询问gcd(x,y),我们找到x的“迷之分解”,然后如果分解的一部分小于等于$\sqrt{n}$ 那就查表,否则那就是一个质数,分类讨论一下就行了。

伪代码:

O(1) 查询gcd

UPD:实际测试了一下随机数据跑得并没有沙茶gcd快。可能是我实现的姿势不够优越(雾

大家可以测试一下跑gcd(5702887,9227465)这个算法比沙茶gcd不知道快到哪里去了

//跑得比谁都快的gcd?
#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <algorithm>
#include <string.h>
#include <vector>
#include <math.h>
#include <time.h>
#include <limits>
#include <set>
#include <map>
using namespace std;
const int N=;
const int sn=sqrt(N);
bool np[N+];
int ps[N+],pn=;
int cs[N+][];
void xs()
{
np[]=cs[][]=cs[][]=cs[][]=;
for(int i=;i<=N;i++)
{
if(!np[i]) {cs[i][]=cs[i][]=; cs[i][]=i; ps[++pn]=i;}
for(int j=;j<=pn&&i*ps[j]<=N;j++)
{
np[i*ps[j]]=;
int cm=cs[i][]*ps[j];
if(cm<cs[i][])
{
cs[i*ps[j]][]=cm;
cs[i*ps[j]][]=cs[i][];
cs[i*ps[j]][]=cs[i][];
}
else if(cm<cs[i][])
{
cs[i*ps[j]][]=cs[i][];
cs[i*ps[j]][]=cm;
cs[i*ps[j]][]=cs[i][];
}
else
{
cs[i*ps[j]][]=cs[i][];
cs[i*ps[j]][]=cs[i][];
cs[i*ps[j]][]=cm;
}
if(i%ps[j]);else break;
}
}
}
int gcdd[sn+][sn+];
void smgcd()
{
for(int i=;i<=sn;i++) gcdd[i][]=gcdd[][i]=i;
for(int i=;i<=sn;i++)
{
for(int j=;j<=i;j++) gcdd[i][j]=gcdd[j][i]=gcdd[i-j][j];
}
}
void pre_gcd() {xs(); smgcd();}
int gcd(int a,int b)
{
if(a>N||b>N)
{
puts("Fuck You\n");
return -;
}
int *x=cs[a],g=;
for(int i=;i<;i++)
{
int d;
if(x[i]<=sn) d=gcdd[x[i]][b%x[i]];
else if(b%x[i]) d=;
else d=x[i];
g*=d; b/=d;
}
return g;
}
int euclid_gcd(int x,int y)
{
while(y)
{
int t=x%y; x=y; y=t;
}
return x;
}
int tmd=-;
void gc()
{
if(tmd==-) tmd=clock();
else
{
printf("Passed: %dms\n",clock()-tmd);
tmd=-;
}
}
int main()
{
int seed=time();
//1kw个随机数测试
int ans;
printf("Euclid gcd...\n");
srand(seed);
gc();
ans=;
for(int i=;i<=;i++)
{
int a=(rand()*+rand())%N+,b=(rand()*+rand())%N+;
ans^=euclid_gcd(a,b);
}
printf("Ans = %d\n",ans);
gc();
printf("New gcd...\n");
srand(seed);
gc();
pre_gcd();
ans=;
for(int i=;i<=;i++)
{
int a=(rand()*+rand())%N+,b=(rand()*+rand())%N+;
ans^=gcd(a,b);
}
printf("Ans = %d\n",ans);
gc();
}