给书配代码-电力经济调度（2）：计及动态约束及节能环保要求的经济调度

Love The Way You Lie 2022-10-31 05:27 93阅读 0赞

电气博文传送门  
[学好电气全靠它，个人电气博文目录（持续更新中…）][Link 1]

本博文系列为我解读电力调度书籍，然后根据自己的理解配上程序。

## 知识点 ##

1.动态模型  
电力经济调度（ED）通常只针对一个时段进行优化，如果计人发电机组的跨时段约束，如爬坡率限制等，它将变为如下的动态模型。  
![2021022209020726.png_pic_center][]  
![20210222090312833.png_pic_center][]

![20210222091022609.png_pic_center][]

2.适应节能环保需求的经济调度扩展模型  
对于火电机组，节能要求主要反映在一段时间的总能耗需控制在一-定范围内， 为此得到扩展模型如下  
![20210222091235754.png_pic_center][]

## 算例 ##

![2021022209192156.png_pic_center][] ![20210222092111333.png_pic_center][]

书上给的答案

![20210222093409113.png_pic_center][]

## 程序实现(语言python） ##

[电力经济调度程序复现（1）：基于拉格朗日及运筹规划方法的经济调度算法][1]  
**1.机组运行参数求解**  
第一台机组运行参数求解（把100-200之间分为4段，每25一段）

import numpy as np
    X0=np.arange(100,201,1)#100-200的整数
    y=4+0.3*x0+0.0007*x0*x0
    #线性规划参数
    print('b',4+0.3*100+0.0007*100*100)
    print('PI11阶段参数',((4+0.3*125+0.0007*125*125)-(4+0.3*100+0.0007*100*100))/(125-100))
    print('PI12阶段参数',((4+0.3*150+0.0007*150*150)-(4+0.3*125+0.0007*125*125))/(150-125))
    print('PI13阶段参数',((4+0.3*175+0.0007*175*175)-(4+0.3*150+0.0007*150*150))/(175-150))
    print('PI14阶段参数',((4+0.3*200+0.0007*200*200)-(4+0.3*175+0.0007*175*175))/(200-175))

和书上给的答案很相近  
![20210222094600949.png_pic_center][]

