洛谷P3216 [HNOI2011]数学作业 题解 矩阵快速幂

题目链接:https://www.luogu.com.cn/problem/P3216

题目大意:定义函数 (f(n)) 为将 (1 sim n) 连起来得到的数,求 (f(n) ext{ mod } m) 的结果。

解题思路:

状态转移方程为:

[f_i = 10^{ lfloor log_{10} i floor } imes f_{i-1} + i ]

转移矩阵为:

[egin{bmatrix} f_n \ n \ 1 end{bmatrix} = egin{bmatrix} 10^k & 1 & 1 \ 0 & 1 & 1 \ 0 & 0 & 1 end{bmatrix} imes egin{bmatrix} f_{n-1} \ n-1 \ 1 end{bmatrix} ]

其中的 (k) 表示 (lfloor log_{10} i floor + 1)(即 (n) 所占的位数)。但是 (k) 是变化的,没有办法直接矩阵乘法。但是可以发现:

  • (n in [1, 9]) 时,(k = 1)
  • (n in [10, 99]) 时,(k = 2)
  • (n in [100, 999]) 时,(k=3)
  • …… ……
  • (n = 10^{18}) 时,(k=19)

所以 (k) 总共有 (19) 种不同的情况,我们只需要针对不同的 (k) 做矩阵快速幂就可以了。

示例代码:

#include <bits/stdc++.h>
using namespace std;
long long n, m;
struct Matrix {
    long long a[3][3];
    Matrix operator * (Matrix b) const {
        Matrix c;
        memset(c.a, 0, sizeof(c.a));
        for (int i = 0; i < 3; i ++)
            for (int j = 0; j < 3; j ++)
                for (int k = 0; k < 3; k ++)
                    c.a[i][j] = (c.a[i][j] + a[i][k] * b.a[k][j] % m) % m;
        return c;
    }
    Matrix operator ^ (long long n) const {
        Matrix b, c;
        memcpy(b.a, a, sizeof(a));
        memset(c.a, 0, sizeof(c.a));
        for (int i = 0; i < 3; i ++) c.a[i][i] = 1;
        while (n) {
            if (n % 2) c = c * b;
            b = b * b;
            n /= 2;
        }
        return c;
    }
} ans, tmp;
long long x[20];
int main() {
    cin >> n >> m;
    long long t1 = 1, t2 = 9;
    for (int i = 0; i < 18; i ++) {
        x[i] = min(t2, n) - t1 + 1;
        if (n < t2) break;
        t1 *= 10, t2 = t2 * 10 + 9;
    }
    if (n == 1000000000000000000LL) x[18] = 1;
    memset(ans.a, 0, sizeof(ans.a));
    for (int i = 0; i < 3; i ++) ans.a[i][i] = 1;
    for (int i = 18; i >= 0; i --) {
        if (!x[i]) continue;
        memset(tmp.a, 0, sizeof(tmp.a));
        tmp.a[0][0] = 1;
        for (int j = 0; j <= i; j ++) tmp.a[0][0] = tmp.a[0][0] * 10 % m;
        tmp.a[0][1] = tmp.a[0][2] = tmp.a[1][1] = tmp.a[1][2] = tmp.a[2][2] = 1;
        ans = ans * (tmp ^ x[i]);
    }
    cout << ans.a[0][2] << endl;
    return 0;
}
原文地址:https://www.cnblogs.com/quanjun/p/13957450.html