(x(2)-k2)-sz*k2; %石油变化率在正负sz%之间 (x(5)-x(2))-sz*x(2); (x(8)-x(5))-sz*x(5); (x(11)-x(8))-sz*x(8); (x(14)-x(11))-sz*x(11);
(x(3)-k3)-qz*k3; %天然气变化率在正负qz%之间 (x(6)-x(3))-qz*x(3); (x(9)-x(6))-qz*x(6); (x(12)-x(9))-qz*x(9); (x(15)-x(12))-qz*x(12);
-(x(1)-k1)-mz*k1; -(x(4)-x(1))-mz*x(1); -(x(7)-x(4))-mz*x(4); -(x(10)-x(7))-mz*x(7); -(x(13)-x(10))-mz*x(10);
-(x(2)-k2)-sz*k2; -(x(5)-x(2))-sz*x(2); -(x(8)-x(5))-sz*x(5); -(x(11)-x(8))-sz*x(8); -(x(14)-x(11))-sz*x(11);
-(x(3)-k3)-qz*k3; -(x(6)-x(3))-qz*x(3); -(x(9)-x(6))-qz*x(6); -(x(12)-x(9))-qz*x(9); -(x(15)-x(12))-qz*x(12); ];
ceq=[];
优化运行函数 x0=[16500.33 4283.73 951.88 16500.33 4283.73 951.88 16500.33 4283.73 951.88 16500.33 4283.73 951.88 16500.33 4283.73 951.88 ];
lb=zeros(1,15);
[x,fval,exitflag,output]=fmincon('objfun',x0,[],[],[],[],lb,[],'confun')
附录四
问题三优化程序： 目标函数
function f=objfun31(x)
19
eng2010=25774; eng2011=28051; eng2012=30008; eng2013=32739; eng2014=35752; eng2015=38701; eng2016=45488.0; eng2017=49884.0; eng2018=54705.0; eng2019=59991.8; eng2020=64789.6;
d1=[42.32 46.87 51.9 57.47 63.65 ]; d2=[0.79 0.8 0.82 0.83 0.85 ];
d3=[286.38 377.19 496.79 654.32 861.8];
f=-( (d1(1).*x(1)+d2(1).*x(2)+d3(1).*x(3))/ eng2011+ (d1(2).*x(4)+d2(2).*x(5)+d3(2).*x(6))/
eng2012+(d1(3).*x(7)+d2(3).*x(8)+d3(3).*x(9))/ eng2013+(d1(4).*x(10)+d2(4).*x(11)+d3(4).*x(12))/
eng2014+(d1(5).*x(13)+d2(5).*x(14)+d3(5).*x(15))/ eng2015);
限制条件
function [c,ceq]=confun31(x) gdp2010=41425.48; gdp2011=49175; gdp2012=58020; gdp2013=68456; gdp2014=80769; gdp2015=95297; eng2010=25774; eng2011=28051; eng2012=30008; eng2013=32739; eng2014=35752; eng2015=38701;
mz=0.01; sz=0.01; qz=0.1;
co2010=54571.39;
k1=50.7; k2=22940.4;
20
k3=62.5;
d1=[42.32 46.87 51.9 57.47 63.65 ]; d2=[0.79 0.8 0.82 0.83 0.85 ];
d3=[286.38 377.19 496.79 654.32 861.8]; c= [
gdp2011-(d1(1).*x(1)+d2(1).*x(2)+d3(1).*x(3)); gdp2012-(d1(2).*x(4)+d2(2).*x(5)+d3(2).*x(6)); gdp2013-(d1(3).*x(7)+d2(3).*x(8)+d3(3).*x(9)); gdp2014-(d1(4).*x(10)+d2(4).*x(11)+d3(4).*x(12)); gdp2015-(d1(5).*x(13)+d2(5).*x(14)+d3(5).*x(15));
(x(1)-k1)-mz*k1; (x(4)-x(1))-mz*x(1); (x(7)-x(4))-mz*x(4); (x(10)-x(7))-mz*x(7); (x(13)-x(10))-mz*x(10); (x(2)-k2)-sz*k2; (x(5)-x(2))-sz*x(2); (x(8)-x(5))-sz*x(5); (x(11)-x(8))-sz*x(8); (x(14)-x(11))-sz*x(11); (x(3)-k3)-qz*k3; (x(6)-x(3))-qz*x(3); (x(9)-x(6))-qz*x(6); (x(12)-x(9))-qz*x(9); (x(15)-x(12))-qz*x(12);
-(x(1)-k1)-mz*k1; -(x(4)-x(1))-mz*x(1); -(x(7)-x(4))-mz*x(4); -(x(10)-x(7))-mz*x(7); -(x(13)-x(10))-mz*x(10);
-(x(2)-k2)-sz*k2; -(x(5)-x(2))-sz*x(2); -(x(8)-x(5))-sz*x(5); -(x(11)-x(8))-sz*x(8); -(x(14)-x(11))-sz*x(11);
-(x(3)-k3)-qz*k3; -(x(6)-x(3))-qz*x(3); -(x(9)-x(6))-qz*x(6);
21
-(x(12)-x(9))-qz*x(9); -(x(15)-x(12))-qz*x(12); ];
ceq=[]; 运行
x0=[50.7 22940.4 62.5 50.7 22940.4 62.5 50.7 22940.4 62.5 50.7
22940.4 62.5 50.7 22940.4 62.5 ]; aeq=[
1 1 1 0 0 0 0 0 0 0 0 0 0 0 0; 0 0 0 1 1 1 0 0 0 0 0 0 0 0 0; 0 0 0 0 0 0 1 1 1 0 0 0 0 0 0; 0 0 0 0 0 0 0 0 0 1 1 1 0 0 0; 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1]; beq=[ 28051; 30008; 32739; 35752; 38701];
lb=zeros(1,15);
[x,fval,exitflag,output]=fmincon('objfun31',x0,[],[],aeq,beq,lb,
[],'confun31')
22