acm数论!.docx

ACM 数论基本模板
1.欧几里得 求最大公约数,最小公倍数
(1)递归的写法:int gcd(int a,int b) {return b?gcd(b,a%b):a;}
(2)辗转相除法:
int gcd(int a,int b)
{
 if(a<b)    return gcd(b,a);
 int r;
 while(b) {r=a%b;a=b;b=r;}
 return a;
}
(3)stein+欧几里得 快速求解大数的最大公约数
i64 stein(i64 a,i64 b)
{ 
 if(a<b)   return stein(b,a);
 if(b==0) return a;
 if((a&1)==0&&(b&1)==0) return 2*stein(a>>1,b>>1);//a and b are even  
 if((a&1)==0) return  stein(a>>1,b);   // only a is  even    
 if((b&1)==0) return  stein(a,b>>1);   // only b is  even     
 return    stein((a+b)>>1,(a-b)>>1);   // a and b are odd
}
最小公倍数: int lcm(int a,int b) {return a/gcd(a,b)*b;}
2.扩展欧几里得 求ax=b (mod m) ax+my=b 如果r=gcd(a,m)且b%r==0,则同余方程有解,其最小解为x*(b/r);
ax+by=c 如r=gcd(a,b),则存在x,y,使xa+yb=r;当x+=b,y-=a后仍然成立        
因为xa+yb+ab-ab=r;==>(x+b)a+(y-a)b=r
int exgcd(int a,int b,int &x,int &y)
{
 if(b==0) {x=1;y=0;return a;}
 int r=exgcd(b,a%b,y,x);
 y-=x*(a/b);
 return r;
}
3.素数判定
(1)试除法:
bool  isprime(int n)
{
 int i;
 for(i=2;i<=(int)sqrt(n*1.0);i++)
     if(n%i==0)    return  false;
 return  true;
}
bool isprime(int n)
{
 if(n==2) return true;
 if(n==1||(n&1)==0)  return false;
 for(int i=3;i*i<=n;i+=2)   if(n%i==0) return fals;
 return true;
}
(2)miller-rabin 算法
bool witness(i64 a,i64 n)
{
  i64 x,d=1,i=ceil(log(n-1.0)/log(2.0))-1;
  for(;i>=0;i--)
  {
    x=d; d=(d*d)%n;
    if(d==1&&x!=1&&x!=n-1) return 1;
    if(((n-1)&(1<<i))>0) d=(d*a)%n;
  }
  return  d==1?0:1;
}
bool miller_rabin(i64 n)
{
  if(n==2)   return 1;
  if(n==1||(n&1)==0) return 0;
  i64 j,a;
  for(j=0;j<50;j++)
  {
   a=rand()*(n-2)/RAND_MAX+1;
   if(witness(a,n)) return 0;
  }
  return 1;
}
另一种写法，更好理解

bool witness(i64 a,i64 n)
{
 int i,j=0;
 i64 m=n-1,x,y;
 while(m%2==0)
 {
  m>>=1;
  j++;
 }
 x=pow(a,m,n);///快速幂取模
 for(i=1;i<=j;i++)
 {
  y=pow(x,2,n);
  if(y==1&&x!=1&&x!=n-1) return true;
  x=y;
 }
 return y==1?false:true;
}
bool miller_rabin(i64 n)
{
 if(n==2) return true;
 if(n==1||n%2==0) return false;
 for(int i=1;i<=10;i++)
 {
  i64 a=rand()%(n-1)+1;
  if(witness(a,n)) return false;
 }
 return true;
}

4.素数筛法   //前17个素数 prime[18]={17,2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59}
bool f[100002];//保存判断是否是素数的结果，p[i]=1 是素数，p[i]=0 则不是素数
int prime[78499];//保存素数prime[0]为素数的个数
void PRIME(int M)
{
 int i,i2,k;
 for(i=0;i<=M;i+=2) f[i]=0;
 for(i=1;i<=M;i+=2) f[i]=1;
 f[1]=0;  f[2]=1;
 for(i=3;i<=(int)sqrt(1.0*M);i+=2)
 if(p[i])
 {
  i2=i+i;  k=i*i;
  while(k<=M) {f[k]=0; k+=i2;}
 }
 prime[1]=2; k=1;
 for(i=3;i<=M;i+2) if(f[i]) prime[++k]=i;
 prime[0]=k;
}
(2)
void PRIME(int M)
{
 int i,j,k;
 prime[1]=2; prime[2]=3;
 for(i=5;i<=M;i+=2)
 {
  for(j=1;prime[j]*prime[j]<=i;j++)
  if(i%prime[j]==0) goto loop;
  prime[++k]=i;
  loop:;
 }
 prime[0]=k;
}
5.整数分解
 (1)
