LightOJ1132---Summing up Powers （矩阵快速幂+二项式定理（简单推理））

最新推荐文章于 2022-11-24 22:43:07 发布

原创最新推荐文章于 2022-11-24 22:43:07 发布 · 600 阅读

1 ·

本内容遵循CC 4.0 BY-SA版权协议

ACM的进程同时被 3 个专栏收录

293 篇文章

订阅专栏

ACM竞赛

292 篇文章

订阅专栏

【数论】--矩阵快速幂

15 篇文章

订阅专栏

本文介绍了一种使用矩阵快速幂解决特定类型数列求和问题的方法。通过引入二项式定理，将复杂递推式转化为矩阵形式，利用快速幂运算高效计算大规模数值。适用于N值较大时的优化计算。

【题目来源】：https://vjudge.net/problem/LightOJ-1132
【题意】
题意呢，就像题面描述的一样。。。
【思路】
首先，这道题因为存在递推式，并且，N非常之大，所以想到矩阵快速幂，接着，整体思路为（假设前n项和为Sn）：
S(n+1)=Sn+(n+1)^k
然后按照一般的步骤（将左边，右边化为后，前两个状态，也即是说，后一状态可由前一状态递推而来），但是呢，遇到了一点阻碍，怎么样才可以将(n+1)^k化为关于n的递推式呢？
进而想到二项式定理：
(a+b)^k=C(0,k)*a^k*b^0+C(1,k)*a^k-1*b^1+…+C(k,k)*a^0*b^k;
而在这里（1的多少次方均为1，乘法不做记录）：
(n+1)^k=C(0,k)*n^k+C(1,k)*n^k-1+…+C(k,k)*n^0;
所以这道题化为矩阵快速幂的递推式就出来了：
S(n+1)=Sn+C(0,k)*n^k+C(1,k)*n^k-1+…+C(k,k)*n^0;
写出系数矩阵：
这里写图片描述
【代码】

#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define LL long long
using namespace std;
LL N,mod;
int K;
struct mat
{
    LL a[55][55];
};
LL z_h[55][55];//预处理k的所有可能性的二项式系数
mat operator*(mat s,mat t)
{
    mat r;
    for(int i=1; i<=K+2; i++)
        for(int j=1; j<=K+2; j++)
        {
            r.a[i][j]=0;
            for(int k=1; k<=K+2; k++)
                if(s.a[i][k]&&t.a[k][j])
                {
                    r.a[i][j]+=s.a[i][k]*t.a[k][j]%mod;
                    r.a[i][j]%=mod;
                }
        }

    return r;
}
void init()
{
    mod=1;
    for(int i=1; i<=32; i++)
    {
        mod*=2;
    }
    for(int i=1; i<=50; i++)
    {
        z_h[i][1]=1;
        int x=i;
        for(int j=2; j<=i+1; j++)
        {
            z_h[i][j]=z_h[i][j-1]*x/(j-1);
            x--;
        }
    }
}
mat pow_mat(mat base)
{
    mat ans;
    memset(ans.a,0,sizeof(ans.a));
    for(int i=1; i<=K+2; i++)
    {
        ans.a[i][i]=1;
    }
    N-=1;
    while(N)
    {
        if(N&1)
        {
            ans=ans*base;
        }
        base=base*base;
        N>>=1;
    }
    return ans;

}
int main()
{
    init();
    int T,cases=1;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%lld%d",&N,&K);
        mat base;
        memset(base.a,0,sizeof(base.a));
        int x=K;
        for(int i=2; i<=K+2; i++)
        {
            int y=1;
            for(int j=i; j<=K+2; j++)
            {
                base.a[i][j]=z_h[x][y++];
            }
            x--;
        }
        base.a[1][1]=1;
        base.a[K+2][K+2]=1;
        for(int i=2; i<=K+2; i++)
        {
            base.a[1][i]=base.a[2][i];
        }
        mat ans=pow_mat(base);
        LL sum=0;
        for(int i=1; i<=K+2; i++)
        {
            sum+=ans.a[1][i]%mod;
            sum%=mod;
        }
        printf("Case %d: %lld\n",cases++,sum);
    }
}