電氣博文傳送門(mén)
學(xué)好電氣全靠它，個(gè)人電氣博文目錄（持續(xù)更新中…）

知識(shí)點(diǎn)

電力經(jīng)濟(jì)調(diào)度(Economic Dispatch, ED)的目標(biāo)是追求某個(gè)研究時(shí)段內(nèi)所有開(kāi)機(jī)機(jī)組(Committed Unit)的運(yùn)行成本最小，約束條件包括電力供需平衡及機(jī)組出力限制，優(yōu)化變量是每臺(tái)機(jī)組的出力。由于ED一般是在機(jī)組組合(Unit Commitment, UC)結(jié)束之后進(jìn)行，因此如下的經(jīng)濟(jì)模型中不會(huì)包含反映機(jī)組開(kāi)、停狀態(tài)的0/1二元整型變量。
經(jīng)濟(jì)模型中個(gè)會(huì)包含反映機(jī)組開(kāi)、停狀態(tài)的0/1 二元整型變量。

算例

程序python實(shí)現(xiàn)

非線性規(guī)劃（scipy.optimize.minimize）
一.背景：
現(xiàn)在項(xiàng)目上有一個(gè)用python 實(shí)現(xiàn)非線性規(guī)劃的需求。非線性規(guī)劃可以簡(jiǎn)單分兩種，目標(biāo)函數(shù)為凸函數(shù) or 非凸函數(shù)。

凸函數(shù)的 非線性規(guī)劃，比如fun=x2+y2+x*y，有很多常用的python庫(kù)來(lái)完成，網(wǎng)上也有很多資料，比如CVXPY

非凸函數(shù)的 非線性規(guī)劃（求極值），從處理方法來(lái)說(shuō)，可以嘗試以下幾種：

1.純數(shù)學(xué)方法，求導(dǎo)求極值；

2.使用神經(jīng)網(wǎng)絡(luò)，深度學(xué)習(xí)來(lái)處理，可參考反向傳播算法中鏈?zhǔn)角髮?dǎo)的過(guò)程；

3.尋找一些python庫(kù)來(lái)做，本文介紹scipy.optimize.minimize的使用方法

二.庫(kù)方法介紹

官方文檔：https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html

來(lái)看下改方法的入?yún)?/p>

scipy.optimize.minimize(fun, x0, args=(), method=None, jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None)

解釋：
fun: 求最小值的目標(biāo)函數(shù)
x0:變量的初始猜測(cè)值，如果有多個(gè)變量，需要給每個(gè)變量一個(gè)初始猜測(cè)值。minimize是局部最優(yōu)的解法，所以
args:常數(shù)值,后面demo會(huì)講解，fun中沒(méi)有數(shù)字，都以變量的形式表示，對(duì)于常數(shù)項(xiàng)，需要在這里給值
method:求極值的方法，官方文檔給了很多種。一般使用默認(rèn)。每種方法我理解是計(jì)算誤差，反向傳播的方式不同而已，這塊有很大理論研究空間
constraints:約束條件，針對(duì)fun中為參數(shù)的部分進(jìn)行約束限制

常規(guī)解法
語(yǔ)言：python

from scipy.optimize import minimize import numpy as np #目標(biāo)函數(shù)即min(FG1+FG2+FG3) def fun(args):a1,a2,a3,b1,b2,b3,c1,c2,c3=argsv=lambda x: (a1+a2*x[0]+a3*x[0]*x[0]+b1+b2*x[1]+b3*x[1]*x[1]+c1+c2*x[2]+c3*x[2]*x[2])return vdef con(args):# 約束條件 分為eq 和ineq#eq表示 函數(shù)結(jié)果等于0 ； ineq 表示 表達(dá)式大于等于0d, x0min, x0max,x1min, x1max,x2min,x2max= args1# (700, 100, 200, 120, 250, 150,300) # d1, x0min, x0max,x1min, x1max,x2min,x2max#約束條件d1-x[0]-x[1]-x[2]轉(zhuǎn)換為x[0]+x[1]+x[2]=700#100<=x[0]<=200,120<=x[1]<=250,150<=x[2]<=300cons = ({'type': 'eq', 'fun': lambda x: d-x[0]-x[1]-x[2]},\{'type': 'ineq', 'fun': lambda x:x[0]-x0min}, \{'type': 'ineq', 'fun': lambda x: -x[0] + x0max}, \{'type': 'ineq', 'fun': lambda x: x[1] - x1min}, \{'type': 'ineq', 'fun': lambda x: -x[1] + x1max}, \{'type': 'ineq', 'fun': lambda x: x[2] - x2min}, \{'type': 'ineq', 'fun': lambda x: -x[2] + x2max})return consif __name__ == "__main__":# 定義常量值args = (4,0.3,0.0007,3,0.32,0.0004,3.5,0.3,0.00045) # a1,a2,a3,b1,b2,b3,c1,c2,c3# 設(shè)置參數(shù)范圍/約束條件args1 = (700, 100, 200, 120, 250, 150,300) # d1, x0min, x0max,x1min, x1max,x2min,x2maxcons = con(args1)# 設(shè)置x初始猜測(cè)值x0 = np.array((150, 250, 20))res = minimize(fun(args), x0, method='SLSQP', constraints=cons)print('代價(jià)',res.fun)print(res.success)print('解',res.x)