void split(int n,int *p,int *t)
{
 int i,s,top=0;
 for(i=1;i<=prime[0];i++)
 {
  s=0;
  while(n%prime[i]==0)  {s++;n/=prime[i];}
  if(s)  {p[++top]=prime[i];t[top]=s;}
  if(n==1) break;
  if(n<prime[i]*prime[i]) {p[++top]=n;t[top]=1;n=1;break;}
 }
 p[0]=t[0]=top;
}
(2)分解1-100000的因子,且由prime[n][]保存n的素因子(prime[n][0]为质因子的个数):
void split(int n)//p[]为素数表
{
 int i,x=n;
 prime[n][0]=0;
 for(i=1;i<=p[0];i++) ////if(x%p[i]==0)严重坑爹的bug
 {
  prime[n][++prime[n][0]]=p[i];
  while(x%p[i]==0) x/=p[i];
  if(x==1) break;
 }
 if(x>1) prime[n][++prime[n][0]]=x;
}
(3)Pollard-pho大数分解

i64 Pollard(i64 n,int c)
{
 i64 i=1,k=2,x=rand()%n,y=x,d;
 srand(time(NULL));
 while(true)
 {
  i++;
  x=(mod_mult(x,x,n)+c)%n;
  d=gcd(y-x,n);
  if(d>1&&d<n) return d;
  if(y==x) return n;
  if(i==k)  { y=x; k<<=1; }
 }
}
6.求因子和与因子个数（包含1和本身）

因子和s是积性函数,即:gcd(a,b)=1==>s(a*b)=s(a)*s(b);
如果p是素数==>s(p^X)=1+p+p^2+...+p^X=(p^(X+1)-1)/(p-1);s(p^2x)=1+p+p^2+...+p^2x=(p^(2x+1)-1)/(p-1);
求因子和:
(1)ans=1+n;
for(i=2;i<=n/2;i++)
if(n%i==0)
{
 if(n/i>i) ans+=i+n/i;
 else  if(n/i==i) ans+=i;
 else  break;
}
(2)另一种递推的写法:
for(i=1;i<=lmax;i++)   num[i]=1+i;
for(i=2;i<=lmax/2;i++)
 for(j=i<<1;j<=lmax;j+=i)  num[j]+=i;
求因子个数:
n=p1^t1*p2^t2*p3^t3***pk^tk; 因子个数为:(t1+1)(t2+1)***(tk+1)
for(ret=i=1;i<=prime[0]&&prime[i]<=(int)sqrt(1.0*n);i++)
if(n%prime[i]==0)
{
 k=0;
 while(n%prime[i]==0) {k++,n/=prime[i];}
 ret*=k+1;
}
if(n>1) ret*=2;
当求n^2的因子个数的时候:n^2=p1^(2*t1)*p2^(2*t2)***pk^(2*tk);  因子个数为:（2*t1+1)(2*t2+1)***(2*tk+1)
for(ret=i=1;i<=prime[0]&&prime[i]<=(int)sqrt(1.0*n);i++)
if(n%prime[i]==0)
{
 k=0;
 while(n%prime[i]==0) {k++,n/=prime[i];}
 ret*=2*k+1;
}
if(n>1) ret*=3;
快速求出一个比较大的区间内的所有因子和:
const int lmax=50005;//求出[1,50005]区间内每一个数的因子和(不包括本身),并用facsum[]数组保存
i64  facsum[lmax];
for(i=0;i<=lmax;i++)   facsum[i]=1;
 for(i=2;i*i<=lmax;i++)
 {
  for(j=i+1;j*i<=lmax;j++)  facsum[i*j]+=i+j;
  facsum[i*i]+=i;
 }
