GLPK (python swiglpk) "Problem has no primal feasible solution" 但可以使用 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] 是布尔值）成立但求解器说 问题没有主要可行的解决方案。但是，如果我运行使用库 CVXOPT 而不是 cvxopt.glpk.ilp 相同的问题实例，求解器会找到最佳解决方案。我已经看到 cvxopt 下面的 c 代码并且做了同样的事情所以我怀疑一些我看不到的小东西..

Answer 1

约束条件

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

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

Answer 2

添加到模型：

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 看起来会更好。

GLPK (python swiglpk) "Problem has no primal feasible solution" 但可以使用 CVXPY

GLPK (python swiglpk) "Problem has no primal feasible solution" but ok with CVXPY

python

glpk