Python – 单纯形返回不同的原始值和对偶值

Python – Simplex returning different values for the primal and the dual

我目前正在使用 Python 中的单纯形算法实现最优传输问题,用于原始和对偶。但是,我没有得到相同的原始值和对偶值。

我认为问题可能来自对偶问题,因为我也使用 Sinkhorn 的算法解决了原始问题,并且它返回了与基于单纯形的原始解决方案相同的值。

但是,我真的找不到问题所在,而且我 运行 不知道问题可能出在哪里。这是我的代码,希望大家清楚!

感谢您花时间阅读,提前感谢您可能提供的任何帮助!!!

import numpy as np
from random import random
import scipy
from scipy import optimize

#descriptions of the primal (resp dual) problem line 42 (resp 48 and 50)

n = 10 #number of points in distribution x_a
m = 20 #number of points in distribution x_b

x_a = [random() for i in range(n)] #n points chosen randomly between 0 and 1
x_a.sort() #sorted list
x_b = [2*random() for j in range(m)] #m points chosen randomly between 0 and 2
x_b.sort()

a = [1/n for i in range(n)] #distribution vector for x_a (uniform)
b = [1/m for i in range(m)] #distribution vector for x_b (uniform)
B = a+b #concatenation of a and b

#cost function (quadratic cost)
def cost(x,y) :
    n = len(x)
    p = len(y)
    tab = [[0 for j in range(p)] for i in range(n)]
    for i in range(n):
        for j in range(p):
            tab[i][j] = (x[i]-y[j])**2
    return tab

#definition of the matrix A
a0 = np.kron([1 for i in range(m)], np.eye(n, n))
a1 = np.kron(np.eye(m, m), [1 for i in range(n)])
A = np.concatenate((a0,a1), axis = 0)

#matrix of cost which we put columnwise into a 1D vector
cost_matrix = cost(x_a,x_b)
cost_column = []
for j in range(0,m):
    for i in range(n):
        cost_column += [cost_matrix[i][j]]

#Primal problem : Minimize cost_column*p under the constraint A*p=B (p decision variable)
proba = scipy.optimize.linprog(cost_column, A_ub=None, b_ub=None, A_eq=A, b_eq=B, bounds=None, method='simplex', callback=None, options=None, x0=None)
objective_result_primal = proba.fun
print(objective_result_primal)

A_transpose = np.transpose(A)
#Dual problem : Maximize B*h under the constraint A_transpose*h<=cost_column
B2 = [-B[i] for i in range(m+n)]#we use the opposite of B to turn the maximization problem into a minimization one (so that we can use the simplex)
#The dual problem becomes : Minimize B2*h under the constraint A_transpose*h<=cost_column
proba = scipy.optimize.linprog(B2, A_ub=A_transpose, b_ub=cost_column, A_eq=None, b_eq=None, bounds=None, method='simplex', callback=None, options=None, x0=None)
objective_result_dual = - proba.fun
print(objective_result_dual)

(等式的)运输问题的原始和对偶看起来像:

在对 linprog 的调用中使用 bounds=None 时,我们告诉求解器使用非负变量。这对于原始问题是正确的,但对于对偶问题则不然。对于对偶问题,您需要使用 bounds=(None,None) 表示变量应该是自由的。现在您应该看到原始问题和对偶问题的最佳 objective 值是相同的。