关于数的一些算法

1.寻找一个数的所有质因数（整数分解）

不要一个个找因数，然后再判断是不是素数，这样太慢了，这是比较高效的算法思路是把某一num的素因数一个个剥离出来：（且从小到大存在vector里）

算法复杂度：最差情况：O（根号n）比如素数100000007， 最优情况：O（logN）比如8，16，32，64，3，9，27，81...（因为num/i这步是根号n的复杂度向logN的转变）

vector<int> Prime_Factors;for(i=2;i*i<=num;i++)              //i*i<num 和判断素数的循环条件 一样作用也相同，当num是素数或者1的时候结束循环。 {    if(num%i==0)//只要num%i==0那么i一定是素数。若i不是素数，他的素因数比自身小，已经从前几次循环中从num中剥离掉了。     {        Prime_Factors.push_back(i);        while(num%i==0)            num/=i;           //核心步骤     ｝}if(num>1) Prime_Factors.push_back(num);//最后剩下的num如果不是1的话那么也是它的素因数。 

也可以打一个根号n以内的素数表复杂度O（根号N）这样就不需要一个一个的i++了，直接用素数表的数当作i即可,一般题目是多组数据所以用素数表来做复杂度会更低。

精练一点的代码有：

map<int ,int> prime_factor(int n){map<int,int> res;for(int i=2;i*i<=n;i++){while(n%i==0){++res[i];n/=i;}}if(n!=1) res[n]=1;return res;}

比如hduoj 的Alexandra and Prime Numbers题

2.寻找一个数n的所有约数

按照约数的性质，一个数n有一半的约数存在于根号n以内，另一半多的约数存在于n-根号n范围内，而且由一个已知的约数m，可以推出另一个约数n/m。所以可以遍历根号n以内的所有约数并且记录下来，再根据现有的推出其他的约数（代码是从小到大储存）复杂度O(根号n)

    cin>>n;    int Divisor[10000];    for(int i=1;i*i<=n;i++)        if(n%i==0) Divisor[p++]= i;    for(int j=p-1;j>=0;j--)        if( Divisor[j]*Divisor[j]==n ) continue;        else Divisor[p++]= n/Divisor[j];

精炼点的代码：

vector<int> divisor(n){vector<int >res;for(int i=1;i*i<=n;i++){if(n%i==0) {res.push_back(i);if(i!=n/i) res.push_back(n/i);}}return res;}

3.寻找a，b两个整数的所有公约数

由数论知识可知：a和b的所有公约数都是他们的最大公约数的约数，所以可以先求出gcd 复杂度O（logN）再根据标题2求出gcd的所有约数复杂度O（根号（gcd(a,b)））

    cin>>a>>b;    int n,p=1;    n=gcd(a,b);    int Divisor[10000];    for(int i=1;i*i<=n;i++)        if(n%i==0) Divisor[p++]= i;    for(int j=p-1;j>=0;j--)        if( Divisor[j]*Divisor[j]==n ) continue;        else Divisor[p++]= n/Divisor[j];

比如 http://blog.csdn.net/kalilili/article/details/42640687（求a和b的第k大公约数）

4.素数表模版

素数代码：（简单易懂）

ll prime[N],primenum;//有primenum个素数 include<math.h>void PRIME(ll Max_Prime){primenum=0;prime[primenum++]=2;for(ll i=3;i<=Max_Prime;i+=2)for(ll j=0;j<primenum;j++)if(i%prime[j]==0)break;else if(prime[j]>sqrt((double)i) || j==primenum-1)//||前面那个条件用于剪枝,||后面是必须的{prime[primenum++]=i;break;}}

线性筛代码：（较难，效率高）

