sg_gen_dispatch

SG_GEN_DISPATCH  main function to run Pg setting
   [PG_SETTING, CORR_PG] = SG_GEN_DISPATCH(PGMAX, PL_SETTING, REFSYS_STAT)

   Input:
       PgMax - matrix (Ng by 2) indicating the generation bus numbers
           and generation capacity at each generation bus
       PL_setting - matrix (NL by 2) indicating the load bus numbers
           and load size at each load bus
       refsys_stat - statistics of reference system

   Output:
       Pg_setting - matrix (Ng by 2) indicating the generation bus numbers and
           generation dispatch

0001 function [Pg_setting] = sg_gen_dispatch(PgMax, PL_setting, refsys_stat)
0002 %SG_GEN_DISPATCH  main function to run Pg setting
0003 %   [PG_SETTING, CORR_PG] = SG_GEN_DISPATCH(PGMAX, PL_SETTING, REFSYS_STAT)
0004 %
0005 %   Input:
0006 %       PgMax - matrix (Ng by 2) indicating the generation bus numbers
0007 %           and generation capacity at each generation bus
0008 %       PL_setting - matrix (NL by 2) indicating the load bus numbers
0009 %           and load size at each load bus
0010 %       refsys_stat - statistics of reference system
0011 %
0012 %   Output:
0013 %       Pg_setting - matrix (Ng by 2) indicating the generation bus numbers and
0014 %           generation dispatch
0015 
0016 %   SynGrid
0017 %   Copyright (c) 2017-2018, Electric Power and Energy Systems (EPES) Research Lab
0018 %   by Hamidreza Sadeghian and Zhifang Wang, Virginia Commonwealth University
0019 %
0020 %   This file is part of SynGrid.
0021 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0022 Alpha_mod = refsys_stat.Alpha_mod;
0023 mu_committed = refsys_stat.mu_committed ;
0024 Tab_2D_Pg = refsys_stat.Tab_2D_Pg ;
0025 Ng=length(PgMax);
0026 Max=max(PgMax(:,2));
0027 Norm_PgMax=[PgMax(:,1),PgMax(:,2)./Max];
0028 Norm_Total_load=sum(PL_setting(:,2))./Max;
0029 
0030 %% Uncommitted  units (10~20)% of total number of generation units
0031 Ng_Uncommitted=round(Ng*(10+rand*(10))/100);
0032 Norm_Uncommitted_Units=[];
0033 
0034 for i=1:Ng_Uncommitted
0035     % uniform selection of uncommitted unites between [0-0.6]
0036     T_Uncomm=0.6.*rand(1,1);
0037     Test=abs(Norm_PgMax(:,2)-T_Uncomm);
0038     [Test,ind]=sort(Test);
0039     Norm_PgMax=Norm_PgMax(ind,:);
0040     Norm_Uncommitted_Units=[Norm_Uncommitted_Units;Norm_PgMax(1,:),0];
0041     Norm_PgMax(1,:)=[];
0042 end
0043 Norm_PgMax1=Norm_PgMax;
0044 
0045 %% Committed units (40~50)% of total number of generation units
0046 Ng_Comm=round(Ng*(40+rand*(10))/100);
0047 [Norm_PgMax,Norm_Committed_Units]=Pg_sellection(Norm_PgMax,Ng_Comm,mu_committed);
0048 % Dispatch factor generation
0049 [Alpha]=Alpha_gen(Norm_Committed_Units,Alpha_mod);
0050 % Assignment of Alpha to committed generations
0051 [Norm_Committed_Units,Corr_pg]=Assignment_Pg(Alpha,Norm_Committed_Units,Tab_2D_Pg);
0052 
0053 %% Full load units
0054 Norm_FullLoad_Units=Norm_PgMax;
0055 Norm_FullLoad_Units(:,3)=1;
0056 %% Load balance
0057 [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units);
0058 Norm_tot_generation=sum(AllDispatch_f(:,4));
0059 Norm_Gen_diff=Norm_tot_generation-Norm_Total_load;
0060 
0061 if Norm_Gen_diff > 0 %Generation should be decreased
0062     df_comm=sum(Norm_Committed_Units(:,4))-Norm_Gen_diff;
0063     if df_comm > 0
0064         dk=1;
0065         while dk==1
0066             Norm_Committed_Units(:,3)=Norm_Committed_Units(:,3).*0.99;
0067             [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units);
0068             Norm_tot_generation=sum(AllDispatch_f(:,4));
0069             Norm_Gen_diff=Norm_tot_generation-Norm_Total_load;
0070             if Norm_Gen_diff <= 0
0071                 dk=0;
0072             else dk=1;
0073             end
0074         end
0075     else
0076         dkk=1;
0077         while dkk==1
0078             [max_1, max_1_ind]=max(Norm_FullLoad_Units(:,3));
0079             Norm_FullLoad_Units(max_1_ind(1),3)=0;
0080             [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units);
0081             Norm_tot_generation=sum(AllDispatch_f(:,4));
0082             Norm_Gen_diff=Norm_tot_generation-Norm_Total_load;
0083             if Norm_Gen_diff <= 0
0084                 dkk=0;
0085             else dkk=1;
0086             end
0087         end
0088     end
0089 end
0090             [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units);
0091             Norm_tot_generation=sum(AllDispatch_f(:,4));
0092             Norm_Gen_diff=Norm_tot_generation-Norm_Total_load;
0093 if Norm_Gen_diff < 0 %Generation should be increase
0094     ik=1;
0095     while ik==1
0096         if length(find(Norm_Committed_Units(:,3)<1)) > .5*length(Norm_Committed_Units(:,1))
0097             Norm_Committed_Units(:,3)=Norm_Committed_Units(:,3).*1.05;
0098             Norm_Committed_Units(Norm_Committed_Units(:,3)>1,3)=1;
0099             [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units);
0100             Norm_tot_generation=sum(AllDispatch_f(:,4));
0101             Norm_Gen_diff=Norm_tot_generation-Norm_Total_load;
0102             if Norm_Gen_diff > 0
0103                 ik = 0;
0104             else ik = 1;
0105             end
0106         else
0107             ikk=1;
0108             while ikk==1
0109                 [min_0, min_0_ind]=min(Norm_Uncommitted_Units(:,3));
0110                 if min_0==1
0111                     [min_0, min_0_ind_ful]=min(Norm_FullLoad_Units(:,3));
0112                     Norm_FullLoad_Units(min_0_ind_ful(1),3)=1;
0113                     if min_0==1
0114                     Norm_Committed_Units(ceil(rand*length(Norm_Committed_Units(:,1))),3)=1;
0115                     end
0116                 end
0117                 Norm_Uncommitted_Units(min_0_ind(1),3)=1;
0118                 [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units);
0119                 Norm_tot_generation=sum(AllDispatch_f(:,4));
0120                 Norm_Gen_diff=Norm_tot_generation-Norm_Total_load;
0121                 if Norm_Gen_diff > 0
0122                     ikk=0;
0123                     ik=0;
0124                 else ikk=1;
0125                 end
0126             end
0127         end
0128 %         ik=0;
0129     end
0130 end
0131 Uncommitted_Units=Norm_Uncommitted_Units;
0132 Uncommitted_Units(:,4)=Norm_Uncommitted_Units(:,4).*Max;
0133 Committed_Units=Norm_Committed_Units;
0134 Committed_Units(:,4)=Norm_Committed_Units(:,4).*Max;
0135 FullLoad_Units=Norm_FullLoad_Units;
0136 FullLoad_Units(:,4)=Norm_FullLoad_Units(:,4).*Max;
0137 GenerationDispatch = [Uncommitted_Units(:,[1,4]);Committed_Units(:,[1,4]);FullLoad_Units(:,[1,4])];
0138 Pg_setting = sortrows(GenerationDispatch,1);
0139 
0140 %%
0141 function [Norm_PgMax,Norm_Units]=Pg_sellection(Norm_PgMax,Nunit,mu_unit)
0142 
0143 Ng_unit_90=round(Nunit.*0.99);
0144 Ng_unit_01=Nunit-Ng_unit_90;
0145 Norm_Units=[];
0146 for ic=1:Ng_unit_90
0147     T_unit=sg_exprnd(mu_unit,1,1);
0148     Test=abs(Norm_PgMax(:,2)-T_unit);
0149     [Test,ind]=sort(Test);
0150     Norm_PgMax=Norm_PgMax(ind,:);
0151     Norm_Units=[Norm_Units;Norm_PgMax(1,:)];
0152     Norm_PgMax(1,:)=[];
0153 end
0154 if Ng_unit_01>0
0155     for icc=1:Ng_unit_01
0156         T_comm_01=0.5 + 0.5.*rand(1,1);
0157         Test=abs(Norm_PgMax(:,2)-T_comm_01);
0158         [Test,ind]=sort(Test);
0159         Norm_PgMax=Norm_PgMax(ind,:);
0160         Norm_Units=[Norm_Units;Norm_PgMax(1,:)];
0161         Norm_PgMax(1,:)=[];
0162     end
0163 end
0164 
0165 %%
0166 function [Alpha]=Alpha_gen(Norm_Committed_Units,Alpha_mod)
0167 %This is a fuction for random dispatch factor (alpha) generation
0168 
0169 Nc=length(Norm_Committed_Units(:,1));
0170 
0171 %% Uniform random generation
0172 if ~Alpha_mod
0173     % for NYISO and ERCOT 0<a<1
0174     Alpha=rand(Nc,1);
0175 else
0176     % for WECC negative alpha
0177     Nc995=round(Nc.*0.995);
0178     Nc005=Nc-Nc995;
0179     Alpha995=rand(Nc995,1);
0180     Alpha005=-rand(Nc005,1);
0181     Alpha=[Alpha995;Alpha005];
0182 end
0183 
0184 %%
0185 function [Assigned_Committed_units,Corr_pg]=Assignment_Pg(Alpha,Norm_Committed_Units,Tab_2D_Pg)
0186 % This is a function to assign the normalized generation capacities to the generation buses based
0187 % on node degrees
0188 %
0189 % input:
0190 %       Alpha - matrix (Nc by 1) indicating the value of dispatch factors
0191 %       Norm_committed_units - matrix (Nc by 2) indicating thecommitted units and their normalized
0192 %       value
0193 %       Tab_2D - matrix (14 by 10) indicating the assignment pattern
0194 % output:
0195 %       Generation_dispatch- matrix (Nc by 3) indicating the generation capacity setting at each generation bus
0196 %             gen bus number, norm_node degree,  norm_Pg_max
0197 
0198 Nc=length(Norm_Committed_Units);
0199 % calculate actual nnumber of each cell in 2D table nased on the total
0200 % number of generations
0201 Tab_2D_Nc=round(Tab_2D_Pg.*Nc);
0202 % check the mismatch comes from "round". Add or subtract from maximum
0203 % number in the matrix
0204 if sum(sum(Tab_2D_Nc))<Nc
0205     [Max_tab,ind_max_tab]=max(Tab_2D_Nc(:));
0206     [I_row, I_col] = ind2sub(size(Tab_2D_Nc),ind_max_tab);
0207     Tab_2D_Nc(I_row,I_col)=Tab_2D_Nc(I_row,I_col)+(Nc-sum(sum(Tab_2D_Nc)));
0208 elseif sum(sum(Tab_2D_Nc))>Nc
0209     [Max_tab,ind_max_tab]=max(Tab_2D_Nc(:));
0210     [I_row, I_col] = ind2sub(size(Tab_2D_Nc),ind_max_tab);
0211     Tab_2D_Nc(I_row,I_col)=Tab_2D_Nc(I_row,I_col)-(sum(sum(Tab_2D_Nc))-Nc);
0212 end
0213 %% calculate target number of buses based on node degree and 2D table
0214 N_a=zeros(1,10); % matrix for total number of buses in each dispatch factor category
0215 N_a = sum(Tab_2D_Nc,1);
0216 
0217 %% sort node degrees
0218 Sort_alpha=sort(Alpha);
0219 %% assign buses to the node degree categories
0220 A1=Sort_alpha(1:N_a(1,1),:);
0221 Sort_alpha(1:N_a(1,1),:)=[];
0222 
0223 A2=Sort_alpha(1:N_a(1,2),:);
0224 Sort_alpha(1:N_a(1,2),:)=[];
0225 
0226 A3=Sort_alpha(1:N_a(1,3),:);
0227 Sort_alpha(1:N_a(1,3),:)=[];
0228 
0229 A4=Sort_alpha(1:N_a(1,4),:);
0230 Sort_alpha(1:N_a(1,4),:)=[];
0231 
0232 A5=Sort_alpha(1:N_a(1,5),:);
0233 Sort_alpha(1:N_a(1,5),:)=[];
0234 
0235 A6=Sort_alpha(1:N_a(1,6),:);
0236 Sort_alpha(1:N_a(1,6),:)=[];
0237 
0238 A7=Sort_alpha(1:N_a(1,7),:);
0239 Sort_alpha(1:N_a(1,7),:)=[];
0240 
0241 A8=Sort_alpha(1:N_a(1,8),:);
0242 Sort_alpha(1:N_a(1,8),:)=[];
0243 
0244 A9=Sort_alpha(1:N_a(1,9),:);
0245 Sort_alpha(1:N_a(1,9),:)=[];
0246 
0247 A10=Sort_alpha(1:N_a(1,10),:);
0248 Sort_alpha(1:N_a(1,10),:)=[];
0249 
0250 %% calculate target number of generations in normalized generation capacity categories
0251 N_C=zeros(1,14); % matrix for total number of buses in each node degree category
0252 N_C = sum(Tab_2D_Nc,2)';
0253 
0254 %% sort generation capacities and assign them to normalized generation capacities groups of 2D table
0255 Sort_r_PgCm=sortrows(Norm_Committed_Units,2); %changed by zfwang, because Tab_2D is now both in ascending order from 1->Ng
0256 
0257 G1=Sort_r_PgCm(1:N_C(1,1),:);
0258 Sort_r_PgCm(1:N_C(1,1),:)=[];
0259 
0260 G2=Sort_r_PgCm(1:N_C(1,2),:);
0261 Sort_r_PgCm(1:N_C(1,2),:)=[];
0262 
0263 G3=Sort_r_PgCm(1:N_C(1,3),:);
0264 Sort_r_PgCm(1:N_C(1,3),:)=[];
0265 
0266 G4=Sort_r_PgCm(1:N_C(1,4),:);
0267 Sort_r_PgCm(1:N_C(1,4),:)=[];
0268 
0269 G5=Sort_r_PgCm(1:N_C(1,5),:);
0270 Sort_r_PgCm(1:N_C(1,5),:)=[];
0271 
0272 G6=Sort_r_PgCm(1:N_C(1,6),:);
0273 Sort_r_PgCm(1:N_C(1,6),:)=[];
0274 
0275 G7=Sort_r_PgCm(1:N_C(1,7),:);
0276 Sort_r_PgCm(1:N_C(1,7),:)=[];
0277 
0278 G8=Sort_r_PgCm(1:N_C(1,8),:);
0279 Sort_r_PgCm(1:N_C(1,8),:)=[];
0280 
0281 G9=Sort_r_PgCm(1:N_C(1,9),:);
0282 Sort_r_PgCm(1:N_C(1,9),:)=[];
0283 
0284 G10=Sort_r_PgCm(1:N_C(1,10),:);
0285 Sort_r_PgCm(1:N_C(1,10),:)=[];
0286 
0287 G11=Sort_r_PgCm(1:N_C(1,11),:);
0288 Sort_r_PgCm(1:N_C(1,11),:)=[];
0289 
0290 G12=Sort_r_PgCm(1:N_C(1,12),:);
0291 Sort_r_PgCm(1:N_C(1,12),:)=[];
0292 
0293 G13=Sort_r_PgCm(1:N_C(1,13),:);
0294 Sort_r_PgCm(1:N_C(1,13),:)=[];
0295 
0296 G14=Sort_r_PgCm(1:N_C(1,14),:);
0297 Sort_r_PgCm(1:N_C(1,14),:)=[];
0298 %% assign grouped generation capacities to the 2D table cells
0299 % 1) Save generation buses and their related node degrees in cell arreys of
0300 % K_cell
0301 % 2) Save generation capacities in cell arreys of G_cell
0302 % 3) Assigne generation capacities to the node degrees based  : For each
0303 % cell arrey of K_cell (category of node degrees) starting from first array of G_cell based on the number of
0304 % each cell in Tab_2D_Ng assigne generation capacities and then remove
0305 % asigned generation capacities from G_cell.
0306 Al_cell=cell(1,10);
0307 Al_cell{1}=A1; Al_cell{2}=A2; Al_cell{3}=A3; Al_cell{4}=A4; Al_cell{5}=A5; Al_cell{6}=A6; Al_cell{7}=A7;
0308 Al_cell{8}=A8; Al_cell{9}=A9; Al_cell{10}=A10;
0309 
0310 Cm_cell=cell(1,14);
0311 Cm_cell{1}=G1; Cm_cell{2}=G2; Cm_cell{3}=G3; Cm_cell{4}=G4; Cm_cell{5}=G5; Cm_cell{6}=G6; Cm_cell{7}=G7;
0312 Cm_cell{8}=G8; Cm_cell{9}=G9; Cm_cell{10}=G10; Cm_cell{11}=G11; Cm_cell{12}=G12; Cm_cell{13}=G13; Cm_cell{14}=G14;
0313 
0314 for aa=1:10
0315     K_num=1;
0316     for gg=14:-1:1
0317         if Tab_2D_Nc(gg,aa)>0
0318             [samp_G,ind_G]=sg_datasample(Cm_cell{gg}(:,1),Tab_2D_Nc(gg,aa),'replace',false);
0319             Al_cell{aa}(K_num:(K_num + Tab_2D_Nc(gg,aa)-1),2:3)=Cm_cell{gg}(ind_G,:)  ;
0320             K_num=K_num+Tab_2D_Nc(gg,aa);
0321             Cm_cell{gg}(ind_G,:)=[];
0322         else
0323         Al_cell{aa}(K_num:(K_num + Tab_2D_Nc(gg,aa)-1),2:3)=Cm_cell{gg}(1:Tab_2D_Nc(gg,aa),:);
0324         end
0325     end
0326 end
0327 % Gen_dispatch is a martrix (Ng by 3) inclouds Alpha, bus numbers of
0328 % commited units and their generation capacity
0329 Gen_dispatch=[Al_cell{1};Al_cell{2};Al_cell{3};Al_cell{4};Al_cell{5};Al_cell{6};Al_cell{7};Al_cell{8};Al_cell{9};Al_cell{10}];
0330 Rearrange_Gen_dispatch=[Gen_dispatch(:,2:3),Gen_dispatch(:,1)];
0331 Assigned_Committed_units=sortrows(Rearrange_Gen_dispatch,1);
0332 Corr_pg=corr(Assigned_Committed_units(:,3),Assigned_Committed_units(:,2));
0333 
0334 %%
0335 function [Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units,AllDispatch_f]=dispatch_finalizing(Norm_Uncommitted_Units,Norm_Committed_Units,Norm_FullLoad_Units)
0336 % this fuction calculate final the generation of each generation by
0337 % multiplying dispatch factor and Pgmax
0338 Norm_Uncommitted_Units(:,4)=Norm_Uncommitted_Units(:,2).*Norm_Uncommitted_Units(:,3);
0339 Norm_Committed_Units(:,4)=Norm_Committed_Units(:,2).*Norm_Committed_Units(:,3);
0340 Norm_FullLoad_Units(:,4)=Norm_FullLoad_Units(:,2).*Norm_FullLoad_Units(:,3);
0341 AllDispatch_f=[Norm_Uncommitted_Units;Norm_Committed_Units;Norm_FullLoad_Units];