批量分解素因数（二）

本文中的程序使用压缩的格式存储M以内(包含M)的所有奇数的分解式。

和前一个的程序不同，这个程序侧重于使用尽可能少的内存空间来保存M以内的所有奇数的分解式.
本文中的程序可计算65759以内的所有奇数的因子表，其因子表占用的内存空间不到200KB，32881*2+65535*2=196,832字节.
和上一篇文章一样，M以内的奇数的分解式需要使用两个表格来存储，地址表addrTable和因子表factorTable.
地址表使用WORD类型，addrTable[i]==0,表示(i*2+3)为素数,否则addrTable[i]表示(i*2+3)分解式的第一项的地址.
factorTable[0]不使用，factorTable的最后一项也不使用,故其存储的因子表的总数不能超过65534.
分解式的每一项，仅用一个WORD类型变量来存储，即每一项仅需2个字节。
当奇数x=2*i+3是素数时，其分解式无需存储在factorTable，这时，置addrTable[i]==0。
当奇数x=2*i+3是合数时，且其素因子p的指数为1时，在这种情况下，置bit15=1，用bit0-bit14存储p，这里p的最大取值是32767，故M的最大取值范围不超过3*32767=98,301.
当奇数x=2*i+3是合数时，且其素因子p的指数大于1，在这种情况下，置bit15=0，用bit0-bit8存储p，用bit9-bit14存储指数e，这里p的最大取值不超过511，故M的最大取值范围不超过511*511=261,121.
更多的细节，请参阅源代码。

下面给出源代码

#include <stdlib.h>
#include <stdio.h>
#include <string.h>

typedef unsigned char BYTE;
typedef unsigned long DWORD;
typedef unsigned short WORD;

//筛出n以内(包括n）的所有素数，并统计n以内所有自然数的非重复素因子的个数
//当函数完成后，若pBuff[i]==0, 表明i是素数，
//若pBuff[i]>0, 表明i是合数，且i的非重复素因子的个数是pBuff[i]

void sift(BYTE *pBuff,DWORD n)
{
	DWORD p,m;
	memset(pBuff,0,sizeof(BYTE)*(n+1));
	p=2;
	while (p<=n)
	{
		for (m=p*2;m<=n;m+=p)	//p素数，m是p的倍数，将p的倍数（不包括p）的非重复素因子的个数增加1
			pBuff[m]++;
		
		p++;	//调到下一个数
		while ( pBuff[p] && p<=n) //略过合数
			p++;
	}
}


// 函数 print_expression_on_zipTable输出n以内所有自然数的分解式
// addrTable是地址表，addrTable[i]==0,表示(i*2+3)是素数，否则addrTable[i]指向奇数(i*2+3)的因子表的地址,
// addrTable[i+1]-addrTable[i]表示奇数(i*2+3)的分解式的项数
// factorTable[addrTable[i]]为奇数(i*2+3)的分解式的第一项
void print_expression_on_zipTable(WORD *addrTable,WORD *factorTable,DWORD max)
{
	DWORD i,x,j,t,c;
	WORD *tab;
	
	for (x=3;x<=max;x+=2)
	{
		i=(x-3)/2;
		
		if ( addrTable[i]==0)
		{
			printf("%u=%u\n",x,x);
			continue;
		}
		
		printf("%u=",x);

		for (t=i+1;addrTable[t]==0;)
			t++; //略过素数
		c=addrTable[t]-addrTable[i];
		tab=factorTable+addrTable[i];
		for(j=0;j<c;j++)
		{
			DWORD p,e;

			if (j>0)
				printf("*");
			
			if ( tab[j] & 0x8000)
			{	p=( tab[j] & 0x7fff);  e=1;}
			else
			{	p=( tab[j] & 0x1ff); e=( (tab[j]>>9) & 0x3f);	}
			printf("(%u^%u)",p,e);
		}
		printf("\n");
	}
}