7.欧拉函数
(1)单独求欧拉函数
int eular(int n)
{
 int ret=1,i;
 for(i=2;i*i<=n;i++)
 if(n%i==0)
 {
  n/=i;  ret*=i-1;
  while(n%i==0) {n/=i;ret*=i;}
 }
 if(n>1) ret*=n-1;
 return  ret;
}
int euler(int x)
{
 int i, res=x,tmp=(int)sqrt(x*1.0)+1;
 for(i=2;i<tmp;i++)
 if(x%i==0)
 {
  res=res/i*(i-1);
  while(x%i== 0) x/=i;
 }
 if(x>1) res=res/x*(x-1);
 return res;
}
int eular(int n)
{
 int ret=n,i;
 for(i=2;i*i<=n;i++)
 if(n%i==0)
 {
  ret=ret/i*(i-1);
  while(n%i==0) n/=i;
 }
 if(n>1) ret=ret/n*(n-1);
 return ret;
}
先素数筛法在用欧拉函数(在此仅写其中的一个)
void eular(int n)
{
 int ret=n,i;
 for(i=1;i<=prime[0]&&prime[i]<=(int)sqrt(1.0*n);i++)
 if(n%prime[i]==0)
 {
  ret=ret/prime[i]*(prime[i]-1);
  while(n%prime[i]==0) n/=prime[i];
 }
 if(n>1) ret=ret/n*(n-1);
 return ret;
}
(2)递推求欧拉函数
const int  lmax=300000;
int PHI(int lmax)
{
 int i,j;
 for(i=1;i<=lmax;i++) phi[i]=i&1?i:i/2;
 for(i=3;i<=lmax;i+=2)
 if(phi[i])
 for(j=i;j<=lmax;j+=i) phi[j]=phi[j]/i*(i-1);
}
(3)同时求出欧拉值和素数
int prime[lmax][25],num[lmax],eular[lmax];//prime[n][i]表示n的第i+1个素数因子,num[n]表示n的因子个数,eular[n]表示n的欧拉值
///////下面这个貌似不对
void eular_prime()
{
 int i,j;
 for(i=1;i<=lmax;i++) eular[i]=i,prime[i][0]=0;
 for(i=2;i<=lmax;i++)
 if(eular[i]==i)
 {
  eular[i]=i-1;
  for(j=i<<1;j<=lmax;j+=i)
  {
   eular[j]=eular[j]/i*(i-1);
   prime[j][++prime[j][0]]=i;//prime[j][0]=num[j]表示j的素数因子个数
  }
 }
}
/////上面这个貌似不对
void eular_prime()//每个数的欧拉函数值及筛选法得到数的素因子num[i]为i的因子个数
{
 eular[1]=1;
 for(int i=2;i<lmax;i++)
 {
  if(eular[i]==0)
  for(int j=i;j<lmax;j+=i)
  {
   if(eular[j]==0)  eular[j]=j;
   eular[j]=eular[j]*(i-1)/i;
   prime[j][num[j]++]=i;
  }
  //eular[i]+=eular[i-1];//进行累加（法里数列长度）
 }
}
void eular_prime()
{
 int i,j;
 eular[1]=1;
 for(i=2;i<=lmax;i++)
 if(eular[i]==0)
 for(j=i;j<=lmax;j+=i)
 {
  if(eular[j]==0) eular[j]=j;
  eular[j]=eular[j]/i*(i-1);
  prime[j][++prime[j][0]]=i;
 }
}
欧拉定理的一个重要应用:A^x mod m=A^(x%phi(m)+phi(m)) mod m  (当x>=phi(m)时)
8.求逆元ax=1 (mod m) x是a的逆元

(1)用扩展欧几里得求
int Inv(int a,int m)
{
  int r,x,y;
  r=exgcd(a,m,x,y);
  if(r==1) return (x%m+m)%m;
  return -1;
}
(2)用快速幂取模求a*a^(p-2)=a^(p-1)=1 (mod p) p必须为素数,a的逆元是a^(p-2);
int pow(int a,int n)//a^n%p (n=p-2)
{//这里的做法会让a的值变化,可令t=a;用t代替a计算
 int r=1;
 while(n)
 {
  if(n&1) r=r*a%p;
  a=a*a%p;
  n>>=1;
 }
 return r;
}
9.快速模乘 a*b%p
int mul(int a,int b)
{
 int r=0;
 while(b)
 {
  if(b&1) r=(r+a)%p;
  a=(a<<1)%p;
  b>>=1;
 }
 return r;
}
10.求解模线性方程组(中国剩余定理)
  x=a1 mod m1
  x=a2 mod m2
  ......
  x=an mod mn  其中,a[],m[]已知,m[i]>0且m[i]与m[j]互质,求x.
