最適化(線形計画)の繰り返し

前提・実現したいこと

pythonで、pvパネル、蓄電池を備えた家での24時間でのコスト最小化を目的としています。
入力変数として電力需要、pv発電量を乱数として与えています。
最適化の繰り返しの処理をしなければ上手く最適化することができたので最適化に関するプログラムは正しいと思います。

ここから、繰り返しごとに変数を変更して、最適化を繰り返して複数の最適値を発生させて、その中から最小の最適値を得るプログラムを書こうとしています。
繰り返し回数は３回としています。

発生している問題・エラーメッセージ

最後以外の目的関数の値が上手く出ません。

### 該当のソースコード
from pulp import LpVariable, LpProblem, LpStatus, lpSum, value
import numpy as np
import random,pulp

# 変数（連続）を宣言
demand=[np.random.uniform (2,5) for i in range(24)]
solar=[np.random.uniform (0,10) for i in range(24)]

# Tarrifs
#########
buy = 30.0  # 買電価格buy price (円/kWh)
sell = 10.0  # 売電価格sell price (円/kWh)

# Battery parameters
####################
 
Energymax = 25  # 容量battery capacity (kWh)
powermax = (1 / 6) * Energymax  # 出力battery max power (kW)
eff=0.8 #効率

# Battery variables

Energy = [LpVariable('Energy_{}'.format(i), 0, None) for i in range(24) ] # battery energy (kWh)、バッテリーの残量
pcharge =[LpVariable('pcharge_{}'.format(i), 0, None) for i in range(24) ]  # battery charge (kW)
pdischarge = [LpVariable('pdischarge_{}'.format(i), 0, None) for i in range(24) ]  # battery discharge (kW)


power = [LpVariable('power_{}'.format(i), None, None) for i in range(24)] # total house power (kW)
costpower = [LpVariable('costpower_{}'.format(i), None, None) for i in range(24)] # power cost


for k in range(3):
 objective=[LpVariable('objective_{}'.format(k), None, None) for k in range(3)] 
# Optimisation problem　目的関数(最小化)
 prb = LpProblem('Battery Operation')
# Objective　目的関数
 prb +=objective[k]
 objective[k]=lpSum(costpower)
# Constraints　制約
 for i in range(24):
  prb += power[i] == pcharge[i] - pdischarge[i] + demand[i] - solar[i]  # total power
  prb += costpower[i] >= buy * power[i] 
  prb += costpower[i] >= sell * power[i]  
  prb += Energy[i] <= Energymax  # battery capacity
  prb += pcharge[i] <= powermax
  prb += pdischarge[i] <= powermax
# Battery charge state constraints
# Batteries must start and finish half charged
  prb += Energy[0] == 0.5*Energymax   # starting energy
  prb += Energy[23] == 0.5*Energymax  # finishing energy
 for i in range(1, 24):
  prb += Energy[i] == Energy[i-1] +  eff*pcharge[i] - pdischarge[i]*(1/eff)  # battery transitions
# Solve problem
  prb.solve()


print('Status {}'.format(LpStatus[prb.status]))
print('Cost {}'.format(value(prb.objective)))
for k in range(3):
    print('objective{} = {}'.format( k, pulp.value(objective[k]) ))
    print('objective{} = {}'.format( k, objective[k]))

試したこと

for文を使って最適化を繰り返そうとしましたが上手くいきませんでした。
結果は
Status Unbounded
Cost 0.0
objective0 = None
objective0 = objective_0
objective1 = None
objective1 = objective_1
objective2 = 373.20738500000004
objective2 = costpower_0 + costpower_1 + costpower_10 + costpower_11 + costpower_12 + costpower_13 + costpower_14 + costpower_15 + costpower_16 + costpower_17 + costpower_18 + costpower_19 + costpower_2 + costpower_20 + costpower_21 + costpower_22 + costpower_23 + costpower_3 + costpower_4 + costpower_5 + costpower_6 + costpower_7 + costpower_8 + costpower_9
のようになってしまい、最後しか計算できませんでした。

補足情報（FW/ツールのバージョンなど）

ここにより詳細な情報を記載してください。

行動規範の内容に同意します

回答2件

ベストアンサー

前回回答は「最後以外も計算できる」を満たす最小限のものでしたが、「期待通りの計算ができる」ものにしたコードを提示します。

