P5435 基於值域預處理的快速 GCD - 數論
阿新 • • 發佈:2021-06-29
題意
給你 \(2n\) 個數 \(a_1,a_2,\dots,a_n;b_1,b_2,\dots,b_n\),問 \(a,b\) 兩兩的 \(\gcd\)。
\(n\le 5000\),值域 \(\le 10^6\),時間限制 \(1\mathrm{s}\)。
題解
用 \(V\) 代表值域。
顯然 \(n^2\log V\) 的做法過不了。
考慮這樣一種神奇的做法:把 \(1\sim V\) 中每個整數 \(x\) 分解成三個正整數之積,且這三個數要不然 \(\le \sqrt{x}\),要不然是質數。
具體的分解方法是,找到 \(x\) 的最小質因子 \(p\) 並獲得 \(\frac{x}{p}\)
證明不會。
預處理 \(\sqrt{V}\times\sqrt{V}\) 的 \(\gcd\) 表,當分解出來的數 \(<\sqrt{V}\) 的時候直接查表,否則它一定是個質數,討論一下就行了。
程式碼
#include <cstdio> #include <cstring> #include <cctype> #include <cmath> #include <algorithm> using namespace std; #define For(Ti,Ta,Tb) for(int Ti=(Ta);Ti<=(Tb);++Ti) #define Dec(Ti,Ta,Tb) for(int Ti=(Ta);Ti>=(Tb);--Ti) template<typename T> void Read(T &x){ x=0;int _f=1; char ch=getchar(); while(!isdigit(ch)) _f=(ch=='-'?-1:_f),ch=getchar(); while(isdigit(ch)) x=x*10+(ch^48),ch=getchar(); x=x*_f; } template<typename T,typename... Args> void Read(T &x,Args& ...others){ Read(x);Read(others...); } typedef long long ll; const int Inf=0x3f3f3f3f,N=5005,V=1e6+5,SV=sqrt(V),Mod=998244353; int n,a[N],b[N]; int vis[V],prime[V],pCnt=0,k[V][3]; void Sieve(int mx){ vis[1]=1,k[1][0]=k[1][1]=k[1][2]=1; For(i,2,mx){ if(!vis[i]) prime[++pCnt]=i,k[i][0]=k[i][1]=1,k[i][2]=i; for(int j=1,*cur=k[i*prime[j]];1LL*i*prime[j]<=mx;++j,cur=k[i*prime[j]]){ vis[i*prime[j]]=1; cur[0]=k[i][0]*prime[j],cur[1]=k[i][1],cur[2]=k[i][2]; if(cur[0]>cur[1]) swap(cur[0],cur[1]); if(cur[1]>cur[2]) swap(cur[1],cur[2]); if(i%prime[j]==0) break; } } } int gcd[SV][SV]; void Pre(int mx){ Sieve(mx); For(i,0,SV-1) gcd[0][i]=gcd[i][0]=i; For(i,1,SV-1){ For(j,1,i){ gcd[i][j]=gcd[j][i]=gcd[j][i%j]; } } } int Gcd(int x,int y){ int res=1; For(i,0,2){ if(k[x][i]<SV){ int temp=gcd[k[x][i]][y%k[x][i]]; res*=temp,y/=temp; }else if(y%k[x][i]==0) res*=k[x][i],y/=k[x][i]; }return res; } int main(){ Read(n); int mx=0; For(i,1,n){ Read(a[i]);mx=max(mx,a[i]); } For(i,1,n){ Read(b[i]);mx=max(mx,b[i]); } Pre(mx); For(i,1,n){ ll cur=i,A=0; For(j,1,n){ A=(A+cur*Gcd(a[i],b[j]))%Mod; cur=cur*i%Mod; }printf("%lld\n",A); } return 0; }