from scipy.optimize import linprog

# Coefficients of the objective function (costs)
c = [25000, 28000, 12000]

# Coefficients for the carbon percentage constraints (left side)
A_ub = [
    [-75, -72, -65],  # 0.73 <= (x_A * 75 + x_B * 72 + x_C * 65) / (x_A + x_B + x_C)
    [75, 72, 65],     # (x_A * 75 + x_B * 72 + x_C * 65) / (x_A + x_B + x_C) <= 0.76
]

# Right-hand side for the carbon percentage constraints
b_ub = [0.73, 0.76]

# Coefficients for the sulfur percentage constraints (left side)
A_ub_sulfur = [
    [-23, -26, -30],  # 0.23 <= (x_A * 23 + x_B * 26 + x_C * 30) / (x_A + x_B + x_C)
    [23, 26, 30],     # (x_A * 23 + x_B * 26 + x_C * 30) / (x_A + x_B + x_C) <= 0.25
]

# Right-hand side for the sulfur percentage constraints
b_ub_sulfur = [0.23, 0.25]

# Combine the two constraints
A_ub_total = A_ub + A_ub_sulfur
b_ub_total = b_ub + b_ub_sulfur

# Boundaries for the variables (non-negativity)
x_bounds = [(0, None), (0, None), (0, None)]

# Solving the problem using linear programming
result = linprog(c, A_ub=A_ub_total, b_ub=b_ub_total, bounds=x_bounds, method='highs')

# Display the result
if result.success:
    print(f"Optimal solution:\nCoal A: {result.x[0]} tons\nCoal B: {result.x[1]} tons\nCoal C: {result.x[2]} tons")
else:
    print("Optimization failed.")
