sample.column_generator のソースコード

# -*- coding: utf-8 -*-
from abc import ABCMeta, abstractmethod

from pysimple import Element, Parameter, Variable, Problem, NuoptStatus, Sum
from pysimple.typing import Table

class ColumnGenerator(metaclass=ABCMeta):
    """列生成法のための抽象基底クラス"""
    iternum = property(lambda self: self.__iternum)
    patternnum = property(lambda self: self.__pn0 + self.__iternum)

    @abstractmethod
    def create_init_pattern(self, *, pattern: dict, b: float | dict=1, c: float | dict=1) -> None:
        """set self.pattern, self.b, self.c"""
        self.pattern = pattern
        self.b = b
        self.c = c
        self.__iternum = 0
        self.__pn0 = len(pattern)

    def create_lambda(self, *, silent: bool=True) -> Table:
        """select_pattern の双対問題"""
        pi = Element(value=self.pattern)
        i = Element(set=pi(1).set)
        b = Parameter(index=pi(1), value=self.b)
        c = Parameter(index=pi(0), value=self.c)
        pattern = Parameter(index=pi, value=self.pattern)
        lmb = Variable(index=pi(1), lb=0, name='λ')
        prb = Problem(type=max)
        prb += Sum(pattern[pi]*lmb[pi(1)], pi(1)) <= c[pi(0)]
        prb += Sum(b[i]*lmb[i])
        prb.solve(silent=silent)
        assert prb.status == NuoptStatus.OPTIMAL
        return lmb.val

    @abstractmethod
    def create_new_pattern(self, lmbval: Table, *, silent: bool=True) -> tuple[Table, Table]:
        """組合せを 1 つ生成する"""
        # z[z[i].val>0].val, prb.objective.val

    def update_pattern(self, zval: Table) -> None:
        """update self.pattern, self.c"""
        self.pattern.update({(self.patternnum, i): v for (i,), v in zval.items()})

    # type は create_new_pattern の変数 type と同じ
    def select_pattern(self, *, vtype=bin, silent: bool=True) -> tuple[Table, Table]:  # 主問題
        """生成されたパターンから b を満たす組合せを選択する"""
        pi = Element(value=self.pattern)
        p = Element(set=pi(0).set)
        b = Parameter(index=pi(1), value=self.b)
        c = Parameter(index=pi(0), value=self.c)
        pattern = Parameter(index=pi, value=self.pattern)  # パターン p にアイテム i を含むか(0-1)
        z = Variable(index=p, type=vtype, lb=0)  # パターン p を選択するか(type=bin), 幾つ選択するか(type=int)
        prb = Problem()
        prb += Sum(pattern[pi]*z[pi(0)], pi(0)) >= b[pi(1)], 'cons'  # アイテム i は何れかのパターンで選ぶ
        prb += Sum(c[p]*z[p])
        prb.solve(silent=silent)
        assert prb.status == NuoptStatus.OPTIMAL
        return z[z[p].val>0].val, prb['cons'].dual

    def solve(self, *, eps: float=0.0, maxiter: int=100) -> Table:
        """列生成法で問題を解く"""
        self.create_init_pattern()  # 初期パターン
        objval = eps + 1e-6  # eps より大きい適当な値
        while objval > eps and self.iternum < maxiter:
            self.__iternum += 1

            # dual を求解
            lmbval = self.create_lambda()                 # 双対問題を使う
#            _, lmbval = self.select_pattern(type=float)  # 主問題の双対変数を使う

            # 新しいパターンの生成
            zval, objval = self.create_new_pattern(lmbval)
#            print('update:', dict(zval.items()))
            self.update_pattern(zval)
            print(f'{self.iternum}: {objval}', flush=True)

        # 実行可能解取得
        zval, _ = self.select_pattern()
        return zval