前回回答からの変更点は２点
(1) 毎回のkループでprb = LpProblem('Battery Operation')によりprbを初期化していたため、最適値が上書きされてしまっていました。最適値をv[k]として記録しておくようにしました。
(2) 最適値を記録しても、毎回のkループで説明変数となるdemand、solar（の中身）が同じです。3回別の最適値を求めたいのですから、説明変数も3セット別のものが必要です。よって、最初の乱数を24×3セット生成するようにして、毎回のkループで別のところを参照するように、ループの中のリストインデックスの取り方も変えました。

なお、蛇足ながら、インデントは空白4つがPythonの掟です。せめて2つにしてほしいです。空白1つのインデントをされているコードは初めてみました。インデントもなおしました。

コードは以下です。

Python
1### 該当のソースコード
2from pulp import LpVariable, LpProblem, LpStatus, lpSum, value
3import numpy as np
4import random,pulp
5
6# 変数（連続）を宣言
7demand=[np.random.uniform (2,5) for i in range(24*3)]
8solar=[np.random.uniform (0,10) for i in range(24*3)]
9
10# Tarrifs
11#########
12buy = 30.0  # 買電価格buy price (円/kWh)
13sell = 10.0  # 売電価格sell price (円/kWh)
14
15# Battery parameters
16####################
17
18Energymax = 25  # 容量battery capacity (kWh)
19powermax = (1 / 6) * Energymax  # 出力battery max power (kW)
20eff=0.8 #効率
21
22# Battery variables
23
24Energy = [LpVariable('Energy_{}'.format(i), 0, None) for i in range(24) ] # battery energy (kWh)、バッテリーの残量
25pcharge =[LpVariable('pcharge_{}'.format(i), 0, None) for i in range(24) ]  # battery charge (kW)
26pdischarge = [LpVariable('pdischarge_{}'.format(i), 0, None) for i in range(24) ]  # battery discharge (kW)
27
28
29power = [LpVariable('power_{}'.format(i), None, None) for i in range(24)] # total house power (kW)
30costpower = [LpVariable('costpower_{}'.format(i), None, None) for i in range(24)] # power cost
31
32objective=[LpVariable('objective_{}'.format(k), None, None) for k in range(3)]
33v = [None]*3
34for k in range(3):
35	# Optimisation problem　目的関数(最小化)
36	prb = LpProblem('Battery Operation')
37	# Objective　目的関数
38	prb +=objective[k]
39	objective[k]=lpSum(costpower)
40	# Constraints　制約
41	for i in range(24):
42		prb += power[i] == pcharge[i] - pdischarge[i] + demand[i*(k+1)] - solar[i*(k+1)]  # total power
43		prb += costpower[i] >= buy * power[i] 
44		prb += costpower[i] >= sell * power[i]  
45		prb += Energy[i] <= Energymax  # battery capacity
46		prb += pcharge[i] <= powermax
47		prb += pdischarge[i] <= powermax
48		# Battery charge state constraints
49		# Batteries must start and finish half charged
50		prb += Energy[0] == 0.5*Energymax   # starting energy
51		prb += Energy[23] == 0.5*Energymax  # finishing energy
52	for i in range(1, 24):
53		prb += Energy[i] == Energy[i-1] +  eff*pcharge[i] - pdischarge[i]*(1/eff)  # battery transitions
54		# Solve problem
55		prb.solve()
56	v[k] = pulp.value(objective[k])
57
58print('Status {}'.format(LpStatus[prb.status]))
59print('Cost {}'.format(value(prb.objective)))
60for k in range(3):
61	print('objective{} = {}'.format(k, v[k]))
62	print('objective{} = {}'.format(k, objective[k]))

投稿2020/10/11 00:29

toast-uz

総合スコア3266

simpkins

2020/10/11 00:55

朝早くに本当にありがとうございます。答えだけでなく、インデントのことまで教えていただきありがとうございます。本当に助かりました。

行動規範の内容に同意します

下記の部分の順番が逆ですね。前のループでせっかく計算したobjectiveを、再定義して消してしまっています。

###現状のコード

Python
1for k in range(3):
2 objective=[LpVariable('objective_{}'.format(k), None, None) for k in range(3)]

###正しいコード

Python
1objective=[LpVariable('objective_{}'.format(k), None, None) for k in range(3)] 
2for k in range(3):

投稿2020/10/10 09:11

toast-uz

総合スコア3266

simpkins

2020/10/10 11:22

回答ありがとうございます。ご指摘いただいた通りに修正したらobjecctiveがNoneになることはなくなったのですが、答えがすべてのループで同じになってしまいました。これはすべてのループで、乱数で与えたdemandとsolarが同じになっていることが原因だと考えているのですが、ループ毎に異なる乱数を発生させるにはどうすればいいでしょうか？