HDU5950-Recursive sequence(矩阵快速幂)

题意：

给出n头母牛，第一头报a，第二头报b，第i头报f[i-2]*2+f[i-1]+i^4,问第n头母牛报数多少。

思路：

就是推公式啊，可惜没经验，白搭。

自己推得时候一直不知道i^4怎么消去，原来是可以加上i^3,i^2,i，再算的。总之相乘的矩阵必须都是常数，未知数一定要放到最后乘。

网上找的：

1 0 0 0 0 0 0        1               1

1 1 0 0 0 0 0        i                i+1

1 2 1 0 0 0 0       i2              (i+1)2

1 3 3 1 0 0 0  *   i3     =       (i+1)3

1 4 6 4 1 0 0       i4               (i+1)4

0 0 0 0 0 0 1       f[i-2]           f[i-1]

0 0 0 0 1 2 1      f[i-1]             f[i]

#include <bits/stdc++.h>
using namespace std;
const long long mod =  2147493647;
struct Matrix
{
    long long a[8][8];
};
Matrix s,B;

Matrix matrixmul(Matrix a,Matrix b)
{
    Matrix c;
    for(int i=1; i<8; i++)
    {
        for(int j=1; j<8; j++)
        {
            c.a[i][j]=0;
            for(int k=1; k<8; k++)
                c.a[i][j]+=(a.a[i][k]*b.a[k][j])%mod;
            c.a[i][j]%=mod;
        }

    }
    return c;
}
Matrix mul(Matrix s,int k)
{
    Matrix ans;
    for(int i=1; i<8; i++)
    {
        for(int j=1; j<8; j++)
            ans.a[i][j]=(i==j)?1:0;
    }
    while(k)
    {
        if(k&1)
            ans=matrixmul(ans,s);
        k>>=1;
        s=matrixmul(s,s);
    }
    return ans;
}
int main()
{
    int n,t,a,b;
    for(cin>>t; t--;)
    {
        cin>>n>>a>>b;
        if(n==1)
        {
            cout<<a%mod<<endl;
            continue;
        }
        if(n==2)
        {
            cout<<b%mod<<endl;
            continue;
        }
        if(n==3)
        {
            cout<<(81+2*a%mod+b%mod)%mod<<endl;
            continue;
        }
        n-=2;
        for(int i=1; i<=7; i++)
            for(int j=1; j<=7; j++)
            {
                s.a[i][j]=0;
                B.a[i][j]=0;
            }
        for(int i=1; i<=5; i++)
        {
            s.a[i][1]=1;
        }
        for(int i=2; i<=5; i++)
        {
            s.a[i][2]=i-1;
        }
        s.a[3][3]=1;
        s.a[4][3]=3;
        s.a[5][3]=6;
        s.a[4][4]=1;
        s.a[5][4]=4;
        s.a[5][5]=1;
        s.a[6][5]=1;
        s.a[6][6]=1;
        s.a[7][6]=1;
        s.a[6][7]=2;
        B.a[1][1]=1;
        B.a[2][1]=3;
        B.a[3][1]=9;
        B.a[4][1]=27;
        B.a[5][1]=81;
        B.a[6][1]=b;
        B.a[7][1]=a;
        s=mul(s,n);
        s=matrixmul(s,B);
        cout<<s.a[6][1]%mod<<endl;

    }
    return 0;
}