[CF gym 101741F]GCD

Description

给出n个数ai，求从中删去至多k个数之后剩余数的gcd的最大值。
n<=1e5,k<=n/2,ai<=1e18

Solution

很有趣的一道题，突破口是k<=n/2
这个启示我们什么呢？我们随机一个数x强制保留，那么我们会有至少1/2的概率找到答案。
那么现在问题变成了从剩下的n-1个数中保留一些数使得这些数和x的gcd最大
那么显然这个gcd一定是x的因数，那么我们枚举每个x的因数判断是否合法
是否合法也就是是否剩下是否有至少n-k-1个数都是这个数的倍数。
那么我们可以把剩余的n-1个数和x做gcd，然后打上tag，每个数的实际的tag就是所有是它倍数的数的tag和
这个东西分解质因数就相当于高维前缀和，直接每一维每一维搞过去就好了
然而似乎是我的Pollard_rho写错了死活WA on 68
(:з」∠)

Code

#include <cstdio>
#include <cstring>
#include <algorithm>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define fd(i,a,b) for(int i=a;i>=b;i--)
using namespace std;

typedef long long ll;
typedef long double db;

const int N=2e5+5,S=6;

int ty,n,k,tot,cnt,Cnt[N],e[N][16],id[N];
ll ans,pri[N],p[N],a[N],c[N],suf[N],d[N];

ll mult(ll x,ll y,ll Mo) {
    x%=Mo;y%=Mo;
    ll tmp=(ll)((db)x*y/Mo+1e-8)*Mo;
    return (x*y-tmp+Mo)%Mo;
}

ll pwr(ll x,ll y,ll Mo) {
    ll z=1;x%=Mo;
    for(;y;y>>=1,x=mult(x,x,Mo))
        if (y&1) z=mult(z,x,Mo);
    return z;
}

ll gcd(ll x,ll y) {
    if (!x) return 1;
    if (x<0) return gcd(-x,y);
    return y?gcd(y,x%y):x;
}

bool check(ll a,ll n,ll x,ll t) {
    ll res=pwr(a,x,n),lst=res;
    fo(i,1,t) {
        res=mult(res,res,n);
        if (res==1&&lst!=1&&lst!=n-1) return 1;
        lst=res;
    }
    if (res!=1) return 1;
    else return 0;
} 

bool Miller_Rabin(ll n) {
    if (n<2) return 0;
    if (n==2) return 1;
    if (!(n&1)) return 0;
    ll x=n-1,t=0;
    while (!(x&1)) x>>=1,t++;
    fo(i,1,S) {
        ll a=rand()%(n-1)+1;
        if (check(a,n,x,t)) return 0;
    }
    return 1;
}

ll Pollard_rho(ll n,ll t) {
    int i=1,k=2;
    ll x0=rand()%n,y=x0;
    while (1) {
        i++;
        x0=(mult(x0,x0,n)+t)%n;
        ll d=gcd(y-x0,n);
        if (d!=1&&d!=n) return d;
        if (y==x0) return n;
        if (i==k) y=x0,k<<=1;
    }  
}

void find(ll n) {
    if (n==1) return;
    if (Miller_Rabin(n)) {
        pri[++tot]=n;
        return;
    }
    ll p=n;
    while (p>=n) p=Pollard_rho(p,rand()%(n-1)+1);
    find(p);find(n/p);
}

int Id() {
    int Id=0;
    fo(i,1,cnt) Id+=id[i]*suf[i+1];
    return Id;
}

void dfs(int x,ll y) {
    if (x>cnt) {
        int now=Id();
        d[now]=y;fo(i,1,cnt) e[now][i]=id[i];
        return;
    }
    fo(i,0,c[x]) {
        id[x]=i;
        dfs(x+1,y);
        y*=p[x];
    }
}

int main() {
    scanf("%d%d",&n,&k);
    fo(i,1,n) scanf("%I64d",&a[i]);
    if (k==0) {
        ll g=a[1];
        fo(i,2,n) g=gcd(a[i],g);
        printf("%I64d\n",g);
        return 0;
    }
    srand(19260817);
    for(int s=40;s;s--) {
        int x=rand()%n+1;
        tot=0;find(a[x]);
        sort(pri+1,pri+tot+1);
        cnt=0;
        fo(i,1,tot) 
            if (i==1||pri[i]!=pri[i-1]) p[++cnt]=pri[i],c[cnt]=1;
            else c[cnt]++;
        suf[cnt+1]=1;
        fd(i,cnt,1) suf[i]=suf[i+1]*(c[i]+1);
        int tmp=suf[1];dfs(1,1);
        fo(i,0,tmp-1) Cnt[i]=0;
        fo(i,1,n) 
            if (i!=x) {
                ll g=gcd(a[i],a[x]);
                fo(j,1,cnt) {
                    id[j]=0;
                    while (!(g%p[j])) g/=p[j],id[j]++;
                }
                Cnt[Id()]++;
            }
        fo(i,1,cnt)
            fd(j,tmp-1,0) {
                if (e[j][i]==c[i]) continue;
                Cnt[j]+=Cnt[j+suf[i+1]];
            }
        fo(j,0,tmp-1) if (Cnt[j]>=n-k-1) ans=max(ans,d[j]);
    }   
    printf("%I64d\n",ans);
    return 0;
}