void build_zipFactorTable(DWORD max,BYTE *pBuff, const WORD *addrTable,WORD *factorTable)
{
	DWORD p,m,t;
	WORD e,item;
	int addrTableSize;
	WORD *writePointTab=NULL;
	
	//-----------------------------------------------------
	// addrTable[i]存储奇数(i*2+3)的分解式，故i的最大值为（max-3)/2
	// 另外，需要额外增加一项来存储max+2分解式的地址，以方便计算max分解式的项数
	// 故需要addrTableSize=(max-3)/2+2
	
	addrTableSize=(max-3)/2+2;	
	writePointTab=(WORD *)malloc(sizeof(WORD)*addrTableSize);
	if ( writePointTab==NULL)
	{
		printf("Alloc memory faield\n");
		return ;
	}
	memcpy(writePointTab,addrTable,sizeof(WORD)*addrTableSize);
	
	p=3;
	while (p<=max)
	{
		for (m=p*3;m<=max;m+=(p*2))	//m为p的奇数倍
		{
			e=1;	
			// m可能含有多个素因子p，如120含有3个因子2，下面的代码求出素因子p的指数
			// m是p的倍数，故p一定整除m，故素因子p的指数至少为1
			t=m/p;	
			while (t%p==0) {t/=p; e++; }
			if ( e==1)
				item= (WORD)(0x8000 | p);
			else
				item= (WORD)((e << 9) | p);
			factorTable[writePointTab[(m-3)/2]]=item;
			writePointTab[(m-3)/2]++;
		}
		
		p+=2;	//前进到下一个奇数
		while ( pBuff[p]>0 && p<=max) 
			p+=2;	//略过合数，使得p总是素数
	}
	
	free(writePointTab); writePointTab=NULL;
}

// 函数batchDecompose批量分解n以内所有奇数为素数乘积的形式，并打印其分解式
void batch_odd_Decompose(DWORD n,const char *fileName)
{
	BYTE *buff=NULL;		//用来得到n以内的所有素数，和所有合数的素因子的个数
	WORD *addrTable=NULL;	//地址表addrTable[i]表示奇数(i*2+3)的分解式在factorTable中的偏移地址
	WORD *factorTable=NULL;
	WORD nextPos;
	DWORD i,x,max;
	DWORD total,addrTableSize;

	if ( (n&1)==0)
		n++;	//前进到下一个奇数

	buff=(BYTE *)malloc(sizeof(BYTE)*(n+1));
	if (buff==NULL )
	{ printf("Alloc memory failed\n");	goto free_memory;}

	sift(buff,n);

	for (total=0,x=3;x<=n;x+=2)
	{
		if ( buff[x]>0)
		{
			//buff[x]为数x的分解式的项数，total为所有奇合数的分解式的项数和
			total += buff[x];	 
			
			// 因子表factorTable的地址使用WORD来表示，故最大容量为65536，factorTable[0]和最后一个元素不使用，
			// 故total必须<=65534
			if (total<=65534)	
				max=x; 				  
			else
			{
				total -= buff[x];
				break;
			}
		}
	}
	printf("max=%u\n",max);

	//addrTableSize[i]存储奇数i*2+3的分解式，故i的最大值为（max-3)/2
	//另外，需要额外增加一项来存储 max+2分解式的地址，故需要addrTableSize=(max-3)/2+2
	addrTableSize=(max-3)/2+2;	
	addrTable	=(WORD *)malloc(sizeof(WORD)*(addrTableSize));
	if (addrTable==NULL )
	{ printf("Alloc memory failed\n");	goto free_memory;}

	nextPos=1;
	for (x=3;x<=max;x+=2)
	{
		i=(x-3)/2;
		
		if ( buff[x]==0)		//x是素数
			addrTable[i]=0;
		else
		{
			addrTable[i]=nextPos;	
			nextPos += buff[x];
		}
	}

	i=(x-3)/2;
	addrTable[i]=nextPos;	//增加额外的一项，方便计算max分解式的项数

	factorTable=(WORD *)malloc(sizeof(WORD)*(total+2));
	if (factorTable==NULL)
	{ printf("Alloc memory failed\n");	goto free_memory; }

	build_zipFactorTable(max,buff,addrTable,factorTable);		//将max以内的奇数的分解式存储到factorTable

	print_expression_on_zipTable(addrTable,factorTable,max);	//打印将n以内所有自然数的分解式

free_memory:
	if (buff!=NULL) {free(buff);buff=NULL;}
	if (addrTable!=NULL) {free(addrTable);addrTable=NULL;}
	if (factorTable!=NULL) {free(factorTable);factorTable=NULL;}
}


int main(int argc, char* argv[])
{
	batch_odd_Decompose(3*32767,NULL);
	return 0;
}