# -*- coding: utf-8 -*-
from typing import Any
from pysimple import Set, Element, Parameter, Variable, IntegerVariable, BinaryVariable, Problem, Options, Condition, Sum, Selection, Printf
class data2010:
r = {(1, 'Lead'): 0.2, (1, 'Zinc'): 0.3, (1, 'Tin'): 0.5, (2, 'Lead'): 0.5, (2, 'Zinc'): 0.4, (2, 'Tin'): 0.1, (3, 'Lead'): 0.3, (3, 'Zinc'): 0.2, (3, 'Tin'): 0.5, (4, 'Lead'): 0.3, (4, 'Zinc'): 0.4, (4, 'Tin'): 0.3, (5, 'Lead'): 0.3, (5, 'Zinc'): 0.3, (5, 'Tin'): 0.4, (6, 'Lead'): 0.6, (6, 'Zinc'): 0.3, (6, 'Tin'): 0.1, (7, 'Lead'): 0.4, (7, 'Zinc'): 0.5, (7, 'Tin'): 0.1, (8, 'Lead'): 0.1, (8, 'Zinc'): 0.3, (8, 'Tin'): 0.6, (9, 'Lead'): 0.1, (9, 'Zinc'): 0.1, (9, 'Tin'): 0.8}
c = {1: 7.3, 2: 6.9, 3: 7.3, 4: 7.5, 5: 7.6, 6: 6.0, 7: 5.8, 8: 4.3, 9: 4.1}
b = {'Lead': 0.3, 'Zinc': 0.3, 'Tin': 0.4}
[ドキュメント]def p2010_mixture(**kwds: Any) -> None:
"""配合問題"""
value = data2010
problem = Problem(name='配合問題', type=min)
# 添字
i = Element(value=value.c.keys()) # 市販の合金
j = Element(value=value.b.keys()) # 構成金属
# パラメータ
r = Parameter(index=(i,j), value=value.r, name='構成比率')
c = Parameter(index=i, value=value.c, name='コスト')
b = Parameter(index=j, value=value.b, name='目標比率')
# 変数
x = Variable(index=i, lb=0, name='混合比率')
print(r)
print(c)
print(b)
# 目的関数
problem += Sum(c[i]*x[i])
# 混合比の制約
problem += Sum(x[i]) == 1
problem += Sum(r[i,j]*x[i], i) == b[j]
# 求解
print(problem)
problem.solve(**kwds)
print('errorCode: ', problem.result.errorCode)
print('optValue: ', problem.result.optValue)
# 出力
print(x.val)
class data2020:
upper = {1: 250, 2: 450}
demand = {'a': 200, 'b': 200, 'c': 200}
cost = {(1, 'a'): 3.4, (1, 'b'): 2.2, (1, 'c'): 2.9, (2, 'a'): 3.4, (2, 'b'): 2.4, (2, 'c'): 2.5}
[ドキュメント]def p2020_transport2(**kwds: Any) -> None:
"""輸送問題"""
value = data2020
problem = Problem(name='輸送問題', type=min)
# 添字
i = Element(value=value.upper.keys()) # 工場
j = Element(value=value.demand.keys()) # 店舗
# パラメータ
upper = Parameter(index=i, value=value.upper, name='供給可能量')
demand = Parameter(index=j, value=value.demand, name='需要量')
cost = Parameter(index=(i,j), value=value.cost, name='輸送コスト')
# 変数
x = Variable(index=(i,j), lb=0, name='輸送量')
print(upper)
print(demand)
print(cost)
# 目的関数
problem += Sum(cost[i,j]*x[i,j]), '総コスト'
problem += Sum(x[i,j], j) <= upper[i], '工場の生産量の制約'
problem += Sum(x[i,j], i) == demand[j], '店舗の需要量の制約'
# 求解
print(problem)
problem.solve(**kwds)
# 出力
print(problem.objective.val)
print(x.val)
class data2030:
use = {('I', 'A'): 2, ('I', 'B'): 5, ('II', 'A'): 7, ('II', 'B'): 3}
out = {('I', 1): 30, ('I', 2): 60, ('I', 3): 80, ('II', 1): 20, ('II', 2): 50, ('II', 3): 90}
upper = {('A', 1): 920, ('A', 2): 750, ('A', 3): 500, ('B', 1): 790, ('B', 2): 600, ('B', 3): 480}
costp = {'I': 75, 'II': 50}
costi = {'I': 8, 'II': 7}
[ドキュメント]def p2030_multiplan2(**kwds: Any) -> None:
"""多期間計画問題"""
value = data2030
problem = Problem(name='多期間計画問題', type=min)
# 集合と添字
i = Element(value=value.costp.keys()) # 製品
jt = Element(value=value.upper.keys())
j = Element(set=jt(0).set) # 原料
t = Element(set=jt(1).set) # 期間
# パラメータ
use = Parameter(index=(i,j), value=value.use, name='原料使用量')
out = Parameter(index=(i,t), value=value.out, name='出荷量')
upper = Parameter(index=(j,t), value=value.upper, name='利用可能量')
costp = Parameter(index=i, value=value.costp, name='生産コスト')
costi = Parameter(index=i, value=value.costi, name='在庫コスト')
# 変数
x = Variable(index=(i,t), lb=0, name='生産量')
y = Variable(index=(i,t), lb=0, name='在庫量')
print(use)
print(out)
print(upper)
print(costp)
print(costi)
# 目的関数
problem += Sum(costp[i]*x[i,t] + costi[i]*y[i,t]), '総コスト'
problem += Sum(use[i,j]*x[i,t], i) <= upper[j,t], '原料の制約'
# 出荷量の制約
t1, t2, t3 = t==1, t==2, t==3
problem += x[i,t1] - y[i,t1] == out[i,t1]
problem += x[i,t2] + y[i,t2-1] - y[i,t2] == out[i,t2]
problem += x[i,t3] + y[i,t3-1] == out[i,t3]
# 求解
print(problem)
problem.solve(**kwds)
# 出力
print(problem.objective.val)
print(x.val)
print(y.val)
class data2040:
inD = {(1, 'capacity'): 5, (1, 'worktime'): 24, (2, 'capacity'): 10, (2, 'worktime'): 12, (3, 'capacity'): 20, (3, 'worktime'): 12, (4, 'capacity'): 20, (4, 'worktime'): 24, (5, 'capacity'): 30, (5, 'worktime'): 12, (6, 'capacity'): 50, (6, 'worktime'): 12}
outD = {(1, 'earnings'): 2, (2, 'earnings'): 6, (3, 'earnings'): 10, (4, 'earnings'): 12, (5, 'earnings'): 12, (6, 'earnings'): 20}
[ドキュメント]def p2040_DEA(**kwds: Any) -> None:
"""包絡分析法(DEA)モデル"""
value = data2040
problem = Problem(name='包絡分析法(DEA)モデル', type=max)
# 添字
ki = Element(value=value.inD.keys())
kj = Element(value=value.outD.keys())
i = Element(set=ki(1).set) # 入出力項目
j = Element(set=kj(1).set) # 出力項目
k = Element(set=ki(0).set) # 店舗
kbar = Element(value=[k.set[0]]) # 対象店舗
# 入力データ
inD = Parameter(index=(k,i), value=value.inD, name='入力データ')
# 出力データ
outD = Parameter(index=(k,j), value=value.outD, name='出力データ')
# 変数
# 入力項目に対する重み
inW = Variable(index=i, lb=0)
# 出力
outW = Variable(index=j, lb=0)
print(inD)
print(outD)
# 目的関数
problem += Sum(outD[kbar,j]*outW[j])
# 制約式
problem += Sum(inD[kbar,i]*inW[i], i) == 1
problem += Sum(outD[k,j]*outW[j], j) <= Sum(inD[k,i]*inW[i], i)
problem.options.method = Options.Method.SIMPLEX # 単体法の利用
# 求解
print(problem)
problem.solve(**kwds)
# 結果出力
print(problem.objective.val)
print(inW.val)
print(outW.val)
class data2050:
capacity = 65
value = {('coffee',): 120, ('PET bottles',): 130, ('banana',): 80, ('apple',): 100, ('rice ball',): 250, ('bread',): 185}
size = {('coffee',): 10, ('PET bottles',): 12, ('banana',): 7, ('apple',): 9, ('rice ball',): 21, ('bread',): 16}
[ドキュメント]def p2050_knapsack2(**kwds: Any) -> None:
"""ナップサック問題"""
data = data2050
problem = Problem(name='ナップサック問題', type=max)
# 添字
i = Element(value=data.value.keys()) # 品物
# パラメータを宣言する
capacity = Parameter(value=data.capacity, name='ナップサックの容量')
value = Parameter(index=i, value=data.value, name='品物の価値')
size = Parameter(index=i, value=data.size, name='品物のサイズ')
# 整数変数を宣言する
quantity = IntegerVariable(index=i, lb=0, name='詰め込む個数')
print(capacity)
print(value)
print(size)
# 総価値の最大化
problem += Sum(value[i]*quantity[i]), '総価値'
# 容量に関する制約
problem += Sum(size[i]*quantity[i]) <= capacity
# 求解し結果を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(quantity.val)
class data2060:
cost = {'佐藤': 200, '鈴木': 280, '高橋': 175, '田中': 560, '渡辺': 205, '伊藤': 245, '山本': 80, '中村': 195, '小林': 265, '斉藤': 190}
deliver = dict.fromkeys([
('佐藤','A'), ('佐藤','B'), ('佐藤','C'),
('鈴木','A'), ('鈴木','D'), ('鈴木','F'),
('高橋','B'), ('高橋','E'),
('田中','C'), ('田中','D'), ('田中','E'), ('田中','F'), ('田中','G'),
('渡辺','A'), ('渡辺','F'),
('伊藤','B'), ('伊藤','D'), ('伊藤','F'),
('山本','D'),
('中村','C'), ('中村','G'),
('小林','C'), ('小林','F'), ('小林','G'),
('斉藤','B'), ('斉藤','E'), ('斉藤','G'),
], 1)
[ドキュメント]def p2060_cover2(**kwds: Any) -> None:
"""集合被覆問題"""
value = data2060
problem = Problem(name='集合被覆問題', type=min)
# 添字
ij = Element(value=value.deliver.keys())
i = Element(value=value.cost.keys()) # 候補者
j = Element(value=sorted(ij(1).set)) # エリア
# パラメータおよび変数の宣言
deliver = Parameter(index=(i,j), value=value.deliver, name='配達可能エリア')
cost = Parameter(index=i, value=value.cost, name='コスト')
x = BinaryVariable(index=i, name='採用')
print(deliver[deliver[i,j] == 1])
print(cost)
# 総コストを最小化する
problem += Sum(cost[i]*x[i]), '総コスト'
# 各エリアに配置する
problem += Sum(deliver[i,j]*x[i], i) >= 1
# 求解し結果を表示する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(x.val)
class data2070:
city = 'ABCDXY'
upper = {('A', 'C'): 5, ('B', 'C'): 4, ('B', 'D'): 3, ('C', 'Y'): 8, ('D', 'Y'): 5, ('X', 'A'): 4, ('X', 'B'): 6}
[ドキュメント]def p2070_maxflow2(**kwds: Any) -> None:
"""最大流問題"""
value = data2070
problem = Problem(name='最大流問題', type=max)
# 添字
City = Set(value=value.city)
i = Element(set=City)
j = Element(set=City)
k = Element(set=City)
# 2点間の輸送量の上限をパラメータとして宣言
upper = Parameter(index=(i,j), value=value.upper, name='輸送量の上限')
print(upper[upper[i,j]>0])
# 変数の宣言
total = Variable(name='製粉工場からレストランへの総輸送量')
f = Variable(index=(i,j), lb=0, ub=upper[i,j], name='輸送量')
# 総輸送量の最大化
problem += total
# 各地点での輸送量
problem += total == Sum(f['X',j])
#kXY = (k!='X') & (k!='Y') # こちらでも良い
kXY = k>['X', 'Y']
problem += Sum(f[kXY,j], j) == Sum(f[j,kXY], j)
problem += Sum(f[i,'Y']) == total
# 輸送量に関する上下限制約(lb, ub で設定済み)
# problem += 0 <= f[i,j]
# problem += f[i,j] <= upper[i,j]
# 求解して解を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(f[f[i,j].val>0].val)
print(total.val)
class data2071:
city = 'ABCDXY'
incidence = {('X','e1'): 1, ('A','e1'): -1, ('X','e2'): 1, ('B','e2'): -1, ('A','e3'): 1, ('C','e3'): -1, ('B','e4'): 1, ('C','e4'): -1, ('B','e5'): 1, ('D','e5'): -1, ('C','e6'): 1, ('Y','e6'): -1, ('D','e7'): 1, ('Y','e7'): -1, ('Y','flow'): 1, ('X','flow'): -1}
upper = {'e1': 4, 'e2': 6, 'e3': 5, 'e4': 4, 'e5': 3, 'e6': 8, 'e7': 5}
upper['flow'] = sum(upper.values())
[ドキュメント]def p2071_maxflow3(**kwds: Any) -> None:
"""最大流問題(Matrix)"""
value = data2071
problem = Problem(name='最大流問題(Matrix)', type=max)
# 添字
v = Element(value=value.city) # 頂点
e = Element(value=value.upper.keys()) # 辺
# パラメータおよび変数の宣言
incidence = Parameter(index=(v,e), value=value.incidence) # 接続行列
upper = Parameter(index=e, value=value.upper) # 辺 e の輸送量上限
x = Variable(index=e, lb=0, ub=upper) # 辺 e の輸送量
# 総輸送量の最大化
problem += x['flow'], '総輸送量'
# 輸送量保存則
problem += Sum(incidence[v,e]*x[e], e) == 0, '輸送量保存則'
# 求解して解を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(x.val)
class data2072:
city = 'ABCDXY'
incidence = {('X','e1'): 1, ('A','e1'): -1, ('X','e2'): 1, ('B','e2'): -1, ('A','e3'): 1, ('C','e3'): -1, ('B','e4'): 1, ('C','e4'): -1, ('B','e5'): 1, ('D','e5'): -1, ('C','e6'): 1, ('Y','e6'): -1, ('D','e7'): 1, ('Y','e7'): -1, ('Y','flow'): 1, ('X','flow'): -1}
upper = {'e1': 4, 'e2': 6, 'e3': 5, 'e4': 4, 'e5': 3, 'e6': 8, 'e7': 5}
upper['flow'] = sum(upper.values())
[ドキュメント]def p2072_mincut(**kwds: Any) -> None:
"""最小カット問題"""
value = data2072
problem = Problem(name='最小カット問題')
# 添字
v = Element(value=value.city) # 頂点
e = Element(value=value.upper.keys()) # 辺
# パラメータおよび変数の宣言
incidence = Parameter(index=(v,e), value=value.incidence) # 接続行列
upper = Parameter(index=e, value=value.upper) # 辺 e の輸送量上限
a = Parameter(index=e, value={'flow': 1})
y = Variable(index=v)
z = Variable(index=e, lb=0)
# カット総量の最小化
problem += Sum(upper[e]*z[e]), 'カット総量'
problem += Sum(incidence[v,e]*y[v], v) + z[e] >= a[e], '双対制約式'
# 求解して解を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(y.val)
print(z.val)
class data2080:
city = 'ABCDXY'
upper = {('A', 'C'): 5, ('B', 'C'): 4, ('B', 'D'): 3, ('C', 'Y'): 8, ('D', 'Y'): 5, ('X', 'A'): 4, ('X', 'B'): 6}
cost = {('A', 'C'): 270, ('B', 'C'): 300, ('B', 'D'): 220, ('C', 'Y'): 190, ('D', 'Y'): 170, ('X', 'A'): 250, ('X', 'B'): 200}
supply = {'X': 8, 'Y': -8}
[ドキュメント]def p2080_mincost2(**kwds: Any) -> None:
"""最小費用流問題"""
value = data2080
problem = Problem(name='最小費用流問題', type=min)
# 添字
City = Set(value=value.city)
i = Element(set=City)
j = Element(set=City)
k = Element(set=City)
# パラメータの宣言
upper = Parameter(index=(i,j), value=value.upper, name='輸送量の上限')
cost = Parameter(index=(i,j), value=value.cost, name='輸送費用')
supply = Parameter(index=k, value=value.supply, name='流入量と流出量の差')
print(upper[upper[i,j]>0])
print(cost[cost[i,k]>0])
print(supply[supply[k]!=0])
# 変数の宣言
f = Variable(index=(i,j), lb=0, ub=upper[i,j], name='輸送量')
# 総費用の最小化
problem += Sum(cost[i,j]*f[i,j]), '総費用'
# 各地点での輸送量
problem += Sum(f[k,j], j)-Sum(f[i,k], i) == supply[k]
# 輸送量に関する上下限制約(lb, ub で設定済み)
# problem += 0 <= f[i,j]
# problem += f[i,j] <= upper[i,j]
# 求解して解を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(f[f[i,j].val>0].val)
class data2090:
city = 'ABCXYZ'
upper = {('A', 'C'): 18, ('A', 'Y'): 3, ('B', 'A'): 3, ('B', 'C'): 10, ('C', 'Y'): 10, ('C', 'Z'): 15, ('X', 'A'): 15, ('X', 'B'): 12}
PWheat_cost = {('A', 'C'): 4, ('A', 'Y'): 11, ('B', 'A'): 6, ('B', 'C'): 11, ('C', 'Y'): 9, ('X', 'A'): 12, ('X', 'B'): 10}
UWheat_cost = {('A', 'C'): 5, ('B', 'A'): 7, ('B', 'C'): 9, ('C', 'Z'): 5, ('X', 'A'): 11, ('X', 'B'): 12}
PWheat_supply = {'X': 8, 'Y': -8}
UWheat_supply = {'X': 13, 'Z': -13}
[ドキュメント]def p2090_multiflow2(**kwds: Any) -> None:
"""多品種流問題"""
value = data2090
problem = Problem(name='多品種流問題', type=min)
# 添字
City = Set(value=value.city)
i = Element(set=City)
j = Element(set=City)
k = Element(set=City)
# パラメータの宣言
upper = Parameter(index=(i,j), value=value.upper) # 輸送量の上限
PWheat_cost = Parameter(index=(i,j), value=value.PWheat_cost) # パン用の小麦粉の輸送費用
UWheat_cost = Parameter(index=(i,j), value=value.UWheat_cost) # うどん用の小麦粉の輸送費用
PWheat_supply = Parameter(index=k, value=value.PWheat_supply)
UWheat_supply = Parameter(index=k, value=value.UWheat_supply)
# 変数の宣言
PWheat = Variable(index=(i,j), lb=0) # パン用の小麦粉の輸送数
UWheat = Variable(index=(i,j), lb=0) # うどん用の小麦粉の輸送数
print(upper[upper[i,j]>0])
print(PWheat_cost[PWheat_cost[i,j]>0])
print(UWheat_cost[UWheat_cost[i,j]>0])
print(PWheat_supply[PWheat_supply[k]!=0])
print(UWheat_supply[UWheat_supply[k]!=0])
# 目的関数の宣言
problem += Sum(PWheat_cost[i,j]*PWheat[i,j]+UWheat_cost[i,j]*UWheat[i,j]) # 総輸送費
# 輸送ルートごとの輸送量の上限
problem += PWheat[i,j]+UWheat[i,j] <= upper[i,j]
# パン用の小麦粉に関する保存則
problem += Sum(PWheat[k,j], j)-Sum(PWheat[i,k], i) == PWheat_supply[k]
# うどん用の小麦粉に関する保存則
problem += Sum(UWheat[k,j], j)-Sum(UWheat[i,k], i) == UWheat_supply[k]
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(PWheat[PWheat[i,j].val>0].val)
print(UWheat[UWheat[i,j].val>0].val)
class data2100:
distance = {('A', 'B'): 7, ('A', 'C'): 2, ('A', 'D'): 13, ('A', 'E'): 19, ('A', 'F'): 17,
('B', 'A'): 7, ('B', 'C'): 8, ('B', 'D'): 6, ('B', 'E'): 12, ('B', 'F'): 10,
('C', 'A'): 2, ('C', 'B'): 8, ('C', 'D'): 14, ('C', 'E'): 20, ('C', 'F'): 18,
('D', 'A'): 13, ('D', 'B'): 6, ('D', 'C'): 14, ('D', 'E'): 7, ('D', 'F'): 4,
('E', 'A'): 19, ('E', 'B'): 12, ('E', 'C'): 20, ('E', 'D'): 7, ('E', 'F'): 3,
('F', 'A'): 17, ('F', 'B'): 10, ('F', 'C'): 18, ('F', 'D'): 4, ('F', 'E'): 3}
population = {'A': 800, 'B': 550, 'C': 780, 'D': 600, 'E': 1020, 'F': 360}
class data2110:
distance = {('A', 'B'): 7, ('A', 'C'): 2, ('A', 'D'): 13, ('A', 'E'): 19, ('A', 'F'): 17,
('B', 'A'): 7, ('B', 'C'): 8, ('B', 'D'): 6, ('B', 'E'): 12, ('B', 'F'): 10,
('C', 'A'): 2, ('C', 'B'): 8, ('C', 'D'): 14, ('C', 'E'): 20, ('C', 'F'): 18,
('D', 'A'): 13, ('D', 'B'): 6, ('D', 'C'): 14, ('D', 'E'): 7, ('D', 'F'): 4,
('E', 'A'): 19, ('E', 'B'): 12, ('E', 'C'): 20, ('E', 'D'): 7, ('E', 'F'): 3,
('F', 'A'): 17, ('F', 'B'): 10, ('F', 'C'): 18, ('F', 'D'): 4, ('F', 'E'): 3}
population = {'A': 800, 'B': 550, 'C': 780, 'D': 600, 'E': 1020, 'F': 360}
[ドキュメント]def p2110_center2(**kwds: Any) -> None:
"""pセンター問題"""
value = data2110
problem = Problem(name='pセンター問題', type=min)
# 添字
City = Set(value=value.population.keys()) # 都市
i = Element(set=City)
j = Element(set=City)
# パラメータの宣言
distance = Parameter(index=(i,j), value=value.distance, name='2都市間の距離')
population = Parameter(index=i, value=value.population, name='都市の人口')
# 変数の宣言
x = BinaryVariable(index=(i,j))
v = Variable()
print(distance)
print(population)
# 目的関数の宣言
problem += v
# 制約条件
problem += v >= population[i]*distance[i,j]*x[i,j]
problem += Sum(x[i,i]) == 2
problem += Sum(x[i,j], j) == 1
ij = i!=j # ij = Condition((i,j), i!=j) に同じ
problem += x[ij] <= x[ij(1,1)]
# 求解し出店する都市を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
Printf('{} 市に出店する.', x[i,i].val==1)
class data2120:
dis = {('A', 'B'): 6, ('A', 'C'): 5, ('A', 'D'): 5,
('B', 'A'): 6, ('B', 'C'): 7, ('B', 'D'): 4,
('C', 'A'): 5, ('C', 'B'): 7, ('C', 'D'): 3,
('D', 'A'): 5, ('D', 'B'): 4, ('D', 'C'): 3}
[ドキュメント]def p2120_TSP2(**kwds: Any) -> None:
"""巡回セールスマン問題"""
value = data2120
problem = Problem(name='巡回セールスマン問題', type=min)
# 添字
city = Element(value=value.dis.keys()) # 都市
c1 = Element(set=city(0).set)
c2 = Element(set=city(0).set)
c_ = c1 != c1.set[0] # 出発地点を除く都市
c1_ = Element(set=c_.set)
c2_ = Element(set=c_.set)
# パラメータ
dis = Parameter(index=(c1,c2), value=value.dis, name='距離')
num = len(c_.set)
# 変数
x = BinaryVariable(index=(c1,c2), name='移動')
y = Variable(index=c_, lb=1, ub=num, name='順番')
print(dis)
print(num)
# 目的関数
c12 = c1!=c2 # c12 = Condition((c1,c2), c1!=c2) に同じ
problem += Sum(x[c12]*dis[c12]), '経路長'
# 各都市から1つの別都市へ移動する
problem += Sum(x[c12], c12(1)) == 1
# 各都市に1つの別都市から移動する
problem += Sum(x[c12], c12(0)) == 1
# サブツアーを排除する
c12_ = c1_!=c2_
problem += y[c12_(0)] - y[c12_(1)] + num*x[c12_] <= num-1
# 順番の上下限制約(lb, ub で設定済み)
# problem += 1 <= y[c_]
# problem += y[c_] <= num
# 求解
print(problem)
problem.solve(**kwds)
# 出力
print(problem.objective.val)
print(x.val)
print(y.val)
[ドキュメント]def p2132_fieldassign2(method: Options.Method=Options.Method.AUTO, **kwds: Any) -> None:
"""基礎的なマス埋め割当問題"""
problem = Problem(name='基礎的なマス埋め割当問題', type=min)
I = Set(value=range(1, 5))
J = Set(value=range(1, 5))
i = Element(set=I)
j = Element(set=J)
x = BinaryVariable(index=(i,j))
problem += Sum(x[i,j], i) == 2
problem += Sum(x[i,j], j) == 2
# 以下は出力用のプログラム
print(problem)
problem.options.method = method
problem.solve(**kwds)
print(x[x[i,j].val>0].val)
Printf('{}: '+','.join(['{:.0f}']*len(J)), i, *(x[i,j_].val for j_, in J))
[ドキュメント]def p2132_fieldassign3(method: Options.Method=Options.Method.WCSP, **kwds: Any) -> None:
"""基礎的なマス埋め割当問題(wcsp)"""
p2132_fieldassign2(method=method, **kwds)
class data2133:
JPPair = [('接客','佐藤'), ('接客','鈴木'), ('接客','山本'), ('接客','渡辺'),
('厨房','安藤'), ('厨房','鈴木'), ('厨房','山本'), ('厨房','渡辺'),
('レジ','安藤'), ('レジ','佐藤'), ('レジ','山本'), ('レジ','渡辺'),
('発注','安藤'), ('発注','佐藤'), ('発注','鈴木'), ('発注','渡辺'),
('ごみ捨て','安藤'), ('ごみ捨て','佐藤'), ('ごみ捨て','鈴木'), ('ごみ捨て','山本'),
('買出','安藤'), ('買出','佐藤'), ('買出','鈴木'), ('買出','山本'), ('買出','渡辺'),
('掃除','安藤'), ('掃除','佐藤'), ('掃除','鈴木'), ('掃除','山本'), ('掃除','渡辺'),
('仕込','安藤'), ('仕込','佐藤'), ('仕込','鈴木'), ('仕込','山本'), ('仕込','渡辺')]
cost = dict(zip(JPPair, [1400, 520, 1410, 450, 1800, 1700, 1050, 2300, 800, 1500, 1500, 600, 500, 600, 1000, 1200, 800, 600, 1200, 800, 600, 600, 800, 700, 1300, 1200, 500, 1200, 500, 1300, 1500, 1000, 1200, 1200, 1500]))
jyukuren = dict(zip(JPPair, [3, -2, 3, -4, 5, 3, -4, 5, 0, 3, 3, -1, -3, -1, 1, 2, 1, -1, 3, 1, -1, -1, 1, 0, 2, 2, -2, 2, -3, 4, 5, -2, 0, 1, 5]))
[ドキュメント]def p2133_jobassign3(method: Options.Method=Options.Method.AUTO, **kwds: Any) -> None:
"""仕事割当問題"""
value = data2133
problem = Problem(name='仕事割当問題', type=min)
jp = Element(value=value.JPPair)
cost = Parameter(index=jp, value=value.cost, name='コスト')
jyukuren = Parameter(index=jp, value=value.jyukuren, name='熟練度')
x = BinaryVariable(index=jp)
print(cost)
print(jyukuren)
# 接客,厨房,レジ,掃除,仕込 は 2 人を割り当てる
jp2 = Condition(jp, jp(0)<'接客 厨房 レジ 掃除 仕込'.split())
problem += Sum(x[jp2], jp2(1)) == 2
# 発注,ごみ捨て,買出 は 1 人を割り当てる
jp1 = Condition(jp, jp(0)<'発注 ごみ捨て 買出'.split())
# problem += Sum(x[jp1], jp1(1)) == 1
problem += Selection(x[jp1], jp1(1))
# 各人に割り振る仕事は,最大で 3 つまでとする
problem += Sum(x[jp], jp(0)) <= 3
# 各仕事のクオリティーは負になってはいけない
problem += Sum(x[jp]*jyukuren[jp], jp(1)) >= 0
# 接客,厨房 は別の人が担当する
ep = Condition(jp, jp(0)<'接客 厨房'.split())
problem += Sum(x[ep], ep(0)) <= 1
problem += Sum(cost[jp]*x[jp]), '総コスト'
# 以下出力用プログラム
print(problem)
problem.options.method = method
problem.options.maxTime = 1
problem.solve(**kwds)
print(problem.objective.val)
print(x[x[jp].val>0].val)
[ドキュメント]def p2133_jobassign4(method: Options.Method=Options.Method.WCSP, **kwds: Any) -> None:
"""仕事割当問題(wcsp)"""
p2133_jobassign3(method=method, **kwds)
class data2140:
factory = [1, 2, 3, 4, 5]
location = ['I', 'II', 'III', 'IV', 'V']
F = {(1, 1): 0, (1, 2): 15, (1, 3): 12, (1, 4): 8, (1, 5): 7, (2, 1): 15, (2, 2): 0, (2, 3): 14, (2, 4): 2, (2, 5): 4, (3, 1): 12, (3, 2): 14, (3, 3): 0, (3, 4): 6, (3, 5): 6, (4, 1): 8, (4, 2): 2, (4, 3): 6, (4, 4): 0, (4, 5): 11, (5, 1): 7, (5, 2): 4, (5, 3): 6, (5, 4): 11, (5, 5): 0}
D = {('I', 'I'): 0, ('I', 'II'): 4, ('I', 'III'): 7, ('I', 'IV'): 6, ('I', 'V'): 10, ('II', 'I'): 4, ('II', 'II'): 0, ('II', 'III'): 2, ('II', 'IV'): 5, ('II', 'V'): 4, ('III', 'I'): 7, ('III', 'II'): 2, ('III', 'III'): 0, ('III', 'IV'): 8, ('III', 'V'): 7, ('IV', 'I'): 6, ('IV', 'II'): 5, ('IV', 'III'): 8, ('IV', 'IV'): 0, ('IV', 'V'): 12, ('V', 'I'): 10, ('V', 'II'): 4, ('V', 'III'): 7, ('V', 'IV'): 12, ('V', 'V'): 0}
[ドキュメント]def p2140_QAP(**kwds: Any) -> None:
"""二次割当問題"""
value = data2140
problem = Problem(name='二次割当問題')
# 添字
Factory = Set(value=value.factory) # 工場の集合
Location = Set(value=value.location) # 地区の集合
i = Element(set=Factory)
j = Element(set=Factory)
m = Element(set=Location)
n = Element(set=Location)
# パラメータおよび変数の宣言
F = Parameter(index=(i,j), value=value.F) # 各工場間のフロー
D = Parameter(index=(m,n), value=value.D) # 各地区間の距離
# 工場 i を地区 m に配置する場合 1 それ以外は 0
x = BinaryVariable(index=(i,m))
problem += Sum(F[i,j]*x[i,m]*D[m,n]*x[j,n])
problem += Sum(x[i,m], m) == 1
problem += Sum(x[i,m], i) == 1
# 求解して解を出力する
print(problem)
problem.options.method = Options.Method.WCSP
problem.options.maxTime = 10
problem.solve(**kwds)
print(problem.objective.val)
print(x[x[i,m].val>0].val)
class data2150:
costS = {1: 10, 2: 20, 3: 30}
costP = {1: 100, 2: 120, 3: 150, 4: 160, 5: 280, 6: 300}
p = {1: 5, 2: 10, 3: 20}
q = {1: 5, 2: 6, 3: 8, 4: 10, 5: 15, 6: 20}
D = {1: 10, 2: 10, 3: 12, 4: 12, 5: 15, 6: 20, 7: 30, 8: 30, 9: 35, 10: 40, 11: 50, 12: 56, 13: 52, 14: 47, 15: 40, 16: 35, 17: 42, 18: 50, 19: 48, 20: 40, 21: 30, 22: 20, 23: 15, 24: 12}
[ドキュメント]def p2150_FPP(**kwds: Any) -> None:
"""設備計画問題"""
value = data2150
problem = Problem(name='設備計画問題', type=min)
i = Element(value=value.p.keys())
t = Element(value=value.D.keys())
j = Element(value=value.q.keys())
costS = Parameter(index=i, value=value.costS) # 各商品に対する単位時間あたりのコスト
costP = Parameter(index=j, value=value.costP) # 各発電機に対する導入コスト
p = Parameter(index=i, value=value.p) # 各商品の単位時間あたりの出力電力
q = Parameter(index=j, value=value.q) # 各発電機の単位時間あたりの出力電力
D = Parameter(index=t, value=value.D) # 時刻ごとの電力需要
print(costS)
print(costP)
print(p)
print(q)
print(D)
u = BinaryVariable(index=(i,t)) # 商品.各時間ごとに,どの商品を利用するか
x = BinaryVariable(index=j) # 発電機.どの発電機を導入するか
problem += Sum(costS[i]*u[i,t]) + Sum(costP[j]*x[j]), '総コスト'
# 各時間における商品の使用は 1 つのみ
problem += Sum(u[i,t], i) == 1
# すべての時刻において,購入したスポット電力と自家発電の出力電力の和が電力需要を上回らないといけない
problem += D[t] - Sum(p[i]*u[i,t], i) <= Sum(q[j]*x[j])
print(problem)
problem.solve(**kwds)
# 結果出力
print(problem.objective.val)
print(x[x[j].val>0].val)
print(u[u[i,t].val>0].val)
class data2160:
y = {0.5: 1.22470, 1.0: 0.87334, 1.5: 0.99577, 2.0: 1.34215, 2.5: 2.12172, 3.0: 3.22933, 3.5: 4.44744, 4.0: 6.13509, 4.5: 8.14697, 5.0: 10.50759}
[ドキュメント]def p2160_leastsquare2(**kwds: Any) -> None:
"""最小二乗問題"""
value = data2160
problem = Problem(name='最小二乗問題', type=min)
# 添字
t = Element(value=value.y.keys())
# パラメータ
y = Parameter(index=t, value=value.y) # 観測値
# 変数
a = Variable() # 2 次の項の係数
b = Variable() # 1 次の項の係数
c = Variable() # 定数項
# 誤差を Expression を用いて表現する
e = a*t*t + b*t + c - y[t]
# 目的関数
problem += Sum(e*e, t), '各観測時刻における誤差の二乗和'
# 求解
print(problem)
problem.solve(**kwds)
# 結果出力
print(problem.objective.val)
print(a.val)
print(b.val)
print(c.val)
class data2170:
r = {(1, 1): 0.0000, (1, 2): 0.100, (2, 1): 0.0162, (2, 2): 0.161, (3, 1): 0.0307, (3, 2): 0.197, (4, 1): 0.0419, (4, 2): 0.192, (5, 1): 0.0485, (5, 2): 0.148, (6, 1): 0.0498, (6, 2): 0.084, (7, 1): 0.0458, (7, 2): 0.026, (8, 1): 0.0368, (8, 2): 0.000, (9, 1): 0.0238, (9, 2): 0.016, (10, 1): 0.0082, (10, 2): 0.068, (11, 1): -0.0082, (11, 2): 0.132, (12, 1): -0.0238, (12, 2): 0.184, (13, 1): -0.0368, (13, 2): 0.200, (14, 1): -0.0458, (14, 2): 0.174, (15, 1): -0.0498, (15, 2): 0.116, (16, 1): -0.0485, (16, 2): 0.052, (17, 1): -0.0419, (17, 2): 0.008, (18, 1): -0.0307, (18, 2): 0.003, (19, 1): -0.0162, (19, 2): 0.039, (20, 1): 0.0000, (20, 2): 0.100}
[ドキュメント]def p2170_portfolio1(**kwds: Any) -> None:
"""ポートフォリオ最適化問題"""
value = data2170
problem = Problem(name='ポートフォリオ最適化問題', type=min)
# 添字
j = Element(value=[j for _, j in value.r]) # 銘柄
t = Element(value=[t for t, _ in value.r]) # サンプル
tt = Element(set=t.set) # エラー回避用
# パラメータ
r = Parameter(index=(t,j), value=value.r) # 収益率
rb = Sum(r[tt,j], tt) / len(tt.set) # 平均収益率
# 変数
x = Variable(index=j, lb=0) # 組入比率
# 式
rpb = Sum(rb*x[j], j) # 期待収益率
dev = Sum(r[t,j]*x[j], j) - rpb # 偏差
# 目的関数
problem += Sum(dev*dev, t) / len(t.set), 'リスク'
# 制約条件
problem += Sum(x[j], j) == 1.0, '組入比率の総和は1'
problem += rpb >= 0.0, 'ポートフォリオが与える期待収益率は0以上'
# 求解
print(problem)
problem.solve(**kwds)
# 結果出力
print(problem.objective.val)
print(rpb.val)
print(x.val)
class data2230:
N = [1, 2, 3, 4, 5, 6]
G = {(1, 1): 0, (1, 2): 2, (1, 3): 3, (1, 4): 0, (1, 5): 1, (1, 6): 3, (2, 1): 2, (2, 2): 0, (2, 3): 1, (2, 4): 2, (2, 5): 0, (2, 6): 0, (3, 1): 3, (3, 2): 1, (3, 3): 0, (3, 4): 2, (3, 5): 0, (3, 6): 3, (4, 1): 0, (4, 2): 2, (4, 3): 2, (4, 4): 0, (4, 5): 1, (4, 6): 4, (5, 1): 1, (5, 2): 0, (5, 3): 0, (5, 4): 1, (5, 5): 0, (5, 6): 2, (6, 1): 3, (6, 2): 0, (6, 3): 3, (6, 4): 4, (6, 5): 2, (6, 6): 0}
[ドキュメント]def p2230_maxcut2(**kwds: Any) -> None:
"""隣接行列(最大カット問題)"""
value = data2230
problem = Problem(name='隣接行列(最大カット問題)', type=max)
# 添字
i = Element(value=value.N)
j = Element(value=value.N)
ij = i<j
# パラメータおよび変数の宣言
G = Parameter(index=(i,j), value=value.G) # 隣接行列
# 0-1 の 2 つのグループのどちらに入るか
x = BinaryVariable(index=i)
# カットされたエッジの情報を x[i], x[j] から得るための変数
y = BinaryVariable(index=ij)
problem += Sum(G[ij]*y[ij]), 'カットされるエッジの重みの総和'
# x[i] x[j] y[ij] の整合性をとる制約
problem += -x[ij(0)] - x[ij(1)] + y[ij] <= 0, '(x[i] x[j] y[ij]) = (0 0 1) の排除'
problem += -x[ij(0)] + x[ij(1)] - y[ij] <= 0, '(x[i] x[j] y[ij]) = (0 1 0) の排除'
problem += x[ij(0)] - x[ij(1)] - y[ij] <= 0, '(x[i] x[j] y[ij]) = (1 0 0) の排除'
problem += x[ij(0)] + x[ij(1)] + y[ij] <= 2, '(x[i] x[j] y[ij]) = (1 1 1) の排除'
# 求解して解を出力する
print(problem)
problem.solve(**kwds)
print(problem.objective.val)
print(x[x[i].val>0].val)
print(y[y[ij].val>0].val)
class data2260:
A = {(1, 1): 1, (1, 2): -1, (1, 3): 1, (2, 1): 1, (2, 2): 1, (2, 3): -1, (3, 1): -1, (3, 2): 1, (3, 3): 1, (4, 1): 1, (4, 2): 1, (4, 3): 1}
[ドキュメント]def p2260_pseudoinverse(**kwds: Any) -> None:
"""ムーア・ペンローズ一般逆行列"""
value = data2260
problem = Problem(name='ムーア・ペンローズ一般逆行列')
# 添字
M = Set(value=[m for m, _ in value.A]) # 行
N = Set(value=[n for _, n in value.A]) # 列
m = Element(set=M)
i = Element(set=N)
j = Element(set=N)
k = Element(set=N)
# パラメータおよび変数の宣言
A = Parameter(index=(m,i), value=value.A) # 定数行列成分
ATA = Sum(A[m,i]*A[m,j], m) # trans(A) * A
B = ATA**2
x = Variable(index=(i,j)) # 変数行列成分
problem += B*x[i,j]*B == B
# 求解して解を出力する
print(problem)
problem.solve(**kwds)
print(x[x[i,j].val>0].val)
mm = Element(set=m.set)
A_plus = Sum(ATA[i,j]*x[j,k].val*A[mm,k], (j,k))
Printf(','.join(['{}']*len(m.set)), *(A_plus[i,m_] for m_, in m.set))
if __name__ == '__main__':
from sys import argv
from . import reidaishu
if argv[1:]:
for arg in argv[1:]:
getattr(reidaishu, arg)()
else:
p2010_mixture() # 配合問題
p2020_transport2() # 輸送問題
p2030_multiplan2() # 多期間計画問題
p2040_DEA() # 包絡分析法(DEA)モデル
p2050_knapsack2() # ナップサック問題
p2060_cover2() # 集合被覆問題
p2070_maxflow2() # 最大流問題
p2071_maxflow3() # 最大流問題(Matrix)
p2072_mincut() # 最小カット問題
p2080_mincost2() # 最小費用流問題
p2090_multiflow2() # 多品種流問題
p2100_median2() # pメディアン問題
p2110_center2() # pセンター問題
p2120_TSP2() # 巡回セールスマン問題
p2132_fieldassign2() # 基礎的なマス埋め割当問題
p2132_fieldassign3() # 基礎的なマス埋め割当問題(wcsp)
p2133_jobassign3() # 仕事割当問題
p2133_jobassign4() # 仕事割当問題(wcsp)
p2140_QAP() # 二次割当問題
p2150_FPP() # 設備計画問題
p2160_leastsquare2() # 最小二乗問題
p2170_portfolio1() # ポートフォリオ最適化問題
p2230_maxcut2() # 隣接行列(最大カット問題)
p2260_pseudoinverse() # ムーア・ペンローズ一般逆行列