Pyomo资产组合优化代码故障求助：集合参数约束定义问题

我有一份资产收益率列表：

Re=[0.5346,0.5064,1.0838,0.7665,0.9463,0.7047,0.6735,0.5294,0.7697,0.7299,0.99,1.0856,0.9052,0.3827,0.3804,1.0271,0.9431,0.538,0.9313,0.9423]

需要最大化目标函数：
$$Re(w)=\sum_{i=1}^{n}w_{i}Re_{i}$$
约束条件：

• 资本全额利用：$\sum_{k=1}^{n}w_{k}v_{k}=1$（$v_k$为二元变量，持有资产k取1，否则0）
• 基数约束：$7 \le \sum_{k=1}^{n}v_{k} \le 10$
• 无卖空约束：$w_{k}\ge 0$（k=1,2,...,n）
• 单资产投资比例上下限：$0.01v_{k} \le w_{k} \le 0.3v_{k}$（k=1,2,...,n）

原代码错误分析 #

1. 目标函数定义错误：Objective的rule不需要接收i参数，原代码会导致绑定失败，因为目标是全局求和，不需要按索引生成多个目标。
2. 全局约束错误添加索引：资本全额利用、基数约束都是全局约束，不需要按资产索引i生成多个重复约束，原代码会生成20条相同约束，引发求解器错误。
3. 基数约束不符合需求：原代码写死等于7，实际需要的是7到10的范围约束。
4. 单资产比例约束错误：原代码返回生成器表达式，无法被Pyomo识别，且未结合二元变量v_k，不符合约束要求。
5. 冗余变量定义：model.q是多余的，直接用数值范围约束即可，不需要定义为决策变量。

修复后的代码 #

import pyomo.environ as pyo
from pyomo.opt import SolverFactory

# 创建模型
model = pyo.ConcreteModel()

# 资产集合
model.i = pyo.Set(initialize=[f'a{j+1}' for j in range(20)])

# 收益率参数
Re_values = [0.5346,0.5064,1.0838,0.7665,0.9463,0.7047,0.6735,0.5294,0.7697,0.7299,0.99,1.0856,0.9052,0.3827,0.3804,1.0271,0.9431,0.538,0.9313,0.9423]
model.Re = pyo.Param(model.i, initialize=dict(zip(model.i, Re_values)))

# 决策变量
model.w = pyo.Var(model.i, within=pyo.NonNegativeReals)  # 投资比例
model.v = pyo.Var(model.i, domain=pyo.Binary)  # 是否持有资产

# 目标函数：最大化组合收益率
def obj_rule(model):
    return sum(model.w[i] * model.Re[i] for i in model.i)
model.obj = pyo.Objective(rule=obj_rule, sense=pyo.maximize)

# 约束1：资本全额利用
def cap_full_rule(model):
    return sum(model.w[i] * model.v[i] for i in model.i) == 1
model.cap_full = pyo.Constraint(rule=cap_full_rule)

# 约束2：基数约束（7-10个资产）
def card_low_rule(model):
    return sum(model.v[i] for i in model.i) >= 7
model.card_low = pyo.Constraint(rule=card_low_rule)

def card_high_rule(model):
    return sum(model.v[i] for i in model.i) <= 10
model.card_high = pyo.Constraint(rule=card_high_rule)

# 约束3：单资产投资比例上下限
def w_lower_rule(model, i):
    return model.w[i] >= 0.01 * model.v[i]
model.w_lower = pyo.Constraint(model.i, rule=w_lower_rule)

def w_upper_rule(model, i):
    return model.w[i] <= 0.3 * model.v[i]
model.w_upper = pyo.Constraint(model.i, rule=w_upper_rule)

# 求解模型
solver = SolverFactory('cplex_direct')
results = solver.solve(model, tee=True)

# 输出结果
print("\n求解状态：", results.solver.status)
print("最优目标值：", pyo.value(model.obj))
print("\n选中的资产及投资比例：")
for i in model.i:
    if pyo.value(model.v[i]) > 0.9:
        print(f"资产{i}: 比例={pyo.value(model.w[i]):.4f}")

内容的提问来源于stack exchange，提问作者Reza