ACM必须的算法+数论.doc

ACM必须的算法欧几里德#includeusing namespace std;int hcf(int a,int b) int r=0; while(b!=0) r=a%b; a=b; b=r; return(a); lcd(int u,int v,int h) /u=a，v=b，h为最小公约数=hcf(a,b); return(u*v/h);int main()int a,b,x,y;cinab;x=hcf(a,b);y=lcd(a,b,x);coutx yendl;return 0;扩展欧几里德#include using namespace std;_int64 ext_euclid(_int64 a,_int64 b, _int64 &x, _int64 &y) int t;_int64 d; if (b=0) x=1;y=0;return a; d=ext_euclid(b,a %b,x,y); t=x; x=y; y=t-a/b*y; return d;void modular_equation(_int64 a,_int64 b,_int64 c)/ax = b(mod n) _int64 d; _int64 x,y; d = ext_euclid(a,b,x,y); if ( c % d != 0 )printf(No answern); else x = (x * c/d) % b ;/ 第一次求出的x ; _int64 t = b/d;x = (x%t + t)%t;printf(%I64dn,x);/最小的正数的值for (int i=0;id;i+)printf(The %dth answer is : %ldn,i+1,(x+i*(b/d)%b); /所有的正数值 /*函数返回值为gcd(a,b),并顺带解出ax+by=gcd(a,b)的一个解x,y, 对于不定方程ax+by=c的通解为： x=x*c/d+b/d*t y=y*c/d+a/d*t 当c%gcd(a,b)!=0时,不定方程无解.*/ 中国剩余定理#include using namespace std; int ext_euclid(int a,int b,int &x,int &y) /求gcd(a,b)=ax+by int t,d; if (b=0) x=1;y=0;return a; d=ext_euclid(b,a %b,x,y); t=x; x=y; y=t-a/b*y; return d;int China(int W,int B,int k) /W为按多少排列，B为剩余个数 WB K为组数 int i; int d,x,y,a=0,m,n=1; for (i = 0;ik;i+) n *= Wi; for (i=0;i0)return a; elsereturn(a+n);int main()int B100,W100; 求int k ; a = 2( mod 5 )cin k ; a = 3( mod 13)for(int i = 0 ; i Wi; 5 2cin Bi; 13 3 输出 42coutChina(W,B,k)endl;return 0;欧拉函数求小于n的所有欧拉数#include using namespace std;int phi1000; /数组中储存每个数的欧拉数void genPhi(int n)/求出比n小的每一个数的欧拉数(n-1的)int i, j, pNum = 0 ;memset(phi, 0, sizeof(phi) ;phi1 = 1 ;for(i = 2; i n; i+)if(!phii)for(j = i; j n; j += i)if(!phij)phij = j;phij = phij / i * (i - 1);求n的欧拉数int eular(int n)int ret=1,i;for (i=2;i*i1)ret*=n-1;return ret; /n的欧拉数行列式计算#include using namespace std;int js(int s100100,int n) int z,j,k,r,total=0; int b100100; /*bNN用于存放，在矩阵sNN中元素s0的余子式*/ if(n2) for(z=0;zn;z+) for(j=0;jn-1;j+) for(k=0;k=z) bjk=sj+1k+1; else bjk=sj+1k; if(z%2=0) r=s0z*js(b,n-1); /*递归调用*/ else r=(-1)*s0z*js(b,n-1); total=total+r; else if(n=2)total=s00*s11-s01*s10; return total; 排列long A(long n,long m) /nm long a=1; while(m!=0) a*=n;n-;m-; return a; 组合long C(long n,long m) /nm long i,c=1; i=m; while(i!=0) c*=n;n-;i-;while(m!=0) c/=m;m-; return c; 大数乘大数#include #include using namespace std;char a1000,b1000,s10000;void mult(char a,char b,char s) /a被乘数，b乘数，s为积 int i,j,k=0,alen,blen,sum=0,res6565=0,flag=0; char result65; alen=strlen(a);blen=strlen(b); for (i=0;ialen;i+) for (j=0;j=0;i-) for (j=blen-1;j=0;j-) sum=sum+resi+blen-j-1j; resultk=sum%10; k=k+1; sum=sum/10; for (i=blen-2;i=0;i-) for (j=0;j=i;j+) sum=sum+resi-jj; resultk=sum%10; k=k+1; sum=sum/10; if (sum!=0) resultk=sum;k=k+1; for (i=0;i=0;i-) si=resultk-1-i; sk=0; while(1) if (strlen(s)!=strlen(a)&s0=0) strcpy(s,s+1); else break; int main()cinab;mult(a,b,s);coutsendl;return 0;大数乘小数#include using namespace std;char a100,t1000;void mult(char c,int m,char t) / c为大数,m=10,t为积 int i,l,k,flag,add=0; char s100; l=strlen(c); for (i=0;il;i+) sl-i-1=ci-0; for (i=0;i=10) si=k%10;add=k/10;flag=1; else si=k;flag=0;add=0; if (flag) l=i+1;si=add; else l=i; for (i=0;iai;mult(a,i,t);couttendl;return 0;大数加法#include #include using namespace std;char a1000,b1000,s10000;void add(char a,char b,char s)/a被加数，b加数，s和 int i,j,k,up,x,y,z,l; char *c; if (strlen(a)strlen(b) l=strlen(a)+2; else l=strlen(b)+2; c=(char *) malloc(l*sizeof(char); i=strlen(a)-1; j=strlen(b)-1; k=0;up=0; while(i=0|j=0)if(i0) x=0; else x=ai;if(j9) up=1;z%=10; else up=0;ck+=z+0;i-;j-; if(up) ck+=1; i=0; ck=0; for(k-=1;k=0;k-) si+=ck; si=0; int main()cinab;add(a,b,s);coutsendl;return 0;大数减法#include using namespace std;char a1000,b1000,s10000;void sub(char a,char b,char s) int i,l2,l1,k; l2=strlen(b);l1=strlen(a); sl1=0;l1-; for (i=l2-1;i=0;i-,l1-) if (al1-bi=0) sl1=al1-bi+0; else sl1=10+al1-bi+0; al1-1=bl1-1-1; k=l1; while(ak=0) sl1=al1;l1-;loop: if (s0=0) l1=strlen(a); for (i=0;iab;sub(a,b,s);coutsendl;return 0;大数阶乘#include #include using namespace std;long a10000;int factorial(int n) /n为所求阶乘的n!的nint i,j,c,m=0,w; a0=1; for(i=1;i=n;i+) c=0; for(j=0;j0) m+;am=c; w=m*4+log10(am)+1;printf(%ld,am); / 输出for(i=m-1;i=0;i-) /printf(%4.4ld,ai);/printf(n);return w; /返回值为阶乘的位数储存方法很巧，每一个ai中存四位，不足四位在前加0补齐大数求余int mod(int B) /A为大数，B为小数int i = 0,r = 0;while( Ai != 0 )r=(r*10+Ai+-0)%B; return r ; /r为余数高精度任意进制转换#include #include using namespace std;char s1000,s21000; / s:原进制数字，用字符串表示，s2:转换结果，用字符串表示long d1,d2; / d1:原进制数，d2:需要转换到的进制数void conversion(char s,char s2,long d1,long d2) long i,j,t,num; char c; num=0; for (i=0;si!=0;i+) if (si=0) t=si-0; else t=si-A+10; num=num*d1+t; i=0; while(1) t=num%d2; if (t=9) s2i=t+0; else s2i=t+A-10; num/=d2; if (num=0) break; i+; for (j=0;jsd1d2;conversion(s,s2,d1,d2);couts2endl;return 0;判断一个数是否素数#include /基本方法，n为所求数，返回1位素数，0为合数#include using namespace std;int comp(int n)int i,flag=1; for (i=2;i=sqrt(n);i+) if (n%i=0) flag=0;break; if (flag=1) return 1; else return 0;素数表int prime(int a,int n) /小于n的素数 int i,j,k,x,num,*b; n+; n/=2; b=(int *)malloc(sizeof(int)*(n+1)*2); a0=2;a1=3;num=2; for(i=1;i=2*n;i+) bi=0; for(i=3;i=n;i+=3) for(j=0;j2;j+) x=2*(i+j)-1; while(bx=0) anum+=x; for(k=x;k=2*n;k+=x) bk=1; return num; /小于n的素数的个数bool flag1000000;void prime(int M) /01表int i , j;int sq = sqrt(double(M);for(i = 0 ;i M ;i +)flagi = true;flag1 = false;flag0 = false;for(i = 2 ;i = sq ;i+)if(flagi)for(j = i*i ;j 2和正整数s,该算法出错概率 至多为2(-s)，因此，增大s可以减小出错概 率，一般取s=50就足够了。#include#include using namespace std;int Witness(int a, int n) int i, d = 1, x; for (i = ceil( log( (float) n - 1 ) / log(2.0) ) - 1; i = 0; i-) x = d; d = (d * d) % n; if ( (d = 1) & (x != 1) & (x != n-1) ) return 1; if ( ( (n - 1) & ( 10 )d = (d * a) % n; return (d = 1 ? 0 : 1); int Miller_Rabin(int n, int s) int j, a; for (j = 0; j x;coutMiller_Rabin(x , 50)endl;return 0;KMP#include #include using namespace std; char t1010,p100; /t为大字符串，p为小字符串 int next100;void sn()int j,k; next0=-1; j=1; do k=nextj-1; while(k!=-1 & pk!=pj) k=nextk; nextj=k+1; j+=1; while(pj!=0); int match(int x) /x是大字符串下标,从头开始为0，j为小字符串下标 int k=x,j=0; if(t0=0)return 0; while(tk!=0) while(j!=-1 & pj!=tk) j=nextj; k+=1;j+=1; if(pj=0)return k-j; /返回p所在的下标 return -1; /搜索失败返回-1 int main()int x=0;sn();int r=match( x );coutr ak-1 ，而 bk= ak + k ak-1 + k-1 = bk-1 ak-1 。所以性质1。成立。 2。任意操作都可将奇异局势变为非奇异局势。 事实上，若只改变奇异局势（ak，bk）的某一个分量，那么另一个分量不可能在其他奇异局势中，所以必然是非奇异局势。如果使（ak，bk）的两个分量同时减少，则由于其差不变，且不可能是其他奇异局势的差，因此也是非奇异局势。 3。采用适当的方法，可以将非奇异局势变为奇异局势。 假设面对的局势是（a,b），若 b = a，则同时从两堆中取走 a 个物体，就变为了奇异局势（0，0）；如果a = ak ，b bk，那么，取走b - bk个物体，即变为奇异局势；如果 a = ak ， b ak ，b= ak + k,则从第一堆中拿走多余的数量a - ak 即可；如果a ak ，b= ak + k,分两种情况，第一种，a=aj （j k）,从第二堆里面拿走 b - bj 即可；第二种，a=bj （j k）,从第二堆里面拿走 b - aj 即可。 从如上性质可知，两个人如果都采用正确操作，那么面对非奇异局势，先拿者必胜；反之，则后拿者取胜。 那么任给一个局势（a，b），怎样判断它是不是奇异局势呢？我们有如下公式： ak =k（1+5）/2，bk= ak + k （k=0，1，2，.,n 方括号表示取整函数)奇妙的是其中出现了黄金分割数（1+5）/2 = 1。618.,因此,由ak，bk组成的矩形近似为黄金矩形，由于2/（1+5）=（5-1）/2，可以先求出j=a（5-1）/2，若a=j（1+5）/2，那么a = aj，bj = aj + j，若不等于，那么a = aj+1，bj+1 = aj+1+ j + 1，若都不是，那么就不是奇异局势。然后再按照上述法则进行，一定会遇到奇异局势。（三）尼姆博奕（Nimm Game）：有三堆各若干个物品，两个人轮流从某一堆取任意多的物品，规定每次至少取一个，多者不限，最后取光者得胜。 这种情况最有意思，它与二进制有密切关系，我们用（a，b，c）表示某种局势，首先（0，0，0）显然是奇异局势，无论谁面对奇异局势，都必然失败。第二种奇异局势是（0，n，n），只要与对手拿走一样多的物品，最后都将导致（0，0，0）。仔细分析一下，（1，2，3）也是奇异局势，无论对手如何拿，接下来都可以变为（0，n，n）的情形。 计算机算法里面有一种叫做按位模2加，也叫做异或的运算，我们用符号（+）表示这种运算。这种运算和一般加法不同的一点是1+1=0。先看（1，2，3）的按位模2加的结果：1 =二进制012 =二进制103 =二进制11 （+）0 =二进制00 （注意不进位） 对于奇异局势（0，n，n）也一样，结果也是0。 任何奇异局势（a，b，c）都有a（+）b（+）c =0。如果我们面对的是一个非奇异局势（a，b，c），要如何变为奇异局势呢？假设 a b(1,8,9)奇异局势 乙:(1,8,9)-(1,8,4) 甲:(1,8,4)-(1,5,4)奇异局势 乙:(1,5,4)-(1,4,4) 甲:(1,4,4)-(0,4,4)奇异局势 乙:(0,4,4)-(0,4,2) 甲:(0.4,2)-(0,2,2)奇异局势 乙:(0,2,2)-(0,2,1) 甲:(0,2,1)-(0,1,1)奇异局势 乙:(0,1,1)-(0,1,0) 甲:(0,1,0)-(0,0,0)奇异局势 甲胜。叉乘法求任意多边形面积语法：result=polygonarea(Point *polygon,int N);参数：*polygon：多变形顶点数组N：多边形顶点数目返回值：多边形面积注意：支持任意多边形，凹、凸皆可多边形顶点输入时按顺时针顺序排列源程序：typedef struct double x,y; Point; double polygonarea(Point *polygon,int N)int i,j;double area = 0;for (i=0;iN;i+) j = (i + 1) % N;area += polygoni.x * polygonj.y;area -= polygoni.y * polygonj.x;area /= 2;return(area 0 ? -area : area);求三角形面积语法：result=area3(float x1,float y1,float x2,float y2,float x3,float y3);参数：x13：三角形3个顶点x坐标y13：三角形3个顶点y坐标返回值：三角形面积注意：需要 math.h源程序：float area3(float x1,float y1,float x2,float y2,float x3,float y3)float a,b,c,p,s;a=sqrt(x1-x2)*(x1-x2)+(y1-y2)*(y1-y2);b=sqrt(x1-x3)*(x1-x3)+(y1-y3)*(y1-y3);c=sqrt(x3-x2)*(x3-x2)+(y3-y2)*(y3-y2);p=(a+b+c)/2;s=sqrt(p*(p-a)*(p-b)*(p-c);return s;两点距离（2D、3D）语法：result=distance_2d(float x1,float x2,float y1,float y2);参数：x/y/z12：各点的x、y、z坐标返回值：两点之间的距离注意：需要 math.h源程序：float distance_2d(float x1,float x2,float y1,float y2) return(sqrt(x1-x2)*(x1-x2)+(y1-y2)*(y1-y2);float distance_3d(float x1,float x2,float y1,float y2,float z1,float z2)return(sqrt(x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)+(z1-z2)*(z1-z2);射向法判断点是否在多边形内部语法：result=insidepolygon(Point *polygon,int N,Point p);参数：*polygon：多边形顶点数组N：多边形顶点个数p：被判断点返回值：0：点在多边形内部；1：点在多边形外部注意：若p点在多边形顶点或者边上，返回值不确定，需另行判断需要 math.h源程序：#define MIN(x,y) (x y ? x : y)typedef struct double x,y; Point;int insidepolygon(Point *polygon,int N,Point p) int counter = 0;int i;double xinters;Point p1,p2;p1 = polygon0;for (i=1;i MIN(p1.y,p2.y) if (p.y = MAX(p1.y,p2.y) if (p.x = MAX(p1.x,p2.x) if (p1.y != p2.y) xinters = (p.y-p1.y)*(p2.x-p1.x)/(p2.y-p1.y)+p1.x;if (p1.x = p2.x | p.x = xinters)counter+;p1 = p2;if (counter % 2 = 0)return(OUTSIDE);elsereturn(INSIDE);判断点是否在线段上语法：result=Pointonline(Point p1,Point p2,Point p);参数：p1、p2：线段的两个端点p：被判断点返回值：0：点在不在线段上；1：点在线段上注意：若p线段端点上返回1需要 math.h源程序：#define MIN(x,y) (x y ? x : y) typedef struct double x,y; Point;int FC(double x1,double x2)if (x1-x2-0.000002) return 1; else return 0;int Pointonline(Point p1,Point p2,Point p)double x1,y1,x2,y2;x1=p.x-p1.x;x2=p2.x-p1.x;y1=p.y-p1.y;y2=p2.y-p1.y;if (FC(x1*y2-x2*y1,0)=0) return 0;if (MIN(p1.x,p2.x)=p.x&p.x=MAX(p1.x,p2.x)&(MIN(p1.y,p2.y)=p.y&p.y=MAX(p1.y,p2.y)return 1; else return 0;判断两线段是否相交语法：result=sectintersect(Point p1,Point p2,Point p3,Point p4);参数：p14：两条线段的四个端点返回值：0：两线段不相交；1：两线段相交；2两线段首尾相接注意：p1!=p2;p3!=p4;源程序：#define MIN(x,y) (x y ? x : y) typedef struct double x,y; Point;int lineintersect(Point p1,Point p2,Point p3,Point p4)Point tp1,tp2,tp3;if (p1.x=p3.x&p1.y=p3.y)|(p1.x=p4.x&p1.y=p4.y)|(p2.x=p3.x&p2.y=p3.y)|(p2.x=p4.x&p2.y=p4.y) return 2;/快速排斥试验if (MIN(p1.x,p2.x)p3.x&p3.xMAX(p1.x,p2.x)&MIN(p1.y,p2.y)p3.yMAX(p1.y,p2.y)|(MIN(p1.x,p2.x)p4.x&p3.xMAX(p1.x,p2.x)&MIN(p1.y,p2.y)p3.y=0) return 1; else return 0;点到线段最短距离语法：result=mindistance(Point p1,Point p2,Point q);参数：p1、p2：线段的两个端点q：判断点返回值：点q到线段p1p2的距离注意：需要 math.h源程序：#define MIN(x,y) (x y ? x : y)typedef struct double x,y; Point;double mindistance(Point p1,Point p2,Point q)int flag=1;double k;Point s;if (p1.x=p2.x) s.x=p1.x;s.y=q.y;flag=0;if (p1.y=p2.y) s.x=q.x;s.