8.3 Ext 與 Tor：函子的計算#

import sympy as sp
from sympy import ZZ
from sympy.matrices.normalforms import smith_normal_form

def ab_group_from_matrix(D):
    """A = Z^m / im D 的結構：(自由秩, [不變因子])。"""
    D = sp.Matrix(D)
    m, k = D.shape
    S = smith_normal_form(D, domain=ZZ)
    invars = [int(S[i, i]) for i in range(min(m, k)) if S[i, i] > 1]
    rank = m - sum(1 for i in range(min(m, k)) if S[i, i] != 0)
    return rank, invars

def ext1_Zn_A(n, A_matrix):
    """Ext^1(Z/n, A)，A = Z^m / im(A_matrix)。
    公式：Ext^1(Z/n, A) = A / nA = coker([D | n·I])。"""
    D = sp.Matrix(A_matrix)
    m, k = D.shape
    I_n = n * sp.eye(m)
    M = D.row_join(I_n)                    # m x (k+m)
    return ab_group_from_matrix(M)

def tor1_Zn_A(n, A_matrix):
    """Tor_1(Z/n, A) = A[n] = {a ∈ A : n·a = 0}。
    在 Z^m/im D 上解 n·x ≡ 0 (mod im D)，即 coker 的 n-撓：
    A[n] ≅ ker(n: A -> A)，用 smith 標準形轉換基底後逐一計算。"""
    D = sp.Matrix(A_matrix)
    m, k = D.shape
    S = smith_normal_form(D, domain=ZZ)
    invars = [int(S[i, i]) for i in range(min(m, k)) if S[i, i] > 1]
    torsion = []
    for d in invars:                        # Z/d 的 n-撓 = Z/gcd(n,d)
        g = sp.gcd(n, d)
        if g > 1:
            torsion.append(int(g))
    torsion.sort()
    # 合併成整除鏈（循環直和的標準形）：gcd 們重新做不變因子合併
    return 0, torsion                       # 自由部分無 n-撓

def tor1_by_smith(n, A_matrix):
    """Tor_1(Z/n, A) 的另一算法：用解消 Z --n--> Z -> Z/n 張量 A 的核。
    等價地：ker( n : Z^m/im D -> Z^m/im D )，
    轉成增廣問題：ker([D | nI] : Z^{k+m} -> Z^m) 對 im(D) 的商——
    實務上直接用不變因子：Z/d_j 的 n-撓 = Z/gcd(n,d_j)。"""
    return tor1_Zn_A(n, A_matrix)

import sympy as sp
# A = Z (自由)：D = 零矩陣（無關係）—— A = Z^1，表現：Z -> Z^1 為 [0]
print("A = Z:")
print("  Ext^1(Z/6, Z) =", ext1_Zn_A(6, sp.Matrix([[0]])))   # 預期 Z/6
print("  Tor_1(Z/6, Z) =", tor1_Zn_A(6, sp.Matrix([[0]])))   # 預期 0

# A = Z/12（由 D = [12] 表現）
print("A = Z/12:")
print("  Ext^1(Z/8, Z/12) =", ext1_Zn_A(8, sp.Matrix([[12]])))  # 預期 Z/4
print("  Tor_1(Z/8, Z/12) =", tor1_Zn_A(8, sp.Matrix([[12]])))  # 預期 Z/4

# A = Z/2 ⊕ Z/6（不變因子 2, 6；表現 D = diag(2,6)）
print("A = Z/2 ⊕ Z/6:")
D = sp.diag(2, 6)
print("  Ext^1(Z/6, A) =", ext1_Zn_A(6, D))   # 預期 Z/2 ⊕ Z/6
print("  Tor_1(Z/6, A) =", tor1_Zn_A(6, D))   # 預期 Z/2 ⊕ Z/6

A = Z:
  Ext^1(Z/6, Z) = (0, [6])
  Tor_1(Z/6, Z) = (0, [])
A = Z/12:
  Ext^1(Z/8, Z/12) = (0, [4])
  Tor_1(Z/8, Z/12) = (0, [4])
A = Z/2 ⊕ Z/6:
  Ext^1(Z/6, A) = (0, [2, 6])
  Tor_1(Z/6, A) = (0, [2, 6])

import sympy as sp

