GLPK（python swiglpk）“问题没有原始可行解决方案”，但可以使用 CVXPY

Question

我正在尝试解决一个简单的优化问题：

max   x+y
s.t.  -x <= -1
       x,y in {0,1}^2

使用以下代码

import swiglpk
import numpy as np

def solve_boolean_lp_swig(obj: np.ndarray, aub: np.ndarray, bub: np.ndarray, minimize: bool) -> tuple:

    """
        Solves following optimization problem
        min/max     obj.dot(x)
        s.t         aub.dot(x) <= bub
                    x \in {0, 1}

        obj : m vector
        aub : nxm matrix
        bub : n vector
    """

    # init problem
    ia = swiglpk.intArray(1+aub.size); ja = swiglpk.intArray(1+aub.size)
    ar = swiglpk.doubleArray(1+aub.size)
    lp = swiglpk.glp_create_prob()

    # set obj to minimize if minimize==True else maximize
    swiglpk.glp_set_obj_dir(lp, swiglpk.GLP_MIN if minimize else swiglpk.GLP_MAX)

    # number of rows and columns as n, m 
    swiglpk.glp_add_rows(lp, int(aub.shape[0]))
    swiglpk.glp_add_cols(lp, int(aub.shape[1]))
    
    # setting row constraints (-inf < x <= bub[i])
    for i, v in enumerate(bub):
        swiglpk.glp_set_row_bnds(lp, i+1, swiglpk.GLP_UP, 0.0, float(v))
    
    # setting column constraints (x in {0, 1})
    for i in range(aub.shape[1]):
        # not sure if this is needed but perhaps for presolving
        swiglpk.glp_set_col_bnds(lp, i+1, swiglpk.GLP_FR, 0.0, 0.0) 
        # setting x in {0,1}
        swiglpk.glp_set_col_kind(lp, i+1, swiglpk.GLP_BV)

    # setting aub 
    for r, (i,j) in enumerate(np.argwhere(aub != 0)):
        ia[r+1] = int(i)+1; ja[r+1] = int(j)+1; ar[r+1] = float(aub[i,j])

    # solver settings
    iocp = swiglpk.glp_iocp()
    swiglpk.glp_init_iocp(iocp)
    iocp.msg_lev = swiglpk.GLP_MSG_ALL
    iocp.presolve = swiglpk.GLP_ON
    iocp.binarize = swiglpk.GLP_ON

    # setting objective
    for i,v in enumerate(obj):
        swiglpk.glp_set_obj_coef(lp, i+1, float(v))
    
    swiglpk.glp_load_matrix(lp, r, ia, ja, ar)
    info = swiglpk.glp_intopt(lp, iocp)
    
    # use later
    #status = swiglpk.glp_mip_status(lp)
    
    x = np.array([swiglpk.glp_mip_col_val(lp, int(i+1)) for i in range(obj.shape[0])])

    # for now, keep it simple. info == 0 means optimal 
    # solution (there are others telling feasible solution)
    return (info == 0), x

以及下面的实例（如上所示）

solve_boolean_lp_swig(
    obj = np.array([ 1, 1]),
    aub = np.array([[-1, 0]]),
    bub = np.array([-1]),
    minimize = False
)

在我看来，x=[1,0] 应该是一个有效的解决方案，因为 dot([-1, 0], x) <= -1（并且 [1,0] 是布尔值）成立，但求解器说 PROBLEM HAS NO PRIMAL FEASIBLE SOLUTION。但是，如果我使用 lib CVXOPT 运行相同的问题实例，而使用 cvxopt.glpk.ilp，求解器会找到最佳解决方案。我已经看到了 cvxopt 下面的 c 代码并且做了同样的事情，所以我怀疑有一些我看不到的小东西..

Answer 1

添加到模型中：

swiglpk.glp_write_lp(lp,None,"xxx.lp")

然后你会立即看到问题所在：

\* Problem: Unknown *\

Maximize
 obj: + z_1 + z_2

Subject To
 r_1: 0 z_1 <= -1

Bounds
 0 <= z_1 <= 1
 0 <= z_2 <= 1

Generals
 z_1
 z_2

End

我注意到r=0，所以加载调用的ne 参数已经是错误的。如果你设置r=1，事情看起来会更好。

Answer 2

约束

 x <= -1
 x,y in {0,1}^2

显然是不可行的。我怀疑你的代码没有反映模型。