println("方程组:")
println("x₁ + 2x₂ + 3x₃ = 1")
println("-x₁ + 3x₂ + 7x₃ = 4")
println("9x₁ + 3x₃ = 7")

# 方法1: 直接求解法
println("\n方法1: 直接矩阵求解")
x_direct = A \ b
println("解: x₁ = ", round(x_direct[1], digits=4), 
        ", x₂ = ", round(x_direct[2], digits=4),
        ", x₃ = ", round(x_direct[3], digits=4))

# 方法2: LU分解法
println("\n方法2: LU分解")
F = lu(A)
x_lu = F \ b
println("解: x₁ = ", round(x_lu[1], digits=4), 
        ", x₂ = ", round(x_lu[2], digits=4),
        ", x₃ = ", round(x_lu[3], digits=4))

# 验证解的正确性
println("\n验证解:")
residuals = A * x_direct - b
println("残差: ", round.(residuals, digits=8))
println("总残差平方和: ", round(sum(residuals.^2), digits=10))

# 创建优化模型
model = Model(GLPK.Optimizer)

# 定义决策变量
@variable(model, x >= 0)  # 产品A数量
@variable(model, y >= 0)  # 产品B数量

# 添加约束
@constraint(model, 3x + 2y <= 1000)  # 原料限制

# 定义目标函数
@objective(model, Max, 40x + 30y)  # 最大化利润

# 求解模型
optimize!(model)

# 打印结果
println("\n求解结果:")
println("状态: ", termination_status(model))
if termination_status(model) == MOI.OPTIMAL
    x_opt = value(x)
    y_opt = value(y)
    profit = objective_value(model)
    
    println("最优生产方案:")
    println("  产品A数量 (x) = ", round(x_opt, digits=4))
    println("  产品B数量 (y) = ", round(y_opt, digits=4))
    println("  最大利润 = ", round(profit, digits=4), " 元")
    
    # 分析可行域
    println("\n可行域分析:")
    println("  x ∈ [0, ", min(1000/3, 333.3333), "]")
    println("  y ∈ [0, ", 500, "]")
    
    # 利润线斜率
    slope = -(40/30)
    println("  利润线斜率: ", round(slope, digits=4))
    
    # 计算约束边界交点
    intercept_x = 1000 / 3  # x截距
    intercept_y = 500       # y截距
    
    # 角点分析
    corners = [(0, 0), (intercept_x, 0), (0, intercept_y), (x_opt, y_opt)]
    corner_profits = [40*c[1] + 30*c[2] for c in corners]
    
    println("\n角点利润比较:")
    println("  (0,0): 利润 = 0")
    println("  (333.333,0): 利润 = ", 40*intercept_x)
    println("  (0,500): 利润 = ", 30*intercept_y)
    println("  (", x_opt, ",", y_opt, "): 利润 = ", profit)
else
    println("未找到最优解!")
end

# 第二题：线性规划问题求解
solve_linear_programming()