梅森素数缺点（梅森数与梅森素数）

在古希腊时代人们就发现了4个完全数：6,28,496,8128，它们之间相差不大，但是直到15世纪才发现了第五个完全数33550336，和前面相差竟然如此遥远，有没有简单法则寻找完全数呢？其实公元前300年欧几里得已经得出了一条定理：若2^p-1是素数（也叫质数），则2^(p-1)*(2^p-1)是完全数。因此，对于一个新的p,只要能判断2^p-1 是素数，就相当于找到了一个新的完全数。

梅森是17世纪欧洲科学界一位独特的中心人物，是法兰西学院的奠基人，被选为“100位在世界科学史上有重要地位的科学家”之一。他在欧几里得、费尔马等人研究的基础上对2^p-1型的数做了大量的计算和验证，于1644年在他的《物理数学随感》一书中断言：在不大于257的素数中，当p=2、3、5、7、13、17、19、31、67、127、257时，2^p-1是素数，其它都是合数。前面的7个数（即p=2、3、5、7、13、17、19）已经被前人所证实，而后面的4个数（p=31、67、127、257）则是梅森自己的推断，由于梅森在科学界有崇高的学术地位，当时的人们对其断言都深信不疑（后来人们发现其断言包含若干错漏）。

由于梅森的工作极大地激发了人们研究2^p-1型素数的热情，使其摆脱了作为“完全数”的附庸地位，是一个转折点和里程碑，为了纪念他，数学界就把2^p-1，其中p是素数的数称为“梅森数”，并以Mp记之，即Mp=2^p-1，若p是素数，则2^p-1叫做梅森数，例如：2^2-1=3,2^3-1=7,2^5-1=31,2^7-1=127,前面4个都是素数；而2^11-1=2047=23*89是合数，2^13-1=8191，2^17-1=131071，2^19-1=52428，是素数；2^23-1=8388607=47*178481，又是合数……，如果2^p-1也是素数，则称其为梅森素数。

分解2^p-1 的结果展示：

C语言分解梅森数程序如下：

//分解梅森数Mp=2^p-1,其中p为素数

#include <stdio.h>

#include <math.h>

#include <time.h>

#define N 20

main () //注：如果(2^p-1)是合数,一定能被(2p*k 1)整除,其中k 为自然数中的某一个