如果使用書(shū)上的解得到的代價(jià)為4+0.3x175+0.0007x175x1 75+3+0.32x250+0.0004x250x250+3.5+0.3x275+0.00045275275=305.96875

書(shū)上解法復(fù)現(xiàn)
1、求解參數(shù)

import matplotlib.pyplot as plt from pylab import * mpl.rcParams['font.sans-serif'] = ['SimHei'] mpl.rcParams['axes.unicode_minus'] = False import numpy as npx0=np.arange(100,201,1)#100-200的整數(shù) y=4+0.3*x0+0.0007*x0*x0 plt.plot(x0,y) plt.title('第一臺(tái)機(jī)組運(yùn)行曲線FG1') plt.show() #線性規(guī)劃參數(shù) print('b',4+0.3*100+0.0007*100*100) print('PI11階段參數(shù)',((4+0.3*125+0.0007*125*125)-(4+0.3*100+0.0007*100*100))/(125-100)) print('PI12階段參數(shù)',((4+0.3*150+0.0007*150*150)-(4+0.3*125+0.0007*125*125))/(150-125)) print('PI13階段參數(shù)',((4+0.3*175+0.0007*175*175)-(4+0.3*150+0.0007*150*150))/(175-150)) print('PI14階段參數(shù)',((4+0.3*200+0.0007*200*200)-(4+0.3*175+0.0007*175*175))/(200-175)) 求解得到的參數(shù)和書(shū)上給的解很接近。

同理得到FG2,

x0=np.arange(120,251,1)#120-250的整數(shù) y=3+0.32*x0+0.0004*x0*x0 plt.plot(x0,y) plt.title('第二臺(tái)機(jī)組運(yùn)行曲線FG2') plt.show()#線性規(guī)劃參數(shù) print('b',3+0.32*120+0.0004*120*120) print('PI21階段參數(shù)',((3+0.32*220+0.0004*220*220)-(4+0.32*100+0.0004*100*100))/(220-100)) print('PI22階段參數(shù)',((3+0.32*250+0.0004*250*250)-(3+0.32*220+0.0004*220*220))/(250-220))

FG3的參數(shù)

x0=np.arange(150,301,1)#150-300的整數(shù) y=3.5+0.3*x0+0.00045*x0*x0 plt.plot(x0,y) plt.title('第三臺(tái)機(jī)組運(yùn)行曲線FG3') plt.show()#線性規(guī)劃參數(shù) print('b',3.5+0.3*150+0.00045*150*150) print('PI31階段參數(shù)',((3.5+0.3*200+0.00045*200*200)-(3.5+0.3*150+0.00045*150*150))/(200-150)) print('PI32階段參數(shù)',((3.5+0.3*250+0.00045*250*250)-(3.5+0.3*200+0.00045*200*200))/(250-200)) print('PI33階段參數(shù)',((3.5+0.3*300+0.00045*300*300)-(3.5+0.3*250+0.00045*250*250))/(300-250))

我們求得的運(yùn)行參數(shù)和書(shū)上給的，只是在位數(shù)上不同。為了統(tǒng)一，我后面采取書(shū)上的參數(shù)。
2、運(yùn)行求解

