本程序参考了3篇电力系统TOP-EI期刊，将目前的热点研究数据驱动结合综合综合能源系统，利用分布鲁棒处理不确定性问题，程序中算例丰富，注释清晰，干货满满，小编非常推荐这个程序，下面对程序代码详细介绍！
参考文献：
1. 中国电机工程学报：计及需求响应柔性调节的分布鲁棒 DG 优化配置
2. 电力系统自动化：基于场景聚类的主动配电网分布鲁棒综合优化
3. 电网技术：基于数据驱动的交直流配电网分布鲁棒优化调度
论文亮点：
1. 通过基于多面体的线性化技巧和McCormick方法，将原始混合整数非线性模型线性化，建立混合整数线性规划模型；然后，充分利用需求响应和新能源出力的历史数据，构建数据驱动的两阶段分布鲁棒新能源优化配置模型，其中第一阶段确定新能源的投资方案，第二阶段模拟投入新能源后的系统运行， 并同时考虑不确定性概率分布置信集合的1-范数和∞-范数约束；
2. 在考虑新能源间相关性的基础上，基于原始数据采用拉丁超立方采样得到初始场景，利用手肘法确定K-means算法的聚类数，从而得到典型场景，场景分布的概率置信区间由1-范数和∞-范数约束。
3. 通过二阶锥松弛和McCormick线性化方法，将原始混合整数非线性模型(MINLP)转化为混合整数二阶锥凸优化模型(MISOCP)；然后，结合风光和负荷的典型历史数据以及决策变量的调节特性，构建数据驱动的两阶段分布鲁棒优化调度模型，并综合1-范数和∞-范数同时约束不确定性概率分布置信集合；
第一阶段目标函数：一阶段目标是微燃机启停成本
第二阶段目标函数：分别是电网购电成本、发电机运行成本、蓄电池老化成本（或者运维）、弃风和弃光成本。
程序结果：
部分程序：
%决策变量p_ch=sdpvar(24,4);%蓄电池充电功率p_dis=sdpvar(24,4);%蓄电池放电功率yita=sdpvar(1,1);%两阶段鲁棒优化中的主问题中目标+yita%热储能h_ch=sdpvar(24,4);h_dis=sdpvar(24,4);uh_ch=binvar(24,4);uh_dis=binvar(24,4);%充放热状态01变量%负荷需求响应变量pdz=sdpvar(24,4);%可转移电负荷pdd=sdpvar(24,4);%可中断电负荷phd=sdpvar(24,4);%可消减热负荷%不确定性变量psp=sdpvar(1,4);%ps+，参见高海淑文献公式54，下同];p0_pv=100*repmat([10,10,6,6],24,1).*[0.001605989 0.000643342 0.002270376 0.0011351880.001911248 0.000946739 0.00205834 0.001029170.000884477 0.000231715 0.002484428 0.001242214p_h=[228,228,228,228;240,240,240,240;240,240,240,240;252,252,252,252;264,264,264,264;276,276,276,276;252,252,252,252;228,228,228,228;216,216,216,216;204,204,204,204;180,180,180,180;144,144,144,144;144,144,144,144;132,132,132,132;108,108,108,108;108,108,108,108;120,120,120,120;132,132,132,132;156,156,156,156;168,168,168,168;192,192,192,192;204,204,204,204;216,216,216,216;228,228,228,228];price=[0.48 0.48 0.48 0.48 0.48 0.48 0.48 0.9 1.35 1.35 1.35 0.9 0.9 0.9 0.9 0.9 0.9 0.9 1.35 1.35 1.35 1.35 1.35 0.48];%电价sell=0.8.*price;hchmax=50;sh=300;yitahch=0.85;yitahdis=0.9; %储热atse=0.1;aie=0.1;bie=0.1; %调整幅度constraints=[];%约束开始load=p_l';%约束lam1=sdpvar(96,1);lam11=sdpvar(96,1);lam12=sdpvar(96,1);lam120=sdpvar(96,1);lam2=sdpvar(192,1);lam21=sdpvar(192,1);lam3=sdpvar(4,1);lam4=sdpvar(96,1);lam41=sdpvar(96,1);lam5=sdpvar(96,1);lam51=sdpvar(96,1);lam6=sdpvar(96,1);%设定对偶变量beta1=sdpvar(96,1);beta11=sdpvar(96,1);beta12=sdpvar(96,1);beta120=sdpvar(96,1);beta2=sdpvar(192,1);beta21=sdpvar(192,1);beta3=sdpvar(4,1);beta4=sdpvar(96,1);beta41=sdpvar(96,1);beta5=sdpvar(96,1);beta51=sdpvar(96,1);beta6=sdpvar(96,1);beta7=binvar(672,1);beta8=binvar(288,1);%大m条件中的01变量x=[];P=[];pw=[];pv=[];for i=1:4x=[x,p_buy(:,i)' p_g(:,i)' p_ch(:,i)' p_dis(:,i)' p_gl(:,i)' p_wt(:,i)' p_pv(:,i)'];%x变量P=[P,ps(i).*price ps(i)*cg.*ones(1,24) ps(i)*ccn.*ones(1,48) zeros(1,24) -ps(i)*cq.*ones(1,48)];pw=[pw,ps(i).*p0_wt(:,i)'];pv=[pv,ps(i).*p0_pv(:,i)'];endx=x';%下面相关的计算参数和参考资料一致%B=repmat([c_wt_om.*ones(1,24) c_pv_om.*ones(1,24) zeros(1,24)]',1,4);Q1=[zeros(24,24) eye(24) zeros(24,120) zeros(24,504);zeros(24,168) zeros(24,24) eye(24) zeros(24,120) zeros(24,336);zeros(24,336) zeros(24,24) eye(24) zeros(24,120) zeros(24,168);zeros(24,504) zeros(24,24) eye(24) zeros(24,120)];lint=zeros(tn);for i=1:tn-1lint(i,i)=-1;lint(i,i+1)=1;endlint(tn,tn)=-1;lint(tn,1)=1;Q12=[zeros(24,24) lint zeros(24,120) zeros(24,504);zeros(24,168) zeros(24,24) lint zeros(24,120) zeros(24,336);zeros(24,336) zeros(24,24) lint zeros(24,120) zeros(24,168);zeros(24,504) zeros(24,24) lint zeros(24,120)];Q2=[zeros(48,48) eye(48) zeros(48,72) zeros(48,504);zeros(48,168) zeros(48,48) eye(48) zeros(48,72) zeros(48,336);zeros(48,336) zeros(48,48) eye(48) zeros(48,72) zeros(48,168);zeros(48,504) zeros(48,48) eye(48) zeros(48,72)];Q3=[zeros(1,48) eta.*ones(1,24) -1/eta.*ones(1,24) zeros(1,72) zeros(1,504);zeros(1,168) zeros(1,48) eta.*ones(1,24) -1/eta.*ones(1,24) zeros(1,72) zeros(1,336);zeros(1,336) zeros(1,48) eta.*ones(1,24) -1/eta.*ones(1,24) zeros(1,72) zeros(1,168);zeros(1,504) zeros(1,48) eta.*ones(1,24) -1/eta.*ones(1,24) zeros(1,72)];% Q4=[eye(48) zeros(48,72) zeros(48,360);zeros(48,120) eye(48) zeros(48,72) zeros(48,240);zeros(48,240) eye(48) zeros(48,72) zeros(48,120);zeros(48,360) eye(48) zeros(48,72)];T5=[p_l(:,1)'+pdz(:,1)'-pdd(:,1)'+p_sell(:,1)' p_l(:,2)'+pdz(:,2)'-pdd(:,2)'+p_sell(:,2)' p_l(:,3)'+pdz(:,3)'-pdd(:,3)'+p_sell(:,3)' p_l(:,4)'+pdz(:,3)'-pdd(:,4)'+p_sell(:,4)']';T6=[p_h(:,1)'-phd(:,1)'+h_ch(:,1)'-h_dis(:,1)' p_h(:,2)'-phd(:,2)'+h_ch(:,2)'-h_dis(:,2)' p_h(:,3)'-phd(:,3)'+h_ch(:,3)'-h_dis(:,3)' p_h(:,4)'-phd(:,4)'+h_ch(:,4)'-h_dis(:,4)']';%增加原始约束constraints=[constraints,Q1*x<=Pgmax];constraints=[constraints,Q1*x>=0];constraints=[constraints,Q6*x==T6];constraints=[constraints,0<=h_ch<=hchmax.*uh_ch, 0<=h_dis<=hchmax.*uh_dis, uh_ch+uh_dis<=1];constraints=[constraints,0.1*sh<=0.5*sh+sum(yitahch*h_ch(1:t,j)-h_dis(1:t,j)/yitahdis)<=0.9*sh];%负荷需求响应constraints=[constraints,-atse.*p_l<=pdz<=atse.*p_l,sum(pdz)==0];constraints=[constraints,0<=pdd<=aie.*p_l,0<=phd<=bie.*p_h];constraints=[constraints,0<=p_sell<=20];%kkt% constraints=[constraints,Q1'*lam1-Q1'*lam11+Q12'*lam12-Q12'*lam120+Q2'*lam2-Q2'*lam21+Q3'*lam3+Q4'*lam4-Q4'*lam41+Q5'*lam5+Q6'*lam6>=-P'];% %constraints=[constraints,0<=B+G'*lam6];constraints=[constraints,y_u(j,i)<=uu_m(j,i),y_u(j,i)<=uu_m(j+1,i),y_u(j,i)>=uu_m(j+1,i)+uu_m(j,i)-1];cdz=0.5;cdd=0.5;chd=0.05; %补偿成本c1=sum(price.*p_buy(:,1)'+cg.*p_g(:,1)'+ccn.*p_ch(:,1)'+ccn.*p_dis(:,1)'+cq.*p0_wt(:,1)'-cq.*p_wt(:,1)'+cq.*p0_pv(:,1)'-cq.*p_pv(:,1)'-sell.*p_sell(:,1)'+cdz.*pdz(:,1)'+cdd.*pdd(:,1)'+chd.*phd(:,1)')+0.13*sum(0.82*p_buy(:,1)'+0.49*p_g(:,1)')+0.011*sum(h_ch(:,1)'+p_dis(:,1)');c2=sum(price.*p_buy(:,2)'+cg.*p_g(:,2)'+ccn.*p_ch(:,2)'+ccn.*p_dis(:,2)'+cq.*p0_wt(:,2)'-cq.*p_wt(:,2)'+cq.*p0_pv(:,2)'-cq.*p_pv(:,2)'-sell.*p_sell(:,2)'+cdz.*pdz(:,2)'+cdd.*pdd(:,2)'+chd.*phd(:,2)')+0.13*sum(0.82*p_buy(:,2)'+0.49*p_g(:,2)')+0.011*sum(h_ch(:,2)'+p_dis(:,2)');c3=sum(price.*p_buy(:,3)'+cg.*p_g(:,3)'+ccn.*p_ch(:,3)'+ccn.*p_dis(:,3)'+cq.*p0_wt(:,3)'-cq.*p_wt(:,3)'+cq.*p0_pv(:,3)'-cq.*p_pv(:,3)'-sell.*p_sell(:,3)'+cdz.*pdz(:,3)'+cdd.*pdd(:,3)'+chd.*phd(:,3)')+0.13*sum(0.82*p_buy(:,3)'+0.49*p_g(:,3)')+0.011*sum(h_ch(:,3)'+p_dis(:,3)');c4=sum(price.*p_buy(:,4)'+cg.*p_g(:,4)'+ccn.*p_ch(:,4)'+ccn.*p_dis(:,4)'+cq.*p0_wt(:,4)'-cq.*p_wt(:,4)'+cq.*p0_pv(:,4)'-cq.*p_pv(:,4)'-sell.*p_sell(:,4)'+cdz.*pdz(:,4)'+cdd.*pdd(:,4)'+chd.*phd(:,4)')+0.13*sum(0.82*p_buy(:,4)'+0.49*p_g(:,4)')+0.011*sum(h_ch(:,4)'+p_dis(:,4)');Cz=[c1,c2,c3,c4]';obj_mp=double(obj_mp);p_ch=double(p_ch);p_dis=double(p_dis);p_buy=double(p_buy);p_g=double(p_g);h_dis=value(h_dis);% c1=sum(price.*p_buy(:,1)'+cg.*p_g(:,1)'+ccn.*p_ch(:,1)'+ccn.*p_dis(:,1)'+cq.*p0_wt(:,1)'-cq.*p_wt(:,1)'+cq.*p0_pv(:,1)'-cq.*p_pv(:,1)'-sell.*p_sell(:,1)'+cdz.*pdz(:,1)'+cdd.*pdd(:,1)'+chd.*phd(:,1)')+0.13*sum(0.82*p_buy(:,1)'+0.49*p_g(:,1)')+0.011*sum(h_ch(:,1)'+p_dis(:,1)');% c2=sum(price.*p_buy(:,2)'+cg.*p_g(:,2)'+ccn.*p_ch(:,2)'+ccn.*p_dis(:,2)'+cq.*p0_wt(:,2)'-cq.*p_wt(:,2)'+cq.*p0_pv(:,2)'-cq.*p_pv(:,2)'-sell.*p_sell(:,2)'+cdz.*pdz(:,2)'+cdd.*pdd(:,2)'+chd.*phd(:,2)')+0.13*sum(0.82*p_buy(:,2)'+0.49*p_g(:,2)')+0.011*sum(h_ch(:,2)'+p_dis(:,2)');% c3=sum(price.*p_buy(:,3)'+cg.*p_g(:,3)'+ccn.*p_ch(:,3)'+ccn.*p_dis(:,3)'+cq.*p0_wt(:,3)'-cq.*p_wt(:,3)'+cq.*p0_pv(:,3)'-cq.*p_pv(:,3)'-sell.*p_sell(:,3)'+cdz.*pdz(:,3)'+cdd.*pdd(:,3)'+chd.*phd(:,3)')+0.13*sum(0.82*p_buy(:,3)'+0.49*p_g(:,3)')+0.011*sum(h_ch(:,3)'+p_dis(:,3)');% c4=sum(price.*p_buy(:,4)'+cg.*p_g(:,4)'+ccn.*p_ch(:,4)'+ccn.*p_dis(:,4)'+cq.*p0_wt(:,4)'-cq.*p_wt(:,4)'+cq.*p0_pv(:,4)'-cq.*p_pv(:,4)'-sell.*p_sell(:,4)'+cdz.*pdz(:,4)'+cdd.*pdd(:,4)'+chd.*phd(:,4)')+0.13*sum(0.82*p_buy(:,4)'+0.49*p_g(:,4)')+0.011*sum(h_ch(:,4)'+p_dis(:,4)');% Cz=[c1,c2,c3,c4]';c1=value(c1);c2=value(c2);c3=value(c3);c4=value(c4);Cz=value(Cz);yita=value(yita);qw=ps(1)*sum(p0_wt(:,1)'-p_wt(:,1)')/sum(p0_wt(:,1))+ps(2)*sum(p0_wt(:,2)'-p_wt(:,2)')/sum(p0_wt(:,2))+ps(3)*sum(p0_wt(:,3)'-p_wt(:,3)')/sum(p0_wt(:,3))+ps(4)*sum(p0_wt(:,4)'-p_wt(:,4)')/sum(p0_wt(:,4));qv=ps(1)*sum(p0_pv(:,1)'-p_pv(:,1)')/sum(p0_pv(:,1))+ps(2)*sum(p0_pv(:,2)'-p_pv(:,2)')/sum(p0_pv(:,2))+ps(3)*sum(p0_pv(:,3)'-p_pv(:,3)')/sum(p0_pv(:,3))+ps(4)*sum(p0_pv(:,4)'-p_pv(:,4)')/sum(p0_pv(:,4));% c4buy=sum(price.*p_buy(:,4)');% a=a1+a2+a3+a4;
部分内容源自网络，侵权联系删除！
欢迎感兴趣的小伙伴点击文末阅读原文获取完整版代码，小编会不定期更新高质量的学习资料、文章和程序代码，为您的科研加油助力！
更多创新性预测及优化程序请点击公众号首页：《预测优化》|《综合能源》链接！
