运筹学课程设计

运筹学案例分析发电机经济调度论文目录摘要、关键字、引言 (4)一、问题背景介绍 (4)二、对问题进行数学建模 (5)2.1问题一：某时刻的经济调度 (5)2.2问题二：动态经济调度 (5)2.3问题三：机组组合问题 (6)三、案例分析 (7)3.1曲线拟合 (7)3.2问题一案例分析 (7)3.3问题二案例分析 (9)3.4问题三案例分析 (13)四、研究结论 (15)五、参考文献 (16)摘要针对发电机组功率调度问题，依据发电总成本最小化的评判思想，分别建立数学模型，给出了在某一时刻的最优解，在一段时间中的动态规划最优解，以及考虑到机组启停条件后的一段时间中的动态规划最优解。

关键字：发电机组 动态规划 最优解引言在解决这三个调度问题中，我们建立了所需要的数学模型，并假设所有的约束条件都为线性表达式，建立了目标函数总成本与各个发电机组功率值的线性关系，从而达到了运用运筹学模型建立此问题的目的，同时运用软件编程给出了最优解。

一、问题背景介绍在电力系统的运行过程中，发电机组的出力调度是一个至关重要的问题。

发电经济调度对于保障电网运营安全、保证所输送电能质量、降低实际运营成本、延长发电机组使用年限具有重要意义。

我们对该问题使用运筹学思想进行了分析。

在实际调度中，由于电力输送的特点及发电机组的功率限制，往往需要考虑以下几种限制约束。

1、供求平衡约束由于电能目前难以实现大量储存，发电与用电必须保持一定的平衡，否则会对输送的电能质量及电网安全造成一定影响。

若以j P 表示用户端的功率需求，i P 表示发电机组可以提供的功率，则可将该约束表示为∑∑=jiP P在实际情况下，用电需求∑j P 在一天24个小时中会不断变化，即需不断进行调度以满足约束。

2、机组爬坡约束实际发电机组的发电功率变化有一定的范围，即机组无法在一个固定的时间段内实现功率的无限制增加或减小。

