第 9 章 组合问题与大数运算 9.1 大数运算 在计算机中, 长整型变量的范围是 -2147483648 至 2147483647, 若用长整型变量做乘法运算, 乘积最多不能超过 10 位数 即便用双精度型 (double) 变量, 也仅能保证 16 位有效数字的精度 在某些需要更高精度的乘法运算的场合, 需要用别的办法来实现乘法运算 1. 列表法 例如, 当计 876 234 时, 把乘数与被乘数照表 9.1 列出 : 表 9.1 8 7 6 5 16 14 12 10 2 24 21 18 15 3 32 28 24 20 4 把表 9.1 中的数按图示斜线分组 ( 横纵坐标和相等的数分为一组 ), 把每组数累加起来所 得的和记在表格下方, 如表 9.2 所示 表 9.2 16 14 12 10 24 21 18 15 32 28 24 20 16 38 65 56 39 20 从最低位的 20 开始, 保留个位数字 0, 把个位以外的数 2 进到前一位 ; 把次低位的 39 加上低位进上来的 2 得 41, 保留个位数字 1, 把 4 进到前一位 ; 以此类推, 直至最高位的 16,16 加上低位进上来的 4 得 20, 保留 0, 把 2 进到最高位, 得乘积大数 2057010, 如表 9.3 所示 表 9.3 16 38 65 56 39 20 2 16+4=20 38+7=45 65+6=71 56+4=60 39+2=41 留 2 留 0 进 2 留 5 进 4 留 1 进 7 留 0 进 6 留 1 进 4 留 0 进 2 2 0 5 7 0 1 0 1
2. 程序实现根据以上思路就可以编写 C 程序了, 再经分析可得 : (1) 一个 m 位的整数与一个 位的整数相乘, 乘积为 m+-1 位或 m+ 位 (2) 程序中, 用三个字符数组分别存储乘数 被乘数与乘积 由第 (1) 步分析得知, 存放乘积的字符数组的长度应不小于存放乘数与被乘数的两个数组的长度之和 (3) 可以把第 (2) 步 计算填表 与第三 四步 累加进位 放在一起完成, 可以节省存储表格 9.2 所需的空间 (4) 程序关键部分是两层循环, 内层循环累计一组数的和, 外层循环处理保留数字与进位 程序代码如下 : #defie MAXLENGTH 1000 #iclude <stdio.h> #iclude <strig.h> void compute(char *a, char *b, char *c); void mai(void) char a[maxlength], b[maxlength], c[maxlength * 2]; puts("iput multiplier :"); gets(a); puts("iput multiplicad :"); gets(b); compute(a, b, c); puts("aswer :"); puts(c); getchar(); void compute(char *a, char *b, char *c) it i, j, m, ; log sum, carry; m=strle(a) - 1; =strle(b) - 1; for (i=m; i >=0; i--) 2
a[i] -='0'; for (i=; i >=0; i--) b[i] -='0'; c[m++2]='\0'; carry=0; for (i=m+; i >=0; i--) /* i 为坐标和 */ sum=carry; if ((j=i - m) < 0) j=0; for ( ; j<=i && j<=; j++) /* j 为纵坐标 */ sum+=a[i-j] * b[j]; /* 累计一组数的和 */ c[i+1]=sum % 10+'0'; /* 算出保留的数字 */ carry=sum / 10; /* 算出进位 */ if ((c[0]=carry+'0')=='0') /* if o carry, */ c[0]='\040'; /* c[0] equals to space */ 用以上算法计算 m 位整数乘以 位整数, 需要先进行 m 次乘法运算, 再进行约 m+ 次加法运算和 m+ 次取模运算 ( 实为整数除法 ) 把这个程序稍加修改, 让它自己产生乘数与被乘数, 然后计算随机的 7200 位整数互乘, 在 VC6.0 下耗时 4.5 s 注意以下事实 :8216547 96785 将两数从个位起, 每 3 位分为节, 列出乘法表, 将斜线 3
间的数字相加, 如表 9.4 所示 8 216 547 96 785 表 9.4 8 216 547 768 20736 52512 96 6250 169560 429395 785 768 20736 52512 6250 169560 429395 768 27016 222072 429385 将表中最后一行进行如下处理 : 从个位数开始, 每一个方格里只保留三位数字, 超出 1000 的部分进位到前一个方格里, 如表 9.5 所示 表 9.5 768 27016 222072 429385 768+27 27016+222 222072+429 取 395 进 429 =795 =27238 =222501 395 795 238 501 395 所以 8216547 96785=795238501395 也就是说在计算生成这个二维表时, 不必 1 位 1 位地乘, 而可以 3 位 3 位地乘 ; 在累加时也是满 1000 进位 这样, 计算 m 位整数乘以 位整数只需要进行 m /9 次乘法运算, 再进行约 (m+)/3 次加法运算和 (m+)/3 次取模运算 总体看来, 效率约是前一种算法的 9 倍 有人可能会想 : 既然能够 3 位 3 位地乘, 为什么不 4 位 4 位甚至 5 位 5 位地乘呢? 那不是可以提高 16 乃至 25 倍效率吗? 原因在于 : 本算法在累加表中斜线间的数字时, 如果用无符号长整数 ( 范围 0~4294967295) 作为累加变量, 在最不利的情况下 ( 两个乘数的所有数字均是 9), 能够累加约 4294967295/(999*999)=4300 次, 也就是能够准确计算任意两个均不超过 12900( 每次累加的结果 值 三位, 故 4300*3=12900) 位的整数相乘 如果 4 位 4 位地乘, 在最不利的情况下, 能够累加约 4294967295/(9999*9999)=43 次, 仅能够确保任意两个不超过 172 位的整数相乘, 没有什么实用价值, 更不要说 5 位了 程序代码如下 : #iclude<coio.h> #iclude<strig.h> #iclude<stdio.h> 4
#iclude<stdlib.h> #iclude<time.h> #defie N 7200 // 作 72xx 位的整数乘法 it max(it,it,it); it iitarray(it a[]); void write(it a[],it l); FILE *fp; void mai() it a[5000]=0,b[5000]=0,k[10001]=0; // 声明存放乘数 被乘数与积的数组 clock_t start, ed; // 声明用于计时的变量 usiged log c,d,e; // 声明作累加用的无符号长整数变量 it i,j,la,lb,ma,mi,p,q,t; // 声明其他变量 radomize(); // 初始化随机数 la=iitarray(a); // 产生被乘数, 并返回其长度 lb=iitarray(b); // 产生乘数, 并返回其长度 if(la<lb) // 如果被乘数长度小于乘数, 则交换被乘数与乘数 p=(lb>la)?lb:la; for (q=0;q<p;q++) // 交换被乘数与乘数 t=a[q],a[q]=b[q],b[q]=t; t=la,la=lb,lb=t; // 交换被乘数的长度与乘数的长度 5
start=clock();// 开始计时 c=d=0; // 清空累加变量, 其中 c 用于累加斜线间的数,d 用作进位标志 for(i=la+lb-2;i>=0;i--) // 累加斜线间的数,i 为横纵坐标之和 c=d; // 将前一位的进位标志存入累加变量 c ma=max(0,i-la+1,i-lb+1); // 求累加的下限 mi=(i>la-1)?(la-1):i; // 求累加的上限 for(j=ma;j<=mi;j++) // 计算出横纵坐标之和为 i 的单元内的数, 并累加到 c 中 c+=(log)a[j]*b[i-j]; d=c/1000; // 求进位标志 if(c>999) c%=1000; // 取 c 的末三位 k[i]=c; // 保存至表示乘积的数组 k e=k[0]+1000*d; // 求出乘积的最高位 ed=clock();// 停止计时 fp=fope("result.txt", "w+"); pritf("\the elapsed time was: %3.4f\", (ed-start)/clk_tck); // 打印消耗的时间 6
fpritf(fp,"%d",a[0]); // 打印被乘数最高位 write(a,la); // 打印被乘数其他位 fpritf(fp,"%d",b[0]); // 打印乘数最高位 write(b,lb); // 打印乘数其他位 fpritf(fp,"%ld",e); // 打印乘积最高位 write(k,la+lb-1); // 打印乘积其他位 fclose(fp); max(it a,it b,it c) it d; d=(a>b)?a:b; retur (d>c)?d:c; it iitarray(it a[]) it q,p,i; q=n+radom(100); if(q%3==0) p=q/3; else p=q/3+1; for(i=0;i<p;i++) a[i]=radom(1000); if(q%3==0) 7
a[0]=100+radom(900); if(q%3==2) a[0]=10+radom(90); if(q%3==1) a[0]=1+radom(9); retur p; void write(it a[],it l) it i; char strig[10]; for(i=1;i<l;i++) itoa(a[i],strig,10); if (strle(strig)==1) fpritf(fp,"00"); if (strle(strig)==2) fpritf(fp,"0"); fpritf(fp,"%s",strig); if((i+1)%25==0) fpritf(fp,"\"); fpritf(fp,"\"); fpritf(fp,"\"); 9.2 平面幻方 对平面幻方的构造, 分为三种情况 :N 为奇数 N 为 4 的倍数 N 为其他偶数 (4+2 的形式 ) 1. 奇数阶幻方 8
奇数阶幻方最经典的填法是罗伯特法 ( 也有人称之为楼梯方 ) 填写方法是: 把 1( 或最小的数 ) 放在第一行正中 ; 按以下规律排列剩下的 2-1 个数 : (1) 每一个数放在前一个数的右上一格 ; (2) 如果这个数所要放的格已经超出了顶行, 那么就把它放在底行, 仍然要放在右一列 ; (3) 如果这个数所要放的格已经超出了最右列那么就把它放在最左列, 仍然要放在上一行 ; (4) 如果这个数所要放的格已经超出了顶行且超出了最右列, 那么就把它放在前一个数的下一行同一列的格内 ; (5) 如果这个数所要放的格已经有数填入, 处理方法同 (4) 2. N 为 4 的倍数时互补 : 如果两个数字的和等于幻方最大数和最小数的和, 即 *+1, 称为互补 先看看 4 阶幻方的填法 : 将数字从左到右 从上到下按顺序填写 : 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 将对角线上的数字, 换成与它互补的数字 这里,*+1=4*4+1=17 把 1 换成 17-1=16; 把 6 换成 17-6=11; 把 11 换成 17-11=6, 换完后就是一个四阶幻方 : 16 2 3 13 5 11 10 8 9 7 6 12 4 14 15 1 对于 =4k 阶幻方, 先把数字按顺序填写 写好后, 按 4*4 把它划分成 k*k 个方阵 因为 是 4 的倍数, 一定能用 4*4 的小方阵分割 然后把每个小方阵的对角线, 像制作 4 阶幻方的方法一样, 对角线上的数字换成互补的数字, 就构成幻方 下面是 8 阶幻方的做法 : (1) 先把数字按顺序填 然后, 按 4*4 把它分割成 2*2 个小方阵, 如表 9.6 所示 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 9
(2) 每个小方阵对角线上的数字, 换成与它互补的数 表 9.6 64 2 3 61 60 6 7 57 9 55 54 12 13 51 50 16 17 47 46 20 21 43 42 24 40 26 27 37 36 30 31 33 32 34 35 29 28 38 39 25 41 23 22 44 45 19 18 48 49 15 14 52 53 11 10 56 8 58 59 5 4 62 63 1 程序代码如下 : #iclude"stdio.h" #iclude"math.h" it a[256][256]; it sum; it check(); void is(it ); mai() it i,j,,k,t,p,x; scaf("%d",&); sum=(*+1)*/2; if(%2==1)// 奇数幻方 is(); k=; if(%4==2)// 单偶数幻方 k=/2; is(k); 10
for(i=0;i<k;i++) for(j=0;j<k;j++) a[j+k]=a[j]+2*k*k; a[i+k][j]=a[j]+3*k*k; a[i+k][j+k]=a[j]+k*k; t=(-2)/4; for(i=0;i<k;i++) for(j=0;j<k;j++) if((j<t)&&(i<t)) p=a[j];a[j]=a[i+k][j];a[i+k][j]=p; if((j<t)&&(i>k-t-1)) p=a[j];a[j]=a[i+k][j];a[i+k][j]=p; if((i>=t&&i<=k-t-1)&&(j>=t&&j<t*2)) p=a[j];a[j]=a[i+k][j];a[i+k][j]=p; if(j>1&&j<=t) 11
p=a[j+k];a[j+k]=a[i+k][j+k];a[i+k][j+k]=p; if(%4==0)// 双偶数幻方 x=1; for(i=0;i<;i++) for(j=0;j<;j++) a[j]=x++; for(i=0;i<;i++) for(j=0;j<;j++) if(i%4==0&&abs(i-j)%4==0) for(k=0;k<4;k++) a[i+k][j+k]=*-a[i+k][j+k]+1; else if(i%4==3&&(i+j)%4==3) for(k=0;k<4;k++) a[i-k][j+k]=*-a[i-k][j+k]+1; if(check()==1) for(i=0;i<;i++) 12
for(j=0;j<;j++) pritf("%5d",a[j]); pritf("\"); retur ; it check(it )// 检验是否是幻方 it i,j,sum1=0,sum2; for(i=0;i<;i++) for(j=0;j<;j++) sum1+=a[j]; if(sum1!=sum) retur 0; sum1=0; for(i=0;i<;i++) for(j=0;j<;j++) sum1+=a[j]; if(sum1!=sum) retur 0; sum1=0; for(sum1=0,sum2=0,i=0,j=0;i<;i++,j++) 13
sum1+=a[j]; sum2+=a[-j-1]; if(sum1!=sum) retur 0; if(sum2!=sum) retur 0; else retur 1; void is(it )// 单偶数幻方的输入 it x,y,m; x=0;y=/2; for(m=1;m<=*;m++) a[x][y]=m; if (m%!=0) x--;y++; if(x<0) x=x+; if(y==) y=-y; else x++; if(x==) x=x-; 9.3 Catala 数 Catala 序列是一个整数序列, 其通项公式是从中取出的 C, 也叫做第 个 Catala 数, 前 几个 Catala 数是 :1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796, 58786, 208012, 742900, 2674440, 9694845, 35357670,, 乍看之下没什么特别的, 但是 Catala 数却是许多计数问题的最终形式 14
1. Catala 数的性质 1 2 (2 )! Catala 数的基本公式 C, 0, 但是却有一些变形和具体的性质 : 1 ( 1)!! 2 2 (1) C, 0 1 这是根据原来的式子推导出来的, 推导过程如下 : 2(2 1) (2) C O 1, C 1 C 2 这个递推式可以从公式 (1) 推得 C 1 2 2 2 2 2 1 1 1 (3) C O 1, C 1 CC i i 0 i 0 2 1 (4) C 1 i 0 i 4 (5) C ~ 3 2 π 这个是 Catala 数的增长趋势 个 +1 和 个 -1 构成 2 项 : a1, a2,, a, 其和满足 a1 a2 a k 0, 0 k 2的序 列个数等于第 个 Catala 数 C 假设不满足条件的序列个数为 U, 那么就有 C 2 U 剩下的工作就是求 U 了, 假设有一个最小的 k, 令 a1 a2 a k 0 由于这里 k 是最小的, 所以必有 a1 a2 ak 1 0, a 1, 并且 k 是一个奇数 此时将前 k 项中的 +1 变为 -1, 将 -1 变为 +1, 那么就得到一 k 个有 (+1) 个 +1 和 (-1) 个 -1 的序列了, 这样的序列个数就是要求的 U, 数值大小为 2 U 1 那么就得到了 2 2 2 C U 1 2. 应用 (1) 对于任意的 k, 前 k 个元素中 -1 的个数小等于 +1 的个数的序列计数, 可以不停地 变换形式, 比如将 -1 看成右括号,+1 看成左括号, 就变成了合法括号表达式的个数 比如 2 个左括号和 2 个右括号组成的合法表达式有 C2 2 种, 是 ()() 和 (()) (2) 既然如上一点都把括号加上去了, 那么顺便就再次转换,+1 个数连乘, 乘法顺序 有 C 种, 比如三个数连乘 a*b*c, 那么等于在式子上加括号, 有 2 种乘法顺序, 分别是 (ab)c 和 a(bc) 貌似对应关系比较模糊, 取 为 3 来看看, 为 3 的时候就是 4 个数相乘了, 设为 abcd, 最初的标号定在 a 上, 对于 为 3 得到合法的括号序列有 5 个, 分别是 :((())),()(()), 15