算法总结—数论

来源：互联网发布：linux zip命令解压编辑：程序博客网时间：2024/05/17 23:25

4、数论

4.1 阶乘最后非0位

//求阶乘最后非零位,复杂度O(nlogn)

//返回该位,n以字符串方式传入

#include <string.h>

#define MAXN 10000

int lastdigit(char* buf){

const int mod[20]={1,1,2,6,4,2,2,4,2,8,4,4,8,4,6,8,8,6,8,2};

int len=strlen(buf),a[MAXN],i,c,ret=1;

if (len==1)

return mod[buf[0]-'0'];

for (i=0;i<len;i++)

a[i]=buf[len-1-i]-'0';

for (;len;len-=!a[len-1]){

ret=ret*mod[a[1]%2*10+a[0]]%5;

for (c=0,i=len-1;i>=0;i--)

c=c*10+a[i],a[i]=c/5,c%=5;

}

return ret+ret%2*5;

}

4.2 模线性方程组

#ifdef WIN32

typedef __int64 i64;

#else

typedef long long i64;

#endif

//扩展Euclid求解gcd(a,b)=ax+by

int ext_gcd(int a,int b,int& x,int& y){

int t,ret;

if (!b){

x=1,y=0;

return a;

}

ret=ext_gcd(b,a%b,x,y);

t=x,x=y,y=t-a/b*y;

return ret;

}

//计算m^a, O(loga), 本身没什么用, 注意这个按位处理的方法 :-P

int exponent(int m,int a){

int ret=1;

for (;a;a>>=1,m*=m)

if (a&1)

ret*=m;

return ret;

}

//计算幂取模a^b mod n, O(logb)

int modular_exponent(int a,int b,int n){ //a^b mod n

int ret=1;

for (;b;b>>=1,a=(int)((i64)a)*a%n)

if (b&1)

ret=(int)((i64)ret)*a%n;

return ret;

}

//求解模线性方程ax=b (mod n)

//返回解的个数,解保存在sol[]中

//要求n>0,解的范围0..n-1

int modular_linear(int a,int b,int n,int* sol){

int d,e,x,y,i;

d=ext_gcd(a,n,x,y);

if (b%d)

return 0;

e=(x*(b/d)%n+n)%n;

for (i=0;i<d;i++)

sol[i]=(e+i*(n/d))%n;

return d;

}

//求解模线性方程组(中国余数定理)

// x = b[0] (mod w[0])

// x = b[1] (mod w[1])

// ...

// x = b[k-1] (mod w[k-1])

//要求w[i]>0,w[i]与w[j]互质,解的范围1..n,n=w[0]*w[1]*...*w[k-1]

int modular_linear_system(int b[],int w[],int k){

int d,x,y,a=0,m,n=1,i;

for (i=0;i<k;i++)

n*=w[i];

for (i=0;i<k;i++){

m=n/w[i];

d=ext_gcd(w[i],m,x,y);

a=(a+y*m*b[i])%n;

}

return (a+n)%n;

}

4.3 素数

//用素数表判定素数,先调用initprime

int plist[10000],pcount=0;

int prime(int n){

int i;

if ((n!=2&&!(n%2))||(n!=3&&!(n%3))||(n!=5&&!(n%5))||(n!=7&&!(n%7)))

return 0;

for (i=0;plist[i]*plist[i]<=n;i++)

if (!(n%plist[i]))

return 0;

return n>1;

}

void initprime(){

int i;

for (plist[pcount++]=2,i=3;i<50000;i++)

if (prime(i))

plist[pcount++]=i;

}

//miller rabin

//判断自然数n是否为素数

//time越高失败概率越低,一般取10到50

#include <stdlib.h>

#ifdef WIN32

typedef __int64 i64;

#else

typedef long long i64;

#endif

int modular_exponent(int a,int b,int n){ //a^b mod n

int ret;

for (;b;b>>=1,a=(int)((i64)a)*a%n)

if (b&1)

ret=(int)((i64)ret)*a%n;

return ret;

}

// Carmicheal number: 561,41041,825265,321197185

int miller_rabin(int n,int time=10){

if (n==1||(n!=2&&!(n%2))||(n!=3&&!(n%3))||(n!=5&&!(n%5))||(n!=7&&!(n%7)))

return 0;

while (time--)

if (modular_exponent(((rand()&0x7fff<<16)+rand()&0x7fff+rand()&0x7fff)%(n-1)+1,n-1,n)!=1)

return 0;

return 1;

}

4.4 欧拉函数

int gcd(int a,int b){

return b?gcd(b,a%b):a;

}

inline int lcm(int a,int b){

return a/gcd(a,b)*b;

}

//求1..n-1中与n互质的数的个数

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;

}

0 0