sg_load

SG_LOAD  main function to run load assignment
   [PL_SETTING, CORR_L] = SG_LOAD(LINK_IDS, BTYPES, PGMAX, LOADING, ...
                                   TAB_2D_LOAD)
   Input:
       link_ids - matrix (m by 2) indicating the branch terminal buses
       Btypes - bus type assignment vector (N by 1) of G/L/C (1/2/3)
       PgMax - matrix (Ng x 2) indicating the generation bus numbers and
           generation capacity at each generation bus
       loading - loading level (see SG_OPTIONS)
       Tab_2D_load - PMF map table to assign loads based on node degree

   Output:
       PL_setting - matrix (NL by 2) indicating the load bus numbers
           and load size at each load bus

0001 function [PL_setting] = sg_load(link_ids, Btypes, PgMax, loading, Tab_2D_load)
0002 %SG_LOAD  main function to run load assignment
0003 %   [PL_SETTING, CORR_L] = SG_LOAD(LINK_IDS, BTYPES, PGMAX, LOADING, ...
0004 %                                   TAB_2D_LOAD)
0005 %   Input:
0006 %       link_ids - matrix (m by 2) indicating the branch terminal buses
0007 %       Btypes - bus type assignment vector (N by 1) of G/L/C (1/2/3)
0008 %       PgMax - matrix (Ng x 2) indicating the generation bus numbers and
0009 %           generation capacity at each generation bus
0010 %       loading - loading level (see SG_OPTIONS)
0011 %       Tab_2D_load - PMF map table to assign loads based on node degree
0012 %
0013 %   Output:
0014 %       PL_setting - matrix (NL by 2) indicating the load bus numbers
0015 %           and load size at each load bus
0016 
0017 %   SynGrid
0018 %   Copyright (c) 2017-2018, Electric Power and Energy Systems (EPES) Research Lab
0019 %   by Hamidreza Sadeghian and Zhifang Wang, Virginia Commonwealth University
0020 %
0021 %   This file is part of SynGrid.
0022 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0023 
0024 Gen_tot = sum(PgMax(:,2));
0025 Load_buses = find(Btypes==2);
0026 N=length(Btypes(:,1));
0027 NL=length(Load_buses);
0028 % indicating the node degree at each load bus
0029 Load_nodedegree=nodedegree(Load_buses,link_ids);
0030 % calculate maximum node degree
0031 maxnode=max(Load_nodedegree(:,2));
0032 % Matrix (Nl by 2) indicating bus number of load buses and related
0033 % normalized node degree
0034 Normalized_Load_nodedegree=[Load_nodedegree(:,1),Load_nodedegree(:,2)/maxnode];
0035 % Total load with respect to the network size
0036 switch loading
0037     case 'D'
0038         Total_Load = 10 ^(-0.2*(log10(N))^2+1.98*log10(N)+0.58);
0039     case 'L'
0040         Total_Load = Gen_tot.*(0.3 + rand*0.1);
0041     case 'M'
0042         Total_Load = Gen_tot.*(0.5 + rand*0.1);
0043     case 'H'
0044         Total_Load = Gen_tot.*(0.7 + rand*0.1);
0045 end
0046 % generated loads based on the empirical distribution of realistic power grids
0047 [Normalized_r_PL,MAX_r_PL]=initialLoad(NL,Total_Load);
0048 % Assignment of loads based on the related node degree
0049 % Norm_Load_setting is a matrix ( Nl by 3 ) indicating bus number,
0050 % related normalized node degree and normilized assigned load setting,  respectively.
0051 [Corr_L,Norm_Load_setting]=Assignment_Pl(Normalized_Load_nodedegree,Normalized_r_PL,Tab_2D_load);
0052 % Output: calculate actual values for load settings
0053 Load_setting = Norm_Load_setting;
0054 Load_setting(:,2)= Load_setting(:,2).*maxnode;
0055 Load_setting(:,3)= Load_setting(:,3).*MAX_r_PL;
0056 PL_setting=[Load_setting(:,1),Load_setting(:,3)];
0057 
0058 %%
0059 function Load_nodedegree=nodedegree(Load_buses,link_ids)
0060 % This function calculates the node degree at each load buses
0061 
0062 Load_nodedegree=[];
0063 L1=length(Load_buses(:,1));
0064 L2=length(link_ids(:,1));
0065 
0066 for Bii=1:length(Load_buses(:,1))
0067     count=0;
0068     for Li=1:length(link_ids(:,1))
0069         if link_ids(Li,1)==Load_buses(Bii,1)
0070             count=count+1;
0071         elseif link_ids(Li,2)==Load_buses(Bii,1)
0072             count=count+1;
0073         end
0074     end
0075     Load_nodedegree=[Load_nodedegree;Load_buses(Bii,1),count];
0076 end
0077 
0078 %%
0079 function [Normalized_r_PL,MAX_r_PL]=initialLoad(NL,Total_Load)
0080 
0081 % This is a function to generate loads based on the distribution of realistic grids
0082 % 99% follow the exponential distribution
0083 % 1% of loads take supper large values
0084 
0085 % output:
0086 %        Normalized_r_PL - matrix (NL by 1) indicating the normalized load size at each load bus
0087 
0088 
0089 mu = Total_Load./NL;
0090 P=sg_exprnd(mu,NL,1);
0091 SG=round(0.01*NL);
0092 P_super=max(P)+(2.*rand(SG,1)*max(P));
0093 [samp,ind]=sg_datasample(P,SG,'replace',false);
0094 P(ind)=[];
0095 r_PL=[P;P_super];
0096 
0097 % check scaling property
0098 if (sum(r_PL)>1.05*Total_Load ||  sum(r_PL)<.90*Total_Load)
0099    r_PL=r_PL.*(Total_Load./sum(r_PL));
0100 end
0101 
0102 MAX_r_PL=max(r_PL);
0103 
0104 Normalized_r_PL=r_PL/MAX_r_PL;
0105 
0106 
0107 function [Corr_coL,Norm_Load_setting]=Assignment_Pl(Normalized_Load_nodedegree,Normalized_r_PL,Tab_2D)
0108 % This is a function to assign the normalized generation capacities to the generation buses besed
0109 % on node degrees
0110 %
0111 % input:
0112 %       Normalized_Generation_nodedegree- matrix (2 by Ng) indicating the normalized value of node degrees
0113 %       Normalized_r_Pg- matrix (1 by Ng) indicating the normalized value of generation capacities
0114 %       Tab_2D - matrix (14 by 14) indicating the assignment pattern
0115 % output:
0116 %       Generation_setting- matrix (Ng by 3) indicating the load setting at each bus
0117 
0118 NL=length(Normalized_r_PL);
0119 % calculate actual nnumber of each cell in 2D table nased on the total
0120 % number of generation capacities
0121 Tab_2D_NL=round(Tab_2D.*NL);
0122 % check the mismatch comes from "round". Add or subtract from maximum
0123 % number in the matrix
0124 if sum(sum(Tab_2D_NL))<NL
0125     [Max_tab,ind_max_tab]=max(Tab_2D_NL(:));
0126     [I_row, I_col] = ind2sub(size(Tab_2D_NL),ind_max_tab);
0127     Tab_2D_NL(I_row,I_col)=Tab_2D_NL(I_row,I_col)+(NL-sum(sum(Tab_2D_NL)));
0128 elseif sum(sum(Tab_2D_NL))>NL
0129     [Max_tab,ind_max_tab]=max(Tab_2D_NL(:));
0130     [I_row, I_col] = ind2sub(size(Tab_2D_NL),ind_max_tab);
0131     Tab_2D_NL(I_row,I_col)=Tab_2D_NL(I_row,I_col)-(sum(sum(Tab_2D_NL))-NL);
0132 end
0133 %% calculate target number of buses based on node degree and 2D table
0134 N_K=zeros(1,14); % matrix for total number of buses in each node degree category
0135 N_K = sum(Tab_2D_NL,1);
0136 
0137 %% sort node degrees
0138 Sort_nodedegree=sortrows(Normalized_Load_nodedegree,2);
0139 %% assign buses to the node degree categories
0140 K1=Sort_nodedegree(1:N_K(1,1),:);
0141 Sort_nodedegree(1:N_K(1,1),:)=[];
0142 
0143 K2=Sort_nodedegree(1:N_K(1,2),:);
0144 Sort_nodedegree(1:N_K(1,2),:)=[];
0145 
0146 K3=Sort_nodedegree(1:N_K(1,3),:);
0147 Sort_nodedegree(1:N_K(1,3),:)=[];
0148 
0149 K4=Sort_nodedegree(1:N_K(1,4),:);
0150 Sort_nodedegree(1:N_K(1,4),:)=[];
0151 
0152 K5=Sort_nodedegree(1:N_K(1,5),:);
0153 Sort_nodedegree(1:N_K(1,5),:)=[];
0154 
0155 K6=Sort_nodedegree(1:N_K(1,6),:);
0156 Sort_nodedegree(1:N_K(1,6),:)=[];
0157 
0158 K7=Sort_nodedegree(1:N_K(1,7),:);
0159 Sort_nodedegree(1:N_K(1,7),:)=[];
0160 
0161 K8=Sort_nodedegree(1:N_K(1,8),:);
0162 Sort_nodedegree(1:N_K(1,8),:)=[];
0163 
0164 K9=Sort_nodedegree(1:N_K(1,9),:);
0165 Sort_nodedegree(1:N_K(1,9),:)=[];
0166 
0167 K10=Sort_nodedegree(1:N_K(1,10),:);
0168 Sort_nodedegree(1:N_K(1,10),:)=[];
0169 
0170 K11=Sort_nodedegree(1:N_K(1,11),:);
0171 Sort_nodedegree(1:N_K(1,11),:)=[];
0172 
0173 K12=Sort_nodedegree(1:N_K(1,12),:);
0174 Sort_nodedegree(1:N_K(1,12),:)=[];
0175 
0176 K13=Sort_nodedegree(1:N_K(1,13),:);
0177 Sort_nodedegree(1:N_K(1,13),:)=[];
0178 
0179 K14=Sort_nodedegree(1:N_K(1,14),:);
0180 Sort_nodedegree(1:N_K(1,14),:)=[];
0181 
0182 %% calculate target number of generations in normalized generation capacity categories
0183 N_L=zeros(1,14); % matrix for total number of buses in each node degree category
0184 N_L = sum(Tab_2D_NL,2)';
0185 
0186 %% sort generation capacities and assign them to normalized generation capacities groups of 2D table
0187 % Sort_r_PL=sort(Normalized_r_PL,'descend');
0188 Sort_r_PL=sort(Normalized_r_PL);
0189 G1=Sort_r_PL(1:N_L(1,1),:);
0190 Sort_r_PL(1:N_L(1,1),:)=[];
0191 
0192 G2=Sort_r_PL(1:N_L(1,2),:);
0193 Sort_r_PL(1:N_L(1,2),:)=[];
0194 
0195 G3=Sort_r_PL(1:N_L(1,3),:);
0196 Sort_r_PL(1:N_L(1,3),:)=[];
0197 
0198 G4=Sort_r_PL(1:N_L(1,4),:);
0199 Sort_r_PL(1:N_L(1,4),:)=[];
0200 
0201 G5=Sort_r_PL(1:N_L(1,5),:);
0202 Sort_r_PL(1:N_L(1,5),:)=[];
0203 
0204 G6=Sort_r_PL(1:N_L(1,6),:);
0205 Sort_r_PL(1:N_L(1,6),:)=[];
0206 
0207 G7=Sort_r_PL(1:N_L(1,7),:);
0208 Sort_r_PL(1:N_L(1,7),:)=[];
0209 
0210 G8=Sort_r_PL(1:N_L(1,8),:);
0211 Sort_r_PL(1:N_L(1,8),:)=[];
0212 
0213 G9=Sort_r_PL(1:N_L(1,9),:);
0214 Sort_r_PL(1:N_L(1,9),:)=[];
0215 
0216 G10=Sort_r_PL(1:N_L(1,10),:);
0217 Sort_r_PL(1:N_L(1,10),:)=[];
0218 
0219 G11=Sort_r_PL(1:N_L(1,11),:);
0220 Sort_r_PL(1:N_L(1,11),:)=[];
0221 
0222 G12=Sort_r_PL(1:N_L(1,12),:);
0223 Sort_r_PL(1:N_L(1,12),:)=[];
0224 
0225 G13=Sort_r_PL(1:N_L(1,13),:);
0226 Sort_r_PL(1:N_L(1,13),:)=[];
0227 
0228 G14=Sort_r_PL(1:N_L(1,14),:);
0229 Sort_r_PL(1:N_L(1,14),:)=[];
0230 %% assign grouped generation capacities to the 2D table cells
0231 % 1) Save generation buses and their related node degrees in cell arreys of
0232 % K_cell
0233 % 2) Save load settings in cell arreys of L_cell
0234 % 3) Assigne generation capacities to the node degrees based  : For each
0235 % cell arrey of K_cell (category of node degrees) starting from first array of L_cell based on the number of
0236 % each cell in Tab_2D_Ng assigne load settings and then remove
0237 % asigned generation capacities from L_cell.
0238 K_cell=cell(1,14);
0239 K_cell{1}=K1; K_cell{2}=K2; K_cell{3}=K3; K_cell{4}=K4; K_cell{5}=K5; K_cell{6}=K6; K_cell{7}=K7;
0240 K_cell{8}=K8; K_cell{9}=K9; K_cell{10}=K10; K_cell{11}=K11; K_cell{12}=K12; K_cell{13}=K13; K_cell{14}=K14;
0241 
0242 L_cell=cell(1,14);
0243 L_cell{1}=G1; L_cell{2}=G2; L_cell{3}=G3; L_cell{4}=G4; L_cell{5}=G5; L_cell{6}=G6; L_cell{7}=G7;
0244 L_cell{8}=G8; L_cell{9}=G9; L_cell{10}=G10; L_cell{11}=G11; L_cell{12}=G12; L_cell{13}=G13; L_cell{14}=G14;
0245 
0246 for kk=1:14
0247     K_num=1;
0248     for gg=14:-1:1
0249         if Tab_2D_NL(gg,kk)>0
0250             [samp_L,ind_L]=sg_datasample(L_cell{gg},Tab_2D_NL(gg,kk),'replace',false);
0251             K_cell{kk}(K_num:(K_num + Tab_2D_NL(gg,kk)-1),3)=samp_L;
0252             K_num=K_num+Tab_2D_NL(gg,kk);
0253             L_cell{gg}(ind_L)=[];
0254         else
0255             K_cell{kk}(K_num:(K_num + Tab_2D_NL(gg,kk)-1),3)=L_cell{gg}(1:Tab_2D_NL(gg,kk));
0256         end
0257     end
0258 end
0259 % Load_st is a martrix (Ng by 3) inclouds bus numbers, normalized node degrees and normalized Load settings
0260 Load_st=[K_cell{1};K_cell{2};K_cell{3};K_cell{4};K_cell{5};K_cell{6};K_cell{7};K_cell{8};K_cell{9};K_cell{10};K_cell{11};K_cell{12};K_cell{13};K_cell{14}];
0261 Norm_Load_setting=sortrows(Load_st,1);
0262 Corr_coL=corr(Norm_Load_setting(:,2),Norm_Load_setting(:,3));