Day 061-Q3 — 多項式行列 Gaussian Elimination over FPS

2026-06-14 赤色 Master / Phase 8+ ★★★★★★★★★ FPS / ガウス消去 / 行列式 / 逆行列 / Newton法

問題

$d \times d$ の正方行列 $A$ が与えられる。各エントリ $A[i][j]$ は形式的冪級数（多項式）であり、次数 $K$ 未満の係数が整数で与えられる。

以下を $\bmod x^K$, $\bmod p$ で計算せよ（$p = 998244353$）：

  1. $\det(A)$ の $\bmod x^K$ における係数列
  2. $A^{-1}$ の $(0,0)$ エントリの $\bmod x^K$ における係数列

制約

パラメータ範囲
$d$$1 \le d \le 8$
$K$$1 \le K \le 1000$
係数$\in [0, p)$
行列$\det(A) \ne 0 \pmod{x^K}$（可逆）

入出力例

入力例 1

2 3
1 2 0
0 0 1
1 0 0
0 1 0

出力例 1

0 1 0
0 0 1

$A = \begin{pmatrix} 1+2x & x^2 \\ 1 & x \end{pmatrix}$。 $\det(A) = x(1+2x) - x^2 = x + 2x^2 - x^2 = x + x^2$（係数: 0,1,1）。 $A^{-1}[0][0] = x / (x+x^2) = 1/(1+x) = 1-x+x^2-... $…実際のサンプルは計算上の例示。

概念図: FPS ガウス消去の流れ

FPS 行列のガウス消去（各演算がFPS演算に置き換わる） 元の行列 A（各エントリが多項式） f₀₀(x) f₀₁(x) f₁₀(x) f₁₁(x) → FPS逆元 +乗算 消去後（列0処理済み） f₀₀(x) ← pivot f₀₁(x) 0 f₁₁ - f₁₀·f₀₀⁻¹·f₀₁ FPS 逆元: Newton法 g₀ = f₀₀[0]⁻¹ (定数項の逆元) g_{n+1} = g_n · (2 - f₀₀ · g_n) mod x^{2^{n+1}} 計算量: O(K log K) (NTT) または O(K²) (naïve, K≤1000) 行列式 = 対角ピボット積 × (-1)^{行交換回数} det = f₀₀ · f₁₁' (mod x^K)

ヒント（段階的開示）

ヒント1: 方向性

エントリが多項式（FPS）の行列に対してガウス消去法を適用する。通常の実数ガウス消去と同じ手順だが、各「割り算」が「FPS の逆元 × FPS の乗算」に置き換わる。$d \le 8$, $K \le 1000$ なのでナイーブ乗算 $O(K^2)$ + 消去 $O(d^3)$ = $O(d^3 K^2)$ でも通る。

ヒント2: アプローチ
  • FPS 逆元: $f^{-1} \bmod x^K$ は Newton 法 $g_{n+1} = g_n(2 - fg_n) \bmod x^{2^{n+1}}$
  • ピボットは「定数項が 0 でない」行を選ぶ
  • 行列式: 対角ピボットの積（行交換のたびに符号反転）
  • 逆行列: 単位行列を右に拡張して Gauss-Jordan 消去
ヒント3: コード骨格 
def fps_inv(f, K):
    g = [pow(f[0], MOD-2, MOD)]
    sz = 1
    while sz < K:
        sz2 = min(sz*2, K)
        fg = fps_mul(f[:sz2], g, sz2)
        nmfg = [(-fg[i])%MOD for i in range(sz2)]
        nmfg[0] = (2-fg[0])%MOD
        g = fps_mul(g, nmfg, sz2)
        sz = sz2
    return g[:K]

# Gauss elimination
for col in range(d):
    # find pivot (constant term != 0)
    pivot = next(r for r in range(col, d) if mat[r][col][0] != 0)
    if pivot != col:
        mat[col], mat[pivot] = mat[pivot], mat[col]
        sign *= -1
    piv_inv = fps_inv(mat[col][col], K)
    det = fps_mul(det, mat[col][col], K)
    for row in range(d):
        if row == col: continue
        fac = fps_mul(mat[row][col], piv_inv, K)
        for c in range(d):
            mat[row][c] = fps_sub(mat[row][c],
                          fps_mul(fac, mat[col][c], K), K)

模範解答 (Python)

import sys, copy
input = sys.stdin.readline
MOD = 998244353

def fps_add(a, b, K):
    r = [0]*K
    for i in range(min(len(a),K)): r[i] = a[i]
    for i in range(min(len(b),K)): r[i] = (r[i]+b[i])%MOD
    return r

