cvxpy中构建矩阵逐列最大值位置二元指示变量的问题咨询

你现有实现的错误主要有三点：

• 逻辑冗余：额外引入maxesi、maxesii两个不必要的二元中间变量，最后关联y的约束y[:,bi] >= maxesi[:,bi] + maxesii[:,bi] - 1仅能约束maxesi和maxesii同时为1时y取1，无法覆盖所有场景，多余的变量和约束会压缩可行域甚至导致问题无解，进而影响其他决策变量的输出。
• 数值精度问题：设置的delta=1e-10远小于大部分整数规划求解器的默认浮点精度阈值（通常为1e-6），会被求解器判定为0，导致约束失效。
• 约束逻辑错误：四个大小判断约束存在重复、符号反向的问题，无法正确关联y的取值和D元素与列最大值的相等关系。

首先明确变量定义：

• D：m行n列的原始决策变量矩阵
• maxes：n维向量，每一位为D对应列的最大值，可通过cp.max(D, axis=0)直接得到
• y：m行n列的二元决策变量，即你需要的最大值位置标识变量

仅需要3组约束即可实现需求，M取值只需大于D的最大可能取值范围即可（比如你的示例中D最大为100，M取101就足够）：

1. 当y[i,j] = 1时，强制D[i,j]等于列最大值：D[i,j] >= maxes[j] - M*(1 - y[i,j])，由于maxes[j]本身是列最大值，天然满足D[i,j] <= maxes[j]，因此该约束会直接让D[i,j] = maxes[j]
2. 当y[i,j] = 0时，强制D[i,j]严格小于列最大值：maxes[j] - D[i,j] >= epsilon - M*y[i,j]，其中epsilon取1e-4即可，避免浮点精度问题
3. 每列至少有一个最大值：cp.sum(y[:,j]) >= 1，避免求解器将某列的y全部设为0

对应的cvxpy向量化实现代码如下，无需循环，效率更高：

import cvxpy as cp

m, n = D.shape
y = cp.Variable((m, n), boolean=True)
M = 101  # 可根据D的实际取值范围调整，大于D的最大可能值即可
epsilon = 1e-4

constraints = [
    # 约束1：y=1时D元素等于列最大值
    D >= maxes - M * (1 - y),
    # 约束2：y=0时D元素严格小于列最大值
    maxes - D >= epsilon - M * y,
    # 约束3：每列至少有一个最大值位置
    cp.sum(y, axis=0) >= 1
]

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