Yet Another Number Sequence

题意：

$F_0 = 0, F_1 = 1, F_n = F_{n-1} + F_{n-2}$

求解$\sum_{i=1}^n{ F_i i^K } \ mod \ 10^9+7$。

解法：

记$S(n,m) = \sum_{i=1}^n { F_i i^m}$

这样有：

$$S(2n,m) = \sum_{i=1}^n{F_i i^m} + \sum_{i=1}^n{F_{i+n} (i+n)^m}$$

记$G = \lgroup \begin{matrix} 1 & 1 \\ 1 & 0 \end{matrix} \rgroup$

有$F_n = G_{1,1}$

将$S(n,m)$的含义换为对应求和的矩阵。

1.考虑从$S(n,m)$推到$S(2n,m)$

$$S(2n,m) = \sum_{i=1}^n{G^n i^m} + G^n \sum_{i=1}^n{G^{i} (i+n)^m}$$

$$S(2n,m) = \sum_{i=1}^n{G^n i^m} + G^n \sum_{i=1}^n{G^{i} \sum_{r=0}^m{i^rm^{m-r}C_m^r}}$$

$$S(2n,m) = \sum_{i=1}^n{G^n i^m} + G^n \sum_{r=0}^m{m^{m-r}C_m^r S(n,r)}$$

这样$O(K)$完成单次转移。

2.考虑从$S(n,m)$推到$S(n+1,m)$

$$S(n+1,m) = S(n,m) + (n+1)^m G^{n+1}$$

这样分治下去，总效率$O(K^2 logn)$

（注意本题中n过大，要先对于n取模再进行乘法）

 #include <iostream>

 #include <cstdio>

 #include <cstring>

 #define LL long long

 #define P 1000000007LL

 using namespace std;

 struct MA

 {

     LL a[][];

     void init()

     {

         memset(a,,sizeof(a));

     }

     MA operator*(const MA &x)const

     {

         MA c;

         for(int i=,j,k;i<;i++)

             for(j=;j<;j++)

             {

                 c.a[i][j]=;

                 for(k=;k<;k++)

                 {

                     c.a[i][j]+=a[i][k]*x.a[k][j]%P;

                     if(c.a[i][j]>=P) c.a[i][j]-=P;

                 }

             }

         return c;

     }

     MA operator*(const LL &x)const

     {

         MA c;

         for(int i=,j;i<;i++)

             for(j=;j<;j++)

                 c.a[i][j] = a[i][j]*x%P;

         return c;

     }

     MA operator+(const MA &x)const

     {

         MA c;

         for(int i=,j;i<;i++)

             for(j=;j<;j++)

                 c.a[i][j] = (a[i][j]+x.a[i][j])%P;

         return c;

     }

     void print()

     {

         puts("Matrix");

         for(int i=,j;i<;i++)

         {

             for(j=;j<;j++) cout<<a[i][j]<<' ';

             cout<<endl;

         }

     }

 };

 MA G,Gn;

 MA S[][];

 LL n,power[],C[][];

 int K,now;

 void solve(LL n)

 {

     if(n==)

     {

         now=;

         Gn = G;

         for(int i=;i<=K;i++) S[now][i] = G;

         return;

     }

     solve(n>>);

     now^=;

     MA tmp;

     power[]=;

     for(int k=;k<=K;k++) power[k] = power[k-]*((n>>1LL)%P)%P;

     for(int k=;k<=K;k++)

     {

         tmp.init();

         for(int r=;r<=k;r++)

             tmp = tmp + ( S[now^][r] * (C[k][r]*power[k-r]%P) );

         S[now][k] = S[now^][k] + (Gn * tmp);

     }

     Gn = Gn*Gn;

     if(n&)

     {

         Gn = Gn*G;

         power[]=1LL;

         for(int k=;k<=K;k++)

         {

             S[now][k] = S[now][k] + (Gn * power[k]);

             power[k+] = power[k]*(n%P)%P;

         }

     }

 }

 int main()

 {

     G.a[][]=; G.a[][]=;

     G.a[][]=; G.a[][]=;

     while(~scanf("%I64d%d",&n,&K))

     {

         C[][]=;

         for(int i=;i<=K;i++)

         {

             C[i][]=;

             for(int j=;j<=i;j++)

             {

                 C[i][j] = C[i-][j-]+C[i-][j];

                 if(C[i][j]>=P) C[i][j] -= P;

             }

         }

         solve(n);

         cout << S[now][K].a[][] << endl;

     }

     return ;

 }

巴特西

Yet Another Number Sequence

最新文章

热门文章