前言
首先感谢所有帮助我卡常的大佬们。
题意
求 \((\sum_{i = 0}^{n} a^i b^{n - i})\mod p\) 的值（\(n \leq 10^{18}\)，\(p\) 不一定为质数）。
分析
看到数据范围，首先想到快速幂 \(+\) 乘法逆元或者矩阵乘法。
但是看到 \(p\) 不为质数，那就只有可能是矩阵乘法了。
设 \(S_n = \sum_{i = 0}^{n} a^i b^{n - i}\)，
则有 \(S_{n + 1} = \sum_{i = 0}^{n} a^ib^{n + 1 - i} + a^{n + 1}b^0 = b \sum_{i = 0}^{n} a^ib^{n - i} + a^{n + 1} = S_{n} \times b + a^{n + 1}\)。
有了递推式，我们就可以设出矩阵。
设向量 \(\begin{bmatrix} S_n & a^n \end{bmatrix}\)，则有如下等式：
\[\begin{bmatrix} S_n & a^n \end{bmatrix} \times \begin{bmatrix} b & 0 \\ a & a \end{bmatrix} = \begin{bmatrix} S_{n + 1} & a^{n + 1} \end{bmatrix}\]矩阵乘法不会的点这里矩阵乘法。
剩下的就交给矩阵快速幂了。
几个注意事项：
- 两个
long long相乘会爆，要用__int128。 - 矩阵乘法要展开循环。否则会超时。
- 一定要快读快写！！
代码示例
#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
using namespace std;
#define LL long long
#define PII pair<int, int>
#define PLL pair<LL, LL>
#define Int __int128
LL T;
Int ans[2], t[2][2];
namespace FastIO {
char buf[1 << 21], *_now = buf, *_end = buf;
#define getchar() (_now == _end && (_end = (_now = buf) + fread(buf, 1, 1 << 20, stdin), _now == _end) ? EOF : *_now ++ )
void read() { return; }
template <typename T, typename ...T2>
void read(T &s, T2 &...oth) {
s = 0; int f = 1; char ch = getchar();
for (; ch < '0' || ch > '9'; ch = getchar()) if (ch == '-') f = ~f + 1;
for (; ch >= '0' && ch <= '9'; ch = getchar()) s = (s << 1) + (s << 3) + (ch ^ 48);
s *= f; read(oth...); return;
}
void write() { return; }
template <typename T, typename ...T2>
void write(T s, T2 ...oth) {
if (s == 0) { putchar('0'), putchar('\n'); write(oth...); return; }
int stk[100], top = 0;
if (s < 0) putchar('-'), s = ~s + 1;
while (s) stk[ ++ top] = s % 10, s /= 10;
do putchar(stk[top -- ] + (1 << 4) + (1 << 5)); while (top);
putchar('\n'); write(oth...); return;
}
}
#define read FastIO::read
#define write FastIO::write
#define Inline __inline__ __attribute__((always_inline))
Inline void mul(Int ans[], Int a[], Int b[][2], LL p)
{
Int tmp[2] = {0};
tmp[0] = (tmp[0] + (Int)a[0] * b[0][0]) % p;
tmp[0] = (tmp[0] + (Int)a[1] * b[1][0]) % p;
tmp[1] = (tmp[1] + (Int)a[0] * b[0][1]) % p;
tmp[1] = (tmp[1] + (Int)a[1] * b[1][1]) % p;
memcpy(ans, tmp, sizeof tmp);
}
Inline void mul(Int ans[][2], Int a[][2], Int b[][2], LL p) {
Int tmp[2][2] = {0};
tmp[0][0] = (tmp[0][0] + (Int)a[0][0] * b[0][0]) % p;
tmp[0][0] = (tmp[0][0] + (Int)a[0][1] * b[1][0]) % p;
tmp[0][1] = (tmp[0][1] + (Int)a[0][0] * b[0][1]) % p;
tmp[0][1] = (tmp[0][1] + (Int)a[0][1] * b[1][1]) % p;
tmp[1][0] = (tmp[1][0] + (Int)a[1][0] * b[0][0]) % p;
tmp[1][0] = (tmp[1][0] + (Int)a[1][1] * b[1][0]) % p;
tmp[1][1] = (tmp[1][1] + (Int)a[1][0] * b[0][1]) % p;
tmp[1][1] = (tmp[1][1] + (Int)a[1][1] * b[1][1]) % p;
memcpy(ans, tmp, sizeof tmp);
}
Inline void solve(LL a, LL b, LL n, LL p) {
ans[0] = 1, ans[1] = 1;
t[0][0] = b, t[0][1] = 0, t[1][0] = a, t[1][1] = a;
while (n) {
if (n & 1) mul(ans, ans, t, p);
mul(t, t, t, p); n >>= 1;
}
write(ans[0]);
}
int main() {
read(T);
while (T--) {
LL n, a, b, p;
read(n, a, b, p);
solve(a, b, n, p);
}
return 0;
}