椭圆曲线质因数分解
这个程序分解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\n",i);
break;
}
}
if(tm!=1 && tm!=99999640000243ll) printf("f:%lld\n",tm);
return 0;
}
int main(){
for(int i=0;i<10;++i) maint();
return 0;
}
跑得还好...