{ void shuchu(int,unsigned [N]); //输出函数原型声明

int i,k,x,lbz,lb=1,lcz=1,lc=1; //循环变量i,k,x;被除数总位数lbz,单元数lb,除数总位数lcz,单元数lc

int p,ss,l,jw,g=0,jr=0; //指数p,试商ss,积的单元数l,进位jw,个数g,进入标志jr

int lb1,lc1,lc2,b5,q,c3; //lb1=lb-1,lc1=lc-1,lc2=lcz*2-1,被除数前5位b5及其算术根q,除数前3位c3

unsigned b[N]={0,1},c[N]={0,1},s[N]={},y[N*2]={},xj; //被除数b,除数c,商s,余数y,新积xj

printf("请输入素数p：");scanf("%u",&p);

float t0=clock();

//求2^p-1:

for(i=1;i<=p;i )

{ jw=0;

for(k=1;k<=lb;k )

{ b[k]=b[k]*2 jw;

if(b[k]<=9999)jw=0; else{jw=1;b[k]-=10000;} //每4位为1个单元

}

if(jw>=1){lb ;b[lb]=1;}

}

b[1]--;

printf("2^%u-1=",p);

shuchu(lb,b); //调用输出函数：输出被分解数b的各单元(lb为单元数,b为数组名)

//开始分解：

printf(" = 1");

p*=2;c[1] =p; c3=c[1];//求(2*p 1)

while (lcz<=lbz)

{ lc2=lcz*2;

if(lc2<lbz||(lc2==lbz&&c3<=q))

{ lbz=lb*4;lb1=lb-1;

if(b[lb]>=1000) {b5=b[lb]*10 b[lb1]/1000;}

else if(b[lb]>=100) {lbz--;b5=b[lb]*100 b[lb1]/100;}

else if(b[lb]>=10) {lbz-=2;b5=b[lb]*1000 b[lb1]/10;}

else {lbz-=3;b5=b[lb]*10000 b[lb1];}

q=sqrt(b5 1); lc1=lc-1;

for(x=1;x<=lb;x ) {y[x]=b[x];} //做除法：

for(i=lb;i>=lc;i--)

{ y[i] =y[i 1]*10000;y[i 1]=0; s[i]=0;

while(y[i]>c[lc])

{ if(y[i]>=429496) ss=y[i]/(c[lc] 1); //2^32=42 9496 7296

else ss=(y[i]*10000 y[i-1])/(c[lc]*10000 c[lc1] 1);

if(ss==0) ss=1;

jw=0;s[i] =ss;

for(k=1;k<=lc1;k ) //lc1=lc-1

{ xj=c[k]*ss jw;

if(xj<=9999)jw=0; else{jw=xj/10000;xj%=10000;}

l=k i-lc;

if(y[l]<xj) {y[l] =10000;y[l 1]--;}

y[l]-=xj;

}

xj=c[lc]*ss jw;

y[i]-=xj;

}

}

while(y[lc]>=c[lc])

{ for(x=lc;x>=1;x--)

{ if(y[x]>c[x]) break;

if(y[x]<c[x]) goto tc;

}

s[lc] ;

for(x=1;x<=lc1;x ) //lc1=lc-1

{ if(y[x]<c[x]) {y[x] =10000;y[x 1]--;}

y[x]-=c[x];

}

y[lc]-=c[lc];

}

tc:

for(x=lc;x>=1;x--)

{ if(y[x]!=0) break;}

if(x!=0) //余数 !=0,求新的除数：

{ c[1] =p; g ;

if(jr!=0)

{ if(gQ051!=0) //因51051=3*7*11*13*17

{ while((g%3==0||c[1]%5==0||g%7==0||g==0||g==0||g==0)==1)

{ g ;c[1] =p; }

}

else {g=1;c[1] =p;}

}

else {if((c[2]*10000 c[1])Q051==0) {jr=1;g=1;c[1] =p;}}

if(c[1]>=10000)

{ c[2] ;c[1]-=10000;

for(x=2;x<=lc;x ) { if(c[x]>=10000){c[x 1] ;c[x]-=10000;} }

if(c[lc 1]>=1) lc ;

lcz=lc*4; lc1=lc-1;

if(c[lc]>=1000) {c3=c[lc]/10;}

else if(c[lc]>=100) {lcz--;c3=c[lc];}

else if(c[lc]>=10) {lcz-=2;c3=c[lc]*10 c[lc1]/1000;}

else {lcz-=3;c3=c[lc]*100 c[lc1]/100;}

}

}

else // 余数 =0,输出因数：

{ printf(" * "); shuchu(lc,c); //调用输出函数：输出因数c的各单元

if(s[lb]==0) lb--;lb1=lb-1;

for(x=lc;x<=lb1;x )

{ if(s[x]>=10000) {s[x 1] ;s[x]-=10000;}

}

for(x=lc;x<=lb;x ) {b[x-lc1]=s[x];} //求新的被除数：

lb-=lc1; //lc1=lc-1

}

}

else //分解完了输出最后一个因数：

{ printf(" * "); shuchu(lb,b); //调用输出函数：输出最后一个因数b的各单元

break;

}

}

printf("\n用时%.3f秒",(clock()-t0)/1000);

}

//输出函数：

void shuchu(int lb,unsigned b[N])

{ int x;

printf("%u",b[lb]); lb--;

for(x=lb;x>=1;x--)

{ if(b[x]>=1000) printf(" %u",b[x]);

else if(b[x]>=100) printf(" 0%u",b[x]);

else if(b[x]>=10) printf(" 00%u",b[x]);

else printf(" 000%u",b[x]);

}

}

,