import numpy as np
from scipy.optimize import linprog

# 定義問(wèn)題：最大化 3x + 5y，約束 x + 2y ≤ 100，x,y ≥ 0且為整數(shù)
def solve_integer_problem():
    # 1. 求解連續(xù)線性規(guī)劃（作為整數(shù)解的參考）
    res = linprog(
        c=[-3, -5],  # 目標(biāo)函數(shù)系數(shù)（取負(fù)實(shí)現(xiàn)最大化）
        A_ub=[[1, 2]],  # 約束系數(shù)矩陣
        b_ub=[100],     # 約束上限
        bounds=[(0, None), (0, None)],  # 變量非負(fù)約束
        method='highs'
    )
    
    if not res.success:
        return "求解失敗：無(wú)法找到可行解"
    
    # 2. 在連續(xù)解附近生成候選整數(shù)解
    x_cont, y_cont = res.x
    candidates = [(x, y) for x in range(max(0, int(x_cont)-2), int(x_cont)+3)
                        for y in range(max(0, int(y_cont)-2), int(y_cont)+3)
                        if x + 2*y <= 100]  # 過(guò)濾滿足約束的解
    
    # 3. 篩選最優(yōu)整數(shù)解
    best_x, best_y = max(candidates, key=lambda xy: 3*xy[0] + 5*xy[1])
    max_val = 3*best_x + 5*best_y
    
    return f"最優(yōu)整數(shù)解：\nx = {best_x}, y = {best_y}\n目標(biāo)函數(shù)最大值 = {max_val}"

# 執(zhí)行并打印結(jié)果
print(solve_integer_problem())

from gurobipy import Model, GRB

# 創(chuàng)建模型
model = Model("SupplyChain")

# 添加整數(shù)變量
x = model.addVar(vtype=GRB.INTEGER, name="x")
y = model.addVar(vtype=GRB.INTEGER, name="y")
# 添加約束條件
model.addConstr(x + 2*y <= 100, "ResourceConstraint")

# 設(shè)置目標(biāo)函數(shù)（最大化）
model.setObjective(3*x + 5*y, GRB.MAXIMIZE)



# 執(zhí)行優(yōu)化
model.optimize()

# 檢查優(yōu)化是否成功
if model.status == GRB.OPTIMAL:
    # 輸出變量的最優(yōu)值
    print(f"最優(yōu)解:")
    print(f"x = {x.x}")
    print(f"y = {y.x}")
    print(f"目標(biāo)函數(shù)值 = {model.objVal}")
else:
    print("未找到最優(yōu)解")

import pulp

# 創(chuàng)建問(wèn)題實(shí)例，指定為最大化問(wèn)題
model = pulp.LpProblem("SupplyChain", pulp.LpMaximize)

# 定義決策變量（非負(fù)整數(shù)）
x = pulp.LpVariable('x', lowBound=0, cat='Integer')  # lowBound=0表示x≥0
y = pulp.LpVariable('y', lowBound=0, cat='Integer')  # cat='Integer'指定為整數(shù)變量

# 設(shè)置目標(biāo)函數(shù)
model += 3 * x + 5 * y, "TotalProfit"  # 第二個(gè)參數(shù)是目標(biāo)函數(shù)的名稱(chēng)

# 添加約束條件
model += x + 2 * y <= 100, "ResourceConstraint"  # 約束名稱(chēng)為"ResourceConstraint"

# 求解模型
status = model.solve(pulp.PULP_CBC_CMD(msg=0))  # msg=0表示不顯示求解器日志

# 輸出求解狀態(tài)
print(f"求解狀態(tài): {pulp.LpStatus[status]}")

# 輸出最優(yōu)解
if pulp.LpStatus[status] == "Optimal":
    print(f"最優(yōu)解:")
    print(f"x = {x.varValue}")
    print(f"y = {y.varValue}")
    print(f"目標(biāo)函數(shù)值 = {pulp.value(model.objective)}")
else:
    print("未找到最優(yōu)解")

from pyomo.environ import ConcreteModel, Var, Objective, Constraint, Integers, maximize, SolverFactory

# 創(chuàng)建具體模型實(shí)例
model = ConcreteModel(name="SupplyChain")

