单向/双向V2G环境下分布式电源与电动汽车充电站联合配置方法(matlab代码)

T = 96;%时段数为1小时
nb = 31;%节点数,根节点为31
nl = 30;%支路数
P1=1e-5.*[500 200 250 180 180 120 240 240 240 240 240 160 160 400 320 300 300 300 300 800 0 80 200 100 60 40 500 200 80 30 0]';
Q1=1e-5.*[300 90 160 90 60 60 80 80 80 120 120 90 100 75 100 125 100 100 100 320 0 40 100 40 25 15 160 120 30 10 0]';
load pc_jm_w;
load pc_jm_wd;
load px_jm_w;
load px_jm_wd;
load pq_jm_w;
load pq_jm_wd;
load pd_jm_w;
load pd_jm_wd;
load pc_sc_w;
load pc_sc_wd;
load px_sc_w;
load px_sc_wd;
load pq_sc_w;
load pq_sc_wd;
load pd_sc_w;
load pd_sc_wd;
load pc_bg_w;
load pc_bg_wd;
load px_bg_w;
load px_bg_wd;
load pq_bg_w;
load pq_bg_wd;
load pd_bg_w;
load pd_bg_wd;
%电动汽车
prl=100;%容量
pev=30*0.25;
cdz=[7 30 60;400 3250 5600;40 325 560;10 10 10];%充电桩数据
cr=106.5;%单位容量增容成本
cl=80;%单位网损费用
plmax=600*1e-4;%线路最大电流
d=0.03;%折现率
num_peak=[1 2 5 1.*ones(1,28)];%测试数据，下为正式数据，运行速度太慢
%num_peak=[30 12 15 11 11 7 14 14 14 14 14 10 10 24 19 18 18 18 18 48 0 5 12 6 4 2 30 12 5 2 0]';%峰值电动汽车停车数量
load arr_jm_w;
load arr_jm_wd;
load arr_sc_w;
load arr_sc_wd;
load arr_bg_w;
load arr_bg_wd;
load stay_jm_w;
load stay_jm_wd;
load stay_sc_w;
load stay_sc_wd;
load stay_bg_w;
load stay_bg_wd;
sty_jd=[2 2 1 3 3 2 1 1 1 1 1 1 1 2 2 2 2 1 1 2 1 3 2 2 3 3 3 1 2 3 1];%节点区域类型
%构建负荷矩阵
for i=1:nb
   if sty_jd(i)==1
       pload(i,1:T)=P1(i).*pc_jm_w;%春 工作日
       pload(i,T+1:2*T)=P1(i).*pc_jm_wd;%春 周末
       pload(i,2*T+1:3*T)=P1(i).*px_jm_w;%夏 工作日
       pload(i,3*T+1:4*T)=P1(i).*px_jm_wd;%夏 周末
       pload(i,4*T+1:5*T)=P1(i).*pq_jm_w;%秋 工作日
       pload(i,5*T+1:6*T)=P1(i).*pq_jm_wd;%秋 周末
       pload(i,6*T+1:7*T)=P1(i).*pd_jm_w;%冬 工作日
       pload(i,7*T+1:8*T)=P1(i).*pd_jm_wd;%冬 周末
   elseif sty_jd(i)==2
       pload(i,1:T)=P1(i).*pc_sc_w;%春 工作日
       pload(i,T+1:2*T)=P1(i).*pc_sc_wd;%春 周末
       pload(i,2*T+1:3*T)=P1(i).*px_sc_w;%夏 工作日
       pload(i,3*T+1:4*T)=P1(i).*px_sc_wd;%夏 周末
       pload(i,4*T+1:5*T)=P1(i).*pq_sc_w;%秋 工作日
       pload(i,5*T+1:6*T)=P1(i).*pq_sc_wd;%秋 周末
       pload(i,6*T+1:7*T)=P1(i).*pd_sc_w;%冬 工作日
       pload(i,7*T+1:8*T)=P1(i).*pd_sc_wd;%冬 周末
   else
       pload(i,1:T)=P1(i).*pc_bg_w;%春 工作日
       pload(i,T+1:2*T)=P1(i).*pc_bg_wd;%春 周末
       pload(i,2*T+1:3*T)=P1(i).*px_bg_w;%夏 工作日
       pload(i,3*T+1:4*T)=P1(i).*px_bg_wd;%夏 周末
       pload(i,4*T+1:5*T)=P1(i).*pq_bg_w;%秋 工作日
       pload(i,5*T+1:6*T)=P1(i).*pq_bg_wd;%秋 周末
       pload(i,6*T+1:7*T)=P1(i).*pd_bg_w;%冬 工作日
       pload(i,7*T+1:8*T)=P1(i).*pd_bg_wd;%冬 周末
   end
