• 椭圆曲线质因数分解


    这个程序分解99999640000243这个数.

    #include <map>
    #include <cmath>
    #include <cstdio>
    #include <vector>
    #include <cstring>
    #include <cstdlib>
    #include <algorithm>
    #define mod 1000000007
    #define maxint 2147483648ll
    #define nc1
    #define rg register
    #define rgv
    #define rg1 rgv
    #define rg2 rgv
    #define rg3 rgv
    #define rg4 rgv
    #define rg5 rg
    #define rg6 rg
    #define sorc
    #define psh(a,l,p) (a)[(l)++]=(p)
    typedef long long ll;
    namespace basicmath{
    	inline ll mp(ll a,ll b,ll p){
    	#ifdef __int128
    		return (__int128)a*(__int128)b%p;
    	#else
    		if(p<3037000000ll) return a*b%p;
    #ifdef nc1
    		if(p<1099511627776ll) return (((a*(b>>20)%p)<<20)+(a*(b&((1<<20)-1))))%p;
    #endif
    		ll d=(ll)(a*(long double)b/p+0.5); ll ret=(a*b-d*p)%p;
    		if (ret<0) ret+=p;
    		return ret;
    	#endif
    	}
    	inline ll fp(int a,int b,int p){ a%=p,b%=p-1; int ans=1; for(;b;b>>=1,a=(ll)a*a%p) if(b&1) ans=(ll)ans*a%p; return ans; }
    	inline ll fp(ll a,ll b,ll p){ a%=p,b%=p-1; ll ans=1ll; for(;b;b>>=1,a=mp(a,a,p)) if(b&1) ans=mp(a,ans,p); return ans; }
    	template<int max=10000007> struct sieve{
    		int q[max+5]; int pr[max/4],pl,rxsiz;
    		inline void generate(rg3 int n=max){
    			if(n<=rxsiz) return; if(n>max) return;
    			rxsiz=n;
    			for(rg1 int i=2;i<=n;++i){
    				if(!q[i]) pr[pl++]=i,q[i]=i;
    				for(rg2 int j=0;j<pl;++j){
    					rg4 int t=i*pr[j];
    					if(t>n) break;
    					q[t]=pr[j]; if(q[i]==pr[j]) break;
    				}
    			}
    		}
    		inline int operator()(int n){
    			#ifdef sorc
    			if(n>=pl) return 0;
    			if(n<0) return 0;
    			#endif
    			return pr[n];
    		}
    		#ifdef sorc
    		int sbcc;
    		inline int operator[](int n){
    			if(n>rxsiz) return sbcc;
    			if(n<0) return sbcc;
    			return q[n];
    		}
    		#else
    		inline int operator[](int n){
    			return q[n];
    		}
    		#endif
    	};
    	namespace primeTesting{
    		const int llrtm=7;
    		ll q[llrtm]={2, 325, 9375, 28178, 450775, 9780504, 1795265022ll};//for numbers less than 2^64 test these seven numbers are enough
    		inline bool witness(ll a,rg5 ll n){
    			a%=n;
    			if(!a) return 0;
    			rg5 int t=0; ll u=n-1; for(;~u&1;u>>=1) ++t;
    			rg5 ll x=fp(a,u,n),xx=0; while(t--){
    				xx=mp(x,x,n); if(xx==1 && x!=1 && x!=n-1) return 1;
    				x=xx;
    			}
    			return xx!=1;
    		}
    		inline bool miller(ll n){
    			if(!n) return 0;
    			if(n==2) return 1;
    			for(int i=0;i<llrtm;++i) if(witness(q[i],n)) return 0;
    			return 1;
    		}
    		template<class T> inline bool miller(ll n,T& sv){//T should have the same API as sieve
    			if(n<=sv.rxsiz && sv[n]==n) return 1;
    			return miller(n);
    		}
    		//the both miller tests are deterministic
    		struct pmt{
    			virtual inline bool operator()(ll n){
    				return 1;
    			}
    		};
    		struct mr_s100k:public pmt{
    			sieve<100000> sv;
    			inline mr_s100k(){
    				sv.generate(100000);
    			}
    			inline bool operator()(ll n){
    				return miller(n,sv);
    			}
    		};
    		struct mr_0:public pmt{
    			inline bool operator()(ll n){
    				return miller(n);
    			}
    		};
    	}
    	using namespace primeTesting;
    	namespace gcd{
    		inline ll euclid(ll a,ll b){return a?euclid(b%a,a):b;}
    		inline ll stein (rg6 ll a,rg6 ll b){
    //			if(a<0) a=-a;
    			rg6 ll ret=1;
    			while(a){
    				if((~a&1)&&(~b&1)) ret<<=1,a>>=1,b>>=1; else
    				if(~a&1) a>>=1; else if(~b&1) b>>=1; else
    				{ if(a<b) std::swap(a,b); a-=b; }
    			}
    			return b*ret;
    		}
    	}
    	namespace exgcd{
    		inline ll euclid(ll a,ll b,ll& x,ll& y){
    			if(b){
    				ll d=euclid(b,a%b,x,y);
    				ll t=x; x=y; y=t-a/b*y;
    				return d;
    			}else{
    				x=1,y=0;
    				return a;
    			}
    		}
    	}
    	namespace qmutils{
    		inline ll qm1(rg ll a,rg ll q){//if a < 0 a+=q
    			return a+((-((unsigned long long)a>>63))&q);
    		}
    		inline ll qm2(rg ll a,rg ll q){//if a>= q a-=q
    			return a-((-((unsigned long long)(q+(~a))>>63))&q);
    		}
    	}
    	using namespace qmutils;
    	namespace ec{
    		using namespace exgcd;
    		struct pi{
    			ll x,y,t;
    			inline pi(ll x=0,ll y=0,ll t=0):x(x),y(y),t(t){}
    			inline pi(const pi& p):x(p.x),y(p.y),t(p.t){}
    		};
    		inline bool operator==(const pi& a,const pi& b){
    			return a.x==b.x && a.y==b.y;
    		}
    		struct cv{//determines a curve
    			ll a,b;
    			ll x,y,p;
    			inline ll cb(ll x){ return mp(mp(x,x,p),x,p); }
    			inline pi getpi(){
    				return pi(x,y,0);
    			}
    			inline void random(ll px){
    				p=px;
    				x=rand()%(p-1)+1; y=rand()%(p-1)+1;
    				a=rand()%(p-1)+1;
    				b=qm1(qm1(mp(y,y,p)-cb(x),p)-mp(x,a,p),p);
    			}
    			inline ll slope(pi a,pi b,ll& f){
    				if(a==b){
    					ll yy=(3*mp(a.x,a.x,p))%p,xx=qm2(a.y+a.y,p);
    					ll xp,yp;
    					f=exgcd::euclid(p,xx,xp,yp);
    					yp=qm2(qm1(yp,p),p);
    					return mp(yy,yp,p);
    				}
    				ll yy=qm1(a.y-b.y,p),xx=qm1(a.x-b.x,p);
    				ll xp,yp;
    				f=exgcd::euclid(p,xx,xp,yp);
    				yp=qm2(qm1(yp,p),p);
    				return mp(yy,yp,p);
    			}
    			inline pi  add (pi a,pi b,ll& f){
    				if(a.t){
    					f=1;
    					return b;
    				}
    				if(b.t){
    					f=1;
    					return a;
    				}
    				if(a.x==b.x && a.y==-b.y){
    					f=1;
    					return pi(0,0,1);
    				}
    				ll sl=slope(a,b,f);
    				ll xp=qm1(qm1(mp(sl,sl,p)-a.x,p)-b.x,p);
    				return pi(
    					// x: sl^2 - xp - xq
    					xp ,
    					// y: yp + sl * (xR - xp)
    					p-qm2(a.y+mp(sl,qm1(xp-a.x,p),p),p) ,
    					0
    				);
    			}
    			inline pi mult(pi a,ll b,ll& f){
    				f=1;
    				pi t=pi(0,0,1);
    				for(;b;b>>=1,a=add(a,a,f)){
    					if(f!=1 && f!=p) return t;
    					if(b&1) t=add(t,a,f);
    					if(f!=1 && f!=p) return t;
    				}
    			}
    		};
    	}
    	#define dgcd gcd::euclid
    	namespace factor{
    		struct factr{
    			ll q[30],p[30],c[30];
    			int tl;
    			virtual inline int operator()(ll a,pmt& q){
    				return 0;
    			}
    		};
    		struct rho:public factr{
    			ll sq[70],sql;
    			private:
    				ll rhox(ll n,ll c,ll u){
    					ll i=1,k=2,y=u,x0=u;
    					while(1){
    						++i; x0=qm2(mp(x0,x0,n)+c,n);
    						ll d=dgcd(y-x0,n);
    						if(d!=1 && d!=n) return d;
    						if(y==x0) return n;
    						if(i==k) y=x0,k<<=1;
    					}
    				}
    				void fact_f(ll a,pmt& q){
    					if(q(a)) psh(sq,sql,a); else
    					{
    						ll x=a;
    						while(x==a){
    							ll c=rand()%a,d=(rand()&31)+3;
    							x=rhox(a,d,c);
    						}
    						fact_f(x,q),fact_f(a/x,q);
    					}
    				}
    				inline void norm(){
    					std::sort(sq,sq+sql);
    					tl=0;
    					for(rg int i=0,j=1;i<sql;i=(j++)){
    						q[tl]=sq[i],p[tl]=1,c[tl]=sq[i];
    						while(sq[j]==sq[i] && j<sql) c[tl]*=sq[j++],++p[tl];
    						++tl;
    					}
    				}
    			public:
    				inline int operator()(ll a,pmt& q){
    					sql=0;
    					fact_f(a,q);
    					norm();
    					return tl;
    				}
    		};
    		template<int nx=1000000> struct rhos:public factr{
    			ll sq[70],sql;
    			sieve<nx> sv;
    			private:
    				ll rhox(ll n,ll c,ll u){
    					ll i=1,k=2,y=u,x0=u;
    					while(1){
    						++i; x0=qm2(mp(x0,x0,n)+c,n);
    						ll d=dgcd(y-x0,n);
    						if(d!=1 && d!=n) return d;
    						if(y==x0) return n;
    						if(i==k) y=x0,k<<=1;
    					}
    				}
    				void fact_f(ll a,pmt& q){
    					if(a<=nx){//directly factor a
    						while(a!=1){
    							psh(sq,sql,sv[a]);
    							a/=sv[a];
    						}
    						return;
    					}
    					if(q(a)) psh(sq,sql,a); else
    					{
    						ll x=a;
    						while(x==a){
    							ll c=rand()%a,d=(rand()&31)+3;
    							x=rhox(a,d,c);
    						}
    						fact_f(x,q),fact_f(a/x,q);
    					}
    				}
    				inline void norm(){
    					std::sort(sq,sq+sql);
    					tl=0;
    					for(rg int i=0,j=1;i<sql;i=(j++)){
    						q[tl]=sq[i],p[tl]=1,c[tl]=sq[i];
    						while(sq[j]==sq[i] && j<sql) c[tl]*=sq[j++],++p[tl];
    						++tl;
    					}
    				}
    			public:
    				rhos(){
    					sv.generate(nx);
    				}
    				inline int operator()(ll a,pmt& q){
    					sql=0;
    					fact_f(a,q);
    //					for(int i=0;i<sql;++i) printf("%d ",sq[i]);
    					norm();
    					return tl;
    				}
    		};
    	}
    	using namespace factor;
    	struct factorsGetter{
    		//get all factors
    		std::vector<ll> fac;
    		void get(ll q,int le,factr& p){
    			if(le==p.tl) fac.push_back(q); else{
    				get(q,le+1,p);
    				for(int i=1;i<=p.p[le];++i) get((q*=p.q[le]),le+1,p);
    			}
    		}
    		int operator()(factr& p){
    			fac.clear();
    			get(1,0,p);
    			return fac.size();
    		}
    	};
    }
    using namespace basicmath::ec;
    int maint(){
    	cv nec;
    	nec.random(99999640000243ll);
    	pi t=nec.getpi();
    	ll tm;
    	pi q=t;
    	for(int i=2;i<=100;++i){
    		q=nec.mult(q,i,tm);
    		if(tm!=1 && tm!=99999640000243ll){
    			printf("%d
    ",i);
    			break;
    		}
    	}
    	if(tm!=1 && tm!=99999640000243ll) printf("f:%lld
    ",tm);
    	return 0;
    }
    int main(){
    	for(int i=0;i<10;++i) maint();
    	return 0;
    }
    

    跑得还好...

  • 相关阅读:
    TCP/IP笔记 一.综述
    Makefile的规则
    u盘安装ubuntu10.04 server.txt
    浅谈数据库技术,磁盘冗余阵列,IP分配,ECC内存,ADO,DAO,JDBC
    cocos2d-js 热更新具体解释(一)
    C#一个托付的样例
    JAVA学习之 异常处理机制
    阿里巴巴校招内推简历筛选方案
    《凑硬币》 动态规划算法入门
    android 读取xml
  • 原文地址:https://www.cnblogs.com/tmzbot/p/5008547.html
Copyright © 2020-2023  润新知