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爲「迷之分解」當且僅當k一、k二、k3爲質數或小於等於$\sqrt{k}$ 。

咱們發現線篩的時候對於一個數x，設x最小的質因子爲p，x/p=g，那麼x的「迷之分解」能夠經過g的「迷之分解」中三個數最小的一個乘上p獲得。

證實彷佛能夠用數學概括法證(然而我證不出來啊

而後對於每兩個小於等於$\sqrt{n}$ 的數咱們能夠打一張gcd表出來。

最後若是咱們要詢問gcd(x,y)，咱們找到x的「迷之分解」，而後若是分解的一部分小於等於$\sqrt{n}$ 那就查表，不然那就是一個質數，分類討論一下就好了。

僞代碼:

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=10000000;
const int sn=sqrt(N);
bool np[N+5];
int ps[N+5],pn=0;
int cs[N+5][3];
void xs()
{
    np[1]=cs[1][0]=cs[1][1]=cs[1][2]=1;
    for(int i=2;i<=N;i++)
    {
        if(!np[i]) {cs[i][0]=cs[i][1]=1; cs[i][2]=i; ps[++pn]=i;}
        for(int j=1;j<=pn&&i*ps[j]<=N;j++)
        {
            np[i*ps[j]]=1;
            int cm=cs[i][0]*ps[j];
            if(cm<cs[i][1])
            {
                cs[i*ps[j]][0]=cm;
                cs[i*ps[j]][1]=cs[i][1];
                cs[i*ps[j]][2]=cs[i][2];
            }
            else if(cm<cs[i][2])
            {
                cs[i*ps[j]][0]=cs[i][1];
                cs[i*ps[j]][1]=cm;
                cs[i*ps[j]][2]=cs[i][2];
            }
            else
            {
                cs[i*ps[j]][0]=cs[i][1];
                cs[i*ps[j]][1]=cs[i][2];
                cs[i*ps[j]][2]=cm;
            }
            if(i%ps[j]);else break;
        }
    }
}
int gcdd[sn+3][sn+3];
void smgcd()
{
    for(int i=0;i<=sn;i++) gcdd[i][0]=gcdd[0][i]=i;
    for(int i=1;i<=sn;i++)
    {
        for(int j=1;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 -1;
    }
    int *x=cs[a],g=1;
    for(int i=0;i<3;i++)
    {
        int d;
        if(x[i]<=sn) d=gcdd[x[i]][b%x[i]];
        else if(b%x[i]) d=1;
        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=-1;
void gc()
{
    if(tmd==-1) tmd=clock();
    else
    {
        printf("Passed: %dms\n",clock()-tmd);
        tmd=-1;
    }
}
int main()
{
    int seed=time(0); 
    //1kw個隨機數測試 
    int ans;
    printf("Euclid gcd...\n");
    srand(seed);
    gc();
    ans=0;
    for(int i=1;i<=10000000;i++)
    {
        int a=(rand()*32768+rand())%N+1,b=(rand()*32768+rand())%N+1;
        ans^=euclid_gcd(a,b);
    }
    printf("Ans = %d\n",ans);
    gc();
    printf("New gcd...\n");
    srand(seed);
    gc();
    pre_gcd();
    ans=0;
    for(int i=1;i<=10000000;i++)
    {
        int a=(rand()*32768+rand())%N+1,b=(rand()*32768+rand())%N+1;
        ans^=gcd(a,b);
    }
    printf("Ans = %d\n",ans);
    gc();
}