end
qload=repmat(Q1,1,8*T);
%根据pre_ev生成工作日和周末的电动汽车矩阵 djz-工作日可调度；djzd-周末可调度；前面加n是不可调度
djz=[65,16,0.217993798750154,16,32;6,5,0.0196210667655841,19,42;80,19,0.520129415310742,19,17;10,28,0.108016694136759,28,50;15,17,0.00457962394732347,16,74;82,28,0.986968274783658,28,15;27,3,0.507848830829537,19,56;73,3,0.661596193082714,19,24;17,29,0.590483177142572,19,42;52,17,0.680065530083361,16,35;24,18,0.866887054672508,19,39;12,14,0.300184401213900,1,39;80,13,0.390175938130607,1,17;14,8,0.0868151008658337,11,41;25,9,0.424858411704626,16,11;48,22,0.243573372680951,27,49;8,12,0.00339412296430741,1,21;1,6,0.142484054820971,19,26;18,4,0.598885610367865,19,79;90,7,0.482671375396414,1,7;51,8,0.0683572204708292,11,42;18,1,0.954678274080450,1,41;68,20,0.552309844507730,11,21;74,7,0.370864718973405,1,23;82,12,0.318019117496735,1,15;87,28,0.591594408907438,28,10;82,14,0.904355478217944,1,3;52,22,0.179301843821789,27,32;19,10,0.403856711230207,11,53;6,17,0.274811404837501,16,23;24,5,0.956416361008024,19,73;79,23,0.175811728275530,25,18;19,0,0.316419513732543,28,67;60,17,0.439037203387669,16,28;49,24,0.762408048788714,25,48;72,20,0.123219518324506,11,25;34,3,0.147849468032519,19,19;65,13,0.694403909684554,1,25;2,17,0.279391965252965,16,91;81,30,0.0569328854352480,19,16;56,21,0.719432752595998,27,41;70,12,0.581582083212866,1,11;6,30,0.284823726860618,19,57;39,12,0.131114707043005,1,42;10,19,0.0109790922908117,19,55;76,7,0.448019713464051,1,21;7,15,0.642315234585299,11,21;81,30,0.846372887693173,19,16;27,23,0.236930384112870,25,70;25,27,0.911067053397794,27,67;70,7,0.576053456321354,1,27;39,31,0.0899988149868883,16,31;50,2,0.725687923545844,16,47;51,26,0.858759034071804,25,46;31,14,0.752227970049942,1,11;11,8,0.524637345396311,11,86;68,10,0.291457127647727,11,29;81,18,0.948108735396022,19,6;57,9,0.827732173448263,16,18;43,12,0.826573979042765,1,54;21,10,0.133810985356126,11,64;55,5,0.147655777151737,19,42;77,14,0.382646229559959,1,20;36,17,0.711656705981267,16,61;32,20,0.974836148002758,11,7;57,13,0.309136426466267,1,25;73,31,0.186571444141369,16,24;20,31,0.802261569764289,16,41;70,15,0.808990267199045,11,27;8,18,0.910187830728141,19,19;42,23,0.0391844866475830,25,55;74,17,0.183842944465775,16,23;50,31,0.854851683090675,16,47;65,13,0.934979086178052,1,32;23,12,0.705077475112789,1,54;73,31,0.962431404366634,16,24;56,16,0.901208092653144,16,41;77,4,0.478474472902863,19,20;36,21,0.169608813454314,27,27;20,6,0.326839648349976,19,77;46,13,0.179231476423447,1,51;40,26,0.615325097130850,25,36;84,24,0.464954283314534,25,13;86,13,0.334329419622220,1,11;87,22,0.377455102070865,27,10;8,6,0.0929889268706780,19,44;4,6,0.454966145002097,19,12];
djzd=[2,23,0.354116466056586,25,75;10,18,0.241084055169022,19,81;82,30,0.488899786160604,19,15;22,17,0.762109709211147,16,33;45,20,0.917336040866845,11,16;69,18,0.433298977209119,19,28;38,6,0.633333580320122,19,59;32,25,0.999477858635892,25,65;13,7,0.0236324666582251,1,58;12,13,0.884076806661962,1,53;36,6,0.440943276524324,19,61;13,15,0.856896327782193,11,4;67,30,0.283267898492137,19,13;66,28,0.610868982383243,28,31;19,23,0.346260716664775,25,40;16,25,0.624923527846277,25,71;77,2,0.950790316505566,16,20;73,23,0.831129575710113,25,15;44,19,0.932183355047051,19,53;86,18,0.582746806361791,19,11;72,5,0.143908328674868,19,25;25,10,0.401790917846999,11,39;38,19,0.166890722666246,19,18;10,10,0.769596929272335,11,22;71,21,0.824078213366380,27,26;29,10,0.523029798441993,11,31;80,25,0.556998261278609,25,17;66,7,0.456425364613112,1,31;52,31,0.755220401677830,16,45;23,16,0.0514360535669078,16,73;58,27,0.988277139258114,27,39;40,0,0.540878322763043,28,20;22,10,0.0959494166075845,11,72;72,17,0.338132316524053,16,25;53,30,0.892833319318002,19,34;53,11,0.622802790659130,11,44;40,23,0.781374000275963,25,35;72,28,0.242603386279670,28,12;22,11,0.287084613146176,11,75;6,18,0.162898900217732,19,80;17,16,0.999329479025719,16,34;5,7,0.397839124102705,1,32;23,29,0.683188784699983,19,74;43,29,0.00583432457237088,19,54;77,7,0.932468675920979,1,20;80,18,0.792581664120226,19,17;22,10,0.584523478528556,11,75;29,12,0.862057298433146,1,59;25,15,0.399075226726792,11,58;77,3,0.821442194126392,19,20;35,13,0.572239218412685,1,62;72,23,0.389128763569227,25,25;92,18,0.849722054218674,19,5;60,18,0.963468465806992,19,8;49,16,0.0901660021689635,16,48;85,14,0.781722612917166,1,12;60,8,0.445656211439575,11,37;20,9,0.483294567734286,16,32;77,31,0.159047554722395,16,20;68,12,0.973704902674984,1,29;62,27,0.401883398008086,27,35;95,17,0.933591915842209,16,2;47,20,0.887636865889604,11,19;39,31,0.402351616935948,16,58;87,31,0.653163281094690,16,10;4,19,0.567144361897996,19,92;72,21,0.523313318412010,27,25;50,2,0.890036233228213,16,32;23,4,0.310922713453606,19,22;63,2,0.275431379860917,16,27;85,14,0.755914120907961,1,12;75,4,0.978563885159978,19,22;6,14,0.325653274948034,1,60;23,18,0.603156320439897,19,58;44,1,0.513814835038277,1,39];
ndjz=[13,21,0.599585548122814,27,5;90,20,0.451945709260386,11,7;92,24,0.00734866110284693,25,5;87,12,0.0248552338448721,1,10;92,6,0.193039815977027,19,5;90,12,0.273216707999963,1,7;60,19,0.172604501647669,19,9;88,20,0.255370297944443,11,9;87,17,0.0329398927498766,16,5;93,4,0.0514482932302053,19,4;42,17,0.712414805789522,16,2;92,17,0.540105832359366,16,5;2,28,0.642741739133104,28,0];
ndjzd=[42,14,0.0492131803564864,1,5;4,27,0.407730811355982,27,3;72,4,0.822394006759020,19,2;95,6,0.827209084946655,19,2;92,17,0.0302701634197687,16,5];
[num_w,~]=size(djz);
[num_wd,~]=size(djzd);
%% 1.设参
mpc = real31;
%pload = mpc.Pload(:,t)*(pload1(t)/sum(mpc.Pload(:,t)))/100;%节点有功负荷
%qload = mpc.Qload(:,1);%节点无功负荷
branch = mpc.branch;
branch(:,3) = branch(:,3)*100/(12.66^2);%求阻抗标幺值
r=real(branch(:,3));
x=imag(branch(:,3));
r=r(1:30);
x=x(1:30);
upstream=zeros(nb,nl);%代表流入节点支路
dnstream=zeros(nb,nl);%代表流出节点支路
for i=1:nl
    upstream(i,i)=1;
