当前位置: 代码迷 >> 综合 >> HDU 4059 The Boss on Mars-矩阵+容斥
  详细解决方案

HDU 4059 The Boss on Mars-矩阵+容斥

热度:18   发布时间:2023-12-19 10:39:31.0

错了29遍,终成正果。。。。。

根据题意,很容易的可以想到容斥。

然后的问题就是如何求

sum(n)=1^4+2^4+3^4+....+n^4;

有三种道路:

很显然:1^4+2^4+3^4+....+n^4=(n^5)/5+(n^4)/2+(n^3)/3-n/30;

则1,用java的大数去敲这个的代码。

2,用c++敲,但是用到分数取模,求逆元。

3,用c++敲,但是不用这个公式,用矩阵去构造sum(n).

我用的是第三种。但是第三种有的缺陷,就是时间复杂度有点高。

接下来的问题就是如何优化时间复杂度。

可以预处理1~200W的sum,然后剩下的再用矩阵去构造。

具体构造方法看代码。。


#include<stdio.h>
#include<string.h>
#include<algorithm>
#include<iostream>
#include<math.h>
#include<vector>
using namespace std;
#define LL __int64
#define MOD 1000000007
#define maxn 2000007
LL yu[maxn];
struct matrix
{LL mat[7][7];matrix(){memset(mat,0,sizeof(mat));}friend matrix operator *(matrix A,matrix B){int i,j,k;matrix C;for(i=1; i<=6; i++){for(j=1; j<=6; j++){for(k=1; k<=6; k++){C.mat[i][j]=(C.mat[i][j]+(A.mat[i][k]*B.mat[k][j])%MOD)%MOD;}C.mat[i][j]=C.mat[i][j]%MOD;}}return C;}
} ONE,AA[30];
matrix powmul(matrix A,LL k)
{matrix B;for(int i=1; i<=6; i++)B.mat[i][i]=1;int l=1;while(k){if(k&1)B=B*AA[l];l++;k>>=1;}return B;
}
vector<int>vec;
void init()
{int i,j;//构造矩阵ONE.mat[1][1]=1;ONE.mat[1][2]=1;ONE.mat[1][3]=4;ONE.mat[1][4]=6;ONE.mat[1][5]=4;ONE.mat[1][6]=1;ONE.mat[2][1]=0;ONE.mat[2][2]=1;ONE.mat[2][3]=4;ONE.mat[2][4]=6;ONE.mat[2][5]=4;ONE.mat[2][6]=1;ONE.mat[3][1]=0;ONE.mat[3][2]=0;ONE.mat[3][3]=1;ONE.mat[3][4]=3;ONE.mat[3][5]=3;ONE.mat[3][6]=1;ONE.mat[4][4]=1;ONE.mat[4][5]=2;ONE.mat[4][6]=1;ONE.mat[5][4]=0;ONE.mat[5][5]=1;ONE.mat[5][6]=1;ONE.mat[6][6]=1;yu[0]=0;yu[1]=1;LL x;for(i=2;i<maxn;i++)//预处理1~200W的sum[x]{x=i%MOD;x=x*i;x=x%MOD;x=x*i;x=x%MOD;x=x*i;x=x%MOD;yu[i]=yu[i-1]+x;yu[i]=yu[i]%MOD;}AA[1]=ONE;for(i=2;i<30;i++)AA[i]=AA[i-1]*AA[i-1];
}matrix A;
LL kan(LL x)//求sum[x]
{if(x<maxn)return yu[x];matrix OP;for(int i=1;i<=6;i++)OP.mat[i][1]=1;A=powmul(ONE,x-1);OP=A*OP;return OP.mat[1][1];
}
LL pan(LL n,LL x)
{if(x>=n)return 0;LL p=1;LL ns=4;while(ns--){p=p*x;p=p%MOD;}ns=(n-1)/x;p=p*kan(ns);p=p%MOD;return p;
}
void dos(LL n)
{LL p=1;LL ans=0;int i,j,leap;LL m=n;vec.clear();for(i=2;i*i<=n;i++)//求因子的过程{if(n%i==0){vec.push_back(i);}while(n%i==0)n=n/i;}if(n!=1)vec.push_back(n);n=m;int t=1<<(vec.size());for(i=1;i<t;i++)//容斥的过程{leap=0;p=1;for(j=0;j<vec.size();j++){if(i&(1<<j)){leap++;p=p*vec[j];}}leap=leap%2;if(leap)ans+=pan(n,p);else ans-=pan(n,p);ans=(ans+MOD)%MOD;}ans=kan(n-1)-ans;ans=(ans+MOD)%MOD;printf("%I64d\n",ans);
}
int main()
{int T;LL n;init();scanf("%d",&T);while(T--){scanf("%I64d",&n);if(n==1){cout<<"0"<<endl;continue;}dos(n);}return 0;
}