from scipy.optimize import minimize import numpy as np #目標(biāo)函數(shù)即min(FG1+FG2+FG3) def fun(args):Z0,a0,a1,a2,a3,Z1,b0,b1,Z2,c0,c1,c2=args#Z0+x[0]*a0+x[1]*a1+x[2]*a2+x[3]*a3為FG1#Z1+x[4]*b0+x[5]*b1為FG2#Z2+x[6]*c0+x[7]*c1+x[8]*c2 為FG3v=lambda x: (Z0+x[0]*a0+x[1]*a1+x[2]*a2+x[3]*a3+Z1+x[4]*b0+x[5]*b1+Z2+x[6]*c0+x[7]*c1+x[8]*c2)return vdef con(args):# 約束條件 分為eq 和ineq#eq表示 函數(shù)結(jié)果等于0 ； ineq 表示 表達(dá)式大于等于0d, x1min, x1max,x21min, x21max,x22min,x22max,x3min,x3max= args1#等式約束條件d-x[0]-x[1]-x[2]-x[3]-x[4]-x[5]-x[6]-x[7]-x[8]=700-100-120-150-x[0]-x[1]-x[2]-x[3]-x[4]-x[5]-x[6]-x[7]-x[8]#不等式約束條件#0<=x[0,1,2,3]<=25 x1min=0,x1max=25 PG1的粵省事#0<=x[4]<=100 x21min=0,x21max=100,0<=x[5]<=30 x22min=0,x22max=30 PG2的約束#0<=x[6,7,8]<=50 x3min=0,x3max=50 PG3的約束cons = ({'type': 'eq', 'fun': lambda x,: d-x[0]-x[1]-x[2]-x[3]-x[4]-x[5]-x[6]-x[7]-x[8]},\{'type': 'ineq', 'fun': lambda x:x[0]-x1min},{'type': 'ineq', 'fun': lambda x: -x[0] + x1max},{'type': 'ineq', 'fun': lambda x: x[1] - x1min},{'type': 'ineq', 'fun': lambda x: -x[1] + x1max},{'type': 'ineq', 'fun': lambda x: x[2] - x1min},{'type': 'ineq', 'fun': lambda x: -x[2] + x1max},{'type': 'ineq', 'fun': lambda x: x[3] - x1min},{'type': 'ineq', 'fun': lambda x: -x[3] + x1max},{'type': 'ineq', 'fun': lambda x: x[4] - x21min},{'type': 'ineq', 'fun': lambda x: -x[4] + x21max},{'type': 'ineq', 'fun': lambda x: x[5] - x22min},{'type': 'ineq', 'fun': lambda x: -x[5] + x22max},{'type': 'ineq', 'fun': lambda x: x[6] - x3min},{'type': 'ineq', 'fun': lambda x: -x[6] + x3max},{'type': 'ineq', 'fun': lambda x: x[7] - x3min},{'type': 'ineq', 'fun': lambda x: -x[7] + x3max},{'type': 'ineq', 'fun': lambda x: x[8] - x3min},{'type': 'ineq', 'fun': lambda x: -x[8] + x3max},)return consif __name__ == "__main__":# 定義常量值args = (41,0.456,0.496,0.524,0.564,47.2,0.456,0.508,58.6,0.458,0.502,0.548) # Z0,a0,a1,a2,a3,,Z1,b0,b1,Z2,c0,c1,c2# 設(shè)置參數(shù)范圍/約束條件args1 = (330, 0,25,0,100,0,30,0,50) # d=700-100-120-150, x1min, x1max,x21min, x21max,x22min,x22max,x3min,x3maxcons = con(args1)# 設(shè)置x初始猜測(cè)值x0 = np.array((10, 10, 10,10,20,20,20,10,10))res = minimize(fun(args), x0, method='SLSQP', constraints=cons)print('代價(jià)',res.fun)print(res.success)print('解',[np.around(i) for i in res.x])print('PG1',res.x[0]+res.x[1]+res.x[2]+res.x[3]+100)print('PG2', res.x[4] + res.x[5] + 120)print('PG3', res.x[6] + res.x[7] + res.x[8]+150)

作者：電氣-余登武

總結(jié)

以上是生活随笔為你收集整理的给书配代码-电力经济调度（1）：基于拉格朗日及运筹规划方法的经济调度算法的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問(wèn)題。

如果覺(jué)得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。