清北学堂p69

希望

直接用栈枚举即可。

#include <cstdio>
#include <cstring>
#include <iostream>
std::string back[1001], forward[1001];
int back_head = 0, back_tail = 0, forward_head = 0, forward_tail = 0;
int main() {
    //freopen("kami.in", "r", stdin);
    //freopen("kami.out", "w", stdout);
    std::string s, now = "http://www.acm.org/";
    while (std::cin >> s) {
        if (s == "BACK") {
            if (back_tail == back_head) {
                printf("Ignored\n");
                continue;
            } else {
                forward_tail++;
                forward[forward_tail] = now;
                now = back[back_tail];
                back_tail--;
            }
        } else if (s == "FORWARD") {
            if (forward_tail == forward_head) {
                printf("Ignored\n");
                continue;
            } else {
                back_tail++;
                back[back_tail] = now;
                now = forward[forward_tail];
                forward_tail--;
            }
        } else if (s == "VISIT") {
            back_tail++;
            back[back_tail] = now;
            std::cin >> now;
            forward_tail = forward_head;
        } else if (s == "QUIT") {
            return 0;
        }
        std::cout << now << std::endl;
    }
    //fclose(stdin);
    //fclose(stdout);
    return 0;
}

矩阵快速幂+皮萨诺周期

·矩阵快速幂

矩阵运算之矩阵乘法

矩阵相乘只有在第一个矩阵的列数和第二个矩阵的行数相同时才有意义。

A P \times M 的矩阵, B M \times Q 的矩阵,设矩阵 C 为矩阵 A B 的乘积,

其中矩阵 C 中的第 i 行第 j 列元素可以表示为:

C_{i,j} = \sum_{k=1}^MA_{i,k}B_{k,j}

如果没看懂上面的式子,没关系。通俗的讲,在矩阵乘法中,结果 C 矩阵的第 i 行第 j 列的数,就是由矩阵 A i M 个数与矩阵 B j M 个数分别相乘再相加得到的。

矩阵乘法满足结合律,不满足一般的交换律。

利用结合律,矩阵乘法可以利用快速幂的思想来优化。

矩阵快速幂

在斐波那契数列当中, F_1 = F_2 = 1 F_i = F_{i - 1} + F_{i - 2}(i \geq 3)

如果有一道题目让你求斐波那契数列第 n 项的值,最简单的方法莫过于直接递推了。但是如果 n 的范围达到了 10^{18} 级别,递推就不行了,稳 TLE。考虑矩阵加速递推。

Fib(n) 表示一个 1 \times 2 的矩阵 \left[ \begin{array}{ccc}F_n & F_{n-1} \end{array}\right] 。我们希望根据 Fib(n-1)=\left[ \begin{array}{ccc}F_{n-1} & F_{n-2} \end{array}\right] 推出 Fib(n)

试推导一个矩阵 \text{base} ,使 Fib(n-1) \times \text{base} = Fib(n) ,即 \left[\begin{array}{ccc}F_{n-1} & F_{n-2}\end{array}\right] \times \text{base} = \left[ \begin{array}{ccc}F_n & F_{n-1} \end{array}\right]

怎么推呢?因为 F_n=F_{n-1}+F_{n-2} ,所以 \text{base} 矩阵第一列应该是 \left[\begin{array}{ccc} 1 \\ 1 \end{array}\right] ,这样在进行矩阵乘法运算的时候才能令 F_{n-1} F_{n-2} 相加,从而得出 F_n 。同理,为了得出 F_{n-1} ,矩阵 \text{base} 的第二列应该为 \left[\begin{array}{ccc} 1 \\ 0 \end{array}\right]

综上所述: \text{base} = \left[\begin{array}{ccc} 1 & 1 \\ 1 & 0 \end{array}\right] 原式化为 \left[\begin{array}{ccc}F_{n-1} & F_{n-2}\end{array}\right] \times \left[\begin{array}{ccc} 1 & 1 \\ 1 & 0 \end{array}\right] = \left[ \begin{array}{ccc}F_n & F_{n-1} \end{array}\right]

转化为代码,应该怎么求呢?

定义初始矩阵 \text{ans} = \left[\begin{array}{ccc}F_2 & F_1\end{array}\right] = \left[\begin{array}{ccc}1 & 1\end{array}\right], \text{base} = \left[\begin{array}{ccc} 1 & 1 \\ 1 & 0 \end{array}\right] 。那么, F_n 就等于 \text{ans} \times \text{base}^{n-2} 这个矩阵的第一行第一列元素,也就是 \left[\begin{array}{ccc}1 & 1\end{array}\right] \times \left[\begin{array}{ccc} 1 & 1 \\ 1 & 0 \end{array}\right]^{n-2} 的第一行第一列元素。

注意,矩阵乘法不满足交换律,所以一定不能写成 \left[\begin{array}{ccc} 1 & 1 \\ 1 & 0 \end{array}\right]^{n-2} \times \left[\begin{array}{ccc}1 & 1\end{array}\right] 的第一行第一列元素。另外,对于 n \leq 2 的情况,直接输出 1 即可,不需要执行矩阵快速幂。

上述内容摘自OI-wiki

因此便可以得到矩阵快速幂的模板。

