洛谷P1397 [NOI2013] 矩阵游戏题解

题目链接：P1397 [NOI2013] 矩阵游戏

题意：

cxy 是个喜欢矩阵的小朋友，有一天她想用电脑生成一个巨大的 \(n\) 行 \(m\) 列的矩阵 (你不用担心她如何存储)。她生成的这个矩阵满足一个神奇的性质：

若用 \(F_{i, j}\) 来表示矩阵中第 \(i\) 行第 \(j\) 列的元素，则 \(F_{i, j}\) 满足下面的递推式：

\[ \begin{cases} F_{1, 1} & = & 1 \\ F_{i, j} & = & a \times F_{i, j-1} + b & j \neq 1 \\ F_{i, 1} & = & c \times F_{i-1, m} + d & i \neq 1 \end{cases} \] 递推式中 \(a, b, c, d\) 都是给定的常数。

现在 cxy 想知道 \(F_{n, m}\) 的值是多少，请你帮助她。由于最终结果可能很大，你只需要输出 \(F_{n, m}\) 除以 \(10^9 + 7\) 的余数。

输入格式：

共一行，包含六个正整数 \(n, m, a, b, c, d~(1\le n,m\le 10^{1000000},~1\le a,b,c,d\le 10^9)\) ，意义如题所述。

输出格式：

输出一行一个整数，表示 \(F_{n, m}\) 除以 \(10^9 + 7\) 的余数。

不容易，推了1个小时，总算推出来这题了。

综合性极强的题目，涉及：求解数列通项、扩展欧拉定理（费马小定理）

首先考虑行的递推式，即寻找 \(F_{i,m}\) 与 \(F_{i,1}\) 的关系

定义 \[ f^0(x) = x \\f^1(x) = ax+b \\f^{n}(x) = a\times f^{n-1}(x) + b \] 注意到这符合数列递推式的性质，并且是最简单的一阶线性递推数列。于是我们求个通项。

其特征方程为 \(x = ax + b \Rightarrow (1-a)x = b\) （不会特征方程法可以戳上面那个链接哦 >w<）

当 \(a=1\) 时，不难发现 \(f^n(x) = x + nb\)
当 \(a\ne 1\) 时，\(x_0 = \frac{b}{1-a}\) 。根据公式不难发现 \[ f^n(x) = a^nx - b\dfrac{a^n-1}{1-a} \]

这里有一丢丢边界的细节，仔细想一下，可以得到 \[ F_{i,m} = f^{m-1} (F_{i,1} ) \] 记 \(A = a^n,~B=- b\dfrac{a^n-1}{1-a}\) ，则 \[ F_{i,m} = A\times F_{i,1} + B \] 然后考虑 \(F_{i,1}\) 和 \(F_{i-1,1}\) 的关系。

首先有 \[ F_{i,1} = c\times F_{i-1,m}+d = cA \times F_{i-1,1} + cB + d \] 可以发现这又是一个一阶线性递推数列的形式（第二维的 \(1\) 是不变的）

记 \(C=cA,~D=cB+d\) 。

定义 \[ g^0 = F_{1,1}=1 \\g^{n} = C\times g^{n-1} + D \] 此时就是和上面一样的了。

我们可以实现一个函数 fun(a,b,k) ，具体看代码&注释即可。

最后考虑 \(n,m\le 10^{1000000}\) 这件事。~~据说可以用十进制快速幂？~~

根据扩展欧拉定理（~~其实这题用费马小定理就可以了，因为这题模数是质数~~） \[ \begin{aligned} a^b \equiv \begin{cases} a^b,&b< \varphi(m)\\\\ a^{b \ \bmod \ \varphi(m)+\varphi(m)},&b\ge \varphi(m) \end{cases} \pmod{m} \end{aligned} \] 然后我们只要边读边模就好了（实际上也可以实现一个函数来搞这个）

时间复杂度 \(O(\lg n + \log p)\)

代码：

#include <iostream>
#include <string>
#include <vector>
#include <algorithm>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cstdarg>
#include <cmath>
#include <iomanip>
#include <random>
using namespace std;
#define int long long
#define INF 0x3f3f3f3f3f3f3f3f
#define N ((int)(1e6+15))
const int mod = 1e9+7, phi_mod = mod-1;


int qpow(int a,int b)
{
    int ans=1,base=a%mod;
    if(b == -1) b = mod-2;
    for(; b; b>>=1)
    {
        if(b&1) ans = ans * base % mod;
        base = base * base % mod;
    } return ans;
}
int getnum(char *s,int p)
{
    int res=0;
    for(char *c = s; *c; ++c)
        res = (res * 10 + (*c - 48)) % p;
    return res;
}
/*
    fun(a,b,k)

    a=1 -> return {1, bk}
    a≠1 -> return {a^{k-1}, -b * \frac { a^{k-1} -1 }{ 1-a } }
    
    ps: return is &A,&B 
*/
void fun(int a,int b,char *k,int &A,int &B)
{
    if(a == 1) { A=1; B = b * (getnum(k,mod) + mod - 1) % mod; return;}
    int t = b * qpow(mod + 1 - a, -1) % mod; // b * (1-a)^-1
    A = qpow(a, (getnum(k, phi_mod) - 1) % phi_mod + phi_mod); // a^{k-1}
    B = (mod + 1 - A) * t % mod; // balabala >w< ~
}
char n[N],m[N];
int a,b,c,d,u,v,A,B;
signed main()
{
    ios::sync_with_stdio(0);
    cin.tie(0);cout.tie(0);
    // freopen("check.in","r",stdin);
    // freopen("check.out","w",stdout);
    cin >> n >> m >> a >> b >> c >> d;
    fun(a,b,m,A,B);
    fun(c * A % mod, (c * B + d) % mod,n,u,v);
    cout << (A * (u + v) + B) % mod << '\n';
    return 0;
}