[BZOJ 3240][NOI 2013]矩阵游戏(数学+乘法逆元)

题目链接

http://www.lydsy.com/JudgeOnline/problem.php?id=3240

思路

考虑在同一行$i$，从最左边的列$1$到最右边的列$m$递推：

$$f_{i,m}=a^{m-1}f_{i,1}+b\frac{1-a^{m-1}}{1-a}(a\neq 1)$$

$$f_{i,m}=1+(m-1)b(a=1)$$

然后考虑在同一列$1$，从行$j$到行$j+1$的递推：

$$f_{j+1,1}=a^{m-1}cf_{i,1}+bc\frac{1-a^{m-1}}{1-a}+d(a\neq 1)$$
$$f_{j+1,1}=c+(m-1)bc+d(a=1)$$
可以发现其实这两个式子和上面的两个式子结构是完全一样的，那么也可以很容易推出 $f_{n+1,1}和f_{1,1}$的关系，这样就能快速求出 $f_{n+1,1}$了，而 $f_{n,m}=\frac{(f_{n+1,1}-d)}c$，由于是在模意义下，因此要用逆元做除法。

然而还有一个问题需要注意，因为$n,m$太大了，直接上快速幂肯定是做不了的，有两种做法：1、十进制拆位快速幂，即预处理$a^1,a^{10},a^{100}...$，然后按位做类似二进制上的快速幂。2、费马小定理，由于$a^{p-1}\equiv1(\mod p)$，因此$a^n=a^{n\mod p-1}$，我们用字符串读入$n和m$，预处理出$n_1=n\mod p$，$m_1=m\mod p$，$n_2=n\mod p-1$，$m_2=m\mod p-1$即可。
个人认为方法2更方便，因此推荐使用方法2

实际上这个题也可以构造矩阵，用十进制拆位的矩阵快速幂来做，这样做还是有点卡，不过可以玩玩卡时过掉此题，归根结底，矩阵的做法还是麻烦些。

代码

#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <algorithm>

#define MAXN 10
#define MOD 1000000007

using namespace std;

typedef long long int LL;

char sn[1100000],sm[1100000];

LL fastPow(LL base,LL pow)
{
    LL ans=1;
    while(pow)
    {
        if(pow&1) ans=ans*base%MOD;
        base=base*base%MOD;
        pow>>=1;
    }
    return ans;
}

inline LL extGCD(LL a,LL b,LL &x,LL &y)
{
    if(!b)
    {
        x=1;
        y=0;
        return a;
    }
    LL gcd=extGCD(b,a%b,x,y);
    LL t=x;
    x=y;
    y=t-(a/b)*y;
    return gcd;
}

inline LL rev(LL t) //求t的逆元
{
    LL x,y;
    extGCD(t,MOD,x,y);
    x=(x%MOD+MOD)%MOD;
    return x;
}

LL n1,m1,n2,m2,a,b,c,d,tmp; //n1=n%MOD,n2=n%(MOD-1)

int main()
{
    scanf("%s%s",sn,sm);
    scanf("%lld%lld%lld%lld",&a,&b,&c,&d);
    int lenn=strlen(sn),lenm=strlen(sm);
    for(int i=0;i<lenn;i++)
        n1=(n1*10+sn[i]-'0')%MOD;
    for(int i=0;i<lenm;i++)
        m1=(m1*10+sm[i]-'0')%MOD;
    for(int i=0;i<lenn;i++)
        n2=(n2*10+sn[i]-'0')%(MOD-1);
    for(int i=0;i<lenm;i++)
        m2=(m2*10+sm[i]-'0')%(MOD-1);
    if(a==1)
    {
        b=((m1-1)*b%MOD*c+d)%MOD;
        a=c;
    }
    else
    {
        b=(b*c%MOD*(fastPow(a,m2-1)-1)%MOD*rev(a-1)%MOD+d)%MOD;
        a=fastPow(a,m2-1)*c%MOD;
    }
    if(a==1)
        tmp=(1+n1*b%MOD)%MOD;
    else
        tmp=(fastPow(a,n2)+(fastPow(a,n2)-1)*b%MOD*rev(a-1)%MOD)%MOD;
    printf("%lld\n",((tmp-d)*rev(c)%MOD+MOD)%MOD);
    return 0;
}