def snf_tracked(A):
    """整數矩陣 A 的 Smith 標準形，含么模變換追蹤。
    回傳 (D, U, V)：D = U * A * V，U、V 皆么模（行列式 ±1）。"""
    if isinstance(A, sp.MatrixBase):
        m, n = A.shape                       # 保留形狀（0 列的矩陣也要記得欄數）
        A = A.tolist()
    else:
        m = len(A)
        n = len(A[0]) if m > 0 else 0
    a = [[int(x) for x in row] for row in A]
    A0 = sp.Matrix(m, n, lambda i, j: a[i][j])   # 原始矩陣（保留形狀）
    U = [[1 if i == j else 0 for j in range(m)] for i in range(m)]
    V = [[1 if i == j else 0 for j in range(n)] for i in range(n)]

    def rowop(E):                       # a <- E a, U <- E U
        nonlocal a, U
        a = [[sum(E[i][k]*a[k][j] for k in range(m)) for j in range(n)] for i in range(m)]
        U = [[sum(E[i][k]*U[k][j] for k in range(m)) for j in range(m)] for i in range(m)]

    def colop(F):                       # a <- a F, V <- V F
        nonlocal a, V
        a = [[sum(a[i][k]*F[k][j] for k in range(n)) for j in range(n)] for i in range(m)]
        V = [[sum(V[i][k]*F[k][j] for k in range(n)) for j in range(n)] for i in range(n)]

    def swapmat(n_, i, j):
        E = [[1 if (r == c and r != i and r != j) else 0 for c in range(n_)] for r in range(n_)]
        E[i][j] = E[j][i] = 1
        return E

    def rowswap(i, j):
        rowop(swapmat(m, i, j))

    def colswap(i, j):
        colop(swapmat(n, i, j))

    def rowcomb(i, j, q):               # row_i <- row_i + q*row_j
        E = sp.eye(m).tolist(); E[i][j] = q
        rowop(E)

    def colcomb(i, j, q):               # col_i <- col_i + q*col_j
        F = sp.eye(n).tolist(); F[j][i] = q
        colop(F)

    for k in range(min(m, n)):
        while True:
            piv = pr = pc = None
            for i in range(k, m):       # 找絕對值最小的非零元當軸心
                for j in range(k, n):
                    if a[i][j] != 0 and (piv is None or abs(a[i][j]) < abs(piv)):
                        piv, pr, pc = a[i][j], i, j
            if piv is None:
                break
            if pr != k: rowswap(k, pr)
            if pc != k: colswap(k, pc)
            piv = a[k][k]; ok = True
            for j in range(k+1, n):     # 消去第 k 列其餘元素
                if a[k][j] % piv != 0: ok = False
                colcomb(j, k, -(a[k][j] // piv))
            for i in range(k+1, m):     # 消去第 k 行其餘元素
                if a[i][k] % piv != 0: ok = False
                rowcomb(i, k, -(a[i][k] // piv))
            if ok:                      # 軸心整除整列整行，收工
                break
        if a[k][k] < 0:                 # 對角元取非負
            for j in range(n): a[k][j] = -a[k][j]
            for j in range(n): V[j][k] = -V[j][k]
    D = sp.Matrix(m, n, lambda i, j: a[i][j])
    Um, Vm = sp.Matrix(U), sp.Matrix(V)
    assert Um.det() in (1, -1) and Vm.det() in (1, -1)
    assert (Um * A0 * Vm - D).is_zero_matrix
    return D, Um, Vm

def kernel_basis(B):
    """B: Z^k -> Z^m 的整數核基底（回傳 k x d 矩陣，每一欄是一個基底向量）。"""
    B = sp.Matrix(B)
    m, k = B.shape
    D, U, V = snf_tracked(B)
    # D y = 0 的整數解由「零對角元所在的座標」生成；x = V y
    free_cols = [j for j in range(k) if j >= min(m, k) or D[j, j] == 0]
    return sp.Matrix([[V[i, j] for j in free_cols] for i in range(k)])

def coker_structure(A):
    """Z^m / im A 的結構：回傳 (自由秩, [不變因子])。"""
    A = sp.Matrix(A)
    m, n = A.shape
    D, U, V = snf_tracked(A)
    invars = [int(D[i, i]) for i in range(min(m, n)) if D[i, i] > 1]
    rank = sum(1 for i in range(min(m, n)) if D[i, i] != 0)
    return m - rank, invars

def homology(d_higher, d_lower):
    """H = ker d_lower / im d_higher（兩者的餘域都是 C_q = Z^m）。
    回傳 (自由秩 beta_q, [撓不變因子])。"""
    d_lower = sp.Matrix(d_lower)
    d_higher = sp.Matrix(d_higher)
    K = kernel_basis(d_lower)              # m x d，核基底（欄為基底向量）
    d = K.shape[1]
    if d == 0:
        return 0, []
    C = sp.zeros(d, d_higher.shape[1])     # 像(d_higher) 在核基底下的座標
    for j in range(d_higher.shape[1]):
        sol, _ = K.gauss_jordan_solve(d_higher[:, j])
        assert all(x.is_Integer for x in sol)   # 座標必為整數
        C[:, j] = sol
    return coker_structure(C)

import sympy as sp

# ---- 前置：本節「共用工具」的 snf_tracked / kernel_basis / coker_structure / homology ----

def ext_via_resolution(n, A_matrix):
    """按定義計算 Ext^1(Z/n, A)：
    解消 0 -> Z --n--> Z -> Z/n -> 0，
    取 Hom(-, A)：Hom(Z, A) = A --Hom(×n)--> Hom(Z, A) = A，
    Ext^1 = ker( Hom(d_1) ) / im( Hom(d_2) ) = ker(×n : A -> A) / 0 = A[n]？錯：
    注意 Hom 反變：Ext^1(Z/n, A) = coker( ×n : A -> A ) = A / nA。
    （ker(×n) 是 Tor_1 一側。）這裡示範用矩陣機器算 coker：
    A = Z^m / im D，×n 誘導 A -> A，其 coker = Z^m / (im D + n Z^m)。"""
    D = sp.Matrix(A_matrix)
    m = D.shape[0]
    M = D.row_join(n * sp.eye(m))
    return coker_structure(M)          # 直接給 (自由秩, 不變因子)

def tor_via_resolution(n, A_matrix):
    """按定義計算 Tor_1(Z/n, A)：
    張量解消得 0 -> Tor_1 -> A --×n--> A -> (Z/n)⊗A -> 0，
    Tor_1 = ker(×n : A -> A) = A[n]。
    A = Z^m/im D：類 x̄ ∈ A 滿足 n x̄ = 0 ⟺ n x ∈ im D ⟺ 存在 y 使 n x = D y，
    即 [D | -nI](x, y) = 0。取此方程的整數核基底，其 y 分量給出
    S = ker(×n : A -> A) 的格基底；再把 im D ∩ S 的生成元用該基底表出，
    商 Z^d / (im D ∩ S) 即 Tor_1。"""
    D = sp.Matrix(A_matrix)
    m, k = D.shape
    M = D.row_join(-n * sp.eye(m))     # [D | -nI] : Z^{k+m} -> Z^m
    K = kernel_basis(M)                # (k+m) x t，整數核基底
    Y = K[k:, :]                       # y 分量：S = {y ∈ Z^m : n y ∈ im D} 的格基底
    d = Y.shape[1]
    if d == 0:
        return 0, []
    C = sp.zeros(d, k)                 # im D 的生成元在 S 的基底下的座標（必為整數）
    for j in range(k):
        sol, _ = Y.gauss_jordan_solve(D[:, j])
        assert all(x.is_Integer for x in sol)
        C[:, j] = sol
    return coker_structure(C)          # Z^d / (im D ∩ S) = Tor_1

# 用 Z/12 驗證兩種途徑一致
D12 = sp.Matrix([[12]])
print("按定義 Ext:", ext_via_resolution(8, D12))   # 預期 (0, [4])
print("按定義 Tor:", tor_via_resolution(8, D12))   # 與 A[8] = Z/4 對照

按定義 Ext: (0, [4])
按定義 Tor: (0, [4])

import sympy as sp

# H_*(RP^2; Z) = (Z, Z/2, 0)；H_1 = Z/2
# Künneth: H_2(X × Y) ≅ (H_0⊗H_2) ⊕ (H_1⊗H_1) ⊕ (H_2⊗H_0) ⊕ Tor_1(H_1(X), H_1(Y))
print(tor1_Zn_A(2, sp.Matrix([[2]])))   # Tor_1(Z/2, Z/2) = Z/2
# → H_2(RP^2 × RP^2) = Z/2 ⊕ Z/2（一項來自 H_1⊗H_1，一項來自 Tor）

(0, [2])