const int maxn = 100000003;int p[6666666], tot;bool vis[maxn];void get_prime(){    int i, j;    for (i = 2;i < maxn; i++) {        if (!vis[i]) p[++tot] = i;        for (j = 1; j <= tot; j++){            if (p[j] * i >= maxn)break;            vis[p[j]*i] = 1;      //根据线性筛的理论知识知,应该有p[j]*(1...i)均标记,但这里只标记了p[j]*i来优化,原因自己思考             if (i % p[j] == 0)break; //如果i不是素数,标记完p[j]*i后就不需再标记p[j...]*i了,又是一个优化.         }    }}

上面的代码虽然优化得较好但是不易码出来，可以用易理解的模版：

int prime[MAX_N];bool is_prime[MAX_N+1];int seive(int n){    int p=0;    for(int i=0;i<=n;i++) is_prime[i]=true;    is_prime[0]=is_prime[1]=false;    for(int i=2;i<=n;i++)        if(is_prime[i])        {            prime[p++]=i;            for(int j=2*i(或i*i 注意防止爆int);j<=n;j+=i) is_prime[j] = false;        }    return p;}

复杂度：n+n/2+n/3+n/4+...+n/n=nlnn 所以复杂度上界是O(nlnn) 上面那个线性筛的模版经过优化更多，复杂度会更优

打表代码：

int main(){      PRIME(10000);      FILE *fp=fopen("素数表","w");      for(ll i=0;i<primenum;i++)      fprintf(fp,"%lld,",prime[i]);      fclose(fp);      return 0;  } 

5.区间筛法

给定整数a,b 求区间[a,b)的素数，当a，b很大的时候不能直接用筛法，所以要空间和时间优化，先筛去[2，根号b)的合数，再向上取整筛去[a,b)的合数。

typedef long long ll;bool is_prime[MAX_L];bool is_prime_small[MAX_SQRT_B];// 对区间[a,b)内的整数执行筛法,is_prime[i-a]=true <==> i是素数void segment_sieve(ll a,ll b){    for(int i=0;(ll)i*i<b;i++) is_prime_small[i]=true;    for(int i=0;i<b-a;i++) is_prime[i]=true;    if(1-a>=0) is_prime[1-a]=false;   //易错因为1不是素数也不是合数，这也是区间筛的一个易错bug         for(int i=2;(ll)i*i<b;i++){        if(is_prime_small[i]){            for(int j=2*i;(ll)j*j<b;j+=i) is_prime_small[j] = false; //筛[2,根号b)            for(int j=max((ll)2,(a+i-1)/i)*i;j<b;j+=i) is_prime[j-a]=false;            // ((a-1)/i+1)*i，向上取整找到a以上的第一个i的倍数，这里易错，有可能a<=i,            //这个时候a以上的第一个i的倍数就是i。这个时候不要把它筛去了，所以用到了max        }    }}

6.扩展欧几里德算法

a 和 b 的最大公约数是 gcd ，那么，我们一定能够找到这样的 x 和 y ，使得: a*x + b*y = gcd 这是一个不定方程，有多解是一定的，但是只要我们找到一组特殊的解 x0 和 y0 那么，我们就可以用 x0 和 y0 表示出整个不定方程的通解：

 x = x0 + (b/gcd)*t       

y = y0 –(a/gcd)*t

假设当前我们要处理的是求出 a 和 b的最大公约数，并求出 x 和 y 使得 a*x + b*y= gcd ，而我们已经求出了下一个状态：b 和 a%b 的最大公约数，并且求出了一组x1 和y1 使得： b*x1 + (a%b)*y1 = gcd ，

由a%b = a - (a/b)*b（这里的 “/” 指的是整除）

得：gcd = b*x1 + (a-(a/b)*b)*y1 = b*x1 + a*y1 – (a/b)*b*y1  = a*y1 + b*(x1 – a/b*y1)

解：x = y1，y = x1 – a/b*y1

int e_gcd(int a,int b,int&x,int&y){if(!a&&!b) return -1;//无最大公约数int d=a;if(b!=0){ans=e_gcd(b,a%b,y,x); //把x,y颠倒，由于引用传参，所以在下一步就省了swap(x,y) y-=(a/b)*x;    //x，y相当已经颠倒了，所以别把x，y写反了 }else{x=1;y=0;}return d;}

一般写成这样即可：

ll e_gcd(ll a,ll b,ll&x,ll&y)  {      ll ans;      if(b==0) ans=a,x=1,y=0;      else ans=e_gcd(b,a%b,y,x),y-=(a/b)*x;      return ans;  }  

相关题目 POJ 1061 http://poj.org/problem?id=1061 HDU 2669 http://acm.hdu.edu.cn/showproblem.php?pid=2669

7.求乘法逆元

 对a∈Z*n，存在x∈Z*n，使得a×x ≡1 (mod n)，则称x为a的乘法逆元。可以看出a，n必须互质才有解。

可以等价于求这样的方程的解x： a*x + n*y = 1，所以乘法逆元有多组解，一般都是求最小的非负整数解，所以可以用扩展欧几里德解出一个特解x0，所以通解是 x = x0 + (n/gcd(a,n))*t   gcd(a,n)=1，所以通解为x = x0 + n*t

所以求最小的非负整数解x即把特解x0模n即可，这里指数学上的模，在计算机上当x0是负数时%n还是负数，所以还要处理一步，当x0负数那么解为 ：x0模上 abs（m）再加上abs(m) ，于是，我们不难写出下面的代码求解:

//ax= 1(mod n)  --->  ax+ny=1long long mod_reverse(long long a,long long n)  {      long long x,y;      long long d=extend_gcd(a,n,x,y);      if(d==1) return (x%n+n)%n;  // 小技巧     else return -1;  }  

到此就结束了，下面稍稍浪一下，由求逆元的解法可以拓展到用扩展欧几里德求不定方程ax≡c(mod n)即ax+ny=c的解，同理如果c不能整除gcd(a,n)那么方程无解，先用扩展欧几里德求a*x + n*y = gcd(a,n)的特解x0，y0,然后原方程的特解就为x1=x0*c/gcd(a,n)，y1=y0*c/gcd(a,n)。所以原不定方程的通解为 x = x1 + (n/gcd(a,n))*t，   y = y1 –(a/gcd(a,n))*t