随机自适应鲁棒方法！基于恶劣场景辨别法的微网随机自适应鲁棒模型程序代码！

前言

随着分布式电源(distributed generator，DG)接入电网比例的不断提高，其出力的随机性和间歇性对电力系统的安全稳定运行构成威胁。微网凭借先进的智能计量技术、协调控制技术以及信息通信技术，为有效解决DG并网提供了新的思路。然而，由于电价的波动性和风光出力的随机性，不平衡惩罚依然存在，微网参与市场竞标仍存在很大的风险。微网往往聚合可再生能源、储能、需求响应(demand response，DR)等多种分布式能源，通过协调内部各机组出力，实现作为一个有机整体参与电网运行与调度，极大减小了分布式电源单独并网对公网造成的冲击，并提高了其市场竞争力。多阶段鲁棒优化调度模型能够在仅知道不确定性因素边界值的基础下，严格保证极端条件下的系统运行可靠性，并通过动态迭代降低保守性，克服了静态鲁棒优化的缺点。

恶劣场景概念及现实意义

从现实意义上来讲，所谓最恶劣场景，就是使得微网实时调度收益最低的场景，例如，当光伏出力波动较大，在本该出力较高的10－16（一般认为）时，由于天气等原因，此时光伏出力值较低，为了满足较高的负荷需求，微网在此时需要增加向电网的购电电量或增开燃气轮机组等其余发电设备，相较于成本较低廉的光伏发电，该举措使得微网的发电成本增加，从而其收益有所降低。

事实上，整体目标函数包括两部分，第一阶段基准场景下的日前优化调度收益最大化，由于基准场景一般给定，所以第一阶段为确定性的优化模型。第二阶段考虑了所有可能出现的光伏场景，设该场景集合为Ω，第二阶段的min-max问题其实质含义为，假设Ω集合中共计有ns种光伏场景，则内层max求解的是每一个光伏场景微网的实时调度收益最大化问题，共计有ns种最优解，分别对应不同场景。而外层的min主要用于辨别这ns种场景下使得微网实时收益最低的场景，若假设当取场景si时，对应的微网实时调度收益在所有场景中最低，则可认为此时的si最恶劣场景，也称恶劣场景。

鲁棒性的定义与意义

鲁棒性是指算法模型在面对各种异常情况时，能够保持良好的性能和稳定性。在实际应用中，数据可能存在噪声、缺失、异常值等问题，算法模型需要具备一定的鲁棒性才能在这些情况下正常工作。提高算法模型的鲁棒性可以提高其在真实环境中的适应性和可靠性，为实际问题的解决提供更好的保障。

自适应鲁棒优化求解思路

目前求解自适应鲁棒优化问题主要有4种思路：1)仿射算法；2)Benders分解算法；3)列约束生成法；4)场景法。

首先，仿射算法利用线性决策规则建立决策变量与不确定参数之间的仿射关系，从而将两阶段问题转化为单阶段优化问题，但是结果较为保守;

Benders分解算法以及列约束生成法在线性化过程中会引入大量的边界参数和整数变量，导致问题规模较大时模型求解复杂;

场景法需要列举大量不确定场景集合，导致计算的效率偏低，求解时间过长。

迭代流程图

程序介绍

针对光伏出力的随机性与间歇性，采用动态鲁棒优化法对其进行处理。构建了考虑电价和光伏出力不确定性的微网两阶段鲁棒优化调度模型，并采用恶劣场景辨别算法将原问题分解为主问题和子问题进行迭代求解。子问题用来辨别最恶劣的光伏出力情景，并通过主问题对该情景下的单层优化模型进行求解，从而极大地削减了所需求解情景数量，提高了模型的计算效率。程序中算例丰富，注释清晰，干货满满，创新性和可扩展性很高，足以撑起一篇高水平论文！下面对程序做简要介绍！

程序适用平台：Matlab+Yalmip+Cplex

参考文献：《基于恶劣场景辨别法的微网随机自适应鲁棒模型》-南方电网技术

复现情况：本程序未复现考虑微网同时参与日前能量市场、实时能量市场以及碳交易市场的情况，感兴趣的小伙伴可以自己尝试添加。

程序结果

部分程序

kmp=length(j);P_MP=sdpvar(1,1);​u_GTon=binvar(1,24);
u_GToff=binvar(1,24);u_GT=binvar(1,24);​C_GT=sdpvar(kmp,24);
P_DA=sdpvar(1,24);S_DA=sdpvar(1,24);​P_RT=sdpvar(kmp,24);
S_RT=sdpvar(kmp,24);g_GT=sdpvar(kmp,24);P_D=sdpvar(kmp,24);​
C_GT>=0,​Pmax>=P_DA>=0,
Pmax>=S_DA>=0,Pmax>=P_RT>=0,​Pmax>=S_RT>=0,
g_GT>=0,P_Dmax>=P_D>=0,​
u_GT(t)-u_GT(t-1)<= u_GTon(t),u_GT(t-1)-u_GT(t)<= u_GToff(t),      ​
​2*u_GToff(t)<=2-u_GT(t)-u_GT(t+1);    
2*u_GTon(t)<= u_GT(t)+u_GT(t+1);​
C_GT(s,t)==a*u_GT(t)+k*g_GT(s,t)+SUC*u_GTon(t)+SUD*u_GToff(t),​
g_GT(s,t) <= u_GT(t)*g_GT_max;     ​
​g_GT(s,t)-g_GT(s,t-1)<=RU;  
​g_GT(s,t)-g_GT(s,t-1)>=-RD; 
S_DA(t)+S_RT(s,t) == P_DA(t)+P_RT(s,t)+P_RES(s,t)+g_GT(s,t);    ​
​P_MP<=sum(x_RT.*(S_RT(s,:)-1.1*P_RT(s,:)))-sum(C_GT(s,:));
Obj_MP=Obj_MP+x_DA(t)*(S_DA(t)-1.1*P_DA(t))-SUC*u_GTon(t)-SUD*u_GToff(t);

部分内容源自网络，侵权联系删除！

欢迎感兴趣的小伙伴关注并私信获取完整版代码，小编会不定期更新高质量的学习资料、文章和程序代码，为您的科研加油助力！