设m1,m2,...,mn是两两互素的正数,则对任意的整数a1,a2,...,an,同余方程组
其解为:X=((M_1*M1*a1)+(M_2*M2*a2)+...+(M_n*Mn*an)) mod m;
其中m=m1*m2*...*mn;  Mi=m/mi; M_i是Mi的逆元
int china(int *a,int *m,int n)
{
   int M=1,ans=0,mi,i,x,y;
   for(i=0;i<n;i++) M*=m[i];
   for(i=0;i<n;i++)
   {
       mi=M/m[i];
       exgcd(m[i],mi,x,y);
       ans=(ans+mi*y*a[i])%M;
   }
   return (ans%M+M)%M;
}
ACM 数论基本模板
(2012-08-08 21:50:02)
转载▼
标签：
数论
模板
杂谈
分类： ACM
转自 ：
ACM 数论基本模板
1.欧几里得 求最大公约数,最小公倍数
(1)递归的写法:int gcd(int a,int b) {return b?gcd(b,a%b):a;}
(2)辗转相除法:
int gcd(int a,int b)
{
 if(a<b)    return gcd(b,a);
 int r;
 while(b) {r=a%b;a=b;b=r;}
 return a;
}
(3)stein+欧几里得 快速求解大数的最大公约数
i64 stein(i64 a,i64 b)
{ 
 if(a<b)   return stein(b,a);
 if(b==0) return a;
 if((a&1)==0&&(b&1)==0) return 2*stein(a>>1,b>>1);//a and b are even  
 if((a&1)==0) return  stein(a>>1,b);   // only a is  even    
 if((b&1)==0) return  stein(a,b>>1);   // only b is  even     
 return    stein((a+b)>>1,(a-b)>>1);   // a and b are odd
}
最小公倍数: int lcm(int a,int b) {return a/gcd(a,b)*b;}
2.扩展欧几里得 求ax=b (mod m) ax+my=b 如果r=gcd(a,m)且b%r==0,则同余方程有解,其最小解为x*(b/r);
ax+by=c 如r=gcd(a,b),则存在x,y,使xa+yb=r;当x+=b,y-=a后仍然成立        
因为xa+yb+ab-ab=r;==>(x+b)a+(y-a)b=r
int exgcd(int a,int b,int &x,int &y)
{
 if(b==0) {x=1;y=0;return a;}
 int r=exgcd(b,a%b,y,x);
 y-=x*(a/b);
 return r;
}
3.素数判定
(1)试除法:
bool  isprime(int n)
{
 int i;
 for(i=2;i<=(int)sqrt(n*1.0);i++)
     if(n%i==0)    return  false;
 return  true;
}
bool isprime(int n)
{
 if(n==2) return true;
 if(n==1||(n&1)==0)  return false;
 for(int i=3;i*i<=n;i+=2)   if(n%i==0) return fals;
 return true;
}
(2)miller-rabin 算法
bool witness(i64 a,i64 n)
{
  i64 x,d=1,i=ceil(log(n-1.0)/log(2.0))-1;
  for(;i>=0;i--)
  {
    x=d; d=(d*d)%n;
    if(d==1&&x!=1&&x!=n-1) return 1;
    if(((n-1)&(1<<i))>0) d=(d*a)%n;
  }
  return  d==1?0:1;
}
bool miller_rabin(i64 n)
{
  if(n==2)   return 1;
  if(n==1||(n&1)==0) return 0;
  i64 j,a;
  for(j=0;j<50;j++)
  {
   a=rand()*(n-2)/RAND_MAX+1;
   if(witness(a,n)) return 0;
  }
  return 1;
}
另一种写法，更好理解

bool witness(i64 a,i64 n)
{
 int i,j=0;
 i64 m=n-1,x,y;
 while(m%2==0)
 {
  m>>=1;
  j++;
 }
 x=pow(a,m,n);///快速幂取模
 for(i=1;i<=j;i++)
 {
  y=pow(x,2,n);
  if(y==1&&x!=1&&x!=n-1) return true;
  x=y;
 }
 return y==1?false:true;
}
bool miller_rabin(i64 n)
{
 if(n==2) return true;
 if(n==1||n%2==0) return false;
 for(int i=1;i<=10;i++)
 {
  i64 a=rand()%(n-1)+1;
  if(witness(a,n)) return false;
 }
 return true;
}

