如何在Pyomo中更高效实现含多索引集合的优化模型？

我需要在Pyomo中实现如下优化模型（最小示例）：

由于使用的数据中IJK、JKL、KLM集合的元素数量庞大，需要高效实现满足索引映射（(i,j,k)属于IJK、(j,k,l)属于JKL且(k,l,m)属于KLM）的x变量求和逻辑。随着集合基数增大，遍历所有集合元素的方式效率低下。

当前实现代码 #

模型定义代码 #

import pyomo.environ as pyo

def model(I, IJK, JKL, KLM):
    model = pyo.ConcreteModel()

    model.I = pyo.Set(initialize=I)
    model.IJK = pyo.Set(initialize=IJK)
    model.JKL = pyo.Set(initialize=JKL)
    model.KLM = pyo.Set(initialize=KLM)

    model.z = pyo.Param(default=1)

    x_list = [(i, j, k, l, m) for (i, j, k) in model.IJK
              for (jj, kk, l) in model.JKL if (jj == j) and (kk == k)
              for (kkk, ll, m) in model.KLM if (kkk == k) and (ll == l)]

    model.x_list = pyo.Set(initialize=x_list)

    model.x = pyo.Var(model.x_list, domain=pyo.NonNegativeReals)

    model.OBJ = pyo.Objective(expr=model.z)

    model.ei = pyo.Constraint(model.I, rule=ei_rule)


def ei_rule(model, i):
    lhs = [model.x[k] for k in model.x_list if k[0] == i]
    if len(lhs) == 0:
        return pyo.Constraint.Skip
    else:
        return sum(lhs) >= 0

随机数据生成代码 #

import numpy as np
import pandas as pd

def create_random_data(n, m):
    I = [f'i{x}' for x in range(1, n + 1)]
    J = [f'j{x}' for x in range(1, m + 1)]
    K = [f'k{x}' for x in range(1, m + 1)]
    L = [f'l{x}' for x in range(1, m + 1)]
    M = [f'm{x}' for x in range(1, m + 1)]

    ijk = pd.DataFrame(np.random.binomial(1, 0.05, size=(len(I)*len(J)*len(K))),
                       index=pd.MultiIndex.from_product(
                           [I, J, K], names=['i', 'j', 'k']),
                       columns=['value']).reset_index()
    jkl = pd.DataFrame(np.random.binomial(1, 0.05, size=(len(J)*len(K)*len(L))),
                       index=pd.MultiIndex.from_product(
                           [J, K, L], names=['j', 'k', 'l']),
                       columns=['value']).reset_index()
    klm = pd.DataFrame(np.random.binomial(1, 0.05, size=(len(K)*len(L)*len(M))),
                       index=pd.MultiIndex.from_product(
                           [K, L, M], names=['k', 'l', 'm']),
                       columns=['value']).reset_index()

    IJK = [tuple(x) for x in ijk.loc[ijk['value'] == 1]
           [['i', 'j', 'k']].to_dict('split')['data']]
    JKL = [tuple(x) for x in jkl.loc[jkl['value'] == 1]
           [['j', 'k', 'l']].to_dict('split')['data']]
    KLM = [tuple(x) for x in klm.loc[klm['value'] == 1]
           [['k', 'l', 'm']].to_dict('split')['data']]
    return I, J, K, L, M, IJK, JKL, KLM

测试代码及结果 #

I, J, K, L, M, IJK, JKL, KLM = create_random_data(300, 20) 
%timeit -r 7 -n 10 model(I, IJK, JKL, KLM)

输出结果：582 ms ± 65 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

问题 #

请问是否存在更高效的Pyomo实现方式？

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