#include <cstdio>
#include <cstring>
#define MOD 1000000007
long long n, a[3], mul[3][3], res[3][3], tmp[3][3], tp[3]; 
void mul_1() {
    memset(tmp, 0, sizeof(tmp)); 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            for(int k = 1; k <= 2; k++)
                tmp[i][j] = (tmp[i][j] + res[i][k] * mul[k][j]) % MOD; 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            res[i][j] = tmp[i][j]; 
}
void mul_2() {
    memset(tmp, 0, sizeof(tmp)); 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            for(int k = 1; k <= 2; k++)
                tmp[i][j] = (tmp[i][j] + mul[i][k] * mul[k][j]) % MOD; 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            mul[i][j] = tmp[i][j]; 
}
void solve() {
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            tp[i] = (tp[i] + res[i][j] * a[j]) % MOD; 
    printf("%lld\n", tp[1]); 
}
int main() {
    scanf("%lld", &n); 
    if(n <= 2) {
        printf("1\n"); 
    } else {
        a[1] = a[2] = 1; 
        for(int i = 1; i <= 2; i++)
            res[i][i] = 1; 
        for(int i = 1; i <= 2; i++)
            for(int j = 1; j <= 2; j++)
                mul[i][j] = 1; 
        mul[2][2] = 0; 
        n -= 2; 
        while (n) {
            if (n & 1) mul_1(); 
            n >>= 1; 
            mul_2(); 
        }
        solve(); 
    }
    return 0; 
}

·皮亚诺周期

虽然我们知道了矩阵快速幂,但是, n=10^{100} ,直接求 f(f(n)) 的话无疑时间空间都得爆。

而斐波那契数列有一个很神奇的性质,皮亚诺周期。

自然数 n 皮萨诺周期(通常记为 π(n) )是指斐波那契数列模 n 后的周期

详见 维基百科-皮亚诺周期

知道了这一个之后,我们就能够求 f(f(n)) 辣。

首先求 f(f(n))\ \text{mod}\ {10^9+7} 的周期,打表程序如下:

#include <cstdio>
#define MOD 1000000007
int main() {
    long long a1 = 1, a2 = 1, cnt = 2, now;
    while (1) {
        cnt++;
        now = (a1 + a2 + MOD) % MOD;
        if (now == 1 && (now + a1 + MOD) % MOD == 1) {
            printf("%lld\n", cnt - 1);
            return 0;
        }
        a2 = a1;
        a1 = now;
    }
    return 0;
}

跑了好长时间,最后输出了 2000000016 这个数字,因此 f(f(n))\ \text{mod}\ {10^9+7} 的周期就是 2000000016

因此, f(n) 每次可以对 2000000016 取模,反正跑出来的答案是一样的。

f(n) 里面的 n 最大值可达 10^{100} ,这个程序仍然没法跑,我们便求 f(n)\ \text{mod}\ 2000000016 的周期,达标程序如下:

#include <cstdio>
#define MOD 2000000016
int main() {
    long long a1 = 1, a2 = 1, cnt = 2, now;
    while (1) {
        cnt++;
        now = (a1 + a2 + MOD) % MOD;
        if (now == 1 && (now + a1 + MOD) % MOD == 1) {
            printf("%lld\n", cnt - 1);
            return 0;
        }
        a2 = a1;
        a1 = now;
    }
    return 0;
}

最后得到 f(n) 的周期为 329616 ,同理, n 每次对 329616 取模就行了。

最终程序:

#include <cstdio>
#include <cstring>
#include <iostream>
long long MOD;
long long a[3], mul[3][3], res[3][3], tmp[3][3], tp[3]; 
void mul_1() {
    memset(tmp, 0, sizeof(tmp)); 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            for(int k = 1; k <= 2; k++)
                tmp[i][j] = (tmp[i][j] + res[i][k] * mul[k][j]) % MOD; 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            res[i][j] = tmp[i][j]; 
}
void mul_2() {
    memset(tmp, 0, sizeof(tmp)); 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            for(int k = 1; k <= 2; k++)
                tmp[i][j] = (tmp[i][j] + mul[i][k] * mul[k][j]) % MOD; 
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            mul[i][j] = tmp[i][j]; 
}
long long solve() {
    for(int i = 1; i <= 2; i++)
        for(int j = 1; j <= 2; j++)
            tp[i] = (tp[i] + res[i][j] * a[j]) % MOD; 
    return tp[1]; 
}
long long pow(long long n, long long mo) {
    memset(a, 0, sizeof(a));
    memset(mul, 0, sizeof(mul));
    memset(res, 0, sizeof(res));
    memset(tmp, 0, sizeof(tmp));
    memset(tp, 0, sizeof(tp));
    if(n <= 2) {
        return 1; 
    } else {
        MOD = mo;
        a[1] = a[2] = 1; 
        for(int i = 1; i <= 2; i++)
            res[i][i] = 1; 
        for(int i = 1; i <= 2; i++)
            for(int j = 1; j <= 2; j++)
                mul[i][j] = 1; 
        mul[2][2] = 0; 
        n -= 2; 
        while (n) {
            if (n & 1) mul_1(); 
            n >>= 1; 
            mul_2(); 
        }
        return solve(); 
    }
    return 0;
}
int main() {
    //freopen("na.in", "r", stdin);
    //freopen("na.out", "w", stdout);
    long long T;
    scanf("%lld", &T);
    while(T--) {
        std::string s;
        std::cin >> s;
        long long n = 0;
        for (int i = 0; i < s.length(); i++) {
            n = (n * 10 + (s[i] - '0') + 329616) % 329616;
        }
        //printf("%lld\n", n);
        long long ans1 = pow(n, 2000000016);
        long long ans = pow(ans1, 1000000007);
        printf("%lld\n", ans);
    }
    //fclose(stdin);
    //fclose(stdout);
    return 0;
}