HDU 6395 Sequence 杜教板子题
题目意思非常明确,就是叫你求第n项,据我们学校一个大佬说他推出了矩阵,但是我是菜鸡,那么肯定是用简单的方法水过啦!我们先p^(1/2)的复杂度处理出i=[i,p]范围内的所有种类的(int)(p/i),然后我们就可以知道种可能的除数的范围,就是分成几块
这里我不太会表达,看代码比较好
1 /** 2 求n/i的所有结果 3 **/ 4 #include<stdio.h> 5 int main( ){ 6 int n; 7 scanf("%d", &n); 8 long long nex; 9 for(long long i=1; i<=n; i=nex+1){ 10 nex=n/(n/i); 11 printf("%I64d, %I64d\n", n/i, nex); 12 /** 13 n/i保存到数组里从小到大排序,你们就知道块是按这个分的了 14 **/ 15 /** 16 看一下输出之类的应该能发现就是类似于 17 int n; 18 vector<int> V; 19 scanf("%d", &n); 20 for(int i=1; i<=n; ++i){ 21 V.push_back((n/i)); 22 } 23 ......数组去重 24 会发现这就是n/i的所有结果 25 26 **/ 27 } 28 }
对这些分出来的块我们判断一下,如果两块之间距离很小,那么就没必要用杜教板子推第n项;直接暴力;
其他的就加入该块的前几个元素去推;
具体看恶心的代码
1 #include<bits/stdc++.h> 2 using namespace std; 3 #define rep(i,a,n) for (int i=a;i<n;i++) 4 #define pb push_back 5 #define SZ(x) ((int)(x).size()) 6 typedef vector<int> VI; 7 typedef long long ll; 8 typedef pair<int,int> PII; 9 const ll mod=1000000007; 10 ll powmod(ll a,ll b) {ll res=1;a%=mod; assert(b>=0); for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;} 11 int _, n; 12 namespace linear_seq { 13 const int N=10010; 14 ll res[N],base[N],_c[N],_md[N]; 15 vector<int> Md; 16 void mul(ll *a,ll *b,int k) { 17 rep(i,0,k+k) _c[i]=0; 18 rep(i,0,k) if (a[i]) rep(j,0,k) _c[i+j]=(_c[i+j]+a[i]*b[j])%mod; 19 for (int i=k+k-1;i>=k;i--) if (_c[i]) 20 rep(j,0,SZ(Md)) _c[i-k+Md[j]]=(_c[i-k+Md[j]]-_c[i]*_md[Md[j]])%mod; 21 rep(i,0,k) a[i]=_c[i]; 22 } 23 int solve(ll n,VI a,VI b) { // a 系数 b 初值 b[n+1]=a[0]*b[n]+... 24 // printf("%d\n",SZ(b)); 25 ll ans=0,pnt=0; 26 int k=SZ(a); 27 assert(SZ(a)==SZ(b)); 28 rep(i,0,k) _md[k-1-i]=-a[i];_md[k]=1; 29 Md.clear(); 30 rep(i,0,k) if (_md[i]!=0) Md.push_back(i); 31 rep(i,0,k) res[i]=base[i]=0; 32 res[0]=1; 33 while ((1ll<<pnt)<=n) pnt++; 34 for (int p=pnt;p>=0;p--) { 35 mul(res,res,k); 36 if ((n>>p)&1) { 37 for (int i=k-1;i>=0;i--) res[i+1]=res[i];res[0]=0; 38 rep(j,0,SZ(Md)) res[Md[j]]=(res[Md[j]]-res[k]*_md[Md[j]])%mod; 39 } 40 } 41 rep(i,0,k) ans=(ans+res[i]*b[i])%mod; 42 if (ans<0) ans+=mod; 43 return ans; 44 } 45 VI BM(VI s) { 46 VI C(1,1),B(1,1); 47 int L=0,m=1,b=1; 48 rep(n,0,SZ(s)) { 49 ll d=0; 50 rep(i,0,L+1) d=(d+(ll)C[i]*s[n-i])%mod; 51 if (d==0) ++m; 52 else if (2*L<=n) { 53 VI T=C; 54 ll c=mod-d*powmod(b,mod-2)%mod; 55 while (SZ(C)<SZ(B)+m) C.pb(0); 56 rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod; 57 L=n+1-L; B=T; b=d; m=1; 58 } else { 59 ll c=mod-d*powmod(b,mod-2)%mod; 60 while (SZ(C)<SZ(B)+m) C.pb(0); 61 rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod; 62 ++m; 63 } 64 } 65 return C; 66 } 67 int gao(VI a,ll n) { 68 VI c=BM(a); 69 c.erase(c.begin()); 70 rep(i,0,SZ(c)) c[i]=(mod-c[i])%mod; 71 return solve(n,c,VI(a.begin(),a.begin()+SZ(c))); 72 } 73 }; 74 ll sav[100007], ans[5]; 75 int main(){ 76 int T; 77 register int i, j; 78 ll c, d, p, num, nex, top, a, b; 79 scanf("%d", &T); 80 while(T--){ 81 vector<int> v; 82 scanf("%I64d%I64d%I64d%I64d%I64d%I64d", &ans[1], &ans[2], &c, &d, &p, &num); 83 a=ans[1]; 84 b=ans[2]; 85 top=0; 86 if(num<51){ 87 if(num<=2){ 88 printf("%I64d\n", ans[num]); 89 continue; 90 } 91 for(i=3; i<=num; ++i){ 92 ans[3]=(c*ans[1]%mod+d*ans[2]%mod+(p/i))%mod; 93 ans[1]=ans[2]; 94 ans[2]=ans[3]; 95 } 96 printf("%I64d\n", ans[2]); 97 continue; 98 } 99 for(i=1; i<=p; i=nex+1){ 100 nex=p/(p/i); 101 sav[top++]=(p/i); 102 } 103 for(i=0; i<=top/2; ++i){ 104 swap(sav[i], sav[top-i-1]); 105 } 106 sav[top]=num; 107 ans[3]=(c*ans[1]%mod+d*ans[2]%mod+(p/3))%mod;///这里处理是因为,在区间[sav[i], sav[i+1]]内sav[i]和[sav[i]+1, sav[i+1]]对p的除数是不同的;所以说第一次处理一下sav[i],就可以让区间变成(sav[i], sav[i+1]]半开半闭区间 108 ans[1]=ans[2]; 109 ans[2]=ans[3]; 110 for(i=2; i<top; ++i){ 111 if(sav[i+1]>=num){ 112 if(num-sav[i]<=32){ 113 for(j=sav[i]+1; j<=num; ++j){ 114 ans[3]=(c*ans[1]%mod+d*ans[2]%mod+(p/j))%mod; 115 ans[1]=ans[2]; 116 ans[2]=ans[3]; 117 } 118 printf("%I64d\n", ans[2]); 119 break; 120 }else{ 121 v.clear(); 122 for(j=sav[i]+1; j<=sav[i]+11; ++j){ 123 ans[3]=(c*ans[1]%mod+d*ans[2]%mod+(p/j))%mod; 124 v.push_back((int)ans[3]); 125 ans[1]=ans[2]; 126 ans[2]=ans[3]; 127 } 128 n=num-sav[i]; 129 printf("%d\n", linear_seq::gao(v, n-1)); 130 break; 131 } 132 }else if(sav[i+1]-sav[i]<=32){ 133 for(j=sav[i]+1; j<=sav[i+1]; ++j){ 134 ans[3]=(c*ans[1]%mod+d*ans[2]%mod+(p/j))%mod; 135 ans[1]=ans[2]; 136 ans[2]=ans[3]; 137 } 138 }else{ 139 v.clear(); 140 for(j=sav[i]+1; j<=sav[i]+11; ++j){ 141 ans[3]=(c*ans[1]%mod+d*ans[2]%mod+(p/j))%mod; 142 v.push_back((int)ans[3]); 143 ans[1]=ans[2]; 144 ans[2]=ans[3]; 145 } 146 n=sav[i+1]-sav[i]; 147 ans[1]=linear_seq::gao(v, n-2);///这里还不是结束,也就是要至少保留连续的两项 148 ans[2]=linear_seq::gao(v, n-1); 149 } 150 } 151 } 152 }