📐 Chuyên đề 1 — Hệ phương trình bậc nhất ba ẩn¶

import numpy as np

# Hệ PT từ Ví dụ 3 (SGK tr.8):
# x + y + z = 2
# 7x + 3y + z = 4
# -5x + 7y - 2z = 5

A = np.array([
    [1,  1,  1],
    [7,  3,  1],
    [-5, 7, -2]
], dtype=float)

b = np.array([2, 4, 5], dtype=float)

x = np.linalg.solve(A, b)
print(f'Nghiệm: x={x[0]:.4f}, y={x[1]:.4f}, z={x[2]:.4f}')

# Kiểm tra
print('Kiểm tra A @ x == b:', np.allclose(A @ x, b))

def gauss_elimination(A, b, verbose=True):
    """Phương pháp khử Gauss với in từng bước."""
    n = len(b)
    # Ma trận tăng cường [A|b]
    M = np.hstack([A.astype(float), b.reshape(-1,1).astype(float)])
    
    if verbose:
        print('=== Ma trận ban đầu [A|b] ===')
        print(np.array2string(M, precision=3, suppress_small=True))
    
    # Khử xuôi (forward elimination)
    for col in range(n):
        # Tìm pivot
        pivot_row = np.argmax(np.abs(M[col:, col])) + col
        M[[col, pivot_row]] = M[[pivot_row, col]]
        
        for row in range(col+1, n):
            if M[col, col] == 0:
                continue
            factor = M[row, col] / M[col, col]
            M[row] -= factor * M[col]
        
        if verbose:
            print(f'\n--- Sau khi khử cột {col+1} ---')
            print(np.array2string(M, precision=3, suppress_small=True))
    
    # Thế ngược (back substitution)
    x = np.zeros(n)
    for i in range(n-1, -1, -1):
        if abs(M[i, i]) < 1e-10:
            print('⚠️ Hệ vô nghiệm hoặc vô số nghiệm!')
            return None
        x[i] = (M[i, -1] - np.dot(M[i, i+1:n], x[i+1:n])) / M[i, i]
    
    return x

A = np.array([[1,1,1],[7,3,1],[-5,7,-2]])
b = np.array([2, 4, 5])
sol = gauss_elimination(A, b)
if sol is not None:
    print(f'\n✅ Nghiệm: x={sol[0]:.3f}, y={sol[1]:.3f}, z={sol[2]:.3f}')

# Bài 1.3a (SGK tr.14): Giải hệ phương trình
# 2x - y - z = 2
# x + y = 3
# x - y + z = 2

# TODO: Điền vào ma trận A và vector b
A_bt = np.array([
    [2, -1, -1],
    [1,  1,  0],
    [1, -1,  1]
], dtype=float)

b_bt = np.array([2, 3, 2], dtype=float)

# Giải và in kết quả
x_bt = np.linalg.solve(A_bt, b_bt)
print(f'Nghiệm bài 1.3a: x={x_bt[0]:.3f}, y={x_bt[1]:.3f}, z={x_bt[2]:.3f}')

# Bài 1.4 (SGK tr.14) - Bài toán lương
# Gọi x, y, z là lương quản lý kho, văn phòng, tài xế (triệu đồng)
# x + y = 164
# x + z = 156
# x - z = 8

A_luong = np.array([[1,1,0],[1,0,1],[1,0,-1]], dtype=float)
b_luong = np.array([164, 156, 8], dtype=float)
sol_luong = np.linalg.solve(A_luong, b_luong)
print(f'Lương quản lý kho: {sol_luong[0]:.0f} triệu đồng/năm')
print(f'Lương quản lý văn phòng: {sol_luong[1]:.0f} triệu đồng/năm')
print(f'Lương tài xế: {sol_luong[2]:.0f} triệu đồng/năm')

# Hệ PT: x + y + z = 16100
#        20x - 3y - 3z = 0
#        10x + 5y + z = 87200

A_rung = np.array([[1,1,1],[20,-3,-3],[10,5,1]], dtype=float)
b_rung = np.array([16100, 0, 87200], dtype=float)
sol_rung = np.linalg.solve(A_rung, b_rung)
print(f'Cây bần:  {sol_rung[0]:.0f} cây — sinh khối {sol_rung[0]*10/1000:.1f} tấn/ha')
print(f'Cây đước: {sol_rung[1]:.0f} cây — sinh khối {sol_rung[1]*5/1000:.2f} tấn/ha')
print(f'Cây mắm:  {sol_rung[2]:.0f} cây — sinh khối {sol_rung[2]*1/1000:.3f} tấn/ha')