NYOJ 301 递推求值【矩阵快速幂取模】

递推求值

时间限制： 1000 ms  |  内存限制： 65535 KB

难度： 4

描述

给你一个递推公式：

f(x)=a*f(x-2)+b*f(x-1)+c

并给你f(1),f(2)的值，请求出f(n)的值，由于f(n)的值可能过大，求出f(n)对1000007取模后的值。

注意:-1对3取模后等于2

输入

第一行是一个整数T,表示测试数据的组数(T<=10000)
随后每行有六个整数，分别表示f(1),f(2),a,b,c,n的值。
其中0<=f(1),f(2)<100,-100<=a,b,c<=100,1<=n<=100000000 (10^9)

输出

输出f(n)对1000007取模后的值

样例输入

2
1 1 1 1 0 5
1 1 -1 -10 -100 3

样例输出

5
999896

来源

经典题目

上传者

张云聪

原题链接：http://acm.nyist.net/JudgeOnline/bestcode.php?pid=301

分析：由于n的值比较大，所以常规方法肯定会超时。根据递推式求第n个表达式的值时，通常用矩阵乘法来做。

本题要构造两个矩阵，其中一个为矩阵A，作为初始矩阵

f2  0   0
f1  0   0
1   0   0
另一个为矩阵B
b   a   c
1   0   0
0   0   1
因为F(2)和F(1)是已知的，当n>=3时，每次都乘以矩阵B，就能推出下一个矩阵。而矩阵的第一行第一列的元素就是所求的结果。

所以利用矩阵快速幂能够快速准确地求出结果。

//B^(n-2)采用矩阵快速幂计算

参考博客：http://www.cnblogs.com/520xiuge/archive/2016/05/08/5469577.html

http://blog.csdn.net/lyhvoyage/article/details/22926265

AC代码：

#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
typedef long long LL;
const int maxn=100+5;
const int mod=1000007;
const int N=3;
struct Mat
{
    LL mat[N][N];
};
Mat mul(Mat a,Mat b)
{
    Mat ans;
    for(int i=0; i<N; i++)
    {
        for(int j=0; j<N; j++)
        {
            ans.mat[i][j]=0;
            for(int k=0; k<N; k++)
            {
                ans.mat[i][j]+=a.mat[i][k]*b.mat[k][j];
                ans.mat[i][j]%=mod;
            }
        }
    }
    return ans;
}
Mat quickPow(Mat a,int n)
{
    Mat ans=
    {
        1,0,0,
        0,1,0,
        0,0,1
    };
    while(n)
    {
        if(n&1)
            ans=mul(ans,a);
        n>>=1;
        a=mul(a,a);
    }
    return ans;
}
int main()
{
    int T;
    cin>>T;
    int a,b,c,n;
    int f[2];
    Mat tmp,ans;
    Mat e=
    {
        0,0,0,
        1,0,0,
        0,0,1
    };
    tmp.mat[2][0]=1;
    while(T--)
    {
        cin>>f[0]>>f[1]>>a>>b>>c>>n;
        if(n==1||n==2)
        {
            cout<<((f[n-1]+mod)%mod)<<endl;
            continue;
        }
        e.mat[0][0]=b;
        e.mat[0][1]=a;
        e.mat[0][2]=c;

        tmp.mat[0][0]=f[1];
        tmp.mat[1][0]=f[0];
        //

        ans=mul(quickPow(e,n-2),tmp);
        cout<<(ans.mat[0][0]+mod)%mod<<endl;
    }
    return 0;
}

NYOJ排行第一代码：

运行号：548458
运行人:张燚

#include<stdio.h>
long long m,n,a,b,c,k;
long long mol=1000007;
long long x1,x2,x3;
long long y1,y2,y3;
long long z1,z2,z3;
void f()
{
    long long h=k;
    if(h==1)
    {
        printf("%lld\n",m);
        return;
    }
    if(h==2)
    {
        printf("%lld\n",n);
        return;
    }
    long long r1=1,r2=0,r3=0;
    long long s1=0,s2=1,s3=0;
    long long t1=0,t2=0,t3=1;
    h-=2;
    while(h>=1)
    {
        int o1,o2,o3,p1,p2,p3,q1,q2,q3;
        if(h%2==1)
        {
            o1=(r1*x1+s1*x2+t1*x3)%mol;
            p1=(r1*y1+s1*y2+t1*y3)%mol;
            q1=(r1*z1+s1*z2+t1*z3)%mol;
            o2=(r2*x1+s2*x2+t2*x3)%mol;
            p2=(r2*y1+s2*y2+t2*y3)%mol;
            q2=(r2*z1+s2*z2+t2*z3)%mol;
            o3=(r3*x1+s3*x2+t3*x3)%mol;
            p3=(r3*y1+s3*y2+t3*y3)%mol;
            q3=(r3*z1+s3*z2+t3*z3)%mol;

            r1=o1,r2=o2,r3=o3;
            s1=p1,s2=p2,s3=p3;
            t1=q1,t2=q2,t3=q3;
        }
        o1=(x1*x1+y1*x2+z1*x3)%mol;
        p1=(x1*y1+y1*y2+z1*y3)%mol;
        q1=(x1*z1+y1*z2+z1*z3)%mol;
        o2=(x2*x1+y2*x2+z2*x3)%mol;
        p2=(x2*y1+y2*y2+z2*y3)%mol;
        q2=(x2*z1+y2*z2+z2*z3)%mol;
        o3=(x3*x1+y3*x2+z3*x3)%mol;
        p3=(x3*y1+y3*y2+z3*y3)%mol;
        q3=(x3*z1+y3*z2+z3*z3)%mol;

        x1=o1,x2=o2,x3=o3;
        y1=p1,y2=p2,y3=p3;
        z1=q1,z2=q2,z3=q3;
        h/=2;
    }
    n=(m*r2+n*s2+t2)%mol;
    n=(n+mol)%mol;
    printf("%lld\n",n);
}
int main()
{
    int N;
    scanf("%d",&N);
    while(N--)
    {
        scanf("%lld%lld%lld%lld%lld%lld",&m,&n,&a,&b,&c,&k);
        x1=0,y1=1,z1=0;
        x2=a,y2=b,z2=c;
        x3=0,y3=0,z3=1;
        f();
    }
}
标程代码：

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<iostream>
using namespace std;
#define CLR(arr,val) memset(arr,val,sizeof(arr))
const int MAX=100;
char str[MAX+2];
struct Matrix;
int mod=1000007;
struct Matrix
{
    Matrix() {}
    void init(int r,int c) //全零矩阵
    {
        row=r;
        col=c;
        CLR(data,0);
    }
    void init(int size)
    {
        row=size;
        col=size;
        CLR(data,0);
        for(int i=0; i!=size; i++) data[i][i]=1;
    }

