题解 P3216 【[HNOI2011]数学作业】

· · 题解

简单的矩阵递推练手题

前言：听完了矩阵递推的网课后我打开了团队练习题单，发现自己并不会矩阵游戏中矩阵满足费马小定理的证明，于是只好来切这道题，但事实上，抛开细节对代码构筑能力的考验，这道题比矩阵游戏简单，甚至可以说在你明白这是道矩阵递推的题后，简直就是水题。\ 一开始先不想矩阵递推，不难推出公式：

long long ans=0; 
for(long long i=1;i<=n;i++)
{
    //ws[i]是i的位数 
    ans=ans*ws[i]+i;//不模的简洁版 
    ans=(ans*ws[i]%mod+i)%mod;//取模的细节版 
}

明显会爆longlong， n 可是小于等于 10^{18} 的啊。\ 思考矩阵递推，发现在矩阵里面搞位数这玩意儿实属不会，看看 n 的范围，想到按位数来搞矩阵快速幂，总共做 18 次，这样就不难处理位数了，而且代码跑的比香港记者还快。 然后想出矩阵如下（以下以位数为 1 时为例子）：

10 & 1 & 0 \\ 0 & 1 & 1 \\ 0 & 0 & 1 \end{vmatrix} \times sum[i-1] \\ i \\ 1 \end{vmatrix} = sum[i]\times10+i \\ i+1 \\ 1 \end{vmatrix} = sum[i] \\ i+1 \\ 1 \end{vmatrix}

又比如位数为3时：

1000 & 1 & 0 \\ 0 & 1 & 1 \\ 0 & 0 & 1 \end{vmatrix} \times sum[i-1] \\ i \\ 1 \end{vmatrix} = i+1 \\ 1 \end{vmatrix} = sum[i] \\ i+1 \\ 1 \end{vmatrix}

举了两个例子后思路就清晰了，我们可以先求1到10，再求10到100，然后是100到1000等等。

struct node
{
    int x,y;//行列数 
    int a[5][5];
    void empty()//清空哦 
    {
        x=y=0;
        for(int i=0;i<=4;i++)
        for(int j=0;j<=4;j++)
        {
            a[i][j]=0;
        }
    }
    void check()//验错用的
    {
        for(int i=1;i<=x;i++)
        {
            for(int j=1;j<=y;j++)
            printf("%lld ",a[i][j]);
            printf("\n");
        }
    }
}q[20],p,p2;//q表示对应位数下使用的3X1初始矩阵,p是用来维护递推关系的柿子，p是单位矩阵0. 
/*
    q[1]：
    0
    1
    1
    q[2]： 
    123456789
    10
    1 
*/ 
void pre()//初始化 
{
    int temp=1;
    for(int i=1;i<=18;i++)
    {
        q[i].x=3,q[i].y=1;
        q[i].a[2][1]=temp;
        temp=temp*10;
        q[i].a[3][1]=1;
    }
    p.x=p.y=p2.x=p2.y=3;
    p.a[1][1]=10;
    p.a[2][2]=p.a[1][2]=p.a[2][3]=p.a[3][3]=1;
    p2.a[1][1]=p2.a[2][2]=p2.a[3][3]=1;
}
void work(int now,int pre,int to)//现在的位数,上一次计算到的位置，这一次该计算到的位置 
{
    //q[now]是递推出来的3X1矩阵，p是左边用于递推的3X3矩阵 
    int temp=q[now-1].a[1][1];//继承 
    q[now].a[1][1]=temp;
    if(now==ws_ans)//如果当前位数和答案一致,可以合法的到达答案了 
    {
        node b=quickpow(p,n-pre);//矩阵快速幂 
        q[now]=b*q[now];//计算出最终的矩阵 
        printf("%lld",q[now].a[1][1]);
        return;
    }
    else q[now]=quickpow(p,to-pre)*q[now];//直接求解下一个位数的第一个矩阵 
    p.a[1][1]=p.a[1][1]*10%mod;//此处记得取模 
    work(now+1,to,(to+1)*10-1);//递归下一位数 
}

然后还要注意取模的问题，不要爆longlong，基本上自己手测几个样例之后就能过了。\ 送一组大样例，输入19191145141919 415411 ，输出246506。\ 原代码如下：（有些变量名与上面不一样，注意区分）

#include<cstdio>
#define int unsigned long long
using namespace std;
int n,m,ws;
struct node
{
    int x,y; 
    int a[5][5];
    void empty()
    {
        x=y=0;
        for(int i=0;i<=4;i++)
        for(int j=0;j<=4;j++)
        {
            a[i][j]=0;
        }
    }
    void check()
    {
        for(int i=1;i<=x;i++)
        {
            for(int j=1;j<=y;j++)
            printf("%lld ",a[i][j]);
            printf("\n");
        }
    }
}q[20],p,p2;
void pre()
{
    int temp=1;
    for(int i=1;i<=18;i++)
    {
        q[i].x=3,q[i].y=1;
        q[i].a[2][1]=temp;
        temp=temp*10;
        q[i].a[3][1]=1;
    }
    p.x=p.y=p2.x=p2.y=3;
    p.a[1][1]=10;
    p.a[2][2]=p.a[1][2]=p.a[2][3]=p.a[3][3]=1;
    p2.a[1][1]=p2.a[2][2]=p2.a[3][3]=1;
}
node operator*(node a,node b)
{
    node o;
    o.empty();
    int x=a.x;
    int y=a.y;
    int z=b.y;
    for(int i=1;i<=x;i++)
    for(int j=1;j<=y;j++)
    for(int k=1;k<=z;k++)
        o.a[i][k]=(o.a[i][k]+(a.a[i][j]%m)*(b.a[j][k]%m))%m;
    o.x=x;
    o.y=z;
    return o;
}
node quickpow(node a,int b)
{
    node ans=p2,temp=a;
    while(b)
    {
        if(b&1)ans=ans*temp;
        temp=temp*temp;
        b/=2;
    }
    return ans;
}
void work(int now,int pre,int to)
{
    int vra=q[now-1].a[1][1];
    q[now].a[1][1]=vra;
    if(now==ws)
    {
        node b=quickpow(p,n-pre);
        q[now]=b*q[now];
        printf("%lld",q[now].a[1][1]);
        return;
    }
    else q[now]=quickpow(p,to-pre)*q[now];
    p.a[1][1]=p.a[1][1]*10%m;
    work(now+1,to,(to+1)*10-1);
} 
signed main()
{
    scanf("%lld%lld",&n,&m);
    pre();
    int temp=n;
    while(temp)
    {
        temp/=10;
        ws++;
    }
    work(1,0,9);
}

撰写不易，求通过。