J.u's的影响力(矩阵快速幂，费马小定理)

（牛客第一场）J.u's的影响力(矩阵快速幂，费马小定理)

第一天的影响力为 $x$ ，第二天的影响力为 $y$ ，从第三天开始，每一天的影响力为前两天的影响力的乘积再乘以 $a$ 的 $b$ 次方。用数学语言描述： $f(i)=f(i-1)*f(i-2)*a^b$ 。输出第n天的影响力，对1e9+7取模。

由于 $n<=10^{12}$ ，所以用快速幂。

先推几个公式，找规律： $f(1)=x\\f(2)=y\\f(3)=f(1)*f(2)*a^b=x*y*a^b\\f(4)=f(2)*f(3)*a^b=x*y^2*a^{2b} \\f(5)=f(3)*f(4)*a^b=x^2*y^3*a^{4b}$
从第三项开始x的幂次为1,1,2,3,5...是斐波那契数列

$y$ 的幂次为 $1,2,3,5...$ 是斐波那契数列

$a$ 的幂次为 $1,2,4,7,12...$ 不是斐波那契数列，但为前两项之和+1

同时还需用到一个数论的结论：费马小定理。

费马小定理： $p$ 为素数，对于任意整数 $x$ 都有 $x^p\equiv x(mod\space p)$ 。

若 $x$ 无法被 $p$ 整除，则有 $x^{p-1}\equiv 1(mod\space p)$ 。

根据费马小定理，我们可以通过矩阵快速幂求得逆元，因为由上式有 $x^{p-2}\equiv x^{-1}(mod\space p)$ ，不过跟这题没关系。

矩阵快速幂求斐波那契数列

由于数据范围很大，如果用递推求斐波那契数列显然会超时。

模板

typedef vector<int> vec;
typedef vector<vector<int>> mat;
typedef long long ll;

const ll mod=1e9+7;

// 首先求出两个矩阵乘积
mat mul(mat &A,mat &B){
    mat C(A.size(),vec(B[0].size()));
    for(int i=0;i<A.size();i++)
        for(int k=0;k<B.size();k++)
            for(int j=0;j<B[0].size();j++)
                C[i][j]=(C[i][j]+A[i][k]*B[k][j])%mod;
    return C;
}

//算快速幂
mat poww(mat A,ll n){
    mat B(A.size(),vec(A.size()));
    for(int i=0;i<A.size();i++)
        B[i][i]=1;
    //相当于朴素快速幂里面res=1; 将B矩阵对角线元素初始化为1
    while(n>0){
        if(n&1) B=mul(B,A);
        A=mul(A,A);
        n>>=1;
    }
    return B;

}

void fib(int n){
    mat A(2,vec(2));
    A[0][0]=A[0][1]=A[1][0]=1;
    A[1][1]=0;
    A=poww(A,n-1);
    printf("%d",A[0][0]);
    //A[0][0] 对应斐波那契数列的第 n 项
}

题解代码

#include <iostream>
#include <cstdio>
#include <vector>
using namespace std;

typedef long long ll;
typedef vector<ll> vec;
typedef vector<vector<ll>> mat;
const ll mod=1e9+6;   //费马小定理，用的模数
const ll mod2=1e9+7;  //答案所用模数

// 首先求出两个矩阵乘积
mat mul(mat &A,mat &B){        //用模数1，矩阵快速幂算的是x,y,a^b的幂，幂的模要用到费马小定理
    mat C(A.size(),vec(B[0].size()));
    for(int i=0;i<A.size();i++)
        for(int k=0;k<B.size();k++)
            for(int j=0;j<B[0].size();j++)
                C[i][j]=(C[i][j]+A[i][k]*B[k][j]%mod)%mod;
    return C;
}

//算快速幂
mat poww(mat A,ll n){
    mat B(A.size(),vec(A.size()));
    for(int i=0;i<A.size();i++)
        B[i][i]=1;
    //相当于朴素快速幂里面res=1; 将B矩阵对角线元素初始化为1
    while(n>0){
        if(n&1) B=mul(B,A);
        A=mul(A,A);
        n>>=1;
    }
    return B;

}
ll poww(ll a,ll n){                        //用模数2，算出最终答案
    // 底数为0时，直接返回0， 否则会卡10% 的数据
    if(a==0) return 0;
    ll res=1;
    while(n>0){
        if(n&1) res=res*a%mod2;
        a=a*a%mod2;
        n>>=1;
    }
    return res;

}


int main(){
    mat A(2,vec(2));
    A[0][0]=A[0][1]=A[1][0]=1;
    A[1][1]=0;
    //矩阵1
    /*
            1 1
            1 0
            用来算斐波那契数，也就是x,y的模数
    */
    mat B(3,vec(3));
    B={{1,1,1},{1,0,0},{0,0,1}};
    //矩阵2
    /*
        1 1 1
        1 0 0
        0 0 1
        用来算a^b的幂次
    */
    ll n,x,y,a,b;
    cin>>n>>x>>y>>a>>b;
    x%=mod2;y%=mod2;
    a%=mod2;                    //分别对mod2取模

    if(n==1){
        printf("%lld\n",x);                    //输出前两项
        return 0;
    }
    if(n==2){
        printf("%lld\n",y);
        return 0;
    }
    if(a==0||x==0||y==0){        //由前面费马小定理幂次同余模p的条件，a,x,y能直接被mod2整除，输出0
        printf("0\n");
        return 0;
    }
    a=poww(a,b);
    A=poww(A,n-2);
    B=poww(B,n-2);
    ll p=A[0][0], q=A[1][0], r=(B[1][0]+B[1][2])%mod;//分别对应y 的幂次，x的幂次，a^b的幂次
    ll ans=(poww(x,q)*poww(y,p)%mod2*poww(a,r))%mod2;  //快速幂求答案
    printf("%lld\n",ans%mod2);
    return 0;
}