    int data[MAX][MAX];
    int row,col;
    const Matrix& pow(int m);
};
Matrix ret,temp,pret;
const Matrix& operator*(const Matrix& m1,const Matrix &m2)
{
    ret.init(m1.row,m2.col);
    for(int i=0; i!=m1.row; i++)
        for(int j=0; j!=m1.col; j++)
            if(m1.data[i][j]) for(int k=0; k!=m2.col; k++)
                    if(m2.data[j][k]) ret.data[i][k]=(ret.data[i][k]+((long long)m1.data[i][j]*m2.data[j][k]))%mod;
    return ret;
}
const Matrix& Matrix::pow(int m)
{
    temp=*this;
    pret.init(row);
    while(m)
    {
        if(m&1)
            pret=temp*pret;
        temp=temp*temp;
        m>>=1;
    }
    return pret;
}
Matrix m1,m2;
int main()
{
    int t,k;
    m1.init(3,3);
    m1.data[1][0]=1;
    m1.data[2][2]=1;
    m2.init(1,3);
    m2.data[0][2]=1;
    scanf("%d",&t);
    while(t--)
    {
        scanf("%d%d%d%d%d%d",&m2.data[0][0],&m2.data[0][1],&m1.data[0][1],&m1.data[1][1],&m1.data[2][1],&k);
        if(k==1) printf("%d\n",m2.data[0][0]);
        else printf("%d\n",((m2*m1.pow(k-2)).data[0][1]+mod)%mod);
    }
}


AC代码3：

#include <stdio.h>
#include <string.h>
#define N 3
#define mod 1000007
#define ll long long
struct Matrix{
    ll mat[N][N];//乘积时(可能接近mod^2)超出int,
};
struct Matrix mul(struct Matrix a,struct Matrix b)
{
    int i,j,k;
    struct Matrix res;
    for(i=0;i<N;i++)
    {
        for(j=0;j<N;j++)
        {
            res.mat[i][j] = 0;
            for(k=0;k<N;k++)
            {
                res.mat[i][j] += a.mat[i][k]*b.mat[k][j];
                res.mat[i][j] %= mod;//必须每次取余
            }
        }
    }
    return res;
}
struct Matrix mul_matrix(struct Matrix b,int n)
{
    struct Matrix res = {
        1,0,0,
        0,1,0,
        0,0,1
    };//单位阵
    while(n)
    {
        if(n&1)
            res = mul(res,b);
        n >>= 1;
        b = mul(b,b);
    }
    return res;
}
int main()
{
    int T,f[2],a,b,c,n;
    struct Matrix tmp,res;
    struct Matrix tmp_matrix = {//中间矩阵
        0,0,0,//每组测试将此三处的值更新为b,a,c
        1,0,0,
        0,0,1
    };
    memset(tmp.mat,0,sizeof(tmp.mat));//tmp值只有(0,0),(1,0)需要填入f(2),f(1)的值,其余不变
    tmp.mat[2][0] = 1;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d%d%d%d%d",&f[0],&f[1],&a,&b,&c,&n);
        if(n == 1 || n == 2)
            printf("%d\n",(f[n-1]+mod)%mod);//注意对负值的处理
        else{
            tmp_matrix.mat[0][0] = b;
            tmp_matrix.mat[0][1] = a;
            tmp_matrix.mat[0][2] = c;
            tmp.mat[0][0] = f[1];
            tmp.mat[1][0] = f[0];
            //实际上最终我们只需进行矩阵和向量(f(2),f(1),1)'的运算,但是因为在进行矩阵快速幂时我们已经定义了矩阵乘积
            //所以不妨借用,最后取(0,0)的值即可
            res = mul(mul_matrix(tmp_matrix,n-2),tmp);
            printf("%lld\n",(res.mat[0][0]+mod)%mod);
        }
    }
    return 0;
}        
尊重原创，转载请注明出处：http://blog.csdn.net/hurmishine