# 定義決策變量（非負(fù)整數(shù)）
model.x = Var(domain=Integers, bounds=(0, None), name="x")  # bounds=(0, None)表示x≥0
model.y = Var(domain=Integers, bounds=(0, None), name="y")  # domain=Integers指定為整數(shù)變量

# 定義目標(biāo)函數(shù)（最大化）
def objective_rule(model):
    return 3 * model.x + 5 * model.y
model.objective = Objective(rule=objective_rule, sense=maximize, name="TotalProfit")

# 定義約束條件
def constraint_rule(model):
    return model.x + 2 * model.y <= 100
model.constraint = Constraint(rule=constraint_rule, name="ResourceConstraint")

# 注意/注意、選擇求解器并求解（這里使用CBC求解器，需要提前安裝）
solver = SolverFactory('cbc')  # 使用開(kāi)源的CBC求解器
results = solver.solve(model)

# 輸出求解狀態(tài)
print(f"求解狀態(tài): {results.solver.status}")
print(f"求解終止條件: {results.solver.termination_condition}")

# 輸出最優(yōu)解
if results.solver.termination_condition == 'optimal':
    print("\n最優(yōu)解:")
    print(f"x = {model.x.value}")
    print(f"y = {model.y.value}")
    print(f"目標(biāo)函數(shù)值 = {model.objective.value}")
else:
    print("未找到最優(yōu)解")

from ortools.linear_solver import pywraplp

# 創(chuàng)建求解器（使用GLOP線性規(guī)劃求解器，OR-Tools內(nèi)置）
solver = pywraplp.Solver.CreateSolver("SCIP")  # SCIP是OR-Tools中用于整數(shù)規(guī)劃的求解器

# 定義整數(shù)變量（非負(fù)）
x = solver.IntVar(0, solver.infinity(), "x")  # IntVar表示整數(shù)變量，第一個(gè)參數(shù)是下界
y = solver.IntVar(0, solver.infinity(), "y")

print(f"變量數(shù)量: {solver.NumVariables()}")

# 添加約束條件 x + 2y ≤ 100
constraint = solver.Constraint(-solver.infinity(), 100)  # 約束范圍(-∞, 100]
constraint.SetCoefficient(x, 1)  # x的系數(shù)為1
constraint.SetCoefficient(y, 2)  # y的系數(shù)為2

print(f"約束數(shù)量: {solver.NumConstraints()}")

# 定義目標(biāo)函數(shù)：最大化3x + 5y
objective = solver.Objective()
objective.SetCoefficient(x, 3)
objective.SetCoefficient(y, 5)
objective.SetMaximization()  # 設(shè)置為最大化問(wèn)題

# 求解問(wèn)題
status = solver.Solve()

# 輸出結(jié)果
if status == pywraplp.Solver.OPTIMAL:
    print("找到最優(yōu)解：")
    print(f"x = {x.solution_value()}")
    print(f"y = {y.solution_value()}")
    print(f"目標(biāo)函數(shù)最大值 = {objective.Value()}")
else:
    print("未找到最優(yōu)解")

# 輸出求解器信息
print(f"\n求解器使用時(shí)間: {solver.WallTime()} 毫秒")

from gekko import GEKKO

# 創(chuàng)建模型
m = GEKKO(remote=False)  # remote=False表示使用本地求解器

# 定義整數(shù)變量（非負(fù)）
x = m.Var(integer=True, lb=0, name='x')  # integer=True指定為整數(shù)變量，lb=0設(shè)置下界為0
y = m.Var(integer=True, lb=0, name='y')

# 設(shè)置目標(biāo)函數(shù)（最大化3x + 5y）
# Gekko默認(rèn)是最小化，所以通過(guò)取負(fù)數(shù)實(shí)現(xiàn)最大化
m.Maximize(3*x + 5*y)

# 添加約束條件
m.Equation(x + 2*y <= 100)

# 求解模型
m.solve(disp=False)  # disp=False不顯示詳細(xì)求解過(guò)程

# 輸出結(jié)果
print("最優(yōu)解:")
print(f"x = {x.value[0]}")
print(f"y = {y.value[0]}")
print(f"目標(biāo)函數(shù)最大值 = {m.options.OBJFCNVAL * -1}")  # 乘以-1還原最大值

