如何在scipy.optimize中动态生成约束?

How to generate constraints dynamically in scipy.optimize?

好吧,我想做的是使用 scipy.optimize.minimize 对以下内容进行建模。

我要优化的是这个函数及其约束: 这里变量 V 是一个变量列表,列表的长度等于 Omega 的大小。 到目前为止我所做的如下:

import numpy as np
from scipy.linalg import norm
from scipy.optimize import minimize

# set of concepts
M = ['linear algebra','seq2seq', 'artificial neural network','pointer networks']

#subchapters
S1=['linear algebra', 'differential calculus','backpropagation']
S2=['linear algebra','seq2seq', 'artificial neural network']
S3=['linear algebra','seq2seq', 'artificial neural network','pointer networks']

#vector representation of the subchapter in the concept space
x=[[1,1,0,0],
   [1,1,1,0],
   [1,1,1,1]]

# set of prerequisite relations among subchapters (entered manually for testing)
Omega = [(1, 2),(2,3),(1,3)]

# number of concepts
m = len(M)

# define of theta and lambda constants (manual)
theta = 2
lamda = 1

# matrix A is a m x m matrix , where m is the number of concepts
# it represents the prerequisite relations among concepts
# A is generated randomly
np.random.seed(43)
#A = np.zeros((m,m), dtype = int)
A = np.random.randint(2, size=(m,m))

# define the slack variable V as an array of size equal to len(Omega)
V = np.empty(len(Omega), dtype=float)

bnds=[]
# bounds -1 and 1 , create the array
# -1 <= a[i][j] <= 1
bnds_a_s_t = [bnds.append((-1, 1)) for _ in range(np.size(A))]
# bounds for the slack variable V, V is positive
bnds_V_i_j = [bnds.append((0,np.inf)) for _ in range(np.size(V))]

#constraints
cons=[]
#equality constraint
#a[s][t] + a[t][s] = 0
def equality_constraint(X):
    A_no_flatten=X[:len(X)-len(Omega)]
    #reconstruct matrix A
    A=np.reshape(A_no_flatten,(m,m))

    for s in range(m):
        for t in range(m):
            r=A[s][t]+A[t][s]
            #r=0 constraint
            con = {'type': 'eq', 'fun': lambda X: r} 
            cons.append(con)

# inequality constraint
#x[i].T @ (C[i][j] * A) @ x[j]
def inequality_constraint(X,x):
    for couple in Omega:
        # get the i and j
        i = couple[0]
        j = couple[1]
        
        #initialize C to 1s
        C = np.full((m,m), 1, dtype = int)
        # take all elements from X except last len(Omega) elements
        A_no_flatten=X[:len(X)-len(Omega)]
        # reconstruct list V
        V=X[-len(Omega):]
        #index for V
        f=0
        #reconstruct matrix A
        A=np.reshape(A_no_flatten,(m,m))
        
        #construct C[i][j]
        for s in range(m):
            for t in range(m):
                if x[i][t]>0 or x[j][s]>0:
                    C[s][t]=0
                else:
                    C[s][t]=1
        
        first= x[i].T
        second = C*A
        third = x[j]
        
        first_sec = first@second
        res=first_sec@third
        ineq_con = {'type': 'ineq', 'fun': lambda X: res -theta +V[f]}
        f+=1
        cons.append(ineq_con)
        
# arguments passed to the function
# here we pass x matrix 
# arguments are passed and used in constraints and in the objective function
# the objective function will minimize A and V which are matrix A and slack variable V
arguments=(x,)

# objective function
def objective(X, x):
    A_no_flatten=X[:len(X)-len(Omega)]
    # reconstruct list V
    V=X[-len(Omega):]
    #reconstruct matrix A
    A=np.reshape(A_no_flatten,(m,m))
    
    # sum of square V
    sum_square=0.0
    for it in V:
        sum_square+=it**2
    # sum of square V * lambda
    sum_square_lambda=sum_square*lamda
    
    return norm(A, 1) + sum_square_lambda

# list of variables to pass to the objective function
#pass the x0.flatten() which is the A + V combined, and then when in objective function we recreate them
# the first one A is all except the last s items where s is the size of V
# and then V is the rest

B = A.flatten()
p0 = np.append(B,V)

# scipy minimize
sol = minimize(objective, x0 = p0, args=arguments, bounds=bnds, constraints=cons)
print(sol.x)

我得到的是:

[-7.73458566e-010  0.00000000e+000  4.99999997e-001  1.00000000e+000
  1.00000000e+000  0.00000000e+000 -5.00000003e-001  1.00000000e+000
  1.00000000e+000  1.00000000e+000  4.99999997e-001 -7.73458566e-010
 -7.73458566e-010  0.00000000e+000  4.99999997e-001 -7.73458566e-010
  6.01347002e-154  1.07176259e-311  0.00000000e+000]

这不符合约束条件,也不是我所期望的

我不知道这样添加约束是否正确,因为我似乎没有调用约束函数,我需要在循环中添加它们,每个函数都依赖于 X这是要最小化的列表。 当我打印 cons 数组时,它是空的,我知道,但是我没有找到另一种添加约束的方法 a[s][t]+a[t][s]=0 和另一种方法,我不知道我的方法是否正确,提前感谢您的帮助,非常感谢。

这不是一个完整的答案,但它可以帮助您入门。正如评论中已经提到的,当传递给 minimize 方法时,您的约束列表 cons 是空的。因此,让我们考虑您的第一个平等约束。有几个问题:

  • 每次调用函数 equality_constraint 时,您都会向列表添加新的约束条件 cons

  • 由于将每个约束 A[s][t] + A[t][s] == 0 作为标量函数传递,因此非常麻烦。相反,您可以使用单个向量值函数:

def constraint1(X):
    A_no_flatten = X[:len(X)-len(Omega)]
    A = np.reshape(A_no_flatten,(m,m))
    return A.flatten() + A.T.flatten()

顾名思义,.flatten() 方法将矩阵展平为向量,而 A.T 只是 A 的转置。现在您可以添加此约束:

cons = []
cons.append({'type': 'eq', 'fun': constraint1})

对另一个约束进行类似的操作。