end
for i=[1:5,7:8,21:25]
    dnstream(i,i+1)=1;
end
dnstream(1,[7 10 12 14 15])=1;
dnstream(10,11)=1;
dnstream(12,13)=1;
dnstream(2,[16 18 20])=1;
dnstream(16,17)=1;
dnstream(18,19)=1;
dnstream(22,[27 28])=1;
dnstream(23,[29 30])=1;
dnstream(31,1)=1;
dnstream(31,21)=1;
Vmax=[1.05*1.05*ones(30,1);1.05*1.05*ones(1,1)];
Vmin=[0.95*0.95*ones(30,1);1.05*1.05*ones(1,1)];
Pgmax=[zeros(30,1);ones(1,1)];
Qgmax=[zeros(30,1);ones(1,1)];
%定义变量
V = sdpvar(nb,8*T);%电压的平方
I = sdpvar(nl,8*T);%电流的平方
P = sdpvar(nl,8*T);%线路有功
Q = sdpvar(nl,8*T);%线路无功
Pg = sdpvar(nb,8*T);%发电机有功
Qg = sdpvar(nb,8*T);%发电机无功
Ns = intvar(nb,1);%充电桩数量
Npv = intvar(8,1);%光伏节点安装数量
Ng = intvar(6,1);%燃气轮机安装数量
evc = binvar(num_w,T);%工作日电动汽车充电状态变量
evd = binvar(num_w,T);%工作日电动汽车放电状态变量
evdc = binvar(num_wd,T);%周末电动汽车充电状态变量
evdd = binvar(num_wd,T);%周末电动汽车放电状态变量