hduoj 4196 ——Remoteland

这是2011 Southwestern Europe Regional Contest Problems的一道题目，数论题。但关键点是程序的常数级优化。

【题目大意】：找一个数D，是一个平方数，且它是由<=n中的若干个数构成。D要尽可能大，求D是多少。答案mod  1000000007

思路倒是不难，我们对n！进行标准分解，如果幂次大于1，则保留偶次幂，奇次幂减去1，幂次小于1的直接break掉即可。比如4！=2²*3则我们就只需保留2²即可。此时结果就是4了。

^{下面简要说明下这样做的正确性：每个正整数由算术基本定理可以知道都可以分解为标准形式，那么n！当然也可以分解为标准形式了。接着我们进行了上述的操作之后我们肯定能保证剩下的就是一个平方数了。}

^{但是这里实际上是有一个问题的，我们考虑题目意思的时候，选取每个数字时是作为一个整体选进去的，也就是选取的时候同时把它的标准分解式选取了进去，但是我们进行上述的调整调整操作的时候却是以单个素因子进行的，那么就会有一个潜在的问题了，我们会不会在选取素因子的时候把一个数的某个素因子选取了进去，却没有把其它的素因子选取进去呢？比如34=2×17，我们肯能会选取了2，却没有选取17，但是对于34只有选与不选，不能只要人家的某个因子，是吧？当时我纠结于此，以为这个题目难度很高，还在想怎么确保选取的时候不会出现这种情况。}

^{不过后来考虑到换一种思路，立即可以证明这种担心是多余的：我们考虑对于每个我们要去掉的素因子，我们把它们写成一行，比如：}

^{23！=219*39*54*73*112*13*17*19*23}

^{进行上述操作之后，将要去掉的就是：}

^{2*3*7*13*17*19*23}

^{注意到上述这些数字必然会在1～23中出现的，那么我们实际上不取这写数字就可以了。这样就可以清晰的说明必然能够选取适当的数使其乘积为上述操作后留下的那个数。}

^{那么怎么证明它是最大的呢？}

^{这个问题即我们找不到一个更大的满足条件的数。依然考虑23！，这个结果最大即为23！，但是显然23！不满足是平方数的要求，对于2这个素因子，我们肯定不能把它的指数加1变为偶数的，那将会大于23！，只能减去1变为偶数。同理对于剩下的素因子我们有相同的考虑，最后实际上就是我们讲到的上述操作了，又由于素因子自己的不同类，所以即证明了这样取的数就是要求的数。}

^{思路理清之后，代码自然不难写，但是这道题目时间卡的很紧，需要各种优化，我的代码思路是，预处理素数表，对n！进行标准分解，对素因子幂进行二分快速幂。重要优化之一就是筛法没有必要筛到10000000，而是5000000，大概能节约几百ms，另一个就是进行二分幂的时候并不是对每个素因子都进行二分幂的，而是把指数相同的算一个进行一次二分幂，这个优化保证了我的代码没有TLE！}

^{PS：网上看到一种做法说mod1000000007要用到乘法逆元，那么乘法逆元不是解决除法取模时用到的吗？现在完全没有除法啊，怎么会用到逆元呢？不解，有知道的大神请指教啊。}

 
 
#include <stdio.h>
#include <math.h>
#define max 5000500//不需要弄到10000000.
bool prime[max]={1,1};
int a[348540];
int main()
{
	int i,j,m,n,k,t;
	long long ret,rett;
	m=sqrt(max+0.0001);
	for(i=2;i<m;++i)
		if(prime[i]==0)
			for(j=i*i;j<max;j+=i)
				prime[j]=1;
	for(i=2,j=0;i<max;++i)
		if(prime[i]==0)
			a[j++]=i;//好经典的筛法代码啊！
	while(scanf("%d",&n),n)
	{
		ret=1;
		rett=1;
		t=0;
		for(i=0;i<348539;++i)
		{
			m=n;j=0;
			while(m)
			{
				j+=m/a[i];
				m=m/a[i];
			}
			if(j>1)
			{
				j=j&1?j-1:j;
				if(j==t)
					rett=rett*a[i]%1000000007;//不是进来一个素因子就进行二分幂的，而是先“存”起来。
				else
				{
					while(t)
					{      
						if(t&1)
							ret=ret*rett%1000000007;
						rett=rett*rett%1000000007;//发现大数取模运算很费时间的。
						t=t>>1;
					}
					t=j;
					rett=a[i];
				}
			}
			else
			{
				while(t)
				{      
					if(t&1)
						ret=ret*rett%1000000007;
					rett=rett*rett%1000000007;
					t=t>>1;
				}
				break;
			}
		}
		printf("%I64d\n",ret);
	}
}