题意:
给k,求i是素数且在1~k内并且2^i-1是合数的情况,并将它分解。
分析:
枚举1至i然后用miller_rabin素数判定和pollard_rho因数分解即可。
代码:
//poj 2191 //sep9 #include <iostream> #include <map> #include <vector> #define gcc 10007 #define max_prime 200000 using namespace std; typedef long long ll; vector<int> primes; vector<bool> is_prime; ll gcd(ll a,ll b) { ll m=1; if(!b) return a; while(m){ m=a%b; a=b; b=m; } return a; } ll multi_mod(ll a,ll b,ll mod) { ll sum=0; while(b){ if(b&1) sum=(sum+a)%mod; a=(a+a)%mod; b>>=1; } return sum; } ll pow_mod(ll a,ll mod) { ll sum=1; while(b) { if(b&1) sum=multi_mod(sum,a,mod); a=multi_mod(a,mod); b>>=1; } return sum; } bool miller_rabin(ll n,int times) { if(n<2) return false; if(n==2) return true; if(!(n&1)) return false; ll q=n-1; int k=0; while(!(q&1)){ ++k; q>>=1; } for(int i=0;i<times;++i){ ll a=rand()%(n-1)+1; ll x=pow_mod(a,q,n); if(x==1) continue; for(int j=0;j<k;++j){ ll buf=multi_mod(x,x,n); if(buf==1&&x!=1&&x!=n-1) return false; x=buf; } if(x!=1) return false; } return true; } ll pollard_rho(ll n) { ll d,i,k,y; x=rand()%(n-1)+1; y=x; i=1; k=2; do { ++i; d=gcd(n+y-x,n); if(d>1&&d<n) return d; if(i==k) y=x,k*=2; x=((multi_mod(x,n)-gcc)%n+n)%n; }while(y!=x); return n; } bool isPrime(ll n) { if(n<=max_prime) return is_prime[n]; return miller_rabin(n,20); } void factorize(ll n,map<ll,int>& factors) { if(isPrime(n)) ++factors[n]; else{ for(int i=0;i<primes.size();++i){ int p=primes[i]; while(n%p==0){ ++factors[p]; n/=p; } } if(n!=1){ if(isPrime(n)) ++factors[n]; else{ ll d=pollard_rho(n); factorize(d,factors); factorize(n/d,factors); } } } } void ini_primes() { is_prime=vector<bool>(max_prime+1,true); is_prime[0]=is_prime[1]=false; for(int i=2;i<=max_prime;++i){ if(is_prime[i]){ primes.push_back(i); for(int j=i*2;j<=max_prime;j+=i) is_prime[j]=false; } } } int main() { ini_primes(); int i,k; scanf("%d",&k); ll p=1; for(i=1;i<=k;++i){ p=2*p; if(is_prime[i]) if(isPrime(p-1)==0){ map<ll,int> factor; factorize(p-1,factor); int first=0; for(map<ll,int>::iterator it=factor.begin();it!=factor.end();++it) while(it->second){ if(first==0){ printf("%I64d ",it->first); first=1; } else printf("* %I64d ",it->first); --it->second; } printf("= %I64d = ( 2 ^ %d ) - 1\n",p-1,i); } } return 0; }