第二台机组运行参数求解（把120-250之间分为2段，(100-220),(220-250）

import numpy as np
    x0=np.arange(120,251,1)#120-250的整数
    y=3+0.32*x0+0.0004*x0*x0
    #线性规划参数
    print('b',3+0.32*120+0.0004*120*120)
    print('PI21阶段参数',((3+0.32*220+0.0004*220*220)-(4+0.32*100+0.0004*100*100))/(220-100))
    print('PI22阶段参数',((3+0.32*250+0.0004*250*250)-(3+0.32*220+0.0004*220*220))/(250-220))

![20210222095051763.png_pic_center][]

第三台机组运行参数求解（分成3段，（200-150），(250-200)，(300-250)

import numpy as np
    x0=np.arange(150,301,1)#150-300的整数
    y=3.5+0.3*x0+0.00045*x0*x0
    plt.plot(x0,y)
    plt.title('第三台机组运行曲线FG3')
    plt.show()
    
    #线性规划参数
    print('b',3.5+0.3*150+0.00045*150*150)
    print('PI31阶段参数',((3.5+0.3*200+0.00045*200*200)-(3.5+0.3*150+0.00045*150*150))/(200-150))
    print('PI32阶段参数',((3.5+0.3*250+0.00045*250*250)-(3.5+0.3*200+0.00045*200*200))/(250-200))
    print('PI33阶段参数',((3.5+0.3*300+0.00045*300*300)-(3.5+0.3*250+0.00045*250*250))/(300-250))

![20210222095407982.png_pic_center][]

**2运行求解**  
我使用非线性规划来求解  
scipy.optimize 中t minimize

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

来看下改方法的入参

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

#!/usr/bin/env python3
    # -*- coding: utf-8 -*-
    # @Author: yudengwu(余登武）
    # @Date : 2021/2/22
    #@email:1344732766@qq.com
    # -*- coding: utf-8 -*-
    from scipy import  optimize as opt
    import numpy as np
    from scipy.optimize import minimize
    # 目标函数
    def objective(x): #=原文中的mincost #式子3-17
        #500*(41+x[0]*0.456+x[1]*0.496+x[2]*0.524+x[3]*0.564)为第一台机组第一小时购煤成本
        #550*(47.2+x[4]*0.456+x[5]*0.508) 为第二台机组第一小时购煤成本
        #600*(58.6+x[6]*0.458+x[7]*0.502+x[8]*0.548)为第三台机组第一小时购煤成本
    
        #500*(41+x[9]*0.456+x[10]*0.496+x[11]*0.524+x[12]*0.564)为第一台机组第二小时购煤成本
        #550*(47.2+x[13]*0.456+x[14]*0.508) 为第二台机组第二小时购煤成本
        #600 * (58.6 + x[15] * 0.458 + x[16] * 0.502 + x[17] * 0.548)为第三台机组第二小时购煤成本
    
        return 500*(41+x[0]*0.456+x[1]*0.496+x[2]*0.524+x[3]*0.564)+550*(47.2+x[4]*0.456+x[5]*0.508)+600*(58.6+x[6]*0.458+x[7]*0.502+x[8]*0.548)+\
            500*(41+x[9]*0.456+x[10]*0.496+x[11]*0.524+x[12]*0.564)+550*(47.2+x[13]*0.456+x[14]*0.508)+600*(58.6+x[15]*0.458+x[16]*0.502+x[17]*0.548)
    
    
    
    #
    # 约束条件
    
    def constraint1(x):#式子3-18
        return 100+x[0]+x[1]+x[2]+x[3]+120+x[4]+x[5]+150+x[6]+x[7]+x[8]-700 #等式约束，第一小时，三台机组 输出和为700MW
    
    def constraint2(x):#式子3-18
        return 100+x[9]+x[10]+x[11]+x[12]+120+x[13]+x[14]+150+x[15]+x[16]+x[17]-500 #等式约束，第二小时，三台机组 输出和为500MW
    
    #边界约束 式子3-19
    
    b1=(0.0,25)#0<=x[0,1,2,3]<=25 ,0<=x[9,10,11,12]<=25,PG1
    b2=(0.0,100)#0<=x[4]<=100 x21min=0,x21max=100, 0<=x[13]<=100 PG2的约束
    b3=(0.0,30)#0<=x[5]<=30 x22min=0,x22max=30 0<=x[14]<=30 PG2的约束
    b4=(0.0,50)##0<=x[6,7,8]<=50 x3min=0,x3max=50 0<=x[15,16,17]<=50 PG3的约束
    
    bnds = (b1,b1,b1,b1,b2,b3,b4,b4,b4,b1,b1,b1,b1,b2,b3,b4,b4,b4)#边界约束 式子3-19
    
    #机组1爬坡约束 式子3-20
    def constraint3(x):
        return x[9]+x[10]+x[11]+x[12]-x[0]-x[1]-x[2]-x[3]+50 #不等式约束，机组1的爬坡功率约束为50MW
    #机组2爬坡约束 式子3-20
    def constraint4(x):
        return x[13]+x[14]-x[4]-x[5]+60 #不等式约束，机组2的爬坡功率约束为60MW
    #机组3爬坡约束式子3-20
    def constraint5(x):
        return x[15]+x[16]+x[17]-x[6]-x[7]-x[8]+150#不等式约束，机组3的爬坡功率约束为150MW
    
    #第一台机组在第一小时消耗煤炭不超过70t 式子3-22
    def constraint6(x):
        return 70-(41+x[0]*0.456+x[1]*0.496+x[2]*0.524+x[3]*0.564) #不等式约束
    
    #所有机组在二小时内总的SO2的排放量不超过8t,式子3-23,3-24
    
    def constraint7(x):
        return 8-0.02*(41+x[0]*0.456+x[1]*0.496+x[2]*0.524+x[3]*0.564)-0.02*(41+x[9]*0.456+x[10]*0.496+x[11]*0.524+x[12]*0.564)-\
            0.015*(47.2+x[4]*0.456+x[5]*0.508)-0.015*(47.2+x[13]*0.456+x[14]*0.508)-\
            0.01*(58.6+x[6]*0.458+x[7]*0.502+x[8]*0.548)-0.01*(58.6+x[15]*0.458+x[16]*0.502+x[17]*0.548) #不等式约束
    #所有机组在二小时内总的SO2的排放量大于0t,式子3-23,3-24
    def constraint8(x):
        return 0.02*(41+x[0]*0.456+x[1]*0.496+x[2]*0.524+x[3]*0.564)+0.02*(41+x[9]*0.456+x[10]*0.496+x[11]*0.524+x[12]*0.564)+\
            0.015*(47.2+x[4]*0.456+x[5]*0.508)+0.015*(47.2+x[13]*0.456+x[14]*0.508)+\
            0.01*(58.6+x[6]*0.458+x[7]*0.502+x[8]*0.548)+0.01*(58.6+x[15]*0.458+x[16]*0.502+x[17]*0.548) #不等式约束
    
    if __name__ =='__main__':
        con1 = { 'type': 'eq', 'fun': constraint1}
        con2 = { 'type': 'eq', 'fun': constraint2}
        con3 = { 'type': 'ineq', 'fun': constraint3}
        con4 = { 'type': 'ineq', 'fun': constraint4}
        con5 = { 'type': 'ineq', 'fun': constraint5}
        con6 = { 'type': 'ineq', 'fun': constraint6}
        con7 = { 'type': 'ineq', 'fun': constraint7}
        con8 = { 'type': 'ineq', 'fun': constraint8}
        cons = ([con1, con2, con3, con4, con5, con6, con7, con8])  # 8个约束条件
    
        # 初始值
        x0 = np.random.uniform(10, 40, 18)
        # 计算
        solution = minimize(objective, x0, method='SLSQP', \
                            bounds=bnds, constraints=cons)
        x = solution.x
    
        print('目标值: ' + str(objective(x)))
        print('答案为')
        print('第一小时PG1', x[0] + x[1] + x[2] + x[3] + 100)
        print('第一小时PG2', x[4] + x[5] + 120)
        print('第一小时PG3', x[6] + x[7] + x[8] + 150)
        print('--------------')
        print('第二小时PG1', x[9] + x[10] + x[11] + x[12] + 100)
        print('第二小时PG2', x[13] + x[14] + 120)
        print('第二小时PG3', x[15] + x[16] + x[17] + 150)

答案和书上的答案一模一样  
![20210222120602976.png_pic_center][]

进阶版论文链接：[电气期刊论文实现：考虑排放目标和输电损耗的经济调度【有代码】][Link 2]  
![在这里插入图片描述][20210222120823535.jpg_pic_center]  
作者：电气-余登武。写作不容易，请点个赞再走。  
![在这里插入图片描述][20210222120908683.gif_pic_center]

[Link 1]: https://blog.csdn.net/KOBEYU652453/article/details/117459596
[2021022209020726.png_pic_center]: /images/20221024/af92e5d5d61d4ef6a0f9d51a5cea817e.png
[20210222090312833.png_pic_center]: /images/20221024/b39563246d6d4ec4ad3dc3dbeb29ae63.png
[20210222091022609.png_pic_center]: /images/20221024/05644e4e80c14a06a55634dbd4d5a13c.png
[20210222091235754.png_pic_center]: /images/20221024/da5e7ef7b2d84b34be5f87d25df09703.png
[2021022209192156.png_pic_center]: https://img-blog.csdnimg.cn/2021022209192156.png#pic_center
[20210222092111333.png_pic_center]: /images/20221024/6f8de4aeaa144a8297ef7104f417d8b4.png
[20210222093409113.png_pic_center]: /images/20221024/6e72da85fbcc497a8323082977efc5df.png
[1]: https://blog.csdn.net/kobeyu652453/article/details/113876312
[20210222094600949.png_pic_center]: /images/20221024/4cfa8f68b20c407092bc246f628abf8f.png
[20210222095051763.png_pic_center]: /images/20221024/045af1c2b7d64ef5b9513227ff7764cf.png
[20210222095407982.png_pic_center]: /images/20221024/69f2b356873f4d60b449d5659aef227c.png
[20210222120602976.png_pic_center]: /images/20221024/69b4e1b18b67423795bf5084d7420519.png
[Link 2]: https://blog.csdn.net/KOBEYU652453/article/details/117386398
[20210222120823535.jpg_pic_center]: /images/20221024/519d8ece1f4d4d1b9e0f6860495b3df7.png
[20210222120908683.gif_pic_center]: /images/20221024/a3f7ce2df5304b5c8ff210a36e490216.png