参考了某大佬的
我们可以根据(s[n-2], a[n-1]^2, a[n-1]*a[n-2], a[n-2]^2) * A = (s[n-1], a[n]^2, a[n]*a[n-1], a[n-1]^2)
能够求出关系矩阵
|1 0 0 0 |
A = |1 x^2 x 1 |
|0 2*x*y y 0 |
|0 y^2 0 0 |
这样就A了!
#include <iostream> #include <cstdio> #include <cstring> #include <algorithm> using namespace std; typedef long long ll; const ll Mod = 10007; const int N = 5; int msize; struct Mat { ll mat[N][N]; }; Mat operator *(Mat a, Mat b) { Mat c; memset(c.mat, 0, sizeof(c.mat)); for(int k = 0; k < msize; ++k) for(int i = 0; i < msize; ++i) if(a.mat[i][k]) for(int j = 0; j < msize; ++j) if(b.mat[k][j]) c.mat[i][j] = (a.mat[i][k] * b.mat[k][j] + c.mat[i][j])%Mod; return c; } Mat operator ^(Mat a, ll k) { Mat c; memset(c.mat,0,sizeof(c.mat)); for(int i = 0; i < msize; ++i) c.mat[i][i]=1; for(; k; k >>= 1) { if(k&1) c = c*a; a = a*a; } return c; } int main() { ll n,x,y; msize = 4; while(~scanf("%I64d%I64d%I64d",&n,&x,&y)) { Mat A; A.mat[0][0] = 1, A.mat[0][1] = 1, A.mat[0][2] = 0, A.mat[0][3] = 0; A.mat[1][0] = 0, A.mat[1][1] = x*x%Mod, A.mat[1][2] = 2*x*y%Mod, A.mat[1][3] = y*y%Mod; A.mat[2][0] = 0, A.mat[2][1] = x, A.mat[2][2] = y, A.mat[2][3] = 0; A.mat[3][0] = 0, A.mat[3][1] = 1, A.mat[3][2] = 0, A.mat[3][3] = 0; A = A^n; printf("%I64d ", (A.mat[0][0] + A.mat[0][1] + A.mat[0][2] + A.mat[0][3])%Mod); } return 0; }