HDU - 6395 Sequence（分块+矩阵快速幂）

题目链接

题意：

$
\begin{cases}
& F_{1} = A \\
& F_{2} = B \\
& F_{n} = C·F_{n-2} + D·F_{n-1}+\lfloor \frac{P}{n}\rfloor
\end{cases}
$

给定递推式，求$F_{n}\mod{1e9+7}$

思路：

一来就会想到用矩阵快速幂，但是其中有一个下取整就肯定不是裸的矩阵快速幂；

有向下取整我们就可以用分块来做，复杂度为$\sqrt{P}log_{2}P$

  1 /*
  2 *  Author: windystreet
  3 *  Date  : 2018-08-14 09:46:10
  4 *  Motto : Think twice, code once.
  5 */
  6 #include<bits/stdc++.h>
  7 
  8 using namespace std;
  9 
 10 #define X first
 11 #define Y second
 12 #define eps  1e-5
 13 #define gcd __gcd
 14 #define pb push_back
 15 #define PI acos(-1.0)
 16 #define lowbit(x) (x)&(-x)
 17 #define bug printf("!!!!!\n");
 18 #define mem(x,y) memset(x,y,sizeof(x))
 19 
 20 typedef long long LL;
 21 typedef long double LD;
 22 typedef pair<int,int> pii;
 23 typedef unsigned long long uLL;
 24 
 25 const int maxn = 1e5+7;
 26 const int INF  = 1<<30;
 27 const int mod  = 1e9+7;
 28 
 29 struct Matrix
 30 {
 31     int n,m;
 32     LL ma[5][5];
 33     Matrix (int x,int y):n(x),m(y){clear();}
 34     void set(int n_,int m_){n = n_,m = m_;}
 35     LL *operator[](int x){return ma[x];}
 36     Matrix operator*(Matrix x){
 37         Matrix res(n,x.m);
 38         for(int i=1;i<=n;i++)
 39             for(int j=1;j<=x.m;j++)
 40                 for(int k=1;k<=m;k++)
 41                     (res[i][j]+=ma[i][k]*x[k][j]%mod+mod)%=mod;
 42         return res;
 43     }
 44     Matrix operator ^(int y){
 45         Matrix x(n,m);
 46         for(int i=1;i<=n;i++)
 47             for(int j=1;j<=m;j++)
 48                 x[i][j]=ma[i][j];
 49         Matrix res(x.n,x.n);
 50         for(int i=1;i<=x.n;i++)
 51             res[i][i]=1;
 52         for(;y;y>>=1,x=x*x)
 53             if(y&1)res=res*x;
 54         return res;
 55     }
 56     void print(){
 57         for(int i=1;i<=n;i++)
 58             for(int j=1;j<=m;j++)
 59                 printf("%lld%c",ma[i][j]," \n"[j==m]);
 60     }
 61     void clear(){mem(ma,0);}
 62 };
 63 int a ,b,c,d ,p,n;
 64 pii f(int x,int y){
 65     if(!x)return make_pair(y,n);
 66     int l = max(y,(p+x+1)/(x+1)),r = min(n,p/x);
 67     return make_pair(l,r);
 68 }
 69 
 70 void solve(){
 71     scanf("%d%d%d%d%d%d",&a,&b,&c,&d,&p,&n);
 72     if(n==1)printf("%d\n",a);
 73     if(n==2)printf("%d\n",b);
 74     else{
 75         Matrix x(1,3),y(3,3);
 76         x[1][1] = b,x[1][2]=a,x[1][3]=1;
 77         for(int i=3;i<=n;){
 78             pii seg = f(p/i,i);
 79             y[1][1]=d,y[1][2] = 1,y[1][3] = 0;
 80             y[2][1]=c,y[2][2] = 0,y[2][3] = 0;
 81             y[3][1]=p/i,y[3][2]=0,y[3][3] = 1;
 82             x = x*(y^(seg.Y - seg.X  + 1));
 83             i = seg.Y+1;
 84         }
 85         printf("%lld\n",x[1][1] );
 86     }
 87     
 88     return;
 89 }
 90 
 91 int main()
 92 {
 93 //    freopen("F:\\in.txt","r",stdin);
 94 //    freopen("out.txt","w",stdout);
 95 //    ios::sync_with_stdio(false);
 96     int t = 1;
 97     scanf("%d",&t);
 98     while(t--){
 99     //    printf("Case %d: ",cas++);
100         solve();
101     }
102     return 0;
103 }

 

posted @ 2018-08-14 11:19  windystreet  阅读(385)  评论(0)  编辑  收藏  举报