Number Sequence（HDU-1005）

---

1 题意

  已知$f_1=f_2=1$，求$f_n=(A{f}_{n-1}+B{f}_{n-2})mod7$，其中$1⩽A,B⩽1000,1⩽n⩽10^8$
  链接：link

---

2 思路

2.1 循环节

  命题：$f_n$序列是循环序列（从某个点开始），且循环节长度小于50。
  证明：设$a,b$是不超过7的非负整数，则二元组$(a,b)$一共有49种可能，所以$(f_i,f_{i+1}),(i⩾1)$的前$50$项一定存在重复元素。
  假设$(f_n,f_{n+1})=(f_m,f_{m+1}),(1⩽n<m⩽50)$。即$f_n=f_m,f_{n+1}=f_{m+1}$，由递推公式$f_n=(A{f}_{n-1}+B{f}_{n-2})mod7$可得出$f_{n+2}=f_{m+2}$，进而$f_{n+3}=f_{m+3},f_{n+4}=f_{m+4},...$这样就可得知该序列从$n$开始以$m-n$为循环长度循环。

  推论1：如果存在$f_n=f_m,f_{n+1}=f_{m+1}$，则该序列从$n$开始以$m-n$为循环长度循环。

  推论2：存在$1⩽n⩽49$，使得$f_n=f_{50},f_{n+1}=f_{51}$。
  证明：已知一定存在$(f_n,f_{n+1})=(f_m,f_{m+1}),(1⩽n<m⩽50)$，并且可以推出$f_{n+2}=f_{m+2},f_{n+3}=f_{m+3},f_{n+4}=f_{m+4},...$那么当$m+k=50$时，可使得$f_{n+k}=f_{50},f_{n+k+1}=f_{51}$，其中$1⩽n+k<50$。

  根据推论1和推论2，可在$[1,50]$中找到一个长度小于$50$的循环节。

2.1.1 时间复杂度分析

1. 计算前$51$项的时间复杂度为$\mathcal{O}(51)$。
2. 在$1⩽n⩽49$寻找$f(n)=f(50),f(n+1)=f(51)$的时间复杂度为$\mathcal{O}(49)$。
3. 计算$f(n)$的时间复杂度为$\mathcal{O}(1)$。

  总的时间复杂度为$\mathcal{O}(1)$。

2.1.2 代码

#include<iostream>
#include<cstdio>
using namespace std;
int a[52]={0,1,1},A,B,n;
int main(){
    while(scanf("%d %d %d",&A,&B,&n)&&A+B+n){
        int step,first;
        for(int i=3;i<=51;i++) a[i]=(A*a[i-1]+B*a[i-2])%7;
        for(int i=1;i<=48;i++){
            if(a[50]==a[i]&&a[51]==a[i+1]){
                first=i;
                break;
            }
        }
        step=50-first;
        if(n<50) printf("%d\n",a[n]);
        else printf("%d\n",a[(n-first)%step+first]);
    }
    return 0;
}

2.2 矩阵快速幂

  关于矩阵快速幂的详细介绍请参考：link
  这道题的矩阵快速幂的公式为：$\left[\begin{array}{c}{f}_n\\ f_{n-1}\end{array}\right]={\left[\begin{array}{cc}A& B\\ 1& 0\end{array}\right]}^{n-2}\left[\begin{array}{c}{f}_2\\ f_1\end{array}\right]$

2.2.1 时间复杂度分析

  矩阵快速幂的时间复杂度为$\mathcal{O}(m^3\log n)$，其中$m$为矩阵大小，$n$为幂指数。
  对于这道题来说，时间复杂度为$\mathcal{O}(\log n)$。

2.2.2 代码

#include<iostream>
#include<cstdio>
using namespace std;
const int MAXN=15;
const int MOD=7;
int A,B,n,sz=2;
struct Matrix{
    int dt[MAXN][MAXN];
    Matrix(int tp=0){
        for(int i=1;i<=sz;i++){
            for(int j=1;j<=sz;j++){
                if(i==j) dt[i][j]=tp;
                else dt[i][j]=0;
            }
        }
    }
    Matrix operator+(const Matrix& a){
        Matrix res;
        for(int i=1;i<=sz;i++){
            for(int j=1;j<=sz;j++){
                res.dt[i][j]=(a.dt[i][j]+dt[i][j])%MOD;
            }
        }
        return res;
    }
    Matrix operator-(const Matrix& a){
        Matrix res;
        for(int i=1;i<=sz;i++){
            for(int j=1;j<=sz;j++){
                res.dt[i][j]=(dt[i][j]-a.dt[i][j])%MOD;
            }
        }
        return res;
    }
    Matrix operator*(const Matrix& a){
        Matrix res;
        for(int i=1;i<=sz;i++){
            for(int j=1;j<=sz;j++){
                for(int k=1;k<=sz;k++){
                    res.dt[i][j]=(res.dt[i][j]+dt[i][k]*a.dt[k][j])%MOD;
                }
            }
        }
        return res;
    }
};
Matrix qpow(const Matrix& a,int b){
    Matrix res(1),rem=a;
    while(b){
        if(b&1) res=res*rem;
        rem=rem*rem;
        b>>=1;
    }
    return res;
}
int main(){
    while(scanf("%d %d %d",&A,&B,&n)&&A+B+n){
        if(n<=2){
            printf("1\n");
            continue;
        }
        Matrix mat;
        mat.dt[1][1]=A;mat.dt[1][2]=B;
        mat.dt[2][1]=1;
        Matrix res=qpow(mat,n-2);
        printf("%d\n",(res.dt[1][1]+res.dt[1][2])%7);
    }
    return 0;
}

2.3 算法比较

  从最坏的时间复杂度来看，循环节更具有优势，因为循环节最多计算$51+49+1=101$次，而矩阵快速幂要计算$2^3\times \log (10^8)=212$次。
  总的来说循环节的时间复杂度只与模有关，而矩阵快速幂的时间复杂度只与n有关。所以当模较大时应当采用矩阵快速幂，当模较小且n较大时可采用循环节的做法。