4.素数筛法   //前17个素数 prime[18]={17,2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59}
bool f[100002];//保存判断是否是素数的结果，p[i]=1 是素数，p[i]=0 则不是素数
int prime[78499];//保存素数prime[0]为素数的个数
void PRIME(int M)
{
 int i,i2,k;
 for(i=0;i<=M;i+=2) f[i]=0;
 for(i=1;i<=M;i+=2) f[i]=1;
 f[1]=0;  f[2]=1;
 for(i=3;i<=(int)sqrt(1.0*M);i+=2)
 if(p[i])
 {
  i2=i+i;  k=i*i;
  while(k<=M) {f[k]=0; k+=i2;}
 }
 prime[1]=2; k=1;
 for(i=3;i<=M;i+2) if(f[i]) prime[++k]=i;
 prime[0]=k;
}
(2)
void PRIME(int M)
{
 int i,j,k;
 prime[1]=2; prime[2]=3;
 for(i=5;i<=M;i+=2)
 {
  for(j=1;prime[j]*prime[j]<=i;j++)
  if(i%prime[j]==0) goto loop;
  prime[++k]=i;
  loop:;
 }
 prime[0]=k;
}
5.整数分解
 (1)
void split(int n,int *p,int *t)
{
 int i,s,top=0;
 for(i=1;i<=prime[0];i++)
 {
  s=0;
  while(n%prime[i]==0)  {s++;n/=prime[i];}
  if(s)  {p[++top]=prime[i];t[top]=s;}
  if(n==1) break;
  if(n<prime[i]*prime[i]) {p[++top]=n;t[top]=1;n=1;break;}
 }
 p[0]=t[0]=top;
}
(2)分解1-100000的因子,且由prime[n][]保存n的素因子(prime[n][0]为质因子的个数):
void split(int n)//p[]为素数表
{
 int i,x=n;
 prime[n][0]=0;
 for(i=1;i<=p[0];i++) ////if(x%p[i]==0)严重坑爹的bug
 {
  prime[n][++prime[n][0]]=p[i];
  while(x%p[i]==0) x/=p[i];
  if(x==1) break;
 }
 if(x>1) prime[n][++prime[n][0]]=x;
}
(3)Pollard-pho大数分解

i64 Pollard(i64 n,int c)
{
 i64 i=1,k=2,x=rand()%n,y=x,d;
 srand(time(NULL));
 while(true)
 {
  i++;
  x=(mod_mult(x,x,n)+c)%n;
  d=gcd(y-x,n);
  if(d>1&&d<n) return d;
  if(y==x) return n;
  if(i==k)  { y=x; k<<=1; }
 }
}
6.求因子和与因子个数（包含1和本身）

因子和s是积性函数,即:gcd(a,b)=1==>s(a*b)=s(a)*s(b);
如果p是素数==>s(p^X)=1+p+p^2+...+p^X=(p^(X+1)-1)/(p-1);s(p^2x)=1+p+p^2+...+p^2x=(p^(2x+1)-1)/(p-1);
求因子和:
(1)ans=1+n;
for(i=2;i<=n/2;i++)
if(n%i==0)
{
 if(n/i>i) ans+=i+n/i;
 else  if(n/i==i) ans+=i;
 else  break;
}
(2)另一种递推的写法:
for(i=1;i<=lmax;i++)   num[i]=1+i;
for(i=2;i<=lmax/2;i++)
 for(j=i<<1;j<=lmax;j+=i)  num[j]+=i;
求因子个数:
n=p1^t1*p2^t2*p3^t3***pk^tk; 因子个数为:(t1+1)(t2+1)***(tk+1)
for(ret=i=1;i<=prime[0]&&prime[i]<=(int)sqrt(1.0*n);i++)
if(n%prime[i]==0)
{
 k=0;
 while(n%prime[i]==0) {k++,n/=prime[i];}
 ret*=k+1;
}
if(n>1) ret*=2;
当求n^2的因子个数的时候:n^2=p1^(2*t1)*p2^(2*t2)***pk^(2*tk);  因子个数为:（2*t1+1)(2*t2+1)***(2*tk+1)
for(ret=i=1;i<=prime[0]&&prime[i]<=(int)sqrt(1.0*n);i++)
if(n%prime[i]==0)
{
 k=0;
 while(n%prime[i]==0) {k++,n/=prime[i];}
 ret*=2*k+1;
}
if(n>1) ret*=3;
快速求出一个比较大的区间内的所有因子和:
const int lmax=50005;//求出[1,50005]区间内每一个数的因子和(不包括本身),并用facsum[]数组保存
i64  facsum[lmax];
for(i=0;i<=lmax;i++)   facsum[i]=1;
 for(i=2;i*i<=lmax;i++)
 {
  for(j=i+1;j*i<=lmax;j++)  facsum[i*j]+=i+j;
  facsum[i*i]+=i;
 }
