2018 Multi-University Training Contest 7 Sequence【矩阵快速幂+分块】

Sequence

Time Limit: 4000/2000 MS (Java/Others) Memory Limit: 262144/262144 K (Java/Others)
Total Submission(s): 466 Accepted Submission(s): 146

Problem Description

Let us define a sequence as below

⎧⎩⎨⎪⎪⎪⎪⎪⎪F1F2Fn===ABC⋅Fn−2+D⋅Fn−1+⌊Pn⌋

Your job is simple, for each task, you should output Fn module 109+7 .

Input

The first line has only one integer T , indicates the number of tasks.

Then, for the next T lines, each line consists of 6 integers, A , B , C , D , P , n .

1≤T≤200≤A,B,C,D≤1091≤P,n≤109

Sample Input

2 3 3 2 1 3 5 3 2 2 2 1 4

Sample Output

36 24

Source

2018 Multi-University Training Contest 7

思路：

这一题不能直接套模板，因为p/n是一个变量，不容易从上一个状态转化过来。

我们可以从p/n下手进行优化。

观察当p=25时，n从1到25，p/n的值为：

25 12 8 6 5 4 3 3 2 2 2 2 1 1 1 1 1 1 1 1 1 1 1 1 1

可以发现前sqrt(p)的部分变化较大，后面的部分是一段一段的区间。

并且可以计算出每一段区间的边界。

当p/n=k时，这一段区间为( p/(k+1),p/k ]。

容易证得后面的部分最多有sqrt(p)个区间（0的区间除外）。

这样就可以对前面sqrt(p)的部分暴力求出，后面的部分分块进行矩阵快速幂。当前块的第一个值由上一块的最后两个值转移过来。这样就可以将p/n看成一个常量来计算。

转移矩阵：

$\bigl(\begin{smallmatrix} d & c & 1\\ 1 & 0 & 0\\ 0 & 0 & 1 \end{smallmatrix}\bigr)*\bigl(\begin{smallmatrix} f(n-1)\\ f(n-2)\\ p/n \end{smallmatrix}\bigr)\Rightarrow \bigl(\begin{smallmatrix} f(n)\\ f(n-1)\\ p/n \end{smallmatrix}\bigr)$

注意：

n<=2时直接输出。

p<4时，由于sqrt(p)<2，导致没办法将b录进去，容易出现问题。

n>p时，p/n==0，直接矩阵快速幂。

代码：

#include<bits/stdc++.h>
using namespace std;
#define mod 1000000007
#define ll long long
#define MAXN 200005
typedef long long LL;
ll f[MAXN];
struct state
{
    LL a[3][3];
};
state multi(state a,state b)
{
    state c={0};
    for(int i=0;i<3;i++)
    for(int j=0;j<3;j++)
    for(int k=0;k<3;k++)
    c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j])%mod;
    return c;
}
state quick_pow(state a,LL n)
{
    state b={0};
    for(int i=0;i<3;i++)
        b.a[i][i]=1;
    while(n)
    {
        if(n&1) b=multi(a,b);
        a=multi(a,a);
        n>>=1;
    }
    return b;
}
int main()
{

	LL n,a,b,c,d,p;
	int t;
	scanf("%d",&t);
	while(t--)
    {
        scanf("%lld%lld%lld%lld%lld%lld",&a,&b,&c,&d,&p,&n);
        int mid=sqrt(p);
        f[1]=a;
        f[2]=b;
        ll pre1,pre2,ans;
        for(int i=3;i<=mid && i<=n;i++)
            f[i]=(c*f[i-2]+d*f[i-1]+p/i)%mod;
        if(n<=mid || n<=2)
        {
            printf("%lld\n",f[n]);
            continue;
        }
        if(p<4)
        {
            f[3]=(c*f[3-2]+d*f[3-1]+p/3)%mod;
            if(n==3) printf("%lld\n",f[3]);
            else
            {
                n-=3;
                pre2=f[2];
                pre1=f[3];
                state a={0};
                state b={0};
                a.a[0][0]=d; a.a[0][1]=c; a.a[0][2]=1;
                a.a[1][0]=1;
                a.a[2][2]=1;
                b.a[0][0]=pre1;
                b.a[1][0]=pre2;
                b.a[2][0]=0;
                ans=multi(quick_pow(a,n),b).a[0][0];
                printf("%lld\n",ans);
            }
            continue;
        }
        n-=mid;
        pre1=f[mid];
        pre2=f[mid-1];
        int L,R;
        state a={0};
        state b={0};
        a.a[0][0]=d; a.a[0][1]=c; a.a[0][2]=1;
        a.a[1][0]=1;
        a.a[2][2]=1;
        for(int i=mid-1;i>=1;i--)
        {
            L=p/(i+1)+1;
            R=p/i;
            if(R-L+1==0) continue;
            b.a[0][0]=pre1;
            b.a[1][0]=pre2;
            b.a[2][0]=i;
            if(n<=R-L+1)
            {
                ans=multi(quick_pow(a,n),b).a[0][0];
                n=0;
                break;
            }
            if(R-L+1==1)
            {
                pre2=pre1;
                pre1=multi(quick_pow(a,R-L+1),b).a[0][0];
            }
            else
            {
                state now=multi(quick_pow(a,R-L+1),b);
                pre2=now.a[1][0];
                pre1=now.a[0][0];
            }
            n-=R-L+1;
        }
        if(n)
        {
            b.a[0][0]=pre1;
            b.a[1][0]=pre2;
            b.a[2][0]=0;
            ans=multi(quick_pow(a,n),b).a[0][0];
        }
        printf("%lld\n",ans);
    }
	return 0;
}