POJ 1091主要是他28%的通过率和1900多个接收吸引了我。
题目大意:
Description
比如当N=2,M=18时,持有卡片(10, 15, 18)的跳蚤,就可以完成任务:他可以先向左跳10个单位长度,然后再连向左跳3次,每次15个单位长度,最后再向右连跳3次,每次18个单位长度。而持有卡片(12, 15, 18)的跳蚤,则怎么也不可能跳到距他左边一个单位长度的地方。
当确定N和M后,显然一共有M^N张不同的卡片。现在的问题是,在这所有的卡片中,有多少张可以完成任务。
Input
Output
分析了一下发现……这个跳蚤的题目果然很变态,虽然原理蛮简单,但是实现起来还是有点困难的。
首先仔细阅读题目,就可分析出来是求使得多元一次丢番图方程有解的系数的个数,即:
对于整数(N,M),求解使得a1*x1+a2*x2+a3*x3+……+aN*xN+M*x(N+1)=1 的方程有解的数组(a1,a2,a3……aN,M)的个数,其中ai<=M,(i=1,2,3……N)。
由数论的知识可以知道,方程有解,则 gcd(a1,a2,a3……aN,M)|1,即这里头的数两两互质。
当然了,如果说要直接去求两两互质的这个数组很定很麻烦,所以再继续分析。
这里就需要用到容斥原理了,gcd(a1,a2,a3……aN,M)|1 的数组的数目,应该是总数组的数目减去所有使得gcd(a1,a2,a3……aN,M)不等于1的数组的数目。由于发现数组中总有个固定值M,那么从M下手。于是在所有gcd不等于1的所有数组中,他们的gcd只可能等于M的约数,所以只要a1,a2,a3……aN是M的某个约数的倍数,那么这个数组的gdc肯定不是1。
M的约数还是比较多的,用容斥原理算起来还是麻烦。在分析,发现只要a1,a2,a3……aN是M的某个素因子的倍数就可以了,因为a1,a2,a3……aN是M的某个约数的倍数就肯定是M的某个素因子的倍数。而且在计算容斥原理的时候,只要找出不同的素因子就可以了。那么这么一分析,假设M的不同的素因子为pi,那么使得 gcd(a1,a2,a3……aN,M)=k*pi,的组合的数目就是 (M/pi)^N,所以最后的结果加上容斥原理就是:
M^N-(M/p1)^N-(M/p2)^N-……+(M/(p1*p2))^N……
果然算式很变态……难点有以下几个:
1.首先是寻找M的不同素数因子。这个用小于sqrt(M)的数试除就可以。亏好M <= 100000000,sqrt(M)<10000不算太大。
2.其次是容斥原理的实现,这个有点难度。我的想法是用一个整型变量 i,不断加1,然后记录其2进制表示中1的位数,就可以找出所有的容斥原理的组合了。(灵感来源于《离散数学及其应用第六版》 5.6.3 生成组合)
3.然后最最变态的是,居然还要你实现大数的幂方算法和加减法...这个我在实现的时候借用了我之前做POJ1001的代码,然后自己写了个大数的加减法代码。不过自认为只是实现了,并没有最优化。
4.还有一点,题目中只有(N <= 15 , M <= 100000000)的描述,还要处理一下N=0,M=0的情况。
5.关于各种算法性能的估计,对于求素数算法,由于M <= 100000000,所以他的不同素因子不会超过10个,因为2*3*5*7*9*11*13*17*19*23>1,000,000,000,单次循环最坏情况下不会超过10000,所以性能不会很差。对于容斥原理,不同因子不会超过10个,所以,容斥原理计数的整型变量不会大于2^10=1024,性能也不会太差。最后对于大数算法,最大位数是 15*10也就是150位,这个性能也不会很差。
所以,这个要实现的东西好多啊。。。还是我自己太水了觉得这么点就嫌多了。。。
代码如下:
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
typedef long long int da;
typedef struct __data{
char *data;
int len;
int sign;
} powDa;
void printPowDa(powDa a){
int len=a.len;
int i,prefix;
char *buffer=a.data;
for (i=len-1;i>=0;i--){
if (buffer[i]!='0'){
break;
}
}
prefix=i;
if (a.sign==-1){printf("-");}
for (i=prefix;i>=0;i--){
printf("%c",buffer[i]);
}
printf("\n");
}
powDa *charTimes(powDa *a,powDa *b){ //乘法函数
char *aChar=(a->data);
char *bChar=(b->data);
int aLen=a->len;
int bLen=b->len;
int ansLen=aLen+bLen;
powDa *ans=(powDa *)calloc(1,sizeof(powDa));
char *ansChar=(char *)calloc(ansLen,sizeof(char));
char *temp=(char *)malloc((bLen+1)*aLen*sizeof(char));
int i,j;
int sum;
for (i=0;i<aLen;i++){ //首先是做各位数的乘法,竖式的中间结果保存在一个矩阵中
sum=0;
for (j=0;j<bLen;j++){
sum+=(bChar[j]-'0')*(aChar[i]-'0');
temp[i*(bLen+1)+j]=sum%10+'0';
sum=sum/10;
}
temp[i*(bLen+1)+j]=sum%10+'0';
}
sum=0; //竖式求和
for (i=0;i<ansLen;i++){ //总是在各种限制的边界条件搞错,花费了很多时间,基本功不扎实啊!!
if (i-aLen+1>0){
j=i-aLen+1;
}
else j=0;
for (;((i-j)>=0)&&(j<=bLen);j++){
sum=sum+(temp[(i-j)*(bLen+1)+j]-'0');
}
ansChar[i]=sum%10+'0';
sum=sum/10;
}
free(temp);
ans->data=ansChar;
ans->len=ansLen;
return ans;
}
powDa *power_cal(powDa *a,int power){ //分治法幂递归
powDa *ans;
powDa *temp,*temp2;
if (power==1){
powDa *ans=(powDa *)calloc(1,sizeof(powDa));
char *ansChar=(char *)calloc(a->len,sizeof(char));
ans->len=a->len;
ans->data=ansChar;
memcpy(ansChar,a->data,(a->len)*sizeof(char));
return ans;
}
else if (power==2){
return charTimes(a,a);
}
if (power%2==1){
temp=power_cal(a,power/2);
temp2=power_cal(temp,2);
ans=charTimes(temp2,a);
free(temp->data);
free(temp);
free(temp2->data);
free(temp2);
}
else {
temp=power_cal(a,power/2);
ans=power_cal(temp,2);
free(temp->data);
free(temp);
}
return ans;
}
//转换函数
powDa daToPowDa(da x){
powDa ans;
char a[20]; //容纳20位的
int i,len;
sprintf(a,"%lld",x);
ans.len=strlen(a);
ans.data=(char *)calloc(ans.len,sizeof(char));
len=ans.len;
for (i=0;i<len;i++)
ans.data[i]=a[len-1-i];
if (x>=0) ans.sign=1;
else ans.sign=-1;
return ans;
}
powDa myPow(da x,da y){
powDa res;
powDa *ans;
powDa powX=daToPowDa(x);
ans=power_cal(&powX,(int)y);
res=(*ans);
if (powX.sign==1) res.sign=1;
else if (y%2==0) res.sign=1;
else res.sign=-1;
free(ans);
free(powX.data);
return res;
}
//大数加法
void powDaAdd(powDa *a,powDa b){
int aLen=a->len;
int bLen=b.len;
int ansLen=aLen>bLen?aLen:bLen;
int sum,i;
char *ans;
char *aChar=a->data;
char *bChar=b.data;
int aSign=a->sign;
int bSign=b.sign;
ansLen+=1;
ans=(char *)calloc(ansLen,sizeof(char));
sum=0;
if ((aSign*bSign)==1){
for (i=0;i<ansLen;i++){
if (i<aLen) sum+=(aChar[i]-'0');
if (i<bLen) sum+=(bChar[i]-'0');
ans[i]=sum%10+'0';
sum/=10;
}
}
else {
for (i=0;i<ansLen;i++){
if (i<aLen) sum+=(aChar[i]-'0');
if (i<bLen) sum-=(bChar[i]-'0');
ans[i]=(10+sum)%10+'0';
if (sum<0) sum=-1;
else sum=0;
}
if (sum==-1){
a->sign=-aSign;
sum=0;
for (i=0;i<ansLen;i++){
if (i<aLen) sum-=(aChar[i]-'0');
if (i<bLen) sum+=(bChar[i]-'0');
ans[i]=((10+sum)%10)+'0';
if (sum<0) sum=-1;
else sum=0;
}
}
}
free(aChar);
a->data=ans;
a->len=ansLen;
}
int solvNum(da n,da *a){ //求解不同质因子,以及质因子数目
int k;
da i,sqrtN;
k=1;
while (n!=1){
sqrtN=(da)sqrt(n);
for (i=2;i<=sqrtN;i++){
if (n%i==0){
if (a[k-1]!=i){
a[k]=i;
k++;
}
n=n/i;
break;
}
}
if ((i==(sqrtN+1))||(sqrtN==1)){ //n是个素数。如果是合数,必有小于sqrtN的因子
if (a[k-1]!=n){
a[k]=n;
k++;
}
break;
}
}
return k-1; //a[0]作为哨兵用了
}
void calNum(da N,da M){
da a[15]; //在100000000之内必然小于30个因子实际上是小于log2(100000000)
//理论上不同的因子个数应该小于10,2*3*5*7*9*11*13*17*19*23>1,000,000,000
int num;
int i,j,k,pow2;
da div,w;
powDa sum,add;
if (M==0) {printf("0\n");return;}
if (N==0) {M==1?printf("1\n"):printf("0\n");return;}
memset(a,0,sizeof(da)*15); //清空
num=solvNum(M,a); //分解质因数
//容斥原理剔除M因子的倍数的
sum=myPow(M,N);
pow2=(1<<num);
for(j=1;j<pow2;j++){
i=0;
div=1;
for (k=0;k<num;k++){
if (((j>>k)&((long long int)1))){
div*=a[k+1];
i++;
}
}
w=M/div;
add=myPow(w,N);
if (i%2==1) add.sign*=-1;
powDaAdd(&sum,add);
free(add.data);
}
printPowDa(sum);
free(sum.data);
}
int main(){
da M;
da N;
scanf("%lld %lld",&N,&M);
calNum(N,M);
return 0;
}
不过呢,后来又看了看其他大神的解法,发现了一种不用容斥原理的解法:
扩展欧拉函数!好像很牛逼的样子
不过其他人代码中让我不解的是1.没有用大数?他们是怎么解决100000000^15这种问题的?2.难道真不用考虑N=0;M=0的情况么