Yukicoder No.426 往復漸化式
$b_{i-1}$を$b_{i+1}$だと思い込み$b_n$には気付いたが$b_0$のtypoだと無視した結果、誤読で時間を潰し、editorialを見て気付いた。
ついでに $ \left( \begin{matrix}
6i & 6i+1 & 6i+2 \
6i+3 & 6i+4 & 6i+5
\end{matrix} \right) $ の$i$は虚数だと思ってた。
solution
segment木。$O(N + Q \log N)$。
漸化式を整理すると以下のようになる。
\[\begin{array}{ccl} a_i & = & A\_{i-1} A\_{i-2} \dots A_0 a_0 \\\\ b_i & = & B\_{i+1} A\_{i+2} \dots B_n b_n \\\\ & & + \; S\_{i+1} a\_{i+1} \\\\ & & + \; B\_{i+1} S\_{i+2} A\_{i+1} a\_{i+1} \\\\ & & + \; B\_{i+1} B\_{i+2} S\_{i+3} A\_{i+2} A\_{i+1} a\_{i+1} \\\\ & & + \; \vdots \\\\ & & + \; B\_{i+1} B\_{i+2} \dots B\_{n-1} S_n A\_{n-1} \dots A\_{i+2} A\_{i+1} a\_{i+1} \end{array}\]ただし
\[S_i = \left( \begin{matrix} 6i & 6i+1 & 6i+2 \\\\ 6i+3 & 6i+4 & 6i+5 \end{matrix} \right)\]以下をsegment木で求めるのは容易である。
- $A_{i-1} A_{i-2} \dots A_0$
- $B_{i+1} A_{i+2} \dots B_n$
残る部分について、$a_{i+1}$をくくり出せば、次を求めればよい。
\[\begin{array}{cc} & S\_{i+1} \\\\ + & B\_{i+1} S\_{i+2} A\_{i+1} \\\\ + & B\_{i+1} B\_{i+2} S\_{i+3} A\_{i+2} A\_{i+1} \\\\ + & \vdots \\\\ + & B\_{i+1} B\_{i+2} \dots B\_{n-2} S_n A\_{n-2} \dots A\_{i+2} A\_{i+1} \\\\ + & B\_{i+1} B\_{i+2} \dots B\_{n-2} B\_{n-1} S_n A\_{n-1} A\_{n-2} \dots A\_{i+2} A\_{i+1} \end{array}\]これを変形すると例えば以下のようになる。
\[\begin{array}{crcl} & & S\_{i+1} & \\\\ + & & B\_{i+1} S\_{i+2} A\_{i+1} & \\\\ + & & B\_{i+1} B\_{i+2} S\_{i+3} A\_{i+2} A\_{i+1} & \\\\ + & B\_{i+1} B\_{i+2} B\_{i+3} & S\_{i+4} & A\_{i+3} A\_{i+2} A\_{i+1} \\\\ + & B\_{i+1} B\_{i+2} B\_{i+3} & B\_{i+4} S\_{i+5} A\_{i+4} & A\_{i+3} A\_{i+2} A\_{i+1} \\\\ + & \vdots & \vdots & \vdots \\\\ + & B\_{i+1} B\_{i+2} B\_{i+3} & B\_{i+4} \dots B\_{n-2} S_n A\_{n-2} \dots A\_{i+4} & A\_{i+3} A\_{i+2} A\_{i+1} \\\\ + & B\_{i+1} B\_{i+2} B\_{i+3} & B\_{i+4} \dots B\_{n-2} B\_{n-1} S_n A\_{n-1} A\_{n-2} \dots A\_{i+4} & A\_{i+3} A\_{i+2} A\_{i+1} \end{array}\]こう分割することにより、$A_i, B_i$単体の区間を用いて、三角形のそれの区間を合成できる。
implementation
segment木を再帰なしで実装するやつをやってみたらほぼ$2$倍速ぐらいになった: recursion $\to$ loop
#include <iostream>
#include <vector>
#include <array>
#include <functional>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
using ll = long long;
using namespace std;
const int mod = 1e9+7;
template <typename T, size_t H, size_t W>
using matrix = array<array<T,W>,H>;
template <typename T, size_t A, size_t B, size_t C>
matrix<T,A,C> operator * (matrix<T,A,B> const & a, matrix<T,B,C> const & b) {
matrix<T,A,C> c = {};
repeat (y,A) {
repeat (z,B) {
repeat (x,C) {
c[y][x] = (c[y][x] + a[y][z] *(ll) b[z][x]) % mod;
}
}
}
return c;
}
template <typename T, size_t H, size_t W>
array<T,H> operator * (matrix<T,H,W> const & a, array<T,W> const & b) {
array<T,H> c = {};
repeat (y,H) {
repeat (z,W) {
c[y] = (c[y] + a[y][z] *(ll) b[z]) % mod;
}
}
return c;
}
template <typename T, size_t H, size_t W>
matrix<T,H,W> operator + (matrix<T,H,W> const & a, matrix<T,H,W> const & b) {
matrix<T,H,W> c;
repeat (y,H) {
repeat (x,W) {
c[y][x] = (a[y][x] + b[y][x]) % mod;
}
}
return c;
}
template <typename T, size_t N>
array<T,N> operator + (array<T,N> const & a, array<T,N> const & b) {
array<T,N> c;
repeat (i,N) {
c[i] = (a[i] + b[i]) % mod;
}
return c;
}
template <typename T, size_t H, size_t W>
matrix<T,H,W> matrix_zero() { return {}; }
template <typename T, size_t N>
matrix<T,N,N> matrix_unit() { matrix<T,N,N> a = {}; repeat (i,N) a[i][i] = 1; return a; }
template <typename T>
struct segment_tree { // on monoid
int n;
vector<T> a;
function<T (T,T)> append; // associative
T unit; // unit
segment_tree() = default;
segment_tree(int a_n, T a_unit, function<T (T,T)> a_append) {
n = 1; while (n < a_n) n *= 2;
a.resize(2*n-1, a_unit);
unit = a_unit;
append = a_append;
}
template <typename F>
void point_update(int i, F func) { // 0-based
a[i+n-1] = func(a[i+n-1]);
for (i = (i+n)/2; i > 0; i /= 2) {
a[i-1] = append(a[2*i-1], a[2*i]);
}
}
T point_get(int i) {
return range_concat(i, i+1);
}
T range_concat(int l, int r) { // 0-based, [l, r)
T lacc = unit, racc = unit;
for (l += n, r += n; l < r; l /= 2, r /= 2) {
if (l % 2 == 1) lacc = append(lacc, a[(l ++) - 1]);
if (r % 2 == 1) racc = append(a[(-- r) - 1], racc);
}
return append(lacc, racc);
}
};
struct range_t {
matrix<int,3,3> A;
matrix<int,2,2> B;
matrix<int,2,3> S;
};
int main() {
int n; cin >> n;
array<int,3> a0; repeat (i,3) cin >> a0[i];
array<int,2> bn; repeat (i,2) cin >> bn[i];
segment_tree<range_t> segtree(n+1, (range_t) { matrix_unit<int,3>(), matrix_unit<int,2>(), matrix_zero<int,2,3>() }, [&](range_t const & l, range_t const & r) {
range_t m;
m.A = r.A * l.A;
m.B = l.B * r.B;
m.S = l.S + l.B * r.S * l.A;
return m;
});
repeat (i,n+1) {
segtree.point_update(i, [&](range_t m) {
repeat (y,2) repeat (x,3) m.S[y][x] = 6*i + 3*y+x;
return m;
});
}
int query; cin >> query;
while (query --) {
string type; int i; cin >> type >> i;
if (type == "a") {
segtree.point_update(i, [](range_t m) {
repeat (y,3) repeat (x,3) cin >> m.A[y][x];
return m;
});
} else if (type == "b") {
segtree.point_update(i, [](range_t m) {
repeat (y,2) repeat (x,2) cin >> m.B[y][x];
return m;
});
} else if (type == "ga") {
range_t m = segtree.range_concat(0, i);
array<int,3> ai = m.A * a0;
cout << ai[0] << ' ' << ai[1] << ' ' << ai[2] << endl;
} else if (type == "gb") {
range_t l = segtree.range_concat(0, i+1);
range_t r = segtree.range_concat(i+1, n+1);
array<int,2> bi = r.B * bn + r.S * l.A * a0;
cout << bi[0] << ' ' << bi[1] << endl;
}
}
return 0;
}