目录
1 电力系统环境经济调度数学模型
2 算例——IEEE10节点
2.1 数据
2.2 Python代码学习
3 一点拓展知识
1 电力系统环境经济调度数学模型
2 算例——IEEE10节点
2.1 数据
我弄成一个表格，方便编程读写：
2.2 Python代码学习
多目标遗传算法NSGAII在电力系统多目标问题有广泛的应用，只要把文中的目标函数和约束条件换了，就搞定啦。
部分代码：
#========导入第三方库=============import numpy as npimport pandas as pdimport matplotlib.pyplot as pltimport matplotlib as mplmpl.rcParams['font.sans-serif'] = ['SimHei'] # 指定默认字体mpl.rcParams['axes.unicode_minus'] = False # 解决保存图像是负号'-'显示为方块的问题import matplotlib;matplotlib.use('TkAgg')from tqdm import tqdm # 进度条设置# ========导入数据=============# ~~~~~~~~~~~~~~~~~~~~~读取文件~~~~~~~~~~~~~~~~~~~~~~file=pd.read_csv('IEEE10.csv')#机组参数#B=pd.read_csv('B10.csv')#B=np.reshape()B=np.array([[0.000049, 0.000014,0.000015,0.000015,0.000016,0.000017,0.000017,0.000018,0.000019,0.00002],[0.000014,0.000045,0.000016,0.000016,0.000017,0.000015,0.000015,0.000016,0.000018,0.000018],[0.000015,0.000016,0.000039,0.00001,0.000012,0.000012,0.000014,0.000014,0.000016,0.000016],[0.000015,0.000016,0.00001,0.00004,0.000014,0.00001,0.000011,0.000012,0.000014,0.000015],[0.000016,0.000017,0.000012,0.000014,0.000035,0.000011,0.000013,0.000013,0.000015,0.000016],[0.000017,0.000015,0.000012,0.00001,0.000011,0.000036,0.000012,0.000012,0.000014,0.000015],[0.000017,0.000015,0.000014,0.000011,0.000013,0.000012,0.000038,0.000016,0.000016,0.000018],[0.000018,0.000016,0.000014,0.000012,0.000013,0.000012,0.000016,0.00004,0.000015,0.000016],[0.000019,0.000018,0.000016,0.000014,0.000015,0.000014,0.000016,0.000015,0.000042,0.000019],[0.00002,0.000018,0.000016,0.000015,0.000016,0.000015,0.000018,0.000016,0.000019,0.000044]])P_max = file['b'] #机组上限P_min = file['c'] # #机组下限# ~~~~~~~~~~~~~~~~~机组特性系数~~~~~~~~~~~~~~~~~~~~~~~~·ai = file['d']bi=file['e']ci=file['f']di=file['g']ei=file['h']# ~~~~~~~排放特性系数~~~~~~~~~~~~~~~~·ali=file['g']beti=file['h']gari=file['i']eti=file['l']dali=file['m']# ======多目标优化算法一次 只能求解单个时刻的解===============class GaMultiobjective(object):def __init__(self, Pload1, P_max, P_min, ai, bi, ci, di, ei, ali, beti, gari, eti, dali):self.Pload1 = Pload1 # 负荷self.P_max = P_max # 机组上限self.P_min = P_min # 机组下限self.ai = aiself.bi = biself.ci = ciself.di = diself.ei = eiself.ali = aliself.beti = betiself.gari = gariself.eti = etiself.dali = dali# ===初始化种群====def Initialize(self):X = np.zeros((self.NP, 10)) # 初始化群体，10代表 个机组出力for n in range(self.NP): # 遍历每一个粒子X[n, 0] = np.random.uniform(10, 55, 1)[0] # G1X[n, 1] = np.random.uniform(20, 80, 1)[0] # G2X[n, 2] = np.random.uniform(47, 120, 1)[0] # G3X[n, 3] = np.random.uniform(20, 130, 1)[0] # G4X[n, 4] = np.random.uniform(50, 160, 1)[0] # G5X[n, 5] = np.random.uniform(70, 240, 1)[0] # G6X[n, 6] = np.random.uniform(60, 300, 1)[0] # G7X[n, 7] = np.random.uniform(70, 340, 1)[0] # G8X[n, 8] = np.random.uniform(130, 470, 1)[0] # G9X[n, 9] = np.random.uniform(150, 470, 1)[0] # G10return X# ==========定义目标函数、和对应的惩罚项=============# ===定义函数1 目标函数1：系统运行成本===def function1(self, X1):"""个体目标函数:param X1: （个体[G1，G2，G3，G4，G5，G6,G7,G8,G9,G10]:return: 函数1目标函数值"""SUMCOST = [] # 存储总的成本for i in range (9): # 遍历每一台机组,Python是从0开始索引，0-9：10台机组cost = self.ci[i] * X1[i] * X1[i] + self.bi[i]* X1[i] +self.ai[i]+\np.abs(self.di[i]*np.sin(ei[i]*(self.P_min[i]-X1[i]))) #考虑阀点效应SUMCOST.append(cost)return np.sum(SUMCOST)# ====定义函数2 总污染排放量====def function2(self, X1):"""个体目标函数:param X1: （个体[G1，G2，G3，G4，G5，G6,G7,G8,G9,G10]:return: 函数2目标函数值"""emission=[] #储存总的污染排放量for i in range(9):e=self.ali[i]+self.beti[i]*X1[i]+self.gari[i]*X1[i]*X1[i]+self.eti[i]*np.exp(self.dali[i]*X1[i])emission.append(e)return np.sum(emission)# ===对应的约束 负荷平衡约束(本例子calc_e（） 不起作用，已通过其它方法解决掉负荷平衡约束）===def calc_e(self, X1):"""函数1 对应的个体惩罚项:param X1: （个体[G1，G2，G3，G4，G5，G6,G7,G8,G9,G10]:return:"""Ploss=0for i in range(9):for j in range(9):Ploss+=X1[i]*B[i][j]*X1[j]# cost=np.abs(X1[0]+X1[1]+X1[2]+X1[3]+X1[4]+X1[5]+X1[6]+X1[7]+X1[8]+X1[9]-np.sum(Ploss)-self.Pload1)# return costreturn max(0,abs(sum(X1)-Ploss-Pload1-0.1))
#========导入第三方库=============
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.sans-serif'] = ['SimHei'] # 指定默认字体
mpl.rcParams['axes.unicode_minus'] = False # 解决保存图像是负号'-'显示为方块的问题
import matplotlib;
matplotlib.use('TkAgg')
from tqdm import tqdm # 进度条设置
# ========导入数据=============
# ~~~~~~~~~~~~~~~~~~~~~读取文件~~~~~~~~~~~~~~~~~~~~~~
file=pd.read_csv('IEEE10.csv')#机组参数
#B=pd.read_csv('B10.csv')
#B=np.reshape()
B=np.array([[0.000049, 0.000014 ,0.000015, 0.000015, 0.000016, 0.000017 ,0.000017, 0.000018 ,0.000019, 0.00002],
[0.000014, 0.000045, 0.000016, 0.000016, 0.000017 ,0.000015, 0.000015 ,0.000016 ,0.000018, 0.000018],
[0.000015, 0.000016, 0.000039, 0.00001 ,0.000012, 0.000012 ,0.000014, 0.000014, 0.000016, 0.000016],
[0.000015, 0.000016, 0.00001, 0.00004, 0.000014 ,0.00001 ,0.000011, 0.000012, 0.000014, 0.000015],
[0.000016, 0.000017, 0.000012, 0.000014, 0.000035, 0.000011, 0.000013, 0.000013, 0.000015, 0.000016],
[0.000017, 0.000015, 0.000012, 0.00001, 0.000011 ,0.000036, 0.000012, 0.000012, 0.000014, 0.000015],
[0.000017 ,0.000015 ,0.000014, 0.000011, 0.000013, 0.000012, 0.000038, 0.000016, 0.000016, 0.000018],
[0.000018 ,0.000016, 0.000014, 0.000012, 0.000013, 0.000012, 0.000016 ,0.00004, 0.000015, 0.000016],
[0.000019 ,0.000018 ,0.000016, 0.000014, 0.000015 ,0.000014 ,0.000016, 0.000015, 0.000042 ,0.000019],
[0.00002, 0.000018 ,0.000016 ,0.000015, 0.000016, 0.000015, 0.000018 ,0.000016, 0.000019 ,0.000044]])
P_max = file['b'] #机组上限
P_min = file['c'] # #机组下限
# ~~~~~~~~~~~~~~~~~机组特性系数~~~~~~~~~~~~~~~~~~~~~~~~·
ai = file['d']
bi=file['e']
ci=file['f']
di=file['g']
ei=file['h']
# ~~~~~~~排放特性系数~~~~~~~~~~~~~~~~·
ali=file['g']
beti=file['h']
gari=file['i']
eti=file['l']
dali=file['m']
# ======多目标优化算法一次 只能求解单个时刻的解===============
class GaMultiobjective(object):
def __init__(self, Pload1, P_max, P_min, ai, bi, ci, di, ei, ali, beti, gari, eti, dali):
self.Pload1 = Pload1 # 负荷
self.P_max = P_max # 机组上限
self.P_min = P_min # 机组下限
self.ai = ai
self.bi = bi
self.ci = ci
self.di = di
self.ei = ei
self.ali = ali
self.beti = beti
self.gari = gari
self.eti = eti
self.dali = dali
# ===初始化种群====
def Initialize(self):
X = np.zeros((self.NP, 10)) # 初始化群体，10代表 个机组出力
for n in range(self.NP): # 遍历每一个粒子
X[n, 0] = np.random.uniform(10, 55, 1)[0] # G1
X[n, 1] = np.random.uniform(20, 80, 1)[0] # G2
X[n, 2] = np.random.uniform(47, 120, 1)[0] # G3
X[n, 3] = np.random.uniform(20, 130, 1)[0] # G4
X[n, 4] = np.random.uniform(50, 160, 1)[0] # G5
X[n, 5] = np.random.uniform(70, 240, 1)[0] # G6
X[n, 6] = np.random.uniform(60, 300, 1)[0] # G7
X[n, 7] = np.random.uniform(70, 340, 1)[0] # G8
X[n, 8] = np.random.uniform(130, 470, 1)[0] # G9
X[n, 9] = np.random.uniform(150, 470, 1)[0] # G10
return X
# ==========定义目标函数、和对应的惩罚项=============
# ===定义函数1 目标函数1：系统运行成本===
def function1(self, X1):
"""
个体目标函数
:param X1: （个体[G1，G2，G3，G4，G5，G6,G7,G8,G9,G10]
:return: 函数1目标函数值
"""
SUMCOST = [] # 存储总的成本
for i in range (9): # 遍历每一台机组,Python是从0开始索引，0-9：10台机组
cost = self.ci[i] * X1[i] * X1[i] + self.bi[i]* X1[i] +self.ai[i]+\
np.abs(self.di[i]*np.sin(ei[i]*(self.P_min[i]-X1[i]))) #考虑阀点效应
SUMCOST.append(cost)
return np.sum(SUMCOST)
# ====定义函数2 总污染排放量====
def function2(self, X1):
"""
个体目标函数
:param X1: （个体[G1，G2，G3，G4，G5，G6,G7,G8,G9,G10]
:return: 函数2目标函数值
"""
emission=[] #储存总的污染排放量
for i in range(9):
e=self.ali[i]+self.beti[i]*X1[i]+self.gari[i]*X1[i]*X1[i]+self.eti[i]*np.exp(self.dali[i]*X1[i])
emission.append(e)
return np.sum(emission)
# ===对应的约束 负荷平衡约束(本例子calc_e（） 不起作用，已通过其它方法解决掉负荷平衡约束）===
def calc_e(self, X1):
"""
函数1 对应的个体惩罚项
:param X1: （个体[G1，G2，G3，G4，G5，G6,G7,G8,G9,G10]
:return:
"""
Ploss=0
for i in range(9):
for j in range(9):
Ploss+=X1[i]*B[i][j]*X1[j]
# cost=np.abs(X1[0]+X1[1]+X1[2]+X1[3]+X1[4]+X1[5]+X1[6]+X1[7]+X1[8]+X1[9]-np.sum(Ploss)-self.Pload1)
# return cost
return max(0,abs(sum(X1)-Ploss-Pload1-0.1))
3 一点拓展知识
现代这种“探索、征服”的心态，从世界地图的演变可以看得一目了然。早在历史进到现代之前，许多文化就已经有了自己的世界地图。当然，当时并没有人真正知道全世界是什么样子，在亚非大陆上的人对美洲一无所知，美洲文化也不知道亚非大陆上的情形。但碰到不熟悉的地区，地图上不是一笔未提，就是画上了想象出来的怪物和奇景。这些地图上并没有空白的空间，让人觉得全世界就在自己的掌握之中。
在15、16世纪，欧洲人的世界地图开始出现大片空白。从这点可以看出科学心态的发展，以及欧洲帝国主义的动机。地图上的空白可以说是在心理及思想上的一大突破，清楚表明欧洲人愿意承认自己对于一大部分的世界还一无所知。
图1 1459年欧洲人的世界地图。可以看到地图上似乎巨细靡遗，就算是当时欧洲人根
本一无所知的南非地区，都有密密麻麻的信息。
1492年，哥伦布从西班牙出发向西航行，希望能找到一条前往东亚的新航线。哥伦布当时相信的仍然是旧的世界地图，以为全世界在地图上一览无遗。哥伦布从旧地图推算，日本应该位于西班牙以西大约7000公里远。但事实上，从西班牙到东亚的距离要超过两万公里，而且中间还隔着个他不知道的美洲大陆。1492年10月12日大约凌晨2点，哥伦布一行人与这片未知大陆有了第一次接触。皮塔号（Pinta）的瞭望手胡安·罗德里格斯·贝尔梅霍(Juan Rodriguez Berme jo)从桅杆上看到了现在的巴哈马群岛，高呼着:“有陆地!有陆地!”
哥伦布当时相信这个小岛就位于东亚海外，属于“Indies”(印度地方，包含今日印度、中南半岛及东印度群岛等地），所以他把当地人称为“Indians”(这正是美国原住民也被称为“印第安人”的原因）。一直到他过世，哥伦布都不认为自己犯了一个大错。不论是对他还是许多当代的人来说，说他发现了一个完全未知的大陆，这根本难以想象。毕竟千百年来，不管是那些伟大的思想家和学者甚至是不可能犯错的《圣经》，都只知道有欧洲、非洲和亚洲。怎么有可能他们全错了呢?难道《圣经》居然漏了大半个世界，只字未提?这种情况，就好像是说在1969年阿波罗11号要前往月球的途中，居然撞到了另一个从来没人看到的月亮。而正因为哥伦布不愿意接受自己的无知，我们可以说他仍然是个中世纪的人，深信着自己已经知道了全世界，所以就算已经有了如此重大的发现，也无法说服他。