from sko.GA import GA
import numpy as np
import matplotlib.pyplot as plt
# 目標(biāo)函數(shù)：最大化3x + 5y（轉(zhuǎn)為最小化問(wèn)題，返回負(fù)值）
def objective_func(x):
    return -(3 * x[0] + 5 * x[1])

# 約束條件：x + 2y ≤ 100（滿足約束時(shí)返回值≥0）
def constraint_func(x):
    return 100 - (x[0] + 2 * x[1])

# 創(chuàng)建遺傳算法實(shí)例
ga = GA(
    func=objective_func,          # 目標(biāo)函數(shù)
    n_dim=2,                      # 變量維度（x和y）
    size_pop=500,                 # 種群大小
    max_iter=2000,                 # 迭代次數(shù)
    lb=[0, 0],                    # 變量下界（x≥0, y≥0）
    ub=[100, 50],                 # 變量上界（根據(jù)約束估算）
    prob_mut=0.005,                #變異系數(shù)，默認(rèn)0.001
    constraint_eq=[constraint_func],  # 約束條件

)

# 執(zhí)行優(yōu)化
best_x, best_y = ga.run()

# 輸出結(jié)果
print("遺傳算法優(yōu)化結(jié)果：")
print(f"最優(yōu)解：x = {int(best_x[0])}, y = {int(best_x[1])}")
print(f"目標(biāo)函數(shù)最大值：{int(-best_y)}")  # 還原為正值

# 繪制迭代曲線

plt.rcParams['font.sans-serif'] = ['SimHei']  # 設(shè)置中文
plt.plot(ga.generation_best_Y)
plt.title("遺傳算法迭代過(guò)程")
plt.xlabel("迭代次數(shù)")
plt.ylabel("最優(yōu)目標(biāo)函數(shù)值（負(fù)值）")
plt.show()

import cvxpy as cp
from cvxpy.error import SolverError  # 導(dǎo)入特定異常類(lèi)

# 定義整數(shù)變量（正確指定變量類(lèi)型）
x = cp.Variable(integer=True, name="x")
y = cp.Variable(integer=True, name="y")

# 定義目標(biāo)函數(shù)：最大化3x + 5y
objective = cp.Maximize(3 * x + 5 * y)

# 定義約束條件：
constraints = [
    x + 2 * y <= 100,  # 資源約束
    x >= 0,             # 非負(fù)約束
    y >= 0              # 非負(fù)約束
]

# 創(chuàng)建優(yōu)化問(wèn)題
problem = cp.Problem(objective, constraints)

# 定義求解器嘗試順序
solvers = [
    ("CBC", cp.CBC),
    ("GLPK_MI", cp.GLPK_MI),
    ("ECOS_BB", cp.ECOS_BB),
    ("SCIP", cp.SCIP)
]

result = None
solver_used = None

# 嘗試多個(gè)求解器
for solver_name, solver in solvers:
    try:
        result = problem.solve(solver=solver, verbose=False)
        solver_used = solver_name
        # 檢查是否獲得最優(yōu)解
        if problem.status == cp.OPTIMAL:
            break
    except SolverError:
        print(f"求解器 {solver_name} 不可用，嘗試下一個(gè)...")
    except Exception as e:
        print(f"求解器 {solver_name} 出錯(cuò): {str(e)}")

# 輸出求解結(jié)果
print(f"\n使用的求解器: {solver_used or '無(wú)可用求解器'}")
print(f"求解狀態(tài): {problem.status}")

if problem.status == cp.OPTIMAL:
    # 確保獲取整數(shù)值（直接取整）
    x_val = int(round(x.value))
    y_val = int(round(y.value))
    
    print("\n最優(yōu)解:")
    print(f"x = {x_val}")
    print(f"y = {y_val}")
    print(f"目標(biāo)函數(shù)最大值 = {round(objective.value)}")  # 使用objective.value更準(zhǔn)確
elif problem.status == cp.INFEASIBLE:
    print("問(wèn)題不可行，無(wú)解")
elif problem.status == cp.UNBOUNDED:
    print("問(wèn)題無(wú)界，解無(wú)限大")
else:
    print("未找到最優(yōu)解，請(qǐng)檢查求解器安裝和問(wèn)題設(shè)置")
    print(f"最后嘗試的結(jié)果: x={x.value}, y={y.value}")