若用wi P 表示机组单位时间内功率变化最大极限，t i ,P 表示第i 台机组在t 时刻的功率，则该约束可表示为wi 1)(t i,t i,P P P ≤-+3、机组启停约束实际发电机组可以实现启停，即在用电需求并不太大时，可以选择关闭某些机组以降低成本；当用电需求激增时，可以选择开启机组以加大出力。

二、对问题进行数学建模：2.1问题一：某时刻的经济调度在某一确定时刻下，发电机组必须全部启动，找到此时刻下各个发电机组功率的最优取值。

1.明确问题 ⑴限制条件 发电机全部启动每时刻发出总功率=负载消耗功率每一台机组满足功率最大值与最小值限制 ⑵目标某一时刻发电的总成本最小 ⑶决策变量每台发电机发出的功率Pi 2.数学模型maxmin )()(min i i i ii i i P P P U P b P a Z ≤≤=+=∑∑负载总量2.2问题二：动态经济调度在某一段时间中，发电机组必须全部启动，找到此段时间中各个发电机组功率在各个时间点中最优取值，实现对所有发电机组的一个动态规划。

1.明确问题 ⑴限制条件 发电机全部启动每时刻发出总功率=负载消耗功率每一台机组满足功率最大值与最小值限制 发电机功率变化的爬坡限制 ⑵目标根据一天的负载量变化实现一天对发电机功率的动态调度，使一天内发电的总成本最小⑶决策变量某一时刻某台发电机发出的功率t i ,P ,（t 为时刻，i 为发电机编号） 2.数学模型wit i t i i t i i t i t i i t t P P P P P U b P a Z Z Y ≤-≤≤=+==+∑∑∑)1(,,max ,min i ,P )(P t )(min 负载总量的取值对于任意一个3.解决方法在这里不妨给出一种近似的优化方法首先，在初始时刻不考虑爬坡限制解出此时刻的最优解，然后根据这一时刻的最优解和爬坡限制，找到下一时刻解得可行数值范围，并求出在此范围内下一时刻的最优解，这样以此类推，找到所有时刻满足爬坡限制的一个近似的优化解，而此方法的初始点可以根据实际选为开始时间点，结束时间点，中间时刻点或者负载功率最大的时间点，等等。

2.3问题三：机组组合问题在某一段时间中，任一发电机组可以开启或者关闭，并且开启与关闭间没有任何时间长度的限制，找到此段时间中各个发电机组功率在各个时间点中最优取值，实现对所有发电机组的一个动态规划。

1.明确问题 ⑴限制条件每台发电机可以根据实际需要开启或者关闭，中间不需要有间隔条件 每时刻发出总功率=负载消耗功率每一台机组满足功率最大值与最小值限制 发电机功率变化的爬坡限制 ⑵目标根据一天的负载量变化实现一天对发电机功率的动态调度，使一天内发电的总成本最小⑶决策变量某一时刻某台发电机发出的功率t i,P （t 为时刻，i 为发电机编号）t i ,S （机组i 在t 时刻关闭，其值为0，否则，其值为1） 2.数学模型)01(''''''0'''t )b '(a Z 'Z minY 0,1,,1,,1,,,max ,min ,i ,i t ,,,t =⨯⨯≤-⨯⨯=≤≤=+=•==+++∑∑∑i t i t i wi t i t i t i t i t i i t i i tti t i t i t i t i S P P P P P P P P P P P U P P S P P ，则令若初始时刻为）（或（负载总量）的取值对于任意一个）（）（）（三、案例分析： ３.１曲线拟合：在曲线拟合的过程中，我们取了能表示曲线大体趋势的六个点，然后用软件进行了拟合，但是，我们必须看到，这个过程不是最优的，因为，运用多项式模型进行拟合，拟合之后是一个7次多项式，而多项式最高次数越高，其对初值变化的敏感性就越高，又因为各个测量数据都不是十分精确的，因此，差别很小的初值可能会带来趋势明显不同的拟合曲线形态，所以，原则上对所取得六个点使用三阶样条插值多项式拟合的方法应该是最好的，其思想就是在每两个点间都插入一个不同的三次多项式，在保证曲线光滑的前提下进行拟合，这样能减小曲线对初值的敏感性，会使拟合情况更准确。

在实际情况下，我们对普通城市功率的日变化情况取了12点，利用matlab 进行了曲线拟合，结果如下：函数为：197.5266.6x -197.6x 47.4x -5.418x 0.3238x -0.009782x x -0.0001179P 234567++++=拟合后选取24小时的点为：(0,200)(1,87)(2,157)(3,267)(4,355)(5,399)(6,402)(7,379)(8,346)(9,318)(10,302)(11,303) (12,318)(13,340)(14,362)(15,382)(16,399)(17,419)(18,450)(19,498)(20,572)(21,656)(22,702)(23,610)(24,200)３.２问题1案例分析：对于问题1，我们分别利用cplex 和matlab 对程序进行了编写： #include<ilcplex/ilocplex.h> #include<iostream> using namespace std; ILOSTLBEGIN int main() {IloEnv env; try{cout<<"Entering..."<<'\n'; IloModel model(env);IloNumVarArray vars(env);vars.add(IloNumVar(env,20.0,100.0));vars.add(IloNumVar(env,50.0,368.0));vars.add(IloNumVar(env,75.0,389));model.add(IloMinimize(env,0.3*vars[0]+0.2*vars[1]+0.25*vars[2]+10.5));model.add(vars[0]+vars[1]+vars[2]==330);IloCplex cplex(model);if(!cplex.solve()){env.error()<<"Failed to optimize LP."<<'\n';throw(-1);}IloNumArray vals(env);env.out()<<"Solution statas="<<cplex.getStatus()<<'\n';env.out()<<"Solution value="<<cplex.getObjValue()<<'\n';cplex.getValues(vals,vars);env.out()<<"Values="<<vals<<'\n';cout<<"Success!"<<'\n';}catch(IloException& e){cerr<<"Concert exception caught:"<<e<<'\n';}catch(...){cerr<<"Unknown exception caught"<<'\n';}env.end();return 0;}运行结果如下图所示：上述cplex程序仅对1个小时进行了最优解的求解，对于全天24小时我们采用了matlab进行了编程。

利用matlab中yalmip工具编写程序如下：P1=sdpvar (1,1);%P1代表1号机组功率P2=sdpvar (1,1);% P2代表2号机组功率P3=sdpvar (1,1);% P3代表3号机组功率%定义每台机组每小时功率PA=P1;B=P2;C=P3;%将功率赋予数组f=10.5+0.3*A+0.2*B+0.25*C;%求最小开销F=set(20<=A)+set(A<=100)+set(50<=B)+set(B<=368)+set(75<=C)+set(C<=389); %P约束F=F+set(A(1)+B(1)+C(1)==330); %任意小时功率P之和约束solvesdp(F,f);%进行线性规划D=[A,B,C];E=reshape(D,1,3)';double(f);double(E);%输出结果运行结果为：P1=[20]P2=[235]P3=[75]Min f=82.25３.３问题2案例分析：对于问题2，同样编写cplex程序如下：#include<ilcplex/ilocplex.h>#include<iostream>using namespace std;ILOSTLBEGINint main(){IloEnv env;try{cout<<"Entering..."<<'\n';IloModel model(env);IloNumVarArray vars(env);vars.add(IloNumVar(env,20.0,100.0));vars.add(IloNumVar(env,50.0,368.0));vars.add(IloNumVar(env,75.0,389));//定义第一个小时的三台机组对应的三个变量vars[0][1][2]vars.add(IloNumVar(env,20.0,100.0));vars.add(IloNumVar(env,50.0,368.0));vars.add(IloNumVar(env,75.0,389));//定义第二个小时的三台机组对应的三个变量vars[3][4][5]vars.add(IloNumVar(env,20.0,100.0));vars.add(IloNumVar(env,50.0,368.0));vars.add(IloNumVar(env,75.0,389));//定义第三个小时的三台机组对应的三个变量vars[6][7][8]vars.add(IloNumVar(env));//定义第一个小时的最优解z1vars.add(IloNumVar(env));//定义第二个小时的最优解z2vars.add(IloNumVar(env));//定义第三个小时的最优解z3model.add(IloMinimize(env,vars[9]+vars[10]+vars[11]));//目标函数的意义：使z1+z2+z3最小model.add(vars[9]==10.5+0.3*vars[0]+0.2*vars[1]+0.25*vars[2]);//第一个小时里的目标函数z1=0.3x1+0.2x2+0.25x3+10.5,在这里相当于约束model.add(vars[10]==10.5+0.3*vars[3]+0.2*vars[4]+0.25*vars[5]);//第二个小时里的目标函数z2=0.3x4+0.2x5+0.25x6+10.5,在这里相当于约束model.add(vars[11]==10.5+0.3*vars[6]+0.2*vars[7]+0.25*vars[8]); //第三个小时里的目标函数z3=0.3x7+0.2x8+0.25x9+10.5,在这里相当于约束model.add(vars[0]+vars[1]+vars[2]==200);//第一个小时的功率需求200与三台机组提供功率平衡约束model.add(vars[3]+vars[4]+vars[5]==160);//第二个小时的功率需求160与三台机组提供功率平衡约束model.add(vars[6]+vars[7]+vars[8]==160);//第三个小时的功率需求160与三台机组提供功率平衡约束model.add(vars[3]-vars[0]<=30);model.add(vars[3]-vars[0]>=-30);//一号机组在第二个小时的爬坡约束model.add(vars[6]-vars[3]<=30);model.add(vars[6]-vars[3]>=-30);//一号机组在第三个小时的爬坡约束model.add(vars[4]-vars[1]<=30);model.add(vars[4]-vars[1]>=-30);//二号机组在第二个小时的爬坡约束model.add(vars[7]-vars[4]<=30);model.add(vars[7]-vars[4]>=-30);//二号机组在第三个小时的爬坡约束model.add(vars[5]-vars[2]<=30);model.add(vars[5]-vars[2]>=-30);//三号机组在第二个小时的爬坡约束model.add(vars[8]-vars[5]<=30);model.add(vars[8]-vars[5]>=-30);//三号机组在第三个小时的爬坡约束IloCplex cplex(model);if(!cplex.solve()){env.error()<<"Failed to optimize LP."<<'\n';throw(-1);}IloNumArray vals(env);env.out()<<"Solution statas="<<cplex.getStatus()<<'\n';env.out()<<"Solution value="<<cplex.getObjValue()<<'\n';cplex.getValues(vals,vars);env.out()<<"Values="<<vals<<'\n';cout<<"Success!"<<'\n';}catch(IloException& e){cerr<<"Concert exception caught:"<<e<<'\n';}catch(...){cerr<<"Unknown exception caught"<<'\n';}env.end();return 0;}运行结果如下图所示：同样的，上述cplex程序仅进行了对3个小时内最优解的求解。