連立の漸化式を行列の掛け算に変形する方法はよく知られていますが、今回はその応用を考えます。
数字のすべての部分列の合計を求める
例えば$n=123$のすべての部分列は$\lbrace 1,2,3,12,23,123\rbrace$となりその合計を$A(n)$と表すと$A(123)=1+2+3+12+23+123=164$となります。
【例題 １】$A(123456789)$を求めよ
桁を左から１つずつ増やしながらに部分列の増分に注目します。ここで、
\begin{align}
&d_i : i番目の数字\\
&B_i : i番目の数字を加えたときに増える部分列の合計\\
&C_i : B_iの部分列の個数\\
\end{align}
| i | $d_i$ | $B_i$ | $C_i$ | $B_i$の漸化式 |
|---|---|---|---|---|
| 0 | 0 | 0 | 0 | - |
| 1 | 1 | $1$ | 1 | $1 + (0) \cdot 10 +1 \cdot 0$ |
| 2 | 2 | $2+12$ | 2 | $2 + (1) \cdot 10 +2 \cdot 1$ |
| 3 | 3 | $3+23+123$ | 3 | $3 + (2+12) \cdot 10 +3 \cdot 2$ |
| 4 | 4 | $4+34+234+1234$ | 4 | $4 + (3+23+123) \cdot 10 +4 \cdot 3$ |
| $\cdots$ | ||||
| i | $d_i$ | $\cdots \cdots$ | $C_i$ | $d_i + B_{i-1} \cdot 10 +d_i \cdot C_{i-1}$ |
ここで得られた漸化式をリストします。
\begin{align}
&A_i = A_{i-1} + B_{i} \ \ \ &\dots (1) \\
&B_i = d_i + B_{i-1} \cdot 10 +d_i \cdot c_{i-1} &\dots (2) \\
&C_i = C_{i-1} + 1 &\dots (3)
\end{align}
ここで$d_i$は$i-1$に依存しないので除くと、$(1)$だけが右辺に$B_i$があるのでこれに$(2)$を代入します。
行列にしやすいように多少順序を変えて、定数項のために$(4)$を加えます
\begin{align}
&A_i = A_{i-1} + \color{red}{10 \cdot B_{i-1} +d_i \cdot C_{i-1} + d_i} &\dots (1)’ \\
&B_i =0 + 10 \cdot B_{i-1} +d_i \cdot C_{i-1}+ d_i &\dots (2)' \\
&C_i =0 + 0 + C_{i-1} + 1 &\dots (3)' \\
&1 = 0 + 0 + 0 + 1 &\dots (4) \
\end{align}
ここから以下の行列の漸化式が得られました。
\begin{align}
\begin{bmatrix}
A_i \\
B_i \\
C_i \\
1 \\
\end{bmatrix}
&=
\begin{bmatrix}
1 & 10 & d_i & d_i \\
0 & 10 & d_i & d_i \\
0 & 0 & 1 & 1\\
0 & 0 & 0 & 1\\
\end{bmatrix}
\cdot
\begin{bmatrix}
A_{i-1} \\
B_{i-1} \\
C_{i-1} \\
1 \\
\end{bmatrix} \\
&= M_i \cdot
\begin{bmatrix}
A_{i-1} \\
B_{i-1} \\
C_{i-1} \\
1 \\
\end{bmatrix}
& \cdots (5)
\end{align}
したがって
\begin{align}
&\begin{bmatrix}
A(123456789) \\
- \\
- \\
- \\
\end{bmatrix}
= (M_9 M_8 \cdots M_2 M_1) \cdot
\begin{bmatrix}
0 \\
0 \\
0 \\
1 \\
\end{bmatrix} \\
&よって \\
&A(123456789) = (M_9 M_8 \cdots M_2 M_1)_{1,4}
\end{align}
これをコードにします。基本的には行列の掛け算なのでnumpyを使えば簡潔に書くことができました。
import numpy as np
def A(dl):
rdl = reversed(dl)
ret = np.identity(4, dtype=int) # identity matrix
for d in rdl:
ret @= np.array([[1,10, d, d],[0, 10, d, d],[0, 0, 1, 1], [0, 0, 0, 1]])
return ret
print(A([1,2,3,4,5,6,7,8,9])[0,3])
# 167657325
元の数字がパターンの繰り返しの時
漸化式を行列で表すと、繰り返しの演算が出て来た場合に「繰り返し２乗法」が使えるので高速化することができます。
【例題 ２】$n=123456789123456789...$と$123456789$が$10^3$回繰り返す数だったとき$A(n)$を$10^9+7$で割った余りを求めよ
$n$は$9 \times 10^3$桁のとてつもなく大きな数になりますが「フィボナッチ数列を行列を使ってO(log(n))で計算する」で作ったmtxPowMod を使って以下の式を計算します。
\begin{align}
&A(n) = ((M_9 M_8 \cdots M_2 M_1))^k_{1,4}
\end{align}
M = 10**9 + 7
print(mtxPowMod(A([1,2,3,4,5,6,7,8,9]),10**3, M)[0,3])
# 864817697
(開発環境：Google Colab)
この考え方はProject Euler, Problem 603: Substring Sums of Prime Concatenationsを解くのに役に立ちます。