首先考虑 DP，设 \(f_n,g_n\) 分别表示把 \(n\) 个硬币移动到下一个圈 / 下下一个圈的最小步数。
若要把所有硬币移动到下一个圈，首先要把面值最大的移动到下一个圈，然后再把剩下的移到这个圈。
为了将面值最大的移动到下一个圈，要先把剩下的先移到下下个圈，然后再把最大的移到下一个圈。
显然这样可以构建起 \(f_n\) 的 DP 转移方程： \[f_n = 2g_{n-1}+1\]
若要把所有硬币移动到下下个圈，首先要把面海最大的移动到下下个圈，然后再把剩下的移到这个圈。
为了将面值最大的移动到下下个圈，要先把剩下的先移到下下个圈，然后再把最大的移到下一个圈。
然后再将剩下的移回原来的圈，再将最大的移动到下下个圈，然后再把剩下的移到下下个圈。
总结起来就是 \[g_n = 2g_{n-1} + f_{n-1} + 2\]
（我也知道很啰嗦不过这些我相信读者自己能推出来）
然后这里正解是直接矩阵快速幂 + 分块预处理……
我选择推通项（
设 \(F,G\) 为 \(f,g\) 的生成函数，则 \[ \begin{aligned} F(x) &= 2xG(x) + \frac1{1-x} - 1 \\ G(x) &= 2xG(x) + xF(x) + \frac2{1-x} - 2 \end{aligned} \]
用第一条式子替换掉第二条式子中的 \(F(x)\)，有 \[ \begin{aligned} G(x) &= 2xG(x) + x(2xG(x) + \frac1{1-x} - 1) + \frac2{1-x} - 2 \\ &= 2x^2G(x) + 2xG(x) + \frac{x+2}{1-x} - x - 2 \\ (1-2x-2x^2)G(x) &= \frac x{1-x}(x+2) \\ G(x) &= \frac x{(1-x)[1-(1+\sqrt3)x][1-(1-\sqrt3)x]}(x+2) \\ G(x) &= \frac16(x+2)\left[\frac{1+\sqrt3}{1-(1+\sqrt3)x} + \frac{1-\sqrt3}{1-(1-\sqrt3)x}\right] \\ g_n &= (1+\sqrt3)^n \left(\frac12+\frac13\sqrt3\right) + (1-\sqrt3)^n \left(\frac12-\frac13\sqrt3\right) - 1 \end{aligned} \] （省去了因式分解部分分式等一系列步骤）
然而，可惜的是！\(3\) 不是模 \(998244353\) 的二次剩余！
不过我们依然可以扩域，即把每个数表示为 \(a + b\sqrt3\) 的形式，然后分块预处理幂。
代码： 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
using namespace std;
const int W = 1 << 20;
const int mod = 998244353;
const int inv[] = {499122177,332748118};
typedef pair<int,int> cp;
inline cp operator+(const cp &a,const cp &b)
{
return cp(add(a.first,b.first),add(a.second,b.second));
}
inline cp operator-(const cp &a,const cp &b)
{
return cp(dec(a.first,b.first),dec(a.second,b.second));
}
inline cp operator*(const cp &a,const cp &b)
{
return cp(((long long)a.first * b.first % mod + 3LL * a.second % mod * b.second % mod) % mod,((long long)a.first * b.second % mod + (long long)a.second * b.first % mod) % mod);
}
cp pw1[2][W + 5],pw2[2][W + 5];
inline cp query1(long long n)
{
return pw1[0][n & (W - 1)] * pw1[1][n >> 20];
}
inline cp query2(long long n)
{
return pw2[0][n & (W - 1)] * pw2[1][n >> 20];
}
inline int query_g(long long n)
{
if(n <= 0)
return 0;
int ret = (query1(n) * cp(inv[0],inv[1]) + query2(n) * cp(inv[0],mod - inv[1])).first;
return ret = dec(ret,1);
}
inline int query_f(long long n)
{
if(n <= 0)
return 0;
int ret = query_g(n - 1);
ret = add(ret,ret);
return ret = add(ret,1);
}
int T,ans1,ans2;
long long n;
int main()
{
freopen("riemannian.in","r",stdin),freopen("riemannian.out","w",stdout);
pw1[0][0] = pw1[1][0] = pw2[0][0] = pw2[1][0] = cp(1,0);
for(register int i = 1;i <= W;++i)
pw1[0][i] = pw1[0][i - 1] * cp(1,1),
pw2[0][i] = pw2[0][i - 1] * cp(1,mod - 1);
for(register int i = 1;i < W;++i)
pw1[1][i] = pw1[1][i - 1] * pw1[0][W],
pw2[1][i] = pw2[1][i - 1] * pw2[0][W];
scanf("%d",&T);
for(;T;--T)
scanf("%lld",&n),ans1 ^= query_f(n),ans2 ^= query_g(n);
printf("%d %d\n",ans1,ans2);
}