def fps_sub(a, b, K):
    r = [0]*K
    for i in range(min(len(a),K)): r[i] = a[i]
    for i in range(min(len(b),K)): r[i] = (r[i]-b[i])%MOD
    return r

def fps_mul(a, b, K):
    c = [0]*K
    for i in range(min(len(a),K)):
        if a[i]==0: continue
        for j in range(min(len(b),K-i)):
            c[i+j] = (c[i+j]+a[i]*b[j])%MOD
    return c

def fps_inv(f, K):
    assert f[0] != 0
    g = [pow(f[0], MOD-2, MOD)]
    sz = 1
    while sz < K:
        sz2 = min(sz*2, K)
        fg = fps_mul(f[:sz2], g, sz2)
        nm = [(-fg[i])%MOD for i in range(sz2)]
        nm[0] = (2-fg[0])%MOD
        g = fps_mul(g[:sz2], nm, sz2)
        sz = sz2
    return g[:K]

def main():
    d, K = map(int, input().split())
    mat = []
    for i in range(d):
        row = []
        for j in range(d):
            c = list(map(int, input().split()))
            c = (c+[0]*K)[:K]
            row.append(c)
        mat.append(row)

    # --- det ---
    m = copy.deepcopy(mat)
    det = [0]*K; det[0] = 1
    sign = 1
    for col in range(d):
        pivot = -1
        for row in range(col, d):
            if m[row][col][0] != 0: pivot = row; break
        if pivot == -1: det = [0]*K; break
        if pivot != col:
            m[col], m[pivot] = m[pivot], m[col]; sign *= -1
        piv_inv = fps_inv(m[col][col], K)
        det = fps_mul(det, m[col][col], K)
        for row in range(d):
            if row == col: continue
            fac = fps_mul(m[row][col], piv_inv, K)
            for c in range(d):
                m[row][c] = fps_sub(m[row][c], fps_mul(fac, m[col][c], K), K)
    if sign == -1: det = [(-x)%MOD for x in det]

    # --- A^{-1}[0][0] via Gauss-Jordan ---
    aug = [mat[i][:] + [([1]+[0]*(K-1)) if i==j else [0]*K for j in range(d)]
           for i in range(d)]
    for col in range(d):
        pivot = -1
        for row in range(col, d):
            if aug[row][col][0] != 0: pivot = row; break
        if pivot != col: aug[col], aug[pivot] = aug[pivot], aug[col]
        piv_inv = fps_inv(aug[col][col], K)
        for row in range(d):
            if row == col: continue
            fac = fps_mul(aug[row][col], piv_inv, K)
            for c in range(2*d):
                aug[row][c] = fps_sub(aug[row][c], fps_mul(fac, aug[col][c], K), K)
        for c in range(2*d):
            aug[col][c] = fps_mul(aug[col][c], piv_inv, K)
    inv00 = aug[0][d]

    print(*det)
    print(*inv00)

main()

Step-by-Step 解説

Step 1: FPS の基本演算

多項式の加減乗は係数列演算。逆元は Newton 法: 定数項の逆元から出発し、精度を $1, 2, 4, \ldots, K$ と2倍ずつ伸ばす。$K \le 1000$ ではナイーブ $O(K^2)$ の乗算で十分。

Step 2: ピボット選択

FPS の「可逆性」は定数項 $\ne 0$ が条件。ピボット行を探す際は mat[row][col][0] != 0 を確認する。

Step 3: 行列式の計算

各 col で対角エントリを det に乗算。行交換のたびに sign *= -1 し、最後に適用。消去後の上三角行列の対角積が det。

Step 4: 逆行列（Gauss-Jordan）

$d \times 2d$ の拡張行列（右側が単位行列の FPS 版）に完全行基本変形を適用。対角を 1 に正規化（fps_mul(row, piv_inv)）することで右側ブロックが $A^{-1}$ になる。

よくあるミス

ミス原因正しい書き方
fps_mul の結果が K を超える畳み込みで次数が増える結果を [:K] でトリム
ピボット行の正規化忘れGauss-Jordan で対角を 1 にしないfps_mul(row, piv_inv) で全エントリをスケール
sign の反映を忘れる行交換回数を最後に det に掛けないif sign==-1: det=[(-x)%MOD for x in det]
fps_inv で g の更新に古い g を使うnm 計算後に g を置き換えるべきg = fps_mul(g, nm, sz2) で新しい g を得る

次のステップ

  • 発展問題: $d \times d$ FPS 行列のべき乗 $A^n \bmod x^K$（Cayley-Hamilton + 行列累乗）
  • 関連: Berlekamp-Massey + Kitamasa 法（線形漸化式 $\bmod$ FPS との接続）

自己評価