xudyh的gcd模板
hdu 5019
1 #include <cstdlib> 2 #include <cctype> 3 #include <cstring> 4 #include <cstdio> 5 #include <cmath> 6 #include <algorithm> 7 #include <vector> 8 #include <string> 9 #include <iostream> 10 #include <map> 11 #include <set> 12 #include <queue> 13 #include <stack> 14 #include <bitset> 15 #include <list> 16 #include <cassert> 17 #include <complex> 18 using namespace std; 19 #define rep(i,a,n) for (int i=a;i<n;i++) 20 #define per(i,a,n) for (int i=n-1;i>=a;i--) 21 #define all(x) (x).begin(),(x).end() 22 //#define fi first 23 #define se second 24 #define SZ(x) ((int)(x).size()) 25 #define TWO(x) (1<<(x)) 26 #define TWOL(x) (1ll<<(x)) 27 #define clr(a) memset(a,0,sizeof(a)) 28 #define POSIN(x,y) (0<=(x)&&(x)<n&&0<=(y)&&(y)<m) 29 typedef vector<int> VI; 30 typedef vector<string> VS; 31 typedef vector<double> VD; 32 typedef long long ll; 33 typedef long double LD; 34 typedef pair<int,int> PII; 35 typedef pair<ll,ll> PLL; 36 typedef vector<ll> VL; 37 typedef vector<PII> VPII; 38 typedef complex<double> CD; 39 const int inf=0x20202020; 40 const ll mod=1000000007; 41 const double eps=1e-9; 42 43 ll powmod(ll a,ll b) //return (a*b)%mod 44 {ll res=1;a%=mod;for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;} 45 ll powmod(ll a,ll b,ll mod) //return (a*b)%mod 46 {ll res=1;a%=mod;for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;} 47 ll gcd(ll a,ll b) //return gcd(a,b) 48 { return b?gcd(b,a%b):a;} 49 // head 50 51 namespace Factor { 52 const int N=1010000; 53 ll C,fac[10010],n,mut,a[1001000]; 54 int T,cnt,i,l,prime[N],p[N],psize,_cnt; 55 ll _e[100],_pr[100]; 56 vector<ll> d; 57 58 inline ll mul(ll a,ll b,ll p) { //return (a*b)%p 59 if (p<=1000000000) return a*b%p; 60 else if (p<=1000000000000ll) return (((a*(b>>20)%p)<<20)+(a*(b&((1<<20)-1))))%p; 61 else { 62 ll d=(ll)floor(a*(long double)b/p+0.5); 63 ll ret=(a*b-d*p)%p; 64 if (ret<0) ret+=p; 65 return ret; 66 } 67 } 68 69 void prime_table(){ //prime[1..tot]: prime[i]=ith prime 70 int i,j,tot,t1; 71 for (i=1;i<=psize;i++) p[i]=i; 72 for (i=2,tot=0;i<=psize;i++){ 73 if (p[i]==i) prime[++tot]=i; 74 for (j=1;j<=tot && (t1=prime[j]*i)<=psize;j++){ 75 p[t1]=prime[j]; 76 if (i%prime[j]==0) break; 77 } 78 } 79 } 80 81 void init(int ps) { //initial 82 psize=ps; 83 prime_table(); 84 } 85 86 ll powl(ll a,ll n,ll p) { //return (a^n)%p 87 ll ans=1; 88 for (;n;n>>=1) { 89 if (n&1) ans=mul(ans,a,p); 90 a=mul(a,a,p); 91 } 92 return ans; 93 } 94 95 bool witness(ll a,ll n) { 96 int t=0; 97 ll u=n-1; 98 for (;~u&1;u>>=1) t++; 99 ll x=powl(a,u,n),_x=0; 100 for (;t;t--) { 101 _x=mul(x,x,n); 102 if (_x==1 && x!=1 && x!=n-1) return 1; 103 x=_x; 104 } 105 return _x!=1; 106 } 107 108 bool miller(ll n) { 109 if (n<2) return 0; 110 if (n<=psize) return p[n]==n; 111 if (~n&1) return 0; 112 for (int j=0;j<=7;j++) if (witness(rand()%(n-1)+1,n)) return 0; 113 return 1; 114 } 115 116 ll gcd(ll a,ll b) { 117 ll ret=1; 118 while (a!=0) { 119 if ((~a&1) && (~b&1)) ret<<=1,a>>=1,b>>=1; 120 else if (~a&1) a>>=1; else if (~b&1) b>>=1; 121 else { 122 if (a<b) swap(a,b); 123 a-=b; 124 } 125 } 126 return ret*b; 127 } 128 129 ll rho(ll n) { 130 for (;;) { 131 ll X=rand()%n,Y,Z,T=1,*lY=a,*lX=lY; 132 int tmp=20; 133 C=rand()%10+3; 134 X=mul(X,X,n)+C;*(lY++)=X;lX++; 135 Y=mul(X,X,n)+C;*(lY++)=Y; 136 for(;X!=Y;) { 137 ll t=X-Y+n; 138 Z=mul(T,t,n); 139 if(Z==0) return gcd(T,n); 140 tmp--; 141 if (tmp==0) { 142 tmp=20; 143 Z=gcd(Z,n); 144 if (Z!=1 && Z!=n) return Z; 145 } 146 T=Z; 147 Y=*(lY++)=mul(Y,Y,n)+C; 148 Y=*(lY++)=mul(Y,Y,n)+C; 149 X=*(lX++); 150 } 151 } 152 } 153 154 void _factor(ll n) { 155 for (int i=0;i<cnt;i++) { 156 if (n%fac[i]==0) n/=fac[i],fac[cnt++]=fac[i];} 157 if (n<=psize) { 158 for (;n!=1;n/=p[n]) fac[cnt++]=p[n]; 159 return; 160 } 161 if (miller(n)) fac[cnt++]=n; 162 else { 163 ll x=rho(n); 164 _factor(x);_factor(n/x); 165 } 166 } 167 168 void dfs(ll x,int dep) { 169 if (dep==_cnt) d.push_back(x); 170 else { 171 dfs(x,dep+1); 172 for (int i=1;i<=_e[dep];i++) dfs(x*=_pr[dep],dep+1); 173 } 174 } 175 176 void norm() { 177 sort(fac,fac+cnt); 178 _cnt=0; 179 rep(i,0,cnt) if (i==0||fac[i]!=fac[i-1]) _pr[_cnt]=fac[i],_e[_cnt++]=1; 180 else _e[_cnt-1]++; 181 } 182 183 vector<ll> getd() { 184 d.clear(); 185 dfs(1,0); 186 return d; 187 } 188 189 vector<ll> factor(ll n) { //return all factors of n cnt:the number of factors 190 cnt=0; 191 _factor(n); 192 norm(); 193 return getd(); 194 } 195 196 vector<PLL> factorG(ll n) { 197 cnt=0; 198 _factor(n); 199 norm(); 200 vector<PLL> d; 201 rep(i,0,_cnt) d.push_back(make_pair(_pr[i],_e[i])); 202 return d; 203 } 204 205 bool is_primitive(ll a,ll p) { 206 //assert(miller(p)); 207 vector<PLL> D=factorG(p-1); 208 rep(i,0,SZ(D)) if (powmod(a,(p-1)/D[i].first,p)==1) return 0; 209 return 1; 210 } 211 } 212 213 ll x,y,k; 214 int _; 215 int main() { 216 Factor::init(200000); 217 for (scanf("%d",&_);_;_--) { 218 scanf("%I64d%I64d%I64d",&x,&y,&k); 219 vector<ll> c=Factor::factor(gcd(x,y)); //c:all factors of gcd(x,y) 220 sort(all(c)); // =all common factors of x and y 221 reverse(all(c)); 222 if (SZ(c)<k) puts("-1"); else printf("%I64d\n",c[k-1]); 223 } 224 }
posted on 2014-11-18 19:24 Pentium.Labs 阅读(688) 评论(0) 编辑 收藏 举报