7.欧拉函数
(1)单独求欧拉函数
int eular(int n)
{
 int ret=1,i;
 for(i=2;i*i<=n;i++)
 if(n%i==0)
 {
  n/=i;  ret*=i-1;
  while(n%i==0) {n/=i;ret*=i;}
 }
 if(n>1) ret*=n-1;
 return  ret;
}
int euler(int x)
{
 int i, res=x,tmp=(int)sqrt(x*1.0)+1;
 for(i=2;i<tmp;i++)
 if(x%i==0)
 {
  res=res/i*(i-1);
  while(x%i== 0) x/=i;
 }
 if(x>1) res=res/x*(x-1);
 return res;
}
int eular(int n)
{
 int ret=n,i;
 for(i=2;i*i<=n;i++)
 if(n%i==0)
 {
  ret=ret/i*(i-1);
  while(n%i==0) n/=i;
 }
 if(n>1) ret=ret/n*(n-1);
 return ret;
}
先素数筛法在用欧拉函数(在此仅写其中的一个)
void eular(int n)
{
 int ret=n,i;
 for(i=1;i<=prime[0]&&prime[i]<=(int)sqrt(1.0*n);i++)
 if(n%prime[i]==0)
 {
  ret=ret/prime[i]*(prime[i]-1);
  while(n%prime[i]==0) n/=prime[i];
 }
 if(n>1) ret=ret/n*(n-1);
 return ret;
}
(2)递推求欧拉函数
const int  lmax=300000;
int PHI(int lmax)
{
 int i,j;
 for(i=1;i<=lmax;i++) phi[i]=i&1?i:i/2;
 for(i=3;i<=lmax;i+=2)
 if(phi[i])
 for(j=i;j<=lmax;j+=i) phi[j]=phi[j]/i*(i-1);
}
(3)同时求出欧拉值和素数
int prime[lmax][25],num[lmax],eular[lmax];//prime[n][i]表示n的第i+1个素数因子,num[n]表示n的因子个数,eular[n]表示n的欧拉值
///////下面这个貌似不对
void eular_prime()
{
 int i,j;
 for(i=1;i<=lmax;i++) eular[i]=i,prime[i][0]=0;
 for(i=2;i<=lmax;i++)
 if(eular[i]==i)
 {
  eular[i]=i-1;
  for(j=i<<1;j<=lmax;j+=i)
  {
   eular[j]=eular[j]/i*(i-1);
   prime[j][++prime[j][0]]=i;//prime[j][0]=num[j]表示j的素数因子个数
  }
 }
}
/////上面这个貌似不对
void eular_prime()//每个数的欧拉函数值及筛选法得到数的素因子num[i]为i的因子个数
{
 eular[1]=1;
 for(int i=2;i<lmax;i++)
 {
  if(eular[i]==0)
  for(int j=i;j<lmax;j+=i)
  {
   if(eular[j]==0)  eular[j]=j;
   eular[j]=eular[j]*(i-1)/i;
   prime[j][num[j]++]=i;
  }
  //eular[i]+=eular[i-1];//进行累加（法里数列长度）
 }
}
void eular_prime()
{
 int i,j;
 eular[1]=1;
 for(i=2;i<=lmax;i++)
 if(eular[i]==0)
 for(j=i;j<=lmax;j+=i)
 {
  if(eular[j]==0) eular[j]=j;
  eular[j]=eular[j]/i*(i-1);
  prime[j][++prime[j][0]]=i;
 }
}
欧拉定理的一个重要应用:A^x mod m=A^(x%phi(m)+phi(m)) mod m  (当x>=phi(m)时)
8.求逆元ax=1 (mod m) x是a的逆元

(1)用扩展欧几里得求
int Inv(int a,int m)
{
  int r,x,y;
  r=exgcd(a,m,x,y);
  if(r==1) return (x%m+m)%m;
  return -1;
}
(2)用快速幂取模求a*a^(p-2)=a^(p-1)=1 (mod p) p必须为素数,a的逆元是a^(p-2);
int pow(int a,int n)//a^n%p (n=p-2)
{//这里的做法会让a的值变化,可令t=a;用t代替a计算
 int r=