HDU 6395 Sequence 矩陣冪,下取整,分段
阿新 • • 發佈:2018-11-02
題意:F[1]=A,F[2]=B, F[i]=C*F[i-2]+D*F[i-1]+[p/i]
1<=P,n<=1e9. 求出第n項的值.
[p/i]隨i不同而變化, 如果沒有這一項或者這一項是一個常數的話,顯然就可以用矩陣冪來推.
[p/i]下取整的值最多隻有2*sqrt(p)項. (i<sqrt(p)顯然結果最多sqrt(p)項, i>=sqrt(p) 結果最大為sqrt(p),則[sqrt(p),...1]).
初始前一段的兩個值為F1,F2,然後用矩陣冪快速求出下一段最後兩個值即可.O(sqrt(p) * log n).
#include <bits/stdc++.h> using namespace std; typedef long long ll; const int mod=1e9+7; int T,A,B,C,D,P,n; struct mat{ ll b[3][3]; mat(){memset(b,0,sizeof(b));} mat operator *(mat &tmp){ mat ans; for(int i=0;i<3;i++) for(int j=0;j<3;j++) for(int k=0;k<3;k++) ans.b[i][j]=(ans.b[i][j]+(b[i][k]*tmp.b[k][j])%mod)%mod; return ans; } }unit,num; ll calc(ll i){ ll val=P/i; ll l=i,r=n; while(l<r){ int mid=r-(r-l)/2; if(P/mid==val) l=mid; else if(P/mid<val) r=mid-1; else l=mid+1; } return l; } mat powmod(mat a,ll n){ mat ans=unit; while(n){ if(n&1) ans=ans*a; n>>=1; a=a*a; } return ans; } int main() { memset(unit.b,0,sizeof(unit.b)); unit.b[0][0]=unit.b[1][1]=unit.b[2][2]=1; cin>>T; while(T--){ cin>>A>>B>>C>>D>>P>>n; ll f1=B,f2=A; if(n<=2){ printf("%d\n",n==1?f2:f1); continue; } // b={{D,1,0},{C,0,0},{0,0,1}};//{fi,fi-1,1}*B^n memset(num.b,0,sizeof(num.b)); num.b[0][0]=D,num.b[0][1]=1; num.b[1][0]=C,num.b[2][2]=1; for(ll i=3;i<=n;){ ll j=calc(i),val=P/i;// [i,i+1...j]結果都相同.都為p/i. // cout<<i<<' '<<j<<' '<<val<<'\n'; num.b[2][0]=val; mat aa=powmod(num,j-i+1); ll F1=((f1*aa.b[0][0])%mod+(f2*aa.b[1][0])%mod+aa.b[2][0])%mod; ll F2=((f1*aa.b[0][1])%mod+(f2*aa.b[1][1])%mod+aa.b[2][1])%mod; i=j+1; f1=F1,f2=F2; } printf("%lld\n",f1); } return 0; }